Probabilistische voorspellingen scoren: CRPS, PIT-kalibratie en DeepAR

Deze blog heeft al gepleit voor het voorspellen van distributies in plaats van punten, en heeft al uitgewerkt wat je met het interval moet doen als je het eenmaal hebt: [conforme voorspelling voor risicobewuste positiebepaling] (/en/blog/post/conformal-prediction-trading) leidt distributievrije intervallen af en de maatregel die deze verbruikt, en het artikel [Temporal Fusion Transformer] (/en/blog/post/temporal-fusion-transformer-trading) verzendt een multi-kwantiel uitvoerlaag. Wat geen van beide behandelt, is het deel dat beslist of iets ervan de moeite waard is om te vertrouwen: hoe je een voorspellende verdeling scoort, en hoe je controleert of de aangegeven onzekerheid eerlijk is.

Dat is waar dit artikel over gaat. Drie dingen specifiek:

1. CRPS — de juiste scoreregel voor probabilistische voorspellingen, en de exacte relatie met het flipperkastverlies dat in het TFT-artikel al wordt gerapporteerd.
2. Het PIT-histogram — een kalibratiediagnostiek die afleest hoe een model verkeerd is gekalibreerd, niet alleen of dat zo is.
3. DeepAR — de familie van autoregressieve bemonsteringsmodellen, die alleen op deze blog verschijnt als benchmarkvoetnoot en nooit wordt uitgelegd.

Eerst één frame, want dat bepaalt welke machines je nodig hebt. Een puntvoorspelling mislukt per regime verschillend: in een trend met lage volatiliteit is de voorwaardelijke verdeling smal en vrijwel symmetrisch en is een puntvoorspelling een goede samenvatting; vóór een gepland evenement is het bimodaal, en het voorwaardelijke gemiddelde bevindt zich precies daar waar de prijs het minst waarschijnlijk zal landen; in een crisis domineert de linkerstaart en onderschat het gemiddelde de keerzijde slecht. Drie routes herstellen de volledige verdeling: parametrisch (voorspel de parameters van een veronderstelde familie: snel, risico op verkeerde specificatie), kwantielgebaseerd (voorspel een vast raster: aannamevrij, discreet) en steekproefgebaseerd (Monte Carlo-paden van autoregressieve bemonstering, uitval of ensembles: flexibel, duur). De laatste twee domineren in de financiële wereld, juist omdat de verdeling tussen deze regimes van vorm verandert.

Kwantielvoorspellingen, in het kort

Het flipperkastverlies en een concreet kwantielraster worden al gegeven in de [TFT-kwantieluitvoerlaag] (/en/blog/post/temporal-fusion-transformer-trading), dus ik zal de formule niet opnieuw formuleren. Het enige dat de moeite waard is om te internaliseren, is de asymmetrie ervan: at $\tau = 0.5$ over- en ondervoorspelling kosten hetzelfde en het verlies vermindert tot MAE, maar dan $\tau = 0.95$ ondervoorspelling wordt bestraft $19\times$ meer dan overvoorspelling. Die verhouding $\tau/(1-\tau)$ is het mechanisme – het is wat de gepaste waarde naar een niveau sleept dat slechts 5% van de waarnemingen zou moeten overschrijden.

De praktische rimpel die het gepubliceerde artikel niet vermeldt, is kwantieloverschrijding: niets in een per-kwantielverlies voorkomt $\hat{q}_{0.25} > \hat{q}_{0.75}$. Met voldoende gegevens en een gedeelde ruggengraat is dit zeldzaam, maar het gebeurt wel bij de extreme kwantielen op dunne monsters, en het verbreekt stilletjes elke stroomafwaartse CRPS- of dekkingsberekening. De goedkope oplossing is een post-hoc soort van de voorspelde kwantielvector; de principiële is een monotone outputparametrisering (predict $q_{\min}$ plus niet-negatieve stappen).

Een minimale implementatie vanaf het begin, handig als u kwantitatieve output wilt zonder een voorspellingsraamwerk te gebruiken:

import torch
import torch.nn as nn

class QuantileRegressionNet(nn.Module):
    """Multi-quantile forecasting network for financial returns."""

    def __init__(self, input_dim: int, hidden_dim: int = 128,
                 quantiles: list[float] = [0.01, 0.05, 0.1, 0.25, 0.5, 0.75, 0.9, 0.95, 0.99]):
        super().__init__()
        self.quantiles = quantiles
        self.backbone = nn.Sequential(
            nn.Linear(input_dim, hidden_dim), nn.ReLU(), nn.Dropout(0.2),
            nn.Linear(hidden_dim, hidden_dim), nn.ReLU(), nn.Dropout(0.2),
        )
        self.heads = nn.ModuleList([nn.Linear(hidden_dim, 1) for _ in quantiles])

    def forward(self, x: torch.Tensor) -> torch.Tensor:
        h = self.backbone(x)
        out = torch.cat([head(h) for head in self.heads], dim=-1)
        return torch.sort(out, dim=-1).values


def pinball_loss(predictions: torch.Tensor, targets: torch.Tensor,
                 quantiles: list[float]) -> torch.Tensor:
    """predictions: (batch, n_quantiles); targets: (batch, 1)."""
    errors = targets - predictions
    tau = torch.tensor(quantiles, device=predictions.device).unsqueeze(0)
    return torch.max(tau * errors, (tau - 1) * errors).mean()

De intervalbreedte wordt gekoppeld aan de positiegrootte, en we werken die mapping volledig uit – inclusief het edge-ratio no-trade filter – in conformele voorspelling voor risicobewuste positiebepaling; de $1/\sigma$ gelijkwaardigheid met targeting op volatiliteit wordt afgeleid uit targeting op volatiliteit met GARCH-voorspellingen.

DeepAR: autoregressieve probabilistische prognoses

DeepAR (Salinas et al., 2020) volgt de op monsters gebaseerde route. In plaats van kwantielen rechtstreeks te voorspellen, gebruikt het een autoregressieve RNN om bij elke stap een waarschijnlijkheid te parametriseren, en trekt vervolgens volledige voorspellende verdelingen door die waarschijnlijkheid naar voren te rollen.

Architectuur. Bij elke stap $t$ het netwerk verbruikt de vorige observatie $z_{t-1}$ (geschaald met een factor per reeks), covariaten $\mathbf{x}_t$en optionele statische functies. Een LSTM draagt ​​de staat:

$\mathbf{h}_t = \text{LSTM}(\mathbf{h}_{t-1}, [z_{t-1}, \mathbf{x}_t])$

Een dicht hoofd komt in kaart $\mathbf{h}_t$ waarschijnlijkheidsparameters $\theta_t = \text{MLP}(\mathbf{h}_t)$ — $(\mu_t, \sigma_t)$ voor Gaussiaans, $(\mu_t, \sigma_t, \nu_t)$ voor Student-t. Trainen maximaliseert $\sum_i \sum_t \log p(z_{i,t} \mid \theta_{i,t})$ over alle series. Bij gevolgtrekking monster je uit $p(\cdot \mid \theta_t)$, voer het monster terug als de volgende invoer en herhaal $S$ tijden te krijgen $S$ trajecten.

Voor de handel zijn twee gevolgen van belang. De waarschijnlijkheid is een modelkeuze, dus dikke staarten zijn iets dat je selecteert in plaats van iets waar je op hoopt — Student-t voor zware staarten, of een Gaussiaans mengsel $\sum_k w_k \mathcal{N}(z_t; \mu_k(\mathbf{h}_t), \sigma_k^2(\mathbf{h}_t))$ voor bimodaal gebeurtenisgedrag. Voorspellingen in meerdere stappen zijn coherent: elk monsterpad is een plausibel traject dat de seriële correlatie behoudt, wat u nodig heeft voor elke horizon langer dan één stap. (Het derde gebruikelijke verkoopargument: één mondiaal model over vele series $N$ per-asset modellen – is hetzelfde data-efficiëntie argument dat wordt aangevoerd in het TFT artikel, waarbij schaling per serie en statische covariaten daar de rol spelen van de statische encoders.)

DeepAR met GluonTS

import pandas as pd
from gluonts.dataset.pandas import PandasDataset
from gluonts.torch.model.deepar import DeepAREstimator
from gluonts.torch.distributions import StudentTOutput
from gluonts.evaluation import make_evaluation_predictions, Evaluator


def prepare_crypto_dataset(returns_df: pd.DataFrame, freq: str = "h"):
    """Wide return frame (index=datetime, columns=assets) -> GluonTS dataset.

    from_long_dataframe wants LONG format: one row per (timestamp, item_id)
    with the target in a column. Passing a DataFrame of dicts does not work.
    """
    long_df = (
        returns_df.stack()
        .rename("target")
        .rename_axis(index=["timestamp", "item_id"])
        .reset_index()
        .dropna(subset=["target"])
    )
    return PandasDataset.from_long_dataframe(
        long_df, target="target", item_id="item_id",
        timestamp="timestamp", freq=freq,
    )


estimator = DeepAREstimator(
    prediction_length=24,    # 24 hours ahead
    context_length=168,      # one week of hourly data
    freq="h",
    num_layers=2,
    hidden_size=64,
    dropout_rate=0.1,
    lr=1e-3,
    batch_size=64,
    trainer_kwargs={"max_epochs": 50},
    distr_output=StudentTOutput(),
)

predictor = estimator.train(training_data=train_dataset)

forecast_it, ts_it = make_evaluation_predictions(
    dataset=test_dataset, predictor=predictor, num_samples=500,
)
forecasts, actuals = list(forecast_it), list(ts_it)

evaluator = Evaluator(quantiles=[0.05, 0.25, 0.5, 0.75, 0.95])
agg_metrics, item_metrics = evaluator(actuals, forecasts)
print(f"mean_wQuantileLoss: {agg_metrics['mean_wQuantileLoss']:.4f}")

num_samples regelt het aantal Monte Carlo-paden. Voor live gebruik is 100-200 meestal voldoende; gebruik voor offline evaluatie 500-1000, omdat de CRPS-schatting op basis van te weinig monsters vertekenend laag is.

Andere routes naar een voorspellende verdeling

Drie alternatieven zijn de moeite waard om te kennen, maar hebben geen eigen secties nodig. MC Dropout (Gal en Ghahramani, 2016) houdt de drop-out actief bij gevolgtrekkingen en neemt het gemiddelde en de variantie over $S$ voorwaartse passen; het is bij benadering een vermomde variatie-inferentie, en het is de snelste manier om onzekerheid toe te passen op een model dat je al hebt getraind. Diepe ensembles (Lakshminarayanan et al., 2017) trainen $M$ kopieën van verschillende zaden en behandel de voorspellende distributie als een mengsel – consistent sterk in benchmarks, op $M\times$ de kosten, die dit uitsluiten voor latentiegevoelige horizonten, maar niet voor 4 uur of dagelijks. Normaliserende stromen (Rasul et al., 2021) leren een omkeerbare kaart van een eenvoudige basisdichtheid tot een willekeurig doel, waarbij multimodaliteit en asymmetrie worden vastgelegd zonder een parametrische familie te kiezen. Alle drie produceren monsters, dus alles in de volgende sectie is ongewijzigd op hen van toepassing.

CRPS: de juiste maatstaf voor een voorspellende verdeling

MSE en MAE scoren puntvoorspellingen. Ze kunnen je niet vertellen of een distributie goed was, omdat ze altijd maar naar één samenvatting ervan kijken. De standaardvervanger is de Continuous Rated Probability Score, en dit is het meest bruikbare in dit artikel.

CRPS is de geïntegreerde kwadratische afstand tussen de voorspelde CDF en de gedegenereerde CDF die alle massa legt op wat er feitelijk is gebeurd:

$\text{CRPS}(F, y) = \int_{-\infty}^{\infty} \left(F(z) - \mathbb{1}[y \leq z]\right)^2 dz$

Drie eigendommen verdienen het zijn plaats:

• Het is een goede scoreregel (Gneiting en Raftery, 2007): alleen in verwachting geminimaliseerd wanneer de voorspelde verdeling gelijk is aan de werkelijke verdeling. Noch opzettelijke overmoed, noch hedging met een kunstmatig brede distributie verbetert uw score. MSE-op-de-mediaan heeft deze eigenschap niet en daarom is CRPS niet optioneel.
• Het generaliseert MAE. Voor een gedegenereerde puntvoorspelling stort het in op een absolute fout, dus CRPS leeft in de eenheden van het doel — bij het voorspellen van logrendementen per uur is een CRPS van 0,004 direct vergelijkbaar met een gemiddelde absolute fout van 40 bps, in plaats van een eenheidloos getal dat je niet kunt controleren op gezond verstand.
• Het beloont scherpte na kalibratie. Tussen twee gelijk gekalibreerde voorspellingen scoort de smallere beter. Dat is ook de reden dat CRPS alleen niet voldoende is: een slechte score vertelt je niet welke van de twee voorwaarden niet is gelukt, daar is het volgende onderdeel voor.

De brug naar kwantielverlies

Dit is de verbinding die de rest van de blog mist. Gegeven een kwantielraster $\mathcal{T}$ in plaats van een volledige CDF, wordt CRPS benaderd door het flipperkastverlies:

$\text{CRPS} \approx \frac{2}{|\mathcal{T}|} \sum_{\tau \in \mathcal{T}} \mathcal{L}_\tau(y, \hat{q}_\tau)$

Lees dat letterlijk: het "kwantielverlies" gerapporteerd in het TFT-artikel en de hier besproken CRPS zijn dezelfde hoeveelheid tot een factor 2, waarbij de benadering strakker wordt naarmate het raster dichter wordt. Het zijn geen concurrerende statistieken en er is geen reden om beide te rapporteren.

Eén waarschuwing voor de eenheid, omdat het bijt. GluonTS's mean_wQuantileLoss is datzelfde gemiddelde flipperkastverlies genormaliseerd door de som van absolute doelwaarden, waardoor het dimensieloos en vergelijkbaar is voor activa van verschillende schaal. CRPS zelf is niet genormaliseerd en blijft in retoureenheden. Niet afdrukken mean_wQuantileLoss onder een label met de tekst "CRPS" - de conceptversie van dit artikel deed precies dat, en het is een gemakkelijke manier om twee getallen te vergelijken die niet op dezelfde schaal liggen.

CRPS uit monsters

Met Monte Carlo-monsters $\{y^{(s)}\}_{s=1}^{S}$ (DeepAR, ensembles, stromen), gebruik de energievorm:

$\text{CRPS}(F, y) = \frac{1}{S} \sum_{s=1}^{S} |y^{(s)} - y| - \frac{1}{2S^2} \sum_{s=1}^{S} \sum_{s'=1}^{S} |y^{(s)} - y^{(s')}|$

De eerste term beloont nauwkeurigheid, de tweede bestraft overspreiding. Naïef geschreven is de tweede term $O(S^2)$ en wordt het knelpunt wanneer u duizenden voorspellingen scoort. Sorteren vouwt het eerst samen: voor statistieken met oplopende volgorde $y_{(1)} \le \dots \le y_{(S)}$, de dubbele som is gelijk $2\sum_k (2k - S - 1)\, y_{(k)}$, zo is het allemaal $O(S \log S)$ gedomineerd door de soort.

import numpy as np

def crps_empirical(samples: np.ndarray, observation: float) -> float:
    """CRPS from Monte Carlo samples, O(n log n) via order statistics."""
    n = len(samples)
    mae = np.mean(np.abs(samples - observation))
    x = np.sort(samples)
    k = np.arange(1, n + 1)
    dispersion = np.sum((2 * k - n - 1) * x) / n**2
    return mae - dispersion


def crps_quantile(quantile_predictions: np.ndarray,
                  quantile_levels: np.ndarray,
                  observation: float) -> float:
    """CRPS approximation from a quantile grid (2x mean pinball loss)."""
    errors = observation - quantile_predictions
    pinball = np.where(errors >= 0,
                       quantile_levels * errors,
                       (quantile_levels - 1) * errors)
    return 2.0 * np.mean(pinball)

Beide vormen moeten het nauw eens zijn over dezelfde voorspellende verdeling; als dat niet het geval is, vermoed dan dat er sprake is van kwantieloverschrijding of van te weinig monsters. In productie, properscoring.crps_ensemble(observation, samples) is een beproefde drop-in.

Kalibratie: is de aangegeven onzekerheid eerlijk?

Een goed CRPS garandeert niet dat de intervallen betekenen wat ze beweren. Een model waarvan het nominale interval van 90% 70% van de uitkomsten vastlegt, is overmoedig, en in een maatregel die de intervalbreedte aangeeft, zal deze u precies op de hoogte stellen wanneer dat niet zou moeten. Het TFT-artikel vermeldt de waarschuwing en schrijft conforme voorspelling voor als oplossing; wat volgt is hoe u de storing daadwerkelijk detecteert en diagnosticeert.

Het PIT-histogram

Voor elke waarneming $y_t$, bereken het kwantiel onder de voorspelde CDF van die stap:

$u_t = \hat{F}_t(y_t)$

Als het model is gekalibreerd, $\{u_t\}$ is uniform aan $[0,1]$. Het diagnostische vermogen is dat de vorm van de afwijking de faalmodus aangeeft:

• U-vormig: overmoedig – er komt te veel massa in de staarten terecht, de verdeling is te smal.
• Bultvormig: weinig zelfvertrouwen – waarnemingen clusteren nabij het midden, de verdeling is te breed.
• Linksscheef: het model voorspelt systematisch te veel.
• Rechts scheef: het model voorspelt systematisch te laag.

Een opmerking om verwarring te voorkomen met copula-modellen voor gezamenlijk risico, die ook gebruik maakt van de waarschijnlijkheidsintegraaltransformatie: daar is het een marginale transformatie, een voorverwerkingsstap die GARCH-EVT-marginaliteiten omzet in pseudo-uniforme observaties, zodat er een copula aan kan worden aangepast. Hier wordt de transformatie toegepast op voorspellingen die buiten de steekproef vallen en is de uniformiteit het resultaat dat wordt getest, en niet de input die wordt geproduceerd. Dezelfde wiskunde, tegengestelde richting van gevolgtrekking.

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import kstest

def pit_calibration_check(forecasts, actuals, n_bins=20):
    """PIT histogram + KS test against uniform.

    forecasts: list of SampleForecast objects (GluonTS)
    """
    pit_values = []
    for forecast, actual in zip(forecasts, actuals):
        samples = forecast.samples          # (n_samples, prediction_length)
        h = forecast.prediction_length
        for t in range(h):
            obs = actual.values[-h + t]
            pit_values.append(np.mean(samples[:, t] <= obs))

    pit_values = np.array(pit_values)
    ks_stat, p_value = kstest(pit_values, "uniform")

    fig, ax = plt.subplots(figsize=(8, 4))
    ax.hist(pit_values, bins=n_bins, density=True, alpha=0.7, edgecolor="black")
    ax.axhline(y=1.0, color="red", linestyle="--", label="Perfect calibration")
    ax.set_xlabel("PIT value")
    ax.set_ylabel("Density")
    ax.set_title(f"PIT Histogram (KS={ks_stat:.3f}, p={p_value:.3f})")
    ax.legend()
    plt.tight_layout()
    return pit_values, ks_stat, p_value

Twee kanttekeningen bij de KS-test. Er wordt uitgegaan van onafhankelijke PIT-waarden, en overlappende voorspellingen over meerdere horizonten zijn sterk autogecorreleerd. Beschouw de p-waarde dus als een ruwe vlag en de histogramvorm als het echte bewijs. En bij voldoende waarnemingen verwerpt de test uniformiteit vanwege een miskalibratie die te klein is om er toe te doen; de effectgrootte is wat een beslissing moet bepalen.

Dekkingscontrole

Een grovere maar directer bruikbare diagnose: doe het $\alpha$% intervallen bevatten $\alpha$% van de resultaten?

import numpy as np
import pandas as pd

def coverage_table(forecasts, actuals, levels=(0.50, 0.80, 0.90, 0.95)):
    """Empirical coverage at multiple nominal levels."""
    results = {}
    for level in levels:
        lower_q = (1 - level) / 2
        upper_q = 1 - lower_q
        covered = total = 0
        for forecast, actual in zip(forecasts, actuals):
            samples = forecast.samples
            h = forecast.prediction_length
            for t in range(h):
                obs = actual.values[-h + t]
                lo = np.quantile(samples[:, t], lower_q)
                hi = np.quantile(samples[:, t], upper_q)
                covered += int(lo <= obs <= hi)
                total += 1
        empirical = covered / total
        results[f"{int(level*100)}% interval"] = {
            "nominal": level, "empirical": empirical,
            "gap": empirical - level,
        }
    return pd.DataFrame(results).T

De vuistregel: gaten onder de 2 procentpunten zijn ruis bij typische steekproefgroottes, en gaten groter dan 3-5 punten zijn een reëel kalibratieprobleem dat tot uiting komt in de positiebepaling. Voer het uit per horizon, niet gebundeld; de dekking neemt bijna altijd af naarmate de horizon zich uitstrekt, en door samenvoegen wordt deze verborgen.

Wanneer kalibratie mislukt

Drie standaardreparaties, in oplopende volgorde van garantie. Temperatuurschaling deelt de voorspelde schaalparameter door een geleerde $T$ aangepast op beschikbare gegevens – één parameter, triviaal goedkoop, corrigeert uniforme over-/onderzekerheid, maar niets vormafhankelijk. Isotone herkalibratie brengt voorspelde kwantielniveaus monotoon in kaart met waargenomen frequenties, waarbij vormvervorming wordt verwerkt, maar een redelijk grote kalibratieset nodig is. Conformele voorspelling omvat elk model en biedt dekkingsgaranties voor eindige steekproeven; het volledige split-conforme algoritme, de exacte volgorde-statistische rangorde en de interpolatie- en klemvallen die een samenvatting van één regel uitnodigt, bevinden zich allemaal in [conformele voorspelling voor handel] (/en/blog/post/conformal-prediction-trading).

Praktische overwegingen

Niet-stationariteit. Kalibratie verandert naarmate de volatiliteitsregimes veranderen, dus een vaste kalibratieset vervalt. Het mechanisme dat precies hiervoor is ontworpen is Adaptive Conformal Inference, dat een online misdekkingsniveau bijwerkt en een dekkingsgarantie op lange termijn biedt, zelfs bij vijandige sequenties - zie de ACI-sectie van [conformal forecasting for trading] (/en/blog/post/conformal-prediction-trading).

Rekeningskosten. DeepAR bij 500 monsterpaden kost ongeveer $500\times$ een enkele voorwaartse pass. Voor intraday-werk geeft u de voorkeur aan kwantielregressie (alle kwantielen in één keer) of een parametrische kop (voorspel $\mu, \sigma$ eenmaal); reserveer autoregressieve bemonstering voor 4 uur en dagelijkse horizonten.

Regimewijzigingen. Koppel de voorspeller aan een regimedetector en houd kalibratieparameters per regime bij — regimedetectie met HMM heeft de werkende implementatie en backtest.

Multivariate prognoses. Marginale uitkeringen per activum zijn niet voldoende voor portefeuillerisico; de gezamenlijke staart is waar het om gaat, en copula-modellen voor gezamenlijk risico kwantificeert hoezeer onafhankelijkheidsaannames dit onderschatten.

Downstream consumenten. VaR en het verwachte tekort komen rechtstreeks voort uit een op steekproeven gebaseerde voorspelling als een kwantiel en een voorwaardelijk staartgemiddelde. De definities en het Monte Carlo-recept staan ​​in copula-modellen voor gezamenlijk risico. Het dimensioneren van Kelly is het enige dat niet gratis is: het afleiden van een Kelly-fractie uit een voorspellingsinterval vereist een extra aanname over de verdeling binnen dat interval, en het conforme artikel pleit expliciet tegen het vastpinnen van een intervalverhouding op $f^*$. Zie het Kelly-criterium voor strategieën voor wat de breuk feitelijk vereist.

Conclusie

De evaluatiestapel voor een probabilistische voorspelling is kort en niet onderhandelbaar: CRPS voor de score, omdat deze correct is en in de eenheden van uw doel leeft; het PIT-histogram voor de diagnose, omdat het de foutmodus benoemt in plaats van er slechts één te markeren; dekking per horizon voor de beslissing, omdat dat het getal is dat overeenkomt met de positiegrootte. Alles wat als ‘kwantielverlies’ wordt gerapporteerd, is CRPS tot een factor 2, dus er is hier sprake van één maatstaf, niet van twee.

Het ongemakkelijke is dat dit allemaal machines zijn om erachter te komen dat een model slechter is dan je had gehoopt. Dat is het punt. Een model dat eerlijk de intervallen vergroot als het het niet weet, is strikt genomen nuttiger dan een model dat zelfverzekerd smal blijft, en de enige manier om ze uit elkaar te houden is door ze op de juiste manier te scoren.

Referenties

• Salinas, D., Flunkert, V., Gasthaus, J., & Januschowski, T. (2020). "DeepAR: probabilistische prognoses met autoregressieve terugkerende netwerken." International Journal of Forecasting, 36(3), 1181-1191.
• Koenker, R. & Bassett, G. (1978). "Regressiekwantielen." Econometrica, 46(1), 33-50.
• Gneiting, T. & Raftery, AE (2007). "Strikt juiste scoreregels, voorspelling en schatting." Journal of the American Statistical Association, 102(477), 359-378.
• Gneiting, T., Balabdaoui, F., & Raftery, AE (2007). "Probabilistische voorspellingen, kalibratie en scherpte." Journal of the Royal Statistical Society: Serie B, 69(2), 243-268.
• Gal, Y. & Ghahramani, Z. (2016). "Uitval als een Bayesiaanse benadering: modelonzekerheid in diep leren vertegenwoordigen." ICML.
• Lakshminarayanan, B., Pritzel, A., en Blundell, C. (2017). "Eenvoudige en schaalbare voorspellende onzekerheidsschatting met behulp van diepe ensembles." NeurIPS.
• Rasul, K., Sheikh, A.-S., Schuster, I., Bergmann, U., & Vollgraf, R. (2021). "Multivariate probabilistische tijdreeksvoorspellingen via geconditioneerde normaliserende stromen." ICLR.
• GluonTS: probabilistische tijdreeksmodellering in Python. https://ts.gluon.ai