from autoemulate.emulators import GaussianProcess
from autoemulate.core.compare import AutoEmulate
from autoemulate.simulations.epidemic import Epidemic
from autoemulate.calibration.history_matching import HistoryMatchingWorkflow
from autoemulate.calibration.history_matching_dashboard import HistoryMatchingDashboard

random_seed = 42

History Matching workflow

The HistoryMatching calibration method can be a useful way to iteratively decide which simulations to run to generate data refine the emulator on. The HistoryMatchingWorkflow implements this iterative sample-predict-refit workflow. Each time it is run:

• parameters are sampled from the not ruled out yet (NROY) space

• an emulator is used in combination with HistoryMatching to score the implausability of the samples

• simulations are run for a subset of the NROY samples

• the emulator is refit given the newly simulated data

In this tutorial, we demonstrate how to implement this simulator in the loop workflow.

1. Simulate data and train an emulator

In this example, we’ll use the Epidemic simulator, which returns the peak infection rate given two input parameters, beta(the transimission rate per day) and gamma (the recovery rate per day).

simulator = Epidemic(log_level="error")
x = simulator.sample_inputs(500, random_seed=random_seed)
y = simulator.forward_batch(x)

For the purposes of this tutorial, we will restrict the model choice to GaussianProcess.

ae = AutoEmulate(
    x,
    y,
    models=[GaussianProcess],
    model_tuning=False,
    log_level="error",
    random_seed=random_seed
)

We can verify that the fitted emulator performs well on both the train and test data.

ae.summarise()

	model_name	x_transforms	y_transforms	params	rmse_test	r2_test	r2_test_std	r2_train	r2_train_std
0	GaussianProcess	[StandardizeTransform()]	[StandardizeTransform()]	{}	0.001501	0.999955	0.000012	0.999964	0.000004

model = ae.best_result().model

2. Calibrate

To instantiate the HistoryMatchingWorkflow object, we need an observed mean and, optionally, variance for each simulator output.

observations = {"infection_rate": (0.3, 0.05)}

We also pass the fitted emulator and the simulator to the HistoryMatchingWorkflow.

hmw = HistoryMatchingWorkflow(
    simulator=simulator,
    emulator=model,
    observations=observations,
    threshold=3.0,
    train_x=x,
    train_y=y, 
    random_seed=random_seed
)

The run method implements the iterative sample-predict-refit workflow:

• sample n_test_samples to test from the not ruled out yet (NROY) space

• use emulator to filter out implausible samples and update the NROY space

• run n_simulations predictions for the sampled parameters using the simulator

• refit the emulator using the simulated data

The HistoryMatchingWorkflow object maintains and updates the internal state each time run() is called so the full workflow can be run over a number of iterations.

test_parameters, impl_scores = hmw.run(n_simulations=20, n_test_samples=100)

3. Visualise results

For visualising results, you can use the HistoryMatchingDashboard which implements a number of interactive plots. To initialize it just pass in outputs of the run() method and information about the simulation.

dashboard = HistoryMatchingDashboard(
    samples=test_parameters,
    impl_scores=impl_scores,
    param_names=simulator.param_names,  
    output_names=simulator.output_names, 
)

The plots can be viewed using the display() method. Below you can see an example image of what the dashboard looks like.

# dashboard.display()