%TVCC_SIMU.M
%Simulates equilibrium dynamics in an economy with an aggregate time-varying collateral constraint. For a description of the underlying economic model, see ``Individual Versus Aggregate Collateral Constraints and the Overborrowing Syndrome,'' by Martin Uribe, April 2006.
%(c) Martin Uribe, April 27, 2006 
function [dm, A, Z, R, LA, H, Y, C, BETTA,Q, Ap, u] = tvcc_simu(zL, zH, DL, DH, aL, aH, a1L, a1H, n, RSTAR, SIG, OMEGA, BETTA1, ALFA, PAIZ, KAPA, Kss, T);

if nargin < 18
T = 5000; %length of simulation
end

rand('state',sum(100*clock)); %sets a different state each time.

Tcut = 1000;%number of periods thrown out in computing the test

%Usable periods
Tleft = T-Tcut;

A0 = 5.10;%initial debt position 

b=0;h=0;aa=0;%initial conditions for construction of den Haan and Marcet (RES, 1994) accuracy test

%initial state
x = rand<0.5; 

sL = x; %low state
sH = 1-x; %high state

for t=1:T

%External debt
A(t,1) = A0;

%Generate new state of nature
x = rand<PAIZ(1,1);

%Update state of nature
sL = sL * x + (1-sL) * (1-x);
sH = sH * x + (1-sH) * (1-x);

%Productivity shock
Z(t,1) = sL * zL + sH * zH;

%Chebyshev parameters
a = sL * aL + sH * aH;
a1 = sL * a1L + sH * a1H;

%Assume collateral constraint does not bind
R(t,1) = RSTAR;


X(t,1:n) = chebyreg2(A(t,1),DL,DH,n);

%Marginal utility of consumption
LA(t,1) = exp(X(t,:)*a) .* R(t,1);

%Price of land
Q(t,1) = exp(X(t,:)*a1);

%Hours
H(t,1) = ( (1-ALFA) * exp(Z(t)) * Kss^(ALFA) ) .^(1/(OMEGA+ALFA-1));

%Output
Y(t,1) = exp(Z(t)) * Kss^ALFA * H(t).^(1-ALFA);

%Consumption
C(t,1) = LA(t).^(-1/SIG) + H(t)^OMEGA/OMEGA;

%Next period's debt
Ap(t,1) = RSTAR .* (A(t) + C(t) - Y(t));

%Make sure next-period debt is within debt grid and satisfies collateral constraint
Ap(t,1) = max(Ap(t),DL);
Ap(t,1) = min([Ap(t) DH KAPA*Q(t)*Kss]);

%Consumption
C(t) = Ap(t,1) / RSTAR - A(t)  + Y(t);

%Marginal utility of consumption
LA(t,1) = (C(t) - H(t)^OMEGA/OMEGA)^(-SIG);

%Discount factor
BETTA(t,1) = (1 + C(t) - H(t)^OMEGA/OMEGA)^(-BETTA1);

%Interest rate
R(t,1) = LA(t) / exp(X(t,:)*a); 

%Update initial debt position
A0 = Ap(t);

%Prepare for computing denHaan-Marcet (RES, 1994) statistic 

%discard first Tcut periods
if t>Tcut;
tt = t-Tcut;

%Error in the bond Euler equaiton
u(tt,1) = LA(t-1)-R(t-1).* BETTA(t-1).* LA(t);

%Error in the price-of-land Euler equation
u1(tt,1) = LA(t-1)*Q(t-1) - BETTA(t-1) * LA(t) * (Q(t) + exp(Z(t)) * ALFA * (H(t)/Kss)^(1-ALFA));

%Instruments (unity and current and past values of debt and productivity shock)
h = [A(t-5:t-1)' 1 Z(t-5:t-1)'];

%Product of instruments and errors
UH = kron([u(tt) u1(tt)],h);

%mean of UH
b = b + UH/Tleft;

%Variance-covariance matrix of UH
aa = aa + UH'*UH/Tleft;
end %if t>Tcut
end %for t=1:T

%denHaan-Marcet (RES, 1994) statistic
dm = Tleft * b*inv(aa)*b';