Asymptotic Theory of Simulated Moments

Why Simulated Moments Have Two Sources of Uncertainty

A simulated method of moments estimator matches moments calculated from observed data to moments generated by a model. Its uncertainty therefore has two possible sources: the observed sample is random, and the model moments are approximated using random simulation draws. The asymptotic distribution must account for both. Treating simulated moments as exact can understate uncertainty, while treating every simulation design as equivalent can give the wrong covariance formula.

This tutorial studies a standard setup in which the estimator minimizes a weighted quadratic discrepancy. The key results are consistency, which says the estimator converges to the true parameter, and asymptotic normality, which describes its error after scaling by the square root of the sample size. The covariance formula depends on how simulation draws are organized, so that design is stated explicitly before the result is used.

Estimator and Assumptions

Let \(\theta_0\) be a \(p\)-dimensional true parameter, and let \(m(Y_i)\) be a \(q\)-dimensional vector of moments calculated from observed unit \(i\), where \(q\ge p\). Assume the observed units \(Y_1,\ldots,Y_n\) are independent and identically distributed, with finite second moments. Write the population observed moment as \(\mu_0=E[m(Y_i)]\), and its sample counterpart as \(\bar m_n=n^{-1}\sum_{i=1}^n m(Y_i)\).

At parameter value \(\theta\), let \(\mu(\theta)\) denote the model-implied moment vector. Suppose it can be approximated by simulating \(R\) independent draws for each observed unit. Denote the simulated vector for unit \(i\), draw \(r\), by \(s_{ir}(\theta)\), and define \(\tilde\mu_{n,R}(\theta)=(nR)^{-1}\sum_{i=1}^n\sum_{r=1}^R s_{ir}(\theta)\). The draws are assumed independent across units and repetitions, independent of the observed data, and to satisfy \(E[s_{ir}(\theta)]=\mu(\theta)\). This setup makes the simulation error average across units as well as across repetitions.

For the results below, assume the parameter space is compact, \(W_n\) converges in probability to a positive-definite matrix \(W\), and the moment functions are sufficiently smooth and uniformly well behaved for laws of large numbers and Taylor expansions to apply. Identification requires that the population discrepancy \(h(\theta)=\mu_0-\mu(\theta)\) equal zero only at \(\theta_0\). For asymptotic normality, also assume \(\mu(\theta)\) is differentiable near \(\theta_0\), with Jacobian \(G=\partial\mu(\theta_0)/\partial\theta'\) of full column rank \(p\).

The covariance matrices of the data moments and a single simulated draw at the truth are, respectively, \(\Omega=\operatorname{Var}(m(Y_i))\) and \(\Psi=\operatorname{Var}(s_{ir}(\theta_0))\). These definitions allow dependence among components within each moment vector. The independence assumptions above make the data and simulation contributions uncorrelated.

Consistency

For fixed \(R\), the discrepancy converges to its population counterpart under the assumptions just stated. The observed sample average converges to \(\mu_0\). The simulated average converges to \(\mu(\theta)\), uniformly over the parameter space under the assumed uniform law of large numbers. Therefore \(h_{n,R}(\theta)\) converges uniformly to \(h(\theta)\), and the criterion converges to \(Q(\theta)=h(\theta)'Wh(\theta)\).

Since \(W\) is positive definite, \(Q(\theta)\ge0\), and \(Q(\theta)=0\) exactly when \(h(\theta)=0\). Identification makes \(\theta_0\) the unique minimizer of \(Q\). Uniform convergence of the criterion and uniqueness of that minimizer imply that any minimizer \(\hat\theta\) converges in probability to \(\theta_0\). In this repeated-draws-per-unit design, simulation noise averages over the \(n\) units even if \(R\) is held fixed.

Asymptotic Normality and the Covariance Formula

At the true parameter, the discrepancy is the difference between the observed moment average and the simulated moment average. The multivariate central limit theorem gives \(\sqrt n(\bar m_n-\mu_0)\Rightarrow N(0,\Omega)\). For the simulation average, each unit contributes the average of \(R\) draws, whose covariance is \(\Psi/R\). Consequently, \(\sqrt n(\tilde\mu_{n,R}(\theta_0)-\mu(\theta_0))\Rightarrow N(0,\Psi/R)\). The two terms are independent under the stated assumptions; subtracting one from the other adds their covariance matrices. Since \(\mu_0=\mu(\theta_0)\), \(\sqrt n h_{n,R}(\theta_0)\Rightarrow N(0,V_R)\), where \(V_R=\Omega+\Psi/R\).

To obtain the estimator’s limit, expand the discrepancy around \(\theta_0\). Since its population derivative is \(-G\), \(h_{n,R}(\hat\theta)=h_{n,R}(\theta_0)-G(\hat\theta-\theta_0)+o_p(n^{-1/2})\). For a just-identified problem, \(q=p\) and \(G\) is nonsingular; matching the discrepancy to zero gives the first-order approximation \(\hat\theta-\theta_0=G^{-1}h_{n,R}(\theta_0)+o_p(n^{-1/2})\). For an overidentified problem, minimizing the quadratic criterion projects the discrepancy onto the column space of \(G\), with projection determined by \(W\).

The projection can be derived from the first-order condition. Ignoring terms of smaller order, minimize \((h_{n,R}(\theta_0)-G\delta)'W(h_{n,R}(\theta_0)-G\delta)\) over \(\delta\). Differentiating with respect to \(\delta\) and setting the result to zero gives \(G'WG\delta=G'Wh_{n,R}(\theta_0)\). Full column rank of \(G\) and positive definiteness of \(W\) ensure \(G'WG\) is invertible. Hence \(\delta=(G'WG)^{-1}G'Wh_{n,R}(\theta_0)\), which gives the asymptotic result.

If simulation error is negligible, the covariance reduces to the usual moment-estimator sandwich with \(V_R\) replaced by \(\Omega\). If \(R\) is finite, the term \(\Psi/R\) remains in the root-\(n\) covariance. Increasing \(R\) reduces that contribution, but does not remove the data-sampling contribution \(\Omega\). When \(q>p\), the choice of \(W\) affects how the noisy moment discrepancies are combined. In this setup the asymptotically efficient choice is \(W=V_R^{-1}\), when \(V_R\) is nonsingular; using \(\Omega^{-1}\) alone generally ignores simulation noise.

When Simulation Is Shared Across the Sample

A different design uses one simulated sample of size \(S\) to approximate the model moment vector, rather than drawing separately for each observed unit. Write the simulation average as \(\tilde\mu_S(\theta)\), with per-draw covariance \(\Psi\). Its simulation error has covariance approximately \(\Psi/S\), while the observed average has covariance \(\Omega/n\). Under independent data and simulation, if \(n/S\to\kappa\) for a finite nonnegative constant \(\kappa\), then \(\sqrt n h_{n,S}(\theta_0)\Rightarrow N(0,\Omega+\kappa\Psi)\). The corresponding estimator covariance replaces \(V_R\) in the theorem by \(\Omega+\kappa\Psi\).

Common Pitfalls and Interpretation

• Ignoring simulation error. Using only \(\Omega\) is incorrect when the simulated moments are random at the relevant asymptotic scale. Include \(\Psi/R\) for independent repeated draws per unit, or the appropriate \(n/S\) contribution for a shared simulation sample.
• Mixing simulation designs. A fixed number of draws for each of \(n\) units and a fixed total simulation sample are not interchangeable. In the first design, simulation noise averages across units; in the second, the shared simulation error can persist as \(n\) increases.
• Confusing root-\(n\) covariance with estimator covariance. The theorem gives the covariance of \(\sqrt n(\hat\theta-\theta_0)\). The approximate covariance of \(\hat\theta\) itself is that matrix divided by \(n\).
• Reading a standard error as a causal statement. The asymptotic distribution describes estimation uncertainty under the model and regularity conditions. It does not by itself establish a causal interpretation for a parameter.

In practice, the formula identifies what must be estimated to quantify uncertainty: the data-moment covariance, the simulation-moment covariance, the derivative matrix, and the weighting matrix. Their estimates should correspond to the actual simulation design and to the same moment definitions used to construct the estimator. The asymptotic approximation is most useful when the sample is large enough for the central limit and linearization arguments to be credible.