Problem solving for steady state

Steady state
1

Im trying to solve for the steady states of the following system, originally with 28 variables. However, I have never been able to obtain the steady state, possibly because I’ve never been near a steady state enough. I tried substitution and reduce the system to 14 variables (now the equations are much longer) but still no luck. I was wondering if there is anything wrong with my model or is it just too complicated to be solved? Would there be any other choice for me other than simplifying the model?

Thank you very much.

Original model:

var c1 c2 d1 d2 la l1y l2y k1 k2 w lambc lamtfp lamer c1p c2p d1p d2p lap l1yp l2yp k1p k2p wp lambcp lamtfpp lamerp r g;
parameters nu phi phiprime chi epsilon eta mu omega alphaa delta pay payprime;

nu=0.5;
chi=1;
omega=1; 
alphaa=2/3; 
delta=0.025;

phi=0.1;
phiprime=0.1; 
epsilon=1.1; 
eta=1.1; 
mu=0.01222;
pay=0.1;
payprime=0.1;

model;
g=1+mu*(la(-1)^omega+lap(-1)^omega);
r=c1p/d2;

c1+c2+(d1+d2)/r+k1(1)+k2(1)+payprime*(k2(1)*g(1)-(1-delta)*k2)^2+phiprime*c1p+wp*l1yp=l1y^alphaa*k1^(1-alphaa)+l1yp^alphaa*k2^(1-alphaa)+(1-delta)*(k1+k2)+phi*d2/r+w*l2y/r;

r*(c1p+c2p)+d1p+d2p+k1p(1)+k2p(1)+pay*(k1p(1)*g(1)-(1-delta)*k1p)^2+phi*d2+w*l2y=l2y^alphaa*k1p^(1-alphaa)+l2yp^alphaa*k2p^(1-alphaa)+(1-delta)*(k1p+k2p)+r*phiprime*c1p+r*wp*l1yp;

(c1^nu+c2^nu+d1^nu+d2^nu)^(-1)*c1^(nu-1)=lambc;
(c1^nu+c2^nu+d1^nu+d2^nu)^(-1)*c2^(nu-1)=lambc;
(c1^nu+c2^nu+d1^nu+d2^nu)^(-1)*d1^(nu-1)=lambc/r;
(c1^nu+c2^nu+d1^nu+d2^nu)^(-1)*d2^(nu-1)+lambc*phi/r=lambc/r+r*lamer/d2;
chi*(la^epsilon+(l1y+l2y)^epsilon)^(eta/epsilon-1)*la^(epsilon-1)=lamtfp*mu*omega*la^(omega-1);
chi*(la^epsilon+(l1y+l2y)^epsilon)^(eta/epsilon-1)*(l1y+l2y)^(epsilon-1)=lambc*alphaa*l1y^(alphaa-1)*k1^(1-alphaa);
chi*(la^epsilon+(l1y+l2y)^epsilon)^(eta/epsilon-1)*(l1y+l2y)^(epsilon-1)=lambc*w/r;
wp=alphaa*l1yp^(alphaa-1)*k2^(1-alphaa);
lambc(1)/g(1)*((1-alphaa)*l1y(1)^alphaa*k1(1)^(-alphaa)+1-delta)=lambc;
lambc(1)/g(1)*((1-alphaa)*l1yp(1)^alphaa*k2(1)^(-alphaa)+1-delta+2*(1-delta)*payprime*(k2(2)*g(2)-(1-delta)*k2(1)))=lambc*(1+2*payprime*(k2(1)*g(1)-(1-delta)*k2));
lambc(1)/g(1)*(alphaa*l1y(1)^alphaa*k1(1)^(1-alphaa)+alphaa*l1yp(1)^alphaa*k2(1)^(1-alphaa)+payprime*(k2(2)*g(2)-(1-delta)*k2(1))^2)=lamtfp;
lambc*(d1+d2-phi*d2-w*l2y)/r/r=lamer;

(c1p^nu+c2p^nu+d1p^nu+d2p^nu)^(-1)*c1p^(nu-1)+lambcp*r*phiprime=r*lambcp-r*lamerp/c1p;
(c1p^nu+c2p^nu+d1p^nu+d2p^nu)^(-1)*c2p^(nu-1)=r*lambcp;
(c1p^nu+c2p^nu+d1p^nu+d2p^nu)^(-1)*d1p^(nu-1)=lambcp;
(c1p^nu+c2p^nu+d1p^nu+d2p^nu)^(-1)*d2p^(nu-1)=lambcp;
chi*(lap^epsilon+(l1yp+l2yp)^epsilon)^(eta/epsilon-1)*lap^(epsilon-1)=lamtfpp*mu*omega*lap^(omega-1);
chi*(lap^epsilon+(l1yp+l2yp)^epsilon)^(eta/epsilon-1)*(l1yp+l2yp)^(epsilon-1)=lambcp*wp*r;
chi*(lap^epsilon+(l1yp+l2yp)^epsilon)^(eta/epsilon-1)*(l1yp+l2yp)^(epsilon-1)=lambcp*alphaa*l2yp^(alphaa-1)*k2p^(1-alphaa);
w=alphaa*l2y^(alphaa-1)*k1p^(1-alphaa);
lambcp(1)/g(1)*((1-alphaa)*l2y(1)^alphaa*k1p(1)^(-alphaa)+1-delta+2*(1-delta)*pay*(k1p(2)*g(2)-(1-delta)*k1p(1)))=lambcp*(1+2*pay*(k1p(1)*g(1)-(1-delta)*k1p));
lambcp(1)/g(1)*((1-alphaa)*l2yp(1)^alphaa*k2p(1)^(-alphaa)+1-delta)=lambcp;
lambcp(1)/g(1)*(alphaa*l2y(1)^alphaa*k1p(1)^(1-alphaa)+alphaa*l2yp(1)^alphaa*k2p(1)^(1-alphaa)+pay*(k1p(2)*g(2)-(1-delta)*k1p(1))^2)=lamtfpp;
lambcp*(phiprime*c1p+wp*l1yp-c1p-c2p)=lamerp;
end;

initval;
c1=2.88926199417712;
c2=2.88926199417712;
d1=11.3887510905584;
d2=1.68808278057854;
la=54.3805316898628;
l1y=15.6912341135446;
l2y=239.42601433478;
k1=0.0831311719057078;
k2=0.135914733485534;
w=0.0359395974529077;
lambc=0.010665194907865;
lamtfp=0.0184754516380136;
lamer=0.026156226574055;
c1p=3.38573459383048;
c2p=4.07134273372575;
d1p=4.78174847361626;
d2p=4.78174847361626;
lap=64.6167176936905;
l1yp=180.548846969043;
l2yp=1.90965361764221;
k1p=0.204270712225458;
k2p=182.076409874019;
wp=0.0514779962433012;
lambcp=0.0101672864916624;
lamtfpp=0.0278131543470114;
lamerp=0.0165379700568582;
r=2.00554404996009;
g=8.70680223167034;

end;

steady (solve_algo=2, maxit=200000);

After substitution:

var c1 d2 la l1y l2y k1 k2 c1p d2p lap l1yp l2yp k1p k2p;
parameters nu phi phiprime chi epsilon eta mu omega alphaa delta pay payprime;

nu=0.5;
chi=1;
omega=1; 
alphaa=2/3; 
delta=0.025;

phi=0.1;
phiprime=0.1; 
epsilon=1.1; 
eta=1.1; 
mu=0.01222;
pay=0.1;
payprime=0.1;
 
model;
2*c1+((c1/((c1p/d2)^(1/(nu-1))))+d2)/(c1p/d2)+k1+k2+payprime*(k2*(1+mu*(la^omega+lap^omega))-(1-delta)*k2)^2+phiprime*c1p+(alphaa*l1yp^(alphaa-1)*k2^(1-alphaa))*l1yp=l1y^alphaa*k1^(1-alphaa)+l1yp^alphaa*k2^(1-alphaa)+(1-delta)*(k1+k2)+phi*d2/(c1p/d2)+(alphaa*l2y^(alphaa-1)*k1p^(1-alphaa))*l2y/(c1p/d2);
 
(c1p/d2)*(c1p+(d2p*(c1p/d2)^(1/(nu-1))))+2*d2p+k1p+k2p+pay*(k1p*(1+mu*(la^omega+lap^omega))-(1-delta)*k1p)^2+phi*d2+(alphaa*l2y^(alphaa-1)*k1p^(1-alphaa))*l2y=l2y^alphaa*k1p^(1-alphaa)+l2yp^alphaa*k2p^(1-alphaa)+(1-delta)*(k1p+k2p)+(c1p/d2)*phiprime*c1p+(c1p/d2)*(alphaa*l1yp^(alphaa-1)*k2^(1-alphaa))*l1yp;

d2^(nu-1)+c1^(nu-1)*phi/(c1p/d2)=c1^(nu-1)/(c1p/d2)+(c1p/d2)*(c1^(nu-1)*((c1/((c1p/d2)^(1/(nu-1))))+d2-phi*d2-(alphaa*l2y^(alphaa-1)*k1p^(1-alphaa))*l2y)/(c1p/d2)/(c1p/d2))/d2;
chi*(la^epsilon+(l1y+l2y)^epsilon)^(eta/epsilon-1)*la^(epsilon-1)=(((2*c1^nu+(c1/((c1p/d2)^(1/(nu-1))))^nu+d2^nu)^(-1)*c1^(nu-1))/(1+mu*(la^omega+lap^omega))*(alphaa*l1y^alphaa*k1^(1-alphaa)+alphaa*l1yp^alphaa*k2^(1-alphaa)+payprime*(k2*(1+mu*(la^omega+lap^omega))-(1-delta)*k2)^2))*mu*omega*la^(omega-1);
chi*(la^epsilon+(l1y+l2y)^epsilon)^(eta/epsilon-1)*(l1y+l2y)^(epsilon-1)=((2*c1^nu+(c1/((c1p/d2)^(1/(nu-1))))^nu+d2^nu)^(-1)*c1^(nu-1))*alphaa*l1y^(alphaa-1)*k1^(1-alphaa);
chi*(la^epsilon+(l1y+l2y)^epsilon)^(eta/epsilon-1)*(l1y+l2y)^(epsilon-1)=((2*c1^nu+(c1/((c1p/d2)^(1/(nu-1))))^nu+d2^nu)^(-1)*c1^(nu-1))*(alphaa*l2y^(alphaa-1)*k1p^(1-alphaa))/(c1p/d2);
((1-alphaa)*l1y^alphaa*k1^(-alphaa)+1-delta)=(1+mu*(la^omega+lap^omega));
((1-alphaa)*l1yp^alphaa*k2^(-alphaa)+1-delta+2*(1-delta)*payprime*(k2*(1+mu*(la^omega+lap^omega))-(1-delta)*k2))=(1+2*payprime*(k2*(1+mu*(la^omega+lap^omega))-(1-delta)*k2))*(1+mu*(la^omega+lap^omega));
 
c1p^(nu-1)+d2p^(nu-1)*(c1p/d2)*phiprime=c1p/d2*d2p^(nu-1)-c1p/d2*(d2p^(nu-1)*(phiprime*c1p+(alphaa*l1yp^(alphaa-1)*k2^(1-alphaa))*l1yp-c1p-(d2p*(c1p/d2)^(1/(nu-1)))))/c1p;
chi*(lap^epsilon+(l1yp+l2yp)^epsilon)^(eta/epsilon-1)*lap^(epsilon-1)=(((c1p^nu+(d2p*(c1p/d2)^(1/(nu-1)))^nu+2*d2p^nu)^(-1)*d2p^(nu-1))/(1+mu*(la^omega+lap^omega))*(alphaa*l2y^alphaa*k1p^(1-alphaa)+alphaa*l2yp^alphaa*k2p^(1-alphaa)+pay*(k1p*(1+mu*(la^omega+lap^omega))-(1-delta)*k1p)^2))*mu*omega*lap^(omega-1);
chi*(lap^epsilon+(l1yp+l2yp)^epsilon)^(eta/epsilon-1)*(l1yp+l2yp)^(epsilon-1)=((c1p^nu+(d2p*(c1p/d2)^(1/(nu-1)))^nu+2*d2p^nu)^(-1)*d2p^(nu-1))*(alphaa*l1yp^(alphaa-1)*k2^(1-alphaa))*c1p/d2;
chi*(lap^epsilon+(l1yp+l2yp)^epsilon)^(eta/epsilon-1)*(l1yp+l2yp)^(epsilon-1)=((c1p^nu+(d2p*(c1p/d2)^(1/(nu-1)))^nu+2*d2p^nu)^(-1)*d2p^(nu-1))*alphaa*l2yp^(alphaa-1)*k2p^(1-alphaa);
((1-alphaa)*l2y^alphaa*k1p^(-alphaa)+1-delta+2*(1-delta)*pay*(k1p*(1+mu*(la^omega+lap^omega))-(1-delta)*k1p))=(1+2*pay*(k1p*(1+mu*(la^omega+lap^omega))-(1-delta)*k1p))*(1+mu*(la^omega+lap^omega));
((1-alphaa)*l2yp^alphaa*k2p^(-alphaa)+1-delta)=(1+mu*(la^omega+lap^omega));
end;
 
initval;
c1=0.5243;
d2=14.6371;
la=16.7470;
l1y=1.8784;
l2y=14.8029;
k1=6.0658;
k2=20.9940;
c1p=6.3620;
d2p=4.6059;
lap=13.6277;
l1yp=17.8523;
l2yp=8.0043;
k1p=22.0770;
k2p=28.1082;
end;
 
steady (solve_algo=2, maxit=200000);
2

You should try to analytically compute the steady state. This does not simply involve plugging one equation into another one. You need to find a handle on how to sequentially compute the endogenous variables.

3

Thanks for the reply. I have tried doing so by writing a separate steady state solving file but so far no luck yet. I will try one more time. However, if the model is too complicated for each variable to be solved explicitly, is it correct to say that my only choice then would be to simplify the model’s assumption and reduce the model to an easier one?

4

I have not yet seen a model where the steady state computation cannot be broken down into one or two equations that can be solved numerically in the end. Take a look at the NK_baseline.mod in the Dynare example folder for one case dealing with a rather complicated model.

1 Like
5

I have tried but am still able to reduce the model to the point where a variable can be solved explicitly in terms of just the parameters as in the example provided by dynare. If I cannot express a variable in terms of just the parameters, does it mean there is no way to solve the model analytically? I have attached a pdf writeup of a model. If it is possible, could you give me some hints on how I could solve this analytically? Thank you very much and I really appreciate your help.
model.pdf (79.1 KB)

6

Have you tried calibrating labor and then solving the steady state. That makes the life often easier.

7

Do you mean setting all the labors as parameters over a range, then randomly select their values within that range, then solve for the other variables and finally check the square of residuals in the equations that haven’t been used? That is what I am doing as of now. However, it is not working yet and I feel like I’m just throwing darts randomly. Would there be another way? Thanks

8

No, I mean the often used approach of setting labor in steady state to 1/3 of the available time and then later on finding the labor disutility parameter that makes this value for labor the steady state value. In this case you can often rather easily solve for the other variables

2 Likes
9

In my model the amount of labor that agent can supply is not fixed. There are also 3 types of labor, i.e. R&D labor, production labor dedicated to product 1 and production labor dedicated to product 2. So do you mean that I should set the total quantity of labor fixed and then divide that by 3 to arrive at the calibration value for each type of labor?

10

Please read my statement carefully: it said you could try to fix labor to a value in steady state and then find the parameters that are compatible with this setting. This is often easier than fixing a parameter and then finding labor that is consistent with it. That often simplifies things, because with labor given in a first step, you can often directly compute capital from the Euler equation.

11

Dr. Pfeifer,

I have a question related to this problem. So in my code, I fix steady-state labor l_ss = 1/3 then solve for bunch of steady state value of variables using it. I know that I need to calibrate the disutility of labor parameter (for example, CHIW). But when I do this, i.e., setting CHIW = lambda_ssw_ss(1-l_ss) from the FOC of labor, the code won’t work and dynare will update my steady-state value of l, which is then different than what I initially set. The code works when I set both l_ss = 1/3 and CHIW = 0.5, for example. Can you tell me what I did wrong? Thank you.

12

Usually this implies you did make a mistake in your computations. If the parameter is correctly updated, the steady state should have the correct value for labor.

13

Hi Prof. Pfeifer,

and allow me @valerio88, thank you both of you!! I really learn a lot from you and the forum!!

Indeed, I’ve seen several papers and lecture notes, which compute SS by starting from N = \frac{1}{3} or 0.4. But my question is can I start from other variables?

My object is to do a Ramsey optimal policy exercise with some extra features. Before this section, the setup of the baseline is fairly common. In the Ramsey part, almost all SS need to be (re)built upon the instrument, so does the N. Now we look backward, in the baseline for the competitive equilibrium, it’s natural to attempt to connect N conditional on other variables as well, rather than direct assigning \frac{1}{3} at the very beginning.

For example, A_t = A_{t-1}^{\rho}\exp(\zeta^A_t), at the SS, A=1.

normalize Y_H=1, and mc =\displaystyle \frac{\epsilon-1}{\epsilon} is given.

K=\displaystyle\frac{\alpha Y_H}{r^k}mc, where r^k=\displaystyle\frac{1}{\beta}-(1-\delta) is known \Rightarrow K.

Y_H=AK^{\alpha}N^{1-\alpha}\Rightarrow N

w=\displaystyle\frac{(1-\alpha)Y_H}{N}mc\Rightarrow w

\chi N^{\phi} = \lambda w \Rightarrow \lambda, if, instead, I have set the coefficient of disutility of work \chi before.

\lambda = C^{-\sigma}\Rightarrow C.

This is the main difference between my code and peers.

However, my Matlab/Dynare reports that:
The steadystate file did not compute the steady state
resid shows only the production function has some resids, while model_diagnostics(M_,options_,oo_) says:
No obvious problems with this mod-file were detected

On the one hand, my calibration experience is very limited, prone to error.
On the other hand, before more “trial and error”, let me paraphrase my doubt:
Is it feasible to calibrate N from other variables?

So sorry for my verbosity.

1 Like
14

Yes. The reason you can fix N in steady state of many models is that this allows you to assign the value of the labor disutility parameter. That is a way of calibrating the model and is usually very easy, because you avoid a nonlinear equation in N. But in your context, it needs to be different. You need to fix a value for the labor disutility parameter \chi. However, that usually means that you need to compute N from a nonlinear equation you end up with after substituting for consumption as a function for labor in the budget constraint.

2 Likes
15

Thank you for instant reply! :handshake: :handshake: :handshake:

I worry my way to calibrate N is theoretically okay, but actually impracticable. Anyway, I plan to code the Ramsey case, and then move back, comparing these two cases.

16

I have checked. Both starting from N=\frac{1}{3} and the way that reconcilating with Ramsey are feasible. Results are identical. Thank you once again!

17

Ok, Thank you. But for future reference: the two ways are not always identical.