Silviase/jgov_v1.4_full_ann huggingface dataset:
The optimal wage $ W_{t}^{\star} $ is the nominal wage that insures that the sum of current and expected discounted marginal profits are zero when the union can only revise nominal wages by using the ad hoc rule (1.41).

It is possible to obtain a recursive expression for multiplier $ \Upsilon_{t} $ that appears in the expression for a union profit. Equation (1.38) can be written equivalently in the form

$$ \Upsilon_{t}^{1-\theta_{s}(1+\psi_{s})}=\int_{0}^{1}\left(\frac{W_{t}(\varsigma)}{W_{t}}\right)^{1-\theta_{s}(1+\psi_{s})}\mathrm{d}\varsigma $$

The wage offered by the union at date t appears under the integral sign. This price has been determined optimally j periods before with probability $ (1 - \xi_{w})\xi_{w}^{j} $ . We can then rewrite the integral as

$$ \Upsilon_{t}^{1-\theta_{s}(1+\psi_{s})}=(1-\xi_{w})\sum_{j=0}^{\infty}\xi_{w}^{j}\left(\frac{\Omega_{t,t-j}W_{t-j}^{\star}}{W_{t}}\right)^{1-\theta_{s}(1+\psi_{s})} $$

where $ W_{t-j}^{\star} $ is the optimal wage at date t-j. Finally, one can interpret the infinite sum as the solution of the following recursive equation:

$$ \begin{aligned}\Upsilon_{t}^{1-\theta_{s}(1+\psi_{s})}=&\big(1-\xi_{w}\big)\left(\frac{W_{t}^{\star}}{W_{t}}\right)^{1-\theta_{s}(1+\psi_{s})}\\&+\xi_{w}\left(\frac{\Omega_{t,t-1}}{W_{t}/W_{t-1}}\right)^{1-\theta_{s}(1+\psi_{s})}\Upsilon_{t-1}^{1-\theta_{s}(1+\psi_{s})}\end{aligned} $$

One can get a more explicit expression for the wage that satisfies equation (1.46). Substituting in this equation the expression for marginal profit (1.45) and dividing by $ W_{t}^{\star} - (1 + \psi_{s}) \theta_{s} $ , one gets

$$ \frac{w_{t}^{\star}}{w_{t}}=\frac{\theta_{s}(1+\psi_{s})}{\theta_{s}(1+\psi_{s})-1}\frac{\mathcal{H}_{1,t}}{\mathcal{H}_{2,t}}+\frac{\psi_{s}}{\theta_{s}(1+\psi_{s})-1}\left(\frac{w_{t}^{\star}}{w_{t}}\right)^{1+(1+\psi_{s})\theta_{s}}\frac{\mathcal{H}_{3,t}}{\mathcal{H}_{2,t}} $$

where $ w_{t}^{\star} $ is the real wage obtained by the union at date t when it can adjust the nominal wage in an optimal manner and $ w_{t} $ the real nominal wage in the economy, with

$$ \mathcal{H}_{1,t}=\mathbb{E}_{t}\sum_{j=0}^{\infty}(\beta\xi_{w})^{j}\lambda_{t+j}w_{t+j}^{m}\left(\frac{\Omega_{t+j}}{\frac{w_{t+j}}{w_{t}}\frac{P_{t+j}}{P_{t}}}\right)^{-(1+\psi_{s})\theta_{s}}\Upsilon_{t+j}^{(1+\psi_{s})\theta_{s}}\mathcal{L}_{t+j} $$
$$ \begin{aligned}\mathcal{H}_{2,t}=\mathbb{E}_{t}\sum_{j=0}^{\infty}(\beta\xi_{w})^{j}\lambda_{t+j}\varepsilon_{l,t+j}w_{t+j}\times\\ \left(\frac{\Omega_{t+j}}{\frac{w_{t+j}}{w_{t}}\frac{P_{t+j}}{P_{t}}}\right)^{1-(1+\psi_{s})\theta_{s}}\Upsilon_{t+j}^{(1+\psi_{s})\theta_{s}}\mathcal{L}_{t+j}\\\mathcal{H}_{3,t}=\mathbb{E}_{t}\sum_{j=0}^{\infty}(\beta\xi_{w})^{j}\lambda_{t+j}\varepsilon_{l,t+j}w_{t+j}\frac{\Omega_{t+j}}{\frac{w_{t+j}}{w_{t}}\frac{P_{t+j}}{P_{t}}}\mathcal{L}_{t+j}\end{aligned} $$

Noticing that $ w_{t+j}/w_{t} $ , the growth factor of the real wage between t and $ t+j $ , can be equivalently written as $ \Pi_{i=1}^{j}\varpi_{t+i} $ ( $ \varpi_{t} $ is the growth factor of the real wage between t and t-1) and that we have

$$ \Omega_{t+j,t}=(1+g)^{j}\left(\prod_{h=1}^{j}\mathcal{E}_{t+h}\right)^{\frac{1}{1-\rho_{x}}}\left(\prod_{h=0}^{j-1}\bar{\pi}_{t+h}\right)^{\gamma_{w}}\left(\prod_{h=0}^{j-1}\pi_{t+h}\right)^{1-\gamma_{w}}, $$

we can finally represent variables $\mathcal{H}_{1,t}$, $\mathcal{H}_{2,t}$ and $\mathcal{H}_{3,t}$ in the recursive form

$$ \begin{aligned}\mathcal{H}_{1,t}=&\lambda_{t}w_{t}^{m}\mathcal{L}_{t}\Upsilon_{t}^{(1+\psi_{s})\theta_{s}}\\&+\beta\xi_{w}\mathbb{E}_{t}\left[\left(\frac{\varpi_{t+1}\pi_{t+1}}{(1+g)\mathcal{E}_{t+1}^{\frac{1}{1-\rho_{x}}}\pi_{t-1}^{\gamma_{w}}\bar{\pi}_{t}^{1-\gamma_{w}}}\right)^{(1+\psi_{s})\theta_{s}}\mathcal{H}_{1,t+1}\right]\end{aligned} $$

$$ \begin{aligned}\mathcal{H}_{2,t}=&\lambda_{t}\varepsilon_{l,t}w_{t}\mathcal{L}_{t}\Upsilon_{t}^{(1+\psi_{s})\theta_{s}}\\&+\beta\xi_{w}\mathbb{E}_{t}\left[\left(\frac{\varpi_{t+1}\pi_{t+1}}{(1+g)\mathcal{E}_{t+1}^{\frac{1}{1-\rho_{x}}}\pi_{t-1}^{\gamma_{w}}\bar{\pi}_{t}^{1-\gamma_{w}}}\right)^{(1+\psi_{s})\theta_{s}-1}\mathcal{H}_{2,t+1}\right]\end{aligned} $$

$$ \begin{aligned}\mathcal{H}_{3,t}&=\lambda_{t}\varepsilon_{l,t}w_{t}\mathcal{L}_{t}\\&\quad+\beta\xi_{w}\mathbb{E}_{t}\left[\left(\frac{(1+g)\mathcal{E}_{t+1}^{\frac{1}{1-\rho_{x}}}\pi_{t-1}^{\gamma_{w}}\bar{\pi}_{t}^{1-\gamma_{w}}}{\varpi_{t+1}\pi_{t+1}}\right)\mathcal{H}_{3,t+1}\right]\end{aligned} $$

Noticing that $ \vartheta_{s,t} \equiv \int_{0}^{1} \frac{W_{t}(\varsigma)}{W_{t}} d\varsigma $ can be written in the recursive form

$$ \vartheta_{s,t}=(1-\xi_{w})\frac{w_{t}^{\star}}{w_{t}}+\xi_{w}\frac{(1+g)\mathcal{E}_{t}^{\frac{1}{1-\rho_{x}}}\pi_{t-1}^{\gamma_{w}}\bar{\pi}_{t}^{1-\gamma_{w}}}{\varpi_{t}\pi_{t}}\vartheta_{s,t-1} $$
we can rewrite equation $ (1.39) $ as

$$ \frac{\psi_{s}\vartheta_{s,t}}{1+\psi_{s}}+\frac{\Upsilon_{t}}{1+\psi_{s}}=1 $$

In the end, wage dynamics are characterized by equations (1.44), (1.51), (1.47), (1.48), (1.50a), (1.50b), (1.50c).

### 1.4 Government and monetary authority

#### 1.4.1 Fiscal policy

We assume that exogenous government expenditures $ G_{t} = g_{t}Y_{t} $ are exactly financed by lump sum taxes:

$$ T_{t}=P_{t}G_{t}. $$

#### 1.4.2 Central Bank

We assume that the behavior of the central bank is adequately described by the following Taylor rule:

$$ R_{t}=\max\left\{1,R_{t-1}^{\rho_{R}}\left[R^{\star}\left(\frac{\pi_{t-1}}{\bar{\pi}_{t}}\right)^{r_{\pi}}\left(\frac{Y_{t}}{\mathcal{Y}_{t}}\right)^{r_{Y}}\right]^{1-\rho_{R}}\varepsilon_{R,t}\right\} $$

where $ \bar{\pi}_{t} $ is the inflation target of the central bank, $ Y_{t} $ is the reference output level that would be attained by an economy without nominal rigidities, $ \log \varepsilon_{R,t} $ is an AR(1) stationary process with zero mean. The max function constrains the nominal interest factor to be greater or equal than one. Equation (1.53) defines two regimes. In the first one, the nominal interest rate is constant (equal to zero), whereas in the second one the nominal interest rate reacts to excess inflation and output gap fluctuations. In the sequel, we consider the case where the steady state of the economy is in the second regime, ie the long run level of the nominal interest rate is assumed to be strictly positive $ ^{2} $ . Nevertheless, during the transitions to this steady state the economy can hit the Zero Lower Bound for the nominal interest rate. Because of this constraint, the model is non differentiable everywhere.
### 1.5 General equilibrium

#### 1.5.1 Price distortion

Prices in the intermediary good sector are heterogeneous. However, it is possible to show that this heterogeneity doesn't hinder aggregation. We know that the firms of this sector choose all the same mix of factor of production in the sense that the ratio of capital demand to labor demand is constant across firms (see equation 1.23). Expressing labor demand of firm $ \iota $ as a function of its demand of physical capital, we can write the production of this firm

$$ y_{t}(\iota)=\left(\frac{A_{t}L_{t}^{d}}{K_{t}^{d}}\right)^{1-\alpha}K_{t}^{d}(\iota) $$

When $ K_{t}^{d} \equiv \int_{0}^{1} K_{t}^{d}(\iota) \mathrm{d}\iota $ is aggregate demand for physical capital and $ y_{t} \equiv \int_{0}^{1} Y_{t}(\iota) \mathrm{d}\iota $ represents the sum of intermediary productions, we can write directly

$$ y_{t}=\left(K_{t}^{d}\right)^{\alpha}\left(A_{t}L_{t}^{d}\right)^{1-\alpha} $$

The sum of intermediary productions is different from $ Y_{t} $ , because aggregation technology isn't linear. Integrating the demand function for good $ \iota $ from the final good producers (1.18) over $ \iota $ , we get

$$ y_{t}=\Delta_{p,t}Y_{t} $$

with

$$ \Delta_{p,t}\equiv\frac{1}{1+\psi_{f}}\int_{0}^{1}\left(\left(\frac{P_{t}(\iota)/P_{t}}{\Theta_{t}}\right)^{-(1+\psi_{f})\theta_{f}}+\psi_{f}\right)\mathrm{d}\iota $$

Price distortion can be written recursively in the following manner:

$$ \Delta_{p,t}=\frac{1}{1+\psi_{f}}\Theta_{t}^{(1+\psi_{f})\theta_{f}}\nabla_{p,t}+\frac{\psi_{f}}{1+\psi_{f}} $$

$$ \nabla_{p,t}=\left(1-\xi_{p}\right)\left(\frac{P_{t}^{\star}}{P_{t}}\right)^{-(1+\psi_{f})\theta_{f}}+\xi_{p}\left(\frac{\bar{\pi}_{t}^{\gamma_{p}}\pi_{t-1}^{1-\gamma_{p}}}{\pi_{t}}\right)^{-(1+\psi_{f})\theta_{f}}\nabla_{p,t-1} $$

where the Lagrange multiplier $\Theta_{t}$ is defined recursively as well. Then, we have

$$ \Delta_{p,t}Y_{t}=\left(K_{t}^{d}\right)^{\alpha}\left(A_{t}L_{t}^{d}\right)^{1-\alpha} $$
#### 1.5.2 Wage distortion

Here, we show how to link aggregate labor supply by the households with aggregated labor supply by the employment agency to the firms of the intermediary good sector. This link is affected by the heterogeneity of wages induced by their nominal rigidity. Integrating labor demand for type $ \varsigma $ by the employment agency (1.37) over $ \varsigma $ , we find directly

$$ L_{t}=\Delta_{w,t}\mathcal{L}_{t} $$

with

$$ \Delta_{w,t}\equiv\frac{1}{1+\psi_{s}}\int_{0}^{1}\left(\left(\frac{W_{t}(\varsigma)/W_{t}}{\Upsilon_{t}}\right)^{-(1+\psi_{s})\theta_{s}}+\psi_{s}\right)\mathrm{d}\varsigma $$

Wage distortion can be written in recursive form:

$$ \Delta_{w,t}=\frac{1}{1+\psi_{s}}\Upsilon_{t}^{(1+\psi_{s})\theta_{s}}\nabla_{w,t}+\frac{\psi_{s}}{1+\psi_{s}} $$

$$ \begin{aligned}\nabla_{w,t}=&\big(1-\xi_{w}\big)\left(\frac{w_{t}^{\star}}{w_{t}}\right)^{-(1+\psi_{s})\theta_{s}}\\&+\xi_{s}\left(\frac{(1+g)\mathcal{E}_{t}^{\frac{1}{1-\rho_{x}}}\bar{\pi}_{t}^{\gamma_{w}}\pi_{t-1}^{1-\gamma_{w}}}{\varpi_{t}\pi_{t}}\right)^{-(1+\psi_{s})\theta_{s}}\nabla_{w,t-1}\end{aligned} $$

where the Lagrange multiplier $ \Upsilon_{t} $ is defined recursively in equation (1.47).

#### 1.5.3 Dividends paid by intermediary good firms

Firms in the intermediary good sector interact in monopolistic competition and make profits that are paid to households in the form of dividends. The sum of nominal profits at date t is

$$ \begin{aligned}\Pi_{t}&=\int_{0}^{1}\Pi_{t}(\iota)\mathrm{d}\iota\\&=P_{t}\left(\varepsilon_{y,t}Y_{t}-r_{t}^{k}K_{t}^{d}-w_{t}L_{t}^{d}\right)\end{aligned} $$

As households own the firms, the profits are paid to them. The repartition of these profits between the households in undetermined in general equilibrium, but we know that

$$ \int_{0}^{1}\mathcal{D}_{1,t}(h)\mathrm{d}h=P_{t}\left(\varepsilon_{y,t}Y_{t}-r_{t}^{k}K_{t}^{d}-w_{t}L_{t}^{d}\right) $$
#### 1.5.4 Dividends paid by the unions

In the same way, we can compute aggregate nominal profit of the unions at date t. These profits as well are paid to the households. We have

$$ \begin{aligned}\mathcal{S}_{t}&=\int_{0}^{1}\mathcal{S}_{t}(\varsigma)\mathrm{d}\varsigma\\&=\varepsilon_{l,t}W_{t}\mathcal{L}_{t}-W_{t}^{m}L_{t}\end{aligned} $$

and, then,

$$ \int_{0}^{1}\mathcal{D}_{2,t}(h)\mathrm{d}h=\varepsilon_{l,t}W_{t}\mathcal{L}_{t}-W_{t}^{m}L_{t} $$

#### 1.5.5 Equilibrium in factor markets and in bond markets

In general equilibrium, labor supply forms the employment agency equals aggregate labor demand by firms of the intermediary good market. In the same way, aggregate supply of physical capital by the households equals aggregate demand from these firms. In formal terms,

$$ \mathcal{L}_{t}\equiv\Delta_{w,t}^{-1}\int_{0}^{1}L_{t}(h)\mathrm{d}h=\int_{0}^{1}L_{t}^{d}(\iota)\mathrm{d}\iota\equiv L_{t}^{d} $$

$$ \widetilde{K}_{t}\equiv\int_{0}^{1}z_{t}(h)K_{t-1}(h)\mathrm{d}h=\int_{0}^{1}K_{t}^{d}(\iota)\mathrm{d}\iota\equiv K_{t}^{d} $$

Finally, aggregate demand for bonds must be zero, as we assume a close economy and no government debt.

$$ \int_{0}^{1}B_{t}(h)\mathrm{d}h=0 $$

#### 1.5.6 Equilibrium on the good market

By summing the budget constraints of the households (1.2) over $ h \in [0,1] $ and by substituting the equilibrium conditions on the bond market, the definition of aggregate dividends and the budget constraint of the government, we get

$$ \begin{aligned}P_{t}G_{t}+P_{t}C_{t}+p_{I,t}P_{t}I_{t}=W_{t}^{m}L_{t}+&P_{t}r_{t}^{K}z_{t}K_{t-1}+P_{t}\left(\varepsilon_{y,t}Y_{t}-r_{t}^{k}K_{t}^{d}-w_{t}L_{t}^{d}\right)\\&+\varepsilon_{lt}W_{t}\mathcal{L}_{t}-W_{t}^{m}L_{t}\end{aligned} $$
After simplification and knowing that the factor markets are in equilibrium, we obtain

$$ G_{t}+C_{t}+p_{I,t}I_{t}=\varepsilon_{y,t}Y_{t}+\left(\varepsilon_{l,t}-1\right)w_{t}\mathcal{L}_{t} $$

or equivalently

$$ G_{t}+C_{t}+p_{I,t}I_{t}=\varepsilon_{y,t}\Delta_{p,t}^{-1}y_{t}+\left(\varepsilon_{l,t}-1\right)\Delta_{w,t}^{-1}w_{t}L_{t} $$

## 2 Extended path

Due to the max function appearing in the Taylor rule, the model is not differentiable everywhere. Consequently, a perturbation approach cannot be considered for solving the model. The reason is that, even if we forget that the Taylor approximation does not apply in this case, by approximating the model around the deterministic steady state the ZLB (or more generally any occasionally binding constraint) would not affect the expectations of the agents.

Our model can generically be represented as

$$ \mathbb{E}_{t}\left[\mathcal{F}\big(y_{t+1},y_{t},y_{t-1},\varepsilon_{t}\big)\right]=0\qquad\forall t $$

where y is a vector of endogenous variables, $ \varepsilon $ is a vector of structural innovations (a Gaussian multivariate white noise), $ E_{t} $ is the conditional expectation operator and F is a non linear function. The standard (numerical or analytical) solution approach is to look for an expression for the decision and transition rules:

$$ y_{t}=\mathcal{G}(y_{t-1},\varepsilon_{t}) $$

where the invariant reduced form G satisfies:

$$ \mathbb{E}_{t}\left[\mathcal{F}\left(\mathcal{G}\left(\mathcal{G}\left(y_{t-1},\varepsilon_{t}\right),\varepsilon_{t+1}\right),\mathcal{G}\left(y_{t-1},\varepsilon_{t}\right),y_{t-1},\varepsilon_{t}\right)\right]=0\qquad\forall t $$

If there is no closed form solution, the function G must be locally (around the deterministic steady state) or globally (over an arbitrary domain for the states variables) approximated. Unfortunately, because of the size of our model a global approximation cannot be considered here.

The extended path approach indirectly characterizes the decision and transition rules G by generating time-series that satisfy the model equations.
Basically, the trick is to move the conditional expectation under function F, that is to replace (2.1) by

$$ \mathcal{F}\left(\mathbb{E}_{t}\left[y_{t+1}\right],y_{t},y_{t-1},\varepsilon_{t}\right)=0\qquad\forall t $$

This is obviously not the same model if the function F is nonlinear, because of the Jensen inequality. The main approximation of the EP method lies here. The magnitude of the errors induced by this approximation will depend on the degree of nonlinearity in the forward terms appearing in the original model. Gagnon (1990) and Love (2009), considering a stochastic growth model, show that the approximation error is reasonable and that the EP approach performs as well (or even better) as a global approximation approach (Galerkin). The advantage of the EP approach over global approximation being that it can solve medium or large scaled models (which is not possible with global approximation methods due to the so-called curse of dimensionality). Technically, the EP method will generate time series for the endogenous variables by calling recursively a perfect foresight model solver (for (2.2)). The algorithm is described in the sequel.

Given initial conditions for the states $ (y_{0}) $ , terminal conditions for the jumping variables $ (y_{N+1} = y^{\star}) $ and a sequence of expected innovations, one can solve a perfect foresight model by solving the following system of non linear equations:

$$ \begin{aligned}\mathcal{F}(y_{2},y_{1},y_{0},\varepsilon_{1})&=0\\\mathcal{F}(y_{3},y_{2},y_{1},\varepsilon_{2})&=0\\\mathcal{F}(y_{4},y_{3},y_{1},\varepsilon_{3})&=0\\&\vdots\\\mathcal{F}(y_{N+1},y_{N},y_{N-1},\varepsilon_{N})&=0\end{aligned} $$

for the equilibrium path $ \{y_{1}, y_{2}, \ldots, y_{N}\} $ . This system of nonlinear equations is solved with a Newton algorithm, taking care of the sparsity of the Jacobian. The sole approximation here concerns the terminal condition: we assume that the endogenous variables reach the steady state at time $ N + 1 $ , while this result is only asymptotic. This size of this approximation can be tailored by evaluating the sensitivity of the solution with respect to the value of the horizon, N. Alternatively, if the steady state of the model is unknown, we could replace the terminal condition by $ y_{N+1} = y_{N} $ . In practice, we may have to set N equal to several hundreds, in order to obtain a reasonable size for the approximation errors. This explains the computational burden.
of this approach (the size of the non linear system of equations to solve grows linearly with N).

We now restrict our attention to the following dynamic problem (at time t):

$$ \mathcal{F}\big(y_{t+j+1},y_{t+j},y_{t+j-1},\varepsilon_{t+j}\big)=0\qquad\forall j=0,1,\ldots,N $$

with no shocks in the future, $ \varepsilon_{j}=0 $ for all j>0. The solution of the corresponding system of non linear equations:

$$ \begin{aligned}\mathcal{F}\big(y_{t+1},y_{t},y_{t-1},\varepsilon_{t}\big)&=0\\\mathcal{F}\big(y_{t+2},y_{t+1},y_{t},0\big)&=0\end{aligned} $$

$$ \begin{aligned}\vdots\\\mathcal{F}(y_{t+N+1},y_{t+N},y_{t+N-1},0)&=0\end{aligned} $$

will be denoted:

$$ \mathcal{Y}_{t,N}\equiv\left\{y_{t},\cdots,y_{t+N}\right\}\equiv\mathcal{Q}_{N,t}(y_{t-1},\varepsilon_{t}) $$

Finally, we further restrict our attention by selecting only the first term of the solution for the equilibrium path :

$$ y_{t}=\mathcal{R}_{N,t}(y_{t-1},\varepsilon_{t}) $$

Note that the mapping between the endogenous variables at time t with the initial condition (endogenous variables at time t - 1) and the time t innovations ( $ \varepsilon_{t} $ ) is not state invariant. If the initial states or the innovations are changed, we have to solve again the perfect foresight model to obtain $ y_{t} $ .

To simulate a time series of T observations we iterate over $ \mathcal{R}_{N,t}(y_{t-1},\varepsilon_{t}) $ . Given an initial conditions for the states and a sequence of unexpected innovations, we have:

$$ \begin{aligned}y_{1}&=\mathcal{H}_{N,1}(y_{0},\varepsilon_{1})\\y_{2}&=\mathcal{H}_{N,2}(y_{1},\varepsilon_{2})\\&\vdots\\y_{T}&=\mathcal{H}_{N,T}(y_{T-1},\varepsilon_{T})\end{aligned} $$

The simulation of a time series with T observations requires the solutions to T perfect foresight models. Obviously, this simulation approach is far more
time consuming than the perturbation (first or second order approximation) approach or the global approach (for small models). In order to reduce the computational burden of this solution method we did not use the traditional terminal condition. Knowing that the approximation error associated to the paths generated by the first order approximation of G are arbitrarily small in the vicinity of the steady state, we replaced the terminal condition $ y_{N+1} = y^{\star} $ by $ y_{N+1} = y_{N+1}^{(1)} $ , where $ y_{N+1}^{(1)} $ is the attained level of the endogenous variables at time $ N+1 $ when we iterate on the first order approximation of G. Using this alternative terminal condition allows us to significantly reduce the value of N. In the same spirit, we also considered the terminal condition $ y_{N+1} = \frac{\partial \mathcal{G}}{\partial y_{N}^{\prime}}(y^{\star})y_{N} $ , but our experiments showed that this strategy was generally slower, due to the propagation of the round off errors.

## 3 Simulated Method of Moments

### 3.1 Intuition and notations

The basic idea is to find the parameters of the model such that the distance between simulated moments (mean, variance, covariance, auto-covariance, skewness, kurtosis, ...) and sample moments is minimized. It can be proved that this estimator provides consistent estimates of the parameters (see chapter 2 in Gourieroux and Monfort (1996)).

Let $ \mathscr{Y}_{T}^{\star}\equiv\{y_{1}^{\star},y_{2}^{\star},\ldots,y_{T}^{\star}\} $ be the sample and $ \mathscr{Y}_{T}^{(s)}(\boldsymbol{\theta})\equiv\{y_{1}^{(s)}(\boldsymbol{\theta}),y_{2}^{(s)}(\boldsymbol{\theta}),\ldots,y_{T}^{(s)}(\boldsymbol{\theta})\} $ be a simulated sample, obtained with the extended path method described in the previous section, for a vector of (estimated) parameters $ \theta $ .

Let $ \mathcal{H}(\mathcal{Y}_{T}^{\star}) $ and $ \mathcal{H}(\mathcal{Y}_{T}^{(s)}(\boldsymbol{\theta})) $ be $ n \times 1 $ vectors of sample and simulated moments, with n > m. Some of the matched moments may require more than one observation to be computed. For instance, if we want to match an order p auto-covariance function, we need to consider p consecutive values of the observed variables. To this end, we define vectors gathering p consecutive sample or simulated variables, $ z_{t}^{\star} \equiv (y_{t}^{\star}, y_{t-1}^{\star}, \ldots, y_{t-p+1}^{\star}) $ and $ z_{t}^{(s)}(\boldsymbol{\theta}) \equiv (y_{t}^{(s)}(\boldsymbol{\theta}), y_{t-1}^{(s)}(\boldsymbol{\theta}), \ldots, y_{t-p+1}^{(s)}(\boldsymbol{\theta})) $ . The sample moments are defined as follows $ ^{3} $

$$ \mathcal{H}(\mathcal{Y}_{T}^{\star})=T^{-1}\sum_{t=1}^{T}h(z_{t}) $$
the same definition applies for the simulated moments. The continuous function $ h : R^{np} \longrightarrow R^m $ defines the set of moments to be matched, $ \mathbb{E}[h(z)] $ . For instance, in the univariate case (m = n = 1), if the we want to match the first order auto-covariance we have $ h(z_t) = (y_t - \bar{y})(y_{t-1} - \bar{y}) $ where $ \bar{y} $ is the arithmetic mean of y.

Finally we define the $ n \times 1 $ moment vector

$$ g_{T}(\boldsymbol{\theta})=\frac{1}{T}\sum_{t=1}^{T}h(z_{t}^{\star})-\frac{1}{S T}\sum_{s=1}^{S}\sum_{t=1}^{T}h(z_{t}^{(s)}(\boldsymbol{\theta})) $$

where S is the number of simulated samples. Note that the seed of the random number routines used to generate simulated samples must be kept constant when the moment vector is evaluated for different values of the vector of parameters $ \theta $ . We will denote $ G(\pmb{\theta}) $ the $ n \times m $ Jacobian matrix associated to the moment vector:

$$ G(\boldsymbol{\theta})=\frac{\partial g}{\partial\theta^{\prime}}(\boldsymbol{\theta}) $$

### 3.2 The simulated moments estimator

The simulated moments estimator of $ \theta $ is defined as:

$$ \hat{\boldsymbol{\theta}}_{T}(W)=\arg\min_{\boldsymbol{\theta}}g_{T}(\boldsymbol{\theta})^{\prime}W g_{T}(\boldsymbol{\theta}) $$

where W is a symmetric positive definite weighting matrix. It can be shown that this estimator is $ \mathcal{O}_{p}(T^{-1/2}) $ and that its asymptotic variance is given by:

$$ \mathbb{V}_{\infty}\left[\sqrt{T}\hat{\boldsymbol{\theta}}_{T}(W)\right]=\kappa(S)\left(G_{0}^{\prime}W G_{0}\right)^{-1}G_{0}^{\prime}W\Omega W G_{0}\left(G_{0}^{\prime}W G_{0}\right)^{-1} $$

where

$$ \kappa(S)=\left(1+\frac{1}{S}\right) $$

is a scale factor monotonically decreasing with the number of simulated samples used to evaluate the moment vector,

$$ G_{0}=G(\pmb{\theta}_{0}) $$

is the Jacobian matrix of the moment vector evaluated at the true value $ \theta_{0} $ , and $ \Omega $ is the long run covariance matrix of the vector of moments:

$$ \Omega=\lim_{T\rightarrow\infty}\mathbb{E}\left[\left(T^{-\frac{1}{2}}\sum_{t=1}^{T}h(z_{t})\right)\left(T^{-\frac{1}{2}}\sum_{t=1}^{T}h(z_{t})\right)^{\prime}\right] $$