diff --git a/src/conditional_filter_proposal.m b/src/conditional_filter_proposal.m index aa82a449a..0ef5aedab 100644 --- a/src/conditional_filter_proposal.m +++ b/src/conditional_filter_proposal.m @@ -1,4 +1,4 @@ -function [ProposalStateVector,Weights] = conditional_filter_proposal(ReducedForm,obs,StateVectors,SampleWeights,Q_lower_triangular_cholesky,H_lower_triangular_cholesky,H,DynareOptions,normconst2) +function [ProposalStateVector,Weights] = conditional_filter_proposal(ReducedForm,obs,StateVectors,SampleWeights,Q_lower_triangular_cholesky,H_lower_triangular_cholesky,H,ParticleOptions,ThreadsOptions,normconst2) % % Computes the proposal for each past particle using Gaussian approximations % for the state errors and the Kalman filter @@ -73,11 +73,11 @@ if isempty(init_flag2) init_flag2 = 1; end -if DynareOptions.particle.proposal_approximation.cubature || DynareOptions.particle.proposal_approximation.montecarlo +if ParticleOptions.proposal_approximation.cubature || ParticleOptions.proposal_approximation.montecarlo [nodes,weights] = spherical_radial_sigma_points(number_of_structural_innovations); weights_c = weights ; -elseif DynareOptions.particle.proposal_approximation.unscented - [nodes,weights,weights_c] = unscented_sigma_points(number_of_structural_innovations,DynareOptions); +elseif ParticleOptions.proposal_approximation.unscented + [nodes,weights,weights_c] = unscented_sigma_points(number_of_structural_innovations,ParticleOptions); else error('Estimation: This approximation for the proposal is not implemented or unknown!') end @@ -85,12 +85,12 @@ end epsilon = Q_lower_triangular_cholesky*(nodes') ; yhat = repmat(StateVectors-state_variables_steady_state,1,size(epsilon,2)) ; -tmp = local_state_space_iteration_2(yhat,epsilon,ghx,ghu,constant,ghxx,ghuu,ghxu,DynareOptions.threads.local_state_space_iteration_2); +tmp = local_state_space_iteration_2(yhat,epsilon,ghx,ghu,constant,ghxx,ghuu,ghxu,ThreadsOptions.local_state_space_iteration_2); PredictedStateMean = tmp(mf0,:)*weights ; PredictedObservedMean = tmp(mf1,:)*weights; -if DynareOptions.particle.proposal_approximation.cubature || DynareOptions.particle.proposal_approximation.montecarlo +if ParticleOptions.proposal_approximation.cubature || ParticleOptions.proposal_approximation.montecarlo PredictedStateMean = sum(PredictedStateMean,2) ; PredictedObservedMean = sum(PredictedObservedMean,2) ; dState = bsxfun(@minus,tmp(mf0,:),PredictedStateMean)'.*sqrt(weights) ; @@ -117,7 +117,7 @@ else end ProposalStateVector = StateVectorVarianceSquareRoot*randn(size(StateVectorVarianceSquareRoot,2),1)+StateVectorMean ; -ypred = measurement_equations(ProposalStateVector,ReducedForm,DynareOptions) ; +ypred = measurement_equations(ProposalStateVector,ReducedForm,ThreadsOptions) ; foo = H_lower_triangular_cholesky \ (obs - ypred) ; likelihood = exp(-0.5*(foo)'*foo)/normconst2 + 1e-99 ; Weights = SampleWeights.*likelihood; diff --git a/src/conditional_particle_filter.m b/src/conditional_particle_filter.m index 3f365f85c..4aef9f392 100644 --- a/src/conditional_particle_filter.m +++ b/src/conditional_particle_filter.m @@ -1,4 +1,4 @@ -function [LIK,lik] = conditional_particle_filter(ReducedForm,Y,start,DynareOptions) +function [LIK,lik] = conditional_particle_filter(ReducedForm,Y,start,ParticleOptions,ThreadsOptions) % % Evaluates the likelihood of a non-linear model with a particle filter % - the proposal is built using the Kalman updating step for each particle. @@ -74,7 +74,7 @@ if isempty(init_flag) number_of_state_variables = length(mf0); number_of_observed_variables = length(mf1); init_flag = 1; - number_of_particles = DynareOptions.particle.number_of_particles ; + number_of_particles = ParticleOptions.number_of_particles ; end % Get covariance matrices @@ -108,14 +108,14 @@ SampleWeights = ones(1,number_of_particles)/number_of_particles ; for t=1:sample_size for i=1:number_of_particles [StateParticles(:,i),SampleWeights(i)] = ... - conditional_filter_proposal(ReducedForm,Y(:,t),StateParticles(:,i),SampleWeights(i),Q_lower_triangular_cholesky,H_lower_triangular_cholesky,H,DynareOptions,normconst2) ; + conditional_filter_proposal(ReducedForm,Y(:,t),StateParticles(:,i),SampleWeights(i),Q_lower_triangular_cholesky,H_lower_triangular_cholesky,H,ParticleOptions,ThreadsOptions,normconst2) ; end SumSampleWeights = sum(SampleWeights) ; lik(t) = log(SumSampleWeights) ; SampleWeights = SampleWeights./SumSampleWeights ; - if (DynareOptions.particle.resampling.status.generic && neff(SampleWeights). -lambda = (DynareOptions.particle.unscented.alpha^2)*(n+DynareOptions.particle.unscented.kappa) - n ; +lambda = (ParticleOptions.unscented.alpha^2)*(n+ParticleOptions.unscented.kappa) - n ; nodes = [ zeros(n,1) ( sqrt(n+lambda).*([ eye(n) -eye(n)]) ) ]' ; W_m = lambda/(n+lambda) ; -W_c = W_m + (1-DynareOptions.particle.unscented.alpha^2+DynareOptions.particle.unscented.beta) ; +W_c = W_m + (1-ParticleOptions.unscented.alpha^2+ParticleOptions.unscented.beta) ; temp = ones(2*n,1)/(2*(n+lambda)) ; W_m = [W_m ; temp] ; W_c = [W_c ; temp] ;