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

Ozzie's suggestions

Browse Source
This commit is contained in:
Sebastian Kosch 2020-07-18 21:34:55 -07:00
parent 9cdcc85d6c
commit f3922761f4
13 changed files with 251 additions and 285 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

163
src/distPlus/distribution/AlgebraicShapeCombination.re
View File

@ -17,25 +17,25 @@ let toDiscretePointMassesFromTriangulars =
let n = s |> XYShape.T.length;
// first, double up the leftmost and rightmost points:
let {xs, ys}: XYShape.T.t = s;
let _ = Js.Array.unshift(xs[0], xs);
let _ = Js.Array.unshift(ys[0], ys);
let _ = Js.Array.push(xs[n - 1], xs);
let _ = Js.Array.push(ys[n - 1], ys);
Js.Array.unshift(xs[0], xs) |> ignore;
Js.Array.unshift(ys[0], ys) |> ignore;
Js.Array.push(xs[n - 1], xs) |> ignore;
Js.Array.push(ys[n - 1], ys) |> ignore;
let n = E.A.length(xs);
// squares and neighbourly products of the xs
let xsSq: array(float) = Belt.Array.makeUninitializedUnsafe(n);
let xsProdN1: array(float) = Belt.Array.makeUninitializedUnsafe(n - 1);
let xsProdN2: array(float) = Belt.Array.makeUninitializedUnsafe(n - 2);
for (i in 0 to n - 1) {
let _ = Belt.Array.set(xsSq, i, xs[i] *. xs[i]);
Belt.Array.set(xsSq, i, xs[i] *. xs[i]) |> ignore;
();
};
for (i in 0 to n - 2) {
let _ = Belt.Array.set(xsProdN1, i, xs[i] *. xs[i + 1]);
Belt.Array.set(xsProdN1, i, xs[i] *. xs[i + 1]) |> ignore;
();
};
for (i in 0 to n - 3) {
let _ = Belt.Array.set(xsProdN2, i, xs[i] *. xs[i + 2]);
Belt.Array.set(xsProdN2, i, xs[i] *. xs[i + 2]) |> ignore;
();
};
// means and variances
@ -45,12 +45,11 @@ let toDiscretePointMassesFromTriangulars =
if (inverse) {
for (i in 1 to n - 2) {
let _ =
Belt.Array.set(
masses,
i - 1,
(xs[i + 1] -. xs[i - 1]) *. ys[i] /. 2.,
);
Belt.Array.set(
masses,
i - 1,
(xs[i + 1] -. xs[i - 1]) *. ys[i] /. 2.,
) |> ignore;
// this only works when the whole triange is either on the left or on the right of zero
let a = xs[i - 1];
@ -71,43 +70,39 @@ let toDiscretePointMassesFromTriangulars =
-. inverseMean
** 2.;
let _ = Belt.Array.set(means, i - 1, inverseMean);
Belt.Array.set(means, i - 1, inverseMean) |> ignore;
let _ = Belt.Array.set(variances, i - 1, inverseVar);
Belt.Array.set(variances, i - 1, inverseVar) |> ignore;
();
};
{n: n - 2, masses, means, variances};
} else {
for (i in 1 to n - 2) {
// area of triangle = width * height / 2
let _ =
Belt.Array.set(
masses,
i - 1,
(xs[i + 1] -. xs[i - 1]) *. ys[i] /. 2.,
);
Belt.Array.set(
masses,
i - 1,
(xs[i + 1] -. xs[i - 1]) *. ys[i] /. 2.,
) |> ignore;
// means of triangle = (a + b + c) / 3
let _ =
Belt.Array.set(means, i - 1, (xs[i - 1] +. xs[i] +. xs[i + 1]) /. 3.);
Belt.Array.set(means, i - 1, (xs[i - 1] +. xs[i] +. xs[i + 1]) /. 3.) |> ignore;
// variance of triangle = (a^2 + b^2 + c^2 - ab - ac - bc) / 18
let _ =
Belt.Array.set(
variances,
i - 1,
(
xsSq[i - 1]
+. xsSq[i]
+. xsSq[i + 1]
-. xsProdN1[i - 1]
-. xsProdN1[i]
-. xsProdN2[i - 1]
)
/. 18.,
);
Belt.Array.set(
variances,
i - 1,
(
xsSq[i - 1]
+. xsSq[i]
+. xsSq[i + 1]
-. xsProdN1[i - 1]
-. xsProdN1[i]
-. xsProdN2[i - 1]
)
/. 18.,
) |> ignore;
();
};
{n: n - 2, masses, means, variances};
@ -163,7 +158,7 @@ let combineShapesContinuousContinuous =
for (i in 0 to t1m.n - 1) {
for (j in 0 to t2m.n - 1) {
let k = i * t2m.n + j;
let _ = Belt.Array.set(masses, k, t1m.masses[i] *. t2m.masses[j]);
Belt.Array.set(masses, k, t1m.masses[i] *. t2m.masses[j]) |> ignore;
let mean = combineMeansFn(t1m.means[i], t2m.means[j]);
let variance =
@ -173,8 +168,8 @@ let combineShapesContinuousContinuous =
t1m.means[i],
t2m.means[j],
);
let _ = Belt.Array.set(means, k, mean);
let _ = Belt.Array.set(variances, k, variance);
Belt.Array.set(means, k, mean) |> ignore;
Belt.Array.set(variances, k, variance) |> ignore;
// update bounds
let minX = mean -. 2. *. sqrt(variance) *. 1.644854;
let maxX = mean +. 2. *. sqrt(variance) *. 1.644854;
@ -194,15 +189,15 @@ let combineShapesContinuousContinuous =
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) { // go through all of the result points
let _ = if (variances[j] > 0. && masses[j] > 0.) {
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])) /. (sqrt(2. *. 3.14159276 *. variances[j]));
let _ = Belt.Array.set(outputYs, i, outputYs[i] +. contribution);
Belt.Array.set(outputYs, i, outputYs[i] +. contribution) |> ignore;
();
};
();
};
} |> ignore;
();
};
@ -219,26 +214,23 @@ let toDiscretePointMassesFromDiscrete = (s: DistTypes.xyShape): pointMassesWithM
let variances: array(float) = Belt.Array.makeUninitializedUnsafe(n);
for (i in 0 to n - 1) {
let _ =
Belt.Array.set(
masses,
i,
ys[i]
);
Belt.Array.set(
masses,
i,
ys[i]
) |> ignore;
let _ =
Belt.Array.set(
means,
i,
xs[i]
);
Belt.Array.set(
means,
i,
xs[i]
) |> ignore;
let _ =
Belt.Array.set(
variances,
i,
0.0
);
Belt.Array.set(
variances,
i,
0.0
) |> ignore;
();
};
@ -246,11 +238,11 @@ let toDiscretePointMassesFromDiscrete = (s: DistTypes.xyShape): pointMassesWithM
};
let combineShapesContinuousDiscrete =
(op: ExpressionTypes.algebraicOperation, s1: DistTypes.xyShape, s2: DistTypes.xyShape)
: array(DistTypes.xyShape) => {
(op: ExpressionTypes.algebraicOperation, continuousShape: DistTypes.xyShape, discreteShape: DistTypes.xyShape)
: DistTypes.xyShape => {
let t1n = s1 |> XYShape.T.length;
let t2n = s2 |> XYShape.T.length;
let t1n = continuousShape |> XYShape.T.length;
let t2n = discreteShape |> XYShape.T.length;
// each x pair is added/subtracted
let fn = Operation.Algebraic.toFn(op);
@ -262,46 +254,45 @@ let combineShapesContinuousDiscrete =
| `Add
| `Subtract => {
for (j in 0 to t2n - 1) {
// for each one of the discrete points
// create a new distribution, as long as the original continuous one
// creates a new continuous shape for each one of the discrete points, and collects them in outXYShapes.
let dxyShape: array((float, float)) =
Belt.Array.makeUninitializedUnsafe(t1n);
for (i in 0 to t1n - 1) {
let _ =
Belt.Array.set(
dxyShape,
i,
(fn(s1.xs[i], s2.xs[j]), s1.ys[i] *. s2.ys[j]),
);
Belt.Array.set(
dxyShape,
i,
(fn(continuousShape.xs[i], discreteShape.xs[j]), continuousShape.ys[i] *. discreteShape.ys[j]),
) |> ignore;
();
};
let _ = Belt.Array.set(outXYShapes, j, dxyShape);
Belt.Array.set(outXYShapes, j, dxyShape) |> ignore;
();
}
}
| `Multiply
| `Divide => {
for (j in 0 to t2n - 1) {
// for each one of the discrete points
// create a new distribution, as long as the original continuous one
// creates a new continuous shape for each one of the discrete points, and collects them in outXYShapes.
let dxyShape: array((float, float)) =
Belt.Array.makeUninitializedUnsafe(t1n);
for (i in 0 to t1n - 1) {
let _ =
Belt.Array.set(
dxyShape,
i,
(fn(s1.xs[i], s2.xs[j]), s1.ys[i] *. s2.ys[j] /. s2.xs[j]),
);
Belt.Array.set(
dxyShape,
i,
(fn(continuousShape.xs[i], discreteShape.xs[j]), continuousShape.ys[i] *. discreteShape.ys[j] /. discreteShape.xs[j]),
) |> ignore;
();
};
let _ = Belt.Array.set(outXYShapes, j, dxyShape);
Belt.Array.set(outXYShapes, j, dxyShape) |> ignore;
();
}
}
};
outXYShapes
|> E.A.fmap(XYShape.T.fromZippedArray);
|> E.A.fmap(XYShape.T.fromZippedArray)
|> E.A.fold_left(
XYShape.PointwiseCombination.combine((+.),
XYShape.XtoY.continuousInterpolator(`Linear, `UseZero)),
XYShape.T.empty);
};

77
src/distPlus/distribution/Continuous.re
View File

@ -3,7 +3,7 @@ open Distributions;
type t = DistTypes.continuousShape;
let getShape = (t: t) => t.xyShape;
let interpolation = (t: t) => t.interpolation;
let make = (interpolation, xyShape, integralSumCache, integralCache): t => {
let make = (~interpolation=`Linear, ~integralSumCache=None, ~integralCache=None, xyShape): t => {
xyShape,
interpolation,
integralSumCache,
@ -19,7 +19,7 @@ let lastY = (t: t) => t |> getShape |> XYShape.T.lastY;
let oShapeMap =
(fn, {xyShape, interpolation, integralSumCache, integralCache}: t)
: option(DistTypes.continuousShape) =>
fn(xyShape) |> E.O.fmap(make(interpolation, _, integralSumCache, integralCache));
fn(xyShape) |> E.O.fmap(make(~interpolation, ~integralSumCache, ~integralCache));
let emptyIntegral: DistTypes.continuousShape = {
xyShape: {xs: [|neg_infinity|], ys: [|0.0|]},
@ -35,7 +35,7 @@ let empty: DistTypes.continuousShape = {
};
let stepwiseToLinear = (t: t): t =>
make(`Linear, XYShape.Range.stepwiseToLinear(t.xyShape), t.integralSumCache, t.integralCache);
make(~integralSumCache=t.integralSumCache, ~integralCache=t.integralCache, XYShape.Range.stepwiseToLinear(t.xyShape));
let combinePointwise =
(
@ -71,16 +71,13 @@ let combinePointwise =
let interpolator = XYShape.XtoY.continuousInterpolator(t1.interpolation, extrapolation);
make(
`Linear,
~integralSumCache=combinedIntegralSum,
XYShape.PointwiseCombination.combine(
(+.),
interpolator,
interpolator,
t1.xyShape,
t2.xyShape,
),
combinedIntegralSum,
None,
);
};
@ -89,7 +86,7 @@ let toLinear = (t: t): option(t) => {
| {interpolation: `Stepwise, xyShape, integralSumCache, integralCache} =>
xyShape
|> XYShape.Range.stepsToContinuous
|> E.O.fmap(make(`Linear, _, integralSumCache, integralCache))
|> E.O.fmap(make(~integralSumCache, ~integralCache))
| {interpolation: `Linear} => Some(t)
};
};
@ -117,13 +114,13 @@ let reduce =
let mapY = (~integralSumCacheFn=_ => None,
~integralCacheFn=_ => None,
fn, t: t) => {
let yMapFn = shapeMap(XYShape.T.mapY(fn));
t
|> yMapFn
|> updateIntegralSumCache(E.O.bind(t.integralSumCache, integralSumCacheFn))
|> updateIntegralCache(E.O.bind(t.integralCache, integralCacheFn));
~fn, t: t) => {
make(
~interpolation=t.interpolation,
~integralSumCache=t.integralSumCache |> E.O.bind(_, integralSumCacheFn),
~integralCache=t.integralCache |> E.O.bind(_, integralCacheFn),
t |> getShape |> XYShape.T.mapY(fn),
);
};
let rec scaleBy = (~scale=1.0, t: t): t => {
@ -131,7 +128,7 @@ let rec scaleBy = (~scale=1.0, t: t): t => {
let scaledIntegralCache = E.O.bind(t.integralCache, v => Some(scaleBy(~scale, v)));
t
|> mapY((r: float) => r *. scale)
|> mapY(~fn=(r: float) => r *. scale)
|> updateIntegralSumCache(scaledIntegralSumCache)
|> updateIntegralCache(scaledIntegralCache)
};
@ -177,24 +174,20 @@ module T =
let truncatedShape =
XYShape.T.fromZippedArray(truncatedZippedPairsWithNewPoints);
make(`Linear, truncatedShape, None, None);
make(truncatedShape)
};
// TODO: This should work with stepwise plots.
let integral = (t) => {
if (t |> getShape |> XYShape.T.isEmpty) {
make(`Linear, {xs: [|neg_infinity|], ys: [|0.0|]}, None, None);
} else {
switch (t.integralCache) {
| Some(cache) => cache
| None =>
t
|> getShape
|> XYShape.Range.integrateWithTriangles
|> E.O.toExt("This should not have happened")
|> make(`Linear, _, None, None)
};
};
let integral = (t) =>
switch (getShape(t) |> XYShape.T.isEmpty, t.integralCache) {
| (true, _) => emptyIntegral
| (false, Some(cache)) => cache
| (false, None) =>
t
|> getShape
|> XYShape.Range.integrateWithTriangles
|> E.O.toExt("This should not have happened")
|> make
};
let downsample = (length, t): t =>
@ -249,23 +242,17 @@ let combineAlgebraicallyWithDiscrete =
t2: DistTypes.discreteShape,
) => {
let t1s = t1 |> getShape;
let t2s = t2.xyShape; // would like to use Discrete.getShape here, but current file structure doesn't allow for that
let t2s = t2.xyShape; // TODO would like to use Discrete.getShape here, but current file structure doesn't allow for that
if (XYShape.T.isEmpty(t1s) || XYShape.T.isEmpty(t2s)) {
empty;
} else {
let shapeArray = AlgebraicShapeCombination.combineShapesContinuousDiscrete(op, t1s, t2s);
let continuousAsLinear = switch (t1.interpolation) {
| `Linear => t1;
| `Stepwise => stepwiseToLinear(t1)
};
let t1Interpolator = XYShape.XtoY.continuousInterpolator(t1.interpolation, `UseZero);
let t2Interpolator = XYShape.XtoY.discreteInterpolator;
let combinedShape =
shapeArray
|> E.A.fold_left(
XYShape.PointwiseCombination.combine((+.),
t1Interpolator,
t2Interpolator),
XYShape.T.empty);
let combinedShape = AlgebraicShapeCombination.combineShapesContinuousDiscrete(op, continuousAsLinear |> getShape, t2s);
let combinedIntegralSum =
Common.combineIntegralSums(
@ -275,7 +262,7 @@ let combineAlgebraicallyWithDiscrete =
);
// TODO: It could make sense to automatically transform the integrals here (shift or scale)
make(t1.interpolation, combinedShape, combinedIntegralSum, None)
make(~interpolation=t1.interpolation, ~integralSumCache=combinedIntegralSum, combinedShape)
};
};
@ -297,6 +284,6 @@ let combineAlgebraically =
t2.integralSumCache,
);
// return a new Continuous distribution
make(`Linear, combinedShape, combinedIntegralSum, None);
make(~integralSumCache=combinedIntegralSum, combinedShape);
};
};

119
src/distPlus/distribution/Discrete.re
View File

@ -2,7 +2,7 @@ open Distributions;
type t = DistTypes.discreteShape;
let make = (xyShape, integralSumCache, integralCache): t => {xyShape, integralSumCache, integralCache};
let make = (~integralSumCache=None, ~integralCache=None, xyShape): t => {xyShape, integralSumCache, integralCache};
let shapeMap = (fn, {xyShape, integralSumCache, integralCache}: t): t => {
xyShape: fn(xyShape),
integralSumCache,
@ -10,9 +10,21 @@ let shapeMap = (fn, {xyShape, integralSumCache, integralCache}: t): t => {
};
let getShape = (t: t) => t.xyShape;
let oShapeMap = (fn, {xyShape, integralSumCache, integralCache}: t): option(t) =>
fn(xyShape) |> E.O.fmap(make(_, integralSumCache, integralCache));
fn(xyShape) |> E.O.fmap(make(~integralSumCache, ~integralCache));
let emptyIntegral: DistTypes.continuousShape = {
xyShape: {xs: [|neg_infinity|], ys: [|0.0|]},
interpolation: `Stepwise,
integralSumCache: Some(0.0),
integralCache: None,
};
let empty: DistTypes.discreteShape = {
xyShape: XYShape.T.empty,
integralSumCache: Some(0.0),
integralCache: Some(emptyIntegral),
};
let empty: t = {xyShape: XYShape.T.empty, integralSumCache: Some(0.0), integralCache: None};
let shapeFn = (fn, t: t) => t |> getShape |> fn;
let lastY = (t: t) => t |> getShape |> XYShape.T.lastY;
@ -37,15 +49,13 @@ let combinePointwise =
// It could be done for pointwise additions, but is that ever needed?
make(
~integralSumCache=combinedIntegralSum,
XYShape.PointwiseCombination.combine(
(+.),
XYShape.XtoY.discreteInterpolator,
XYShape.XtoY.discreteInterpolator,
t1.xyShape,
t2.xyShape,
),
combinedIntegralSum,
None,
);
};
@ -100,27 +110,26 @@ let combineAlgebraically =
let combinedShape = XYShape.T.fromZippedArray(rxys);
make(combinedShape, combinedIntegralSum, None);
make(~integralSumCache=combinedIntegralSum, combinedShape);
};
let mapY = (~integralSumCacheFn=_ => None,
~integralCacheFn=_ => None,
fn, t: t) => {
let yMapFn = shapeMap(XYShape.T.mapY(fn));
t
|> yMapFn
|> updateIntegralSumCache(E.O.bind(t.integralSumCache, integralSumCacheFn))
|> updateIntegralCache(E.O.bind(t.integralCache, integralCacheFn));
~fn, t: t) => {
make(
~integralSumCache=t.integralSumCache |> E.O.bind(_, integralSumCacheFn),
~integralCache=t.integralCache |> E.O.bind(_, integralCacheFn),
t |> getShape |> XYShape.T.mapY(fn),
);
};
let scaleBy = (~scale=1.0, t: t): t => {
let scaledIntegralSumCache = E.O.bind(t.integralSumCache, v => Some(scale *. v));
let scaledIntegralCache = E.O.bind(t.integralCache, v => Some(Continuous.scaleBy(~scale, v)));
let scaledIntegralSumCache = t.integralSumCache |> E.O.fmap((*.)(scale));
let scaledIntegralCache = t.integralCache |> E.O.fmap(Continuous.scaleBy(~scale));
t
|> mapY((r: float) => r *. scale)
|> mapY(~fn=(r: float) => r *. scale)
|> updateIntegralSumCache(scaledIntegralSumCache)
|> updateIntegralCache(scaledIntegralCache)
};
@ -130,29 +139,21 @@ module T =
type t = DistTypes.discreteShape;
type integral = DistTypes.continuousShape;
let integral = (t) =>
if (t |> getShape |> XYShape.T.isEmpty) {
Continuous.make(
`Stepwise,
{xs: [|neg_infinity|], ys: [|0.0|]},
None,
None,
);
} else {
switch (t.integralCache) {
| Some(c) => c
| None => {
let ts = getShape(t);
// The first xy of this integral should always be the zero, to ensure nice plotting
let firstX = ts |> XYShape.T.minX;
let prependedZeroPoint: XYShape.T.t = {xs: [|firstX -. epsilon_float|], ys: [|0.|]};
let integralShape =
ts
|> XYShape.T.concat(prependedZeroPoint)
|> XYShape.T.accumulateYs((+.));
switch (getShape(t) |> XYShape.T.isEmpty, t.integralCache) {
| (true, _) => emptyIntegral
| (false, Some(c)) => c
| (false, None) => {
let ts = getShape(t);
// The first xy of this integral should always be the zero, to ensure nice plotting
let firstX = ts |> XYShape.T.minX;
let prependedZeroPoint: XYShape.T.t = {xs: [|firstX -. epsilon_float|], ys: [|0.|]};
let integralShape =
ts
|> XYShape.T.concat(prependedZeroPoint)
|> XYShape.T.accumulateYs((+.));
Continuous.make(`Stepwise, integralShape, None, None);
}
};
Continuous.make(~interpolation=`Stepwise, integralShape);
}
};
let integralEndY = (t: t) =>
@ -179,17 +180,15 @@ module T =
let currentLength = t |> getShape |> XYShape.T.length;
if (i < currentLength && i >= 1 && currentLength > 1) {
let clippedShape =
t
|> getShape
|> XYShape.T.zip
|> XYShape.Zipped.sortByY
|> Belt.Array.reverse
|> Belt.Array.slice(_, ~offset=0, ~len=i)
|> XYShape.Zipped.sortByX
|> XYShape.T.fromZippedArray;
make(clippedShape, None, None); // if someone needs the sum, they'll have to recompute it
t
|> getShape
|> XYShape.T.zip
|> XYShape.Zipped.sortByY
|> Belt.Array.reverse
|> Belt.Array.slice(_, ~offset=0, ~len=i)
|> XYShape.Zipped.sortByX
|> XYShape.T.fromZippedArray
|> make;
} else {
t;
};
@ -197,17 +196,15 @@ module T =
let truncate =
(leftCutoff: option(float), rightCutoff: option(float), t: t): t => {
let truncatedShape =
t
|> getShape
|> XYShape.T.zip
|> XYShape.Zipped.filterByX(x =>
x >= E.O.default(neg_infinity, leftCutoff)
&& x <= E.O.default(infinity, rightCutoff)
)
|> XYShape.T.fromZippedArray;
make(truncatedShape, None, None);
t
|> getShape
|> XYShape.T.zip
|> XYShape.Zipped.filterByX(x =>
x >= E.O.default(neg_infinity, leftCutoff)
&& x <= E.O.default(infinity, rightCutoff)
)
|> XYShape.T.fromZippedArray
|> make;
};
let xToY = (f, t) =>

4
src/distPlus/distribution/DistPlus.re
View File

@ -94,11 +94,11 @@ module T =
(
~integralSumCacheFn=previousIntegralSum => None,
~integralCacheFn=previousIntegralCache => None,
fn,
~fn,
{shape, _} as t: t,
)
: t =>
Shape.T.mapY(~integralSumCacheFn, fn, shape)
Shape.T.mapY(~integralSumCacheFn, ~fn, shape)
|> updateShape(_, t);
// get the total of everything

7
src/distPlus/distribution/DistTypes.re
View File

@ -14,11 +14,11 @@ type xyShape = {
ys: array(float),
};
type interpolation = [
type interpolationStrategy = [
| `Stepwise
| `Linear
];
type extrapolation = [
type extrapolationStrategy = [
| `UseZero
| `UseOutermostPoints
];
@ -27,14 +27,13 @@ type interpolator = (xyShape, int, float) => float;
type continuousShape = {
xyShape,
interpolation: interpolation,
interpolation: interpolationStrategy,
integralSumCache: option(float),
integralCache: option(continuousShape),
};
type discreteShape = {
xyShape,
/* interpolation is always `Discrete */
integralSumCache: option(float),
integralCache: option(continuousShape),
};

2
src/distPlus/distribution/Distributions.re
View File

@ -4,7 +4,7 @@ module type dist = {
let minX: t => float;
let maxX: t => float;
let mapY:
(~integralSumCacheFn: float => option(float)=?, ~integralCacheFn: DistTypes.continuousShape => option(DistTypes.continuousShape)=?, float => float, t) => t;
(~integralSumCacheFn: float => option(float)=?, ~integralCacheFn: DistTypes.continuousShape => option(DistTypes.continuousShape)=?, ~fn: float => float, t) => t;
let xToY: (float, t) => DistTypes.mixedPoint;
let toShape: t => DistTypes.shape;
let toContinuous: t => option(DistTypes.continuousShape);

20
src/distPlus/distribution/Mixed.re
View File

@ -1,7 +1,7 @@
open Distributions;
type t = DistTypes.mixedShape;
let make = (~continuous, ~discrete, integralSumCache, integralCache): t => {continuous, discrete, integralSumCache, integralCache};
let make = (~integralSumCache=None, ~integralCache=None, ~continuous, ~discrete): t => {continuous, discrete, integralSumCache, integralCache};
let totalLength = (t: t): int => {
let continuousLength =
@ -16,7 +16,7 @@ let scaleBy = (~scale=1.0, t: t): t => {
let scaledContinuous = Continuous.scaleBy(~scale, t.continuous);
let scaledIntegralCache = E.O.bind(t.integralCache, v => Some(Continuous.scaleBy(~scale, v)));
let scaledIntegralSumCache = E.O.bind(t.integralSumCache, s => Some(s *. scale));
make(~discrete=scaledDiscrete, ~continuous=scaledContinuous, scaledIntegralSumCache, scaledIntegralCache);
make(~discrete=scaledDiscrete, ~continuous=scaledContinuous, ~integralSumCache=scaledIntegralSumCache, ~integralCache=scaledIntegralCache);
};
let toContinuous = ({continuous}: t) => Some(continuous);
@ -54,7 +54,7 @@ module T =
let truncatedDiscrete =
Discrete.T.truncate(leftCutoff, rightCutoff, discrete);
make(~discrete=truncatedDiscrete, ~continuous=truncatedContinuous, None, None);
make(~integralSumCache=None, ~integralCache=None, ~discrete=truncatedDiscrete, ~continuous=truncatedContinuous);
};
let normalize = (t: t): t => {
@ -82,7 +82,7 @@ module T =
|> Discrete.scaleBy(~scale=newDiscreteSum /. discreteIntegralSum)
|> Discrete.updateIntegralSumCache(Some(newDiscreteSum));
make(~continuous=normalizedContinuous, ~discrete=normalizedDiscrete, Some(1.0), None);
make(~integralSumCache=Some(1.0), ~integralCache=None, ~continuous=normalizedContinuous, ~discrete=normalizedDiscrete);
};
let xToY = (x, t: t) => {
@ -144,16 +144,12 @@ module T =
let discreteIntegral = Continuous.stepwiseToLinear(Discrete.T.Integral.get(t.discrete));
Continuous.make(
`Linear,
XYShape.PointwiseCombination.combine(
(+.),
XYShape.XtoY.continuousInterpolator(`Linear, `UseOutermostPoints),
XYShape.XtoY.continuousInterpolator(`Linear, `UseOutermostPoints),
Continuous.getShape(continuousIntegral),
Continuous.getShape(discreteIntegral),
),
None,
None,
);
};
};
@ -177,19 +173,19 @@ module T =
(
~integralSumCacheFn=previousIntegralSum => None,
~integralCacheFn=previousIntegral => None,
fn,
~fn,
t: t,
)
: t => {
let yMappedDiscrete: DistTypes.discreteShape =
t.discrete
|> Discrete.T.mapY(fn)
|> Discrete.T.mapY(~fn)
|> Discrete.updateIntegralSumCache(E.O.bind(t.discrete.integralSumCache, integralSumCacheFn))
|> Discrete.updateIntegralCache(E.O.bind(t.discrete.integralCache, integralCacheFn));
let yMappedContinuous: DistTypes.continuousShape =
t.continuous
|> Continuous.T.mapY(fn)
|> Continuous.T.mapY(~fn)
|> Continuous.updateIntegralSumCache(E.O.bind(t.continuous.integralSumCache, integralSumCacheFn))
|> Continuous.updateIntegralCache(E.O.bind(t.continuous.integralCache, integralCacheFn));
@ -332,5 +328,5 @@ let combinePointwise = (~integralSumCachesFn = (_, _) => None, ~integralCachesFn
t2.integralCache,
);
make(~discrete=reducedDiscrete, ~continuous=reducedContinuous, combinedIntegralSum, combinedIntegral);
make(~integralSumCache=combinedIntegralSum, ~integralCache=combinedIntegral, ~discrete=reducedDiscrete, ~continuous=reducedContinuous);
};

8
src/distPlus/distribution/MixedShapeBuilder.re
View File

@ -9,8 +9,8 @@ type assumptions = {
};
let buildSimple = (~continuous: option(DistTypes.continuousShape), ~discrete: option(DistTypes.discreteShape)): option(DistTypes.shape) => {
let continuous = continuous |> E.O.default(Continuous.make(`Linear, {xs: [||], ys: [||]}, Some(0.0), None));
let discrete = discrete |> E.O.default(Discrete.make({xs: [||], ys: [||]}, Some(0.0), None));
let continuous = continuous |> E.O.default(Continuous.make(~integralSumCache=Some(0.0), {xs: [||], ys: [||]}));
let discrete = discrete |> E.O.default(Discrete.make(~integralSumCache=Some(0.0), {xs: [||], ys: [||]}));
let cLength =
continuous
|> Continuous.getShape
@ -24,10 +24,10 @@ let buildSimple = (~continuous: option(DistTypes.continuousShape), ~discrete: op
| (_, _) =>
let mixedDist =
Mixed.make(
~integralSumCache=None,
~integralCache=None,
~continuous,
~discrete,
None,
None,
);
Some(Mixed(mixedDist));
};

12
src/distPlus/distribution/Shape.re
View File

@ -19,8 +19,8 @@ let fmap = ((fn1, fn2, fn3), t: t): t =>
let toMixed =
mapToAll((
m => m,
d => Mixed.make(~discrete=d, ~continuous=Continuous.empty, d.integralSumCache, d.integralCache),
c => Mixed.make(~discrete=Discrete.empty, ~continuous=c, c.integralSumCache, c.integralCache),
d => Mixed.make(~integralSumCache=d.integralSumCache, ~integralCache=d.integralCache, ~discrete=d, ~continuous=Continuous.empty),
c => Mixed.make(~integralSumCache=c.integralSumCache, ~integralCache=c.integralCache, ~discrete=Discrete.empty, ~continuous=c),
));
let combineAlgebraically =
@ -176,11 +176,11 @@ module T =
));
};
let maxX = mapToAll((Mixed.T.maxX, Discrete.T.maxX, Continuous.T.maxX));
let mapY = (~integralSumCacheFn=previousIntegralSum => None, ~integralCacheFn=previousIntegral=>None, fn) =>
let mapY = (~integralSumCacheFn=previousIntegralSum => None, ~integralCacheFn=previousIntegral=>None, ~fn) =>
fmap((
Mixed.T.mapY(~integralSumCacheFn, ~integralCacheFn, fn),
Discrete.T.mapY(~integralSumCacheFn, ~integralCacheFn, fn),
Continuous.T.mapY(~integralSumCacheFn, ~integralCacheFn, fn),
Mixed.T.mapY(~integralSumCacheFn, ~integralCacheFn, ~fn),
Discrete.T.mapY(~integralSumCacheFn, ~integralCacheFn, ~fn),
Continuous.T.mapY(~integralSumCacheFn, ~integralCacheFn, ~fn),
));
let mean = (t: t): float =>

114
src/distPlus/distribution/XYShape.re
View File

@ -146,56 +146,52 @@ module XtoY = {
/* Returns a between-points-interpolating function that can be used with PointwiseCombination.combine.
Interpolation can either be stepwise (using the value on the left) or linear. Extrapolation can be `UseZero or `UseOutermostPoints. */
let continuousInterpolator = (interpolation: DistTypes.interpolation, extrapolation: DistTypes.extrapolation): interpolator => {
switch (interpolation) {
| `Linear => switch (extrapolation) {
| `UseZero => (t: T.t, leftIndex: int, x: float) => {
if (leftIndex < 0) {
0.0
} else if (leftIndex >= T.length(t) - 1) {
0.0
} else {
let x1 = t.xs[leftIndex];
let x2 = t.xs[leftIndex + 1];
let y1 = t.ys[leftIndex];
let y2 = t.ys[leftIndex + 1];
let fraction = (x -. x1) /. (x2 -. x1);
y1 *. (1. -. fraction) +. y2 *. fraction;
};
}
| `UseOutermostPoints => (t: T.t, leftIndex: int, x: float) => {
if (leftIndex < 0) {
t.ys[0];
} else if (leftIndex >= T.length(t) - 1) {
T.lastY(t);
} else {
let x1 = t.xs[leftIndex];
let x2 = t.xs[leftIndex + 1];
let y1 = t.ys[leftIndex];
let y2 = t.ys[leftIndex + 1];
let fraction = (x -. x1) /. (x2 -. x1);
y1 *. (1. -. fraction) +. y2 *. fraction;
};
let continuousInterpolator = (interpolation: DistTypes.interpolationStrategy, extrapolation: DistTypes.extrapolationStrategy): interpolator => {
switch (interpolation, extrapolation) {
| (`Linear, `UseZero) => (t: T.t, leftIndex: int, x: float) => {
if (leftIndex < 0) {
0.0
} else if (leftIndex >= T.length(t) - 1) {
0.0
} else {
let x1 = t.xs[leftIndex];
let x2 = t.xs[leftIndex + 1];
let y1 = t.ys[leftIndex];
let y2 = t.ys[leftIndex + 1];
let fraction = (x -. x1) /. (x2 -. x1);
y1 *. (1. -. fraction) +. y2 *. fraction;
};
}
| (`Linear, `UseOutermostPoints) => (t: T.t, leftIndex: int, x: float) => {
if (leftIndex < 0) {
t.ys[0];
} else if (leftIndex >= T.length(t) - 1) {
t.ys[T.length(t) - 1]
} else {
let x1 = t.xs[leftIndex];
let x2 = t.xs[leftIndex + 1];
let y1 = t.ys[leftIndex];
let y2 = t.ys[leftIndex + 1];
let fraction = (x -. x1) /. (x2 -. x1);
y1 *. (1. -. fraction) +. y2 *. fraction;
};
}
| (`Stepwise, `UseZero) => (t: T.t, leftIndex: int, x: float) => {
if (leftIndex < 0) {
0.0
} else if (leftIndex >= T.length(t) - 1) {
0.0
} else {
t.ys[leftIndex];
}
}
| `Stepwise => switch (extrapolation) {
| `UseZero => (t: T.t, leftIndex: int, x: float) => {
if (leftIndex < 0) {
0.0
} else if (leftIndex >= T.length(t) - 1) {
0.0
} else {
t.ys[leftIndex];
}
}
| `UseOutermostPoints => (t: T.t, leftIndex: int, x: float) => {
if (leftIndex < 0) {
t.ys[0];
} else if (leftIndex >= T.length(t) - 1) {
T.lastY(t);
} else {
t.ys[leftIndex];
}
| (`Stepwise, `UseOutermostPoints) => (t: T.t, leftIndex: int, x: float) => {
if (leftIndex < 0) {
t.ys[0];
} else if (leftIndex >= T.length(t) - 1) {
t.ys[T.length(t) - 1]
} else {
t.ys[leftIndex];
}
}
}
@ -246,7 +242,7 @@ module PointwiseCombination = {
// and interpolates the value on the other side, thus accumulating xs and ys.
// This is written in raw JS because this can still be a bottleneck, and using refs for the i and j indices is quite painful.
function(fn, t1Interpolator, t2Interpolator, t1, t2) {
function(fn, interpolator, t1, t2) {
let t1n = t1.xs.length;
let t2n = t2.xs.length;
let outX = [];
@ -263,7 +259,7 @@ module PointwiseCombination = {
x = t1.xs[i];
ya = t1.ys[i];
yb = t2Interpolator(t2, j, x);
yb = interpolator(t2, j, x);
} else if (i == t1n - 1 && j < t2n - 1 ||
t1.xs[i+1] > t2.xs[j+1]) { // if b has to catch up to a, or if a is already done
j++;
@ -271,7 +267,7 @@ module PointwiseCombination = {
x = t2.xs[j];
yb = t2.ys[j];
ya = t1Interpolator(t1, i, x);
ya = interpolator(t1, i, x);
} else if (i < t1n - 1 && j < t2n && t1.xs[i+1] === t2.xs[j+1]) { // if they happen to be equal, move both ahead
i++;
j++;
@ -384,16 +380,16 @@ module Range = {
let newXs: array(float) = Belt.Array.makeUninitializedUnsafe(2 * length);
let newYs: array(float) = Belt.Array.makeUninitializedUnsafe(2 * length);
let _ = Belt.Array.set(newXs, 0, xs[0] -. epsilon_float);
let _ = Belt.Array.set(newYs, 0, 0.);
let _ = Belt.Array.set(newXs, 1, xs[0]);
let _ = Belt.Array.set(newYs, 1, ys[0]);
Belt.Array.set(newXs, 0, xs[0] -. epsilon_float) |> ignore;
Belt.Array.set(newYs, 0, 0.) |> ignore;
Belt.Array.set(newXs, 1, xs[0]) |> ignore;
Belt.Array.set(newYs, 1, ys[0]) |> ignore;
for (i in 1 to E.A.length(xs) - 1) {
let _ = Belt.Array.set(newXs, i * 2, xs[i] -. epsilon_float);
let _ = Belt.Array.set(newYs, i * 2, ys[i-1]);
let _ = Belt.Array.set(newXs, i * 2 + 1, xs[i]);
let _ = Belt.Array.set(newYs, i * 2 + 1, ys[i]);
Belt.Array.set(newXs, i * 2, xs[i] -. epsilon_float) |> ignore;
Belt.Array.set(newYs, i * 2, ys[i-1]) |> ignore;
Belt.Array.set(newXs, i * 2 + 1, xs[i]) |> ignore;
Belt.Array.set(newYs, i * 2 + 1, ys[i]) |> ignore;
();
};

2
src/distPlus/expressionTree/ExpressionTreeEvaluator.re
View File

@ -112,7 +112,7 @@ module VerticalScaling = {
Shape.T.mapY(
~integralSumCacheFn=integralSumCacheFn(sm),
~integralCacheFn=integralCacheFn(sm),
fn(sm),
~fn=fn(sm),
rs,
),
),

4
src/distPlus/renderers/samplesRenderer/Samples.re
View File

@ -120,7 +120,7 @@ module T = {
|> E.FloatFloatMap.fmap(r => r /. length)
|> E.FloatFloatMap.toArray
|> XYShape.T.fromZippedArray
|> Discrete.make(_, None, None);
|> Discrete.make;
let pdf =
continuousPart |> E.A.length > 5
@ -150,7 +150,7 @@ module T = {
~outputXYPoints=samplingInputs.outputXYPoints,
formatUnitWidth(usedUnitWidth),
)
|> Continuous.make(`Linear, _, None, None)
|> Continuous.make
|> (r => Some((r, foo)));
}
: None;

4
src/distPlus/symbolic/SymbolicDist.re
View File

@ -317,13 +317,13 @@ module T = {
switch (d) {
| `Float(v) =>
Discrete(
Discrete.make({xs: [|v|], ys: [|1.0|]}, Some(1.0), None),
Discrete.make(~integralSumCache=Some(1.0), {xs: [|v|], ys: [|1.0|]}),
)
| _ =>
let xs = interpolateXs(~xSelection=`ByWeight, d, sampleCount);
let ys = xs |> E.A.fmap(x => pdf(x, d));
Continuous(
Continuous.make(`Linear, {xs, ys}, Some(1.0), None),
Continuous.make(~integralSumCache=Some(1.0), {xs, ys}),
);
};
};