Estimation of parameters
Most of the code of this section is in the submodule Estimation
Settings
BASEforHANK.Parsing.EstimationSettings
— TypeEstimationSettings()
Collect settings for the estimation of the model parameters in a struct
.
Use package Parameters
to provide initial values. Input and output file names are stored in the fields mode_start_file
, data_file
, save_mode_file
and save_posterior_file
.
The estimation settings (globally instantiated as e_set
in BASEforHANK.jl
) manage the following areas of the estimation:
- the match between data and model:
data_file
is a path to a .csv-file that contains quarterly observations of several variables (in columns), named inobserved_vars_input
. Missing data is denoted byNaN
. If some column-names do not align with the model variables,data_rename
is used. Level variables that should correspond to growth rates in the data can be selected ingrowth_rate_select
. Measurement errors will be added for observables inmeas_error_input
- estimation of variances:
shock_names
contain the aggregate shocks in the model, whose variances are estimated.me_treatment
defines how measurement errors are treated: for:fixed
, their variances are fixed by the data-variance, otherwise they are estimated either with:bounded
uniform priors, or:unbounded
priors (seeBASEforHANK.Estimation.measurement_error()
). For the latter case, the priors are set inmeas_error_distr
- numerical parameters: the maximum number of iterations to find the mode of the likelihood (see
BASEforHANK.Estimation.mode_finding()
) is set inmax_iter_mode
.ndraws
,burnin
andmhscale
are parameters for the Random-Walk Metropolis Hastings algorithm (seeBASEforHANK.rwmh()
) - estimation flags: whether to estimate the model is set in
estimate_model
.compute_hessian
determines whether the Hessian is computed after mode finding or set to an identity matrix (seeBASEforHANK.Estimation.mode_finding()
).multi_chain_init
sets whether multiple chains in the RWMH (seeBASEforHANK.Estimation.multi_chain_init()
andBASEforHANK.rwmh()
) are started from an overdispersed posterior mode. All flags are set tofalse
by default.
Mode finding
BASEforHANK.Estimation.mode_finding
— Functionmode_finding(XSSaggr, A, B, indexes, indexes_aggr, distrSS, compressionIndexes, m_par, n_par, e_set)
Given definition of observed variables and their transformation (level or growth rate) from e_set
, load the data, construct the observation equation, and maximize likeli()
(the log-likelihood) using the package Optim
.
Save estimation results to e_set.save_mode_file
.
Returns
par_final
: parameter vector that maximizes the likelihoodhessian_final
: Hessian of the log-likelihood atpar_final
posterior_mode
: log-likelihood atpar_final
meas_error
,meas_error_std
: returns frommeasurement_error()
parnames
: names of estimated parameters (including measurement error variances)Data
,Data_missing
: data frome_set.data_file
; marker for missing dataH_sel
: selector matrix for states/controls that are observedpriors
: priors of parameters (including measurement error variances)smoother_output
: output from the Kalman smoother
The main computations are the construction of the likelihood of the model parameters and its maximization. We get the model parameters that are to be estimated, together with their priors, from m_par
(in addition to measurement error variances, see Settings).
The likelihood function
The function BASEforHANK.Estimation.likeli()
computes the log-likelihood of the model parameters par
in the following steps:
call
BASEforHANK.PerturbationSolution.LinearSolution_estim()
to derive the linear state-space representation of the model givenpar
. Differently fromBASEforHANK.PerturbationSolution.LinearSolution()
, differentiate only the system of aggregate equilibrium conditions with respect to aggregate variables, i.e.BASEforHANK.PerturbationSolution.Fsys_agg()
. This is sufficient, as the estimated parameters do not enter in the heterogeneous agent part of the equilibrium system [BBL]. Then, update the derivativesA
andB
of the full model for aggregate variables and conditions, and compute the observation and state transition equations as inBASEforHANK.PerturbationSolution.LinearSolution()
delete rows of the observation equation that correspond to unobserved controls (the selector matrix
H_sel
is constructed inBASEforHANK.Estimation.mode_finding()
). Then, feed the linear state-space system, the data, and the variances of the structural and measurement shocks, into the Kalman filter (seeBASEforHANK.Estimation.kalman_filter()
), which computes the log-likelihood
We find the maximizer of the likelihood function, as well as its Hessian at the maximum, with the package Optim
. Note that in order to obtain the Hessian, you need to set e_set.compute_hessian = true
.
Called functions
BASEforHANK.Estimation.likeli
— Functionlikeli(par, Data, Data_missing, H_sel, priors, meas_error, meas_error_std, sr, lr, m_par, e_set; smoother=false)
Compute the likelihood of Data
, given model-parameters par
and prior priors
(maximize to find MLE of par
).
Solve model with LinearSolution_estim()
, compute likelihood with kalman_filter()
or with kalman_filter_smoother()
(if smoother==True
).
Returns
if smoother==False
:
log_like
,prior_like
,post_like
,alarm
: log-likelihoods (post
is the sum ofprior
and computed likelihood);alarm
indicates error when solving model withLinearSolution_estim
, sets log-likelihood to-9.e15
if smoother==True
:
smoother_output
: returns fromkalman_filter_smoother()
likeli(par, sr, lr, er, m_par, e_set; smoother=false)
Compute the likelihood of er.Data
, given model-parameters par
and prior er.priors
(maximize to find MLE of par
).
Solve model with LinearSolution_estim()
, compute likelihood with kalman_filter()
or with kalman_filter_smoother()
(if smoother==True
).
Returns
if smoother==False
:
log_like
,prior_like
,post_like
,alarm
: log-likelihoods (post
is the sum ofprior
and computed likelihood);alarm
indicates error when solving model withLinearSolution_estim
, sets log-likelihood to-9.e15
if smoother==True
:
smoother_output
: returns fromkalman_filter_smoother()
State2Control
,LOM
: state-to-control and state transition matrizzes
BASEforHANK.PerturbationSolution.LinearSolution_estim
— FunctionLinearSolution_estim(sr, m_par, A, B,;estim)
Calculate the linearized solution to the non-linear difference equations defined by function Fsys
, while only differentiating with respect to the aggregate part of the model, Fsys_agg()
.
The partials of the Jacobian belonging to the heterogeneous agent part of the model are taken from the full-model derivatives provided as arguments, A
and B
(computed by LinearSolution()
).
Arguments
sr
: steady-state structure (variable values, indexes, numerical parameters, ...)A
,B
: derivative ofFsys()
with respect to argumentsX
[B
] andXPrime
[A
]m_par
: model parameters
Returns
as in LinearSolution()
BASEforHANK.Estimation.kalman_filter
— Functionkalman_filter(H,LOM,Data,D_miss,SCov,MCov,e_set)
Compute likelihood of Data
, applying the Kalman filter to the state-space represenation (H
,LOM
) of the model.
Arguments
H::Array{Float64,2}
: observation equationLOM::Array{Float64,2}
: law of motion for statesData::Array{Union{Missing,Float64},2}
,D_miss::BitArray{2}
: data (time $\times$ variable); marker for missing dataSCov::Array{Float64,2}
: covariance of structural shocksMCov::Array{Float64,2}
: covariance of measurement error
Returns
- log-likelihood
BASEforHANK.Estimation.measurement_error
— Functionmeasurement_error(Data,observed_vars,e_set)
Build measurement error.
Arguments
Data
: matrix of observables [nobs * nvar]observed_vars
: vector of observed variable names [nvar * 1]e_set::EstimationSettings
Returns
meas_error
: ordered dictionary of measurement errors linked to observablesmeas_error_prior
: corresponding priors for measurement errorsmeas_error_std
: standard deviations of observables with measurement error
Bayesian estimation
BASEforHANK.sample_posterior
— Functionmcmc_estimation(mr,er;file=e_set.save_posterior_file)
Sample posterior of parameter vector with rwmh()
, take sample mean as parameter estimate, and save all results in file
.
Arguments
sr::SteadyResults
mr::LinearResults
er::EstimResults
We use a Monte Carlo Markov Chain method, specifically the Random-Walk Metropolis Hastings (BASEforHANK.Estimation.rwmh()
) algorithm, to sample from the posterior probability distribution of the parameter vector. The acceptance rate of the algorithm can be adjusted via setting EstimationSettings.mhscale
. To obtain the posterior likelihood of each draw, we call BASEforHANK.Estimation.likeli()
, which evaluates the priors at par
(BASEforHANK.Estimation.prioreval()
) and returns the log-posterior as a sum of the log-prior and the log-likelihood.
Given the draws from the posterior, we can analyze the probabilities of the parameters using the package MCMCChains
. We take the average over the draws as our Bayesian estimate of the parameter vector, par_final
. To obtain an estimate of the underlying state over the data sample period, we call BASEforHANK.Estimation.likeli()
with par_final
and keyword smoother=true
(this calls the Kalman smoother BASEforHANK.Estimation.kalman_filter_smoother()
). The result is stored in smoother_output
, and saved with the other results in e_set.save_posterior_file
.
Called functions
BASEforHANK.Estimation.rwmh
— Functionrwmh(xhat, Σ, sr, lr, er, m_par, e_set)
Sample the posterior of the parameter vector using the Random-Walk Metropolis Hastings algorithm.
Returns
draws::Array{Float64,2}
:e_set.ndraws + e_set.burnin
sampled parameter vectors (row vectors)posterior
: vector of posteriors for the respective drawsaccept_rate
: acceptance rate
BASEforHANK.Estimation.multi_chain_init
— Functionmulti_chain_init(xhat, Σ, sr, lr, er, m_par, e_set)
Draw overdispersed initial values for multi-chain RWMH.
Returns
init_draw
: overdispersed starting value for chaininit_draw
: Bool variable indicating whether search was succesful
BASEforHANK.Estimation.prioreval
— Functionprioreval(par,priors)
Evaluate prior PDF at the parameters given in par
.
Arguments
par
: vector of parameters [npar*1]priors
: vector of prior distributions [npar*1]
Returns
log_priorval
: log prior density [scalar]alarm
: indicator that is 1 if there is a violation of the prior bounds [scalar]
BASEforHANK.Estimation.kalman_filter_smoother
— Functionkalman_filter_smoother(H, LOM, Data, D_nomiss, SCov, MCov, e_set)
Compute likelihood and estimate of underlying states given the full observed Data
by applying the Kalman smoother to the state-space representation (H
,LOM
) of the model.
Arguments
H::Array{Float64,2}
: observation equationLOM::Array{Float64,2}
: law of motion for statesData::Array{Union{Missing,Float64},2}
,D_nomiss::BitArray{2}
: data (time $\times$ variable); marker for existent dataSCov::Array{Float64,2}
: covariance of structural shocksMCov::Array{Float64,2}
: covariance of measurement error
Returns
log_lik
: log-likelihoodxhat_tgt
,xhat_tgT
: estimate of underlying states from forward iteration [xhat_tgt
] and backward iteration [xhat_tgT
]Sigma_tgt
,Sigma_tgT
: estimate of covariance matrices from forward iteration [Sigma_tgt
] and backward iteration [Sigma_tgT
]s
,m
: smoothed structural shocks [s
] and measurement errors [m
]
- BBLFor details, see the paper Shocks, Frictions, and Inequality in US Business Cycles, American Economic Review, forthcoming.