Multifidelity Monte Carlo Estimator for High-Dimensional Quantities of Interest

1 Introduction

• Prior vector approaches either retain tractable closed-form design while treating components separately or exploit cross-output correlations without a closed-form design.The introduction also notes that general ACV and MLBLUE flexibility requires estimated covariance information, which can be costly when pilot data are scarce.Source: §1 · pp. 2–3
• Full-matrix correction can reduce population error at least as much as scalar or componentwise correction, but estimating its covariance matrices can be unstable when pilot samples are scarce relative to output dimension.The paper therefore frames the practical choice as a balance between population optimality and pilot-estimation stability.Source: §1 · p. 3
• The method replaces scalar control-variate coefficients with matrices and minimizes weighted mean-squared error, with an effective correlation extending scalar MFMC’s closed-form design.The effective correlation replaces squared Pearson correlation in the allocation rule, cost–correlation ordering condition, and model selection.Source: §1 · pp. 3–4
• Block-diagonal, diagonal, and scalar coefficient classes form a hierarchy trading variance reduction for lower pilot-data demand, with exact MSE increases over the full class expressed through effective-correlation gaps.Closed-form coefficients and effective correlations are available for compatible weightings.Source: §1 · pp. 3–4
• For a scalar QoI obtained by linearly projecting the state, correcting the full state before projection is at least as accurate as scalar MFMC on projected outputs with population-optimal coefficients.Source: §1 · p. 4
• Experiments on Gaussian, advection–diffusion ROM, and two-phase gas-injection problems show that the suitable coefficient class depends on pilot size relative to the covariance quantities it must estimate.The pilot-based algorithm also allows pilot reuse and an exact reduction to regression on known low-dimensional ROM coordinates.Source: §1 · p. 4

2.1 Scalar-valued MFMC estimator

• Scalar MFMC estimates the high-fidelity model mean under a fixed computational budget using lower-cost surrogate models correlated with the high-fidelity output.Source: §2.1 · pp. 4–5
• Nested sampling evaluates each model on an initial subset of a common i.i.d. sample pool, so cheaper models receive more evaluations than expensive ones.Source: §2.1 · p. 5
• The estimator is unbiased for any control-variate coefficients and nested sample sizes, and its mean-squared error therefore equals its variance.Source: §2.1 · p. 5
• Each correction’s contribution to variance depends on its coefficient and the sample-size gap between adjacent levels, while the design jointly chooses coefficients and sample sizes under the cost constraint.Source: §2.1 · p. 5
• Under the stated correlation ordering and cost–correlation condition, scalar MFMC has closed-form optimal coefficients and sample sizes determined by model costs and squared correlations.For vector outputs, the effective correlation plays the role of squared Pearson correlation in the inherited design.Source: §2.1 · p. 6

3.1 Problem setup and the estimator

• The target is the mean of a vector-valued high-fidelity output, with lower-fidelity models evaluated at lower cost under a fixed computational budget.Source: §3.1 · p. 6
• The estimator retains scalar MFMC’s nested common-sample structure but replaces each scalar coefficient with a matrix control variate.Source: §3.1 · pp. 6–7
• Matrix coefficients allow each high-fidelity component to use linear combinations of low-fidelity components, exploiting cross-component correlations for variance reduction.Source: §3.1 · p. 7
• A user-chosen symmetric positive-semidefinite weighting matrix encodes the relative importance of output components and defines how estimation errors are measured.Examples include equal component weights, diagonal componentwise weights, and weights for linearly mapped quantities of interest.Source: §3.1 · p. 7

3.2 Weighted mean-squared error of the estimator

• The estimator is unbiased for any coefficient matrices and nested sample sizes, so its weighted mean-squared error equals the trace of its covariance weighted by M.Source: §3.2 · p. 8
• Sample means computed from nested shared samples have cross-covariance equal to the model cross-covariance divided by the larger sample size.The correction increments across levels are mutually uncorrelated, while each can covary with the high-fidelity sample mean to cancel part of its sampling error.Source: §3.2 · pp. 8–9
• The weighted MSE consists of the plain high-fidelity Monte Carlo error plus level-specific products of a sample-size factor and a coefficient-dependent factor.A correction helps when its coefficient-dependent term Φ_i(A_i) is negative.Source: §3.2 · p. 11
• The MSE decomposition decouples the design into minimizing each level’s coefficient-dependent term and then choosing nested sample sizes under the computational budget.Source: §3.2 · p. 11

3.3.1 The full matrix class

• The full-matrix class has a closed-form global minimizer, unique when M ≻ 0.Its coefficient uses a d × d covariance matrix, which can be unstable or infeasible to estimate from limited pilot data in high dimensions.Source: §3.3 · p. 11, §3.3.1 · p. 12

3.3.2 The effective correlation

• Effective correlation measures the fraction of weighted high-fidelity variance a model can cancel, generalizing squared Pearson correlation and increasing with the admissible class.It lies between 0 and 1, and the full class has at least as much effective correlation as any contained structural class.Source: §3.3.2 · pp. 12–13

3.3.3 Closed forms for the reduced classes

• Block-diagonal coefficients exploit correlations within each block but ignore those across blocks; only within-block covariance matrices must be estimated and inverted.The partition can group components from the same physical field, spatial region, or coupled subsystem.Source: §3.3.3 · pp. 13–14
• Diagonal coefficients correct each component using only its matching low-fidelity component, while scalar coefficients share one coefficient across all components.The diagonal class ignores cross-component covariance and needs no matrix inversion; the scalar class is the cheapest structure to estimate, using one coefficient per fidelity level.Source: §3.3.3 · pp. 14–15
• The full class requires npilot > d, the block class requires npilot > maxb db, and the scalar and diagonal classes require no inversion.Although a more expressive class is never worse at the population level, its additional covariance estimates increase pilot-data demand; closed forms for reduced classes depend on compatible weightings.Source: §3.3.3 · pp. 13, 15, Table 1 · p. 16

3.4 Optimal sample allocation

• Each ξ_i is the weighted high-fidelity variance share explained by model i but not model i + 1, and this share determines the sampling cost at level i.Source: §3.4 · p. 16

Then the constrained allocation problem

• Under one of the four structural classes and cost–correlation condition (3.10), Theorem 3.18 minimizes weighted MSE over A₂,…,Aₖ ∈ S and sample sizes 0 < m₁ ≤ ··· ≤ mₖ with budget Σwᵢmᵢ = p, using optimal coefficients.The resulting optimal MSE is given by (3.11), with ξᵢ = γ₁,ᵢ − γ₁,ᵢ₊₁.Source: Then the constrained allocation problem · p. 18
• Given 1 = γ1,1 > γ1,2 > ··· > γ1,k > γ1,k+1 = 0 and cost–correlation condition (3.10), weighted MSE E_M(μ̂_MF) is minimized over A2,…,Ak ∈ S and sample sizes 0 < m1 ≤ ··· ≤ mk with budget Σ w_i m_i = p.Source: Then the constrained allocation problem · p. 18
• Replacing squared Pearson correlations with effective correlations preserves scalar MFMC’s allocation ratios, cost–correlation ordering, and model-selection logic.For d = 1, the vector-valued and scalar results coincide.Source: Then the constrained allocation problem · p. 18

• At a fixed nested allocation, the MSE increase from a reduced class is quantified by the gaps between its effective correlations and those of the full class.The loss at level i is proportional to the effective-correlation gap and the allocation-dependent δ_i.Source: §3.5 · p. 19
• The class comparison follows from nested feasible sets for optimized allocations and from subtracting the class-specific MSE expressions at fixed allocation.Source: §3.5 · p. 19
• With each class using its own optimal allocation, the full coefficient class is never worse than a reduced class.Source: §3.5 · p. 19

• Correcting the full state and applying the linear map afterward is equivalent to the weighted vector problem with M = B⊤B.The mapped-output coefficient optimization is equivalent to minimizing the weighted MSE over the state-level coefficients and allocation.Source: §3.6 · p. 20
• For a scalar linear projection, full-state correction is at least as accurate as scalar MFMC on projected outputs, whether each estimator uses its own optimal allocation or a fixed nested allocation.The comparison uses effective correlations evaluated with M = cc⊤.Source: §3.6 · pp. 20–21
• The full-state control variate can use correlation information discarded by projection: γ_scalar_1,i measures correlation between projected outputs, while γ_full_1,i measures how well the full low-fidelity state predicts the scalar high-fidelity quantity.Source: §3.6 · p. 21
• The full-state method’s practical price is that its optimal coefficients involve covariance matrices that must be estimated from pilot data.Source: §3.6 · p. 21

3.7 Practical implementation

• The pilot evaluates all models on common samples to estimate the covariances, coefficients, and effective correlations that determine production allocation.The implementation selects a model subset satisfying the effective-correlation ordering and cost–correlation condition, including the HF-only option, by predicted MSE.Source: §3.7 · p. 21
• Production evaluations use new samples independent of the pilot, with each model evaluated on the allocated prefix of the production sample set.Source: §3.7 · p. 21

Pilot reuse.

• Pilot samples can be reused in the final estimator at no additional cost, but reuse generally makes the estimator biased.The reused pilot means contribute to correction terms while also being used to estimate the coefficients and allocation, producing a generally nonzero covariance term.Source: Pilot reuse · pp. 21–22

4.1 Experiment 1: Gaussian test problem

• With 78 pilots for d=10, all four estimators substantially reduce MSE, and block and full cross-component classes outperform scalar and diagonal classes.The full covariance is estimated accurately in this setting, and allocations closely track the optimal design.Source: Case 1: baseline parameter setting (d… · p. 23, §4.1 · p. 23
• For d=100 with 78 pilots, the full-matrix estimator's MSE is 44 times HF-only MC; five 20-component blocks avoid the instability and achieve the lowest MSE.The pilot-estimated full covariance has rank at most 77, while each block covariance is estimated from 78 samples.Source: Case 2: increasing the dimension to d… · pp. 23, 25, Figure 1 · p. 25
• Adding a cheaper third model lowers every class's predicted ratio; diagonal and block empirical MSEs fall, while scalar and full MSEs are nearly unchanged.The block class slightly beats the full class even though the full class has larger population correlation at both levels.Source: Case 3: a third fidelity level (k = 3,… · pp. 25–26

Case 4: scalar QoI, full-state versus scalar control variates

• For q=c⊤µ, projected positive and negative correlations cancel (γscalar=0.715); correcting components before averaging reaches γfull=0.962 and reduces MSE by a factor of 2.9.The full control vector reweights each component and corrects its sign before averaging.Source: Case 4: scalar QoI, full-state versus… · p. 26

4.2.1 Model problem and the two fidelities

• The advection–diffusion experiment pairs a 75×75-grid HF state with a rank-10 POD–Galerkin ROM whose mean relative state error is about 0.6%.The total sampling budget is 200 HF-solve equivalents, with pilot sizes of 49 and 99.Source: §4.2.1 · pp. 26–27, §4.2 · p. 26

4.2.2 Exact reduction to the ROM coordinates

• The optimal correction regresses the HF output on reduced coordinates q; only G ∈ Rd×r and Cq ∈ Rr×r are estimated, and an r×r system is solved (10×10 rather than 5329×5329).Thus, the full class needs npilot > r rather than npilot > d.Source: §4.2.2 · p. 28

4.2.3 Case 1: full-field quantity of interest

• Every estimator reduces MSE by a factor of 16–27; the full class trails componentwise classes at 49 pilots but becomes best at 99, despite one-third less production budget.The scalar and diagonal classes deteriorate with the larger pilot because their few parameters were already well estimated and the additional pilot samples consume production budget.Source: Table 5 · p. 28, §4.2.3 · p. 28

4.2.4 Case 2: a mapped scalar quantity of interest

• For the regional mean, full-state correction beats scalar MFMC on projected outputs at both pilot sizes, by factors of 2.2 at 49 and 5.1 at 99.The projected-output correlation is 0.928, compared with full-state effective correlation 0.997.Source: §4.2.4 · p. 29, Figure 2 · p. 29

4.3.2 Estimator implementation

• The coarse model approximately resolves the sample-dependent plume front, while many saturation components are usually zero, making pilot-estimated low-fidelity covariances ill-conditioned.The block, full, and full-state estimators use a truncated pseudo-inverse that discards eigenvalues below 10^-2 times the largest.Source: §4.3 · p. 30, §4.3.2 · pp. 30–31

4.3.1 Model problem and the two fidelities

• The two-phase experiment compares a 30 × 30 high-fidelity saturation model with a 15 × 15 coarse-grid model costing 0.185 HF-solve equivalents.The random input is a conditioned horizontal log-permeability field represented by 437 leading Karhunen–Loève modes.Source: §4.3 · p. 29, §4.3.1 · p. 30

4.3.3 Case 1: full saturation field

• At 25 pilots, the scalar and diagonal classes reduce full-field MSE by 16% and 26%, while the full class is worse than HF-only MC; at 50 pilots, the block class has the lowest MSE.With a truncated pseudo-inverse, the full class estimates about 156 active coarse cells from only 25–50 pilot samples; ridge regularization instead makes both block and full classes fail.Source: §4.3.3 · p. 32, Table 8 · p. 32

4.3.4 Case 2: mean gas saturation in a band across the plume

• At both pilot sizes, the full-state control variate lowers band-mean MSE more than projected scalar MFMC, by factors of 1.22 and 1.34 against a predicted 1.35.The projected-output squared correlation is 0.938, compared with a full-state effective correlation of 0.978.Source: §4.3.4 · p. 33, Figure 4 · p. 33, Table 9 · p. 35

4.4 Pilot reuse

• Pilot reuse lowers empirical MSE in 34 of 36 comparisons; the two exceptions rise by less than 1%, within 0.05 standard errors.The reduction grows with the pilot budget share: 38–52% in Experiment 1, 22–51% in Experiment 2, and 6–21% for the Experiment 3 field.Source: §4.4 · p. 33, Figure 5 · p. 34

5 Conclusion

• The full-matrix coefficient and effective correlation preserve the closed-form MFMC allocation, ordering condition, and model selection, while block-diagonal, diagonal, and scalar classes form a nested hierarchy.The hierarchy's population variance reduction does not decrease as the coefficient class becomes more expressive.Source: §5 · p. 35
• The best coefficient class depends on pilot size relative to the covariance quantities it estimates: full classes can lose their advantage or become unstable, while reduced classes remain robust.For a reduced-order model, the full class reduces to regression on r reduced coordinates, so pilot demand scales with r rather than output dimension d.Source: §5 · p. 35
• For a linearly projected scalar QoI, correcting the full state before projection is at least as accurate with accurate covariance information, with larger gains when projection discards correlation.The full-state regression requires the same regularization of the low-fidelity covariance as the full class.Source: §5 · pp. 35–36