History Matching

from autoemulate.emulators import GaussianProcess
from autoemulate.core.compare import AutoEmulate
from autoemulate.simulations.projectile import Projectile
from autoemulate.calibration.history_matching import HistoryMatching

History Matching#

In this tutorial we demonstrate the use of History Matching to determine which points in the input space are plausible given a set of observations.

Performing History Matching requires:

  • a fit emulator that predicts uncertainty (e.g., Gaussian Process) and

  • an observation associated with the simulator output.

The emulator is used to efficiently estimate the simulator output, accounting for all uncertainties. The emulated output is then compared with the observation and parameter inputs that are unlikely to produce the observation can then be “ruled out” as implausible, reducing the input space.

1. Simulate data and train an emulator#

In this example we use the Projectile simulator. It returns distance travelled by a projectile given a drag coefficient c and initial velocity v0.

simulator = Projectile(log_level="error")
x = simulator.sample_inputs(1000)
y = simulator.forward_batch(x)

Since we need an uncertainty aware emulator, we restrict AutoEmulate to only train a Gaussian Process.

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

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()] {'mean_module_fn': <function poly_mean at 0x7f... 26.113054 0.999989 0.000003 0.999991 0.000001 
model = ae.best_result().model

2. Calibrate#

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

#  observed data: (mean, variance)
hm = HistoryMatching(observations={"distance": (2500, 100)})

We also need predictions for a set of query points using the trained emulator. These must have uncertainty estimates.

x_new = simulator.sample_inputs(10)
mean, variance = model.predict_mean_and_variance(x_new)

Primary use of the HistoryMatching class is the calculate_implausibility method, which returns the implausibility metric (a number of standard deviations between the emulator mean and the observation) for the queried points.

implausability = hm.calculate_implausibility(mean, variance)
implausability

tensor([[ 39.9126],
        [ 41.7190],
        [ 12.4470],
        [ 30.7504],
        [ 18.1594],
        [ 36.1814],
        [337.8336],
        [149.1406],
        [ 31.9233],
        [ 42.9756]])

The get_nroy method can be used to filter the queried points given the implausability scores and only retain those that have not been ruled out yet (NROY).

hm.get_nroy(implausability, x_new)

tensor([], size=(0, 2))
Contents