Monte Carlo code

Discussion of State Space and Dynamic Stochastic General Equilibrium Models
Post Reply
jocildo_47
Posts: 7
Joined: Sun Sep 06, 2015 3:53 pm

Monte Carlo code

Post by jocildo_47 »

Please someone can give me one code to replicate the exercise Monte Carlo section 2.4 of Carvalho and Dam (2010) "Estimating the cross-sectional distribution of price stickiness from aggregate data, "National Bank of Denmark, WP 64
TomDoan
Posts: 7814
Joined: Wed Nov 01, 2006 4:36 pm

Re: Monte Carlo code

Post by TomDoan »

That's going to require using DSGE to produce a state-space model---once you have that in state-space model form, you can use DLM with TYPE=SIMULATE to create artificial data. With zeta=1, there's a simple closed form expression for the price variable given the (easily simulated) m and yn variables, but not with zeta<>1.

This sets up and does one simulation of the model. The Monte Carlo experiment requires estimating the model from the data (Y and P I assume) for each set of simulated data. Figuring out a reasonable way to do the estimation is probably the tricky part. We have examples of doing ML estimation with models that require DSGE solution as part of the function evaluation. However, it can often be hard to come up with good guess values when you have just two observables and a dozen deep parameters.

Code: Select all

*
* DSGE model from Carvalho and Dam(2010), "The Cross-Sectional Distribution of
* Price Stickiness Implied by Aggregate Data", FRB New York, Staff Report 419.
*
dec series y yn p m
dec series xk1 xk2 xk3 xk4
dec series pk1 pk2 pk3 pk4
*
dec real rho0 rho1 rho2
dec real delta0 delta1 delta2
dec real zeta beta
dec vect omega(4)
*
frml meq             = m  - (rho0+rho1*m{1}+rho2*m{2})
frml yneq            = yn - (delta0+delta1*yn{1}+delta2*yn{2})
frml(identity) mid   = m  - (p+y)
frml(identity) pdef  = p  - (omega(1)*pk1+omega(2)*pk2+omega(3)*pk3+omega(4)*pk4)
frml(identity) pkeq1 = pk1 - (1.0/1.0)*xk1
frml(identity) pkeq2 = pk2 - (1.0/2.0)*(xk2+xk2{1})
frml(identity) pkeq3 = pk3 - (1.0/3.0)*(xk3+xk3{1}+xk3{2})
frml(identity) pkeq4 = pk4 - (1.0/4.0)*(xk4+xk4{1}+xk4{2}+xk4{3})
*
* This defines a new variable for the summand on the right side of the xk's, to
* simplify writing the formulas.
*
frml(identity) zdef  = z   - (zeta*(m-yn)+(1-zeta)*p)
*
frml(identity) xkeq1 = xk1 - (1-beta)/(1-beta^1)*z
frml(identity) xkeq2 = xk2 - (1-beta)/(1-beta^2)*(z+beta*z{-1})
frml(identity) xkeq3 = xk3 - (1-beta)/(1-beta^2)*(z+beta*z{-1}+beta^2*z{-2})
frml(identity) xkeq4 = xk4 - (1-beta)/(1-beta^3)*(z+beta*z{-1}+beta^2*z{-2}+beta^3*z{-3})
*
compute zeta=.10,omega=||.40,.10,.10,.40||
compute rho0=0.0,rho1=1.43,rho2=-0.45,sigm=.005
compute delta0=0.0,delta1=.35,delta2=.15,sign=.05
compute beta=.99
*
group dsge meq yneq mid pdef zdef pkeq1 pkeq2 pkeq3 pkeq4 xkeq1 xkeq2 xkeq3 xkeq4
dsge(model=dsge,a=adsge,f=fdsge,z=zdsge) m yn y p z pk1 pk2 pk3 pk4 xk1 xk2 xk3 xk4
*
* Simulate 1000 data points of the model with the input parameters.
*
dlm(a=adsge,f=fdsge,z=zdsge,type=simulate,sw=%diag(||sigm^2,sign^2||),yhat=yhat,presample=ergodic) 1 1000 xstates
set psim 1 1000 = xstates(t)(4)
set ysim 1 1000 = xstates(t)(3)
jocildo_47
Posts: 7
Joined: Sun Sep 06, 2015 3:53 pm

Re: Monte Carlo code

Post by jocildo_47 »

Dear Tom, thank you for the code that you sent me. The model worked well, and helped me advance my research. I will return to you in due course with new doubts. My best regards. Jocildo
jocildo_47
Posts: 7
Joined: Sun Sep 06, 2015 3:53 pm

Re: Monte Carlo code

Post by jocildo_47 »

Dear Tom. Please, see, in the attached model, if there is any mistake in the form of imposing the restrictions on omegas, because, after the simulations, are resulting values that do not respect the restrictions when we use pmethod option, and zero variances of these parameters when we estimate without the pmethod option. Thank you very much. Jocildo

Code: Select all

dec series y yn p m z
dec series xk1 xk2 xk3 xk4
dec series pk1 pk2 pk3 pk4
*
dec real rho0 rho1 rho2
dec real delta0 delta1 delta2
dec real omega1 omega2 omega3 omega4
dec real zeta beta
*
frml meq             = m  - (rho0+rho1*m{1}+rho2*m{2})
frml yneq            = yn - (delta0+delta1*yn{1}+delta2*yn{2})
frml(identity) mid   = m  - (p+y)
frml(identity) pdef  = p  - (omega1*pk1+omega2*pk2+omega3*pk3+omega4*pk4)
frml(identity) pkeq1 = pk1 - (1.0/1.0)*xk1
frml(identity) pkeq2 = pk2 - (1.0/2.0)*(xk2+xk2{1})
frml(identity) pkeq3 = pk3 - (1.0/3.0)*(xk3+xk3{1}+xk3{2})
frml(identity) pkeq4 = pk4 - (1.0/4.0)*(xk4+xk4{1}+xk4{2}+xk4{3})
*
* This defines a new variable for the summand on the right side of the xk's, to
* simplify writing the formulas.
*
frml(identity) zdef  = z   - (zeta*(m-yn)+(1-zeta)*p)
*
frml(identity) xkeq1 = xk1 - (1-beta)/(1-beta^1)*z
frml(identity) xkeq2 = xk2 - (1-beta)/(1-beta^2)*(z+beta*z{-1})
frml(identity) xkeq3 = xk3 - (1-beta)/(1-beta^2)*(z+beta*z{-1}+beta^2*z{-2})
frml(identity) xkeq4 = xk4 - (1-beta)/(1-beta^3)*(z+beta*z{-1}+beta^2*z{-2}+beta^3*z{-3})
*
* Replicar 75 vezes
*
compute ndraws = 75
declare rect param(ndraws,13)
*
do i = 1, ndraws
*
compute zeta=.10,omega1=.40,omega2=.10,omega3=.10,omega4=.40
compute rho0=0.0,rho1=1.43,rho2=-0.45,sigm=.005
compute delta0=0.0,delta1=.35,delta2=.15,sign=.05
compute beta=.99
*
group dsge meq yneq mid pdef zdef pkeq1 pkeq2 pkeq3 pkeq4 xkeq1 xkeq2 xkeq3 xkeq4
dsge(model=dsge,a=adsge,f=fdsge,z=zdsge) m yn y p z pk1 pk2 pk3 pk4 xk1 xk2 xk3 xk4
*
* Simulate 1000 data points of the model with the input parameters.
*
dlm(a=adsge,f=fdsge,z=zdsge,type=simulate,sw=%diag(||sigm^2,sign^2||),yhat=yhat,presample=ergodic) 1 1000 xstates
set psim 1 1000 = xstates(t)(4)
set ysim 1 1000 = xstates(t)(3)
set ynsim 1 1000 = xstates(t)(2)
set msim 1 1000 = xstates(t)(1)
*
compute cdsge = %identity(2)~~%zeros(%rows(adsge)-2,2)
*
nonlin(parmset=base) zeta omega1 omega2 omega3 omega4 rho0 rho1 rho2 sigm sign delta0 delta1 delta2
nonlin(parmset=constraints) omega1>=0.0 omega2>=0.0 omega3>=0.0 omega4>=0.0 omega1+omega2+omega3+omega4==1.0
dlm(a=adsge,f=fdsge,z=zdsge,c=cdsge,y=||msim,ysim||,sw=%diag(||sigm^2,sign^2||),presample=ergodic,$
pmethod=simplex,piters=5,$
method=bfgs,parmset=base+constraints,noprint)
*
do j = 1, 13
compute param(i,j)= %beta(j)
end do j
*
end do i
*
dis param
TomDoan
Posts: 7814
Joined: Wed Nov 01, 2006 4:36 pm

Re: Monte Carlo code

Post by TomDoan »

SIMPLEX and GENETIC do not (and really cannot due to how they work) satisfy the >=, <= or == constraints. Those types of constraints are done with Lagrange multiplier methods as part of the BFGS algorithm. If you need more information, read the part about the algorithm in Section 4.4 of the User's Guide.
jocildo_47
Posts: 7
Joined: Sun Sep 06, 2015 3:53 pm

Re: Monte Carlo code

Post by jocildo_47 »

Please, Someone could identify why the attempt to resolve the model bellow, with dsge, is returning the following message: ## FO17. Mismatch between unstable roots and expectational terms.

dec series Y In C R K W L Z P
dec real sigma phi alpha beta delta rho

compute sigma = 2
compute phi = 1.5
compute alpha = 0.35
compute beta = 0.985
compute delta = 0.025
compute Zbar = 1.0
compute rho = 0.95

Frml(identity) f1 = C^sigma*L^phi - W/P
Frml(identity) f2 = beta*((1-delta)+R{-1}/P{-1}) - (C{-1}/C)^sigma
Frml(identity) f3 = (1-delta)*K + In - K{-1}
Frml(identity) f4 = Z*K^alpha*L^(1-alpha) - Y
Frml(identity) f5 = alpha^(y/(R/P)) - K
Frml(identity) f6 = (1-alpha)^(Y/(W/P)) - L
Frml(identity) f7 = ((1/Z)*(W/(1-alpha)^(1-alpha)*(R/alpha)^alpha)) - P
Frml(identity) f8 = C + In - Y
Frml F9 = (1-rho)*log(Zbar) + rho*log(Z{1}) - log(Z)

group model f1 f2 f3 f4 f5 f6 f7 f8 f9
dsge(model=model,expand=loglinear,a=adsge,z=zdsge,f=fdsge,solveby=svd,steadstate=ss,cutoff=1.001) Y In C R K W L Z P
dis adsge
TomDoan
Posts: 7814
Joined: Wed Nov 01, 2006 4:36 pm

Re: Monte Carlo code

Post by TomDoan »

Redate your capital stock equation. As written, that will be interpreted as an expected value of K(t+1). See:

https://estima.com/docs/RATS%2010%20Use ... f#page=363
Post Reply