Guvenen (2009): Asset Pricing with Heterogeneous IES

Guvenen (2009) constructs a two-agent model to account for several salient features of asset pricing moments, such as high risk premium, low and relatively smooth interest rate, and countercyclical movements in risk premium and Sharpe ratio. The two key ingredients of his model are limited stock market participation and heterogeneity in the elasticity of inter-temporal substitution in consumption (EIS).

When all endogenous states are pre-determined, such as capital and bond share in the current example, the expectation of future variables can be parameterized and integrated before solving the system of equilibrium equations, which greatly speeds up the computation. We illustrate how to use parameterized expectation with the toolbox in this example.

The Model

There are two types of infinitely-lived agents: stockholders (\(h\)) with measure \(\mu\), and non-stockholders (\(n\)) with measure \(1-\mu\). Agents have Epstein-Zin utility functions

\[U_{i,t}=\left\{ \left(1-\beta\right)c_{i,t}^{1-\rho^{i}}+\beta\left[\mathbb{E}_{t}\left(U_{i,t+1}^{1-\alpha}\right)\right]^{\frac{1-\rho^{i}}{1-\alpha}}\right\} ^{1/\left(1-\rho^{i}\right)}.\label{eq:utility Guv}\]

for \(i=h,n\). Most importantly, \(\rho^{h}<\rho^{n}\), i.e., the non-stockholders have lower EIS which is inversely proportional to \(\rho^{i}\), and thus they have higher desire for consumption smoothness over time. Each agent has one unit of labor endowment.

Stockholders can trade stock \(s_{t}\) and bond \(b_{h,t}\) at prices \(P_{t}^{s}\) and \(P_{t}^{f}\) respectively. Their budget constraint is

\[c_{h,t}+P_{t}^{f}b_{h,t+1}+P_{t}^{s}s_{t+1}\leq b_{h,t}+s_{t}\left(P_{t}^{s}+D_{t}\right)+W_{t},\]

where \(W_t\) is the labor income and borrowing constraint is

\[b_{h,t+1}\geq-\underline{B},\]

and in calibration \(\underline{B}\) is set to six times of the average monthly wage rate. The non-stockholders have the same constraints. In addition, they are restricted from trading stocks.

A representative firm produces the consumption good using capital \(K_{t}\) and labor \(L_{t}\) based on a Cobb-Douglas production function:

\[Y_{t}=Z_{t}K_{t}^{\theta}L_{t}^{1-\theta},\]

and the technology evolves according to an AR(1) process:

\[\ln Z_{t+1}=\phi\ln Z_{t}+\varepsilon_{t+1},\text{ }\varepsilon\overset{i.i.d.}{\sim}N\left(0,\sigma_{\varepsilon}^{2}\right).\]

The firm maximizes its value \(P_{t}^{s}\) expressed as the sum of its future dividends \(\left\{ D_{t+j}\right\} _{j=1}^{\infty}\) discounted by the shareholders’ marginal rate of substitution process:

\[P_{t}^{s}=\max_{\left\{ I_{t+j},L_{t+j}\right\} }\mathbb{E}_{t}\left[\sum_{j=1}^{\infty}\beta^{j}\frac{\Lambda_{h,t+j}}{\Lambda_{h,t}}D_{t+j}\right].\]

The firm accumulates capital subject to a concave adjustment cost function in investment:

\[K_{t+1}=\left(1-\delta\right)K_{t}+\varPhi\left(\frac{I_{t}}{K_{t}}\right)K_{t}.\]

Each period, the firm sells one-period bonds at price \(P_{t}^{f}\). The bond supply is constant and equals to \(\chi\) fraction of its average capital stock \(\bar{K}\). Thus dividend \(D_{t}\) can be written as

\[D_{t}=Z_{t}K_{t}^{\theta}L_{t}^{1-\theta}-W_{t}L_{t}-I_{t}-\left(1-P_{t}^{f}\right)\chi\bar{K}.\]

A sequential competitive equilibrium is given by sequences of allocations \(\{c_{i,t},b_{i,t+1},s_{t+1},I_{t},K_{t+1},L_{t}\}\) \(i=h,n\) and prices \(\left\{ P_{t}^{s},P_{t}^{f},W_{t}\right\}\) such that (i) given the price sequences, \(\left\{ c_{i,t},b_{i,t+1},s_{t+1}\right\}\) \(i=h,n\) solve the stockholders’ and non-stockholders’ optimization problems; (ii) Given the wage sequence \(\left\{ W_{t}\right\}\) and the law of motion for capital, \(\left\{L_{t},I_{t}\right\}\) are optimal for the representative firm; (iii) all markets clear:

\[\begin{split}\mu b_{h,t+1}+\left(1-\mu\right)b_{n,t+1} & =\chi\bar{K},\label{eq:MC bond Guv}\\ \mu s_{t+1} & =1,\label{eq:MC stock Guv}\\ L_{t} & =1,\nonumber \\ \mu c_{h,t}+\left(1-\mu\right)c_{n,t}+I_{t} & =Y_{t}.\nonumber\end{split}\]

We compute a recursive equilibrium using \(\left\{ K_{t},B_{t}^{n},Z_{t}\right\}\) as the aggregate state variables, where \(B_{t}^{n}=\left(1-\mu\right)b_{n,t}\) is total bond holding by the non-stockholders. We have \(8\) unknowns: {\(c_{h,t}\), \(c_{n,t}\), \(I_{t}\), \(B_{t+1}^{n}\), \(\lambda_{h,t}\), \(\lambda_{n,t}\), \(P_{t}^{s}\), \(P_{t}^{f}\)}, and the \(8\) equations used to solve them are:

  • Euler equations for bond holding:

\[P_{t}^{f}=\beta\left(1+\lambda_{i,t}\right)\mathbb{E}_{t}\left(\frac{\Lambda_{i,t+1}}{\Lambda_{i,t}}\right),\text{ }\forall i=h,n.\]

  • Euler equations for the stockholders’ demand of equity:

\[P_{t}^{s} =\beta\mathbb{E}_{t}\left[\frac{\Lambda_{h,t+1}}{\Lambda_{h,t}}\left(P_{t+1}^{s}+D_{t+1}\right)\right].\]

  • Slackness condition of borrowing limit:

\[\lambda_{i,t}\left(b_{i,t+1}+\underline{B}\right)=0,\text{ }\forall i=h,n.\]

  • The budget constraints (imposing \(s_{t+1}=1/\mu\)):

\[c_{h,t}+P_{t}^{f}b_{h,t+1}+\frac{P_{t}^{s}}{\mu}=P_{t}^{s}+D_{t}+\frac{\chi\bar{K}-B_{t}^{n}}{\mu}+W_{t},\]
\[c_{n,t}+P_{t}^{f}b_{n,t+1}=\frac{B_{t}^{n}}{1-\mu}+W_{t}.\]

  • Firm’s optimal capital accumulation \(K_{t+1}\):

\[q_{t} =\beta\mathbb{E}_{t}\left\{ \frac{\Lambda_{h,t+1}}{\Lambda_{h,t}}\left[\theta Z_{t}K_{t}^{\theta-1}-\frac{I_{t+1}}{K_{t+1}}+q_{t+1}\left(1-\delta+\varPhi\left(\frac{I_{t+1}}{K_{t+1}}\right)\right)\right]\right\} ,\label{eq:capital price Guv}\]

in which capital price \(q_{t}\) is the Lagrangian multiplier on the capital formation equation and satisfies

\[q_{t}\varPhi^{\prime}\left(\frac{I_{t}}{K_{t}}\right)=1.\label{eq:optimal inv Guv}\]

The auxiliary variables can be determined by the utility function, market clearing conditions, and the following two equations:

\[W_{t}= \left(1-\theta\right)Z_{t}\left(\frac{K_{t}}{L_{t}}\right)^{\theta},\]
\[\beta\frac{\Lambda_{i,t+1}}{\Lambda_{i,t}} =\beta^{\frac{1-\alpha}{1-\rho^{i}}}\left(\frac{c_{i,t+1}}{c_{i,t}}\right)^{-\rho^{i}}\left[\frac{\frac{U_{i,t+1}}{c_{i,t}}}{\left[\left(\frac{U_{i,t}}{c_{i,t}}\right)^{1-\rho^{i}}-(1-\beta)\right]^{1/\left(1-\rho^{i}\right)}}\right]^{\rho^{i}-\alpha}.\]

In period \(t\), the 6 future variables in use: \(c_{h,t+1}\), \(c_{n,t+1}\), \(P_{t+1}^{s}+D_{t+1}\), \(I_{t+1}/K_{t+1}\), \(U_{h,t+1}\) and \(U_{n,t+1}\) are functions of \(\left\{ K_{t+1},B_{t+1}^{n},Z_{t+1}\right\}\) and are solved from the previous iteration. Similar to Guvenen (2009), the initial guesses for these functions are obtained by solving a version of the model with no leverage.

The gmod File

The model can be solved with the following guvenen2009.gmod file

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
% Parameters
parameters beta alpha rhoh rhon theta delta mu xsi chi a1 a2 Kss Bbar bn_shr_lb bn_shr_ub varianceScale;

beta = 0.9966;      % discount factor
alpha = 6;          % risk aversion
rhoh = 1/.3;        % inv IES for stockholders
rhon = 1/.1;        % inv IES for non-stockholders
theta = .3;         % capital share
delta = .0066;      % depreciation rate
mu = .2;            % participation rate
xsi = .4;           % adjustment cost coefficient
chi = .005;         % leverage ratio
a1 = (((delta^(1/xsi))*xsi)/(xsi-1));
a2 = (delta/(1-xsi));
Kss = ((1/beta-1+delta)/theta)^(1/(theta-1));
Bbar = -0.1*(1-theta)*Kss^theta; %borrowing constraint
varianceScale = 1e4;

TolEq = 1e-4;
INTERP_ORDER = 4; EXTRAP_ORDER = 4;
PrintFreq = 100;
SaveFreq = inf;

% Shocks
var_shock Z;
shock_num = 15;
phi_z = 0.984;  % productivity AR(1)
mu_z = 0;
sigma_e = 0.015/(1+phi_z^2+phi_z^4).^0.5;
[z,shock_trans,~]=tauchen(shock_num,mu_z,phi_z,sigma_e,2);
Z = exp(z);

% States
var_state K bn_shr;
K_pts = 10;
K = exp(linspace(log(.84*Kss),log(1.2*Kss),K_pts));

bn_shr_lb = (1-mu)*Bbar/(chi*Kss);
bn_shr_ub = (chi*Kss - mu*Bbar)/(chi*Kss);
b_pts = 30;
bn_shr = linspace(bn_shr_lb,bn_shr_ub,b_pts);

% Last period
var_policy_init c_h c_n;

inbound_init c_h 1e-6 100;
inbound_init c_n 1e-6 100;

var_aux_init Y W vh vn vhpow vnpow Ps Pf Div Eulerstock Eulerbondh Eulerbondn Inv dIdK Eulerf;

model_init;  
  Y = Z*(K^theta);
  W = (1-theta)*Z*(K^theta);
  resid1 = 1 - (W + (bn_shr*chi*Kss/(1-mu)))/c_n; % c_n: individual consumption
  resid2 = 1 - (W + (Div/mu) + ((1-bn_shr)*chi*Kss/mu))/c_h; % c_h: individual consumption  
  vh = ((1-beta)*(c_h^(1-rhoh)))^(1/(1-rhoh));
  vn = ((1-beta)*(c_n^(1-rhon)))^(1/(1-rhon));
  vhpow = vh^(1-alpha);
  vnpow = vn^(1-alpha);
  Pf = 0;
  Ps = 0;
  Div = Y - W - (1-Pf)*chi*Kss;  % investment is zero          

  Eulerstock = (vh^(rhoh-alpha))*(c_h^-rhoh)*(Ps + Div);
  Eulerbondh = (vh^(rhoh-alpha))*(c_h^-rhoh);
  Eulerbondn = (vn^(rhon-alpha))*(c_n^-rhon);

  Inv = 0;     
  Knext = 0;                                             
  dIdK = (Inv/K) - (1/a1)*(xsi/(xsi-1))*(Inv/(K*((1/a1)*((Knext/K)-(1-delta)-a2))))*(Knext/K);
  Eulerf = (vh^(rhoh-alpha))*(c_h^-rhoh)*(theta*Z*(K^(theta-1)) - dIdK);
  
  equations;
     resid1;
     resid2;
  end;
end;

var_interp EEulerstock_interp EEulerbondh_interp EEulerbondn_interp EEulerf_interp Evh_interp Evn_interp EPD_interp EPD_square_interp;
initial EEulerstock_interp shock_trans*reshape(Eulerstock,shock_num,[]);
initial EEulerbondh_interp shock_trans*reshape(Eulerbondh,shock_num,[]);
initial EEulerbondn_interp shock_trans*reshape(Eulerbondn,shock_num,[]);
initial EEulerf_interp shock_trans*reshape(Eulerf,shock_num,[]);
initial Evh_interp shock_trans*reshape(vhpow,shock_num,[]);
initial Evn_interp shock_trans*reshape(vnpow,shock_num,[]);
initial EPD_interp shock_trans*reshape(Div,shock_num,[]);
initial EPD_square_interp shock_trans*reshape(Div.^2,shock_num,[]) / varianceScale;

EEulerstock_interp = shock_trans*Eulerstock;
EEulerbondh_interp = shock_trans*Eulerbondh;
EEulerbondn_interp = shock_trans*Eulerbondn;
EEulerf_interp = shock_trans*Eulerf;
Evh_interp = shock_trans*vhpow;
Evn_interp = shock_trans*vnpow;
EPD_interp = shock_trans*(Ps+Div);
EPD_square_interp = shock_trans*(Ps+Div).^2 / varianceScale;

% Endogenous variables, bounds, and initial values
var_policy c_h c_n Ps Pf Inv bn_shr_next lambdah lambdan;

inbound c_h 1e-3 10;
inbound c_n 1e-3 10;
inbound Ps 1e-3 500;
inbound Pf 1e-3 10;
inbound Inv 1e-9 10;
inbound bn_shr_next bn_shr_lb bn_shr_ub;
inbound lambdah 0 2;
inbound lambdan 0 2;

% Other equilibrium variables
var_aux Y W b_h b_n Div dIdKp Eulerstock Eulerbondh Eulerbondn dIdK Eulerf vhpow vnpow omega PDratio Rs R_ep vh vn Knext std_ExcessR SharpeRatio;

model;
  Y = Z*(K^theta);              % output
  W = (1-theta)*Z*(K^theta);    % Wage = F_l
  Div = Y - W - Inv - (1-Pf)*chi*Kss;            % dividends

  Knext = (1-delta)*K + (a1*((Inv/K)^((xsi-1)/xsi))+a2)*K;
  dIdKp = (1/a1)*(xsi/(xsi-1))*(Inv/(K*((1/a1)*((Knext/K)-(1-delta)-a2))));

  b_h = (1-bn_shr)*chi*Kss/mu;
  b_n = bn_shr*chi*Kss/(1-mu);
    
  [EEulerstock_future,EEulerbondh_future,EEulerbondn_future,EEulerf_future,Evh_future,Evn_future,EPD_future,EPD_square_future] = GDSGE_INTERP_VEC(shock,Knext,bn_shr_next);
  EPD_square_future = EPD_square_future*varianceScale;
  
  vh =  ((1-beta)*(c_h^(1-rhoh)) + beta*(Evh_future^((1-rhoh)/(1-alpha))))^(1/(1-rhoh));
  vn =  ((1-beta)*(c_n^(1-rhon)) + beta*(Evn_future^((1-rhon)/(1-alpha))))^(1/(1-rhon));

  Eulerstock = (vh^(rhoh-alpha))*(c_h^-rhoh)*(Ps + Div);
  Eulerbondh = (vh^(rhoh-alpha))*(c_h^-rhoh);
  Eulerbondn = (vn^(rhon-alpha))*(c_n^-rhon);
  
  dIdK = (Inv/K) - (1/a1)*(xsi/(xsi-1))*(Inv/(K*((1/a1)*((Knext/K)-(1-delta)-a2))))*(Knext/K);
  Eulerf = (vh^(rhoh-alpha))*(c_h^-rhoh)*(theta*Z*(K^(theta-1)) - dIdK);

  vhpow = vh^(1-alpha);
  vnpow = vn^(1-alpha);

  omega = (Ps+Div+ mu*b_h)/(Ps+Div+chi*Kss);
  PDratio = Ps/Div;
  Rs = EPD_future/Ps;
  R_ep = Rs - 1/Pf;
  % The following inline implements 
  % std_ExcessR = (GDSGE_EXPECT{(PD_future'/Ps - Rs)^2})^0.5;
  std_ExcessR = (EPD_square_future/(Ps^2) + Rs^2 - 2*EPD_future*Rs/Ps)^0.5;
  SharpeRatio = R_ep/std_ExcessR;

  % Equations:
  err_bdgt_h = 1 - (W + (Div/mu) + b_h - Pf*(chi*Kss*(1-bn_shr_next)/mu))/c_h; % these are individual consumptions
  err_bdgt_n = 1 - (W + b_n - Pf*(bn_shr_next*chi*Kss/(1-mu)))/c_n;
  foc_stock = 1 - (beta*EEulerstock_future*(Evh_future^((alpha-rhoh)/(1-alpha))))/((c_h^(-rhoh))*Ps);
  foc_bondh = 1 - (beta*EEulerbondh_future*(Evh_future^((alpha-rhoh)/(1-alpha))) + lambdah)/((c_h^(-rhoh))*Pf);
  foc_bondn = 1 - (beta*EEulerbondn_future*(Evn_future^((alpha-rhon)/(1-alpha))) + lambdan)/((c_n^-rhon)*Pf);
  foc_f = 1 - (beta*EEulerf_future*(Evh_future^((alpha-rhoh)/(1-alpha))))/((c_h^(-rhoh))*dIdKp);
  
  slack_bn = lambdan*(bn_shr_next - bn_shr_lb);    %mun_lw*bn_shr_next;
  slack_bh = lambdah*(bn_shr_ub - bn_shr_next);    %mun_up*(1-bn_shr_next);

  equations;
    err_bdgt_h;
    err_bdgt_n;
    foc_stock;
    foc_bondh;
    foc_bondn;
    foc_f;
    slack_bn;
    slack_bh;
  end;

end;

simulate;
  num_periods = 10000;
  num_samples = 100;
  
  initial K Kss;
  initial bn_shr 0.5;  
  initial shock 2;
  
  var_simu Y c_h c_n Inv Ps Div Pf bn_shr_next Knext omega PDratio Rs R_ep SharpeRatio std_ExcessR;

  K' = Knext;
  bn_shr' = bn_shr_next;
end;

The implementation uses the exact parameterization and specification (number of grid points, domains of endogenous state variables, number of discretized exogenous shocks, the cubic spline approximation, etc.) in Guvenen (2009). However, the algorithm differs from his. His is based on the algorithm in Krusell and Smith (1998): starting from a conjectured law of motion for state-variables and pricing functions, he solves the agents’ Bellman equation and the agents’ policy functions using standard value function iterations. Then he uses these policy functions and temporary market clearing conditions to obtain a new law of motions and new pricing functions. The current algorithm, like other examples, utilizes the fact that agents’ optimization is concave, replaces the agent’s optimization with first order and complementarity conditions, and solves the agents’ optimization problem together with the market clearing conditions.

Since both state variables \(K_t\) and \(B_t^n\) are predetermined, the expectation of future variables (e.g., the marginal utility entering the bond’s Euler equation, or the price and dividend entering the equity’s Euler equation) can be calculated before-hand conditional on current exogenous states. For example, for the bond’s Euler equation

\[P_{t}^{f}=\beta\left(1+\lambda_{i,t}\right)\mathbb{E}_{t}\left(\frac{\Lambda_{i,t+1}}{\Lambda_{i,t}}\right),\text{ }\forall i=h,n.\]

we can define function

\[\tilde{\Lambda}(z,K',B') \equiv E[\Lambda(z',K',B')|z]\]

and approximate \(\tilde{\Lambda}(z,K',B')\) instead of \(\Lambda(z',K',B')\). This is made possible by that \((K,B)\) are predetermined in that the transition of the two states do not depend on \(z'\) (this is not true when the endogenous state follows an implicit state transition function that depends on the realization of \(z'\), such as wealth share and consumption share). When a projection method or spline is used to approximate \(\tilde{\Lambda}(z,K,B)\), the integration becomes quite simple, since projection is a linear operator and therefore the projection of an expectation equals to the expectation of the projection.

Once \(\tilde{\Lambda}(z,K',B')\) is approximated beforehand, then one only needs to evaluate \(\tilde{\Lambda}\) once for the current exogenous shock (\(z\)) in each evaluation of the system of equations. If instead one is to work with \(\Lambda(z',K',B')\), \(\Lambda\) needs to be evaluated for each \(z'\) and integrated using the transition matrix. This is why integrating \(\Lambda(z',K',B')\) beforehand is more computation-efficient, especially when the size of exogenous shocks is large.

In particular, the parameterized expectation is implemented by declaring the expectation of relevant variables in var_interp in

79
var_interp EEulerstock_interp EEulerbondh_interp EEulerbondn_interp EEulerf_interp Evh_interp Evn_interp EPD_interp EPD_square_interp;

specifying the updating of var_interp by integrating the corresponding variables using the transition matrix

89
90
91
92
93
94
95
96
EEulerstock_interp = shock_trans*Eulerstock;
EEulerbondh_interp = shock_trans*Eulerbondh;
EEulerbondn_interp = shock_trans*Eulerbondn;
EEulerf_interp = shock_trans*Eulerf;
Evh_interp = shock_trans*vhpow;
Evn_interp = shock_trans*vnpow;
EPD_interp = shock_trans*(Ps+Div);
EPD_square_interp = shock_trans*(Ps+Div).^2 / varianceScale;

and evaluating the approximation at the current shock only when evaluating the system of equations

124
  [EEulerstock_future,EEulerbondh_future,EEulerbondn_future,EEulerf_future,Evh_future,Evn_future,EPD_future,EPD_square_future] = GDSGE_INTERP_VEC(shock,Knext,bn_shr_next);

where the GDSGE_INTERP_VEC takes the first argument shock, a reserved keyword referring to the index of the current exogenous shock.

After processed by the online compiler, we can run the returned file to run policy iterations:

>> options = struct;
options.MaxIter = 500;
IterRslt = iter_guvenen2009(options);

Iter:100, Metric:0.512539, maxF:8.65033e-09
Elapsed time is 8.272149 seconds.
Iter:200, Metric:0.285722, maxF:9.39891e-09
Elapsed time is 7.549686 seconds.
Iter:300, Metric:0.174423, maxF:9.01752e-09
Elapsed time is 7.250628 seconds.
Iter:400, Metric:0.113139, maxF:8.89712e-09
Elapsed time is 6.938417 seconds.
Iter:500, Metric:0.0763611, maxF:8.61954e-09
Elapsed time is 6.762700 seconds.

We first solve a simple version of the problem as a warm-up. This simple problem has a lower aggregate net supply of bond (\(\chi=0.005\)) and a tighter borrowing constraint (\(\bar{B}\)). This ensures that starting from a last period solution the equilibrium exists, which is also noted in Guvenen (2009).

Starting with the warm-up solution, we now change the parameters, increase \(\chi\) to \(0.15\), its value in the paper, and specify a much more relaxed borrowing constraint (\(\bar{B}\)). This can be done very conveniently with the toolbox by specifying the new parameters into a structure, and the warm-up solution as a field named WarmUp into the structure.

>> beta = 0.9967;      % discount factor
theta = .3;         % capital share
delta = .0066;      % depreciation rate
mu = .2;            % participation rate
Kss = ((1/beta-1+delta)/theta)^(1/(theta-1));
chi = 0.15;
Bbar = -6*(1-theta)*Kss^theta;
bn_shr_lb = (1-mu)*Bbar/(chi*Kss);
bn_shr_ub = (chi*Kss - mu*Bbar)/(chi*Kss);
b_pts = 30;

options = struct;
options.chi = chi;
options.Bbar = Bbar; %borrowing constraint
options.bn_shr_lb = bn_shr_lb;
options.bn_shr_ub = bn_shr_ub;
options.bn_shr = linspace(bn_shr_lb,bn_shr_ub,b_pts);

options.WarmUp = IterRslt;
options.SkipModelInit=1;

We can now call the routine with the option structure as an argument

>> IterRslt = iter_guvenen2009(options);

Iter:600, Metric:0.0605821, maxF:9.85466e-09
Elapsed time is 6.361178 seconds.

...

Iter:2490, Metric:9.99635e-05, maxF:9.29708e-09
Elapsed time is 5.205234 seconds.

The policy iteration continues from where the warm-up solution ends, toward convergence.

We can inspect the policy functions at convergence by calling in MATLAB:

figure; clf;
subplot(1,2,1);
surf(KMesh,BMesh,squeeze(var_aux.R_ep(z1,:,:))*12);
title('Equity Premium');
xlabel('$K$','Interpreter','Latex');
ylabel('$B_n$','Interpreter','Latex');
set(gca,'FontSize',15);

subplot(1,2,2);
surf(KMesh,BMesh,squeeze((1./var_policy.Pf(z1,:,:))-1)*12);
title('Interest Rate');
xlabel('$K$','Interpreter','Latex');
ylabel('$B_n$','Interpreter','Latex');
set(gca,'FontSize',15);
../../_images/policy_guv.png

As shown, given capital, as the bond holding of non-stock holders increases, the demand for bond increases and for stock decreases, which results in a lower interest rate, lower stock price and lower equity premium.

We can then simulate the economy

>> SimuRslt = simulate_guvenen2009(IterRslt);
Periods: 1000
shock       K  bn_shr       Y     c_h     c_n     Inv      Ps     Div      Pfbn_shr_next   Knext   omega PDratio      Rs    R_epSharpeRatio
    8   126.1  0.3264   4.267   5.103   3.003  0.8448   111.7  0.3939  0.9978  0.3265   126.1  0.9521   283.7   1.0040.002244 0.04976
Elapsed time is 28.379607 seconds.
Periods: 2000
shock       K  bn_shr       Y     c_h     c_n     Inv      Ps     Div      Pfbn_shr_next   Knext   omega PDratio      Rs    R_epSharpeRatio
    9   127.3  0.3652   4.341   5.197    3.04   0.869   119.3  0.4121  0.9989  0.3655   127.3  0.9493   289.5   1.003 0.00224 0.05007
Elapsed time is 25.144899 seconds.
Periods: 3000
shock       K  bn_shr       Y     c_h     c_n     Inv      Ps     Div      Pfbn_shr_next   Knext   omega PDratio      Rs    R_epSharpeRatio
    11   135.2  0.5826   4.547   5.381   3.158  0.9444   136.5  0.4421   1.001  0.5829   135.3   0.928   308.7   1.0010.002332 0.05232
Elapsed time is 14.973382 seconds.

and inspect the distribution of state variables in the ergodic set

startPeriods = 500;
Kss = IterRslt.params.Kss;
K = SimuRslt.K(:,startPeriods:end);
omega = SimuRslt.omega(:,startPeriods:end);

figure; clf;
subplot(1,2,1);
histogram(K/Kss,40,'Normalization','probability');%,'BinLimits',[0.5,0.75]);
title('Distribution of $K$','Interpreter','Latex');
xlabel('$K/K^{ss}$','Interpreter','Latex');
ylabel('Fractions');
set(gca,'FontSize',15);

subplot(1,2,2);
histogram(omega,40,'Normalization','probability');%,'BinLimits',[0.5,0.75]);
title('Wealth Share Distribution','Interpreter','Latex');
xlabel('Wealth of Stockholders','Interpreter','Latex');
ylabel('Fractions');
set(gca,'FontSize',15);
../../_images/dist_guv.png

What’s Next?

This example illustrates how to work with parameterized expectation with the toolbox, which greatly reduces computation cost for models that can be solved with predetermined endogenous state variables. The finite-agent asset pricing model is an important class in the Macro and Macro-finance literature. Since the framework of the toolbox solves allocations and prices jointly respecting the optimality and market clearing conditions, it is able to deal with setups which are traditionally regarded as ill-conditioned, such as models with very high risk aversion levels (see example Barro et al. (2017)) or models with equilibrium asset returns that are close to collinear such as Cao, Evans, and Luo (2020).

Or you can directly proceed to Toolbox API.