function [Y,X] = simu_2nd(gx, hx, gxx, hxx, gss, hss, eta, sig, x0, e) % %[Y,X] = simu_2nd(gx, hx, gxx, hxx, gss, hss, eta, sig, x0, e) simulates times series from the model: % %xp^i = hx^i_a x_a + 1/2 hxx^i_ab x_a x_b + 1/2 hss^i sig^2 + sig * eta^i_c ep_c % %y^i = gx^i_a x_a + 1/2 gxx^i_ab x_a x_b + 1/2 gss^i sig^2 % %where % hx is nx by nx % gx is ny by nx % hxx is nx by nx by nx % gxx is ny by nx by nx % hss is nx by 1 % gss is ny by 1 % eta is nx by ne % sig is a parameter scaling the innovation to the exogenous variables % e is T by ne random shock % T is the length of the simulation % % x0 is the initial condition for X % %Output: vectors X (T by nx) and Y (T by ny) containing time series for x and y and e (T by 1) containing % % (c) Stephanie Schmitt-Grohe and Martin Uribe, January 22, 2002 %Number of periods T=size(e,1)+1; %Number of shocks ne = size(eta,2); %Columns of X and Y nx=size(hx,1); ny=size(gx,1); x0=x0(:); X(1,1:nx)=x0'; %Decompose x into its first- and second-order components xf=x0; xs=0*x0; for t=1:T %Construct y(t) for i=1:ny y(i,1) = gx(i,:) * (xf+xs) + 1/2 * xf' * squeeze(gxx(i,:,:)) * xf + 1/2 * gss(i,1)*sig^2; end %i Y(t,1:ny) = y'; %Remember xf(t) and xs(t) when entering t+1 xoldf = xf; xolds = xs; %Construct x(t+1) if t