Correct the auxiliary filter part.
parent
3c18bf1f6b
commit
df842aa542
|
@ -43,7 +43,7 @@ persistent start_param sample_size number_of_observed_variables number_of_struct
|
||||||
% Set seed for randn().
|
% Set seed for randn().
|
||||||
set_dynare_seed('default') ;
|
set_dynare_seed('default') ;
|
||||||
pruning = DynareOptions.particle.pruning;
|
pruning = DynareOptions.particle.pruning;
|
||||||
second_resample = 1 ;
|
second_resample = 0 ;
|
||||||
variance_update = 1 ;
|
variance_update = 1 ;
|
||||||
|
|
||||||
% initialization of state particles
|
% initialization of state particles
|
||||||
|
@ -122,7 +122,7 @@ for t=1:sample_size
|
||||||
chol_sigma_bar = chol(h_square*sigma_bar)' ;
|
chol_sigma_bar = chol(h_square*sigma_bar)' ;
|
||||||
end
|
end
|
||||||
% Prediction (without shocks)
|
% Prediction (without shocks)
|
||||||
wtilde = zeros(1,number_of_particles) ;
|
tau_tilde = zeros(1,number_of_particles) ;
|
||||||
for i=1:number_of_particles
|
for i=1:number_of_particles
|
||||||
% model resolution
|
% model resolution
|
||||||
[ys,trend_coeff,exit_flag,info,Model,DynareOptions,BayesInfo,DynareResults,ReducedForm] = ...
|
[ys,trend_coeff,exit_flag,info,Model,DynareOptions,BayesInfo,DynareResults,ReducedForm] = ...
|
||||||
|
@ -145,22 +145,23 @@ for t=1:sample_size
|
||||||
tmp = local_state_space_iteration_2(yhat,zeros(number_of_structural_innovations,1),ghx,ghu,constant,ghxx,ghuu,ghxu,DynareOptions.threads.local_state_space_iteration_2);
|
tmp = local_state_space_iteration_2(yhat,zeros(number_of_structural_innovations,1),ghx,ghu,constant,ghxx,ghuu,ghxu,DynareOptions.threads.local_state_space_iteration_2);
|
||||||
end
|
end
|
||||||
PredictionError = bsxfun(@minus,Y(t,:)',tmp(mf1,:));
|
PredictionError = bsxfun(@minus,Y(t,:)',tmp(mf1,:));
|
||||||
wtilde(i) = exp(-.5*(const_lik+log(det(ReducedForm.H))+sum(PredictionError.*(ReducedForm.H\PredictionError),1))) ;
|
% Replace Gaussian density with a Student density with 3 degrees of
|
||||||
|
% freedom for fat tails.
|
||||||
|
z = sum(PredictionError.*(ReducedForm.H\PredictionError),1) ;
|
||||||
|
tau_tilde(i) = weights(i).*(tpdf(z,3*ones(size(z)))+1e-99) ;
|
||||||
|
%tau_tilde(i) = weights(i).*exp(-.5*(const_lik+log(det(ReducedForm.H))+sum(PredictionError.*(ReducedForm.H\PredictionError),1))) ;
|
||||||
end
|
end
|
||||||
% unormalized weights and observation likelihood contribution
|
|
||||||
tau_tilde = weights.*wtilde ;
|
|
||||||
sum_tau_tilde = sum(tau_tilde) ;
|
|
||||||
% particles selection
|
% particles selection
|
||||||
tau_tilde = tau_tilde/sum_tau_tilde ;
|
tau_tilde = tau_tilde/sum(tau_tilde) ;
|
||||||
indx = resample(0,tau_tilde',DynareOptions.particle);
|
indx = resample(0,tau_tilde',DynareOptions.particle);
|
||||||
StateVectors = StateVectors(:,indx) ;
|
StateVectors = StateVectors(:,indx) ;
|
||||||
if pruning
|
if pruning
|
||||||
StateVectors_ = StateVectors_(:,indx) ;
|
StateVectors_ = StateVectors_(:,indx) ;
|
||||||
end
|
end
|
||||||
xparam = bsxfun(@plus,(1-small_a).*m_bar,small_a.*xparam(:,indx)) ;
|
xparam = bsxfun(@plus,(1-small_a).*m_bar,small_a.*xparam(:,indx)) ;
|
||||||
wtilde = wtilde(indx) ;
|
w_stage1 = weights(indx)./tau_tilde(indx) ;
|
||||||
% draw in the new distributions
|
% draw in the new distributions
|
||||||
lnw = zeros(1,number_of_particles) ;
|
wtilde = zeros(1,number_of_particles) ;
|
||||||
i = 1 ;
|
i = 1 ;
|
||||||
while i<=number_of_particles
|
while i<=number_of_particles
|
||||||
candidate = xparam(:,i) + chol_sigma_bar*randn(number_of_parameters,1) ;
|
candidate = xparam(:,i) + chol_sigma_bar*randn(number_of_parameters,1) ;
|
||||||
|
@ -191,18 +192,16 @@ for t=1:sample_size
|
||||||
end
|
end
|
||||||
StateVectors(:,i) = tmp(mf0,:) ;
|
StateVectors(:,i) = tmp(mf0,:) ;
|
||||||
PredictionError = bsxfun(@minus,Y(t,:)',tmp(mf1,:));
|
PredictionError = bsxfun(@minus,Y(t,:)',tmp(mf1,:));
|
||||||
lnw(i) = exp(-.5*(const_lik+log(det(ReducedForm.H))+sum(PredictionError.*(ReducedForm.H\PredictionError),1)));
|
wtilde(i) = w_stage1(i)*exp(-.5*(const_lik+log(det(ReducedForm.H))+sum(PredictionError.*(ReducedForm.H\PredictionError),1)));
|
||||||
i = i+1 ;
|
i = i+1 ;
|
||||||
end
|
end
|
||||||
end
|
end
|
||||||
% importance ratio
|
|
||||||
wtilde = lnw./wtilde ;
|
|
||||||
% normalization
|
% normalization
|
||||||
weights = wtilde/sum(wtilde);
|
weights = wtilde/sum(wtilde);
|
||||||
if (variance_update==1) && (neff(weights)<DynareOptions.particle.resampling.threshold*sample_size)
|
if (variance_update==1) && (neff(weights)<DynareOptions.particle.resampling.threshold*sample_size)
|
||||||
variance_update = 0 ;
|
variance_update = 0 ;
|
||||||
end
|
end
|
||||||
% final resampling (advised)
|
% final resampling (not advised)
|
||||||
if second_resample==1
|
if second_resample==1
|
||||||
indx = resample(0,weights,DynareOptions.particle);
|
indx = resample(0,weights,DynareOptions.particle);
|
||||||
StateVectors = StateVectors(:,indx) ;
|
StateVectors = StateVectors(:,indx) ;
|
||||||
|
@ -234,16 +233,16 @@ for t=1:sample_size
|
||||||
pass2=1 ;
|
pass2=1 ;
|
||||||
pass3=1 ;
|
pass3=1 ;
|
||||||
for j=1:number_of_particles
|
for j=1:number_of_particles
|
||||||
if cumulated_weights(j)>0.025 && pass1==1
|
if cumulated_weights(j)>=0.025 && pass1==1
|
||||||
lb95_xparam(i,t) = (temp(j-1,1)+temp(j,1))/2 ;
|
lb95_xparam(i,t) = temp(j,1) ;
|
||||||
pass1 = 2 ;
|
pass1 = 2 ;
|
||||||
end
|
end
|
||||||
if cumulated_weights(j)>0.5 && pass2==1
|
if cumulated_weights(j)>=0.5 && pass2==1
|
||||||
median_xparam(i,t) = (temp(j-1,1)+temp(j,1))/2 ;
|
median_xparam(i,t) = temp(j,1) ;
|
||||||
pass2 = 2 ;
|
pass2 = 2 ;
|
||||||
end
|
end
|
||||||
if cumulated_weights(j)>0.975 && pass3==1
|
if cumulated_weights(j)>=0.975 && pass3==1
|
||||||
ub95_xparam(i,t) = (temp(j-1,1)+temp(j,1))/2 ;
|
ub95_xparam(i,t) = temp(j,1) ;
|
||||||
pass3 = 2 ;
|
pass3 = 2 ;
|
||||||
end
|
end
|
||||||
end
|
end
|
||||||
|
@ -368,4 +367,4 @@ for plt = 1:nbplt,
|
||||||
fprintf(fidTeX,' \n');
|
fprintf(fidTeX,' \n');
|
||||||
end
|
end
|
||||||
end
|
end
|
||||||
|
|
||||||
|
|
Loading…
Reference in New Issue