Probabilistic and count-aware modeling

dashboards/grafana/promforecast-capacity.json

@@ -453,6 +453,42 @@
           }
         ]
       }
+    },
+    {
+      "type": "timeseries",
+      "id": 10,
+      "title": "Forecast fan chart ($id quantiles)",
+      "description": "Predicted quantiles across the forecast horizon for the selected series. Enable per query with emission.quantiles (e.g. [0.1, 0.5, 0.9, 0.99]); requires the remote_write sink for the forward curve. The central forecast is overlaid for reference.",
+      "datasource": {
+        "type": "prometheus",
+        "uid": "${datasource}"
+      },
+      "gridPos": {
+        "h": 9,
+        "w": 24,
+        "x": 0,
+        "y": 57
+      },
+      "targets": [
+        {
+          "expr": "{__name__=~\".+_forecast_quantile\",group=~\"$group\",id=~\"$id\"}",
+          "legendFormat": "q={{q}} {{instance}}"
+        },
+        {
+          "expr": "{__name__=~\".+_forecast\",group=~\"$group\",id=~\"$id\"}",
+          "legendFormat": "central {{instance}}"
+        }
+      ],
+      "fieldConfig": {
+        "defaults": {
+          "custom": {
+            "fillOpacity": 10,
+            "lineWidth": 1,
+            "showPoints": "never"
+          }
+        },
+        "overrides": []
+      }
     }
   ]
 }

docker-compose.dev.yml

@@ -26,7 +26,7 @@ services:
       - "--retentionPeriod=30d"
       - "--futureRetention=2d"
     ports:
-      - "8428:8428"
+      - "8429:8428"
     volumes:
       - vm-data:/storage
     healthcheck:
@@ -48,7 +48,7 @@ services:
     image: promforecast-http-demo:dev
     container_name: pf-http-demo
     ports:
-      - "8000:8000"
+      - "8001:8000"
     restart: unless-stopped
 
   http-demo-loadgen:
@@ -96,7 +96,7 @@ services:
     container_name: pf-redis
     command: ["redis-server", "--save", "", "--appendonly", "no"]
     ports:
-      - "6379:6379"
+      - "6378:6379"
     healthcheck:
       test: ["CMD", "redis-cli", "ping"]
       interval: 5s
@@ -113,7 +113,7 @@ services:
       - "--storage.tsdb.retention.time=2h"
       - "--web.enable-lifecycle"
     ports:
-      - "9090:9090"
+      - "9080:9090"
     volumes:
       - ./docker/prometheus.yml:/etc/prometheus/prometheus.yml:ro
       - prom-data:/prometheus
@@ -129,7 +129,7 @@ services:
     container_name: pf-forecaster
     command: ["--config", "/etc/promforecast/config.yaml"]
     ports:
-      - "9091:9091"
+      - "9081:9091"
     volumes:
       # Dev-tuned config: short lookback / step / min_points so the first fit
       # happens minutes after `make dev-up` instead of needing days of data.
@@ -156,7 +156,7 @@ services:
       GF_AUTH_ANONYMOUS_ORG_ROLE: Viewer
       GF_USERS_DEFAULT_THEME: dark
     ports:
-      - "3000:3000"
+      - "3001:3000"
     volumes:
       - ./docker/grafana/provisioning:/etc/grafana/provisioning:ro
       - ./dashboards/grafana:/etc/grafana/dashboards:ro

docs/README.md

@@ -31,6 +31,7 @@ Everything the forecaster emits beyond raw forecast curves — derived
 metrics, calibration, explainability, and ad-hoc what-if endpoints.
 
 - [Capacity planning](outputs/capacity.md) — damped trend, time-to-threshold ETAs, growth rate, will-breach predictive boolean, queue burn-down (drain) forecasting
+- [Arbitrary-quantile forecasts](quantile.md) — emit `<id>_forecast_quantile{q}` at any quantile for risk-based alerting and fan charts
 - [SLO error-budget forecasting](outputs/slo.md) — burn rate, exhaustion projection, multi-window alerting
 - [Maintenance window recommendations](outputs/maintenance-windows.md) — top-N lowest-load windows for scheduling reboots, rollouts, and cert rotations
 - [Conformal calibration](outputs/conformal.md) — empirical-coverage band correction

docs/operations/cardinality.md

@@ -40,6 +40,7 @@ Notation:
 | `<id>_forecast_contribution_ratio`          | `S × N_regressors`                     | Only with regressors |
 | `forecast_band_calibration_offset`          | `S × M × L`                            | `accuracy.conformal.enabled: true` |
 | `<id>_forecast_component`                   | `S × M × 4`                            | `emission.decompose: true` (trend/seasonal/level/residual) |
+| `<id>_forecast_quantile`                    | `S × M × N_quantiles`                  | `emission.quantiles` set; add `S × M × N_quantiles` more for the per-quantile `forecast_band_calibration_offset{q}` rows when `accuracy.conformal.enabled: true` (and one extra hold-out fit) |
 | `forecast_drift_score`                      | `S × M`                                | `emission.drift.enabled: true` |
 | Threshold ETA / will_breach                 | `S × M × N_thresholds × (1 + L)`       | Only configured thresholds |
 | Drain ETA (`le_drain`)                      | `S × M × N_drain_thresholds × (1 + L)` | Charged separately from the standard ETA |

docs/quantile.md

@@ -0,0 +1,67 @@
+# Arbitrary-quantile forecasts
+
+The standard forecast emits a central value plus symmetric confidence bands (`<id>_forecast_lower` / `<id>_forecast_upper`).  
+For risk-based alerting (“P[disk full within 24h] > 0.05”) and fan charts, you want predictions at *arbitrary* quantiles instead.
+
+`emission.quantiles` emits a parallel family from the **same single fit**:
+
+```
+<id>_forecast_quantile{q="0.95", id, group, model, horizon}
+```
+
+## Configuring
+
+```yaml
+groups:
+  - name: capacity
+    queries:
+      - id: node_filesystem_avail_bytes
+        promql: node_filesystem_avail_bytes
+        emission:
+          quantiles: [0.1, 0.25, 0.5, 0.75, 0.9, 0.99]
+```
+
+Values must be strictly between 0 and 1. Duplicates are rejected at load time. The existing band family remains for backward compatibility.
+
+## How a quantile is derived
+
+Each requested quantile is read from the appropriate band edge of the level that brackets it:
+
+| Quantile `q` | Reads from |
+|--------------|------------|
+| `q = 0.5`    | central forecast (`yhat`) |
+| `q < 0.5`    | lower edge of level `round((1 − 2q)·100)` |
+| `q > 0.5`    | upper edge of level `round((2q − 1)·100)` |
+
+Quantiles beyond the configured confidence levels are clamped to the outermost band.
+
+## Conformal calibration per quantile
+
+When `accuracy.conformal.enabled: true` is also set, each quantile receives its own empirical-coverage correction. The offset is exposed alongside band-level offsets in the existing gauge, discriminated by the `q` label:
+
+```
+forecast_band_calibration_offset{id, group, model, q="0.95"}
+```
+
+This adds one extra hold-out fit per emitted model per refresh.
+
+## Fan chart & risk-based alerting
+
+The bundled capacity dashboard includes a “Forecast fan chart” panel that overlays all emitted quantiles plus the central forecast.
+
+Risk-based alerts become natural:
+
+```promql
+# Optimistic 10th-percentile still predicts disk fill inside horizon
+node_filesystem_avail_bytes_forecast_quantile{q="0.1", horizon="24h"} < 0
+```
+
+```promql
+# 99th-percentile request rate exceeds capacity
+http_requests_rate_forecast_quantile{q="0.99"} > 5000
+```
+
+## Cardinality
+
+`N_quantiles` rows per (series, model) per opted-in query, plus matching calibration-offset rows when conformal is enabled.  
+`promforecast validate` counts the family.

docs/reference/models.md

@@ -1,8 +1,10 @@
 # Models
 
 promforecast ships a small set of statistical forecasters from the
-[`statsforecast`](https://github.com/Nixtla/statsforecast) library, plus a
-`SeasonalNaive` baseline.
+[`statsforecast`](https://github.com/Nixtla/statsforecast) library, a
+`SeasonalNaive` baseline, and a family of in-tree count-distribution
+models (`Poisson`, `NegativeBinomial`, `INGARCH`) for integer-valued
+metrics.
 
 Configure models via the `models:` array (under `defaults` or per group).  
 - **Explicit mode** (`accuracy.auto_select: false`): every chosen model runs and emits forecasts side-by-side.  
@@ -25,6 +27,9 @@ External models can be added via the `promforecast.models` entry-point group (se
 | **CrostonSBA**      | Syntetos-Boylan bias-corrected Croston (recommended default).                               | Sparse / intermittent demand                  | Continuous series with strong seasonality        | Negligible                       |
 | **ADIDA**           | Aggregate-Disaggregate Intermittent Demand Approach.                                        | Very sparse signals where Croston struggles   | Continuous series with strong seasonality        | Negligible                       |
 | **IMAPA**           | Intermittent Multiple Aggregation Prediction Algorithm.                                     | Noisy intermittent signals                    | Continuous series with strong seasonality        | Negligible                       |
+| **Poisson**         | Seasonal-rate Poisson model; Poisson quantile bands (non-negative, asymmetric).             | Integer counts where variance ≈ mean          | Overdispersed counts; continuous series          | Negligible                       |
+| **NegativeBinomial**| Seasonal-rate model with negative-binomial bands; overdispersion from method of moments.    | Integer counts where variance > mean (common) | Continuous series                                | Negligible                       |
+| **INGARCH**         | Integer-valued GARCH(1,1) — Poisson autoregression `λ_t = ω + α·y_{t-1} + β·λ_{t-1}` by MLE. | Autocorrelated / bursty count series          | Series with no count autocorrelation             | Sub-second per series            |
 
 ### When to choose damped
 
@@ -85,6 +90,61 @@ emission:
 
 Applies to both `/metrics` and `remote_write`.
 
+These safeguards are a *post-processing* fix for a continuous model fit to
+count data — they correct the symptom (negative / fractional output) but
+not the mis-specified band shape. For a proper fit, use a
+count-distribution model.
+
+## Count-distribution models
+
+`Poisson`, `NegativeBinomial`, and `INGARCH` model the integer nature of
+the data directly, so the predicted bands respect the distribution: they
+cannot extend below zero and they are asymmetric where the count
+distribution is. The seasonal rate is the per-cycle-position sample mean
+(the maximum-likelihood rate of a saturated seasonal Poisson); the band
+families differ in how they spread around it.
+
+**Decision rule** — pick by the relationship between the series'
+variance and mean (and whether it is autocorrelated):
+
+| Series characteristic                              | Model              |
+|----------------------------------------------------|--------------------|
+| variance ≈ mean (pure Poisson arrivals)            | `Poisson`          |
+| variance > mean (overdispersed — the common case)  | `NegativeBinomial` |
+| counts show autocorrelation (bursty / self-exciting) | `INGARCH`        |
+
+Turn them on by adding them to `models:` and declaring the count profile on
+the query so auto-select prefers them while keeping the continuous models
+as a safety net:
+
+```yaml
+defaults:
+  models: [AutoARIMA, Poisson, NegativeBinomial, INGARCH]
+  accuracy:
+    auto_select: true            # MASE-scored for count (MAPE is undefined on zeros)
+
+groups:
+  - name: errors
+    queries:
+      - id: app_error_count
+        promql: sum by (instance) (rate(app_errors_total[5m]))
+        data_profile: count       # per-query: prefer the count family, MASE scoring
+```
+
+When `data_profile: count` is set, auto-select scores candidates by MASE
+(MAPE is undefined wherever the actual is zero — common for error counts)
+and prefers the registered count wrappers, falling back to the continuous
+models if none score finitely. The [count-aware clip-and-round
+safeguards](#count-aware-emission-safeguards) remain on by default under
+`data_profile: count` as a fallback for any continuous model that still
+wins the horse race, but they are no longer the primary handling.
+
+Count wrappers ignore `regressors`; `INGARCH` falls back to the
+seasonal-rate Poisson projection when the series is too short for a stable
+three-parameter fit or the optimiser does not converge to a stationary
+solution, so a deployment can never end up with a diverging count
+forecast.
+
 ## Cold-start seasonality borrowing
 
 See [`../inputs/cold-start.md`](../inputs/cold-start.md). Short version: when history < `min_points`, set `cold_start.borrow_from.promql` to synthesise a `SeasonalNaive`-style forecast from peer siblings. Emits with `model="ColdStart"` and `cold_start="true"`. Quality score is capped (default 0.5).

forecaster/pyproject.toml

@@ -37,6 +37,11 @@ dependencies = [
     #     errors at runtime" UX when the extras aren't installed.
     "ruptures>=1.1",
     "statsmodels>=0.14",
+    # ``scipy``: count-distribution quantile functions (``scipy.stats``)
+    # and the INGARCH maximum-likelihood fit (``scipy.optimize``). Already
+    # present transitively via statsforecast/statsmodels; pinned directly
+    # because ``promforecast.models.count_distribution`` imports it.
+    "scipy>=1.13",
 ]
 
 [project.optional-dependencies]
@@ -86,6 +91,12 @@ croston            = "promforecast.models.intermittent:CrostonModel"
 croston_sba        = "promforecast.models.intermittent:CrostonSBAModel"
 adida              = "promforecast.models.intermittent:ADIDAModel"
 imapa              = "promforecast.models.intermittent:IMAPAModel"
+# Count-distribution models for integer-valued metrics (request counts,
+# queue depths, error counts). Bands respect the count distribution
+# (non-negative, asymmetric) rather than relying on clip-and-round.
+poisson            = "promforecast.models.count_distribution:PoissonModel"
+negative_binomial  = "promforecast.models.count_distribution:NegativeBinomialModel"
+ingarch            = "promforecast.models.count_distribution:INGARCHModel"
 
 [project.urls]
 Homepage      = "https://github.com/esops-dev/promforecast"
@@ -170,6 +181,16 @@ ignore_missing_imports = true
 module = ["ruptures", "statsmodels.tsa.seasonal"]
 ignore_missing_imports = true
 
+[[tool.mypy.overrides]]
+# ``scipy`` ships partial type information but not for the
+# ``scipy.stats`` distribution objects or ``scipy.optimize.minimize``
+# result attributes the count-distribution models reach into. Those
+# imports are confined to ``promforecast.models.count_distribution``
+# behind a narrow typed surface, so silencing missing imports here keeps
+# mypy --strict clean without chasing incomplete third-party stubs.
+module = ["scipy", "scipy.*"]
+ignore_missing_imports = true
+
 [[tool.mypy.overrides]]
 # The ``holidays`` package (optional extras dependency for the holiday
 # regressor) ships partial stubs

forecaster/src/promforecast/config.py

@@ -688,7 +688,7 @@ class BandCoverageConfig(BaseModel):
 class QueryEmissionConfig(BaseModel):
     """Per-query emission toggles.
 
-    Today six families of settings live here:
+    Several families of settings live here:
 
     * :class:`GrowthRateEmissionConfig` — forecast growth-rate metric.
     * ``clip_non_negative`` and ``round_to_integer`` — count-aware
@@ -699,6 +699,10 @@ class QueryEmissionConfig(BaseModel):
       gauges.
     * ``conditional`` — parallel ``<id>_forecast_conditional*`` family
       that re-runs the fit with regressors held at their current value.
+    * ``decompose`` — ``<id>_forecast_component`` trend/seasonal/level/
+      residual split of the predicted curve.
+    * ``quantiles`` — parallel ``<id>_forecast_quantile{q}`` family at
+      arbitrary quantiles, derived from the same fit.
 
     Every family is kept on the per-query :class:`EmissionConfig` block so
     operators opt in *per metric*: a config that mixes a count-valued
@@ -735,6 +739,19 @@ class QueryEmissionConfig(BaseModel):
     # per (series, model) per horizon step. See `models.py` for which models
     # actually decompose; non-decomposable models bump a skip counter.
     decompose: bool = False
+    # Arbitrary-quantile emission. Empty (default) emits nothing extra; a
+    # list like ``[0.1, 0.5, 0.9, 0.99]`` emits a parallel
+    # ``<id>_forecast_quantile{q="0.9"}`` family alongside the
+    # ``_forecast_lower``/``_forecast_upper`` band family (which stays for
+    # backward compatibility). Each quantile is derived from the same
+    # single fit — no re-fit per quantile — by mapping it to the
+    # statsforecast confidence level that brackets it; ``0.5`` is the
+    # central forecast. Values must be strictly between 0 and 1. When
+    # ``accuracy.conformal.enabled`` is on, a per-quantile coverage
+    # calibration widens each quantile and the offset is exposed as
+    # ``forecast_band_calibration_offset{q="..."}``. See
+    # ``docs/quantile.md``.
+    quantiles: list[float] = Field(default_factory=list)
 
 
 class ColdStartBorrowConfig(BaseModel):
@@ -2077,10 +2094,38 @@ def _validate_query_overrides(*, group_name: str, query: QueryConfig) -> None:
                 )
     if query.refresh_interval is not None and query.refresh_interval.total_seconds() <= 0:
         raise ValueError(f"{label}.refresh_interval must be a positive duration")
+    _validate_quantiles(label=label, query=query)
     _validate_slo(label=label, query=query)
     _validate_recommended_windows(label=label, query=query)
 
 
+def _validate_quantiles(*, label: str, query: QueryConfig) -> None:
+    """Semantic checks for ``emission.quantiles``.
+
+    Pydantic accepts any list of floats; the meaningful constraints
+    (strictly inside ``(0, 1)``, no duplicates) live here so a bad
+    quantile fails at load time with a precise label rather than emitting
+    a degenerate ``q="0"`` / ``q="1"`` row (the 0th/100th percentile maps
+    to an infinite band edge) every refresh.
+    """
+    quantiles = query.emission.quantiles
+    if not quantiles:
+        return
+    seen: set[float] = set()
+    for q in quantiles:
+        if not (0.0 < q < 1.0):
+            raise ValueError(
+                f"{label}.emission.quantiles: every quantile must be strictly "
+                f"between 0 and 1 (got {q}); the 0th/100th percentile maps to an "
+                "infinite band edge"
+            )
+        if q in seen:
+            raise ValueError(
+                f"{label}.emission.quantiles: duplicate quantile {q}; each value must be unique"
+            )
+        seen.add(q)
+
+
 def _validate_slo(*, label: str, query: QueryConfig) -> None:
     """Semantic checks for the SLO block of a query.
 
@@ -2194,6 +2239,8 @@ def _validate_events_profile(*, label: str, query: QueryConfig) -> None:
         rejected.append("emission.decompose")
     if query.emission.growth_rate.enabled:
         rejected.append("emission.growth_rate")
+    if query.emission.quantiles:
+        rejected.append("emission.quantiles")
     if rejected:
         raise ValueError(
             f"{label}: data_profile: events does not support "