Estimation of parameters

Note

Most of the code of this section is in the submodule Estimation

Settings

BASEforHANK.Parsing.EstimationSettingsType

EstimationSettings()

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.

source

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 in observed_vars_input. Missing data is denoted by NaN. 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 in growth_rate_select. Measurement errors will be added for observables in meas_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 (see BASEforHANK.Estimation.measurement_error()). For the latter case, the priors are set in meas_error_distr
  • numerical parameters: the maximum number of iterations to find the mode of the likelihood (see BASEforHANK.Estimation.mode_finding()) is set in max_iter_mode. ndraws, burnin and mhscale are parameters for the Random-Walk Metropolis Hastings algorithm (see BASEforHANK.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 (see BASEforHANK.Estimation.mode_finding()). multi_chain_init sets whether multiple chains in the RWMH (see BASEforHANK.Estimation.multi_chain_init() and BASEforHANK.rwmh()) are started from an overdispersed posterior mode. All flags are set to false by default.

Mode finding

BASEforHANK.Estimation.mode_findingFunction
mode_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 likelihood
  • hessian_final: Hessian of the log-likelihood at par_final
  • posterior_mode: log-likelihood at par_final
  • meas_error,meas_error_std: returns from measurement_error()
  • parnames: names of estimated parameters (including measurement error variances)
  • Data,Data_missing: data from e_set.data_file; marker for missing data
  • H_sel: selector matrix for states/controls that are observed
  • priors: priors of parameters (including measurement error variances)
  • smoother_output: output from the Kalman smoother
source

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:

  1. call BASEforHANK.PerturbationSolution.LinearSolution_estim() to derive the linear state-space representation of the model given par. Differently from BASEforHANK.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 derivatives A and B of the full model for aggregate variables and conditions, and compute the observation and state transition equations as in BASEforHANK.PerturbationSolution.LinearSolution()

  2. delete rows of the observation equation that correspond to unobserved controls (the selector matrix H_sel is constructed in BASEforHANK.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 (see BASEforHANK.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.likeliFunction
likeli(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 of prior and computed likelihood); alarm indicates error when solving model with LinearSolution_estim, sets log-likelihood to -9.e15

if smoother==True:

source
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 of prior and computed likelihood); alarm indicates error when solving model with LinearSolution_estim, sets log-likelihood to -9.e15

if smoother==True:

  • smoother_output: returns from kalman_filter_smoother()
  • State2Control,LOM: state-to-control and state transition matrizzes
source
BASEforHANK.PerturbationSolution.LinearSolution_estimFunction
LinearSolution_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 of Fsys() with respect to arguments X [B] and XPrime [A]
  • m_par: model parameters

Returns

as in LinearSolution()

source
BASEforHANK.Estimation.kalman_filterFunction
kalman_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 equation
  • LOM::Array{Float64,2}: law of motion for states
  • Data::Array{Union{Missing,Float64},2},D_miss::BitArray{2}: data (time $\times$ variable); marker for missing data
  • SCov::Array{Float64,2}: covariance of structural shocks
  • MCov::Array{Float64,2}: covariance of measurement error

Returns

  • log-likelihood
source
BASEforHANK.Estimation.measurement_errorFunction
measurement_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 observables
  • meas_error_prior: corresponding priors for measurement errors
  • meas_error_std: standard deviations of observables with measurement error
source

Bayesian estimation

BASEforHANK.sample_posteriorFunction
mcmc_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
source

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.rwmhFunction
rwmh(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 draws
  • accept_rate: acceptance rate
source
BASEforHANK.Estimation.multi_chain_initFunction
multi_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 chain
  • init_draw: Bool variable indicating whether search was succesful
source
BASEforHANK.Estimation.priorevalFunction
prioreval(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]
source
BASEforHANK.Estimation.kalman_filter_smootherFunction
kalman_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 equation
  • LOM::Array{Float64,2}: law of motion for states
  • Data::Array{Union{Missing,Float64},2},D_nomiss::BitArray{2}: data (time $\times$ variable); marker for existent data
  • SCov::Array{Float64,2}: covariance of structural shocks
  • MCov::Array{Float64,2}: covariance of measurement error

Returns

  • log_lik: log-likelihood
  • xhat_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]
source