正弦回归

在本笔记本中，我们展示了如何使用Fortuna在正弦回归任务中获得校准的预测不确定性估计。

!pip install -q aws-fortuna

生成数据

我们首先生成数据。这些是一些遵循正弦关系的输入和目标变量，受到噪声的干扰。

import numpy as np
import matplotlib.pyplot as plt


def make_sinusoidal(n_data: int, noise: float = 0.1, seed: int = 0):
    np.random.seed(seed)
    w = np.arange(30)[:, None] / 30
    b = 2 * np.pi * np.arange(30)[:, None] / 30

    x = np.random.normal(size=(n_data,))
    y = np.cos(w * x + b).sum(0) + noise * np.random.normal(size=(n_data,))
    return x[:, None], y[:, None]


train_data = make_sinusoidal(n_data=500, seed=0)
val_data = make_sinusoidal(n_data=500, seed=1)
test_data = make_sinusoidal(n_data=500, seed=2)

fig, axes = plt.subplots(1, 3, figsize=(10, 1))
axes[0].scatter(*train_data, s=1, label="training data", c="C0")
axes[0].legend()
axes[1].scatter(*val_data, s=1, label="validation data", c="C1")
axes[1].legend()
axes[2].scatter(*test_data, s=1, label="test data", c="C2")
axes[2].legend()

<matplotlib.legend.Legend at 0x7f6319421ed0>

将数据转换为兼容的数据加载器

Fortuna 帮助你将数组元组转换为兼容的数据加载器。

from fortuna.data import DataLoader

train_data_loader = DataLoader.from_array_data(
    train_data, batch_size=128, shuffle=True, prefetch=True
)
val_data_loader = DataLoader.from_array_data(val_data, batch_size=128, prefetch=True)
test_data_loader = DataLoader.from_array_data(test_data, batch_size=128, prefetch=True)

构建一个概率回归器

让我们构建一个概率回归器。这是一个包含多个可配置属性的接口对象，即model, likelihood_log_variance_model, prior, posterior_approximator, output_calibrator。在这个例子中，我们使用MLP模型来处理似然的均值和对数方差，使用深度集成后验近似器，以及默认的温度缩放输出校准器。

from fortuna.prob_model import ProbRegressor, DeepEnsemblePosteriorApproximator
from fortuna.model import MLP
import flax.linen as nn

output_dim = 1
prob_model = ProbRegressor(
    model=MLP(output_dim=output_dim, activations=(nn.tanh, nn.tanh)),
    likelihood_log_variance_model=MLP(
        output_dim=output_dim, activations=(nn.tanh, nn.tanh)
    ),
    posterior_approximator=DeepEnsemblePosteriorApproximator(),
)

训练概率模型：后验拟合和校准

我们现在可以训练概率模型。这包括拟合后验分布和校准概率模型。

from fortuna.prob_model import FitConfig, FitMonitor, CalibConfig, CalibMonitor
from fortuna.metric.regression import rmse

status = prob_model.train(
    train_data_loader=train_data_loader,
    val_data_loader=val_data_loader,
    calib_data_loader=val_data_loader,
    fit_config=FitConfig(
        monitor=FitMonitor(early_stopping_patience=2, metrics=(rmse,))
    ),
    calib_config=CalibConfig(monitor=CalibMonitor(early_stopping_patience=2)),
)

Epoch: 100 | loss: 1801.01233 | rmse: 0.34807: 100%|██████████| 100/100 [00:03<00:00, 29.28it/s]
Epoch: 100 | loss: 1750.73816 | rmse: 0.29804: 100%|██████████| 100/100 [00:01<00:00, 80.62it/s]
Epoch: 100 | loss: 1840.17017 | rmse: 0.37487: 100%|██████████| 100/100 [00:01<00:00, 80.49it/s]
Epoch: 100 | loss: 1830.52966 | rmse: 0.35051: 100%|██████████| 100/100 [00:01<00:00, 80.18it/s]
Epoch: 100 | loss: 1798.47156 | rmse: 0.32457: 100%|██████████| 100/100 [00:01<00:00, 81.05it/s]
Epoch: 15 | loss: -183.82822:  14%|█▍        | 14/100 [00:01<00:06, 12.79it/s]

plt.figure(figsize=(6, 3))
for i, s in enumerate(status["fit_status"]):
    plt.plot(s["loss"], label=f"run #{i+1}")
plt.legend()
plt.title("loss decays", fontsize=12)

Text(0.5, 1.0, 'loss decays')

估计预测统计

我们现在可以通过调用概率回归器的predictive属性和感兴趣的方法来计算一些预测统计量。大多数预测统计量，例如均值或方差，需要一个输入数据点的加载器。你可以通过调用数据加载器的方法to_inputs_loader轻松获得这个加载器。

test_log_probs = prob_model.predictive.log_prob(data_loader=test_data_loader)
test_inputs_loader = test_data_loader.to_inputs_loader()
test_means = prob_model.predictive.mean(inputs_loader=test_inputs_loader)
test_aleatoric_variances = prob_model.predictive.aleatoric_variance(
    inputs_loader=test_inputs_loader
)
test_epistemic_variances = prob_model.predictive.epistemic_variance(
    inputs_loader=test_inputs_loader
)
test_variances = prob_model.predictive.variance(
    inputs_loader=test_inputs_loader,
    aleatoric_variances=test_aleatoric_variances,
    epistemic_variances=test_epistemic_variances,
)
test_stds = prob_model.predictive.std(
    inputs_loader=test_inputs_loader, variances=test_variances
)
test_cred_intervals = prob_model.predictive.credible_interval(
    inputs_loader=test_inputs_loader
)

from fortuna.data import InputsLoader

mesh = np.linspace(-4, 4.2)
mesh_loader = InputsLoader.from_array_inputs(mesh)
mesh_mean = prob_model.predictive.mean(mesh_loader)
mesh_std = prob_model.predictive.std(mesh_loader)
plt.figure(figsize=(6, 3))
plt.scatter(*test_data, s=1, color="C0", label="test data")
plt.plot(mesh, mesh_mean, color="C1", label="mean")
plt.fill_between(
    mesh,
    (mesh_mean - 2 * mesh_std).squeeze(1),
    (mesh_mean + 2 * mesh_std).squeeze(1),
    alpha=0.3,
    color="C0",
    label=f"+/- {2} std",
)
plt.legend(loc="lower right")

<matplotlib.legend.Legend at 0x7f62bbd19490>

计算指标

根据预测统计数据，我们计算指标以评估概率模型对数据的拟合程度，以及不确定性估计的校准情况。

from fortuna.metric.regression import (
    root_mean_squared_error,
    prediction_interval_coverage_probability,
)

test_targets = test_data_loader.to_array_targets()
rmse = root_mean_squared_error(preds=test_means, targets=test_targets)
cred_picp = prediction_interval_coverage_probability(
    lower_bounds=test_cred_intervals[:, 0],
    upper_bounds=test_cred_intervals[:, 1],
    targets=test_targets,
)
print(f"Test RMSE: {rmse}")
print(f"PICP for 95% credible intervals of test inputs: {cred_picp}")

Test RMSE: 0.5550875067710876
PICP for 95% credible intervals of test inputs: 0.9300000667572021

保形区间

PICP指标显示，上述95%可信区间并未完全校准。保形预测方法提供了一种校正它们并改善其校准的方法。

from fortuna.conformal import QuantileConformalRegressor

val_inputs_loader = val_data_loader.to_inputs_loader()
val_cred_intervals = prob_model.predictive.credible_interval(
    inputs_loader=val_inputs_loader
)
test_conformal_intervals = QuantileConformalRegressor().conformal_interval(
    val_lower_bounds=val_cred_intervals[:, 0],
    val_upper_bounds=val_cred_intervals[:, 1],
    test_lower_bounds=test_cred_intervals[:, 0],
    test_upper_bounds=test_cred_intervals[:, 1],
    val_targets=val_data_loader.to_array_targets(),
    error=0.05,
)
conformal_picp = prediction_interval_coverage_probability(
    lower_bounds=test_conformal_intervals[:, 0],
    upper_bounds=test_conformal_intervals[:, 1],
    targets=test_targets,
)
print(f"PICP for 95% conformal intervals of test inputs: {conformal_picp}")

PICP for 95% conformal intervals of test inputs: 0.9340000152587891

另一种可能性是从一维不确定性统计量开始获得保形区间，例如标准差。

from fortuna.conformal import OneDimensionalUncertaintyConformalRegressor

val_means = prob_model.predictive.mean(inputs_loader=val_inputs_loader)
val_stds = prob_model.predictive.std(inputs_loader=val_inputs_loader)

test_conformal_intervals2 = (
    OneDimensionalUncertaintyConformalRegressor().conformal_interval(
        val_preds=val_means,
        val_uncertainties=val_stds,
        test_preds=test_means,
        test_uncertainties=test_stds,
        val_targets=val_data_loader.to_array_targets(),
        error=0.05,
    )
)
conformal_picp2 = prediction_interval_coverage_probability(
    lower_bounds=test_conformal_intervals2[:, 0],
    upper_bounds=test_conformal_intervals2[:, 1],
    targets=test_targets,
)
print(f"PICP for 95% conformal intervals of test inputs: {conformal_picp2}")

PICP for 95% conformal intervals of test inputs: 0.9260000586509705

如果我们有模型输出作为起点怎么办？

如果您已经训练了一个模型并获得了模型输出，您仍然可以使用Fortuna来校准它们，并估计不确定性。仅出于教育目的，让我们将预测均值和对数方差的连接作为模型输出，并假装这些是由其他框架生成的。此外，我们存储了验证和测试目标变量的数组，并假设这些也是给定的。

import numpy as np

calib_outputs = np.concatenate(
    (val_means, np.log(prob_model.predictive.variance(val_inputs_loader))), axis=-1
)
test_outputs = np.concatenate((test_means, np.log(test_variances)), axis=-1)

calib_targets = val_data_loader.to_array_targets()
test_targets = test_data_loader.to_array_targets()

我们现在调用一个校准分类器，使用默认的温度缩放输出校准器，并校准模型输出。

from fortuna.output_calib_model import OutputCalibRegressor, Config, Monitor

calib_model = OutputCalibRegressor()
calib_status = calib_model.calibrate(
    calib_outputs=calib_outputs,
    calib_targets=calib_targets,
    config=Config(monitor=Monitor(early_stopping_patience=2)),
)

Epoch: 100 | loss: -2.3487: 100%|██████████| 100/100 [00:00<00:00, 250.92it/s]

与上述类似，我们现在可以计算预测统计量。

test_log_probs = calib_model.predictive.log_prob(
    outputs=test_outputs, targets=test_targets
)
test_means = calib_model.predictive.mean(outputs=test_outputs)
test_variances = calib_model.predictive.variance(outputs=test_outputs)
test_stds = calib_model.predictive.std(outputs=test_outputs, variances=test_variances)
test_cred_intervals = calib_model.predictive.credible_interval(outputs=test_outputs)

然后可以计算指标和保形区间，完全如上所述。