求解资产组合：两资产 HANK 模型
在这个例子中，我们求解了 Auclert, Bardóczy, Rognlie, and Straub (2021, 附录 B.3 和 附录 E.1) 中的两资产 HANK 模型。模型中，每个家庭选择最优的流动性与非流动性的资产组合，并面临着一个凸的资产组合调整成本。此外，新凯恩斯粘性价格企业决定中间产品的价格，面临着 Rotemberg-类型的价格调整成本，并且根据二次型的资本调整成本作出动态的投资决策。工会通过最大化其成员的平均福利来控制工资水平。
这个例子展示了该工具箱解决具有多维内生状态变量、控制变量以及大型均衡方程系统的复杂投资组合问题的能力。
我们遵循原文中的符号、参数化和校准程序。关于均衡方程组的详细推导，建议读者阅读原文。
模型
家庭
家庭具有异质性的马尔可夫劳动效率 (\(e\))、初始流动资产持有量 (\(b\)) 和非流动性资产持有量 (\(a\))，并求解贝尔曼方程：
\[\begin{split}
V_t\left(e_{it}, b_{it-1}, a_{it-1}\right) =\max_{c_{it}, b_{it}, a_{it}}\left\{\frac{c_{it}^{1-\sigma}}{1-\sigma} - \varphi \frac{N_{t}^{1+\nu}}{1+\nu} +\beta \mathbb{E}_t V_{t+1}\left(e_{it+1}, b_{it}, a_{it}\right)\right\}
\\
\quad s.t. \quad c_{it} + a_{it} + b_{it} = \left(1-\tau_{t}\right)w_t e_{it} N_{t} + \left(1+r_{t}^{a}\right) a_{it-1}
\\
+\left(1+r_t^b\right) b_{it-1}-\Phi_{t}\left(a_{it}, a_{it-1}\right)
\\
a_{it} \geq 0, \quad b_{it} \geq \underline{b}
\end{split}\]
其中，\(\tau_t\) 表示统一税率，\(w_t\) 为真实工资，\(N_t\) 为劳动时间，\(r_t^b\) 和 \(r_t^{a}\) 分别为流动性资产和非流动性资产的真实利率，\(\Phi_t\) 为调整成本函数，定义为：
\[
\Phi_{t}\left(a_{it}, a_{it-1}\right)=\frac{\chi_1}{\chi_2}\left|\frac{a_{it}-\left(1+r_t^a\right) a_{it-1}}{\left(1+r_t^a\right) a_{it-1}+\chi_0}\right|^{\chi_2}\left[\left(1+r_t^a\right) a_{it-1}+\chi_0\right]
\]
其中 \(\chi_0 > 0\)，\(\chi_1 > 0\)，并且 \(\chi_2 > 1\)。
家庭选择消费 (\(c_{it}\)) 以及资产组合 \(\left(b_{it}, a_{it}\right)\)，使其满足借贷约束 \(a_{it} \geq 0\) 和 \(b_{it} \geq \underline{b}\)，并提供一致的且由工会决定的劳动 (\(N_{t}\)) 以满足企业的劳动需求。
以 \(\Gamma_t\) 表示家庭状态分布的度量，以 \(g_{c,t}, g_{b,t}, g_{a,t}\) 表示策略函数。
金融中介
家庭将流动性资产和非流动性资产存入代表性金融中介机构。金融中介机构以单位流动性转换成本 \(\omega\) 将流动性的存款 (\(\int{b_{it}di}\)) 投资于非流动性的政府债券之中，从而实现 \(r_t^b = r_t - \omega\)，其中 \(r_t\) 表示政府债券在时期 \(t\) 的实际利率。同时，金融中介机构还将非流动性存款 (\(\int{a_{it}di}\)) 投资于由政府债券和公司股权组成的共同基金之中，这意味着：
\[
1 + r_t^a = \Theta_{pt-1} \frac{D_t + p_t}{p_{t-1}} + (1 - \Theta_{pt-1}) \left(1 + r_t\right)
\]
其中，共同基金股权份额 \(\Theta_{pt-1} = \frac{p_{t-1}}{p_{t-1} + B^{g}_{t-1} - B^{h}_{t-1}}\) 由时期 \(t-1\) 的股权价格 (\(p_{t-1}\))，总体政府债券 \(\left(B^{g}_{t-1}\right)\), 以及总体流动性资产持有量 \(\left(B^{h}_{t-1} = \int{b_{it-1}di}\right)\) 决定。在完美预期均衡中，无套利条件意味着时期 \(t \geq 1\) 未来回报的均等化。
\[
1 + r_t^a = \frac{D_t + p_t}{p_{t-1}} = 1 + r_t, \forall t \geq 1
\]
然而，\(t = 0\) 的资产回报可能会因 \(t = 0\) 发生的意外冲击而发生改变。
企业
中间品厂商使用柯布-道格拉斯生产技术，利用资本和劳动力生产差异化的产品。这些中间产品由竞争性的最终品厂商使用具有替代弹性 \(\mu_{p}/\left(\mu_p - 1\right)\) 的 CES 聚合器进行加总。中间品厂商在中间产品市场上实行垄断竞争，使得动态投资决策服从于一个二次型的调整成本，并面临对数二次型的价格调整成本。
中间品厂商 \(j\) 的需求、技术和调整成本可通过以下公式定义：
产品需求函数：
\[
y_{jt} = \left(\frac{p_{jt}}{P_{t}}\right)^{-\frac{\mu_{p}}{\mu_p - 1}} Y_{t}
\]
生产函数：
\[
y_{jt} = z_t k_{jt-1}^\alpha n_{jt}^{1-\alpha}
\]
资本调整成本函数：
\[
\zeta \left(k_{jt},k_{jt-1}\right) k_{jt-1} = \frac{1}{2 \delta \epsilon_I}\left(\frac{k_{jt}}{k_{jt-1}} - 1\right)^2 k_{jt-1}
\]
价格调整成本函数：
\[
\psi^p_t \left(p_{jt},p_{jt-1}\right)=\frac{\mu_p}{\mu_p-1} \frac{1}{2 \kappa_p}\left[\log \left(\frac{p_{jt}}{p_{jt-1}}\right)\right]^2 Y_t
\]
资本积累方程：
\[
i_{jt} = k_{jt} - \left(1 - \delta\right) k_{jt-1} + \zeta \left(k_{jt},k_{jt-1}\right) k_{jt-1}
\]
股利方程：
\[
d_{jt} = y_{jt} - w_{t} n_{jt} - i_{jt} - \psi^p_t \left(p_{jt},p_{jt-1}\right)
\]
在对称均衡中，企业的最优解可由以下方程表示：
价格菲利普斯曲线：
\[
\log(1+\pi_t) = \kappa_p \left(mc_t - \frac{1}{\mu_p} \right) + \frac{1}{1+r_{t+1}} \frac{Y_{t+1}}{Y_t} \log(1+\pi_{t+1})
\]
投资方程：
\[
Q_t = 1 + \frac{1}{\delta \epsilon_I}\left(\frac{K_t-K_{t-1}}{K_{t-1}}\right)
\]
估值方程：
\[
(1+r_{t+1})Q_{t} = \alpha \frac{Y_{t+1}}{K_t} mc_{t+1} - \left[\frac{K_{t+1}}{K_t} - (1-\delta) + \frac{1}{2\delta \epsilon_I}\left(\frac{K_{t+1} - K_t}{K_t}\right)^2\right] + \frac{K_{t+1}}{K_t}Q_{t+1}
\]
其中 \(mc_t\) 是如下函数的简写：
边际成本函数：
\[
mc_t = \frac{w_t }{(1-\alpha)\frac{Y_t}{N_t}}
\]
工会
竞争性雇佣者加总了由垄断竞争工会连续统提供的差异化劳动力。每个工会 \(k\) 面临一个具有替代弹性 \(\mu_{w}/\left(\mu_{w} - 1\right)\) 的劳动需求函数，并设定名义工资率 (\(W_{kt}\)) 以最大化其成员的平均效用，并服从效用中的二次型调整成本。
\[
\psi_t^w \left(W_{kt},W_{kt-1}\right)=\frac{\mu_w}{\mu_w - 1} \frac{1}{2 \kappa_w}\left[\log \left(\frac{W_{kt}}{W_{kt-1}}\right)\right]^2
\]
在对称均衡中，工会的最优解会得出工资菲利普斯曲线：
\[
\log \left(1+\pi_t^w\right)=\kappa_w\left[\varphi N_t^{1+v}-\frac{\left(1-\tau_t\right) w_t N_t}{\mu_w} \int{e_{it}c_{it}^{-\sigma}di}\right]+\beta \log \left(1+\pi_{t+1}^w\right)
\]
其中 \(1 + \pi_t^w = W_t /W_{t-1} = (1 + \pi_t) w_t/w_{t-1}\) 是时期 \(t\) 的工资通胀率。
货币政策与财政政策
货币当局实施泰勒规则，将稳定状态目标通胀率设定为零:
\[
i_t = r_{\ast} + \phi \pi_t + \phi_y \left(Y_t - Y_{\ast}\right) + m_t
\]
其中 \(m_t\) 是货币政策冲击，\(r_{\ast}\) and \(Y_{\ast}\) 分别为稳态的真实利率与 GDP。
名义利率、真实利率与通胀率通过费雪方程相互关联：
\[
1+r_t=\frac{1+i_{t-1}}{1+\pi_{t}}
\]
政府通过实施预算平衡来维持不变的债务 (\(B^g\)) 实际价值，并通过税收 \(\tau_t\) 来筹集外部购买支出 (\(G_t\)) 和债券利息支付所需的资金。
\[
\tau_t =\frac{r_{t}B^{g} + G_t}{w_t N_t}
\]
市场出清条件
资本市场出清：
\[
p_{t} + B^{g}=\int{\left(a_{it} + b_{it}\right)di}
\]
商品市场出清：
\[
Y_t = \int{c_{it}di} + \int{\Phi_{t}\left(a_{it}, a_{it-1}\right)di} + \omega \int{b_{it}di} + I_t + \psi^p \left(P_{t},P_{t-1}\right) + G_t
\]
注意，通过使用 \(N_t\) 来表示企业问题中的劳动力需求和工资菲利普斯曲线中的劳动力供给，我们隐含地强加了劳动力市场出清条件。
均衡方程组：转移路径
给定家庭状态的初始分布 \(\Gamma_0\)，完美预期均衡是：(1) 家庭状态分布 \(\{\Gamma_t\}_{t=0}^{\infty}\)；(2) 家庭值函数与策略函数 \(\{V_t,g_{c,t},g_{b,t},g_{a,t}\}_{t=0}^{\infty}\)；以及 (3) 总体数量与价格 \(\{Y_t,K_t,N_t,I_t,\psi^{p}_{t},D_t,B_t^h,w_t,mc_t,p_t,Q_t,\pi_t,\pi_{t}^{w},i_t,r_t,r_t^a,r_t^b,\tau_t\}_{t=0}^{\infty}\) 的序列，使其满足：
\[\begin{split}
Y_t = z_t K_{t-1}^\alpha N_{t}^{1-\alpha}
\\
I_t = K_{t} - \left(1 - \delta\right) K_{t-1} + \zeta \left(K_{t},K_{t-1}\right) K_{t-1}
\\
\psi^p_t =\frac{\mu_p}{\mu_p-1} \frac{1}{2 \kappa_p}\left[\log \left(1 + \pi_t\right)\right]^2 Y_t
\\
D_t = Y_t - w_t N_t - I_t - \psi^{p}_{t}
\\
w_t = (1-\alpha)\frac{Y_t}{N_t} mc_t
\\
\log(1+\pi_t) = \kappa_p \left(mc_t - \frac{1}{\mu_p} \right) + \frac{1}{1+r_{t+1}} \frac{Y_{t+1}}{Y_t} \log(1+\pi_{t+1})
\\
Q_t = 1 + \frac{1}{\delta \epsilon_I}\left(\frac{K_t-K_{t-1}}{K_{t-1}}\right)
\\
(1+r_{t+1})Q_{t} = \alpha \frac{Y_{t+1}}{K_t} mc_{t+1} - \left[\frac{K_{t+1}}{K_t} - (1-\delta) + \frac{1}{2\delta \epsilon_I}\left(\frac{K_{t+1} - K_t}{K_t}\right)^2\right] + \frac{K_{t+1}}{K_t}Q_{t+1}
\\
1 + \pi_t^w = (1 + \pi_t) w_t/w_{t-1}
\\
\log \left(1+\pi_t^w\right)=\kappa_w\left[\varphi N_t^{1+v}-\frac{\left(1-\tau_t\right) w_t N_t}{\mu_w} \int{e g_{c,t}^{-\sigma}d \Gamma_t}\right]+\beta \log \left(1+\pi_{t+1}^w\right)
\\
i_t = r_{\ast} + \phi \pi_t + \phi_y \left(Y_t - Y_{\ast}\right) + m_t
\\
1+r_t=\frac{1+i_{t-1}}{1+\pi_{t}}
\\
r_t^b = r_t - \omega
\\
\frac{D_{t+1} + p_{t+1}}{p_{t}} = 1 + r_{t+1}
\\
\Theta_{pt-1} = \frac{p_{t-1}}{p_{t-1} + B^{g} - B^{h}_{t-1}}
\\
1 + r_t^a = \Theta_{pt-1} \frac{D_t + p_t}{p_{t-1}} + (1 - \Theta_{pt-1}) \left(1 + r_t\right)
\\
\tau_t =\frac{r_{t}B^{g} + G_t}{w_t N_t}
\\
B^{h}_{t} = \int{g_{b,t} d \Gamma_t}
\\
p_{t} + B^{g}=\int{\left(g_{a,t} + g_{b,t}\right)d \Gamma_t}
\end{split}\]
同时，商品市场出清条件隐含在瓦尔拉斯定律之中。该系统可以进一步简化。请参阅下面定义简化方程组的 hmod 文件。
均衡方程组：稳态均衡
在用常数 \(\left(z_t, \pi_t, G_t\right) = \left(z_{\ast}, 0, G_{\ast}\right) \) 表示的稳态均衡中，均衡方程组可简化为：
\[\begin{split}
Y_{\ast} = z_{\ast} K_{\ast}^\alpha N_{\ast}^{1-\alpha}
\\
I_{\ast} = \delta K_{\ast}
\\
\psi^p_{\ast} = 0
\\
D_{\ast} = Y_{\ast} - w_{\ast} N_{\ast} - \delta K_{\ast}
\\
w_{\ast} = (1-\alpha)\frac{Y_{\ast}}{N_{\ast}} mc_{\ast}
\\
mc_{\ast} = \frac{1}{\mu_p}
\\
Q_{\ast} = 1
\\
r_{\ast} + \delta = \alpha \frac{Y_{\ast}}{K_{\ast}} mc_{\ast}
\\
\pi_{\ast}^w = 0
\\
\varphi N_{\ast}^{v} = \frac{\left(1-\tau_{\ast}\right) w_{\ast}}{\mu_w} \int{e g_{c,{\ast}}^{-\sigma}d \Gamma_{\ast}}
\\
i_{\ast} = r_{\ast}
\\
i_{\ast} = r_{\ast}
\\
r_{\ast}^b = r_{\ast} - \omega
\\
D_{\ast} = r_{\ast} p_{\ast}
\\
\Theta_{p{\ast}} = \frac{p_{\ast}}{p_{\ast} + B^{g} - B^{h}_{\ast}}
\\
r_{\ast}^a = r_{\ast}
\\
\tau_{\ast} =\frac{r_{\ast}B^{g} + G_{\ast}}{w_{\ast} N_{\ast}}
\\
B^{h}_{\ast} = \int{g_{b,{\ast}} d \Gamma_{\ast}}
\\
p_{\ast} + B^{g}=\int{\left(g_{a,{\ast}} + g_{b,{\ast}}\right)d \Gamma_{\ast}}
\end{split}\]
hmod 文件
模型可以由如下的 hmod 文件表示：
% The fixed (non-calibrated) parameter values are from Auclert, Bardóczy, Rognlie, and Straub (2021, Appendix B.3)
parameters beta sigma varphi nu chi0 chi1 omega ra r w tau N;
beta = 0.976; % Calibrated Parameter: HANS final calibrated value is 0.977599682193789
sigma = 2.0; % Fixed Parameter
% Note: a typo of varphi in Table B.III: varphi = 2.073 in the table,
% but the actual varphi value in the SSJ code is around 1.713
varphi = 1.713; % Calibrated Parameter: HANS final calibrated value is 1.681157908161363
nu = 1.0; % Fixed Parameter
chi0 = 0.25; % Fixed Parameter
chi1 = 6.416; % Calibrated Parameter: HANS final calibrated value is 6.518415164643450
omega = 0.005; % Fixed Parameter
var_shock e;
shock_trans = [
0.966289 0.033422 0.000289
0.016711 0.966578 0.016711
0.000289 0.033422 0.966289
]; % need this to be the full transition matrix
e = [0.18315644, 0.67277917, 2.47128522];
% State Variable grid
var_state b a;
% For replication, the following grids were used in SSJ
% b = 10.^(linspace(log10(0.25), log10(50+0.25), 50)) - 0.25;
% a = 10.^(linspace(log10(0.25), log10(4000+0.25), 70)) - 0.25;
b = 10.^(linspace(log10(0.25), log10(40+0.25), 50)) - 0.25;
a = 10.^(linspace(log10(0.25), log10(120+0.25), 70)) - 0.25;
% Precompute cash on hand, as it does not change over vfi
var_pre_vfi coh;
coh = (1+ra)*a + (1+r-omega)*b + (1-tau)*w*e*N;
% Declare three choice variables in vfi and provide initial guess
var_policy bp ap c;
initial bp 0.0;
initial ap (1+ra)*a;
initial c (1+r-omega)*b + (1-tau)*w*e*N;
% Compute additional variables from vfi for later use:
% (a) chi is used later in Goods Market Clearing Condition
% (b) uce is used later in the Wage Phillips Curve
var_aux chi uce;
vfi;
% Define constraint equations
chi = (chi1/2)*(ap-(1+ra)*a)^2/((1+ra)*a+chi0);
%c + ap + bp == (1+ra)*a + (1+r-omega)*b + (1-tau)*w*e*N - chi;
c + ap + bp == coh - chi;
% Define the Bellman Equation: note two endogenous state variables (bp,ap)
u = c^(1.0 - sigma)/(1.0 - sigma) - varphi*N^(1.0+nu)/(1.0+nu);
Tv = u + beta*EXPECT(v(bp,ap));
% Define Bound Constraints
bp >= 0.0;
ap >= 0.0;
c >= 1e-8;
% post evaluation after value function iteration
uce = e*(c^(-sigma)); % efficiency unit weighted marginal utility used in Wage Phillips Curve
end;
% Endogenous Aggregate Variables
var_agg r pii Y K N w p d Bh;
% Fixed Parameters
delta = 0.02; % depreciation
epsI = 4.0; % capital adjustment cost parameter
kappap = 0.1; % price adjustment cost parameter
muw = 1.1; % wage market power parameter
kappaw = 0.1; % wage adjustment cost parameter
% Calibration Targets: Six Variables
r = 0.0125; % calibration target
Y = 1.0; % calibration target
K = 10.0; % calibration target
N = 1.00; % calibration target
p = 11.2; % calibration target
Bh = 1.04; % calibration target
% Parameters implied by Calibration Targets
d = 0.14; % d = p*r = 11.2*0.0125 = 0.14
ra = r; % steady state no-arbitrage condition
Bg = 2.8; % Debt/GDP ratio
G = 0.2; % G/GDP ratio
% Calibrated Parameters: 3 preference parameters above and 3 production parameters below
Z = 0.468; % HANKS Final Calibrated Parameters: Z = 0.467789814531232
alpha = 0.33; % HANKS Final Calibrated Parameters: alpha = 0.329949238578680
mup = 1.015; % HANKS Final Calibrated Parameters: mup = 1.015228426395939;
% Parameters implied by Calibration Targets and Calibrated Parameters
mc = 0.9850; % mu = 1/mup
w = 0.66; % w = (1-alpha)*Y/N * mc
% Taylor Rule Parameters
pii = 0.0000; % inflation
phi = 1.5; % Taylor rule
phi_y = 0.0; % Taylor rule
pi_star = 0.0; % Taylor rule
rstar = 0.0125; % Taylor rule
Y_star = 1.0; % Taylor rule
var_agg_shock m_shock;
m_shock = 0.0; % placeholder for monetary shocks
% Calibration Block: Six Calibration Targets (r,Y,K,N,p,Bh) with Six Parameters (beta,varphi,chi1,Z,alpha,mup),
% among which four parameters (varphi,Z,alpha,mup) can be solved analytically.
% Hence, calibration involves an equation solver with two variables and two equilibrium conditions
model_cali(beta, chi1);
% The calibration of mc follows from the following analytical conditions:
% d = rp and d = Y- w N - delta K implies
% rp = Y- w N - delta K so that r (p - K) = Y- w N - ( r + delta) K
% Now use ( r + delta) = mc*alpha*Y/K and w = mc * (1 - alpha) * Y / N to get
% r (p - K) = Y - mc Y so that mc = 1 - r (p - K) / Y;
mc = 1 - r * (p - K) / Y;
mup = 1 / mc; % steady state Price Phillips Curve
alpha = (r + delta) * K / Y / mc; % stead state Valuation Equation: r + delta = mc*alpha*Y/K
Z = Y/(K^alpha*N^(1-alpha)); % production function: Y = Z K^alpha N^(1-alpha)
w = mc * (1 - alpha) * Y / N; % marginal cost function: mc = w/((1-alpha)*Y/N)
ra = r; % no-arbitrage condition
tau = (r*Bg + G)/(w*N); % government budget
% equilibrium conditions
ap + bp == p + Bg; % asset demand = asset supply
Bh == bp; % liquid asset supply = liquid asset demand
% Bound Constraints for Calibrated Parameters
beta <= 1 - 1e-10;
beta >= 1e-10;
chi1 >= 1e-10;
% post evaluation
varphi =((1-tau)*w*uce/muw)/(N.^nu); % Wage Phillips Curve in Steady State
Y_star = Y; % update Y_star in Taylor Rule
end;
% Stationary Equilibrium: An equation solver with two variables and two equilibrium conditions
model_ss(r, Y);
ra = r; % no-arbitrage condition
mc = 1/mup; % steady state Price Phillips Curve
% Infer K from the stead state Valuation Equation: r + delta = mc*alpha*Y/K;
K = mc*alpha*Y/(r + delta);
% Infer N from the production function: Y = z K^alpha N^(1-alpha)
N = (Y/Z/(K^alpha))^(1/(1-alpha));
% Infer w from the marginal cost function: mc = w/((1-alpha)*Y/N)
w = mc*(1-alpha)*Y/N;
% Infer d from I = delta*K and zero steady state price adjustment cost
d = Y - w*N - delta*K;
% Infer p from asset pricing equation
p = d/r;
% Infer tau from government budget
tau = (r*Bg + G)/(w*N);
% equilibrium conditions
ap + bp == p + Bg; % asset demand = asset supply
N == ((1-tau)*w*uce/(varphi*muw)).^(1/nu); % Wage Phillips Curve
% Bound Constraints for Equilibrium Variables
r <= (1/beta - 1 - 1e-10);
Y >= 0.5;
% post evaluation
rstar = r; % update rstar in Taylor Rule
Y_star = Y; % update Y_star in Taylor Rule
pii = pi_star; % exogenously fixed steady state inflation
Bh = bp; % liquid asset
end;
% Transition Path Equilibrium Block: Nine Equations with Nine Unknowns
% The Equation System can be further simplified by utilizing the analytical structure.
% We keep it in the current form for ease of presentation.
model;
ii = rstar + phi*pii(-1) + phi_y*(Y(-1) - Y_star) + m_shock(-1); % Taylor Rule: recall ii was determined from last period
pshare = p(-1)/(p(-1)+Bg-Bh(-1)); % equity share of mutual fund portfolio
ra = pshare*(d+p)/p(-1) + (1- pshare)*(1+r) - 1; % the return of illiquid asset
I = K - (1-delta)*K(-1) + ((K/K(-1) - 1).^2)*K(-1)/(2*delta*epsI); % capital accumulation equation
psip = ((log(1+pii))^2)*Y*mup/((mup - 1)*2*kappap); % price adjustment cost
Q = 1 + (K/K(-1) - 1)/(delta*epsI); % Tobin's Q from the investment equation
tau = (r*Bg + G)/(w*N); % government budget constraint
% equilibrium conditions: 9 equations with 9 unknowns
r == (1.0 + ii)/(1.0 + pii) - 1.0; % Fisher Equation
ap + bp == p + Bg; % asset demand = asset supply
log(1 + pii) + log(w/w(-1)) == kappaw*(varphi*N^(1+nu) - (1-tau)*w*N*uce/muw) + beta*(log(1+pii(+1)) + log(w(+1)/w)); % Wage Phillips Curve
log(1 + pii) == kappap*(w*N/((1-alpha)*Y) - 1/mup) + log(1+pii(+1))*Y(+1)/Y/(1+r(+1)); % Price Phillips Curve
% valuation equation
(1 + r(+1))*(1 + (K/K(-1) - 1)/(delta*epsI)) == alpha*w(+1)*N(+1)/((1-alpha)*K) - (K(+1)/K - (1-delta) + (K(+1)/K - 1).^2/(2*delta*epsI)) + (1 + (K(+1)/K - 1)/(delta*epsI))*K(+1)/K;
Y == Z*(K(-1)^alpha)*(N^(1-alpha)); % production function
d == Y - w*N - I - psip; % dividend equation
0 == (d(+1) + p(+1)) - p*(1.0 + r(+1)); % no-arbitrage equation
bp == Bh; % liquid asset supply = demand
% Bound Constraints for Unknowns
r <= (1/beta - 1 - 1e-10);
pii >= beta - 1.0 + 1e-10; % 1+r <= 1/beta implies 1+pii >=beta*(1+ii)>= beta
N >= 0.5; % tau <=1 implies N > = (r*Bg + G)/w = 0.3556 when r = 0.0125
% post evaluation
goods_mkt = c + I + G + psip + chi + omega*bp - Y; % residual = aggregate demand - aggregate supply;
end;
使用工具箱
解析 .hmod 文件
假设 hank2_ssj.hmod 已经位于工作文件夹之中。在 MATLAB 中调用：
clear_hans
hans hank2_ssj
如上所示，HANS 解析了 hmod 文件，编译了 mex 文件，并且返回了一个描述模型的结构体。
第一行中的 clear_hans 清空了之前生成的所有缓存文件。
校准与稳态均衡
调用生成的代码 solve_cali 来校准模型参数，使其满足稳态均衡的数据目标。通过构造，校准的结果也返回了校准参数对应的稳态均衡。
tic;
cali_rslt = solve_cali;
fprintf('Time used for Calibration:\n');
toc;
我们可以查看 cali_rslt 中的数值：
cali_rslt
在计算出所有的参数及其对应的稳态均衡之后，我们继续进行转移路径的计算。然而，作为一项额外的验证措施，我们可以通过将返回的 cal_rslt 结构传递给 solve_ss 函数来求解稳态均衡：
tic;
ss_rslt = solve_ss(cali_rslt);
fprintf('Time used for solving stationary equilibrium:\n');
toc;
正如预期的那样，初始条件已经满足了稳态均衡的所有条件。因此，计算在初始迭代之后立即停止。
所有的稳态均衡的信息都在结构 ss_rslt 中返回。
ss_rslt
我们可以检查均衡变量的值并验证其满足了所有校准目标。
ss_rslt.var_agg
为了直观展现结果，我们可以将给定的个体冲击下的策略函数绘制在一个三维图形中：
% Pick an idiosyncratic shock and the corresponding (b,a) mesh in 2-D
shock_idx = 2;
b_mesh_2d = squeeze(ss_rslt.vfi_rslt.var_state_mesh.b_mesh(shock_idx,:,:));
a_mesh_2d = squeeze(ss_rslt.vfi_rslt.var_state_mesh.a_mesh(shock_idx,:,:));
% Plot the VFI Policy in 3-D Mesh
surf(b_mesh_2d,a_mesh_2d,squeeze(ss_rslt.vfi_rslt.var_policy.bp(shock_idx,:,:)))
xlabel('b');
ylabel('a');
title('Policy functions for b');
surf(b_mesh_2d,a_mesh_2d,squeeze(ss_rslt.vfi_rslt.var_policy.ap(shock_idx,:,:)))
xlabel('b');
ylabel('a');
title('Policy functions for a');
surf(b_mesh_2d,a_mesh_2d,squeeze(ss_rslt.vfi_rslt.var_policy.c(shock_idx,:,:)));
xlabel('b');
ylabel('a');
title('Policy functions for c');
我们还可以绘制（b, a）上的联合稳态分布，以及b或a的边际分布。
% Extract (b,a) mesh in 2-D
b_mesh_2d = squeeze(ss_rslt.vfi_rslt.var_state_mesh.b_mesh(1,:,:));
a_mesh_2d = squeeze(ss_rslt.vfi_rslt.var_state_mesh.a_mesh(1,:,:));
bgrid = ss_rslt.vfi_rslt.var_state.b(:); % b grid
agrid = ss_rslt.vfi_rslt.var_state.a(:); % a grid
% Plot the histogram of the Stationary Distribution in 3-D Mesh
dist_2d = squeeze(sum(ss_rslt.dist)); % histogram on (b,a) grid
dist_b = sum(dist_2d,2); % histogram on b grid
dist_a = reshape(sum(dist_2d,1),[],1); % histogram on a grid
surf(b_mesh_2d,a_mesh_2d,dist_2d);
xlabel('b');
ylabel('a');
title('Histogram of Stationary Distribution in (b,a)');
plot(bgrid,dist_b);
xlabel('b');
ylabel('prob');
title('Histogram of Stationary Distribution in b');
plot(agrid,dist_a);
xlabel('a');
ylabel('prob');
title('Histogram of Stationary Distribution in a');
求解暂时性冲击发生后的转移路径
求解线性化的转移路径
我们首先构造针对泰勒规则的货币政策冲击序列：
T = 300; % Transition periods
rho_r = 0.6; % persistence of monetary shock
sig_r = -0.01/4; % a negative interest rate shock of 25 basis points
m_shock_t = sig_r * rho_r.^(0:(T-1)); % the path of interest rate shock
figure;
plot(m_shock_t)
我们随后将货币冲击与求解出的 solved ss_rslt 传递到 solve_trans_linear 中，以求解线性化的转移路径。
options = struct;
options.T = T;
options.m_shock_t = m_shock_t;
tic;
trans_linear_rslt = solve_trans_linear(ss_rslt, options);
fprintf('Time used for solving the linearized transition path:\n');toc;
trans_linear_rslt
我们可以画出转移路径 (注意 K 为下一期资本) 。
figure;
plot(trans_linear_rslt.var_agg_t.K)
title('K');
xlabel('time');
求解非线性转移路径
现在，我们将 solved ss_rslt 传递到 solve_nonlinear_trans (init_ss, final_ss, options) 之中，以求解非线性转移路径。
注意，我们将 ss_rslt 传递到了 init_ss 和 final_ss 之中。
tic;
trans_nonlinear_rslt = solve_trans_nonlinear(ss_rslt, ss_rslt, options);
fprintf('Time used for solving the nonlinear transition path:\n');toc;
trans_nonlinear_rslt
我们可以比较非线性与线性的转移路径结果。
figure;
plot(trans_linear_rslt.var_agg_t.K);
hold on;
plot(trans_nonlinear_rslt.var_agg_t.K, '--');
xlabel('time');
title('K');
legend({'Linear','Nonlinear'});
我们也可以查看转移路径中定义的其他均衡变量。
figure;
plot(trans_nonlinear_rslt.var_agg_t.w);
title('w');
xlabel('time');
