fix(rollups): rebase regression x-axis to a fixed origin for exact accumulators

The per-bucket OLS accumulators stored x on the raw epoch-day axis (x ~ 20540),
so sum_xx accumulated past ~1e10 and the square shed ~10 of a double's
significant digits before the value was ever stored. The sub-day x detail was
lost at accumulation time, which no read-side identity can recover; an
ill-conditioned window's composed slope/r2 then drifted from the live REGR_*
probe past the parity gauge.

Rebase x to a fixed origin (epoch-days of 2020-01-01) at write time so the
squared terms stay O(1e7) and exact. Slope, r2, and population sd are invariant
under the affine x-shift, so the live raw-epoch probe still parity-matches the
cross-bucket compose; composeRegression needs no formula change. The per-bucket
slope/r2 columns stay on the raw axis (shift-invariant).

Migration 0193 recomputes sum_x/sum_xy/sum_xx from measurements on the rebased
basis, one set-based update per granularity. Recompute-from-source is idempotent
and rerun-safe; it touches only rows that already carry the accumulators and
leaves pre-migration NULLs for the boot re-fold.

Un-quarantine the DST-boundary parity test: it now matches live REGR_* within
the existing 1e-9 relative bound, with the WEIGHT case still bit-identical.

prisma/migrations/0193_v1210_rollup_x_rescale/migration.sql

@@ -0,0 +1,148 @@
+-- v1.21.0 — rebase the regression accumulators to a fixed recent x-origin.
+--
+-- Migration 0190 stored the OLS accumulators `sum_x / sum_xy / sum_xx` over a
+-- RAW epoch-day x-axis (`EXTRACT(EPOCH FROM measured_at) / 86400.0`). Raw
+-- epoch-days sit at x ≈ 20 540, so `x²` lands near 4.2e8 and the per-bucket
+-- `SUM(x²)` accumulates past ~1e10. A double carries ~15-16 significant decimal
+-- digits; a value ~1e10 has ~10 integer digits, leaving only ~5-6 fractional
+-- digits — so the SUB-DAY x detail (the time-of-day fraction) is already lost
+-- when `sum_xx` is squared and summed, BEFORE the value is ever stored. No
+-- read-side identity (mean-centering included) can recover bits lost at
+-- accumulation time. That precision floor is why the windowed slope / r² on a
+-- near-flat, ill-conditioned window drifted from the live REGR_* probe past the
+-- 1e-9 parity gauge (the rollup-regression-parity DST case).
+--
+-- The fix is a WINDOW-LOCAL X-RESCALE: accumulate x RELATIVE to a fixed origin
+-- so the squared terms stay small and exact. The origin is 2020-01-01
+-- (epoch-day 18262); rebased x = epoch_days − 18262 stays in the low thousands
+-- for any realistic reading, so `x²` ≤ ~1e7 and `sum_xx` never sheds precision.
+--
+--   sum_x'  = Σ (epoch_days − 18262)
+--   sum_xy' = Σ (epoch_days − 18262) · value
+--   sum_xx' = Σ (epoch_days − 18262)²
+--
+-- Slope / r² / population-sd are INVARIANT under an affine x-shift (Sxx, Sxy,
+-- Syy are unchanged when a constant is subtracted from every x), so the
+-- cross-bucket compose yields the SAME regression as before — only the
+-- unreported intercept would move. The live REGR_* probe stays on raw
+-- epoch-days and still parity-matches because it, too, is shift-invariant.
+-- `sum_yy / sum_value / count / mean / sd / slope / r2` are y-only or
+-- shift-invariant and are NOT touched.
+--
+-- This recomputes the three rebased accumulators DIRECTLY FROM `measurements`,
+-- grouped exactly as the writer groups (per user / type / granularity-bucket /
+-- source over non-deleted rows). Recompute-from-source makes the migration
+-- naturally IDEMPOTENT and self-correcting: a rerun re-derives the same rebased
+-- sums regardless of the column's prior state (raw 0190 basis, partially
+-- rebased, or already rebased). It only writes rows that already carry the 0190
+-- accumulators (`sum_xy IS NOT NULL`), so NULL pre-migration rows stay NULL for
+-- the boot re-fold to fill on the new basis. Set-based, one statement per
+-- granularity; the `(user_id, type, measured_at)` index bounds the scan.
+--
+-- Reversibility (down): re-fold the accumulators on the raw epoch-day basis,
+-- i.e. drop the `− 18262` from each SUM. There is no column shape change to
+-- revert; the legacy basis is recoverable by re-running 0190's populator.
+
+-- DAY buckets.
+UPDATE "measurement_rollups" r
+SET "sum_x"  = s."sum_x",
+    "sum_xy" = s."sum_xy",
+    "sum_xx" = s."sum_xx"
+FROM (
+  SELECT
+    m."user_id" AS user_id,
+    m."type"    AS type,
+    date_trunc('day', m."measured_at") AS bucket_start,
+    m."source"  AS source,
+    SUM(EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - 18262) AS sum_x,
+    SUM((EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - 18262) * m."value") AS sum_xy,
+    SUM((EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - 18262)
+        * (EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - 18262)) AS sum_xx
+  FROM measurements m
+  WHERE m."deleted_at" IS NULL
+  GROUP BY m."user_id", m."type", date_trunc('day', m."measured_at"), m."source"
+) s
+WHERE r."granularity" = 'DAY'
+  AND r."sum_xy" IS NOT NULL
+  AND r."user_id"      = s.user_id
+  AND r."type"         = s.type
+  AND r."bucket_start" = s.bucket_start
+  AND r."source"       = s.source;
+
+-- WEEK buckets.
+UPDATE "measurement_rollups" r
+SET "sum_x"  = s."sum_x",
+    "sum_xy" = s."sum_xy",
+    "sum_xx" = s."sum_xx"
+FROM (
+  SELECT
+    m."user_id" AS user_id,
+    m."type"    AS type,
+    date_trunc('week', m."measured_at") AS bucket_start,
+    m."source"  AS source,
+    SUM(EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - 18262) AS sum_x,
+    SUM((EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - 18262) * m."value") AS sum_xy,
+    SUM((EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - 18262)
+        * (EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - 18262)) AS sum_xx
+  FROM measurements m
+  WHERE m."deleted_at" IS NULL
+  GROUP BY m."user_id", m."type", date_trunc('week', m."measured_at"), m."source"
+) s
+WHERE r."granularity" = 'WEEK'
+  AND r."sum_xy" IS NOT NULL
+  AND r."user_id"      = s.user_id
+  AND r."type"         = s.type
+  AND r."bucket_start" = s.bucket_start
+  AND r."source"       = s.source;
+
+-- MONTH buckets.
+UPDATE "measurement_rollups" r
+SET "sum_x"  = s."sum_x",
+    "sum_xy" = s."sum_xy",
+    "sum_xx" = s."sum_xx"
+FROM (
+  SELECT
+    m."user_id" AS user_id,
+    m."type"    AS type,
+    date_trunc('month', m."measured_at") AS bucket_start,
+    m."source"  AS source,
+    SUM(EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - 18262) AS sum_x,
+    SUM((EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - 18262) * m."value") AS sum_xy,
+    SUM((EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - 18262)
+        * (EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - 18262)) AS sum_xx
+  FROM measurements m
+  WHERE m."deleted_at" IS NULL
+  GROUP BY m."user_id", m."type", date_trunc('month', m."measured_at"), m."source"
+) s
+WHERE r."granularity" = 'MONTH'
+  AND r."sum_xy" IS NOT NULL
+  AND r."user_id"      = s.user_id
+  AND r."type"         = s.type
+  AND r."bucket_start" = s.bucket_start
+  AND r."source"       = s.source;
+
+-- YEAR buckets.
+UPDATE "measurement_rollups" r
+SET "sum_x"  = s."sum_x",
+    "sum_xy" = s."sum_xy",
+    "sum_xx" = s."sum_xx"
+FROM (
+  SELECT
+    m."user_id" AS user_id,
+    m."type"    AS type,
+    date_trunc('year', m."measured_at") AS bucket_start,
+    m."source"  AS source,
+    SUM(EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - 18262) AS sum_x,
+    SUM((EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - 18262) * m."value") AS sum_xy,
+    SUM((EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - 18262)
+        * (EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - 18262)) AS sum_xx
+  FROM measurements m
+  WHERE m."deleted_at" IS NULL
+  GROUP BY m."user_id", m."type", date_trunc('year', m."measured_at"), m."source"
+) s
+WHERE r."granularity" = 'YEAR'
+  AND r."sum_xy" IS NOT NULL
+  AND r."user_id"      = s.user_id
+  AND r."type"         = s.type
+  AND r."bucket_start" = s.bucket_start
+  AND r."source"       = s.source;

prisma/schema.prisma

@@ -1313,9 +1313,12 @@ model MeasurementRollup {
   /// `n = count` and `Σy = mean·count` these close the additive set so a
   /// windowed slope / r² / sd composes from summed accumulators across DAY
   /// buckets with a bit-identical result to live REGR_SLOPE / REGR_R2 /
-  /// STDDEV_POP. The x-axis is epoch days (`EXTRACT(EPOCH FROM
-  /// measured_at) / 86400.0`), stored un-rebased so the cross-bucket sum
-  /// reconstructs the live regression. Nullable to keep migration 0190
+  /// STDDEV_POP. The x-axis is epoch days REBASED to 2020-01-01
+  /// (`EXTRACT(EPOCH FROM measured_at) / 86400.0 − 18262`; v1.21.0 migration
+  /// 0193) so `sum_xx` stays O(1e7) instead of O(1e10) and never sheds
+  /// precision squaring a ~20 540 raw epoch-day. Slope / r² / sd are
+  /// invariant under the x-shift, so the live (raw-epoch) REGR_* probe still
+  /// parity-matches the cross-bucket compose. Nullable to keep migration 0190
   /// additive — the boot backfill fills NULLs for pre-migration rows; new
   /// rows always carry them (they ride the same aggregate as mean / sum).
   sumX        Float?            @map("sum_x")

src/lib/rollups/measurement-read.ts

@@ -53,6 +53,34 @@ export interface DailyMeanRow {
   maxValue: number;
 }
 
+/**
+ * v1.21.0 — regression x-origin (epoch-days of 2020-01-01 UTC).
+ *
+ * The regression accumulators (`sum_x / sum_xy / sum_xx`) store x RELATIVE to
+ * this fixed origin: `x = EXTRACT(EPOCH FROM measured_at)/86400.0 − ORIGIN`.
+ * Raw epoch-days sit at x ≈ 20 540, so `x²` lands near 4.2e8 and `Σxx`
+ * accumulates past ~1e10 — squaring shaves ~10 of double's ~15-16 significant
+ * decimal digits BEFORE the value is ever stored, so the sub-day fractional x
+ * (the time-of-day) is already lossy at write time. No read-side identity can
+ * recover bits lost at accumulation. Rebasing to a recent origin keeps x in the
+ * low thousands (x² ≤ ~1e7 for any realistic window), so the squared terms stay
+ * comfortably exact.
+ *
+ * Slope / r² / population-sd are INVARIANT under an affine x-shift (subtracting
+ * a constant from every x leaves Sxx, Sxy, Syy unchanged), so the composed
+ * result is identical to the un-shifted basis — only the unreported intercept
+ * would move. The writer SQL and migration 0193 share this exact literal; the
+ * live REGR_* probe stays on raw epoch-days and matches because slope / r² / sd
+ * do not depend on the origin.
+ *
+ * 2020-01-01 = 18 262 days after the Unix epoch. It MUST stay a stable
+ * compile-time constant — changing it would orphan every stored accumulator
+ * until a full re-fold. Pick a value comfortably before any plausible reading.
+ */
+export const REGRESSION_X_ORIGIN_DAYS = Math.floor(
+  Date.UTC(2020, 0, 1) / 86_400_000,
+); // 18262
+
 /**
  * v1.20.0 F6 — per-bucket OLS regression accumulators (epoch-day x-axis).
  * Mirrors the four columns migration 0190 added to `measurement_rollups`
@@ -112,14 +140,18 @@ export interface ComposedRegression {
  *
  * This is algebraically identical to the determinant form (multiply numerator
  * and denominator by n) and to Postgres' REGR_* / STDDEV_POP, but numerically
- * stable: on the epoch-day x-axis (x ≈ 20 540) the determinant form computes
- * `n·Σxx − Σx²` as a difference of two ~1e10 magnitudes and sheds ~10
- * significant digits to catastrophic cancellation before the divide. The
- * centered form subtracts `Σx²/n` from `Σxx` at the same scale per term, so
- * the cancellation is bounded by the true x-variance, not the absolute x
- * magnitude. Σy is read from the exact stored `sumValue` accumulator when
- * present (falling back to `mean·count`), removing the AVG→multiply ULP
- * residual on the y side too.
+ * stable: the centered form subtracts `Σx²/n` from `Σxx` at the same scale per
+ * term, so the cancellation is bounded by the true x-variance, not the absolute
+ * x magnitude.
+ *
+ * v1.21.0 — the accumulators are stored REBASED to `REGRESSION_X_ORIGIN_DAYS`
+ * (x = epoch_days − origin), so `Σxx` stays O(1e7) instead of O(1e10) and the
+ * squared terms never shed precision at write time. Slope / r² / sd are
+ * invariant under the x-shift, so this composes the SAME regression the live
+ * (un-rebased) REGR_* probe folds — the rebase only restores the bits the raw
+ * epoch-day square would have lost. Σy is read from the exact stored `sumValue`
+ * accumulator when present (falling back to `mean·count`), removing the
+ * AVG→multiply ULP residual on the y side too.
  *
  * The caller MUST collapse overlapping sources to the canonical source
  * BEFORE handing rows here — the accumulators are per-source, so summing a

src/lib/rollups/measurement-rollups.ts

@@ -34,6 +34,7 @@ import { isP2002 } from "@/lib/prisma-errors";
 import {
   collapseRollupRowsBySource,
   loadUserSourcePriority,
+  REGRESSION_X_ORIGIN_DAYS,
 } from "@/lib/rollups/measurement-read";
 import { startOfUtcDay } from "@/lib/tz/start-of-utc-day";
 import type {
@@ -59,6 +60,30 @@ export const ROLLUP_RECOMPUTE_QUEUE = "rollup-recompute";
 /** Worker concurrency cap — recomputes are read-heavy on Postgres. */
 export const ROLLUP_RECOMPUTE_CONCURRENCY = 2;
 
+/**
+ * v1.21.0 — rebased epoch-day x for the regression accumulators.
+ *
+ * `sum_x / sum_xy / sum_xx` store x RELATIVE to `REGRESSION_X_ORIGIN_DAYS`
+ * (epoch-days of 2020-01-01) so the squared terms stay O(1e7) instead of the
+ * O(1e10) the raw epoch-day (x ≈ 20 540) basis produced — squaring a ~20 540 x
+ * sheds ~10 of double's significant decimal digits BEFORE the sum is stored,
+ * losing the sub-day x detail at write time. Slope / r² / sd are invariant under
+ * the x-shift, so the live REGR_* probe (which stays on raw epoch-days) still
+ * matches the cross-bucket compose. The `slope` / `r2` per-bucket columns also
+ * stay on raw epoch-days (shift-invariant; their stored values are unchanged).
+ *
+ * `REGRESSION_X_ORIGIN_DAYS` is a trusted compile-time integer; assert it stays
+ * one at SQL-build time before splicing it in (no user-input surface). Built
+ * lazily inside `runRollupAggregate` rather than at module init so a partial
+ * test mock of `measurement-read` can't trip the assertion before the SQL runs.
+ */
+function rebasedX(): string {
+  if (!Number.isInteger(REGRESSION_X_ORIGIN_DAYS)) {
+    throw new Error("REGRESSION_X_ORIGIN_DAYS must be an integer");
+  }
+  return `(EXTRACT(EPOCH FROM m."measured_at") / 86400.0 - ${REGRESSION_X_ORIGIN_DAYS})`;
+}
+
 /**
  * Payload `boss.send` carries onto the `rollup-recompute` queue. The
  * worker hands the values straight to `recomputeUserRollups`.
@@ -368,6 +393,10 @@ async function runRollupAggregate(input: {
   }
   const dateTrunc = `date_trunc('${truncUnit}', m."measured_at")`;
 
+  // v1.21.0 — rebased epoch-day x for the regression accumulators (see the
+  // `rebasedX` doc). Built per-call so the assertion fires with the real module.
+  const REBASED_X = rebasedX();
+
   if (input.types && input.types.length > 0) {
     // Build the type-list literal carefully. The enum values come
     // from `MeasurementType` (Prisma-generated TS enum); we whitelist
@@ -401,20 +430,18 @@ async function runRollupAggregate(input: {
           m."value",
           EXTRACT(EPOCH FROM m."measured_at") / 86400.0
         )::double precision                                       AS r2,
-        -- v1.20.0 F6 — regression accumulators (epoch-day x-axis). Summed
-        -- per (type, source, bucket) here so the read tier composes the
-        -- windowed slope / r2 / sd across DAY buckets without a live scan.
-        -- value is NOT NULL on measurements, so every row contributes; the
-        -- sums are over the same rows REGR_/STDDEV_POP fold above.
-        SUM(EXTRACT(EPOCH FROM m."measured_at") / 86400.0)::double precision
-                                                                  AS sum_x,
-        SUM(
-          (EXTRACT(EPOCH FROM m."measured_at") / 86400.0) * m."value"
-        )::double precision                                       AS sum_xy,
-        SUM(
-          (EXTRACT(EPOCH FROM m."measured_at") / 86400.0)
-          * (EXTRACT(EPOCH FROM m."measured_at") / 86400.0)
-        )::double precision                                       AS sum_xx,
+        -- v1.20.0 F6 — regression accumulators. Summed per (type, source,
+        -- bucket) here so the read tier composes the windowed slope / r2 / sd
+        -- across DAY buckets without a live scan. value is NOT NULL on
+        -- measurements, so every row contributes; the sums are over the same
+        -- rows REGR_/STDDEV_POP fold above.
+        -- v1.21.0 — x is REBASED to REGRESSION_X_ORIGIN_DAYS so sum_xx stays
+        -- O(1e7) and never sheds precision squaring a ~20 540 epoch-day. Slope
+        -- / r² / sd are shift-invariant, so the live (raw-epoch) REGR_* probe
+        -- still parity-matches the cross-bucket compose.
+        SUM(${REBASED_X})::double precision                       AS sum_x,
+        SUM(${REBASED_X} * m."value")::double precision           AS sum_xy,
+        SUM(${REBASED_X} * ${REBASED_X})::double precision        AS sum_xx,
         SUM(m."value" * m."value")::double precision              AS sum_yy
       FROM measurements m
       WHERE m."user_id" = $1
@@ -451,18 +478,14 @@ async function runRollupAggregate(input: {
         m."value",
         EXTRACT(EPOCH FROM m."measured_at") / 86400.0
       )::double precision                                       AS r2,
-      -- v1.20.0 F6 — regression accumulators (epoch-day x-axis). See the
-      -- typed-list query above for the rationale; the read tier composes
-      -- windowed slope / r² / sd from these summed across DAY buckets.
-      SUM(EXTRACT(EPOCH FROM m."measured_at") / 86400.0)::double precision
-                                                                AS sum_x,
-      SUM(
-        (EXTRACT(EPOCH FROM m."measured_at") / 86400.0) * m."value"
-      )::double precision                                       AS sum_xy,
-      SUM(
-        (EXTRACT(EPOCH FROM m."measured_at") / 86400.0)
-        * (EXTRACT(EPOCH FROM m."measured_at") / 86400.0)
-      )::double precision                                       AS sum_xx,
+      -- v1.20.0 F6 — regression accumulators. See the typed-list query above
+      -- for the rationale; the read tier composes windowed slope / r² / sd from
+      -- these summed across DAY buckets.
+      -- v1.21.0 — x is REBASED to REGRESSION_X_ORIGIN_DAYS (see above) so the
+      -- squared terms stay exact; slope / r² / sd are shift-invariant.
+      SUM(${REBASED_X})::double precision                       AS sum_x,
+      SUM(${REBASED_X} * m."value")::double precision           AS sum_xy,
+      SUM(${REBASED_X} * ${REBASED_X})::double precision        AS sum_xx,
       SUM(m."value" * m."value")::double precision              AS sum_yy
     FROM measurements m
     WHERE m."user_id" = $1

tests/integration/rollup-regression-parity.integration.test.ts

@@ -287,30 +287,25 @@ describe("rollup regression accumulators — live parity (v1.20.0 F6)", () => {
     expect(blended.slope).not.toBeCloseTo(live.slope!, 6);
   });
 
-  // TODO(v1.20.1, still quarantined): the DST residual is a test-harness
-  // conditioning artifact, NOT a product bug — production math is covered
-  // bit-identically by the WEIGHT case above (which folds the same raw rows
-  // live REGR reads and passes at the tight 1e-9 relative bound).
+  // v1.21.0 — un-quarantined by the regression x-rescale.
   //
-  // Root cause (per the v1.20.1 perf-parity verification): the original
-  // quarantine attributed the residual to catastrophic cancellation in the
-  // determinant form `n·Σxx − Σx²` on the ~1e10 epoch-day x-axis.
-  // `composeRegression` now evaluates the mean-centered (corrected-sum)
-  // identities (`Sxx = Σxx − Σx²/n`, …), which subtract the x mean at matched
-  // scale and shrink this residual by ~30× (was failing at ~1e-6-class noise;
-  // now ~9e-8). It is still above the 1e-9 relative gauge for THIS deliberately
-  // near-flat 5-point seed: the remaining gap is the difference between the
-  // composed closed form and Postgres' OWN internal REGR accumulator algorithm
-  // on an ill-conditioned (near-zero x-variance) window, not a divergence in
-  // the stored data. Both answers round to the same value at the read tier's
-  // 2–3 dp.
+  // The DST residual was never a data divergence: it was precision lost at
+  // ACCUMULATION time. Migration 0190 stored `sum_xx = Σ(epoch_days²)` over the
+  // raw epoch-day x-axis (x ≈ 20 540), so `sum_xx` accumulated past ~1e10 and
+  // the square shed ~10 of double's ~15-16 significant digits BEFORE the value
+  // was stored — the sub-day x detail was already gone, and no read-side
+  // identity (mean-centering included) could recover it. The earlier
+  // mean-centered fix shrank the residual to ~9e-8 but could not cross 1e-9 on
+  // this deliberately ill-conditioned 5-point near-flat seed.
   //
-  // The real close is the migration-based x-rescale noted for backlog: store /
-  // compose with x offset to a window-local origin (Σx → O(window) instead of
-  // ~20 540) so the products never reach 1e10. Until that lands the test stays
-  // quarantined rather than masked behind a loosened tolerance. See
-  // .planning/v1200-backlog.md.
-  it.skip("DST boundary: the epoch-day x-axis is UTC, so a spring-forward reading parity-matches", async () => {
+  // The write path now REBASES x to `REGRESSION_X_ORIGIN_DAYS` (2020-01-01), so
+  // `sum_xx` stays O(1e7) and the squared terms stay exact. Slope / r² / sd are
+  // invariant under the x-shift, so the live (raw-epoch) REGR_* probe this case
+  // compares against still parity-matches the cross-bucket compose — now within
+  // the tight 1e-9 relative bound. The WEIGHT case above stays bit-identical
+  // (it was already exact; the rebase only restores precision on the
+  // ill-conditioned tail).
+  it("DST boundary: the epoch-day x-axis is UTC, so a spring-forward reading parity-matches", async () => {
     const prisma = getPrismaClient();
     const user = await seedUser(prisma);