quantified.uncertainty/squiggle
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

feat: Audit SymbolicDist.res

Browse Source 
- Fix buggy lognormal multiplication code
- Add precision to 90% confidence intervals code
- Simplified lognormal code
- Added sources for many of the manipulations
This commit is contained in:
NunoSempere 2022-04-14 16:03:54 -04:00
parent c43b373681
commit ec9c67f090
1 changed files with 15 additions and 11 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

22
packages/squiggle-lang/src/rescript/Distributions/SymbolicDist/SymbolicDist.res
View File

@ -11,7 +11,7 @@ module Normal = {
let from90PercentCI = (low, high) => { let from90PercentCI = (low, high) => {
let mean = E.A.Floats.mean([low, high]) let mean = E.A.Floats.mean([low, high])
let stdev = (high -. low) /. (2. *. 1.644854) let stdev = (high -. low) /. (2. *. 1.6448536269514722)
#Normal({mean: mean, stdev: stdev}) #Normal({mean: mean, stdev: stdev})
} }
let inv = (p, t: t) => Jstat.Normal.inv(p, t.mean, t.stdev) let inv = (p, t: t) => Jstat.Normal.inv(p, t.mean, t.stdev)
@ -21,12 +21,12 @@ module Normal = {
let add = (n1: t, n2: t) => { let add = (n1: t, n2: t) => {
let mean = n1.mean +. n2.mean let mean = n1.mean +. n2.mean
let stdev = sqrt(n1.stdev ** 2. +. n2.stdev ** 2.) let stdev = Js.Math.sqrt(n1.stdev ** 2. +. n2.stdev ** 2.)
#Normal({mean: mean, stdev: stdev}) #Normal({mean: mean, stdev: stdev})
} }
let subtract = (n1: t, n2: t) => { let subtract = (n1: t, n2: t) => {
let mean = n1.mean -. n2.mean let mean = n1.mean -. n2.mean
let stdev = sqrt(n1.stdev ** 2. +. n2.stdev ** 2.) let stdev = Js.Math.sqrt(n1.stdev ** 2. +. n2.stdev ** 2.)
#Normal({mean: mean, stdev: stdev}) #Normal({mean: mean, stdev: stdev})
} }
@ -115,19 +115,22 @@ module Lognormal = {
let mean = (t: t) => Ok(Jstat.Lognormal.mean(t.mu, t.sigma)) let mean = (t: t) => Ok(Jstat.Lognormal.mean(t.mu, t.sigma))
let sample = (t: t) => Jstat.Lognormal.sample(t.mu, t.sigma) let sample = (t: t) => Jstat.Lognormal.sample(t.mu, t.sigma)
let toString = ({mu, sigma}: t) => j`Lognormal($mu,$sigma)` let toString = ({mu, sigma}: t) => j`Lognormal($mu,$sigma)`
let from90PercentCI = (low, high) => { let from90PercentCI = (low, high) => {
let logLow = Js.Math.log(low) let logLow = Js.Math.log(low)
let logHigh = Js.Math.log(high) let logHigh = Js.Math.log(high)
let mu = E.A.Floats.mean([logLow, logHigh]) let mu = E.A.Floats.mean([logLow, logHigh])
let sigma = (logHigh -. logLow) /. (2.0 *. 1.645) let sigma = (logHigh -. logLow) /. (2.0 *. 1.6448536269514722)
#Lognormal({mu: mu, sigma: sigma}) #Lognormal({mu: mu, sigma: sigma})
} }
let fromMeanAndStdev = (mean, stdev) => { let fromMeanAndStdev = (mean, stdev) => {
// https://math.stackexchange.com/questions/2501783/parameters-of-a-lognormal-distribution
// https://wikiless.org/wiki/Log-normal_distribution?lang=en#Generation_and_parameters
if stdev > 0.0 { if stdev > 0.0 {
let variance = Js.Math.pow_float(~base=stdev, ~exp=2.0) let variance = stdev ** 2.
let meanSquared = Js.Math.pow_float(~base=mean, ~exp=2.0) let meanSquared = mean ** 2.
let mu = Js.Math.log(mean) -. 0.5 *. Js.Math.log(variance /. meanSquared +. 1.0) let mu = 2*Js.Math.log(mean) -. 0.5 *. Js.Math.log(variance +. meanSquared)
let sigma = Js.Math.pow_float(~base=Js.Math.log(variance /. meanSquared +. 1.0), ~exp=0.5) let sigma = Js.Math.sqrt(Js.Math.log((variance /. meanSquared) +. 1) )
Ok(#Lognormal({mu: mu, sigma: sigma})) Ok(#Lognormal({mu: mu, sigma: sigma}))
} else { } else {
Error("Lognormal standard deviation must be larger than 0") Error("Lognormal standard deviation must be larger than 0")
@ -135,8 +138,9 @@ module Lognormal = {
} }
let multiply = (l1, l2) => { let multiply = (l1, l2) => {
// https://wikiless.org/wiki/Log-normal_distribution?lang=en#Multiplication_and_division_of_independent,_log-normal_random_variables
let mu = l1.mu +. l2.mu let mu = l1.mu +. l2.mu
let sigma = l1.sigma +. l2.sigma let sigma = Math.sqrt(l1.sigma ** 2. +. l2.sigma ** 2.)
#Lognormal({mu: mu, sigma: sigma}) #Lognormal({mu: mu, sigma: sigma})
} }
let divide = (l1, l2) => { let divide = (l1, l2) => {