diff --git a/packages/squiggle-lang/__tests__/ReducerInterface/ReducerInterface_Distribution_test.res b/packages/squiggle-lang/__tests__/ReducerInterface/ReducerInterface_Distribution_test.res
index 605797b9..a4569a88 100644
--- a/packages/squiggle-lang/__tests__/ReducerInterface/ReducerInterface_Distribution_test.res
+++ b/packages/squiggle-lang/__tests__/ReducerInterface/ReducerInterface_Distribution_test.res
@@ -54,6 +54,7 @@ describe("eval on distribution functions", () => {
   describe("subtract", () => {
     testEval("10 - normal(5, 1)", "Ok(Normal(5,1))")
     testEval("normal(5, 1) - 10", "Ok(Normal(-5,1))")
+    testEval("mean(1 - toPointSet(normal(5, 2)))", "Ok(-4.002309896304692)")
   })
   describe("multiply", () => {
     testEval("normal(10, 2) * 2", "Ok(Normal(20,4))")
diff --git a/packages/squiggle-lang/src/rescript/Distributions/PointSetDist/AlgebraicShapeCombination.res b/packages/squiggle-lang/src/rescript/Distributions/PointSetDist/AlgebraicShapeCombination.res
index 673f53f7..cf563d91 100644
--- a/packages/squiggle-lang/src/rescript/Distributions/PointSetDist/AlgebraicShapeCombination.res
+++ b/packages/squiggle-lang/src/rescript/Distributions/PointSetDist/AlgebraicShapeCombination.res
@@ -209,9 +209,13 @@ let combineShapesContinuousDiscrete = (
       // 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 {
+        // When this operation is flipped (like 1 - normal(5, 2)) then the
+        // x axis coordinates would all come out the wrong order. So we need
+        // to fill them out in the opposite direction
+        let index = flip ? t1n - 1 - i : i
         Belt.Array.set(
           dxyShape,
-          i,
+          index,
           (
             fn(continuousShape.xs[i], discreteShape.xs[j]),
             continuousShape.ys[i] *. discreteShape.ys[j],