dsge: Dynamic Stochastic General Equilibrium Models

Description

Specify, solve, and estimate dynamic stochastic general equilibrium (DSGE) models by maximum likelihood and Bayesian methods. Supports both linear models via an equation-based formula interface and nonlinear models via string-based equations with first-order perturbation (linearization around deterministic steady state). Solution uses the method of undetermined coefficients (Klein, 2000 doi:10.1016/S0165-1889(99)00045-7). Likelihood evaluated via the Kalman filter. Bayesian estimation uses adaptive Random-Walk Metropolis-Hastings with prior specification. Additional tools include Kalman smoothing, historical shock decomposition, local identification diagnostics, parameter sensitivity analysis, second-order perturbation, occasionally binding constraints, impulse-response functions, forecasting, and robust standard errors.

Author(s)

Maintainer: Mustapha Wasseja Mohammed muswaseja@gmail.com

Compute Equation Hessians via Central Differences

Description

Compute Equation Hessians via Central Differences

Usage

.compute_equation_hessians(fn, x0, n_eq, n_vars, eps = 1e-05)

Forward simulate a linear state-space system

Description

Forward simulate a linear state-space system

Usage

.forward_simulate(H, G, M, x0, shock_path, horizon)

Parse a constraint from a string

Description

Parse a constraint from a string

Usage

.parse_constraint_string(s)

Expectation Operator (Alias for lead)

Description

A user-friendly alias for lead(x, 1). Represents the one-period-ahead model-consistent expectation of variable x.

Usage

E(x)

Arguments

Details

E(x) is equivalent to lead(x, 1). The parser translates E(x) to lead(x, 1) internally.

Value

This function is not meant to be called directly; it always throws an error. It is recognized as a syntactic marker by the equation parser inside dsge_model().

See Also

lead(), dsge_model()

Estimate a DSGE Model by Bayesian Methods

Description

Estimates the parameters of a DSGE model using Random-Walk Metropolis-Hastings (RWMH) with adaptive proposal covariance. Supports both linear models (dsge_model) and nonlinear models (dsgenl_model). For nonlinear models, the steady state is re-solved and the model re-linearized at each candidate parameter vector.

Usage

bayes_dsge(
  model,
  data,
  priors,
  chains = 2L,
  iter = 5000L,
  warmup = floor(iter/2),
  thin = 1L,
  proposal_scale = 0.1,
  demean = TRUE,
  seed = NULL
)

Arguments

Details

The sampler operates in unconstrained space with appropriate transformations (log for positive parameters, logit for unit-interval parameters) and Jacobian corrections. The proposal covariance is adapted during warmup to target approximately 25\

Chains are initialized by drawing from the prior. If any draw yields a non-finite log-posterior, the starting values are jittered until a valid point is found.

For nonlinear models (dsgenl_model), each posterior evaluation involves: (1) solving the deterministic steady state at the candidate parameters, (2) computing a first-order linearization, (3) solving the resulting linear system, and (4) evaluating the Kalman filter likelihood. If any stage fails (e.g., steady-state non-convergence, Blanchard-Kahn violation), the proposal is safely rejected. This is computationally more expensive than linear Bayesian estimation.

Value

An object of class "dsge_bayes" containing posterior draws and diagnostics. For nonlinear models, the result also includes solve_failures, the number of parameter draws where steady-state, linearization, or solution failed.

Examples

m <- dsge_model(
  obs(y ~ z),
  state(z ~ rho * z),
  start = list(rho = 0.5)
)

set.seed(42)
z <- numeric(200); for (i in 2:200) z[i] <- 0.8 * z[i-1] + rnorm(1)
dat <- data.frame(y = z)

fit <- bayes_dsge(m, data = dat,
                  priors = list(rho = prior("beta", shape1 = 2, shape2 = 2)),
                  chains = 2, iter = 2000, seed = 1)
summary(fit)

Check Local Identification of DSGE Parameters

Description

Assesses local identification by computing the Jacobian of the mapping from structural parameters to model-implied autocovariance moments. Uses an SVD decomposition to detect rank deficiency (non-identification) and near-collinearity (weak identification).

Usage

check_identification(x, ...)

## S3 method for class 'dsge_fit'
check_identification(x, n_lags = 4L, tol = 1e-06, ...)

## S3 method for class 'dsge_bayes'
check_identification(x, n_lags = 4L, tol = 1e-06, ...)

Arguments

Details

The identification check constructs the moment vector m(\theta) = \mathrm{vec}(\Gamma(0), \Gamma(1), \ldots, \Gamma(K)) where \Gamma(k) is the autocovariance of observables at lag k, implied by the state-space solution. The Jacobian J = \partial m / \partial \theta is computed numerically.

A parameter is locally identified if the Jacobian has full column rank. If the rank is deficient, some linear combination of parameters cannot be distinguished from the data.

Per-parameter identification strength is measured by the norm of the corresponding Jacobian column: parameters with small column norms have little influence on the moments and may be weakly identified.

The condition number of J flags near-collinearity: a large condition number indicates that some parameter combinations are hard to distinguish.

Value

An object of class "dsge_identification" containing:

jacobian

The Jacobian matrix (n_moments x n_params).

svd

SVD decomposition of the Jacobian.

rank

Numerical rank of the Jacobian.

identified

Logical: are all parameters locally identified?

singular_values

Vector of singular values.

strength

Per-parameter identification strength (norm of corresponding Jacobian column).

condition_number

Condition number of the Jacobian.

param_names

Character vector of parameter names.

summary

Data frame with per-parameter diagnostics.

Examples

m <- dsge_model(
  obs(y ~ z),
  state(z ~ rho * z),
  start = list(rho = 0.5)
)
set.seed(1)
z <- numeric(100); for (i in 2:100) z[i] <- 0.8*z[i-1]+rnorm(1)
fit <- estimate(m, data = data.frame(y = z))
id <- check_identification(fit)
print(id)

Define a Linear DSGE Model

Description

Constructs a linear DSGE model object from a set of equations. Each equation is wrapped in obs(), unobs(), or state() to indicate the role of the left-hand-side variable.

Usage

dsge_model(..., fixed = list(), start = list())

Arguments

Details

A linear DSGE model is specified as a system of equations in which variables enter linearly but parameters may enter nonlinearly.

Use lead(x) or E(x) within formulas to denote the one-period-ahead model-consistent expectation of variable x.

The number of exogenous state variables (those with shocks) must equal the number of observed control variables.

Value

An object of class "dsge_model" containing:

equations

List of parsed equation objects.

variables

List with elements observed, unobserved, exo_state, endo_state.

parameters

Character vector of all parameter names.

free_parameters

Character vector of free (estimable) parameter names.

fixed

Named list of fixed parameter values.

start

Named list of starting values.

Examples

# Simple New Keynesian model
nk <- dsge_model(
  obs(p   ~ beta * lead(p) + kappa * x),
  unobs(x ~ lead(x) - (r - lead(p) - g)),
  obs(r   ~ psi * p + u),
  state(u ~ rhou * u),
  state(g ~ rhog * g),
  fixed = list(beta = 0.96),
  start = list(kappa = 0.1, psi = 1.5, rhou = 0.7, rhog = 0.9)
)

Define a Nonlinear DSGE Model

Description

Constructs a nonlinear DSGE model from string-based equations. Each equation is a character string of the form "LHS = RHS". State equations are detected automatically when the left-hand side is of the form VAR(+1).

Usage

dsgenl_model(
  ...,
  observed = character(0),
  unobserved = character(0),
  exo_state,
  endo_state = character(0),
  fixed = list(),
  start = list(),
  ss_guess = NULL,
  ss_function = NULL
)

Arguments

Details

Equations are written in standard mathematical notation with = as the equality sign and VAR(+1) for leads. For example:

"1/C = beta / C(+1) * (alpha * K^(alpha-1) + 1 - delta)"

State equations must have exactly one lead variable on the left-hand side (e.g., "K(+1) = K^alpha - C + (1 - delta) * K"). Control equations are all remaining equations.

Control equations must be provided in the order matching the controls vector (observed first, then unobserved).

Value

An object of class "dsgenl_model".

Examples

# Simple RBC model
rbc <- dsgenl_model(
  "1/C = beta / C(+1) * (alpha * exp(Z) * K^(alpha-1) + 1 - delta)",
  "K(+1) = exp(Z) * K^alpha - C + (1 - delta) * K",
  "Z(+1) = rho * Z",
  observed = "C",
  endo_state = "K",
  exo_state = "Z",
  fixed = list(alpha = 0.33, beta = 0.99, delta = 0.025),
  start = list(rho = 0.9)
)

Estimate a Linear DSGE Model by Maximum Likelihood

Description

Estimates the parameters of a linear DSGE model by maximizing the log-likelihood computed via the Kalman filter.

Usage

estimate(
  model,
  data,
  start = NULL,
  fixed = NULL,
  method = "BFGS",
  control = list(),
  shock_start = NULL,
  demean = TRUE,
  hessian = TRUE
)

Arguments

Details

The estimator optimizes over the structural parameters and the log standard deviations of the shocks. Shock standard deviations are parameterized in log-space to ensure positivity.

If the optimizer encounters parameter values for which the model is not saddle-path stable, the log-likelihood is set to -Inf.

Value

An object of class "dsge_fit".

Examples

# Define a simple AR(1) model
m <- dsge_model(
  obs(y ~ z),
  state(z ~ rho * z),
  start = list(rho = 0.5)
)

# Simulate some data
set.seed(42)
e <- rnorm(200)
z <- numeric(200)
for (i in 2:200) z[i] <- 0.8 * z[i-1] + e[i]
dat <- data.frame(y = z)

fit <- estimate(m, data = dat)
summary(fit)

Fitted values from a DSGE model

Description

Returns filtered fitted values of observed variables (using all information up to and including time t).

Usage

## S3 method for class 'dsge_fit'
fitted(object, ...)

Arguments

Value

A matrix of fitted values with columns named by observed variables.

Forecast from a DSGE Model

Description

Generic function for forecasting from estimated DSGE models.

Usage

forecast(object, ...)

Arguments

Value

A forecast object.

See Also

forecast.dsge_fit()

Forecast from a Fitted DSGE Model

Description

Produces dynamic multi-step forecasts from a fitted DSGE model. Forecasts are generated by iterating the state-space solution forward from the last filtered state.

Usage

## S3 method for class 'dsge_fit'
forecast(object, horizon = 12L, ...)

Arguments

Value

An object of class "dsge_forecast" containing:

forecasts

Data frame with columns period, variable, value.

horizon

The forecast horizon.

states

Matrix of forecasted state vectors.

Geweke convergence diagnostic

Description

Computes Geweke's (1992) convergence diagnostic, which compares the means of the first and last portions of each chain using a z-test.

Usage

geweke_test(object, ...)

Arguments

Value

An object of class "dsge_geweke" with z-scores and p-values for each parameter and chain.

Compute Impulse-Response Functions

Description

Computes the impulse-response functions (IRFs) from a fitted or solved DSGE model. An IRF traces the dynamic response of control and state variables to a one-standard-deviation shock.

Usage

## S3 method for class 'dsge_bayes'
irf(
  x,
  periods = 20L,
  impulse = NULL,
  response = NULL,
  se = TRUE,
  level = 0.95,
  n_draws = 200L,
  ...
)

irf(
  x,
  periods = 20L,
  impulse = NULL,
  response = NULL,
  se = TRUE,
  level = 0.95,
  ...
)

Arguments

Value

An object of class "dsge_irf" containing a data frame with columns: period, impulse, response, value, and optionally se, lower, upper.

Methods (by class)

• irf(dsge_bayes): Compute posterior IRFs from a Bayesian DSGE fit. Returns pointwise posterior median and credible bands.

Examples

m <- dsge_model(
  obs(y ~ z),
  state(z ~ rho * z),
  start = list(rho = 0.5)
)
sol <- solve_dsge(m, params = c(rho = 0.8))
irfs <- irf(sol, periods = 10)

Generalized IRFs Using Second-Order Approximation

Description

Computes impulse-response functions using the second-order solution. These differ from first-order IRFs because responses depend on the initial state and shock sign.

Usage

irf_2nd_order(sol, shock, size = 1, periods = 40L, initial = NULL)

Arguments

Value

A data frame with columns: period, variable, response, type.

Forward Lead Operator for DSGE Equations

Description

Marks a variable as a one-period-ahead model-consistent expectation in a DSGE equation formula. This is the core primitive for forward-looking variables.

Usage

lead(x, k = 1L)

Arguments

Details

lead(x) in a DSGE equation represents the expectation of variable x one period ahead, conditional on the model. In the literature, this corresponds to E_t[x_{t+1}].

This function is not meant to be called directly. It is recognized by the equation parser inside dsge_model().

Value

This function is not meant to be called directly; it always throws an error. It is recognized as a syntactic marker by the equation parser inside dsge_model().

See Also

E() for a user-friendly alias, dsge_model()

Linearize a Nonlinear DSGE Model

Description

Computes a first-order Taylor expansion of the nonlinear model around its deterministic steady state and returns the structural matrices in the canonical linear form.

Usage

linearize(model, steady_state, params = NULL)

Arguments

Value

A list of structural matrices (A0, A1, A2, A3, A4, B0, B1, B2, B3, C, D) plus the steady-state values. A4 captures lead-state coefficients in control equations (often zero).

Marginal likelihood estimation

Description

Estimates the log marginal likelihood using the modified harmonic mean estimator (Geweke, 1999). This provides a practical Bayesian model comparison tool via Bayes factors: BF12 = exp(logML1 - logML2).

Usage

marginal_likelihood(object, ...)

Arguments

Details

The harmonic mean estimator is known to be numerically unstable in some cases. The modified version (with truncation parameter tau) reduces this instability. Results should be interpreted with caution and compared across models only when both use similar MCMC settings.

Value

An object of class "dsge_marginal_likelihood".

MCMC diagnostic summary

Description

Comprehensive MCMC diagnostic summary combining ESS, R-hat, Geweke, and acceptance rate information.

Usage

mcmc_diagnostics(object, ...)

Arguments

Value

An object of class "dsge_mcmc_summary".

Model-implied covariance and correlation matrices

Description

Computes the unconditional (model-implied) covariance and correlation matrices of observable variables from a solved or estimated DSGE model. These are the theoretical second moments implied by the model at the given parameter values.

Usage

model_covariance(x, variables = NULL, n_lags = 0L, ...)

Arguments

Value

An object of class "dsge_covariance" containing:

covariance

Covariance matrix of selected variables.

correlation

Correlation matrix of selected variables.

std_dev

Standard deviations (square root of diagonal).

autocovariances

List of lagged autocovariance matrices (empty if n_lags = 0).

variables

Variable names.

n_lags

Number of autocovariance lags computed.

Examples

mod <- dsge_model(
  obs(pi ~ beta * lead(pi) + kappa * x),
  unobs(x ~ lead(x) - (r - lead(pi) - g)),
  obs(r ~ psi * pi + u),
  state(u ~ rhou * u),
  state(g ~ rhog * g),
  fixed = list(beta = 0.99),
  start = list(kappa = 0.1, psi = 1.5, rhou = 0.5, rhog = 0.5)
)
p <- list(kappa = 0.1, psi = 1.5, rhou = 0.5, rhog = 0.5)
s <- c(u = 0.5, g = 0.5)
sol <- solve_dsge(mod, params = p, shock_sd = s)
model_covariance(sol)

Create an Occasionally Binding Constraint

Description

Specifies an inequality constraint for use with simulate_occbin.

Usage

obc_constraint(variable, type = ">=", bound = 0, shock = NULL)

Arguments

Value

An obc_constraint object.

Examples

# Zero lower bound on nominal interest rate
obc_constraint("r", ">=", 0)

# Upper bound on debt ratio
obc_constraint("b", "<=", 0.6)

Define an Observed Control Variable Equation

Description

Wraps a formula to mark it as an equation for an observed control variable in a linear DSGE model.

Usage

obs(formula)

Arguments

Value

A list with class "dsge_equation" containing the parsed equation and its type.

See Also

unobs(), state(), dsge_model()

Parameter Sensitivity Analysis for DSGE Models

Description

Evaluates the sensitivity of key model outputs to one-at-a-time parameter perturbations. For each free parameter, the model is re-solved at theta +/- delta, and changes in the log-likelihood, impulse responses, steady state, and policy matrix are recorded.

Usage

parameter_sensitivity(x, ...)

## S3 method for class 'dsge_fit'
parameter_sensitivity(
  x,
  what = c("loglik", "irf"),
  delta = 0.01,
  irf_horizon = 20L,
  ...
)

## S3 method for class 'dsge_bayes'
parameter_sensitivity(
  x,
  what = c("loglik", "irf"),
  delta = 0.01,
  irf_horizon = 20L,
  ...
)

Arguments

Details

For each parameter \theta_j, the model is solved at \theta_j (1 + \delta) and \theta_j (1 - \delta). The numerical derivative is approximated as a central difference. Elasticities are reported as (\theta_j / f) \cdot (df / d\theta_j), representing the percentage change in the output for a 1 percent change in the parameter.

Value

An object of class "dsge_sensitivity" containing:

loglik

Data frame of log-likelihood sensitivities (if requested).

irf

Data frame of IRF sensitivities (if requested).

steady_state

Data frame of steady-state sensitivities (if requested).

policy

Data frame of policy matrix sensitivities (if requested).

param_names

Character vector of parameter names.

delta

Perturbation fraction used.

Examples

m <- dsge_model(
  obs(y ~ z),
  state(z ~ rho * z),
  start = list(rho = 0.5)
)
set.seed(1)
z <- numeric(100); for (i in 2:100) z[i] <- 0.8*z[i-1]+rnorm(1)
fit <- estimate(m, data = data.frame(y = z))
sa <- parameter_sensitivity(fit)
print(sa)

Perfect Foresight / Deterministic Transition Paths

Description

Simulate deterministic transition paths for DSGE models under perfect foresight. Supports temporary shocks, permanent shocks, and initial condition experiments using the linearized solution.

Usage

perfect_foresight(
  x,
  shocks = NULL,
  initial = NULL,
  horizon = 40L,
  params = NULL,
  shock_sd = NULL,
  in_sd = FALSE
)

Arguments

Details

The deterministic transition path is computed using the linearized state-space representation:

x_{t+1} = H x_t + M \varepsilon_{t+1}

y_t = G x_t

where x_t are state deviations from steady state, y_t are control deviations, and \varepsilon_t are deterministic shocks.

This uses the first-order linearized solution, so results are approximate for large shocks. For small to moderate shocks, the linearized paths are accurate.

Value

An object of class "dsge_perfect_foresight" containing:

states

Matrix (horizon x n_states) of state deviations from SS

controls

Matrix (horizon x n_controls) of control deviations

state_levels

Matrix of state levels (SS + deviation), if SS available

control_levels

Matrix of control levels, if SS available

steady_state

Named numeric vector of steady-state values

shock_path

Matrix (horizon x n_shocks) of applied shocks

initial

Named vector of initial state deviations

horizon

Integer horizon

state_names

Character vector of state names

control_names

Character vector of control names

shock_names

Character vector of shock names

H

State transition matrix used

G

Policy matrix used

M

Shock impact matrix used

Examples

# Simple AR(1) model
mod <- dsge_model(
  obs(p ~ x),
  state(x ~ rho * x),
  start = list(rho = 0.9)
)
sol <- solve_dsge(mod, params = list(rho = 0.9), shock_sd = c(x = 0.01))

# One-time shock at period 1
pf <- perfect_foresight(sol, shocks = list(x = 0.01), horizon = 40)
plot(pf)

# Displaced initial condition
pf2 <- perfect_foresight(sol, initial = c(x = 0.05), horizon = 40)
plot(pf2)

Plot Bayesian DSGE Results

Description

Produces diagnostic plots for posterior draws from a Bayesian DSGE fit.

Usage

## S3 method for class 'dsge_bayes'
plot(
  x,
  type = c("trace", "density", "prior_posterior", "running_mean", "acf", "pairs", "all",
    "irf"),
  pars = NULL,
  ...
)

Arguments

Details

All plot types except pairs and irf handle any number of parameters by paginating across multiple plot pages (up to 4 parameters per page). In interactive sessions, devAskNewPage() is used to prompt between pages.

For the "pairs" plot, at most 1000 draws are used to keep the plot readable. The correlation matrix is printed to the console.

Forecast plotting is not currently supported for Bayesian fits. Use irf() for posterior impulse-response analysis.

Value

Invisibly returns the dsge_bayes object x. Called for the side effect of producing diagnostic plots on the active graphics device.

Examples

m <- dsge_model(
  obs(y ~ z),
  state(z ~ rho * z),
  start = list(rho = 0.5)
)
set.seed(42)
z <- numeric(200); for (i in 2:200) z[i] <- 0.8 * z[i-1] + rnorm(1)
fit <- bayes_dsge(m, data = data.frame(y = z),
                  priors = list(rho = prior("beta", shape1 = 2, shape2 = 2)),
                  chains = 2, iter = 2000, seed = 1)
plot(fit, type = "trace")
plot(fit, type = "density")
plot(fit, type = "prior_posterior")
plot(fit, type = "running_mean")
plot(fit, type = "acf")
plot(fit, type = "all")
# Parameter selection
plot(fit, type = "trace", pars = "rho")

Plot Historical Shock Decomposition

Description

Creates a stacked bar chart showing the contribution of each structural shock to the observed variables over time.

Usage

## S3 method for class 'dsge_decomposition'
plot(x, which = NULL, ...)

Arguments

Value

No return value, called for the side effect of producing stacked bar charts of the historical shock decomposition on the active graphics device.

Plot DSGE Forecasts

Description

Plots forecast paths for observed variables.

Usage

## S3 method for class 'dsge_forecast'
plot(x, ...)

Arguments

Value

No return value, called for the side effect of producing forecast path plots on the active graphics device.

Plot Impulse-Response Functions

Description

Creates a multi-panel plot of impulse-response functions with optional confidence bands.

Usage

## S3 method for class 'dsge_irf'
plot(x, impulse = NULL, response = NULL, ci = TRUE, ...)

Arguments

Value

No return value, called for the side effect of producing a multi-panel impulse-response plot on the active graphics device.

Plot OccBin Simulation Results

Description

Plot OccBin Simulation Results

Usage

## S3 method for class 'dsge_occbin'
plot(x, vars = NULL, compare = TRUE, shade = TRUE, max_panels = 9L, ...)

Arguments

Value

No return value, called for the side effect of producing OccBin simulation plots on the active graphics device. Constrained and unconstrained paths are shown, with shaded regions where constraints bind.

Plot Perfect Foresight Transition Paths

Description

Plot the deterministic transition paths from a perfect_foresight result.

Usage

## S3 method for class 'dsge_perfect_foresight'
plot(x, vars = NULL, type = "deviation", max_panels = 9L, ...)

Arguments

Value

No return value, called for the side effect of producing transition path plots on the active graphics device.

Plot Smoothed States

Description

Plot Smoothed States

Usage

## S3 method for class 'dsge_smoothed'
plot(x, which = NULL, type = c("states", "fit"), ...)

Arguments

Value

No return value, called for the side effect of producing smoothed state or fit plots on the active graphics device.

Extract Policy Matrix

Description

Returns the policy matrix G from a fitted or solved DSGE model. The policy matrix maps state variables to control variables: y_t = G x_t.

Usage

policy_matrix(x, se = TRUE, level = 0.95)

Arguments

Value

If se = FALSE, returns the G matrix. If se = TRUE, returns a list with matrix, se, lower, upper, and a data frame table.

Posterior predictive check

Description

Simulates data from the posterior predictive distribution and compares summary statistics (variance, autocorrelation) with the observed data. This helps assess whether the estimated model can reproduce key features of the data.

Usage

posterior_predictive(object, ...)

Arguments

Value

An object of class "dsge_ppc" with posterior predictive distributions and p-values for each statistic.

Predict Method for DSGE Models

Description

Computes one-step-ahead predictions or filtered state estimates from a fitted DSGE model.

Usage

## S3 method for class 'dsge_fit'
predict(
  object,
  type = c("observed", "state"),
  method = c("onestep", "filter"),
  newdata = NULL,
  ...
)

Arguments

Value

A matrix of predictions or state estimates.

Prediction accuracy measures for a fitted DSGE model

Description

Computes root mean squared error (RMSE), mean absolute error (MAE), and mean error (bias) of one-step-ahead predictions.

Usage

prediction_accuracy(object, ...)

Arguments

Value

An object of class "dsge_prediction_accuracy" containing:

rmse

Named numeric vector of RMSE by variable.

mae

Named numeric vector of MAE by variable.

bias

Named numeric vector of mean error by variable.

n

Number of observations.

variables

Variable names.

Examples

  m <- dsge_model(
    obs(y ~ z),
    state(z ~ rho * z),
    start = list(rho = 0.5)
  )
  set.seed(1)
  z <- numeric(100); for (i in 2:100) z[i] <- 0.8 * z[i-1] + rnorm(1)
  fit <- estimate(m, data = data.frame(y = z))
  acc <- prediction_accuracy(fit)
  print(acc)

Prediction intervals for DSGE models

Description

Computes point predictions and prediction intervals using the one-step-ahead innovation variance from the Kalman filter.

Usage

prediction_interval(object, level = 0.95, ...)

Arguments

Value

An object of class "dsge_prediction_interval" with components fit, lower, upper, se, level, and variables.

Examples

  m <- dsge_model(
    obs(y ~ z),
    state(z ~ rho * z),
    start = list(rho = 0.5)
  )
  set.seed(1)
  z <- numeric(100); for (i in 2:100) z[i] <- 0.8 * z[i-1] + rnorm(1)
  fit <- estimate(m, data = data.frame(y = z))
  pi <- prediction_interval(fit, level = 0.95)
  print(pi)

Specify a Prior Distribution

Description

Creates a prior distribution object for use in Bayesian DSGE estimation.

Usage

prior(distribution, ...)

Arguments

Details

Distribution parameterizations:

normal

mean, sd

beta

shape1, shape2 (alpha, beta parameters)

gamma

shape, rate

uniform

min, max

inv_gamma

shape, scale — density: p(x) \propto x^{-(shape+1)} \exp(-scale/x)

Value

An object of class "dsge_prior".

Examples

prior("normal", mean = 0, sd = 1)
prior("beta", shape1 = 2, shape2 = 2)
prior("inv_gamma", shape = 0.01, scale = 0.01)

Prior-Posterior Update Diagnostics

Description

Computes diagnostics measuring how informative the data was for each parameter, by comparing the posterior distribution to the prior.

Usage

prior_posterior_update(x, ...)

## S3 method for class 'dsge_bayes'
prior_posterior_update(x, ...)

Arguments

Details

For each estimated parameter, the following diagnostics are computed:

sd_ratio

Ratio of posterior SD to prior SD. Values near 1 indicate the data was uninformative (posterior tracks the prior). Values much less than 1 indicate strong data information.

mean_shift

Absolute difference between posterior mean and prior mean, measured in units of prior SD. Large shifts indicate the data substantially updated beliefs.

update

Classification: "strong" if sd_ratio < 0.5 or mean_shift > 2; "moderate" if sd_ratio < 0.8 or mean_shift > 1; "weak" otherwise (posterior closely resembles the prior).

The SD ratio is the primary indicator of data informativeness. A parameter with sd_ratio close to 1 and small mean_shift is effectively determined by the prior, not the data.

Value

An object of class "dsge_prior_posterior" containing:

summary

Data frame with per-parameter diagnostics including prior mean/sd, posterior mean/sd, SD ratio, mean shift, and update classification.

param_names

Character vector of parameter names.

Examples

m <- dsge_model(
  obs(y ~ z),
  state(z ~ rho * z),
  start = list(rho = 0.5)
)
set.seed(1)
z <- numeric(100); for (i in 2:100) z[i] <- 0.8*z[i-1]+rnorm(1)
fit <- bayes_dsge(m, data = data.frame(y = z),
  priors = list(rho = prior("beta", shape1 = 2, shape2 = 2)),
  chains = 2, iter = 2000, seed = 1)
pp <- prior_posterior_update(fit)
print(pp)

Residuals from a fitted DSGE model

Description

Returns one-step-ahead prediction errors.

Usage

## S3 method for class 'dsge_fit'
residuals(object, ...)

Arguments

Value

A matrix of prediction errors.

Robust (sandwich) variance-covariance matrix

Description

Computes the sandwich (Huber-White) variance-covariance matrix for ML-estimated DSGE model parameters. This provides standard errors that are robust to model misspecification.

Usage

robust_vcov(object, ...)

Arguments

Details

The sandwich estimator is: V_robust = inv(H) B inv(H), where H is the Hessian of the negative log-likelihood and B is the outer product of the per-observation score vectors.

Value

An object of class "dsge_robust_vcov" containing:

vcov

Robust variance-covariance matrix.

se

Robust standard errors.

se_conventional

Conventional (Hessian-based) standard errors for comparison.

param_names

Parameter names.

Examples

  m <- dsge_model(
    obs(y ~ z),
    state(z ~ rho * z),
    start = list(rho = 0.5)
  )
  set.seed(1)
  z <- numeric(100); for (i in 2:100) z[i] <- 0.8 * z[i-1] + rnorm(1)
  fit <- estimate(m, data = data.frame(y = z))
  rv <- robust_vcov(fit)
  print(rv)

Historical Shock Decomposition

Description

Decomposes the observed variables into the contributions of each structural shock. At each time t, the observed deviation from steady state is written as a sum of contributions from current and past shocks plus the initial condition contribution.

Usage

shock_decomposition(x, ...)

## S3 method for class 'dsge_fit'
shock_decomposition(x, ...)

## S3 method for class 'dsge_bayes'
shock_decomposition(x, ...)

Arguments

Details

The state-space solution gives:

x_t = H^t x_0 + \sum_{j=1}^{t} H^{t-j} M \varepsilon_j

The historical decomposition partitions the observed variables y_t = Z x_t into the contribution of each structural shock \varepsilon_j^{(k)} accumulated through the propagation mechanism. The sum of all contributions (including the initial condition term) reproduces the smoothed observables exactly.

Value

An object of class "dsge_decomposition" containing:

decomposition

A 3D array with dimensions [T, n_obs, n_shocks + 1]. The last slice contains the initial condition contribution.

obs_names

Character vector of observed variable names.

shock_names

Character vector of shock names (plus "initial").

observed

T x n_obs matrix of observed data (deviations).

Examples

m <- dsge_model(
  obs(y ~ z),
  state(z ~ rho * z),
  start = list(rho = 0.5)
)
set.seed(1)
e <- rnorm(100)
z <- numeric(100); for (i in 2:100) z[i] <- 0.8*z[i-1]+e[i]
fit <- estimate(m, data = data.frame(y = z))
hd <- shock_decomposition(fit)
plot(hd)

Simulate Using Second-Order Approximation (Pruned)

Description

Simulates sample paths using the pruned second-order approximation following Kim, Kim, Schaumburg, and Sims (2008).

Usage

simulate_2nd_order(sol, n = 200L, n_burn = 100L, seed = NULL)

Arguments

Value

A list with states, controls, state_levels, control_levels matrices (n x n_vars).

Simulate with Occasionally Binding Constraints

Description

Computes deterministic transition paths under piecewise-linear occasionally binding constraints using an iterative shadow-shock method.

Usage

simulate_occbin(
  x,
  constraints,
  shocks = NULL,
  initial = NULL,
  horizon = 40L,
  max_iter = 100L,
  tol = 1e-08,
  in_sd = FALSE
)

Arguments

Details

The algorithm iteratively:

1. Simulates the path with current shadow shocks

2. Identifies periods where constraints are violated

3. Computes shadow shocks to enforce constraints at binding periods

4. Removes shadow shocks at non-binding periods

5. Repeats until the binding regime stabilises

This captures the feedback effect of constraint enforcement on future dynamics through the state transition.

Value

An object of class "dsge_occbin" containing:

states

Matrix of state deviations (constrained path)

controls

Matrix of control deviations (constrained path)

states_unc

Matrix of state deviations (unconstrained path)

controls_unc

Matrix of control deviations (unconstrained path)

binding

Logical matrix (horizon x n_constraints) of binding indicators

shadow_shocks

Matrix of shadow shocks applied

n_iter

Number of OccBin iterations to convergence

converged

Logical: did the algorithm converge?

constraints

List of constraint objects

steady_state

Steady state values if available

horizon

Simulation horizon

state_names

State variable names

control_names

Control variable names

Examples

# Simple NK model with ZLB
nk <- dsge_model(
  obs(pi ~ beta * lead(pi) + kappa * x),
  unobs(x ~ lead(x) - (r - lead(pi) - g)),
  obs(r ~ psi * pi + u),
  state(u ~ rhou * u),
  state(g ~ rhog * g),
  fixed = list(beta = 0.99, kappa = 0.1, psi = 1.5),
  start = list(rhou = 0.5, rhog = 0.5)
)
sol <- solve_dsge(nk, params = list(rhou = 0.5, rhog = 0.5),
                  shock_sd = c(u = 0.5, g = 0.5))
obc <- simulate_occbin(sol,
  constraints = list("r >= 0"),
  shocks = list(g = -0.05),
  horizon = 40)
plot(obc)

Extract Smoothed Structural Shocks

Description

Recovers the structural shocks from the smoothed states using the state transition equation: \hat{\varepsilon}_{t+1} = M^+ (x_{t+1|T} - H x_{t|T}) where M^+ is the Moore-Penrose pseudo-inverse of M.

Usage

smooth_shocks(x, ...)

## S3 method for class 'dsge_fit'
smooth_shocks(x, ...)

## S3 method for class 'dsge_bayes'
smooth_shocks(x, ...)

Arguments

Details

From the state transition x_{t+1} = H x_t + M \varepsilon_{t+1}, the smoothed innovation is x_{t+1|T} - H x_{t|T}. The structural shocks are recovered by projecting onto M: \hat{\varepsilon}_{t+1} = (M'M)^{-1} M' (x_{t+1|T} - H x_{t|T}).

Value

An object of class "dsge_smoothed_shocks" containing:

shocks

(T-1) x n_shocks matrix of smoothed structural shocks.

shock_names

Character vector of shock names.

Smoothed State Estimates from an Estimated DSGE Model

Description

Computes the Rauch-Tung-Striebel (RTS) smoother to produce optimal state estimates using all available observations. Compared to the filtered states (which only use past data), smoothed states also incorporate future observations.

Usage

smooth_states(x, ...)

## S3 method for class 'dsge_fit'
smooth_states(x, ...)

## S3 method for class 'dsge_bayes'
smooth_states(x, ...)

Arguments

Details

The smoother uses the state-space representation:

x_{t+1} = H x_t + M \varepsilon_{t+1}

y_t = Z x_t

where Z = D \cdot G. The smoothed states are the expectation of the state vector conditional on all observations: x_{t|T} = E[x_t | y_1, \ldots, y_T].

For Bayesian models, the smoother is evaluated at the posterior mean.

Value

An object of class "dsge_smoothed" containing:

smoothed_states

T x n_s matrix of smoothed state estimates.

filtered_states

T x n_s matrix of filtered state estimates.

smoothed_obs

T x n_obs matrix of smoothed observable fits.

residuals

T x n_obs matrix of observation residuals.

state_names

Character vector of state variable names.

obs_names

Character vector of observed variable names.

steady_state

Steady-state values (if available).

Examples

m <- dsge_model(
  obs(y ~ z),
  state(z ~ rho * z),
  start = list(rho = 0.5)
)
set.seed(1)
e <- rnorm(100)
z <- numeric(100); for (i in 2:100) z[i] <- 0.8 * z[i-1] + e[i]
fit <- estimate(m, data = data.frame(y = z))
sm <- smooth_states(fit)

Solve Second-Order Perturbation

Description

Computes the second-order approximation of a nonlinear DSGE model around its deterministic steady state.

Usage

solve_2nd_order(model, params, shock_sd, tol = 1e-06)

Arguments

Details

The second-order terms capture nonlinear effects including:

• Asymmetric responses to positive vs negative shocks

• Risk/precautionary effects (g_ss, h_ss corrections)

• State-dependent dynamics (quadratic policy terms)

Uses the method of Schmitt-Grohe and Uribe (2004).

Value

A dsge_solution object with additional second-order fields: order, g_xx, h_xx, g_ss, h_ss.

Solve a Linear or Linearized DSGE Model

Description

Computes the state-space solution of a DSGE model using the Klein (2000) method. Accepts both linear models (dsge_model) and nonlinear models (dsgenl_model). For nonlinear models, the steady state is computed and the model is linearized automatically.

Usage

solve_dsge(model, params = NULL, shock_sd = NULL, tol = 1e-06, order = 1L)

Arguments

Details

The method forms a companion system from the structural matrices and solves via undetermined coefficients iteration. Saddle-path stability requires that all eigenvalues of H have modulus less than 1.

For nonlinear models, the solver first computes the deterministic steady state, then linearizes the model via first-order Taylor expansion, and finally solves the resulting linear system.

Value

An object of class "dsge_solution" containing:

G

Policy matrix (n_controls x n_states).

H

State transition matrix (n_states x n_states).

M

Shock coefficient matrix (n_states x n_shocks).

D

Observation selection matrix.

eigenvalues

Complex vector of eigenvalues.

stable

Logical: is the system saddle-path stable?

n_stable

Number of stable eigenvalues.

params

The parameter values used.

model

Reference to the model object.

Check Stability of DSGE Model

Description

Checks saddle-path stability of the model by examining eigenvalues. A model is stable when the number of eigenvalues with modulus less than 1 equals the number of state variables.

Usage

stability(x)

Arguments

Value

A list of class "dsge_stability" with:

stable

Logical: is the system saddle-path stable?

eigenvalues

Complex eigenvalue vector.

moduli

Moduli of eigenvalues.

classification

Character vector: "stable" or "unstable" for each.

n_stable

Number of stable eigenvalues.

n_states

Number of state variables (required stable count).

Define a State Variable Equation

Description

Wraps a formula to mark it as an equation for a state variable in a linear DSGE model. State equations describe the evolution of state variables one period ahead.

Usage

state(formula, shock = TRUE)

Arguments

Value

A list with class "dsge_equation" containing the parsed equation and its type.

See Also

obs(), unobs(), dsge_model()

Solve for the Deterministic Steady State

Description

Finds the steady-state values of all model variables by solving the system of nonlinear equations with all leads set equal to current values.

Usage

steady_state(model, ...)

## S3 method for class 'dsgenl_model'
steady_state(
  model,
  params = NULL,
  guess = NULL,
  maxiter = 200L,
  tol = 1e-10,
  ...
)

Arguments

Value

An object of class "dsgenl_steady_state" with:

values

Named numeric vector of steady-state values.

residuals

Equation residuals at the solution.

params

Parameter values used.

converged

Logical: did the solver converge?

iterations

Number of iterations used.

Summary of Perfect Foresight Transition

Description

Summary of Perfect Foresight Transition

Usage

## S3 method for class 'dsge_perfect_foresight'
summary(object, ...)

Arguments

Value

Invisibly returns the dsge_perfect_foresight object. Called for the side effect of printing impact effects, peak deviations, and convergence diagnostics to the console.

Extract State Transition Matrix

Description

Returns the transition matrix H from a fitted or solved DSGE model. The transition matrix describes state evolution: x_{t+1} = H x_t + M \varepsilon_{t+1}.

Usage

transition_matrix(x, se = TRUE, level = 0.95)

Arguments

Value

Same structure as policy_matrix().

Define an Unobserved Control Variable Equation

Description

Wraps a formula to mark it as an equation for an unobserved control variable in a linear DSGE model.

Usage

unobs(formula)

Arguments

Value

A list with class "dsge_equation" containing the parsed equation and its type.

See Also

obs(), state(), dsge_model()

Robust vcov via vcov generic

Description

When type = "robust" is passed to vcov(), returns the sandwich variance-covariance matrix.

Usage

## S3 method for class 'dsge_fit'
vcov(object, type = c("conventional", "robust"), ...)

Arguments

Value

Variance-covariance matrix.