From 8d09cf9bebee1dec094cfcfb86b23665f742d51e Mon Sep 17 00:00:00 2001 From: Sebastian Kosch Date: Mon, 6 Jul 2020 23:43:25 -0700 Subject: [PATCH] Fix rendering of uniforms; add normalization constant in convolution code --- src/components/DistBuilder.re | 3 +- .../distribution/AlgebraicShapeCombination.re | 16 ++++----- src/distPlus/distribution/Distributions.re | 33 ++++++++++--------- .../expressionTree/ExpressionTreeEvaluator.re | 7 ++-- src/distPlus/symbolic/SymbolicDist.re | 4 +-- 5 files changed, 34 insertions(+), 29 deletions(-) diff --git a/src/components/DistBuilder.re b/src/components/DistBuilder.re index f087a6bf..6ab65759 100644 --- a/src/components/DistBuilder.re +++ b/src/components/DistBuilder.re @@ -146,7 +146,6 @@ module DemoDist = { (), ); let response = DistPlusRenderer.run(inputs); - Js.log(response); switch (RenderTypes.DistPlusRenderer.Outputs.distplus(response)) { | Some(distPlus) => | _ => @@ -172,7 +171,7 @@ let make = () => { ~onSubmit=({state}) => {None}, ~initialState={ //guesstimatorString: "mm(normal(-10, 2), uniform(18, 25), lognormal({mean: 10, stdev: 8}), triangular(31,40,50))", - guesstimatorString: "truncate(normal(100, 60), 50, 100)", + guesstimatorString: "uniform(1,2) * uniform(2, 3)", // , triangular(30, 40, 60) domainType: "Complete", xPoint: "50.0", xPoint2: "60.0", diff --git a/src/distPlus/distribution/AlgebraicShapeCombination.re b/src/distPlus/distribution/AlgebraicShapeCombination.re index 5fadda1c..49fc0d49 100644 --- a/src/distPlus/distribution/AlgebraicShapeCombination.re +++ b/src/distPlus/distribution/AlgebraicShapeCombination.re @@ -139,8 +139,8 @@ let combineShapesContinuousContinuous = // converts the variances and means of the two inputs into the variance of the output let combineVariancesFn = switch (op) { - | `Add => ((v1, v2, m1, m2) => v1 +. v2) - | `Subtract => ((v1, v2, m1, m2) => v1 +. v2) + | `Add => ((v1, v2, _, _) => v1 +. v2) + | `Subtract => ((v1, v2, _, _) => v1 +. v2) | `Multiply => ( (v1, v2, m1, m2) => v1 *. v2 +. v1 *. m2 ** 2. +. v2 *. m1 ** 2. ) @@ -176,8 +176,8 @@ let combineShapesContinuousContinuous = let _ = Belt.Array.set(means, k, mean); let _ = Belt.Array.set(variances, k, variance); // update bounds - let minX = mean -. variance *. 1.644854; - let maxX = mean +. variance *. 1.644854; + let minX = mean -. 2. *. sqrt(variance) *. 1.644854; + let maxX = mean +. 2. *. sqrt(variance) *. 1.644854; if (minX < outputMinX^) { outputMinX := minX; }; @@ -193,11 +193,11 @@ let combineShapesContinuousContinuous = let outputXs: array(float) = E.A.Floats.range(outputMinX^, outputMaxX^, nOut); let outputYs: array(float) = Belt.Array.make(nOut, 0.0); // now, for each of the outputYs, accumulate from a Gaussian kernel over each input point. - for (j in 0 to E.A.length(masses) - 1) { - let _ = if (variances[j] > 0.) { - for (i in 0 to E.A.length(outputXs) - 1) { + for (j in 0 to E.A.length(masses) - 1) { // go through all of the result points + let _ = if (variances[j] > 0. && masses[j] > 0.) { + for (i in 0 to E.A.length(outputXs) - 1) { // go through all of the target points let dx = outputXs[i] -. means[j]; - let contribution = masses[j] *. exp(-. (dx ** 2.) /. (2. *. variances[j])); + let contribution = masses[j] *. exp(-. (dx ** 2.) /. (2. *. variances[j])) /. (sqrt(2. *. 3.14159276 *. variances[j])); let _ = Belt.Array.set(outputYs, i, outputYs[i] +. contribution); (); }; diff --git a/src/distPlus/distribution/Distributions.re b/src/distPlus/distribution/Distributions.re index f71afccb..963508ab 100644 --- a/src/distPlus/distribution/Distributions.re +++ b/src/distPlus/distribution/Distributions.re @@ -112,7 +112,7 @@ module Continuous = { t2.knownIntegralSum, ); - make( + let res = make( `Linear, XYShape.PointwiseCombination.combine( ~xsSelection=ALL_XS, @@ -123,6 +123,7 @@ module Continuous = { ), combinedIntegralSum, ); + res }; let toLinear = (t: t): option(t) => { @@ -219,7 +220,9 @@ module Continuous = { let integral = (~cache, t) => if (t |> getShape |> XYShape.T.length > 0) { switch (cache) { - | Some(cache) => cache + | Some(cache) => { + cache; + } | None => t |> getShape @@ -243,7 +246,7 @@ module Continuous = { t.knownIntegralSum |> E.O.default(t |> integral(~cache) |> lastY); let integralXtoY = (~cache, f, t: t) => t |> integral(~cache) |> shapeFn(XYShape.XtoY.linear(f)); - let integralYtoX = (~cache, f, t: t) => + let integralYtoX = (~cache, f, t: t) => t |> integral(~cache) |> shapeFn(XYShape.YtoX.linear(f)); let toContinuous = t => Some(t); let toDiscrete = _ => None; @@ -254,7 +257,7 @@ module Continuous = { |> updateKnownIntegralSum(Some(1.0)); }; - let normalizedToContinuous = t => Some(t); // TODO: this should be normalized + let normalizedToContinuous = t => Some(t |> normalize); let normalizedToDiscrete = _ => None; let mean = (t: t) => { @@ -1059,18 +1062,17 @@ module Shape = { Continuous.T.normalizedToContinuous, )); let minX = mapToAll((Mixed.T.minX, Discrete.T.minX, Continuous.T.minX)); - let integral = (~cache) => { + let integral = (~cache) => mapToAll(( - Mixed.T.Integral.get(~cache), - Discrete.T.Integral.get(~cache), - Continuous.T.Integral.get(~cache), + Mixed.T.Integral.get(~cache=None), + Discrete.T.Integral.get(~cache=None), + Continuous.T.Integral.get(~cache=None), )); - }; let integralEndY = (~cache) => mapToAll(( - Mixed.T.Integral.sum(~cache), + Mixed.T.Integral.sum(~cache=None), Discrete.T.Integral.sum(~cache), - Continuous.T.Integral.sum(~cache), + Continuous.T.Integral.sum(~cache=None), )); let integralXtoY = (~cache, f) => { mapToAll(( @@ -1178,7 +1180,6 @@ module DistPlus = { let normalize = (t: t): t => { let normalizedShape = t |> toShape |> Shape.T.normalize; - t |> updateShape(normalizedShape); // TODO: also adjust for domainIncludedProbabilityMass here. }; @@ -1190,7 +1191,6 @@ module DistPlus = { t |> updateShape(truncatedShape); }; - // TODO: replace this with let normalizedToContinuous = (t: t) => { t |> toShape @@ -1236,8 +1236,10 @@ module DistPlus = { : t => Shape.T.mapY(~knownIntegralSumFn, fn, shape) |> updateShape(_, t); - let integralEndY = (~cache as _, t: t) => + // get the total of everything + let integralEndY = (~cache as _, t: t) => { Shape.T.Integral.sum(~cache=Some(t.integralCache), toShape(t)); + } // TODO: Fix this below, obviously. Adjust for limits let integralXtoY = (~cache as _, f, t: t) => { @@ -1247,8 +1249,9 @@ module DistPlus = { // TODO: This part is broken when there is a limit, if this is supposed to be taken into account. let integralYtoX = (~cache as _, f, t: t) => { - Shape.T.Integral.yToX(~cache=Some(t.integralCache), f, toShape(t)); + Shape.T.Integral.yToX(~cache=None, f, toShape(t)); }; + let mean = (t: t) => { Shape.T.mean(t.shape); }; diff --git a/src/distPlus/expressionTree/ExpressionTreeEvaluator.re b/src/distPlus/expressionTree/ExpressionTreeEvaluator.re index 23853db8..d5cd807f 100644 --- a/src/distPlus/expressionTree/ExpressionTreeEvaluator.re +++ b/src/distPlus/expressionTree/ExpressionTreeEvaluator.re @@ -168,8 +168,11 @@ module Normalize = { let rec operationToLeaf = (toLeaf, renderParams, t: node): result(node, string) => { switch (t) { - | `RenderedDist(s) => - Ok(`RenderedDist(Distributions.Shape.T.normalize(s))) + | `RenderedDist(s) => { + Js.log2("normed", Distributions.Shape.T.normalize(s)); + //Js.log2("normed integral", Distributions.Shape.T.normalize(s))); + Ok(`RenderedDist(Distributions.Shape.T.normalize(s))); + } | `SymbolicDist(_) => Ok(t) | _ => t diff --git a/src/distPlus/symbolic/SymbolicDist.re b/src/distPlus/symbolic/SymbolicDist.re index 3b4fb065..d4cbbf4f 100644 --- a/src/distPlus/symbolic/SymbolicDist.re +++ b/src/distPlus/symbolic/SymbolicDist.re @@ -258,11 +258,11 @@ module T = { (~xSelection: [ | `Linear | `ByWeight]=`Linear, dist: symbolicDist, n) => { switch (xSelection, dist) { | (`Linear, _) => E.A.Floats.range(min(dist), max(dist), n) - /* | (`ByWeight, `Uniform(n)) => + | (`ByWeight, `Uniform(n)) => // In `ByWeight mode, uniform distributions get special treatment because we need two x's // on either side for proper rendering (just left and right of the discontinuities). let dx = 0.00001 *. (n.high -. n.low); - [|n.low -. dx, n.low +. dx, n.high -. dx, n.high +. dx|]; */ + [|n.low -. dx, n.low +. dx, n.high -. dx, n.high +. dx|]; | (`ByWeight, _) => let ys = E.A.Floats.range(minCdfValue, maxCdfValue, n); ys |> E.A.fmap(y => inv(y, dist));