Approximate Control Variate Monte Carlo

This tutorial builds upon Two model Approximate Control Variate Monte Carlo, Multi-level Monte Carlo, and Multi-fidelity Monte Carlo. MLMC and MFMC are two approaches which can utilize an ensemble of models of vary cost and accuracy to efficiently estimate the expectation of the highest fidelity model. In this tutorial we introduce a general framework for ACVMC when using two or more models to estimate vector-valued statistics.

The approximate control variate estimator of a vector-valued statistic Q0RS that uses M lower fidelity models is

QACV(Z0,Z1,Z1,,ZM,ZM)=Q0(Z0)+[η1,1η1,SMη2,1η2,SMηS,1ηS,SM][Q1(Z1)Q1(Z1)Q2(Z2)Q2(Z2)QM(ZM)QM(ZM)]RS

or in more compact notation

QACV(ZACV)=Q0(Z0)+ηΔ(ZΔ),Δ(ZΔ)=[Δ1(Z1,Z1)ΔM(ZM,ZM)]RSM,ηRS×SM,

where the entries of η are called control variate weights, ZΔ={Z1,Z1,,ZM,ZM}, and ZACV={Z0,ZΔ}.

Here η=[η1,,ηM]T, Δ=[Δ1,,ΔM]T, and Zα, Zα are sample sets that may or may not be disjoint. Specifying the exact nature of these sets, including their cardinality, can be used to design different ACV estimators which will discuss later.

This estimator is constructed by evaluating each model at two sets of samples Zα={z(n)}n=1Nα and Zα={z(n)}n=1Nα where some samples may be shared between sets such that in some cases ZαZβ.

For any η(ZACV), the covariance of the ACV estimator is

V[QACV]=V[Q0]+ηCov[Δ,Δ]η+ηCov[Δ,Q0]+Cov[Δ,Q0]η,RS×S

The control variate weights that minimize the determinant of the ACV estimator covariance are

η=Cov[Δ,Δ]1Cov[Δ,Q0]

The ACV estimator covariance using the optimal control variate weights is

Cov[QACV,QACV](ZACV)=Cov[Q0,Q0]Cov[Q0,Δ]Cov[Δ,Δ]1Cov[Q0,Δ]

Computing the ACV estimator covariance requires computing Cov[Q0,Δ] andCov[Δ,Δ].

First

Cov[Δ,Δ]=[Cov[Δ1,Δ1]Cov[Δ1,Δ2]Cov[Δ1,ΔM]Cov[Δ2,Δ1]Cov[Δ1,Δ2]Cov[ΔM,Δ1]Cov[ΔM,ΔM]]

and second

Cov[Q0,Δ]=[Cov[Q0,Δ1],,Cov[Q0,ΔM]]

where

Cov[Q0,Δα]=Cov[Q0(Z0),Δα(Zα)]Cov[Q0(Z0),Δα(Zα)]]

When estimating the statistics in Monte Carlo Quadrature: Beyond Mean Estimation the above coavariances involving Δ can be computed using the formula in Delta-Based Covariance Formulas For Approximate Control Variates.

The variance of an ACV estimator is dependent on how samples are allocated to the sets Zα,Zα, which we call the sample allocation A. Specifically, A specifies: the number of samples in the sets Zα,Zα,α, denoted by Nα and Nα, respectively; the number of samples in the intersections of pairs of sets, that is Nαβ=|ZαZβ|, Nαβ=|ZαZβ|, Nαβ=|ZαZβ|; and the number of samples in the union of pairs of sets Nαβ=|ZαZβ| and similarly Nαβ, Nαβ. Thus, finding the best ACV estimator can be theoretically be found by solving the constrained non-linear optimization problem

minAADet[Cov[QACV,QACV]](A)s.t.C(c,A~)Cmax,

Here, A~ is the set of all possible sample allocations, and given the computational costs of evaluating each model once c=[c0,c1,,cM], the constraint ensures that the computational cost of computing the ACV estimator

C(c,A)=α=0MNααcα

is smaller than a user-specified computational budget Cmax.

Unfortunately, to date, no method has been devised to solve the above optimization problem for all possible allocations A~. Consequently, all existing ACV methods restrict the optimization space to AA~. These restricted search spaces are formulated in terms of a set of M+1 independent sample partitions Pm that contain pm samples drawn indepedently from the PDF of the random variable z. Each subset Zα is then assumed to consist of a combinations of these indepedent partitions. ACV methods encode the relationship between the sets Zα and Pm via allocation matrices A. For example, the allocation matrix used by MFMC, which is an ACV method, when applied to three models is given by

A=[01111111000111110000011100000001]

An entry of one indicates in the ith row of the jth column indicates that the ith independent partition Pi is used in the corresponding set Zj if j is odd or Zj if j is even. The first column will always only contain zeros because the set Z0 is never used by ACV estimators.

Once an allocation matrix is specified, the optimal sample allocation can be obtained by optimizing the number of samples in each partition p=[p0,,pM], that is

minpDet[Cov[QACV,QACV]](p;A)s.t.C(c,p;A)Cmax,

The following tutorials introduce different ACV methods and their allocation matrices that have been introduced in the literature.

The following compares the estimator variances of three different ACV estimators. No one estimator type is best for all problems. Consequently for any given problem all possible estimators should be constructed. This requires estimates of the model covariance using a pilot study see Pilot Studies

import numpy as np
import matplotlib.pyplot as plt

from pyapprox.util.visualization import mathrm_labels
from pyapprox.benchmarks import setup_benchmark
from pyapprox.multifidelity.factory import (
    get_estimator, compare_estimator_variances, multioutput_stats)
from pyapprox.multifidelity.visualize import (
    plot_estimator_variance_reductions)

np.random.seed(1)
benchmark = setup_benchmark("polynomial_ensemble")
model = benchmark.fun
cov = model.get_covariance_matrix()
nmodels = cov.shape[0]
target_costs = np.array([1e2], dtype=int)
costs = np.asarray([10**-ii for ii in range(nmodels)])
model_labels = [r'$f_0$', r'$f_1$', r'$f_2$', r'$f_3$', r'$f_4$']

stat = multioutput_stats["mean"](benchmark.nqoi)
stat.set_pilot_quantities(cov)
estimators = [
    get_estimator("mlmc", stat, costs),
    get_estimator("mfmc", stat, costs),
    get_estimator("gmf", stat, costs,
                  recursion_index=np.zeros(nmodels-1, dtype=int))]
est_labels = mathrm_labels(["MLMC", "MFMC", "ACVMF"])
optimized_estimators = compare_estimator_variances(
    target_costs, estimators)

axs = [plt.subplots(1, 1, figsize=(8, 6))[1]]

# get estimators for target cost = 100
ests_100 = [ests[0] for ests in optimized_estimators]
_ = plot_estimator_variance_reductions(
    ests_100, est_labels, axs[0])
plot many model acv

Video

Click on the image below to view a video tutorial on approximate control variate Monte Carlo quadrature

../../_images/acv-thumbnail.png

Total running time of the script: ( 0 minutes 0.345 seconds)

Gallery generated by Sphinx-Gallery