Move integralCache into shapes; clean up interpolator functions for PointwiseCombinations

src/components/DistBuilder.re

@@ -147,7 +147,10 @@ module DemoDist = {
              );
            let response = DistPlusRenderer.run(inputs);
            switch (RenderTypes.DistPlusRenderer.Outputs.distplus(response)) {
-           | Some(distPlus) => <DistPlusPlot distPlus />
+           | Some(distPlus) => {
+             let normalizedDistPlus = DistPlus.T.normalize(distPlus);
+             <DistPlusPlot distPlus={normalizedDistPlus} />;
+           }
            | _ =>
              "Correct Guesstimator string input to show a distribution."
              |> R.ste

src/components/DistBuilder2.re

@@ -44,7 +44,7 @@ module DemoDist = {
             DistPlus.make(
               ~shape=
                 Continuous(
-                  Continuous.make(`Linear, {xs, ys}, None),
+                  Continuous.make(`Linear, {xs, ys}, None, None),
                 ),
               ~domain=Complete,
               ~unit=UnspecifiedDistribution,

src/components/Drawer.re

@@ -177,7 +177,8 @@ module Convert = {
     let continuousShape: Types.continuousShape = {
       xyShape,
       interpolation: `Linear,
-      knownIntegralSum: None,
+      integralSumCache: None,
+      integralCache: None,
     };
     
     let integral = XYShape.Analysis.integrateContinuousShape(continuousShape);
@@ -189,7 +190,8 @@ module Convert = {
         ys,
       },
       interpolation: `Linear,
-      knownIntegralSum: Some(1.0),
+      integralSumCache: Some(1.0),
+      integralCache: None,
     };
     continuousShape;
   };
@@ -673,7 +675,7 @@ module State = {
           pdf,
         );
 
-      let cdf = Continuous.T.integral(~cache=None, _pdf);
+      let cdf = Continuous.T.integral(_pdf);
       let xs = [||];
       let ys = [||];
       for (i in 1 to 999) {

src/components/charts/DistPlusPlot.re

@@ -51,13 +51,13 @@ let table = (distPlus, x) => {
           </td>
           <td className="px-4 py-2 border ">
             {distPlus
-             |> DistPlus.T.Integral.xToY(~cache=None, x)
+             |> DistPlus.T.Integral.xToY(x)
              |> E.Float.with2DigitsPrecision
              |> ReasonReact.string}
           </td>
           <td className="px-4 py-2 border ">
             {distPlus
-             |> DistPlus.T.Integral.sum(~cache=None)
+             |> DistPlus.T.Integral.sum
              |> E.Float.with2DigitsPrecision
              |> ReasonReact.string}
           </td>
@@ -70,15 +70,9 @@ let table = (distPlus, x) => {
           <td className="px-4 py-2">
             {"Continuous Total" |> ReasonReact.string}
           </td>
-          <td className="px-4 py-2">
-            {"Scaled Continuous Total" |> ReasonReact.string}
-          </td>
           <td className="px-4 py-2">
             {"Discrete Total" |> ReasonReact.string}
           </td>
-          <td className="px-4 py-2">
-            {"Scaled Discrete Total" |> ReasonReact.string}
-          </td>
         </tr>
       </thead>
       <tbody>
@@ -87,17 +81,7 @@ let table = (distPlus, x) => {
             {distPlus
              |> DistPlus.T.toContinuous
              |> E.O.fmap(
-                  Continuous.T.Integral.sum(~cache=None),
+                  Continuous.T.Integral.sum
-                )
-             |> E.O.fmap(E.Float.with2DigitsPrecision)
-             |> E.O.default("")
-             |> ReasonReact.string}
-          </td>
-          <td className="px-4 py-2 border ">
-            {distPlus
-             |> DistPlus.T.normalizedToContinuous
-             |> E.O.fmap(
-                  Continuous.T.Integral.sum(~cache=None),
                 )
              |> E.O.fmap(E.Float.with2DigitsPrecision)
              |> E.O.default("")
@@ -106,15 +90,7 @@ let table = (distPlus, x) => {
           <td className="px-4 py-2 border ">
             {distPlus
              |> DistPlus.T.toDiscrete
-             |> E.O.fmap(Discrete.T.Integral.sum(~cache=None))
+             |> E.O.fmap(Discrete.T.Integral.sum)
-             |> E.O.fmap(E.Float.with2DigitsPrecision)
-             |> E.O.default("")
-             |> ReasonReact.string}
-          </td>
-          <td className="px-4 py-2 border ">
-            {distPlus
-             |> DistPlus.T.normalizedToDiscrete
-             |> E.O.fmap(Discrete.T.Integral.sum(~cache=None))
              |> E.O.fmap(E.Float.with2DigitsPrecision)
              |> E.O.default("")
              |> ReasonReact.string}
@@ -143,42 +119,42 @@ let percentiles = distPlus => {
         <tr>
           <td className="px-4 py-2 border">
             {distPlus
-             |> DistPlus.T.Integral.yToX(~cache=None, 0.01)
+             |> DistPlus.T.Integral.yToX(0.01)
              |> showFloat}
           </td>
           <td className="px-4 py-2 border">
             {distPlus
-             |> DistPlus.T.Integral.yToX(~cache=None, 0.05)
+             |> DistPlus.T.Integral.yToX(0.05)
              |> showFloat}
           </td>
           <td className="px-4 py-2 border">
             {distPlus
-             |> DistPlus.T.Integral.yToX(~cache=None, 0.25)
+             |> DistPlus.T.Integral.yToX(0.25)
              |> showFloat}
           </td>
           <td className="px-4 py-2 border">
             {distPlus
-             |> DistPlus.T.Integral.yToX(~cache=None, 0.5)
+             |> DistPlus.T.Integral.yToX(0.5)
              |> showFloat}
           </td>
           <td className="px-4 py-2 border">
             {distPlus
-             |> DistPlus.T.Integral.yToX(~cache=None, 0.75)
+             |> DistPlus.T.Integral.yToX(0.75)
              |> showFloat}
           </td>
           <td className="px-4 py-2 border">
             {distPlus
-             |> DistPlus.T.Integral.yToX(~cache=None, 0.95)
+             |> DistPlus.T.Integral.yToX(0.95)
              |> showFloat}
           </td>
           <td className="px-4 py-2 border">
             {distPlus
-             |> DistPlus.T.Integral.yToX(~cache=None, 0.99)
+             |> DistPlus.T.Integral.yToX(0.99)
              |> showFloat}
           </td>
           <td className="px-4 py-2 border">
             {distPlus
-             |> DistPlus.T.Integral.yToX(~cache=None, 0.99999)
+             |> DistPlus.T.Integral.yToX(0.99999)
              |> showFloat}
           </td>
         </tr>
@@ -225,10 +201,11 @@ module DistPlusChart = {
   [@react.component]
   let make = (~distPlus: DistTypes.distPlus, ~config: chartConfig, ~onHover) => {
     open DistPlus;
-    let discrete = distPlus |> T.normalizedToDiscrete |> E.O.fmap(Discrete.getShape);
+
+    let discrete = distPlus |> T.toDiscrete |> E.O.fmap(Discrete.getShape);
     let continuous =
       distPlus
-      |> T.normalizedToContinuous
+      |> T.toContinuous
       |> E.O.fmap(Continuous.getShape);
     let range = T.xTotalRange(distPlus);
 
@@ -236,7 +213,7 @@ module DistPlusChart = {
     // let minX =
     //   switch (
     //     distPlus
-    //     |> DistPlus.T.Integral.yToX(~cache=None, 0.0001),
+    //     |> DistPlus.T.Integral.yToX(0.0001),
     //     range,
     //   ) {
     //   | (min, Some(range)) => Some(min -. range *. 0.001)
@@ -244,11 +221,11 @@ module DistPlusChart = {
     //   };
 
     let minX = {
-      distPlus |> DistPlus.T.Integral.yToX(~cache=None, 0.00001);
+      distPlus |> DistPlus.T.Integral.yToX(0.00001);
     };
 
     let maxX = {
-      distPlus |> DistPlus.T.Integral.yToX(~cache=None, 0.99999);
+      distPlus |> DistPlus.T.Integral.yToX(0.99999);
     };
 
     let timeScale = distPlus.unit |> DistTypes.DistributionUnit.toJson;
@@ -283,11 +260,11 @@ module IntegralChart = {
       |> Continuous.toLinear
       |> E.O.fmap(Continuous.getShape);
     let minX = {
-      distPlus |> DistPlus.T.Integral.yToX(~cache=None, 0.00001);
+      distPlus |> DistPlus.T.Integral.yToX(0.00001);
     };
 
     let maxX = {
-      distPlus |> DistPlus.T.Integral.yToX(~cache=None, 0.99999);
+      distPlus |> DistPlus.T.Integral.yToX(0.99999);
     };
     let timeScale = distPlus.unit |> DistTypes.DistributionUnit.toJson;
     <DistributionPlot
@@ -334,6 +311,7 @@ let make = (~distPlus: DistTypes.distPlus) => {
   let (x, setX) = React.useState(() => 0.);
   let (state, dispatch) =
     React.useReducer(DistPlusPlotReducer.reducer, DistPlusPlotReducer.init);
+
   <div>
     {state.distributions
      |> E.L.fmapi((index, config) =>

src/distPlus/distribution/AlgebraicShapeCombination.re

@@ -247,7 +247,7 @@ let toDiscretePointMassesFromDiscrete = (s: DistTypes.xyShape): pointMassesWithM
 
 let combineShapesContinuousDiscrete =
     (op: ExpressionTypes.algebraicOperation, s1: DistTypes.xyShape, s2: DistTypes.xyShape)
-    : DistTypes.xyShape => {
+    : array(DistTypes.xyShape) => {
 
   let t1n = s1 |> XYShape.T.length;
   let t2n = s2 |> XYShape.T.length;
@@ -303,10 +303,5 @@ let combineShapesContinuousDiscrete =
   };
 
   outXYShapes
-  |> E.A.fmap(XYShape.T.fromZippedArray)
+  |> E.A.fmap(XYShape.T.fromZippedArray);
-  |> E.A.fold_left(
-      XYShape.PointwiseCombination.combine((+.),
-                                           XYShape.XtoY.linearBetweenPointsExtrapolateZero,
-                                           XYShape.XtoY.linearBetweenPointsExtrapolateZero),
-      XYShape.T.empty);
 };

src/distPlus/distribution/Continuous.re

@@ -3,30 +3,45 @@ open Distributions;
 type t = DistTypes.continuousShape;
 let getShape = (t: t) => t.xyShape;
 let interpolation = (t: t) => t.interpolation;
-let make = (interpolation, xyShape, knownIntegralSum): t => {
+let make = (interpolation, xyShape, integralSumCache, integralCache): t => {
   xyShape,
   interpolation,
-  knownIntegralSum,
+  integralSumCache,
+  integralCache,
 };
-let shapeMap = (fn, {xyShape, interpolation, knownIntegralSum}: t): t => {
+let shapeMap = (fn, {xyShape, interpolation, integralSumCache, integralCache}: t): t => {
   xyShape: fn(xyShape),
   interpolation,
-  knownIntegralSum,
+  integralSumCache,
+  integralCache,
 };
 let lastY = (t: t) => t |> getShape |> XYShape.T.lastY;
 let oShapeMap =
-    (fn, {xyShape, interpolation, knownIntegralSum}: t)
+    (fn, {xyShape, interpolation, integralSumCache, integralCache}: t)
     : option(DistTypes.continuousShape) =>
-  fn(xyShape) |> E.O.fmap(make(interpolation, _, knownIntegralSum));
+  fn(xyShape) |> E.O.fmap(make(interpolation, _, integralSumCache, integralCache));
 
+let emptyIntegral: DistTypes.continuousShape = {
+  xyShape: {xs: [|neg_infinity|], ys: [|0.0|]},
+  interpolation: `Linear,
+  integralSumCache: Some(0.0),
+  integralCache: None,
+};
 let empty: DistTypes.continuousShape = {
   xyShape: XYShape.T.empty,
   interpolation: `Linear,
-  knownIntegralSum: Some(0.0),
+  integralSumCache: Some(0.0),
+  integralCache: Some(emptyIntegral),
 };
+
+let stepwiseToLinear = (t: t): t =>
+  make(`Linear, XYShape.Range.stepwiseToLinear(t.xyShape), t.integralSumCache, t.integralCache);
+
 let combinePointwise =
     (
-      ~knownIntegralSumsFn,
+      ~integralSumCachesFn=(_, _) => None,
+      ~integralCachesFn: (t, t) => option(t) =(_, _) => None,
+      ~extrapolation=`UseZero,
       fn: (float, float) => float,
       t1: DistTypes.continuousShape,
       t2: DistTypes.continuousShape,
@@ -36,61 +51,89 @@ let combinePointwise =
   // can just sum them up. Otherwise, all bets are off.
   let combinedIntegralSum =
     Common.combineIntegralSums(
-      knownIntegralSumsFn,
+      integralSumCachesFn,
-      t1.knownIntegralSum,
+      t1.integralSumCache,
-      t2.knownIntegralSum,
+      t2.integralSumCache,
     );
 
+  // TODO: does it ever make sense to pointwise combine the integrals here?
+  // It could be done for pointwise additions, but is that ever needed?
+
+  // If combining stepwise and linear, we must convert the stepwise to linear first,
+  // i.e. add a point at the bottom of each step
+  let (t1, t2) = switch (t1.interpolation, t2.interpolation) {
+  | (`Linear, `Linear) => (t1, t2);
+  | (`Stepwise, `Stepwise) => (t1, t2);
+  | (`Linear, `Stepwise) => (t1, stepwiseToLinear(t2));
+  | (`Stepwise, `Linear) => (stepwiseToLinear(t1), t2);
+  };
+
+  let interpolator = XYShape.XtoY.continuousInterpolator(t1.interpolation, extrapolation);
+
   make(
     `Linear,
     XYShape.PointwiseCombination.combine(
       (+.),
-      XYShape.XtoY.linearBetweenPointsExtrapolateZero,
+      interpolator,
-      XYShape.XtoY.linearBetweenPointsExtrapolateZero,
+      interpolator,
       t1.xyShape,
       t2.xyShape,
     ),
     combinedIntegralSum,
+    None,
   );
 };
 
 let toLinear = (t: t): option(t) => {
   switch (t) {
-  | {interpolation: `Stepwise, xyShape, knownIntegralSum} =>
+  | {interpolation: `Stepwise, xyShape, integralSumCache, integralCache} =>
     xyShape
     |> XYShape.Range.stepsToContinuous
-    |> E.O.fmap(make(`Linear, _, knownIntegralSum))
+    |> E.O.fmap(make(`Linear, _, integralSumCache, integralCache))
   | {interpolation: `Linear} => Some(t)
   };
 };
 let shapeFn = (fn, t: t) => t |> getShape |> fn;
-let updateKnownIntegralSum = (knownIntegralSum, t: t): t => {
+
+let updateIntegralSumCache = (integralSumCache, t: t): t => {
   ...t,
-  knownIntegralSum,
+  integralSumCache,
+};
+
+let updateIntegralCache = (integralCache, t: t): t => {
+  ...t,
+  integralCache,
 };
 
 let reduce =
     (
-      ~knownIntegralSumsFn: (float, float) => option(float)=(_, _) => None,
+      ~integralSumCachesFn: (float, float) => option(float)=(_, _) => None,
+      ~integralCachesFn: (t, t) => option(t)=(_, _) => None,
       fn,
       continuousShapes,
     ) =>
   continuousShapes
-  |> E.A.fold_left(combinePointwise(~knownIntegralSumsFn, fn), empty);
+  |> E.A.fold_left(combinePointwise(~integralSumCachesFn, ~integralCachesFn, fn), empty);
 
-let mapY = (~knownIntegralSumFn=_ => None, fn, t: t) => {
+let mapY = (~integralSumCacheFn=_ => None,
-  let u = E.O.bind(_, knownIntegralSumFn);
+            ~integralCacheFn=_ => None,
+            fn, t: t) => {
   let yMapFn = shapeMap(XYShape.T.mapY(fn));
 
-  t |> yMapFn |> updateKnownIntegralSum(u(t.knownIntegralSum));
+  t
+  |> yMapFn
+  |> updateIntegralSumCache(E.O.bind(t.integralSumCache, integralSumCacheFn))
+  |> updateIntegralCache(E.O.bind(t.integralCache, integralCacheFn));
 };
 
-let scaleBy = (~scale=1.0, t: t): t => {
+let rec 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(scaleBy(~scale, v)));
+
   t
   |> mapY((r: float) => r *. scale)
-  |> updateKnownIntegralSum(
+  |> updateIntegralSumCache(scaledIntegralSumCache)
-       E.O.bind(t.knownIntegralSum, v => Some(scale *. v)),
+  |> updateIntegralCache(scaledIntegralCache)
-     );
 };
 
 module T =
@@ -135,46 +178,48 @@ module T =
       let truncatedShape =
         XYShape.T.fromZippedArray(truncatedZippedPairsWithNewPoints);
 
-      make(`Linear, truncatedShape, None);
+      make(`Linear, truncatedShape, None, None);
     };
 
     // TODO: This should work with stepwise plots.
-    let integral = (~cache, t) =>
+    let integral = (t) => {
-      if (t |> getShape |> XYShape.T.length > 0) {
+      if (t |> getShape |> XYShape.T.isEmpty) {
-        switch (cache) {
+        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)
+          |> make(`Linear, _, None, None)
         };
-      } else {
+      };
-        make(`Linear, {xs: [|neg_infinity|], ys: [|0.0|]}, None);
       };
 
-    let downsample = (~cache=None, length, t): t =>
+    let downsample = (length, t): t =>
       t
       |> shapeMap(
            XYShape.XsConversion.proportionByProbabilityMass(
              length,
-             integral(~cache, t).xyShape,
+             integral(t).xyShape,
            ),
          );
-    let integralEndY = (~cache, t: t) =>
+    let integralEndY = (t: t) =>
-      t.knownIntegralSum |> E.O.default(t |> integral(~cache) |> lastY);
+      t.integralSumCache |> E.O.default(t |> integral |> lastY);
-    let integralXtoY = (~cache, f, t: t) =>
+    let integralXtoY = (f, t: t) =>
-      t |> integral(~cache) |> shapeFn(XYShape.XtoY.linear(f));
+      t |> integral |> shapeFn(XYShape.XtoY.linear(f));
-    let integralYtoX = (~cache, f, t: t) =>
+    let integralYtoX = (f, t: t) =>
-      t |> integral(~cache) |> shapeFn(XYShape.YtoX.linear(f));
+      t |> integral |> shapeFn(XYShape.YtoX.linear(f));
     let toContinuous = t => Some(t);
     let toDiscrete = _ => None;
 
     let normalize = (t: t): t => {
       t
-      |> scaleBy(~scale=1. /. integralEndY(~cache=None, t))
+      |> updateIntegralCache(Some(integral(t)))
-      |> updateKnownIntegralSum(Some(1.0));
+      |> scaleBy(~scale=1. /. integralEndY(t))
+      |> updateIntegralSumCache(Some(1.0));
     };
 
     let normalizedToContinuous = t => Some(t |> normalize);
@@ -203,7 +248,6 @@ module T =
    each discrete data point, and then adds them all together. */
 let combineAlgebraicallyWithDiscrete =
     (
-      ~downsample=false,
       op: ExpressionTypes.algebraicOperation,
       t1: t,
       t2: DistTypes.discreteShape,
@@ -211,27 +255,35 @@ let combineAlgebraicallyWithDiscrete =
   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 t1n = t1s |> XYShape.T.length;
+  if (XYShape.T.isEmpty(t1s) || XYShape.T.isEmpty(t2s)) {
-  let t2n = t2s |> XYShape.T.length;
+    empty;
+  } else {
+    let shapeArray = AlgebraicShapeCombination.combineShapesContinuousDiscrete(op, t1s, t2s);
 
-  if (t1n > 0 && t2n > 0) {
+    let t1Interpolator = XYShape.XtoY.continuousInterpolator(t1.interpolation, `UseZero);
-    let combinedShape = AlgebraicShapeCombination.combineShapesContinuousDiscrete(op, t1s, t2s);
+    let t2Interpolator = XYShape.XtoY.discreteInterpolator;
+
+    let combinedShape =
+      shapeArray
+      |> E.A.fold_left(
+          XYShape.PointwiseCombination.combine((+.),
+                                              t1Interpolator,
+                                              t2Interpolator),
+          XYShape.T.empty);
 
     let combinedIntegralSum =
       Common.combineIntegralSums(
         (a, b) => Some(a *. b),
-        t1.knownIntegralSum,
+        t1.integralSumCache,
-        t2.knownIntegralSum,
+        t2.integralSumCache,
       );
 
-    make(`Linear, combinedShape, combinedIntegralSum);
+    // TODO: It could make sense to automatically transform the integrals here (shift or scale)
-  } else {
+    make(t1.interpolation, combinedShape, combinedIntegralSum, None)
-    empty;
+  };
-  }
 };
 
-let combineAlgebraically =
+let combineAlgebraically = (op: ExpressionTypes.algebraicOperation, t1: t, t2: t) => {
-    (~downsample=false, op: ExpressionTypes.algebraicOperation, t1: t, t2: t) => {
   let s1 = t1 |> getShape;
   let s2 = t2 |> getShape;
   let t1n = s1 |> XYShape.T.length;
@@ -244,10 +296,10 @@ let combineAlgebraically =
     let combinedIntegralSum =
       Common.combineIntegralSums(
         (a, b) => Some(a *. b),
-        t1.knownIntegralSum,
+        t1.integralSumCache,
-        t2.knownIntegralSum,
+        t2.integralSumCache,
       );
     // return a new Continuous distribution
-    make(`Linear, combinedShape, combinedIntegralSum);
+    make(`Linear, combinedShape, combinedIntegralSum, None);
   };
 };

src/distPlus/distribution/Discrete.re

@@ -2,23 +2,25 @@ open Distributions;
 
 type t = DistTypes.discreteShape;
 
-let make = (xyShape, knownIntegralSum): t => {xyShape, knownIntegralSum};
+let make = (xyShape, integralSumCache, integralCache): t => {xyShape, integralSumCache, integralCache};
-let shapeMap = (fn, {xyShape, knownIntegralSum}: t): t => {
+let shapeMap = (fn, {xyShape, integralSumCache, integralCache}: t): t => {
   xyShape: fn(xyShape),
-  knownIntegralSum,
+  integralSumCache,
+  integralCache
 };
 let getShape = (t: t) => t.xyShape;
-let oShapeMap = (fn, {xyShape, knownIntegralSum}: t): option(t) =>
+let oShapeMap = (fn, {xyShape, integralSumCache, integralCache}: t): option(t) =>
-  fn(xyShape) |> E.O.fmap(make(_, knownIntegralSum));
+  fn(xyShape) |> E.O.fmap(make(_, integralSumCache, integralCache));
 
-let empty: t = {xyShape: XYShape.T.empty, knownIntegralSum: Some(0.0)};
+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;
 
 let combinePointwise =
     (
-      ~knownIntegralSumsFn,
+      ~integralSumCachesFn = (_, _) => None,
+      ~integralCachesFn: (DistTypes.continuousShape, DistTypes.continuousShape) => option(DistTypes.continuousShape) = (_, _) => None,
       fn,
       t1: DistTypes.discreteShape,
       t2: DistTypes.discreteShape,
@@ -26,38 +28,49 @@ let combinePointwise =
     : DistTypes.discreteShape => {
   let combinedIntegralSum =
     Common.combineIntegralSums(
-      knownIntegralSumsFn,
+      integralSumCachesFn,
-      t1.knownIntegralSum,
+      t1.integralSumCache,
-      t2.knownIntegralSum,
+      t2.integralSumCache,
     );
 
+  // TODO: does it ever make sense to pointwise combine the integrals here?
+  // It could be done for pointwise additions, but is that ever needed?
+
   make(
     XYShape.PointwiseCombination.combine(
       (+.),
-      XYShape.XtoY.assumeZeroBetweenPoints,
+      XYShape.XtoY.discreteInterpolator,
-      XYShape.XtoY.assumeZeroBetweenPoints,
+      XYShape.XtoY.discreteInterpolator,
       t1.xyShape,
       t2.xyShape,
     ),
     combinedIntegralSum,
+    None,
   );
 };
 
 let reduce =
-    (~knownIntegralSumsFn=(_, _) => None, fn, discreteShapes)
+  (~integralSumCachesFn=(_, _) => None,
-    : DistTypes.discreteShape =>
+   ~integralCachesFn=(_, _) => None,
+   fn, discreteShapes)
+  : DistTypes.discreteShape =>
   discreteShapes
-  |> E.A.fold_left(combinePointwise(~knownIntegralSumsFn, fn), empty);
+  |> E.A.fold_left(combinePointwise(~integralSumCachesFn, ~integralCachesFn, fn), empty);
 
-let updateKnownIntegralSum = (knownIntegralSum, t: t): t => {
+let updateIntegralSumCache = (integralSumCache, t: t): t => {
   ...t,
-  knownIntegralSum,
+  integralSumCache,
+};
+
+let updateIntegralCache = (integralCache, t: t): t => {
+  ...t,
+  integralCache,
 };
 
 /* This multiples all of the data points together and creates a new discrete distribution from the results.
    Data points at the same xs get added together. It may be a good idea to downsample t1 and t2 before and/or the result after. */
 let combineAlgebraically =
-    (op: ExpressionTypes.algebraicOperation, t1: t, t2: t) => {
+    (op: ExpressionTypes.algebraicOperation, t1: t, t2: t): t => {
   let t1s = t1 |> getShape;
   let t2s = t2 |> getShape;
   let t1n = t1s |> XYShape.T.length;
@@ -66,8 +79,8 @@ let combineAlgebraically =
   let combinedIntegralSum =
     Common.combineIntegralSums(
       (s1, s2) => Some(s1 *. s2),
-      t1.knownIntegralSum,
+      t1.integralSumCache,
-      t2.knownIntegralSum,
+      t2.integralSumCache,
     );
 
   let fn = Operation.Algebraic.toFn(op);
@@ -87,31 +100,45 @@ let combineAlgebraically =
 
   let combinedShape = XYShape.T.fromZippedArray(rxys);
 
-  make(combinedShape, combinedIntegralSum);
+  make(combinedShape, combinedIntegralSum, None);
 };
 
-let mapY = (~knownIntegralSumFn=previousKnownIntegralSum => None, fn, t: t) => {
+let mapY = (~integralSumCacheFn=_ => None,
-  let u = E.O.bind(_, knownIntegralSumFn);
+            ~integralCacheFn=_ => None,
+            fn, t: t) => {
   let yMapFn = shapeMap(XYShape.T.mapY(fn));
 
-  t |> yMapFn |> updateKnownIntegralSum(u(t.knownIntegralSum));
+  t
+  |> yMapFn
+  |> updateIntegralSumCache(E.O.bind(t.integralSumCache, integralSumCacheFn))
+  |> updateIntegralCache(E.O.bind(t.integralCache, integralCacheFn));
 };
 
+
 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)));
+
   t
   |> mapY((r: float) => r *. scale)
-  |> updateKnownIntegralSum(
+  |> updateIntegralSumCache(scaledIntegralSumCache)
-       E.O.bind(t.knownIntegralSum, v => Some(scale *. v)),
+  |> updateIntegralCache(scaledIntegralCache)
-     );
 };
 
 module T =
   Dist({
     type t = DistTypes.discreteShape;
     type integral = DistTypes.continuousShape;
-    let integral = (~cache, t) =>
+    let integral = (t) =>
-      if (t |> getShape |> XYShape.T.length > 0) {
+      if (t |> getShape |> XYShape.T.isEmpty) {
-        switch (cache) {
+        Continuous.make(
+          `Stepwise,
+          {xs: [|neg_infinity|], ys: [|0.0|]},
+          None,
+          None,
+        );
+      } else {
+        switch (t.integralCache) {
         | Some(c) => c
         | None => {
             let ts = getShape(t);
@@ -123,20 +150,14 @@ module T =
               |> XYShape.T.concat(prependedZeroPoint)
               |> XYShape.T.accumulateYs((+.));
 
-            Continuous.make(`Stepwise, integralShape, None);
+            Continuous.make(`Stepwise, integralShape, None, None);
           }
         };
-      } else {
-        Continuous.make(
-          `Stepwise,
-          {xs: [|neg_infinity|], ys: [|0.0|]},
-          None,
-        );
       };
 
-    let integralEndY = (~cache, t: t) =>
+    let integralEndY = (t: t) =>
-      t.knownIntegralSum
+      t.integralSumCache
-      |> E.O.default(t |> integral(~cache) |> Continuous.lastY);
+      |> E.O.default(t |> integral |> Continuous.lastY);
     let minX = shapeFn(XYShape.T.minX);
     let maxX = shapeFn(XYShape.T.maxX);
     let toDiscreteProbabilityMassFraction = _ => 1.0;
@@ -147,14 +168,14 @@ module T =
 
     let normalize = (t: t): t => {
       t
-      |> scaleBy(~scale=1. /. integralEndY(~cache=None, t))
+      |> scaleBy(~scale=1. /. integralEndY(t))
-      |> updateKnownIntegralSum(Some(1.0));
+      |> updateIntegralSumCache(Some(1.0));
     };
 
     let normalizedToContinuous = _ => None;
     let normalizedToDiscrete = t => Some(t); // TODO: this should be normalized!
 
-    let downsample = (~cache=None, i, t: t): t => {
+    let downsample = (i, t: t): t => {
       // It's not clear how to downsample a set of discrete points in a meaningful way.
       // The best we can do is to clip off the smallest values.
       let currentLength = t |> getShape |> XYShape.T.length;
@@ -170,7 +191,7 @@ module T =
           |> XYShape.Zipped.sortByX
           |> XYShape.T.fromZippedArray;
 
-        make(clippedShape, None); // if someone needs the sum, they'll have to recompute it
+        make(clippedShape, None, None); // if someone needs the sum, they'll have to recompute it
       } else {
         t;
       };
@@ -188,7 +209,7 @@ module T =
            )
         |> XYShape.T.fromZippedArray;
 
-      make(truncatedShape, None);
+      make(truncatedShape, None, None);
     };
 
     let xToY = (f, t) =>
@@ -198,11 +219,11 @@ module T =
       |> E.O.default(0.0)
       |> DistTypes.MixedPoint.makeDiscrete;
 
-    let integralXtoY = (~cache, f, t) =>
+    let integralXtoY = (f, t) =>
-      t |> integral(~cache) |> Continuous.getShape |> XYShape.XtoY.linear(f);
+      t |> integral |> Continuous.getShape |> XYShape.XtoY.linear(f);
 
-    let integralYtoX = (~cache, f, t) =>
+    let integralYtoX = (f, t) =>
-      t |> integral(~cache) |> Continuous.getShape |> XYShape.YtoX.linear(f);
+      t |> integral |> Continuous.getShape |> XYShape.YtoX.linear(f);
 
     let mean = (t: t): float => {
       let s = getShape(t);

src/distPlus/distribution/DistPlus.re

@@ -2,8 +2,7 @@ open DistTypes;
 
 type t = DistTypes.distPlus;
 
-let shapeIntegral = shape =>
+let shapeIntegral = shape => Shape.T.Integral.get(shape);
-  Shape.T.Integral.get(~cache=None, shape);
 let make =
     (
       ~shape,
@@ -34,6 +33,7 @@ let update =
 };
 
 let updateShape = (shape, t) => {
+  Js.log("Updating the shape, recalculating the integral");
   let integralCache = shapeIntegral(shape);
   update(~shape, ~integralCache, t);
 };
@@ -59,7 +59,6 @@ module T =
     let normalize = (t: t): t => {
       let normalizedShape = t |> toShape |> Shape.T.normalize;
       t |> updateShape(normalizedShape);
-      // TODO: also adjust for domainIncludedProbabilityMass here.
     };
 
     let truncate = (leftCutoff, rightCutoff, t: t): t => {
@@ -71,6 +70,7 @@ module T =
       t |> updateShape(truncatedShape);
     };
 
+    // TODO: is this still needed?
     let normalizedToContinuous = (t: t) => {
       t
       |> toShape
@@ -82,6 +82,7 @@ module T =
          );
     };
 
+    // TODO: is this still needed?
     let normalizedToDiscrete = (t: t) => {
       t
       |> toShape
@@ -105,34 +106,33 @@ module T =
       shapeFn(Shape.T.toDiscreteProbabilityMassFraction);
 
     // This bit is kind of awkward, could probably use rethinking.
-    let integral = (~cache, t: t) =>
+    let integral = (t: t) =>
       updateShape(Continuous(t.integralCache), t);
 
-    let downsample = (~cache=None, i, t): t =>
+    let downsample = (i, t): t =>
       updateShape(t |> toShape |> Shape.T.downsample(i), t);
     // todo: adjust for limit, maybe?
     let mapY =
         (
-          ~knownIntegralSumFn=previousIntegralSum => None,
+          ~integralSumCacheFn=previousIntegralSum => None,
+          ~integralCacheFn=previousIntegralCache => None,
           fn,
           {shape, _} as t: t,
         )
         : t =>
-      Shape.T.mapY(~knownIntegralSumFn, fn, shape)
+      Shape.T.mapY(~integralSumCacheFn, fn, shape)
       |> updateShape(_, t);
 
     // get the total of everything
-    let integralEndY = (~cache as _, t: t) => {
+    let integralEndY = (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) => {
+    let integralXtoY = (f, t: t) => {
       Shape.T.Integral.xToY(
-        ~cache=Some(t.integralCache),
         f,
         toShape(t),
       )
@@ -140,8 +140,8 @@ module T =
     };
 
     // 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) => {
+    let integralYtoX = (f, t: t) => {
-      Shape.T.Integral.yToX(~cache=None, f, toShape(t));
+      Shape.T.Integral.yToX(f, toShape(t));
     };
 
     let mean = (t: t) => {

src/distPlus/distribution/DistPlusTime.re

@@ -23,6 +23,6 @@
     include DistPlus.T.Integral;
     let xToY = (f: TimeTypes.timeInVector, t: t) => {
       timeInVectorToX(f, t)
-      |> E.O.fmap(x => DistPlus.T.Integral.xToY(~cache=None, x, t));
+      |> E.O.fmap(x => DistPlus.T.Integral.xToY(x, t));
     };
-  };
+  };

src/distPlus/distribution/DistTypes.re

@@ -14,20 +14,36 @@ type xyShape = {
   ys: array(float),
 };
 
+type interpolation = [
+    | `Stepwise
+    | `Linear
+];
+type extrapolation = [
+    | `UseZero
+    | `UseOutermostPoints
+];
+
+type interpolator = (xyShape, int, float) => float;
+
 type continuousShape = {
   xyShape,
-  interpolation: [ | `Stepwise | `Linear],
+  interpolation: interpolation,
-  knownIntegralSum: option(float),
+  integralSumCache: option(float),
+  integralCache: option(continuousShape),
 };
 
 type discreteShape = {
   xyShape,
-  knownIntegralSum: option(float),
+  /* interpolation is always `Discrete */
+  integralSumCache: option(float),
+  integralCache: option(continuousShape),
 };
 
 type mixedShape = {
   continuous: continuousShape,
   discrete: discreteShape,
+  integralSumCache: option(float),
+  integralCache: option(continuousShape),
 };
 
 type shapeMonad('a, 'b, 'c) =

src/distPlus/distribution/Distributions.re

@@ -4,7 +4,7 @@ module type dist = {
   let minX: t => float;
   let maxX: t => float;
   let mapY:
-    (~knownIntegralSumFn: float => option(float)=?, float => float, t) => t;
+    (~integralSumCacheFn: float => option(float)=?, ~integralCacheFn: DistTypes.continuousShape => option(DistTypes.continuousShape)=?, float => float, t) => t;
   let xToY: (float, t) => DistTypes.mixedPoint;
   let toShape: t => DistTypes.shape;
   let toContinuous: t => option(DistTypes.continuousShape);
@@ -13,13 +13,13 @@ module type dist = {
   let normalizedToContinuous: t => option(DistTypes.continuousShape);
   let normalizedToDiscrete: t => option(DistTypes.discreteShape);
   let toDiscreteProbabilityMassFraction: t => float;
-  let downsample: (~cache: option(integral)=?, int, t) => t;
+  let downsample: (int, t) => t;
   let truncate: (option(float), option(float), t) => t;
 
-  let integral: (~cache: option(integral), t) => integral;
+  let integral: (t) => integral;
-  let integralEndY: (~cache: option(integral), t) => float;
+  let integralEndY: (t) => float;
-  let integralXtoY: (~cache: option(integral), float, t) => float;
+  let integralXtoY: (float, t) => float;
-  let integralYtoX: (~cache: option(integral), float, t) => float;
+  let integralYtoX: (float, t) => float;
 
   let mean: t => float;
   let variance: t => float;
@@ -59,10 +59,23 @@ module Common = {
   let combineIntegralSums =
       (
         combineFn: (float, float) => option(float),
-        t1KnownIntegralSum: option(float),
+        t1IntegralSumCache: option(float),
-        t2KnownIntegralSum: option(float),
+        t2IntegralSumCache: option(float),
       ) => {
-    switch (t1KnownIntegralSum, t2KnownIntegralSum) {
+    switch (t1IntegralSumCache, t2IntegralSumCache) {
+    | (None, _)
+    | (_, None) => None
+    | (Some(s1), Some(s2)) => combineFn(s1, s2)
+    };
+  };
+
+  let combineIntegrals =
+      (
+        combineFn: (DistTypes.continuousShape, DistTypes.continuousShape) => option(DistTypes.continuousShape),
+        t1IntegralCache: option(DistTypes.continuousShape),
+        t2IntegralCache: option(DistTypes.continuousShape),
+      ) => {
+    switch (t1IntegralCache, t2IntegralCache) {
     | (None, _)
     | (_, None) => None
     | (Some(s1), Some(s2)) => combineFn(s1, s2)

src/distPlus/distribution/Mixed.re

@@ -1,7 +1,7 @@
 open Distributions;
 
 type t = DistTypes.mixedShape;
-let make = (~continuous, ~discrete): t => {continuous, discrete};
+let make = (~continuous, ~discrete, integralSumCache, integralCache): t => {continuous, discrete, integralSumCache, integralCache};
 
 let totalLength = (t: t): int => {
   let continuousLength =
@@ -11,31 +11,17 @@ let totalLength = (t: t): int => {
   continuousLength + discreteLength;
 };
 
-let scaleBy = (~scale=1.0, {discrete, continuous}: t): t => {
+let scaleBy = (~scale=1.0, t: t): t => {
-  let scaledDiscrete = Discrete.scaleBy(~scale, discrete);
+  let scaledDiscrete = Discrete.scaleBy(~scale, t.discrete);
-  let scaledContinuous = Continuous.scaleBy(~scale, continuous);
+  let scaledContinuous = Continuous.scaleBy(~scale, t.continuous);
-  make(~discrete=scaledDiscrete, ~continuous=scaledContinuous);
+  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);
 };
 
 let toContinuous = ({continuous}: t) => Some(continuous);
 let toDiscrete = ({discrete}: t) => Some(discrete);
 
-let combinePointwise = (~knownIntegralSumsFn, fn, t1: t, t2: t) => {
-  let reducedDiscrete =
-    [|t1, t2|]
-    |> E.A.fmap(toDiscrete)
-    |> E.A.O.concatSomes
-    |> Discrete.reduce(~knownIntegralSumsFn, fn);
-
-  let reducedContinuous =
-    [|t1, t2|]
-    |> E.A.fmap(toContinuous)
-    |> E.A.O.concatSomes
-    |> Continuous.reduce(~knownIntegralSumsFn, fn);
-
-  make(~discrete=reducedDiscrete, ~continuous=reducedContinuous);
-};
-
 module T =
   Dist({
     type t = DistTypes.mixedShape;
@@ -61,30 +47,35 @@ module T =
       let truncatedDiscrete =
         Discrete.T.truncate(leftCutoff, rightCutoff, discrete);
 
-      make(~discrete=truncatedDiscrete, ~continuous=truncatedContinuous);
+      make(~discrete=truncatedDiscrete, ~continuous=truncatedContinuous, None, None);
     };
 
     let normalize = (t: t): t => {
+      let continuousIntegral = Continuous.T.Integral.get(t.continuous);
+      let discreteIntegral = Discrete.T.Integral.get(t.discrete);
+
+      let continuous = t.continuous |> Continuous.updateIntegralCache(Some(continuousIntegral));
+      let discrete = t.discrete |> Discrete.updateIntegralCache(Some(discreteIntegral));
 
       let continuousIntegralSum =
-        Continuous.T.Integral.sum(~cache=None, t.continuous);
+        Continuous.T.Integral.sum(continuous);
       let discreteIntegralSum =
-        Discrete.T.Integral.sum(~cache=None, t.discrete);
+        Discrete.T.Integral.sum(discrete);
       let totalIntegralSum = continuousIntegralSum +. discreteIntegralSum;
 
       let newContinuousSum = continuousIntegralSum /. totalIntegralSum;
       let newDiscreteSum = discreteIntegralSum /. totalIntegralSum;
 
       let normalizedContinuous =
-        t.continuous
+        continuous
         |> Continuous.scaleBy(~scale=newContinuousSum /. continuousIntegralSum)
-        |> Continuous.updateKnownIntegralSum(Some(newContinuousSum));
+        |> Continuous.updateIntegralSumCache(Some(newContinuousSum));
       let normalizedDiscrete =
-        t.discrete
+        discrete
         |> Discrete.scaleBy(~scale=newDiscreteSum /. discreteIntegralSum)
-        |> Discrete.updateKnownIntegralSum(Some(newDiscreteSum));
+        |> Discrete.updateIntegralSumCache(Some(newDiscreteSum));
 
-      make(~continuous=normalizedContinuous, ~discrete=normalizedDiscrete);
+      make(~continuous=normalizedContinuous, ~discrete=normalizedDiscrete, Some(1.0), None);
     };
 
     let xToY = (x, t: t) => {
@@ -98,23 +89,22 @@ module T =
 
     let toDiscreteProbabilityMassFraction = ({discrete, continuous}: t) => {
       let discreteIntegralSum =
-        Discrete.T.Integral.sum(~cache=None, discrete);
+        Discrete.T.Integral.sum(discrete);
       let continuousIntegralSum =
-        Continuous.T.Integral.sum(~cache=None, continuous);
+        Continuous.T.Integral.sum(continuous);
       let totalIntegralSum = discreteIntegralSum +. continuousIntegralSum;
 
       discreteIntegralSum /. totalIntegralSum;
     };
 
-    let downsample = (~cache=None, count, {discrete, continuous}: t): t => {
+    let downsample = (count, t: t): t => {
       // We will need to distribute the new xs fairly between the discrete and continuous shapes.
       // The easiest way to do this is to simply go by the previous probability masses.
 
-      // The cache really isn't helpful here, because we would need two separate caches
       let discreteIntegralSum =
-        Discrete.T.Integral.sum(~cache=None, discrete);
+        Discrete.T.Integral.sum(t.discrete);
       let continuousIntegralSum =
-        Continuous.T.Integral.sum(~cache=None, continuous);
+        Continuous.T.Integral.sum(t.continuous);
       let totalIntegralSum = discreteIntegralSum +. continuousIntegralSum;
 
       // TODO: figure out what to do when the totalIntegralSum is zero.
@@ -124,7 +114,7 @@ module T =
           int_of_float(
             float_of_int(count) *. (discreteIntegralSum /. totalIntegralSum),
           ),
-          discrete,
+          t.discrete,
         );
 
       let downsampledContinuous =
@@ -132,10 +122,10 @@ module T =
           int_of_float(
             float_of_int(count) *. (continuousIntegralSum /. totalIntegralSum),
           ),
-          continuous,
+          t.continuous,
         );
 
-      {discrete: downsampledDiscrete, continuous: downsampledContinuous};
+      {...t, discrete: downsampledDiscrete, continuous: downsampledContinuous};
     };
 
     let normalizedToContinuous = (t: t) => Some(normalize(t).continuous);
@@ -143,66 +133,69 @@ module T =
     let normalizedToDiscrete = ({discrete} as t: t) =>
       Some(normalize(t).discrete);
 
-    let integral = (~cache, {continuous, discrete}: t) => {
+    let integral = (t: t) => {
-      switch (cache) {
+      switch (t.integralCache) {
       | Some(cache) => cache
       | None =>
         // note: if the underlying shapes aren't normalized, then these integrals won't be either -- but that's the way it should be.
-        let continuousIntegral =
+        let continuousIntegral = Continuous.T.Integral.get(t.continuous);
-          Continuous.T.Integral.get(~cache=None, continuous);
+        let discreteIntegral = Continuous.stepwiseToLinear(Discrete.T.Integral.get(t.discrete));
-        let discreteIntegral = Discrete.T.Integral.get(~cache=None, discrete);
 
         Continuous.make(
           `Linear,
           XYShape.PointwiseCombination.combine(
             (+.),
-            XYShape.XtoY.linearBetweenPointsExtrapolateFlat,
+            XYShape.XtoY.continuousInterpolator(`Linear, `UseOutermostPoints),
-            XYShape.XtoY.stepwiseBetweenPointsExtrapolateFlat,
+            XYShape.XtoY.continuousInterpolator(`Linear, `UseOutermostPoints),
             Continuous.getShape(continuousIntegral),
             Continuous.getShape(discreteIntegral),
           ),
           None,
+          None,
         );
       };
     };
 
-    let integralEndY = (~cache, t: t) => {
+    let integralEndY = (t: t) => {
-      integral(~cache, t) |> Continuous.lastY;
+      t |> integral |> Continuous.lastY;
     };
 
-    let integralXtoY = (~cache, f, t) => {
+    let integralXtoY = (f, t) => {
-      t |> integral(~cache) |> Continuous.getShape |> XYShape.XtoY.linear(f);
+      t |> integral |> Continuous.getShape |> XYShape.XtoY.linear(f);
     };
 
-    let integralYtoX = (~cache, f, t) => {
+    let integralYtoX = (f, t) => {
-      t |> integral(~cache) |> Continuous.getShape |> XYShape.YtoX.linear(f);
+      t |> integral |> Continuous.getShape |> XYShape.YtoX.linear(f);
     };
 
     // This pipes all ys (continuous and discrete) through fn.
-    // If mapY is a linear operation, we might be able to update the knownIntegralSums as well;
+    // If mapY is a linear operation, we might be able to update the integralSumCaches as well;
     // if not, they'll be set to None.
     let mapY =
         (
-          ~knownIntegralSumFn=previousIntegralSum => None,
+          ~integralSumCacheFn=previousIntegralSum => None,
+          ~integralCacheFn=previousIntegral => None,
           fn,
-          {discrete, continuous}: t,
+          t: t,
         )
         : t => {
-      let u = E.O.bind(_, knownIntegralSumFn);
+      let yMappedDiscrete: DistTypes.discreteShape =
-
+        t.discrete
-      let yMappedDiscrete =
-        discrete
         |> Discrete.T.mapY(fn)
-        |> Discrete.updateKnownIntegralSum(u(discrete.knownIntegralSum));
+        |> Discrete.updateIntegralSumCache(E.O.bind(t.discrete.integralSumCache, integralSumCacheFn))
+        |> Discrete.updateIntegralCache(E.O.bind(t.discrete.integralCache, integralCacheFn));
 
-      let yMappedContinuous =
+      let yMappedContinuous: DistTypes.continuousShape =
-        continuous
+        t.continuous
         |> Continuous.T.mapY(fn)
-        |> Continuous.updateKnownIntegralSum(u(continuous.knownIntegralSum));
+        |> Continuous.updateIntegralSumCache(E.O.bind(t.continuous.integralSumCache, integralSumCacheFn))
+        |> Continuous.updateIntegralCache(E.O.bind(t.continuous.integralCache, integralCacheFn));
 
       {
         discrete: yMappedDiscrete,
-        continuous: Continuous.T.mapY(fn, continuous),
+        continuous: yMappedContinuous,
+        integralSumCache: E.O.bind(t.integralSumCache, integralSumCacheFn),
+        integralCache: E.O.bind(t.integralCache, integralCacheFn),
       };
     };
 
@@ -211,10 +204,8 @@ module T =
       let continuousMean = Continuous.T.mean(continuous);
 
       // the combined mean is the weighted sum of the two:
-      let discreteIntegralSum =
+      let discreteIntegralSum = Discrete.T.Integral.sum(discrete);
-        Discrete.T.Integral.sum(~cache=None, discrete);
+      let continuousIntegralSum = Continuous.T.Integral.sum(continuous);
-      let continuousIntegralSum =
-        Continuous.T.Integral.sum(~cache=None, continuous);
       let totalIntegralSum = discreteIntegralSum +. continuousIntegralSum;
 
       (
@@ -228,10 +219,8 @@ module T =
 
     let variance = ({discrete, continuous} as t: t): float => {
       // the combined mean is the weighted sum of the two:
-      let discreteIntegralSum =
+      let discreteIntegralSum = Discrete.T.Integral.sum(discrete);
-        Discrete.T.Integral.sum(~cache=None, discrete);
+      let continuousIntegralSum = Continuous.T.Integral.sum(continuous);
-      let continuousIntegralSum =
-        Continuous.T.Integral.sum(~cache=None, continuous);
       let totalIntegralSum = discreteIntegralSum +. continuousIntegralSum;
 
       let getMeanOfSquares = ({discrete, continuous}: t) => {
@@ -260,7 +249,7 @@ module T =
   });
 
 let combineAlgebraically =
-    (~downsample=false, op: ExpressionTypes.algebraicOperation, t1: t, t2: t)
+    (op: ExpressionTypes.algebraicOperation, t1: t, t2: t)
     : t => {
   // Discrete convolution can cause a huge increase in the number of samples,
   // so we'll first downsample.
@@ -268,43 +257,79 @@ let combineAlgebraically =
   // An alternative (to be explored in the future) may be to first perform the full convolution and then to downsample the result;
   // to use non-uniform fast Fourier transforms (for addition only), add web workers or gpu.js, etc. ...
 
-  let downsampleIfTooLarge = (t: t) => {
+  // TODO: figure out when to downsample strategically. Could be an evaluationParam?
+  /*let downsampleIfTooLarge = (t: t) => {
     let sqtl = sqrt(float_of_int(totalLength(t)));
     sqtl > 10. && downsample ? T.downsample(int_of_float(sqtl), t) : t;
   };
 
   let t1d = downsampleIfTooLarge(t1);
   let t2d = downsampleIfTooLarge(t2);
+  */
 
   // continuous (*) continuous => continuous, but also
   // discrete (*) continuous => continuous (and vice versa). We have to take care of all combos and then combine them:
   let ccConvResult =
     Continuous.combineAlgebraically(
-      ~downsample=false,
       op,
-      t1d.continuous,
+      t1.continuous,
-      t2d.continuous,
+      t2.continuous,
     );
   let dcConvResult =
     Continuous.combineAlgebraicallyWithDiscrete(
-      ~downsample=false,
       op,
-      t2d.continuous,
+      t2.continuous,
-      t1d.discrete,
+      t1.discrete,
     );
   let cdConvResult =
     Continuous.combineAlgebraicallyWithDiscrete(
-      ~downsample=false,
       op,
-      t1d.continuous,
+      t1.continuous,
-      t2d.discrete,
+      t2.discrete,
     );
   let continuousConvResult =
     Continuous.reduce((+.), [|ccConvResult, dcConvResult, cdConvResult|]);
 
   // ... finally, discrete (*) discrete => discrete, obviously:
   let discreteConvResult =
-    Discrete.combineAlgebraically(op, t1d.discrete, t2d.discrete);
+    Discrete.combineAlgebraically(op, t1.discrete, t2.discrete);
 
-  {discrete: discreteConvResult, continuous: continuousConvResult};
+  let combinedIntegralSum =
+    Common.combineIntegralSums(
+      (a, b) => Some(a *. b),
+      t1.integralSumCache,
+      t2.integralSumCache,
+    );
+
+  {discrete: discreteConvResult, continuous: continuousConvResult, integralSumCache: combinedIntegralSum, integralCache: None};
+};
+
+let combinePointwise = (~integralSumCachesFn = (_, _) => None, ~integralCachesFn = (_, _) => None, fn, t1: t, t2: t): t => {
+  let reducedDiscrete =
+    [|t1, t2|]
+    |> E.A.fmap(toDiscrete)
+    |> E.A.O.concatSomes
+    |> Discrete.reduce(~integralSumCachesFn, ~integralCachesFn, fn);
+
+  let reducedContinuous =
+    [|t1, t2|]
+    |> E.A.fmap(toContinuous)
+    |> E.A.O.concatSomes
+    |> Continuous.reduce(~integralSumCachesFn, ~integralCachesFn, fn);
+
+  let combinedIntegralSum =
+    Common.combineIntegralSums(
+      integralSumCachesFn,
+      t1.integralSumCache,
+      t2.integralSumCache,
+    );
+
+  let combinedIntegral =
+    Common.combineIntegrals(
+      integralCachesFn,
+      t1.integralCache,
+      t2.integralCache,
+    );
+
+  make(~discrete=reducedDiscrete, ~continuous=reducedContinuous, combinedIntegralSum, combinedIntegral);
 };

src/distPlus/distribution/MixedShapeBuilder.re

@@ -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)));
+  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)));
+  let discrete = discrete |> E.O.default(Discrete.make({xs: [||], ys: [||]}, Some(0.0), None));
   let cLength =
     continuous
     |> Continuous.getShape
@@ -25,7 +25,9 @@ let buildSimple = (~continuous: option(DistTypes.continuousShape), ~discrete: op
     let mixedDist =
       Mixed.make(
         ~continuous,
-        ~discrete
+        ~discrete,
+        None,
+        None,
       );
     Some(Mixed(mixedDist));
   };

src/distPlus/distribution/Shape.re

@@ -15,11 +15,12 @@ let fmap = ((fn1, fn2, fn3), t: t): t =>
   | Continuous(m) => Continuous(fn3(m))
   };
 
+
 let toMixed =
   mapToAll((
     m => m,
-    d => Mixed.make(~discrete=d, ~continuous=Continuous.empty),
+    d => Mixed.make(~discrete=d, ~continuous=Continuous.empty, d.integralSumCache, d.integralCache),
-    c => Mixed.make(~discrete=Discrete.empty, ~continuous=c),
+    c => Mixed.make(~discrete=Discrete.empty, ~continuous=c, c.integralSumCache, c.integralCache),
   ));
 
 let combineAlgebraically =
@@ -27,19 +28,18 @@ let combineAlgebraically =
   switch (t1, t2) {
   | (Continuous(m1), Continuous(m2)) =>
     DistTypes.Continuous(
-      Continuous.combineAlgebraically(~downsample=true, op, m1, m2),
+      Continuous.combineAlgebraically(op, m1, m2),
     )
   | (Continuous(m1), Discrete(m2))
   | (Discrete(m2), Continuous(m1)) =>
     DistTypes.Continuous(
-      Continuous.combineAlgebraicallyWithDiscrete(~downsample=false, op, m1, m2),
+      Continuous.combineAlgebraicallyWithDiscrete(op, m1, m2),
     )
   | (Discrete(m1), Discrete(m2)) =>
     DistTypes.Discrete(Discrete.combineAlgebraically(op, m1, m2))
   | (m1, m2) =>
     DistTypes.Mixed(
       Mixed.combineAlgebraically(
-        ~downsample=true,
         op,
         toMixed(m1),
         toMixed(m2),
@@ -49,20 +49,25 @@ let combineAlgebraically =
 };
 
 let combinePointwise =
-    (~knownIntegralSumsFn=(_, _) => None, fn, t1: t, t2: t) =>
+  (~integralSumCachesFn: (float, float) => option(float) = (_, _) => None,
+   ~integralCachesFn: (DistTypes.continuousShape, DistTypes.continuousShape) => option(DistTypes.continuousShape) = (_, _) => None,
+   fn,
+   t1: t,
+   t2: t) =>
   switch (t1, t2) {
   | (Continuous(m1), Continuous(m2)) =>
     DistTypes.Continuous(
-      Continuous.combinePointwise(~knownIntegralSumsFn, fn, m1, m2),
+      Continuous.combinePointwise(~integralSumCachesFn, ~integralCachesFn, fn, m1, m2),
     )
   | (Discrete(m1), Discrete(m2)) =>
     DistTypes.Discrete(
-      Discrete.combinePointwise(~knownIntegralSumsFn, fn, m1, m2),
+      Discrete.combinePointwise(~integralSumCachesFn, ~integralCachesFn, fn, m1, m2),
     )
   | (m1, m2) =>
     DistTypes.Mixed(
       Mixed.combinePointwise(
-        ~knownIntegralSumsFn,
+        ~integralSumCachesFn,
+        ~integralCachesFn,
         fn,
         toMixed(m1),
         toMixed(m2),
@@ -87,7 +92,7 @@ module T =
     let toContinuous = t => None;
     let toDiscrete = t => None;
 
-    let downsample = (~cache=None, i, t) =>
+    let downsample = (i, t) =>
       fmap(
         (
           Mixed.T.downsample(i),
@@ -108,8 +113,14 @@ module T =
       );
 
     let toDiscreteProbabilityMassFraction = t => 0.0;
+
     let normalize =
-      fmap((Mixed.T.normalize, Discrete.T.normalize, Continuous.T.normalize));
+      fmap((
+        Mixed.T.normalize,
+        Discrete.T.normalize,
+        Continuous.T.normalize
+      ));
+
     let toContinuous =
       mapToAll((
         Mixed.T.toContinuous,
@@ -143,38 +154,38 @@ module T =
         Continuous.T.normalizedToContinuous,
       ));
     let minX = mapToAll((Mixed.T.minX, Discrete.T.minX, Continuous.T.minX));
-    let integral = (~cache) =>
+    let integral =
       mapToAll((
-        Mixed.T.Integral.get(~cache=None),
+        Mixed.T.Integral.get,
-        Discrete.T.Integral.get(~cache=None),
+        Discrete.T.Integral.get,
-        Continuous.T.Integral.get(~cache=None),
+        Continuous.T.Integral.get,
       ));
-    let integralEndY = (~cache) =>
+    let integralEndY =
       mapToAll((
-        Mixed.T.Integral.sum(~cache=None),
+        Mixed.T.Integral.sum,
-        Discrete.T.Integral.sum(~cache),
+        Discrete.T.Integral.sum,
-        Continuous.T.Integral.sum(~cache=None),
+        Continuous.T.Integral.sum,
       ));
-    let integralXtoY = (~cache, f) => {
+    let integralXtoY = (f) => {
       mapToAll((
-        Mixed.T.Integral.xToY(~cache, f),
+        Mixed.T.Integral.xToY(f),
-        Discrete.T.Integral.xToY(~cache, f),
+        Discrete.T.Integral.xToY(f),
-        Continuous.T.Integral.xToY(~cache, f),
+        Continuous.T.Integral.xToY(f),
       ));
     };
-    let integralYtoX = (~cache, f) => {
+    let integralYtoX = (f) => {
       mapToAll((
-        Mixed.T.Integral.yToX(~cache, f),
+        Mixed.T.Integral.yToX(f),
-        Discrete.T.Integral.yToX(~cache, f),
+        Discrete.T.Integral.yToX(f),
-        Continuous.T.Integral.yToX(~cache, f),
+        Continuous.T.Integral.yToX(f),
       ));
     };
     let maxX = mapToAll((Mixed.T.maxX, Discrete.T.maxX, Continuous.T.maxX));
-    let mapY = (~knownIntegralSumFn=previousIntegralSum => None, fn) =>
+    let mapY = (~integralSumCacheFn=previousIntegralSum => None, ~integralCacheFn=previousIntegral=>None, fn) =>
       fmap((
-        Mixed.T.mapY(~knownIntegralSumFn, fn),
+        Mixed.T.mapY(~integralSumCacheFn, ~integralCacheFn, fn),
-        Discrete.T.mapY(~knownIntegralSumFn, fn),
+        Discrete.T.mapY(~integralSumCacheFn, ~integralCacheFn, fn),
-        Continuous.T.mapY(~knownIntegralSumFn, fn),
+        Continuous.T.mapY(~integralSumCacheFn, ~integralCacheFn, fn),
       ));
 
     let mean = (t: t): float =>
@@ -197,8 +208,8 @@ let pdf = (f: float, t: t) => {
   mixedPoint.continuous +. mixedPoint.discrete;
 };
 
-let inv = T.Integral.yToX(~cache=None);
+let inv = T.Integral.yToX;
-let cdf = T.Integral.xToY(~cache=None);
+let cdf = T.Integral.xToY;
 
 let sample = (t: t): float => {
   // this can go, already taken care of in Ozzie's sampling branch

src/distPlus/distribution/XYShape.re

@@ -19,6 +19,7 @@ module T = {
   let ys = (t: t) => t.ys;
   let length = (t: t) => E.A.length(t.xs);
   let empty = {xs: [||], ys: [||]};
+  let isEmpty = (t: t) => length(t) == 0;
   let minX = (t: t) => t |> xs |> E.A.Sorted.min |> extImp;
   let maxX = (t: t) => t |> xs |> E.A.Sorted.max |> extImp;
   let firstY = (t: t) => t |> ys |> E.A.first |> extImp;
@@ -143,60 +144,66 @@ module XtoY = {
     n;
   };
 
-  /* The extrapolators below can be used e.g. with PointwiseCombination.combine */
+  /* Returns a between-points-interpolating function that can be used with PointwiseCombination.combine.
-  let linearBetweenPointsExtrapolateFlat = (t: T.t, leftIndex: int, x: float) => {
+     Interpolation can either be stepwise (using the value on the left) or linear. Extrapolation can be `UseZero or `UseOutermostPoints. */
-    if (leftIndex < 0) {
+  let continuousInterpolator = (interpolation: DistTypes.interpolation, extrapolation: DistTypes.extrapolation): interpolator => {
-      t.ys[0];
+    switch (interpolation) {
-    } else if (leftIndex >= T.length(t) - 1) {
+    | `Linear => switch (extrapolation) {
-      T.lastY(t);
+      | `UseZero => (t: T.t, leftIndex: int, x: float) => {
-    } else {
+        if (leftIndex < 0) {
-      let x1 = t.xs[leftIndex];
+          0.0
-      let x2 = t.xs[leftIndex + 1];
+        } else if (leftIndex >= T.length(t) - 1) {
-      let y1 = t.ys[leftIndex];
+          0.0
-      let y2 = t.ys[leftIndex + 1];
+        } else {
-      let fraction = (x -. x1) /. (x2 -. x1);
+          let x1 = t.xs[leftIndex];
-      y1 *. (1. -. fraction) +. y2 *. fraction;
+          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;
+        };
+      }
+    }
+    | `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];
+        }
+      }
+    }
     }
   };
 
-  let linearBetweenPointsExtrapolateZero = (t: T.t, leftIndex: int, x: float) => {
+  /* Returns a between-points-interpolating function that can be used with PointwiseCombination.combine.
-    if (leftIndex < 0) {
+     For discrete distributions, the probability density between points is zero, so we just return zero here. */
-      0.0
+  let discreteInterpolator: interpolator = (t: T.t, leftIndex: int, x: float) => 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;
-    }
-  };
-
-  let stepwiseBetweenPointsExtrapolateFlat = (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];
-    }
-  };
-
-  let stepwiseBetweenPointsExtrapolateZero = (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];
-    }
-  };
-
-  let assumeZeroBetweenPoints = (t: T.t, leftIndex: int, x: float) => {
-    0.0;
-  };
 };
 
 module XsConversion = {
@@ -371,6 +378,11 @@ module Range = {
 
   let derivative = mapYsBasedOnRanges(delta_y_over_delta_x);
 
+  let stepwiseToLinear = t => {
+    // adds points at the bottom of each step.
+    t;
+  };
+
   // TODO: It would be nicer if this the diff didn't change the first element, and also maybe if there were a more elegant way of doing this.
   let stepsToContinuous = t => {
     let diff = T.xTotalRange(t) |> (r => r *. 0.00001);

src/distPlus/expressionTree/ExpressionTreeEvaluator.re

@@ -61,7 +61,8 @@ module VerticalScaling = {
       (evaluationParams: evaluationParams, scaleOp, t, scaleBy) => {
     // scaleBy has to be a single float, otherwise we'll return an error.
     let fn = Operation.Scale.toFn(scaleOp);
-    let knownIntegralSumFn = Operation.Scale.toKnownIntegralSumFn(scaleOp);
+    let integralSumCacheFn = Operation.Scale.toIntegralSumCacheFn(scaleOp);
+    let integralCacheFn = Operation.Scale.toIntegralCacheFn(scaleOp);
     let renderedShape = render(evaluationParams, t);
 
     switch (renderedShape, scaleBy) {
@@ -69,7 +70,8 @@ module VerticalScaling = {
       Ok(
         `RenderedDist(
           Shape.T.mapY(
-            ~knownIntegralSumFn=knownIntegralSumFn(sm),
+            ~integralSumCacheFn=integralSumCacheFn(sm),
+            ~integralCacheFn=integralCacheFn(sm),
             fn(sm),
             rs,
           ),
@@ -82,13 +84,14 @@ module VerticalScaling = {
 };
 
 module PointwiseCombination = {
-  let pointwiseAdd = (evaluationParams: evaluationParams, t1, t2) => {
+  let pointwiseAdd = (evaluationParams: evaluationParams, t1: t, t2: t) => {
     switch (render(evaluationParams, t1), render(evaluationParams, t2)) {
     | (Ok(`RenderedDist(rs1)), Ok(`RenderedDist(rs2))) =>
       Ok(
         `RenderedDist(
           Shape.combinePointwise(
-            ~knownIntegralSumsFn=(a, b) => Some(a +. b),
+            ~integralSumCachesFn=(a, b) => Some(a +. b),
+            ~integralCachesFn=(a, b) => Some(Continuous.combinePointwise(~extrapolation=`UseOutermostPoints, (+.), a, b)),
             (+.),
             rs1,
             rs2,
@@ -101,7 +104,7 @@ module PointwiseCombination = {
     };
   };
 
-  let pointwiseMultiply = (evaluationParams: evaluationParams, t1, t2) => {
+  let pointwiseMultiply = (evaluationParams: evaluationParams, t1: t, t2: t) => {
     // TODO: construct a function that we can easily sample from, to construct
     // a RenderedDist. Use the xMin and xMax of the rendered shapes to tell the sampling function where to look.
     Error(
@@ -110,7 +113,7 @@ module PointwiseCombination = {
   };
 
   let operationToLeaf =
-      (evaluationParams: evaluationParams, pointwiseOp, t1, t2) => {
+      (evaluationParams: evaluationParams, pointwiseOp: pointwiseOperation, t1: t, t2: t) => {
     switch (pointwiseOp) {
     | `Add => pointwiseAdd(evaluationParams, t1, t2)
     | `Multiply => pointwiseMultiply(evaluationParams, t1, t2)
@@ -227,7 +230,6 @@ let toLeaf =
       node: t,
     )
     : result(t, string) => {
-  Js.log2("EVALUATION PARAMS", evaluationParams);
   switch (node) {
   // Leaf nodes just stay leaf nodes
   | `SymbolicDist(_)

src/distPlus/expressionTree/Operation.re

@@ -64,11 +64,17 @@ module Scale = {
     | `Log => {j|verticalLog($value, $scaleBy) |j}
     };
 
-  let toKnownIntegralSumFn =
+  let toIntegralSumCacheFn =
     fun
     | `Multiply => ((a, b) => Some(a *. b))
     | `Exponentiate => ((_, _) => None)
     | `Log => ((_, _) => None);
+
+  let toIntegralCacheFn =
+    fun
+    | `Multiply => ((a, b) => None) // TODO: this could probably just be multiplied out (using Continuous.scaleBy)
+    | `Exponentiate => ((_, _) => None)
+    | `Log => ((_, _) => None);
 };
 
 module T = {

src/distPlus/renderers/samplesRenderer/Guesstimator.re

@@ -7,7 +7,7 @@ type discrete = {
 let jsToDistDiscrete = (d: discrete): DistTypes.discreteShape => {xyShape: {
   xs: xsGet(d),
   ys: ysGet(d),
-}, knownIntegralSum: None};
+}, integralSumCache: None, integralCache: None};
 
 [@bs.module "./GuesstimatorLibrary.js"]
 external stringToSamples: (string, int) => array(float) = "stringToSamples";

src/distPlus/renderers/samplesRenderer/Samples.re

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

src/distPlus/symbolic/SymbolicDist.re

@@ -299,13 +299,13 @@ module T = {
     switch (d) {
     | `Float(v) =>
       Discrete(
-        Discrete.make({xs: [|v|], ys: [|1.0|]}, Some(1.0)),
+        Discrete.make({xs: [|v|], ys: [|1.0|]}, Some(1.0), None),
       )
     | _ =>
       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)),
+        Continuous.make(`Linear, {xs, ys}, Some(1.0), None),
       );
     };
 };

src/distPlus/utility/E.re

@@ -464,4 +464,4 @@ module JsArray = {
          Rationale.Option.toExn("Warning: This should not have happened"),
        );
   let filter = Js.Array.filter;
-};
+};