Investigate alternative approaches for pre- and post-fit model uncertainty calculation (e.g. for data/model plots and tables). Approximate methods could be particularly useful for very complex models where computational cost is a concern.

iminuit 2.5.0 introduced a utility for error propagation via a numerically computed Jacobi matrix. A notebook shows its use. Applied to pyhf, it looks like this:

result, result_obj = pyhf.infer.mle.fit(data, model, return_result_obj=True)

from iminuit.util import propagate
y, ycov = propagate(
    lambda p: model.expected_data(p, include_auxdata=False),
    result_obj.minuit.values,
    result_obj.minuit.covariance,
)

where y are the post-fit model predictions per bin and np.diag(ycov)** 0.5 the uncertainty per bin.

The bootstrap version looks similar:

result, result_obj = pyhf.infer.mle.fit(data, model, return_result_obj=True)

import numpy as np
rng = np.random.default_rng(1)
par_b = rng.multivariate_normal(
    result_obj.minuit.values, result_obj.minuit.covariance, size=50000
)
y_b = [model.expected_data(p, include_auxdata=False) for p in par_b]
yerr_boot = np.std(y_b, axis=0)

where yerr_boot are the uncertainties per bin. This should probably make use of batching in model.expected_data.

The results of the two methods are similar to the model_utils.calculate_stdev implementation, which agrees well with TRExFitter, example below:

• iminuit.util.propagate (0.02 sec):

[[1.58724387, 5.67483153, 4.58648218, 2.45736349, 2.01580335, 1.08836720],
 [1.23555849, 2.11819107, 0.84599747]]

• bootstrap (12.7 sec for 50000 samples):

[[1.52698188 6.18734267 4.6615714  2.43653558 2.02972604 1.19337253], 
 [1.18446726 2.30379686 0.92395196]]

• model_utils.calculate_stdev (1 sec 0.03 sec after perf: vectorize yield uncertainty calculation #316, calculates additional things):

[[1.51192818, 5.84456980, 4.44760681, 2.37522912, 2.01413137, 1.15397110],
 [1.13756210, 2.11783607, 0.78566074]]

• TRExFitter reference (completely independent, including fit):

[[1.50978849, 5.85530619, 4.46335616, 2.37452751, 2.01563069, 1.14129006],
 [1.13406873, 2.11857512, 0.78459717]]