cleaning diffuse algorithms
git-svn-id: https://www.dynare.org/svn/dynare/dynare_v4@1922 ac1d8469-bf42-47a9-8791-bf33cf982152time-shift
parent
bfe76b3c4b
commit
2c8a9d2f60
|
@ -151,27 +151,19 @@ function [fval,cost_flag,ys,trend_coeff,info] = DsgeLikelihood(xparam1,gend,data
|
|||
%------------------------------------------------------------------------------
|
||||
% 3. Initial condition of the Kalman filter
|
||||
%------------------------------------------------------------------------------
|
||||
kalman_algo = options_.kalman_algo;
|
||||
if options_.lik_init == 1 % Kalman filter
|
||||
Pstar = lyapunov_symm(T,R*Q*R',options_.qz_criterium);
|
||||
Pinf = [];
|
||||
kalman_algo = 1;
|
||||
Pstar = lyapunov_symm(T,R*Q*R',options_.qz_criterium);
|
||||
Pinf = [];
|
||||
elseif options_.lik_init == 2 % Old Diffuse Kalman filter
|
||||
Pstar = 10*eye(np);
|
||||
Pinf = [];
|
||||
kalman_algo = 1;
|
||||
Pstar = 10*eye(np);
|
||||
Pinf = [];
|
||||
elseif options_.lik_init == 3 % Diffuse Kalman filter
|
||||
if options_.kalman_algo < 4
|
||||
Pstar = zeros(np,np);
|
||||
ivs = bayestopt_.restrict_var_list_stationary;
|
||||
R1 = R(ivs,:);
|
||||
Pstar(ivs,ivs) = lyapunov_symm(T(ivs,ivs),R1*Q*R1',options_.qz_criterium);
|
||||
% Pinf = bayestopt_.Pinf;
|
||||
% by M. Ratto
|
||||
RR=T(:,bayestopt_.restrict_var_list_nonstationary);
|
||||
i=find(abs(RR)>1.e-10);
|
||||
R0=zeros(size(RR));
|
||||
R0(i)=sign(RR(i));
|
||||
Pinf=R0*R0';
|
||||
% by M. Ratto
|
||||
else
|
||||
if kalman_algo == 1
|
||||
kalman_algo == 3
|
||||
end
|
||||
[QT,ST] = schur(T);
|
||||
e1 = abs(ordeig(ST)) > 2-options_.qz_criterium;
|
||||
[QT,ST] = ordschur(QT,ST,e1);
|
||||
|
@ -183,104 +175,131 @@ function [fval,cost_flag,ys,trend_coeff,info] = DsgeLikelihood(xparam1,gend,data
|
|||
Pstar = zeros(np,np);
|
||||
B = QT'*R*Q*R'*QT;
|
||||
for i=np:-1:nk+2
|
||||
if ST(i,i-1) == 0
|
||||
if i == np
|
||||
c = zeros(np-nk,1);
|
||||
else
|
||||
c = ST(nk1:i,:)*(Pstar(:,i+1:end)*ST(i,i+1:end)')+...
|
||||
ST(i,i)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i);
|
||||
end
|
||||
q = eye(i-nk)-ST(nk1:i,nk1:i)*ST(i,i);
|
||||
Pstar(nk1:i,i) = q\(B(nk1:i,i)+c);
|
||||
Pstar(i,nk1:i-1) = Pstar(nk1:i-1,i)';
|
||||
else
|
||||
if i == np
|
||||
c = zeros(np-nk,1);
|
||||
c1 = zeros(np-nk,1);
|
||||
else
|
||||
c = ST(nk1:i,:)*(Pstar(:,i+1:end)*ST(i,i+1:end)')+...
|
||||
ST(i,i)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i)+...
|
||||
ST(i,i-1)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i-1);
|
||||
c1 = ST(nk1:i,:)*(Pstar(:,i+1:end)*ST(i-1,i+1:end)')+...
|
||||
ST(i-1,i-1)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i-1)+...
|
||||
ST(i-1,i)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i);
|
||||
end
|
||||
q = [eye(i-nk)-ST(nk1:i,nk1:i)*ST(i,i) -ST(nk1:i,nk1:i)*ST(i,i-1);...
|
||||
-ST(nk1:i,nk1:i)*ST(i-1,i) eye(i-nk)-ST(nk1:i,nk1:i)*ST(i-1,i-1)];
|
||||
z = q\[B(nk1:i,i)+c;B(nk1:i,i-1)+c1];
|
||||
Pstar(nk1:i,i) = z(1:(i-nk));
|
||||
Pstar(nk1:i,i-1) = z(i-nk+1:end);
|
||||
Pstar(i,nk1:i-1) = Pstar(nk1:i-1,i)';
|
||||
Pstar(i-1,nk1:i-2) = Pstar(nk1:i-2,i-1)';
|
||||
i = i - 1;
|
||||
end
|
||||
if ST(i,i-1) == 0
|
||||
if i == np
|
||||
c = zeros(np-nk,1);
|
||||
else
|
||||
c = ST(nk1:i,:)*(Pstar(:,i+1:end)*ST(i,i+1:end)')+...
|
||||
ST(i,i)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i);
|
||||
end
|
||||
q = eye(i-nk)-ST(nk1:i,nk1:i)*ST(i,i);
|
||||
Pstar(nk1:i,i) = q\(B(nk1:i,i)+c);
|
||||
Pstar(i,nk1:i-1) = Pstar(nk1:i-1,i)';
|
||||
else
|
||||
if i == np
|
||||
c = zeros(np-nk,1);
|
||||
c1 = zeros(np-nk,1);
|
||||
else
|
||||
c = ST(nk1:i,:)*(Pstar(:,i+1:end)*ST(i,i+1:end)')+...
|
||||
ST(i,i)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i)+...
|
||||
ST(i,i-1)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i-1);
|
||||
c1 = ST(nk1:i,:)*(Pstar(:,i+1:end)*ST(i-1,i+1:end)')+...
|
||||
ST(i-1,i-1)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i-1)+...
|
||||
ST(i-1,i)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i);
|
||||
end
|
||||
q = [eye(i-nk)-ST(nk1:i,nk1:i)*ST(i,i) -ST(nk1:i,nk1:i)*ST(i,i-1);...
|
||||
-ST(nk1:i,nk1:i)*ST(i-1,i) eye(i-nk)-ST(nk1:i,nk1:i)*ST(i-1,i-1)];
|
||||
z = q\[B(nk1:i,i)+c;B(nk1:i,i-1)+c1];
|
||||
Pstar(nk1:i,i) = z(1:(i-nk));
|
||||
Pstar(nk1:i,i-1) = z(i-nk+1:end);
|
||||
Pstar(i,nk1:i-1) = Pstar(nk1:i-1,i)';
|
||||
Pstar(i-1,nk1:i-2) = Pstar(nk1:i-2,i-1)';
|
||||
i = i - 1;
|
||||
end
|
||||
end
|
||||
if i == nk+2
|
||||
c = ST(nk+1,:)*(Pstar(:,nk+2:end)*ST(nk1,nk+2:end)')+ST(nk1,nk1)*ST(nk1,nk+2:end)*Pstar(nk+2:end,nk1);
|
||||
Pstar(nk1,nk1)=(B(nk1,nk1)+c)/(1-ST(nk1,nk1)*ST(nk1,nk1));
|
||||
c = ST(nk+1,:)*(Pstar(:,nk+2:end)*ST(nk1,nk+2:end)')+ST(nk1,nk1)*ST(nk1,nk+2:end)*Pstar(nk+2:end,nk1);
|
||||
Pstar(nk1,nk1)=(B(nk1,nk1)+c)/(1-ST(nk1,nk1)*ST(nk1,nk1));
|
||||
end
|
||||
|
||||
Z = QT(mf,:);
|
||||
R1 = QT'*R;
|
||||
end
|
||||
[u,s,v]=svd(Z*ST(:,1:nk),0);
|
||||
k = find(abs(diag(s)) < 1e-8);
|
||||
if length(k) > 0
|
||||
[junk,k1] = max(abs(v(:,k)));
|
||||
dd =ones(nk,1);
|
||||
dd(k1) = zeros(length(k1),1);
|
||||
Pinf(1:nk,1:nk) = diag(dd);
|
||||
end
|
||||
|
||||
end
|
||||
%------------------------------------------------------------------------------
|
||||
% 4. Likelihood evaluation
|
||||
%------------------------------------------------------------------------------
|
||||
if any(any(H ~= 0)) % should be replaced by a flag
|
||||
if options_.kalman_algo == 1
|
||||
LIK = DiffuseLikelihoodH1(T,R,Q,H,Pinf,Pstar,data,trend,start);
|
||||
if isinf(LIK) & ~estim_params_.ncn %% The univariate approach considered here doesn't
|
||||
%% apply when H has some off-diagonal elements.
|
||||
LIK = DiffuseLikelihoodH3(T,R,Q,H,Pinf,Pstar,data,trend,start);
|
||||
elseif isinf(LIK) & estim_params_.ncn
|
||||
LIK = DiffuseLikelihoodH3corr(T,R,Q,H,Pinf,Pstar,data,trend,start);
|
||||
end
|
||||
elseif options_.kalman_algo == 3
|
||||
if ~estim_params_.ncn %% The univariate approach considered here doesn't
|
||||
%% apply when H has some off-diagonal elements.
|
||||
LIK = DiffuseLikelihoodH3(T,R,Q,H,Pinf,Pstar,data,trend,start);
|
||||
else
|
||||
LIK = DiffuseLikelihoodH3corr(T,R,Q,H,Pinf,Pstar,data,trend,start);
|
||||
end
|
||||
end
|
||||
if kalman_algo == 1
|
||||
LIK = DiffuseLikelihoodH1(T,R,Q,H,Pinf,Pstar,data,trend,start);
|
||||
if isinf(LIK)
|
||||
kalman_algo = 2;
|
||||
if ~estim_params_.ncn
|
||||
%% The univariate approach considered here doesn't
|
||||
%% apply when H has some off-diagonal elements.
|
||||
LIK = DiffuseLikelihoodH3(T,R,Q,H,Pinf,Pstar,data,trend,start);
|
||||
else
|
||||
LIK = DiffuseLikelihoodH3corr(T,R,Q,H,Pinf,Pstar,data,trend,start);
|
||||
end
|
||||
end
|
||||
elseif kalman_algo == 2
|
||||
if ~estim_params_.ncn
|
||||
%% The univariate approach considered here doesn't
|
||||
%% apply when H has some off-diagonal elements.
|
||||
LIK = DiffuseLikelihoodH3(T,R,Q,H,Pinf,Pstar,data,trend,start);
|
||||
else
|
||||
LIK = DiffuseLikelihoodH3corr(T,R,Q,H,Pinf,Pstar,data,trend,start);
|
||||
end
|
||||
elseif kalman_algo == 3
|
||||
data1 = data - trend;
|
||||
LIK = DiffuseLikelihoodH1_Z(ST,Z,R1,Q,H,Pinf,Pstar,data1,trend,start);
|
||||
if isinf(LIK)
|
||||
kalman_algo = 4;
|
||||
if ~estim_params_.ncn
|
||||
%% The univariate approach considered here doesn't
|
||||
%% apply when H has some off-diagonal elements.
|
||||
LIK = DiffuseLikelihoodH3_Z(ST,Z,R1,Q,H,Pinf,Pstar,data1,trend,start);
|
||||
else
|
||||
LIK = DiffuseLikelihoodH3corr_Z(ST,Z,R1,Q,H,Pinf,Pstar,data1,trend,start);
|
||||
end
|
||||
end
|
||||
elseif kalman_algo == 4
|
||||
data1 = data - trend;
|
||||
if ~estim_params_.ncn
|
||||
%% The univariate approach considered here doesn't
|
||||
%% apply when H has some off-diagonal elements.
|
||||
LIK = DiffuseLikelihoodH3_Z(ST,Z,R1,Q,H,Pinf,Pstar,data1,trend,start);
|
||||
else
|
||||
LIK = DiffuseLikelihoodH3corr_Z(ST,Z,R1,Q,H,Pinf,Pstar,data1,trend,start);
|
||||
end
|
||||
end
|
||||
else
|
||||
if options_.kalman_algo == 1
|
||||
%nv = size(bayestopt_.Z,1);
|
||||
%LIK = kalman_filter(bayestopt_.Z,zeros(nv,nv),T,R,Q,data,zeros(size(T,1),1),Pstar,'u');
|
||||
%tic
|
||||
LIK = DiffuseLikelihood1(T,R,Q,Pinf,Pstar,data,trend,start);
|
||||
%toc
|
||||
%tic
|
||||
%LIK1 = Diffuse_Likelihood1(T,R,Q,Pinf,Pstar,data,trend,start);
|
||||
%toc
|
||||
%LIK1-LIK
|
||||
%if abs(LIK1-LIK)>0.0000000001
|
||||
% disp(['LIK1 and LIK are not equal! ' num2str(abs(LIK1-LIK))])
|
||||
%end
|
||||
if isinf(LIK)
|
||||
LIK = DiffuseLikelihood3(T,R,Q,Pinf,Pstar,data,trend,start);
|
||||
end
|
||||
elseif options_.kalman_algo == 3
|
||||
LIK = DiffuseLikelihood3(T,R,Q,Pinf,Pstar,data,trend,start);
|
||||
elseif options_.kalman_algo == 4
|
||||
data1 = data - trend;
|
||||
LIK = DiffuseLikelihood1_Z(ST,Z,R1,Q,Pinf,Pstar,data1,start);
|
||||
if isinf(LIK)
|
||||
LIK = DiffuseLikelihood3_Z(ST,Z,R1,Q,Pinf,Pstar,data1,start);
|
||||
end
|
||||
elseif options_.kalman_algo == 5
|
||||
data1 = data - trend;
|
||||
LIK = DiffuseLikelihood3_Z(ST,Z,R1,Q,Pinf,Pstar,data1,start);
|
||||
end
|
||||
if kalman_algo == 1
|
||||
LIK = DiffuseLikelihood1(T,R,Q,Pinf,Pstar,data,trend,start);
|
||||
if isinf(LIK)
|
||||
kalman_algo = 2
|
||||
LIK = DiffuseLikelihood3(T,R,Q,Pinf,Pstar,data,trend,start);
|
||||
end
|
||||
elseif kalman_algo == 2
|
||||
LIK = DiffuseLikelihood3(T,R,Q,Pinf,Pstar,data,trend,start);
|
||||
elseif options_.kalman_algo == 3
|
||||
data1 = data - trend;
|
||||
LIK = DiffuseLikelihood1_Z(ST,Z,R1,Q,Pinf,Pstar,data1,start);
|
||||
if isinf(LIK)
|
||||
kalman_algo = 4
|
||||
LIK = DiffuseLikelihood3_Z(ST,Z,R1,Q,Pinf,Pstar,data1,start);
|
||||
end
|
||||
elseif kalman_algo == 4
|
||||
data1 = data - trend;
|
||||
LIK = DiffuseLikelihood3_Z(ST,Z,R1,Q,Pinf,Pstar,data1,start);
|
||||
end
|
||||
end
|
||||
if imag(LIK) ~= 0
|
||||
likelihood = bayestopt_.penalty;
|
||||
likelihood = bayestopt_.penalty;
|
||||
else
|
||||
likelihood = LIK;
|
||||
likelihood = LIK;
|
||||
end
|
||||
% ------------------------------------------------------------------------------
|
||||
% Adds prior if necessary
|
||||
% ------------------------------------------------------------------------------
|
||||
lnprior = priordens(xparam1,bayestopt_.pshape,bayestopt_.p1,bayestopt_.p2,bayestopt_.p3,bayestopt_.p4);
|
||||
fval = (likelihood-lnprior);
|
||||
fval = (likelihood-lnprior);
|
||||
options_.kalman_algo = kalman_algo;
|
|
@ -92,27 +92,17 @@ function [alphahat,etahat,epsilonhat,ahat,SteadyState,trend_coeff,aK,T,R,P,PK,d,
|
|||
Q = M_.Sigma_e;
|
||||
H = M_.H;
|
||||
|
||||
kalman_algo = options_.kalman_algo;
|
||||
if options_.lik_init == 1 % Kalman filter
|
||||
Pstar = lyapunov_symm(T,R*Q*transpose(R),options_.qz_criterium);
|
||||
Pinf = [];
|
||||
kalman_algo = 1;
|
||||
Pstar = lyapunov_symm(T,R*Q*transpose(R),options_.qz_criterium);
|
||||
Pinf = [];
|
||||
elseif options_.lik_init == 2 % Old Diffuse Kalman filter
|
||||
Pstar = 10*eye(np);
|
||||
Pinf = [];
|
||||
kalman_algo = 1;
|
||||
Pstar = 10*eye(np);
|
||||
Pinf = [];
|
||||
elseif options_.lik_init == 3 % Diffuse Kalman filter
|
||||
if options_.kalman_algo < 4
|
||||
Pstar = zeros(np,np);
|
||||
ivs = bayestopt_.restrict_var_list_stationary;
|
||||
R1 = R(ivs,:);
|
||||
Pstar(ivs,ivs) = lyapunov_symm(T(ivs,ivs),R1*Q*R1',options_.qz_criterium);
|
||||
% Pinf = bayestopt_.Pinf;
|
||||
% by M. Ratto
|
||||
RR=T(:,bayestopt_.restrict_var_list_nonstationary);
|
||||
i=find(abs(RR)>1.e-10);
|
||||
R0=zeros(size(RR));
|
||||
R0(i)=sign(RR(i));
|
||||
Pinf=R0*R0';
|
||||
% by M. Ratto
|
||||
else
|
||||
kalman_algo = 3;
|
||||
[QT,ST] = schur(T);
|
||||
e1 = abs(ordeig(ST)) > 2-options_.qz_criterium;
|
||||
[QT,ST] = ordschur(QT,ST,e1);
|
||||
|
@ -124,99 +114,140 @@ function [alphahat,etahat,epsilonhat,ahat,SteadyState,trend_coeff,aK,T,R,P,PK,d,
|
|||
Pstar = zeros(np,np);
|
||||
B = QT'*R*Q*R'*QT;
|
||||
for i=np:-1:nk+2
|
||||
if ST(i,i-1) == 0
|
||||
if i == np
|
||||
c = zeros(np-nk,1);
|
||||
else
|
||||
c = ST(nk1:i,:)*(Pstar(:,i+1:end)*ST(i,i+1:end)')+...
|
||||
ST(i,i)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i);
|
||||
end
|
||||
q = eye(i-nk)-ST(nk1:i,nk1:i)*ST(i,i);
|
||||
Pstar(nk1:i,i) = q\(B(nk1:i,i)+c);
|
||||
Pstar(i,nk1:i-1) = Pstar(nk1:i-1,i)';
|
||||
else
|
||||
if i == np
|
||||
c = zeros(np-nk,1);
|
||||
c1 = zeros(np-nk,1);
|
||||
else
|
||||
c = ST(nk1:i,:)*(Pstar(:,i+1:end)*ST(i,i+1:end)')+...
|
||||
ST(i,i)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i)+...
|
||||
ST(i,i-1)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i-1);
|
||||
c1 = ST(nk1:i,:)*(Pstar(:,i+1:end)*ST(i-1,i+1:end)')+...
|
||||
ST(i-1,i-1)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i-1)+...
|
||||
ST(i-1,i)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i);
|
||||
end
|
||||
q = [eye(i-nk)-ST(nk1:i,nk1:i)*ST(i,i) -ST(nk1:i,nk1:i)*ST(i,i-1);...
|
||||
-ST(nk1:i,nk1:i)*ST(i-1,i) eye(i-nk)-ST(nk1:i,nk1:i)*ST(i-1,i-1)];
|
||||
z = q\[B(nk1:i,i)+c;B(nk1:i,i-1)+c1];
|
||||
Pstar(nk1:i,i) = z(1:(i-nk));
|
||||
Pstar(nk1:i,i-1) = z(i-nk+1:end);
|
||||
Pstar(i,nk1:i-1) = Pstar(nk1:i-1,i)';
|
||||
Pstar(i-1,nk1:i-2) = Pstar(nk1:i-2,i-1)';
|
||||
i = i - 1;
|
||||
end
|
||||
if ST(i,i-1) == 0
|
||||
if i == np
|
||||
c = zeros(np-nk,1);
|
||||
else
|
||||
c = ST(nk1:i,:)*(Pstar(:,i+1:end)*ST(i,i+1:end)')+...
|
||||
ST(i,i)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i);
|
||||
end
|
||||
q = eye(i-nk)-ST(nk1:i,nk1:i)*ST(i,i);
|
||||
Pstar(nk1:i,i) = q\(B(nk1:i,i)+c);
|
||||
Pstar(i,nk1:i-1) = Pstar(nk1:i-1,i)';
|
||||
else
|
||||
if i == np
|
||||
c = zeros(np-nk,1);
|
||||
c1 = zeros(np-nk,1);
|
||||
else
|
||||
c = ST(nk1:i,:)*(Pstar(:,i+1:end)*ST(i,i+1:end)')+...
|
||||
ST(i,i)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i)+...
|
||||
ST(i,i-1)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i-1);
|
||||
c1 = ST(nk1:i,:)*(Pstar(:,i+1:end)*ST(i-1,i+1:end)')+...
|
||||
ST(i-1,i-1)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i-1)+...
|
||||
ST(i-1,i)*ST(nk1:i,i+1:end)*Pstar(i+1:end,i);
|
||||
end
|
||||
q = [eye(i-nk)-ST(nk1:i,nk1:i)*ST(i,i) -ST(nk1:i,nk1:i)*ST(i,i-1);...
|
||||
-ST(nk1:i,nk1:i)*ST(i-1,i) eye(i-nk)-ST(nk1:i,nk1:i)*ST(i-1,i-1)];
|
||||
z = q\[B(nk1:i,i)+c;B(nk1:i,i-1)+c1];
|
||||
Pstar(nk1:i,i) = z(1:(i-nk));
|
||||
Pstar(nk1:i,i-1) = z(i-nk+1:end);
|
||||
Pstar(i,nk1:i-1) = Pstar(nk1:i-1,i)';
|
||||
Pstar(i-1,nk1:i-2) = Pstar(nk1:i-2,i-1)';
|
||||
i = i - 1;
|
||||
end
|
||||
end
|
||||
if i == nk+2
|
||||
c = ST(nk+1,:)*(Pstar(:,nk+2:end)*ST(nk1,nk+2:end)')+ST(nk1,nk1)*ST(nk1,nk+2:end)*Pstar(nk+2:end,nk1);
|
||||
Pstar(nk1,nk1)=(B(nk1,nk1)+c)/(1-ST(nk1,nk1)*ST(nk1,nk1));
|
||||
c = ST(nk+1,:)*(Pstar(:,nk+2:end)*ST(nk1,nk+2:end)')+ST(nk1,nk1)*ST(nk1,nk+2:end)*Pstar(nk+2:end,nk1);
|
||||
Pstar(nk1,nk1)=(B(nk1,nk1)+c)/(1-ST(nk1,nk1)*ST(nk1,nk1));
|
||||
end
|
||||
|
||||
Z = QT(mf,:);
|
||||
R1 = QT'*R;
|
||||
end
|
||||
end
|
||||
% -----------------------------------------------------------------------------
|
||||
% 4. Kalman smoother
|
||||
% -----------------------------------------------------------------------------
|
||||
if any(any(H ~= 0)) % should be replaced by a flag
|
||||
if options_.kalman_algo == 1
|
||||
[alphahat,epsilonhat,etahat,ahat,aK] = DiffuseKalmanSmootherH1(T,R,Q,H,Pinf,Pstar,Y,trend,nobs,np,smpl,mf);
|
||||
if all(alphahat(:)==0)
|
||||
[alphahat,epsilonhat,etahat,ahat,aK] = DiffuseKalmanSmootherH3(T,R,Q,H,Pinf,Pstar,Y,trend,nobs,np,smpl,mf);
|
||||
end
|
||||
elseif options_.kalman_algo == 3
|
||||
[alphahat,epsilonhat,etahat,ahat,aK] = DiffuseKalmanSmootherH3(T,R,Q,H,Pinf,Pstar,Y,trend,nobs,np,smpl,mf);
|
||||
end
|
||||
else
|
||||
if options_.kalman_algo == 1
|
||||
[alphahat,etahat,ahat,aK] = DiffuseKalmanSmoother1(T,R,Q,Pinf,Pstar,Y,trend,nobs,np,smpl,mf);
|
||||
if all(alphahat(:)==0)
|
||||
options_.kalman_algo = 3;
|
||||
[alphahat,etahat,ahat,aK] = DiffuseKalmanSmoother3(T,R,Q,Pinf,Pstar,Y,trend,nobs,np,smpl,mf);
|
||||
end
|
||||
elseif options_.kalman_algo == 3
|
||||
[alphahat,etahat,ahat,aK] = DiffuseKalmanSmoother3(T,R,Q,Pinf,Pstar,Y,trend,nobs,np,smpl,mf);
|
||||
elseif options_.kalman_algo == 4 | options_.kalman_algo == 5
|
||||
data1 = Y - trend;
|
||||
if options_.kalman_algo == 4
|
||||
[alphahat,etahat,ahat,P,aK,PK,d,decomp] = DiffuseKalmanSmoother1_Z(ST, ...
|
||||
Z,R1,Q,Pinf,Pstar,data1,nobs,np,smpl);
|
||||
if any(any(H ~= 0)) % should be replaced by a flag
|
||||
if kalman_algo == 1
|
||||
[alphahat,epsilonhat,etahat,ahat,aK] = ...
|
||||
DiffuseKalmanSmootherH1(T,R,Q,H,Pinf,Pstar,Y,trend,nobs,np,smpl,mf);
|
||||
if all(alphahat(:)==0)
|
||||
options_.kalman_algo = 5;
|
||||
[alphahat,etahat,ahat,P,aK,PK,d,decomp] = DiffuseKalmanSmoother3_Z(ST, ...
|
||||
kalman_algo = 2;
|
||||
if ~estim_params.ncn
|
||||
[alphahat,epsilonhat,etahat,ahat,aK] = ...
|
||||
DiffuseKalmanSmootherH3(T,R,Q,H,Pinf,Pstar,Y,trend,nobs,np,smpl,mf);
|
||||
else
|
||||
[alphahat,epsilonhat,etahat,ahat,aK] = ...
|
||||
DiffuseKalmanSmootherH3corr(T,R,Q,H,Pinf,Pstar,Y,trend, ...
|
||||
nobs,np,smpl,mf);
|
||||
end
|
||||
end
|
||||
elseif options_.kalman_algo == 2
|
||||
if ~estim_params.ncn
|
||||
[alphahat,epsilonhat,etahat,ahat,aK] = ...
|
||||
DiffuseKalmanSmootherH3(T,R,Q,H,Pinf,Pstar,Y,trend,nobs,np,smpl,mf);
|
||||
else
|
||||
[alphahat,epsilonhat,etahat,ahat,aK] = ...
|
||||
DiffuseKalmanSmootherH3corr(T,R,Q,H,Pinf,Pstar,Y,trend, ...
|
||||
nobs,np,smpl,mf);
|
||||
end
|
||||
elseif kalman_algo == 3
|
||||
data1 = Y - trend;
|
||||
[alphahat,epsilonhat,etahat,ahat,aK] = ...
|
||||
DiffuseKalmanSmootherH1_Z(ST,Z,R1,Q,H,Pinf,Pstar,data1,nobs,np,smpl);
|
||||
if all(alphahat(:)==0)
|
||||
kalman_algo = 4;
|
||||
if ~estim_params.ncn
|
||||
[alphahat,epsilonhat,etahat,ahat,aK] = ...
|
||||
DiffuseKalmanSmootherH3_Z(ST,Z,R1,Q,H,Pinf,Pstar,data1,nobs,np,smpl);
|
||||
else
|
||||
[alphahat,epsilonhat,etahat,ahat,aK] = ...
|
||||
DiffuseKalmanSmootherH3corr_Z(ST,Z,R1,Q,H,Pinf,Pstar,data1, ...
|
||||
nobs,np,smpl);
|
||||
end
|
||||
end
|
||||
elseif options_.kalman_algo == 4
|
||||
data1 = Y - trend;
|
||||
if ~estim_params.ncn
|
||||
[alphahat,epsilonhat,etahat,ahat,aK] = ...
|
||||
DiffuseKalmanSmootherH3_Z(ST,Z,R1,Q,H,Pinf,Pstar,data1, ...
|
||||
nobs,np,smpl);
|
||||
else
|
||||
[alphahat,epsilonhat,etahat,ahat,aK] = ...
|
||||
DiffuseKalmanSmootherH3corr_Z(ST,Z,R1,Q,H,Pinf,Pstar,data1, ...
|
||||
nobs,np,smpl);
|
||||
end
|
||||
end
|
||||
else
|
||||
if kalman_algo == 1
|
||||
[alphahat,etahat,ahat,aK] = DiffuseKalmanSmoother1(T,R,Q,Pinf,Pstar,Y,trend,nobs,np,smpl,mf);
|
||||
if all(alphahat(:)==0)
|
||||
kalman_algo = 2;
|
||||
[alphahat,etahat,ahat,aK] = DiffuseKalmanSmoother3(T,R,Q,Pinf,Pstar,Y,trend,nobs,np,smpl,mf);
|
||||
end
|
||||
elseif kalman_algo == 2
|
||||
[alphahat,etahat,ahat,aK] = DiffuseKalmanSmoother3(T,R,Q,Pinf,Pstar,Y,trend,nobs,np,smpl,mf);
|
||||
elseif kalman_algo == 3 | kalman_algo == 4
|
||||
data1 = Y - trend;
|
||||
if options_.kalman_algo == 3
|
||||
[alphahat,etahat,ahat,P,aK,PK,d,decomp] = DiffuseKalmanSmoother1_Z(ST, ...
|
||||
Z,R1,Q,Pinf,Pstar,data1,nobs,np,smpl);
|
||||
if all(alphahat(:)==0)
|
||||
options_.kalman_algo = 4;
|
||||
[alphahat,etahat,ahat,P,aK,PK,d,decomp] = DiffuseKalmanSmoother3_Z(ST, ...
|
||||
Z,R1,Q,Pinf,Pstar, ...
|
||||
data1,nobs,np,smpl);
|
||||
end
|
||||
else
|
||||
[alphahat,etahat,ahat,P,aK,PK,d,decomp] = DiffuseKalmanSmoother3_Z(ST, ...
|
||||
end
|
||||
else
|
||||
[alphahat,etahat,ahat,P,aK,PK,d,decomp] = DiffuseKalmanSmoother3_Z(ST, ...
|
||||
Z,R1,Q,Pinf,Pstar,data1,nobs,np,smpl);
|
||||
end
|
||||
alphahat = QT*alphahat;
|
||||
ahat = QT*ahat;
|
||||
if options_.nk > 0
|
||||
nk = options_.nk;
|
||||
for jnk=1:nk
|
||||
aK(jnk,:,:) = QT*squeeze(aK(jnk,:,:));
|
||||
for i=1:size(PK,4)
|
||||
PK(jnk,:,:,i) = QT*squeeze(PK(jnk,:,:,i))*QT';
|
||||
end
|
||||
for i=1:size(decomp,4)
|
||||
decomp(jnk,:,:,i) = QT*squeeze(decomp(jnk,:,:,i));
|
||||
end
|
||||
end
|
||||
for i=1:size(P,4)
|
||||
P(:,:,i) = QT*squeeze(P(:,:,i))*QT';
|
||||
end
|
||||
end
|
||||
end
|
||||
alphahat = QT*alphahat;
|
||||
ahat = QT*ahat;
|
||||
if options_.nk > 0
|
||||
nk = options_.nk;
|
||||
for jnk=1:nk
|
||||
aK(jnk,:,:) = QT*squeeze(aK(jnk,:,:));
|
||||
for i=1:size(PK,4)
|
||||
PK(jnk,:,:,i) = QT*squeeze(PK(jnk,:,:,i))*QT';
|
||||
end
|
||||
for i=1:size(decomp,4)
|
||||
decomp(jnk,:,:,i) = QT*squeeze(decomp(jnk,:,:,i));
|
||||
end
|
||||
end
|
||||
for i=1:size(P,4)
|
||||
P(:,:,i) = QT*squeeze(P(:,:,i))*QT';
|
||||
end
|
||||
end
|
||||
end
|
||||
end
|
||||
|
|
|
@ -37,6 +37,12 @@ if options_.order > 1
|
|||
options_.order = 1;
|
||||
end
|
||||
|
||||
if (~isempty(options_.unit_root_vars) || options_.diffuse_filter == 1)
|
||||
if options_.lik_init == 1
|
||||
options_.lik_init = 3
|
||||
end
|
||||
end
|
||||
|
||||
if options_.prefilter == 1
|
||||
options_.noconstant = 1;
|
||||
end
|
||||
|
|
|
@ -104,49 +104,50 @@ function global_initialization()
|
|||
options_.ramsey_policy = 0;
|
||||
|
||||
% estimation
|
||||
options_.load_mh_file = 0;
|
||||
options_.MaxNumberOfBytes = 1e6;
|
||||
options_.MaximumNumberOfMegaBytes = 111;
|
||||
options_.PosteriorSampleSize = 1000;
|
||||
options_.bayesian_irf = 0;
|
||||
options_.bayesian_th_moments = 0;
|
||||
options_.cutoff = 1e-12;
|
||||
options_.diffuse_d = [];
|
||||
options_.diffuse_filter = 0;
|
||||
options_.filter_step_ahead = [];
|
||||
options_.filtered_vars = 0;
|
||||
options_.first_obs = 1;
|
||||
options_.prefilter = 0;
|
||||
options_.presample = 0;
|
||||
options_.kalman_algo = 0;
|
||||
options_.kalman_tol = 1e-12;
|
||||
options_.lik_algo = 1;
|
||||
options_.lik_init = 1;
|
||||
options_.mh_replic = 20000;
|
||||
options_.load_mh_file = 0;
|
||||
options_.load_mh_file = 0;
|
||||
options_.logdata = 0;
|
||||
options_.loglinear = 0;
|
||||
options_.markowitz = 0.5;
|
||||
options_.mh_conf_sig = 0.90;
|
||||
options_.mh_drop = 0.5;
|
||||
options_.mh_jscale = 0.2;
|
||||
options_.mh_init_scale = 2*options_.mh_jscale;
|
||||
options_.mode_file = '';
|
||||
options_.mode_compute = 4;
|
||||
options_.mode_check = 0;
|
||||
options_.prior_trunc = 1e-10;
|
||||
options_.mh_conf_sig = 0.90;
|
||||
options_.mh_mode = 1;
|
||||
options_.mh_nblck = 2;
|
||||
options_.load_mh_file = 0;
|
||||
options_.mh_recover = 0;
|
||||
options_.nodiagnostic = 0;
|
||||
options_.loglinear = 0;
|
||||
options_.unit_root_vars = [];
|
||||
options_.bayesian_irf = 0;
|
||||
options_.bayesian_th_moments = 0;
|
||||
options_.smoother = 0;
|
||||
options_.mh_replic = 20000;
|
||||
options_.mode_check = 0;
|
||||
options_.mode_compute = 4;
|
||||
options_.mode_file = '';
|
||||
options_.moments_varendo = 0;
|
||||
options_.filtered_vars = 0;
|
||||
options_.kalman_algo = 1;
|
||||
options_.kalman_tol = 1e-12;
|
||||
options_.posterior_mode_estimation = 1;
|
||||
options_.MaxNumberOfBytes = 1e6;
|
||||
options_.filter_step_ahead = [];
|
||||
options_.diffuse_d = [];
|
||||
options_.logdata = 0;
|
||||
options_.use_mh_covariance_matrix = 0;
|
||||
options_.noconstant = 0;
|
||||
options_.markowitz = 0.5;
|
||||
options_.nodiagnostic = 0;
|
||||
options_.posterior_mode_estimation = 1;
|
||||
options_.prefilter = 0;
|
||||
options_.presample = 0;
|
||||
options_.prior_trunc = 1e-10;
|
||||
options_.simulation_method = 0;
|
||||
options_.cutoff = 1e-12;
|
||||
options_.smoother = 0;
|
||||
options_.student_degrees_of_freedom = 3;
|
||||
options_.subdraws = [];
|
||||
options_.PosteriorSampleSize = 1000;
|
||||
options_.MaximumNumberOfMegaBytes = 111;
|
||||
options_.unit_root_vars = [];
|
||||
options_.use_mh_covariance_matrix = 0;
|
||||
|
||||
% Misc
|
||||
options_.conf_sig = 0.6;
|
||||
|
|
Loading…
Reference in New Issue