Forecasting and Trend Analysis Questions

-- Step 1: compute historical proportions
WITH hist AS (
  SELECT product_id, GREATEST(SUM(revenue), 0) AS h
  FROM sales WHERE period = 'last_year' GROUP BY product_id
),
base AS (
  SELECT p.product_id,
         COALESCE(h,0) AS h,
         COALESCE(override_value, NULL) AS override,
         COALESCE(min_target, 0) AS min_t,
         COALESCE(max_target, 1e12) AS max_t
  FROM products p LEFT JOIN hist h ON p.product_id = h.product_id
)
SELECT * -- further iterative logic handled in app or stored proc
FROM base;

import numpy as np

def ensemble_with_bootstrap(forecasts, val_residuals, n_boot=2000, alpha=0.10, random_state=None):
    """
    forecasts: list or array of shape (M, H) or (N, H) where M = models, H = horizons
    val_residuals: list of arrays, each array shape (T_val, H) or (T_val,) if residuals pooled across horizons
    n_boot: number of bootstrap samples
    alpha: 1 - coverage (0.10 -> 90% PI)
    Returns: dict with 'point' (H,), 'lower' (H,), 'upper' (H,), 'weights' (M,)
    """
    rng = np.random.default_rng(random_state)
    F = np.asarray(forecasts)
    M, H = F.shape
    # compute RMSE per model across validation residuals
    rmses = np.empty(M)
    for i in range(M):
        r = np.asarray(val_residuals[i])
        if r.ndim == 1:
            rmses[i] = np.sqrt(np.mean(r**2))
        else:
            rmses[i] = np.sqrt(np.mean(r**2))
    # guard against zero RMSE
    eps = 1e-8
    inv = 1.0 / (rmses + eps)
    weights = inv / np.sum(inv)
    # point forecast: weighted average across models
    point = weights.reshape(-1, 1) * F
    point = point.sum(axis=0)  # shape (H,)
    # prepare residual pool per model: keep same horizon shape if available
    # For bootstrap, sample a model's residual row then add to ensemble forecast using model weights
    # Precompute residual matrices as (T_i, H)
    residual_mats = []
    for i in range(M):
        r = np.asarray(val_residuals[i])
        if r.ndim == 1:
            # treat as same residual applied to all horizons
            residual_mats.append(r.reshape(-1, 1) * np.ones((1, H)))
        else:
            if r.shape[1] != H:
                raise ValueError("Validation residual horizon mismatch")
            residual_mats.append(r)
    # bootstrap ensemble residual by sampling one residual vector per model (with replacement) each bootstrap
    sims = np.empty((n_boot, H))
    for b in range(n_boot):
        sampled = np.empty((M, H))
        for i in range(M):
            ri = residual_mats[i]
            idx = rng.integers(0, ri.shape[0])
            sampled[i] = ri[idx]
        # ensemble residual = weighted sum of sampled model residuals
        ens_resid = np.dot(weights, sampled)
        sims[b] = point + ens_resid
    lower = np.percentile(sims, 100 * (alpha / 2), axis=0)
    upper = np.percentile(sims, 100 * (1 - alpha / 2), axis=0)
    return {"point": point, "lower": lower, "upper": upper, "weights": weights}