| Title: | Quarterly Projection Models for Monetary Policy Analysis |
|---|---|
| Description: | An end-to-end implementation of the semi-structural quarterly projection models used in central-bank forecasting and policy analysis systems: model declaration with model-consistent expectations, a generalized Schur solver with Blanchard-Kahn diagnostics following Klein (2000) <doi:10.1016/S0165-1889(99)00045-7>, Kalman filtering and smoothing for latent states such as the output gap and the neutral rate, historical shock decompositions, conditional forecasts that distinguish announced from unanticipated policy paths, an auditable judgment ledger, forecast rounds with revision decompositions, Bayesian estimation with identification diagnostics following Iskrev (2010) <doi:10.1016/j.jmoneco.2009.12.007>, and reporting. The canonical small open economy model of Berg, Karam and Laxton (2006) <doi:10.5089/9781451863413.001> ships as a calibrated template, with extension blocks for disaggregated food inflation and managed exchange rates. |
| Authors: | Mustapha Mohammed [aut, cre] |
| Maintainer: | Mustapha Mohammed <[email protected]> |
| License: | MIT + file LICENSE |
| Version: | 1.1.0 |
| Built: | 2026-09-29 16:30:19 UTC |
| Source: | https://github.com/cran/qpmR |
Apply an extension block to a model
add_block(model, block)add_block(model, block)
model |
A |
block |
A |
The extended qpm_model, with the block recorded in
meta$blocks.
m <- add_block(qpm_template("bkl"), block_food_cpi(weight = 0.4)) mm <- add_block(qpm_template("bkl"), block_food_cpi(weight = 0.4)) m
Central-bank forecasts are never raw model output: the desk knows
about the announced electricity tariff, the tax change, the one-off
the model cannot see. add_judgment() makes that adjustment a
first-class, logged operation: you state the change you want (in
percentage points, relative to the current forecast), qpmR back-solves
the structural shocks needed to support it while keeping the whole
forecast model-consistent, records who imposed it and why, and flags
judgment that requires implausibly large shocks.
add_judgment( fc, ..., author = "desk", rationale = "", anticipated = NULL, instruments = NULL )add_judgment( fc, ..., author = "desk", rationale = "", anticipated = NULL, instruments = NULL )
fc |
A |
... |
Named adjustments: one argument per variable, each a named
vector of additions (percentage points, relative to the current
forecast) by period, e.g. |
author |
Who is imposing the judgment (logged). |
rationale |
Why (logged; make it meaningful – the ledger is the audit trail read back before the policy meeting). |
anticipated |
Expectation mode for the re-solve; defaults to the forecast's current mode, or unanticipated. |
instruments |
Shocks allowed to move; defaults to the forecast's current instruments, or all shocks. |
Judgment entries are stored as absolute targets, so the ledger is
replayable; the full set of conditions and judgment is re-solved
jointly each time. Inspect the ledger with judgment_log().
The adjusted qpm_forecast with the entry appended to its
judgment ledger.
sol <- qpm_solve(qpm_template("bkl")) fc <- qpm_forecast(sol, horizon = 8) fc2 <- add_judgment(fc, pi = c(h2 = 0.5), author = "desk", rationale = "announced electricity tariff increase") judgment_log(fc2)sol <- qpm_solve(qpm_template("bkl")) fc <- qpm_forecast(sol, horizon = 8) fc2 <- add_judgment(fc, pi = c(h2 = 0.5), author = "desk", rationale = "announced electricity tariff increase") judgment_log(fc2)
Recalibrate a model at an estimate's point values
apply_estimate(est, type = "mean")apply_estimate(est, type = "mean")
est |
A |
type |
Passed to |
The estimation model recalibrated at the chosen point estimate (structural parameters and shock sds).
# see ?qpm_estimate# see ?qpm_estimate
Splits headline inflation into food and core. Food is 30-50 percent of the consumption basket across most of sub-Saharan Africa and South Asia, and a single-inflation model is unusable there: supply shocks to food dominate headline, but monetary policy should look through the relative-price component. Practically every technical-assistance engagement rebuilds this split by hand.
block_food_cpi( weight = 0.35, persistence = 0.5, demand = 0.2, passthrough = 0.25, correction = 0.1, sd = 3 )block_food_cpi( weight = 0.35, persistence = 0.5, demand = 0.2, passthrough = 0.25, correction = 0.1, sd = 3 )
weight |
Food share of the CPI basket ( |
persistence |
Food inflation persistence ( |
demand |
Output-gap coefficient in food inflation ( |
passthrough |
Real-exchange-rate coefficient in food inflation
( |
correction |
Error-correction speed on the relative food price
( |
sd |
Standard deviation of the food supply shock. |
The block replaces headline inflation with an identity and adds:
pi_core — the Phillips curve, now for core inflation (it keeps the
template's b1, b2, b3 and the eps_pi shock);
pi_food — food inflation: its own persistence, expectations of
headline, demand, a stronger exchange-rate pass-through than core,
and error correction on the relative food price;
rp_food — the relative food price gap, which accumulates the
food-core inflation differential and mean-reverts through f4.
Headline is pi = w_food * pi_food + (1 - w_food) * pi_core, so
everything downstream (the 4-quarter average, the Fisher equation, the
policy rule) continues to use headline. To target core instead,
replace the policy rule with another block.
A qpm_block().
m <- add_block(qpm_template("bkl"), block_food_cpi(weight = 0.45)) sol <- qpm_solve(m) plot(irf(sol, shock = "eps_pifood"), vars = c("pi_food", "pi", "pi_core", "i"))m <- add_block(qpm_template("bkl"), block_food_cpi(weight = 0.45)) sol <- qpm_solve(m) plot(irf(sol, shock = "eps_pifood"), vars = c("pi_food", "pi", "pi_core", "i"))
Turns the template's floating exchange rate into a managed one. A
leaning-against-the-wind rule responds to the real exchange rate gap,
and intervention enters the UIP block directly, so the same model
spans a continuum of regimes: intensity = 0 is a free float,
moderate values a managed float, and large values approach a peg.
Program countries — where reserves, not just the policy rate, are the
operative instrument — live in the middle of that range.
block_fx_intervention(intensity = 1, persistence = 0.6, sd = 1)block_fx_intervention(intensity = 1, persistence = 0.6, sd = 1)
intensity |
Scales both the intervention response to the RER gap
and its effect on the exchange rate. |
persistence |
Persistence of intervention ( |
sd |
Standard deviation of the discretionary intervention shock. |
fx_int is intervention intensity, positive meaning sales of foreign
exchange in support of the domestic currency (which appreciates the
real exchange rate, lowering q).
A qpm_block().
float <- qpm_solve(qpm_template("bkl")) managed <- qpm_solve(add_block(qpm_template("bkl"), block_fx_intervention(intensity = 1))) # the same risk-premium shock moves the exchange rate less under managementfloat <- qpm_solve(qpm_template("bkl")) managed <- qpm_solve(add_block(qpm_template("bkl"), block_fx_intervention(intensity = 1))) # the same risk-premium shock moves the exchange rate less under management
Produces the chart set a forecast round is discussed from: the forecast fans, the filtered latent states, the historical shock decomposition of inflation, and the model's monetary transmission. Output is a multi-page PDF by default, or numbered PNGs.
chart_pack( round, file = NULL, charts = NULL, vars = NULL, width = 9, height = 6.5 )chart_pack( round, file = NULL, charts = NULL, vars = NULL, width = 9, height = 6.5 )
round |
A |
file |
Output file. A |
charts |
Which charts to include, from |
vars |
Variables for the forecast page. Default: the model's
headline set, whichever of |
width, height
|
Page size in inches. |
The output path (invisibly), or NULL when drawing to the
current device.
m <- qpm_template("bkl") obs <- simulate(qpm_solve(m), nsim = 40, seed = 1, burn = 20) obs$period <- next_quarters("2016-Q1", 40) r <- qpm_round("demo", m, obs[, c("period", "pi", "i", "q")], horizon = 8) pdf_path <- file.path(tempdir(), "pack.pdf") chart_pack(r, pdf_path)m <- qpm_template("bkl") obs <- simulate(qpm_solve(m), nsim = 40, seed = 1, burn = 20) obs$period <- next_quarters("2016-Q1", 40) r <- qpm_round("demo", m, obs[, c("period", "pi", "i", "q")], horizon = 8) pdf_path <- file.path(tempdir(), "pack.pdf") chart_pack(r, pdf_path)
Decomposes the forecast revision between two rounds – "inflation for 2027-Q1 is 0.4pp higher than we said in June: why?" – into the contributions of new data (outturns), data revisions, calibration changes, changed conditions, and judgment, by re-running the full pipeline swapping one ingredient at a time. The contributions telescope, so they sum to the total revision exactly; the final step is verified against the new round's archived forecast.
compare_rounds(old, new, variables = NULL, store = "rounds") ## S3 method for class 'qpm_revision' plot(x, variable = NULL, ...)compare_rounds(old, new, variables = NULL, store = "rounds") ## S3 method for class 'qpm_revision' plot(x, variable = NULL, ...)
old, new
|
|
variables |
Variables to decompose (default: all). |
store |
Round store used when |
x |
A |
variable |
Variable to plot. |
... |
Unused. |
Requirements: both rounds must use the same model structure and the
same observables, and their data must carry "YYYY-Qq" period labels
so the calendars align. Comparison covers the calendar quarters both
rounds forecast.
An object of class qpm_revision: a data frame with one row
per variable and overlap period, columns old, new, total, and
the five contributions. Print shows waterfall tables; plot()
draws stacked contribution bars across the overlap.
m <- qpm_template("bkl") obs <- simulate(qpm_solve(m), nsim = 44, seed = 1, burn = 20) obs$period <- next_quarters("2015-Q4", 44) rA <- qpm_round("June", m, obs[1:40, c("period", "pi", "i", "q")], horizon = 8) rB <- qpm_round("September", m, obs[, c("period", "pi", "i", "q")], horizon = 8) rB <- add_judgment(rB, pi = stats::setNames(0.3, rB$forecast$periods[2]), author = "desk", rationale = "tariff") rev <- compare_rounds(rA, rB) revm <- qpm_template("bkl") obs <- simulate(qpm_solve(m), nsim = 44, seed = 1, burn = 20) obs$period <- next_quarters("2015-Q4", 44) rA <- qpm_round("June", m, obs[1:40, c("period", "pi", "i", "q")], horizon = 8) rB <- qpm_round("September", m, obs[, c("period", "pi", "i", "q")], horizon = 8) rB <- add_judgment(rB, pi = stats::setNames(0.3, rB$forecast$periods[2]), author = "desk", rationale = "tariff") rev <- compare_rounds(rA, rB) rev
Quarterly data for Czechia, 1996Q1 onward, in the units and sign
conventions of the qpm_template() model, ready for
qpm_filter(). Czechia is the canonical FPAS economy: a small open
inflation targeter (since 1998) with three decades of clean data, a
large disinflation, a currency-crisis start, the GFC, a floor episode,
COVID, and the 2022 inflation wave.
czechiaczechia
A data frame with one row per quarter and columns:
Quarter label, "1996-Q1" style.
CPI inflation, QoQ annualised percent, seasonally adjusted
with stats::stl() on the quarterly log index.
CPI inflation, year on year percent (no adjustment needed).
3-month PRIBOR, percent p.a., quarterly average – a proxy for the CNB policy rate.
Real CZK/EUR exchange rate, 100 times log, CPI-based, normalized so the 2015 average is zero; an increase is a real depreciation of the koruna.
Real GDP growth, QoQ annualised percent, from seasonally and calendar adjusted chain-linked volumes.
3-month EURIBOR, percent p.a., quarterly average.
Euro-area HICP inflation, year on year percent.
Missing values are genuine ragged edges (e.g. the euro exists only
from 1999, so q, istar start later); qpm_filter() handles them.
Observe pi4 or pi, not both – they are linked by an identity
and the filter will report the collinearity.
Compiled by data-raw/czechia.R from FRED series
CLVMNACSCAB1GQCZ (Eurostat national accounts), CZECPIALLMINMEI,
IR3TIB01CZM156N, IR3TIB01EZM156N (OECD Main Economic
Indicators), CP0000EZ19M086NEST (Eurostat HICP), and the ECB
reference exchange rate EXR/Q.CZK.EUR.SP00.A. Retrieved 2026-08-22.
head(czechia) m <- qpm_calibrate(qpm_template("bkl", trends = "rw"), pi_tar = 2, istar_ss = 2, pistar_ss = 2, prem_ss = 1) cz <- czechia[czechia$period >= "1999", c("period", "pi4", "i", "q", "dy_obs", "istar", "pistar")] fit <- qpm_filter(m, cz) fithead(czechia) m <- qpm_calibrate(qpm_template("bkl", trends = "rw"), pi_tar = 2, istar_ss = 2, pistar_ss = 2, prem_ss = 1) cz <- czechia[czechia$period >= "1999", c("period", "pi4", "i", "q", "dy_obs", "istar", "pistar")] fit <- qpm_filter(m, cz) fit
E() marks model-consistent expectations inside eqs() declarations,
e.g. E(pi[+1]). It is parsed by qpm_model() and never evaluated as
an ordinary function.
E(x)E(x)
x |
A variable reference such as |
No return value; calling E() outside of model equations is
an error by design.
Generalized eigenvalues of a solved model
eigen_table(x)eigen_table(x)
x |
A |
A data frame with the moduli of the generalized eigenvalues of the model companion pencil, sorted ascending, with stability flags.
head(eigen_table(qpm_solve(qpm_template("bkl"))))head(eigen_table(qpm_solve(qpm_template("bkl"))))
Each equation is a two-sided formula. Lags and leads use index
notation: x[-1] is the first lag, x[+1] the first lead. Wrap
expectations in E(): E(pi[+1]) is the model-consistent expectation
of next quarter's inflation. Longer lags and leads (e.g. E(pi4[+4]))
are handled automatically via auxiliary state variables.
eqs(...)eqs(...)
... |
Two-sided formulas. |
A list of formulas.
eqs( pi ~ b1 * pi[-1] + (1 - b1) * E(pi[+1]) + b2 * y_gap + eps_pi )eqs( pi ~ b1 * pi[-1] + (1 - b1) * E(pi[+1]) + b2 * y_gap + eps_pi )
Splits the forecast error variance of each variable at each horizon into the contributions of the structural shocks — "how much of inflation uncertainty two years out is the cost-push shock?". For
the h-step forecast error variance is
and the share of shock k is the same sum with only column k of
Q active.
Shares sum to one across shocks for every variable and horizon.
fevd(x, ...) ## S3 method for class 'qpm_solution' fevd(x, horizon = 24, vars = NULL, shocks = NULL, ...) ## S3 method for class 'qpm_fevd' plot(x, var = NULL, drop_zero = TRUE, ...)fevd(x, ...) ## S3 method for class 'qpm_solution' fevd(x, horizon = 24, vars = NULL, shocks = NULL, ...) ## S3 method for class 'qpm_fevd' plot(x, var = NULL, drop_zero = TRUE, ...)
x |
A |
... |
Passed to methods. |
horizon |
Largest forecast horizon. |
vars |
Variables to include; default all declared variables. |
shocks |
Shocks to include; default all. |
var |
Variable to plot. |
drop_zero |
Omit shocks that never contribute. |
Shares remain well defined when the model has unit roots, even though the variances themselves grow without bound.
A data frame of class qpm_fevd in long format with columns
variable, shock, horizon, share and variance.
print() shows the dominant shocks; plot() draws stacked shares.
sol <- qpm_solve(qpm_template("bkl")) fv <- fevd(sol, horizon = 20) fv plot(fv, var = "pi")sol <- qpm_solve(qpm_template("bkl")) fv <- fevd(sol, horizon = 20) fv plot(fv, var = "pi")
Impulse response functions
irf(x, ...) ## S3 method for class 'qpm_solution' irf(x, shock = NULL, horizon = 24, size = NULL, vars = NULL, ...)irf(x, ...) ## S3 method for class 'qpm_solution' irf(x, shock = NULL, horizon = 24, size = NULL, vars = NULL, ...)
x |
A |
... |
Passed to methods. |
shock |
Shock name(s); default all shocks. |
horizon |
Number of quarters after impact. |
size |
Shock size(s). Default is one standard deviation (from the
model's |
vars |
Variables to include; default all declared (non-auxiliary) variables. |
A data frame of class qpm_irf in long format with columns
shock, variable, horizon, value (deviations from steady
state). Printing shows peak effects; plot() draws panel charts.
sol <- qpm_solve(qpm_template("bkl")) ir <- irf(sol, shock = "eps_i", horizon = 16) ir plot(ir, vars = c("pi", "y_gap", "i", "q"))sol <- qpm_solve(qpm_template("bkl")) ir <- irf(sol, shock = "eps_i", horizon = 16) ir plot(ir, vars = c("pi", "y_gap", "i", "q"))
Print a forecast's judgment ledger
judgment_log(fc)judgment_log(fc)
fc |
A |
The ledger data frame, invisibly.
sol <- qpm_solve(qpm_template("bkl")) fc <- add_judgment(qpm_forecast(sol, horizon = 8), pi = c(h2 = 0.5), author = "desk", rationale = "tariff") judgment_log(fc)sol <- qpm_solve(qpm_template("bkl")) fc <- add_judgment(qpm_forecast(sol, horizon = 8), pi = c(h2 = 0.5), author = "desk", rationale = "tariff") judgment_log(fc)
For a qpm_filtration this is the Kalman-filter log-likelihood of
the data under the calibrated model, with zero degrees of freedom
(nothing was estimated). For a qpm_estimate it is the
log-likelihood at the posterior mode (or at the maximum for
method = "mle"), with degrees of freedom equal to the number of
estimated parameters — so stats::AIC() and stats::BIC() work.
## S3 method for class 'qpm_filtration' logLik(object, ...) ## S3 method for class 'qpm_estimate' logLik(object, ...)## S3 method for class 'qpm_filtration' logLik(object, ...) ## S3 method for class 'qpm_estimate' logLik(object, ...)
object |
A |
... |
Unused. |
An object of class logLik.
sol <- qpm_solve(qpm_template("bkl")) obs <- simulate(sol, nsim = 40, seed = 1, burn = 20) fit <- qpm_filter(sol, obs[, c("period", "pi", "i", "q")]) logLik(fit)sol <- qpm_solve(qpm_template("bkl")) obs <- simulate(sol, nsim = 40, seed = 1, burn = 20) fit <- qpm_filter(sol, obs[, c("period", "pi", "i", "q")]) logLik(fit)
Computes the log marginal likelihood of a Bayesian
estimate – the quantity whose differences across models on the same
data are log Bayes factors. Two estimators are reported:
marginal_likelihood(est, taus = seq(0.1, 0.9, by = 0.2))marginal_likelihood(est, taus = seq(0.1, 0.9, by = 0.2))
est |
A |
taus |
Truncation probabilities for the modified harmonic mean. |
Modified harmonic mean (Geweke 1999) from the posterior draws, computed for a range of truncation probabilities; a small spread across truncations indicates a reliable estimate.
Laplace approximation at the posterior mode in transformed space (parametrization-invariant because the Jacobian is included), with a numerically differenced Hessian.
Truncated priors created with truncate() are renormalized
numerically at construction, so they contribute proper densities
here.
An object of class qpm_logml: $logml (harmonic-mean
estimate, averaged over taus), $by_tau, $laplace, and the
spread across truncations.
Geweke, J. (1999). Using simulation methods for Bayesian econometric models. Econometric Reviews, 18(1), 1-73.
m <- qpm_model(variables = vars(x = "x"), shocks = shocks(e), equations = eqs(x ~ rho * x[-1] + e), params = list(rho = 0.5)) obs <- simulate(qpm_solve(qpm_calibrate(m, rho = 0.8)), nsim = 120, seed = 1) est <- qpm_estimate(m, obs, priors(rho = beta(0.5, 0.2)), iter = 300, chains = 2, seed = 2, verbose = FALSE) marginal_likelihood(est)m <- qpm_model(variables = vars(x = "x"), shocks = shocks(e), equations = eqs(x ~ rho * x[-1] + e), params = list(rho = 0.5)) obs <- simulate(qpm_solve(qpm_calibrate(m, rho = 0.8)), nsim = 120, seed = 1) est <- qpm_estimate(m, obs, priors(rho = beta(0.5, 0.2)), iter = 300, chains = 2, seed = 2, verbose = FALSE) marginal_likelihood(est)
The standard calibration check: does the model reproduce the
volatilities and persistence actually observed? Reports the
population standard deviation and autocorrelations implied by the
solved model — from the stationary covariance
and —
next to the same statistics computed from data, plus the shock that
accounts for most of each variable's unconditional variance.
model_properties(x, data = NULL, vars = NULL, lags = c(1, 4))model_properties(x, data = NULL, vars = NULL, lags = c(1, 4))
x |
A |
data |
Optional data frame of observations in levels (columns
named for model variables, an optional |
vars |
Variables to report; default all declared variables. |
lags |
Autocorrelation orders to report. |
Population moments exist only for stationary models. When the model
has unit roots (random-walk trends) they are undefined, and the
function reports that rather than returning nonsense; use
fevd() and qpm_filter() diagnostics instead.
An object of class qpm_properties: a data frame with the
model and (optionally) data moments.
sol <- qpm_solve(qpm_template("bkl")) model_properties(sol, vars = c("y_gap", "pi", "i", "q")) # against simulated data obs <- simulate(sol, nsim = 200, seed = 5, burn = 50) model_properties(sol, data = obs, vars = c("y_gap", "pi", "i"))sol <- qpm_solve(qpm_template("bkl")) model_properties(sol, vars = c("y_gap", "pi", "i", "q")) # against simulated data obs <- simulate(sol, nsim = 200, seed = 5, burn = 50) model_properties(sol, data = obs, vars = c("y_gap", "pi", "i"))
Generate consecutive quarter labels
next_quarters(last, H)next_quarters(last, H)
last |
The starting quarter, |
H |
Number of quarters to generate. |
Character vector of H quarter labels.
next_quarters("2026-Q3", 4)next_quarters("2026-Q3", 4)
Number of observations
## S3 method for class 'qpm_filtration' nobs(object, ...) ## S3 method for class 'qpm_estimate' nobs(object, ...)## S3 method for class 'qpm_filtration' nobs(object, ...) ## S3 method for class 'qpm_estimate' nobs(object, ...)
object |
A |
... |
Unused. |
Number of time periods used.
sol <- qpm_solve(qpm_template("bkl")) obs <- simulate(sol, nsim = 40, seed = 1, burn = 20) nobs(qpm_filter(sol, obs[, c("period", "pi", "i")]))sol <- qpm_solve(qpm_template("bkl")) obs <- simulate(sol, nsim = 40, seed = 1, burn = 20) nobs(qpm_filter(sol, obs[, c("period", "pi", "i")]))
Produces a forecast whose fan bands integrate over the posterior:
for each sampled parameter draw the model is re-solved, the data
re-filtered (so the initial state reflects that draw), and a path
sampled from the resulting Gaussian predictive; pointwise band
quantiles are then taken across draws. The bands therefore combine
future-shock uncertainty with parameter uncertainty — they are wider
than the fixed-parameter fans of qpm_forecast().
posterior_forecast(est, horizon = 12, ndraws = 200, bands = c(0.5, 0.7, 0.9))posterior_forecast(est, horizon = 12, ndraws = 200, bands = c(0.5, 0.7, 0.9))
est |
A |
horizon |
Forecast horizon in quarters. |
ndraws |
Number of posterior draws to propagate. |
bands |
Coverage levels for the fan. |
A qpm_forecast-compatible object (posterior-mean path,
quantile bands, smoothed history from the posterior-mean model);
$posterior records the number of draws used.
m <- qpm_model(variables = vars(x = "x"), shocks = shocks(e), equations = eqs(x ~ rho * x[-1] + e), params = list(rho = 0.5)) obs <- simulate(qpm_solve(qpm_calibrate(m, rho = 0.8)), nsim = 100, seed = 1) est <- qpm_estimate(m, obs, priors(rho = beta(0.5, 0.2)), iter = 300, chains = 1, seed = 2, verbose = FALSE) fc <- posterior_forecast(est, horizon = 8, ndraws = 15) plot(fc, vars = "x")m <- qpm_model(variables = vars(x = "x"), shocks = shocks(e), equations = eqs(x ~ rho * x[-1] + e), params = list(rho = 0.5)) obs <- simulate(qpm_solve(qpm_calibrate(m, rho = 0.8)), nsim = 100, seed = 1) est <- qpm_estimate(m, obs, priors(rho = beta(0.5, 0.2)), iter = 300, chains = 1, seed = 2, verbose = FALSE) fc <- posterior_forecast(est, horizon = 8, ndraws = 15) plot(fc, vars = "x")
One named argument per estimated quantity: a structural parameter name or a shock name (meaning that shock's standard deviation). Everything without a prior stays calibrated.
priors(...)priors(...)
... |
Named prior declarations, e.g.
|
Available distributions (usable only inside priors(), so base R's
beta() and gamma() functions are never masked):
normal(mean, sd)
beta(mean, sd) — on (0, 1), mean/sd parametrization
gamma(mean, sd) — on (0, Inf)
invgamma(mean, sd) — on (0, Inf); the usual choice for shock sds
uniform(min, max)
truncate(d, lower, upper) — restrict any of the above; the
normalizing constant is dropped (harmless for modes and MCMC)
An object of class qpm_priors.
p <- priors( b1 = beta(0.7, 0.1), b2 = gamma(0.25, 0.1), c2 = truncate(normal(1.5, 0.25), lower = 1), eps_pi = invgamma(1, 0.3) ) pp <- priors( b1 = beta(0.7, 0.1), b2 = gamma(0.25, 0.1), c2 = truncate(normal(1.5, 0.25), lower = 1), eps_pi = invgamma(1, 0.3) ) p
A block is a documented, reusable bundle of variables, shocks,
parameters and equations that adapts a template to a country. Blocks
are how a technical-assistance engagement customizes the canonical
model without forking it: add_block() applies one, qpm_diff()
shows exactly what changed, and the result is an ordinary
qpm_model.
qpm_block( name, variables = NULL, shocks = NULL, params = list(), sigma = NULL, equations = NULL, description = "" )qpm_block( name, variables = NULL, shocks = NULL, params = list(), sigma = NULL, equations = NULL, description = "" )
name |
Block name, used in printing and in the model's history. |
variables |
New variables, from |
shocks |
New shocks, from |
params |
Named list of new parameters (optional). |
sigma |
Named vector of standard deviations for the new shocks (optional; default 1). |
equations |
Equations, from |
description |
One-line description of what the block does. |
Equations whose left-hand side names a variable that already exists replace that variable's equation; equations for newly declared variables are appended. Left-hand sides must be bare variable names.
An object of class qpm_block.
block_food_cpi(), block_fx_intervention()
b <- qpm_block( name = "risk premium shifter", params = list(prem_extra = 0), equations = eqs(prem ~ rho_prem * prem[-1] + (1 - rho_prem) * (prem_ss + prem_extra) + eps_prem) ) m <- add_block(qpm_template("bkl"), b)b <- qpm_block( name = "risk premium shifter", params = list(prem_extra = 0), equations = eqs(prem ~ rho_prem * prem[-1] + (1 - rho_prem) * (prem_ss + prem_extra) + eps_prem) ) m <- add_block(qpm_template("bkl"), b)
Update a model's calibration
qpm_calibrate(model, ..., sigma = NULL)qpm_calibrate(model, ..., sigma = NULL)
model |
A |
... |
Named parameter values to update. Every name must already exist in the model (typo protection). |
sigma |
Optional named vector of shock standard deviations to update. |
The updated qpm_model.
m <- qpm_template("bkl") m <- qpm_calibrate(m, b2 = 0.3, c2 = 2)m <- qpm_template("bkl") m <- qpm_calibrate(m, b2 = 0.3, c2 = 2)
qpm_diff() compares model structure — what a country team changed
in the code. This compares model behaviour: the transmission of a
given shock and the moments each specification implies. When a
calibration is revised, both questions matter, and the second is the
one a policy audience asks.
qpm_compare_models(models, shock = NULL, vars = NULL, horizon = 20) ## S3 method for class 'qpm_model_comparison' plot(x, vars = NULL, ...)qpm_compare_models(models, shock = NULL, vars = NULL, horizon = 20) ## S3 method for class 'qpm_model_comparison' plot(x, vars = NULL, ...)
models |
A named list of |
shock |
Shock whose impulse responses are compared. Default: the first shock common to every model. |
vars |
Variables to compare. Default: those common to all models, capped at the usual headline set. |
horizon |
Impulse-response horizon. |
x |
A |
... |
Unused. |
An object of class qpm_model_comparison holding the
impulse responses and, for stationary models, the implied moments.
base <- qpm_template("bkl") flat <- qpm_calibrate(base, b2 = 0.05) # a much flatter Phillips curve cmp <- qpm_compare_models(list(baseline = base, flat = flat), shock = "eps_y") cmp plot(cmp, vars = c("pi", "i"))base <- qpm_template("bkl") flat <- qpm_calibrate(base, b2 = 0.05) # a much flatter Phillips curve cmp <- qpm_compare_models(list(baseline = base, flat = flat), shock = "eps_y") cmp plot(cmp, vars = c("pi", "i"))
Imposes hard conditions on future values of model variables and finds
the minimum-norm structural shocks (in standard-deviation units,
optionally restricted to instruments) that deliver them. This is how
a policy question becomes a forecast: "what if the policy rate is held
at 3.5 for four quarters?" or "what paths are consistent with
inflation back at target by 2027?".
qpm_condition(fc, ..., anticipated = FALSE, instruments = NULL)qpm_condition(fc, ..., anticipated = FALSE, instruments = NULL)
fc |
A |
... |
Named conditions: one argument per variable, each a named
vector of levels by period, e.g.
|
anticipated |
Logical; announced-at-start ( |
instruments |
Character vector of shocks allowed to move; default all shocks. |
The anticipated switch is the economics: TRUE means the whole
conditioned path is announced at the start of the forecast, so
expectations react immediately (forward-looking variables move before
the conditioning bites); FALSE means the implied shocks arrive as
period-by-period surprises. The two produce materially different
paths in any forward-looking model, and practitioners frequently do
not know which one their tools assume.
Fan bands are recomputed as the Gaussian conditional distribution of
the path given the conditions (using all shocks), so conditioned
points have (near) zero width. The mean path uses only the chosen
instruments; when instruments is restricted, mean and bands answer
slightly different questions – see Antolin-Diaz, Petrella and
Rubio-Ramirez (2021) for the full treatment, which is on the qpmR
roadmap.
The conditioned qpm_forecast, with $shocks_implied (raw
units), $shocks_implied_std (standard deviations), and the
conditions recorded. Printing summarizes the implied shocks and
flags any larger than two standard deviations.
sol <- qpm_solve(qpm_template("bkl")) fc <- qpm_forecast(sol, horizon = 8) hold <- qpm_condition(fc, i = c(h1 = 9.5, h2 = 9.5, h3 = 9.5, h4 = 9.5), anticipated = TRUE, instruments = "eps_i") holdsol <- qpm_solve(qpm_template("bkl")) fc <- qpm_forecast(sol, horizon = 8) hold <- qpm_condition(fc, i = c(h1 = 9.5, h2 = 9.5, h3 = 9.5, h4 = 9.5), anticipated = TRUE, instruments = "eps_i") hold
Rewrites history with some shocks switched off or scaled: "what if the central bank had simply followed its rule through 2022?" is the path implied by setting the policy shocks to zero over that window and re-running the model from the same starting point with all other shocks unchanged.
qpm_counterfactual(fit, shocks, periods = NULL, factor = 0, label = NULL) ## S3 method for class 'qpm_counterfactual' plot(x, vars = NULL, ...)qpm_counterfactual(fit, shocks, periods = NULL, factor = 0, label = NULL) ## S3 method for class 'qpm_counterfactual' plot(x, vars = NULL, ...)
fit |
A |
shocks |
Shocks to modify. |
periods |
Periods over which to modify them (labels as in the data, or integer indices). Default: the whole sample. |
factor |
Multiplier applied to the selected shocks; |
label |
Optional name for the scenario. |
x |
A |
vars |
Variables to plot. |
... |
Unused. |
This differs from qpm_decompose(), which attributes the history that
happened; here the history is replayed under a different assumption.
The counterfactual is only as good as the model's invariance to the
intervention — a Lucas-critique caveat that applies to every exercise
of this kind and is worth stating in any write-up.
An object of class qpm_counterfactual holding the actual and
counterfactual paths and their difference.
sol <- qpm_solve(qpm_template("bkl")) obs <- simulate(sol, nsim = 60, seed = 4, burn = 20) obs$period <- next_quarters("2010-Q4", 60) fit <- qpm_filter(sol, obs[, c("period", "pi", "i", "q")]) cf <- qpm_counterfactual(fit, shocks = "eps_i", label = "no policy surprises") cf plot(cf, vars = c("pi", "i", "y_gap"))sol <- qpm_solve(qpm_template("bkl")) obs <- simulate(sol, nsim = 60, seed = 4, burn = 20) obs$period <- next_quarters("2010-Q4", 60) fit <- qpm_filter(sol, obs[, c("period", "pi", "i", "q")]) cf <- qpm_counterfactual(fit, shocks = "eps_i", label = "no policy surprises") cf plot(cf, vars = c("pi", "i", "y_gap"))
Splits the smoothed history of every variable into the additive
contributions of each structural shock plus the carry-over of the
pre-sample initial state: with smoothed shocks ,
The contributions sum exactly to the smoothed state (deviations from steady state); this identity is verified internally.
qpm_decompose(fit, vars = NULL) ## S3 method for class 'qpm_decomposition' plot(x, var = NULL, drop_zero = TRUE, periods = NULL, ...)qpm_decompose(fit, vars = NULL) ## S3 method for class 'qpm_decomposition' plot(x, var = NULL, drop_zero = TRUE, periods = NULL, ...)
fit |
A |
vars |
Variables to keep (default: all declared variables). |
x |
A |
var |
Variable to plot. |
drop_zero |
Drop components that never contribute. |
periods |
Optional integer window of period indices to display
(e.g. |
... |
Unused. |
A long data frame of class qpm_decomposition with columns
period, variable, component (shock names plus "initial"),
and value (contribution, in deviations from steady state).
plot() draws a stacked-bar decomposition with the smoothed total
overlaid.
sol <- qpm_solve(qpm_template("bkl")) obs <- simulate(sol, nsim = 60, seed = 3, burn = 20) fit <- qpm_filter(sol, obs[, c("period", "pi", "i", "q")]) dec <- qpm_decompose(fit) plot(dec, var = "pi")sol <- qpm_solve(qpm_template("bkl")) obs <- simulate(sol, nsim = 60, seed = 3, burn = 20) fit <- qpm_filter(sol, obs[, c("period", "pi", "i", "q")]) dec <- qpm_decompose(fit) plot(dec, var = "pi")
Shows what a customization actually changed: variables, shocks and parameters added or removed, equations changed, and recalibrations. This is how a model review works when country teams adapt a template.
qpm_diff(old, new, tol = 1e-10)qpm_diff(old, new, tol = 1e-10)
old, new
|
|
tol |
Relative tolerance for calling a parameter changed. |
An object of class qpm_model_diff, printed as a report.
qpm_diff(qpm_template("bkl"), add_block(qpm_template("bkl"), block_food_cpi()))qpm_diff(qpm_template("bkl"), add_block(qpm_template("bkl"), block_food_cpi()))
Turns annual series into quarterly ones consistent with the annual totals. Many of the economies these models are built for publish national accounts annually, so a quarterly projection model has to start by constructing quarterly GDP — usually from an indicator such as industrial production, imports or credit.
qpm_disaggregate( annual, indicator = NULL, frequency = 4, method = c("denton", "chow-lin"), conversion = c("sum", "average"), rho = NULL ) ## S3 method for class 'qpm_disaggregation' plot(x, ...)qpm_disaggregate( annual, indicator = NULL, frequency = 4, method = c("denton", "chow-lin"), conversion = c("sum", "average"), rho = NULL ) ## S3 method for class 'qpm_disaggregation' plot(x, ...)
annual |
Numeric vector of low-frequency values. |
indicator |
Optional numeric vector of high-frequency indicator
values, length |
frequency |
Periods per low-frequency observation (4 for annual-to-quarterly). |
method |
|
conversion |
|
rho |
AR(1) coefficient for |
x |
A |
... |
Unused. |
Two standard methods:
"denton" — Denton-Cholette proportional first differences.
Minimises the squared change in the ratio of the quarterly series to
the indicator (or, without an indicator, in the series itself),
subject to matching the annual figures. Purely a smoothing method:
no regression, no parameters.
"chow-lin" — generalised least squares on the indicator with
AR(1) quarterly residuals, distributing the annual residual across
quarters. Uses the indicator's regression relationship, so it is the
better choice when the indicator genuinely tracks the target.
Both enforce the aggregation constraint exactly: "sum" for flows
(annual GDP is the sum of quarters), "average" for stocks and index
levels.
An object of class qpm_disaggregation: the high-frequency
series, the method used and the fitted parameters.
Denton, F. T. (1971); Chow, G. C. and Lin, A. (1971).
# annual GDP with a quarterly indicator set.seed(1) q_true <- cumsum(rnorm(40, 0.5)) + 100 annual <- colSums(matrix(q_true, nrow = 4)) ind <- q_true + rnorm(40, 0, 1) d <- qpm_disaggregate(annual, ind, method = "chow-lin") d plot(d)# annual GDP with a quarterly indicator set.seed(1) q_true <- cumsum(rnorm(40, 0.5)) + 100 annual <- colSums(matrix(q_true, nrow = 4)) ind <- q_true + rnorm(40, 0, 1) d <- qpm_disaggregate(annual, ind, method = "chow-lin") d plot(d)
Estimates any subset of structural parameters and shock standard
deviations from data, using the Kalman-filter likelihood of the solved
model. method = "bayes" finds the posterior mode (in transformed,
unconstrained space), seeds an adaptive random-walk Metropolis sampler
with the inverse Hessian, and returns posterior draws with split
R-hat and effective-sample-size diagnostics. method = "mle" uses the
same machinery with flat priors on the declared supports.
qpm_estimate( model, data, priors, observables = NULL, measurement_error = 0, kappa = 1e+06, method = c("bayes", "mle"), iter = 4000, burn = floor(iter/2), chains = 2, thin = 1, seed = NULL, verbose = TRUE ) ## S3 method for class 'qpm_estimate' coef(object, type = c("mean", "mode", "median"), ...)qpm_estimate( model, data, priors, observables = NULL, measurement_error = 0, kappa = 1e+06, method = c("bayes", "mle"), iter = 4000, burn = floor(iter/2), chains = 2, thin = 1, seed = NULL, verbose = TRUE ) ## S3 method for class 'qpm_estimate' coef(object, type = c("mean", "mode", "median"), ...)
model |
A |
data |
Data frame in levels (as for |
priors |
A |
observables, measurement_error, kappa
|
Passed to the filter. |
method |
|
iter |
MCMC iterations per chain (including burn-in). |
burn |
Burn-in iterations (default |
chains |
Number of chains (run sequentially). |
thin |
Keep every |
seed |
Optional RNG seed. The previous state of the random number generator is restored on exit, so a seeded call does not disturb the caller's random stream. |
verbose |
Report progress with |
object |
A |
type |
Point estimate: posterior |
... |
Unused. |
Draws that violate Blanchard-Kahn (indeterminacy or explosiveness) receive zero posterior weight — the usual truncation of the prior to the determinacy region. Everything without a prior stays calibrated at its current value.
An object of class qpm_estimate: posterior draws (natural
units), the mode, acceptance rate, split R-hat and effective
sample sizes, and the originating model/data. Use coef() to
extract point estimates and apply_estimate() to recalibrate the
model.
m <- qpm_model(variables = vars(x = "x"), shocks = shocks(e), equations = eqs(x ~ rho * x[-1] + e), params = list(rho = 0.5)) obs <- simulate(qpm_solve(qpm_calibrate(m, rho = 0.8)), nsim = 100, seed = 1) est <- qpm_estimate(m, obs, priors(rho = beta(0.5, 0.2), e = invgamma(1, 0.3)), iter = 300, chains = 2, seed = 2, verbose = FALSE) estm <- qpm_model(variables = vars(x = "x"), shocks = shocks(e), equations = eqs(x ~ rho * x[-1] + e), params = list(rho = 0.5)) obs <- simulate(qpm_solve(qpm_calibrate(m, rho = 0.8)), nsim = 100, seed = 1) est <- qpm_estimate(m, obs, priors(rho = beta(0.5, 0.2), e = invgamma(1, 0.3)), iter = 300, chains = 2, seed = 2, verbose = FALSE) est
Runs the Kalman filter and RTS smoother over the solved model, jointly inferring every latent variable — output gap, neutral rate, equilibrium exchange rate, trend processes — and the historical structural shocks from whatever subset of variables you actually observe. Missing values (ragged edges, gappy series) are handled naturally.
qpm_filter(x, data, observables = NULL, measurement_error = 0, kappa = 1e+06)qpm_filter(x, data, observables = NULL, measurement_error = 0, kappa = 1e+06)
x |
A |
data |
A data frame in levels (model units). Columns whose names
match declared variables are used as observables; an optional
|
observables |
Optional character vector restricting which columns are treated as observed. |
measurement_error |
Measurement-error standard deviation(s): scalar or named vector over observables. Defaults to 0. |
kappa |
Diffuse-prior variance scale used when the model has
unit-root (random-walk) trends; see |
The filter is initialized at the model's stationary distribution (Lyapunov covariance), which is exact for the stationary models qpmR currently supports.
An object of class qpm_filtration: smoothed states in levels
($states), their standard errors ($se), smoothed structural
shocks ($shocks), the log-likelihood ($loglik), innovation
diagnostics ($diag), and the full expanded-state matrix
($states_dev). Feed it to qpm_decompose() for historical shock
decompositions or to qpm_forecast() to forecast from the smoothed
current state.
sol <- qpm_solve(qpm_template("bkl")) obs <- simulate(sol, nsim = 60, seed = 3, burn = 20) fit <- qpm_filter(sol, obs[, c("period", "pi", "i", "q")]) fit plot(fit, vars = c("y_gap", "r_bar"))sol <- qpm_solve(qpm_template("bkl")) obs <- simulate(sol, nsim = 60, seed = 3, burn = 20) fit <- qpm_filter(sol, obs[, c("period", "pi", "i", "q")]) fit plot(fit, vars = c("y_gap", "r_bar"))
Iterates the solved model forward from an initial state and computes analytic forecast uncertainty from the shock variances,
giving Gaussian fan bands around the mean path. The result can then be
conditioned on assumed paths with qpm_condition(), shifted by shock
scenarios with qpm_scenario(), or adjusted with logged judgment via
add_judgment().
qpm_forecast( object, from = NULL, horizon = 12, bands = c(0.5, 0.7, 0.9), sigma = NULL )qpm_forecast( object, from = NULL, horizon = 12, bands = c(0.5, 0.7, 0.9), sigma = NULL )
object |
A |
from |
Initial state: a |
horizon |
Forecast horizon in quarters. |
bands |
Coverage levels for the fan, e.g. |
sigma |
Optional named vector of shock standard deviations. |
An object of class qpm_forecast: a list with paths (long
data frame: variable, h, mean, and lo_*/hi_* per band, in
levels), forecast-period labels in $periods, plus the machinery
needed for conditioning.
sol <- qpm_solve(qpm_template("bkl")) histq <- simulate(sol, nsim = 40, seed = 7, burn = 20) fc <- qpm_forecast(sol, from = histq, horizon = 12) fc plot(fc, vars = c("pi", "i", "y_gap", "q"))sol <- qpm_solve(qpm_template("bkl")) histq <- simulate(sol, nsim = 40, seed = 7, burn = 20) fc <- qpm_forecast(sol, from = histq, horizon = 12) fc plot(fc, vars = c("pi", "i", "y_gap", "q"))
Checks, before any estimation is run, whether the chosen parameters can be told apart by the data. Two Jacobians are analysed numerically at the current calibration, in the spirit of Iskrev (2010):
qpm_identify(model, params = NULL, observables = NULL, lags = 3, h = 1e-05)qpm_identify(model, params = NULL, observables = NULL, lags = 3, h = 1e-05)
model |
A |
params |
Parameters to check: a character vector of structural
parameter and/or shock names, or a |
observables |
Observed variables the moment analysis conditions on. Default: all declared variables. |
lags |
Autocovariance lags in the moment vector. |
h |
Relative step for the central differences. |
solution level: derivatives of the solved transition, shock loading, and observable steady state with respect to the parameters. Rank deficiency here means some parameter movements do not change the model's solution at all.
moment level (stationary models only): derivatives of the
observables' first and second moments (means, and autocovariances
up to lags). Rank deficiency here means some parameter movements
are observationally equivalent in population.
The report names parameters with (numerically) no effect, parameter combinations spanning any null space, and near-collinear pairs of Jacobian columns (correlation above 0.995) that are only jointly identified.
An object of class qpm_identification with the ranks,
singular values, and flagged parameters; printed as a verdict list.
Iskrev, N. (2010). Local identification in DSGE models. Journal of Monetary Economics, 57(2), 189-202.
qpm_identify(qpm_template("bkl"), params = c("b1", "b2", "b3", "c1", "c2"), observables = c("pi", "i", "q", "dy_obs"))qpm_identify(qpm_template("bkl"), params = c("b1", "b2", "b3", "c1", "c2"), observables = c("pi", "i", "q", "dy_obs"))
Runs static checks (unused shocks/parameters, bare leads outside
E(), calibrations outside plausible ranges when the template
documents them) and then attempts to solve the model, reporting the
Blanchard-Kahn outcome. Unit algebra on declared units is on the
roadmap.
qpm_lint(model)qpm_lint(model)
model |
A |
An object of class qpm_lint (a data frame of checks with
status ok, note, warn, or fail), printed with markers.
qpm_lint(qpm_template("bkl"))qpm_lint(qpm_template("bkl"))
Declares a linear semi-structural model with model-consistent
expectations. Equations are parsed immediately, so unknown symbols and
malformed lag/lead references are caught at construction time.
Coefficients are extracted from the current calibration when the model
is solved with qpm_solve().
qpm_model( name = "QPM model", variables, shocks, equations, params = list(), sigma = NULL, meta = list() )qpm_model( name = "QPM model", variables, shocks, equations, params = list(), sigma = NULL, meta = list() )
name |
Model name used in printed output. |
variables |
Variable declarations from |
shocks |
Shock declarations from |
equations |
Equations from |
params |
Named list/vector of parameter values. |
sigma |
Named vector of shock standard deviations. Defaults to 1 for every shock. |
meta |
Optional list of metadata (e.g. plausible parameter
|
An object of class qpm_model.
m <- qpm_model( name = "AR(1) toy", variables = vars(x = "A persistent process"), shocks = shocks(eps_x), equations = eqs(x ~ rho * x[-1] + eps_x), params = list(rho = 0.8) ) mm <- qpm_model( name = "AR(1) toy", variables = vars(x = "A persistent process"), shocks = shocks(eps_x), equations = eqs(x ~ rho * x[-1] + eps_x), params = list(rho = 0.8) ) m
Turns a forecast round into the document a policy meeting is run from: an executive summary with the numbers filled in, the forecast table and fan charts, the filtered gaps, the shock decomposition, the judgment ledger, an optional revision decomposition against the previous round, and a reproducibility appendix.
qpm_report( round, file, compare_to = NULL, store = "rounds", render = TRUE, engine = c("auto", "quarto", "rmarkdown"), quiet = TRUE )qpm_report( round, file, compare_to = NULL, store = "rounds", render = TRUE, engine = c("auto", "quarto", "rmarkdown"), quiet = TRUE )
round |
A |
file |
Output path; required, so that nothing is written unless a
location is named. The extension chooses the format
( |
compare_to |
Optional previous |
store |
Round store used to resolve |
render |
Render the document, or only write the source. |
engine |
|
quiet |
Suppress rendering progress output. |
The report is always written as a self-contained .Rmd source, on
the principle that institutions replace the template's text, not
its plumbing: edit the file, re-render, keep the analysis. Rendering
additionally needs pandoc (via Quarto or RStudio/Positron); where
none is available — air-gapped forecasting machines, bare CI
runners — qpm_report() reports that and returns the .Rmd
unrendered rather than failing.
The path actually produced (the rendered document, or the
.Rmd when rendering was not possible), invisibly.
m <- qpm_template("bkl") obs <- simulate(qpm_solve(m), nsim = 40, seed = 1, burn = 20) obs$period <- next_quarters("2016-Q1", 40) r <- qpm_round("demo", m, obs[, c("period", "pi", "i", "q")], horizon = 8) src <- qpm_report(r, file.path(tempdir(), "mpr.Rmd"), render = FALSE)m <- qpm_template("bkl") obs <- simulate(qpm_solve(m), nsim = 40, seed = 1, burn = 20) obs$period <- next_quarters("2016-Q1", 40) r <- qpm_round("demo", m, obs[, c("period", "pi", "i", "q")], horizon = 8) src <- qpm_report(r, file.path(tempdir(), "mpr.Rmd"), render = FALSE)
Central banks rarely believe their fan charts are symmetric: the
published judgement is usually "risks to inflation are tilted to the
upside". qpm_risk() applies that judgement, replacing the Gaussian
bands of a forecast with two-piece normal bands whose mode stays on
the model's projection while the mean shifts by the stated skew.
Total uncertainty is preserved: the variance implied by the model at
each horizon is held fixed, so a skew redistributes risk rather than
adding it.
qpm_risk(fc, ..., author = "MPC", rationale = "")qpm_risk(fc, ..., author = "MPC", rationale = "")
fc |
A |
... |
Named skews: one argument per variable, each a named vector
of mean minus mode by period, e.g. |
author, rationale
|
Recorded with the risk assessment, as for judgement. |
The skew is stated in the variable's own units as mean minus mode —
a value of 0.3 on inflation means the risks are worth 0.3
percentage points to the upside. Skews may be given for any subset of
variables and horizons; anything unstated keeps symmetric bands.
Unlike add_judgment(), which moves the projection itself and
back-solves the shocks that support it, this changes only the shape of
the distribution around an unchanged central path.
The forecast with skewed bands. $risk records the skews and
$paths gains a mode column alongside mean.
sol <- qpm_solve(qpm_template("bkl")) fc <- qpm_forecast(sol, horizon = 8) risky <- qpm_risk(fc, pi = 0.4, author = "MPC", rationale = "energy prices tilted to the upside") risky plot(risky, vars = c("pi", "i"))sol <- qpm_solve(qpm_template("bkl")) fc <- qpm_forecast(sol, horizon = 8) risky <- qpm_risk(fc, pi = 0.4, author = "MPC", rationale = "energy prices tilted to the upside") risky plot(risky, vars = c("pi", "i"))
A forecast round binds everything that produced a forecast – the
model and its calibration, the data vintage, the filtration, the
forecast with its conditions and judgment – into one object that can
be saved, reloaded, re-run, and compared with later rounds via
compare_rounds(). This is the object an institution archives: next
quarter, "why did the forecast move?" is answered from the rounds, not
from memory.
qpm_round( name, model, data, observables = NULL, horizon = 12, bands = c(0.5, 0.7, 0.9), measurement_error = 0, kappa = 1e+06 )qpm_round( name, model, data, observables = NULL, horizon = 12, bands = c(0.5, 0.7, 0.9), measurement_error = 0, kappa = 1e+06 )
name |
Round name, e.g. |
model |
A |
data |
Data frame in levels with a |
observables |
Columns treated as observed (default: all columns matching declared variables). |
horizon |
Forecast horizon in quarters. |
bands |
Fan coverage levels. |
measurement_error, kappa
|
Passed to |
qpm_round() runs the standard pipeline (solve, filter, baseline
forecast). Conditions, scenarios, and judgment are then applied to the
round directly: qpm_condition(), qpm_scenario(), and
add_judgment() all accept a round and update its forecast.
An object of class qpm_round with elements model, data,
solution, fit, and forecast.
m <- qpm_template("bkl") sol <- qpm_solve(m) obs <- simulate(sol, nsim = 40, seed = 1, burn = 20) obs$period <- next_quarters("2016-Q1", 40) r <- qpm_round("test round", m, obs[, c("period", "pi", "i", "q")], horizon = 8) rm <- qpm_template("bkl") sol <- qpm_solve(m) obs <- simulate(sol, nsim = 40, seed = 1, burn = 20) obs$period <- next_quarters("2016-Q1", 40) r <- qpm_round("test round", m, obs[, c("period", "pi", "i", "q")], horizon = 8) r
Answers the question a policy committee actually asks — what if we responded differently? — by re-solving the model over a grid of rule parameters and scoring each one by the unconditional loss
computed from the model's stationary covariance rather than by simulation, so it is exact. Tracing the resulting variance pairs gives the inflation-output variability frontier (the Taylor curve).
qpm_rule_eval(model, grid, loss = c(pi = 1, y_gap = 0.5), diff_loss = NULL) ## S3 method for class 'qpm_rule_eval' plot(x, xvar = NULL, yvar = NULL, ...)qpm_rule_eval(model, grid, loss = c(pi = 1, y_gap = 0.5), diff_loss = NULL) ## S3 method for class 'qpm_rule_eval' plot(x, xvar = NULL, yvar = NULL, ...)
model |
A |
grid |
A data frame of parameter values, one row per rule and one
column per parameter (e.g. from |
loss |
Named weights on the variances of levels, e.g.
|
diff_loss |
Named weights on the variances of first differences,
e.g. |
x |
A |
xvar, yvar
|
Axes of the frontier. Either a variable name, whose
level variance is used, or a scored column name directly — so
|
... |
Unused. |
Rules that violate Blanchard-Kahn are reported as such rather than dropped: a policy response too weak to deliver determinacy is a finding, not a missing row.
A data frame of class qpm_rule_eval: the grid, the variance
of each targeted variable, the loss, and the Blanchard-Kahn outcome.
m <- qpm_template("bkl") grid <- expand.grid(c2 = c(1.2, 1.5, 2, 3), c3 = c(0, 0.5, 1)) ev <- qpm_rule_eval(m, grid, loss = c(pi = 1, y_gap = 0.5), diff_loss = c(i = 0.5)) ev plot(ev)m <- qpm_template("bkl") grid <- expand.grid(c2 = c(1.2, 1.5, 2, 3), c3 = c(0, 0.5, 1)) ev <- qpm_rule_eval(m, grid, loss = c(pi = 1, y_gap = 0.5), diff_loss = c(i = 0.5)) ev plot(ev)
Adds a specified path of structural shocks to a forecast – "what if
oil pushes the exchange rate 10 percent weaker next quarter?" –
without any inversion. Under anticipated = TRUE the shock path is
announced at the start of the forecast and expectations react ahead
of it.
qpm_scenario(fc, shocks, anticipated = FALSE, label = NULL)qpm_scenario(fc, shocks, anticipated = FALSE, label = NULL)
fc |
A |
shocks |
Named list: one entry per shock, each a named vector of
shock sizes (raw equation units) by period, e.g.
|
anticipated |
Announced at start ( |
label |
Optional scenario label for printing. |
The shifted qpm_forecast (bands unchanged: a deterministic
scenario shifts the mean, not the uncertainty).
sol <- qpm_solve(qpm_template("bkl")) fc <- qpm_forecast(sol, horizon = 8) dep <- qpm_scenario(fc, shocks = list(eps_q = c(h1 = 3)), label = "10pct depreciation")sol <- qpm_solve(qpm_template("bkl")) fc <- qpm_forecast(sol, horizon = 8) dep <- qpm_scenario(fc, shocks = list(eps_q = c(h1 = 3)), label = "10pct depreciation")
Reduces the model to first-order form (adding auxiliary states for lags/leads beyond one quarter), computes the steady state, and solves for the unique stable rational-expectations solution
via the generalized Schur (QZ) decomposition (Klein 2000), with full Blanchard-Kahn diagnostics.
qpm_solve(model, tol = 1e-07)qpm_solve(model, tol = 1e-07)
model |
A |
tol |
Numerical tolerance for the solution residual check. |
An object of class qpm_solution with elements P, Q
(transition and impact matrices over the expanded state vector),
ss (steady state), and an eigenvalue table (see eigen_table()).
Klein, P. (2000). Using the generalized Schur form to solve a multivariate linear rational expectations model. Journal of Economic Dynamics and Control, 24(10), 1405-1423.
sol <- qpm_solve(qpm_template("bkl")) solsol <- qpm_solve(qpm_template("bkl")) sol
"bkl" is the canonical semi-structural small-open-economy quarterly
projection model in the tradition of Berg, Karam and Laxton (2006, IMF
WP/06/80-81): an IS curve, a hybrid Phillips curve, a forward-looking
inflation-targeting policy rule, and a dampened (hybrid) UIP block,
plus equilibrium-trend and foreign processes, and an observation block
for real GDP growth (dy_obs = dy_bar + 4 * (y_gap - y_gap[-1])), so
the model can be filtered on actual national-accounts data without
modelling the level of potential output. The default calibration is
illustrative, for a higher-inflation emerging economy ("Meridia"); it
is not any actual country. See czechia for a real dataset and a
matching recalibration example.
qpm_template( name = c("bkl", "bkl_food", "managed_fx"), trends = c("stationary", "rw") )qpm_template( name = c("bkl", "bkl_food", "managed_fx"), trends = c("stationary", "rw") )
name |
Template name: |
trends |
|
trends selects the equilibrium processes:
"stationary" (default): AR(1) trends anchored at steady-state
parameters; the model is fully stationary.
"rw": driftless random walks for the equilibrium real exchange
rate and potential growth (whose levels are pure normalizations);
the neutral rate stays anchored by real interest parity
(rstar + prem), since a free random walk there would leave
steady-state gaps indeterminate. The model then has unit roots:
free trend levels are normalized to minimum norm in the steady
state, and qpm_filter() switches to diffuse initialization
automatically. This is the configuration for real data, where
trends drift.
Conventions: gaps in percentage points; inflation QoQ annualised;
interest rates in percent per annum; q is 100 times the log real
exchange rate, an increase is a real depreciation; dy_obs is QoQ
annualised real GDP growth.
Country-shaped variants are shipped as shortcuts for the canonical
model plus an extension block (see add_block()):
"bkl_food" — headline CPI split into food and core
(block_food_cpi()), the configuration for economies where food is
a large share of the basket.
"managed_fx" — a leaning-against-the-wind intervention rule
entering the UIP block (block_fx_intervention()).
A calibrated qpm_model.
Berg, A., Karam, P., and Laxton, D. (2006). A Practical Model-Based Approach to Monetary Policy Analysis - Overview. IMF Working Paper 06/80; and the companion How-To guide, IMF WP 06/81.
m <- qpm_template("bkl") summary(m) m_rw <- qpm_template("bkl", trends = "rw")m <- qpm_template("bkl") summary(m) m_rw <- qpm_template("bkl", trends = "rw")
qpmR ships a compiled (C++) Kalman filter and an equivalent reference implementation in R. The compiled one is used by default because estimation runs it once per posterior draw; the R one is kept because the two agreeing to machine precision is what makes the compiled path trustworthy, and it is useful when debugging.
qpm_use_cpp()qpm_use_cpp()
Set options(qpmR.use_cpp = FALSE) to force the R implementation.
TRUE if the compiled filter will be used.
qpm_use_cpp()qpm_use_cpp()
residuals() returns the filter's one-step-ahead prediction errors
for the observed series (type = "innovation"), the same divided by
their standard deviations ("standardized", which is what the
outlier flags use), or the smoothed structural shocks
("shock"). fitted() returns the one-step-ahead predictions of the
observables, so that observed = fitted + innovation.
## S3 method for class 'qpm_filtration' residuals(object, type = c("innovation", "standardized", "shock"), ...) ## S3 method for class 'qpm_filtration' fitted(object, ...)## S3 method for class 'qpm_filtration' residuals(object, type = c("innovation", "standardized", "shock"), ...) ## S3 method for class 'qpm_filtration' fitted(object, ...)
object |
A |
type |
Which residuals to return. |
... |
Unused. |
A data frame with a period column and one column per series.
sol <- qpm_solve(qpm_template("bkl")) obs <- simulate(sol, nsim = 40, seed = 1, burn = 20) fit <- qpm_filter(sol, obs[, c("period", "pi", "i", "q")]) head(residuals(fit, "standardized")) head(fitted(fit))sol <- qpm_solve(qpm_template("bkl")) obs <- simulate(sol, nsim = 40, seed = 1, burn = 20) fit <- qpm_filter(sol, obs[, c("period", "pi", "i", "q")]) head(residuals(fit, "standardized")) head(fitted(fit))
Print a forecast's balance-of-risks assessment
risk_log(fc)risk_log(fc)
fc |
A |
The risk table, invisibly.
sol <- qpm_solve(qpm_template("bkl")) fc <- qpm_risk(qpm_forecast(sol, horizon = 8), pi = 0.4, rationale = "upside energy risk") risk_log(fc)sol <- qpm_solve(qpm_template("bkl")) fc <- qpm_risk(qpm_forecast(sol, horizon = 8), pi = 0.4, rationale = "upside energy risk") risk_log(fc)
A round store is a plain directory: each round lives in
store/<slug>/ as a self-contained round.rds plus human-readable
sidecars (forecast.csv, data.csv, calibration.csv, and
judgment.csv when judgment was applied) so that a round can be
audited without R.
save_round(round, store, overwrite = FALSE) load_round(name, store = "rounds") list_rounds(store = "rounds")save_round(round, store, overwrite = FALSE) load_round(name, store = "rounds") list_rounds(store = "rounds")
round |
A |
store |
Directory of the round store. |
overwrite |
Allow replacing an existing round of the same name. |
name |
Round name (or its slug), or a direct path to a round
directory or |
save_round() returns the round directory invisibly;
load_round() returns the qpm_round; list_rounds() returns a
data frame of the store's contents.
m <- qpm_template("bkl") obs <- simulate(qpm_solve(m), nsim = 40, seed = 1, burn = 20) obs$period <- next_quarters("2016-Q1", 40) r <- qpm_round("demo", m, obs[, c("period", "pi", "i", "q")], horizon = 8) store <- file.path(tempdir(), "rounds") save_round(r, store) list_rounds(store) r2 <- load_round("demo", store)m <- qpm_template("bkl") obs <- simulate(qpm_solve(m), nsim = 40, seed = 1, burn = 20) obs$period <- next_quarters("2016-Q1", 40) r <- qpm_round("demo", m, obs[, c("period", "pi", "i", "q")], horizon = 8) store <- file.path(tempdir(), "rounds") save_round(r, store) list_rounds(store) r2 <- load_round("demo", store)
Accepts bare names or character strings.
shocks(...)shocks(...)
... |
Shock names, e.g. |
A character vector of shock names.
shocks(eps_y, eps_pi, eps_i)shocks(eps_y, eps_pi, eps_i)
Draws Gaussian shocks with the model's standard deviations and iterates
the solved transition. Returns variables in levels (steady state plus
simulated deviations). The final full state is stored as an attribute
so a simulation can be handed straight to qpm_forecast().
## S3 method for class 'qpm_solution' simulate(object, nsim = 40, seed = NULL, sigma = NULL, burn = 0, ...)## S3 method for class 'qpm_solution' simulate(object, nsim = 40, seed = NULL, sigma = NULL, burn = 0, ...)
object |
A |
nsim |
Number of quarters to simulate. |
seed |
Optional seed for reproducibility. The previous RNG state is restored on exit. |
sigma |
Optional named vector of shock standard deviations overriding the model's. |
burn |
Burn-in quarters discarded from the start. |
... |
Unused. |
A data frame of class qpm_sim (period plus one column per
declared variable), with attributes state (final expanded state,
deviations) and shocks (the drawn shocks).
sol <- qpm_solve(qpm_template("bkl")) hist <- simulate(sol, nsim = 60, seed = 42, burn = 20) head(hist)sol <- qpm_solve(qpm_template("bkl")) hist <- simulate(sol, nsim = 60, seed = 42, burn = 20) head(hist)
Exposes the exact matrices qpmR itself uses for filtering, so other estimators and filters can be built on top of a solved model. The representation (in deviations from steady state) is
with d the steady state of the observables and P1 the stationary
(Lyapunov) covariance used to initialize the filter.
state_space(solution, observables = NULL, measurement_error = 0, kappa = 1e+06)state_space(solution, observables = NULL, measurement_error = 0, kappa = 1e+06)
solution |
A |
observables |
Character vector of observed variables (a subset of the declared variables). Default: all declared variables. |
measurement_error |
Measurement-error standard deviation(s): a scalar recycled over observables, or a named vector. |
kappa |
Diffuse-prior variance scale for unit-root directions (only used when the model has unit roots). |
For stationary models P1 is the exact stationary covariance. When
the model has unit roots (random-walk trends), an approximate diffuse
initialization is used: P1 solves the Lyapunov equation for the
slightly damped transition sqrt(1 - 1/kappa) * T, which reproduces
the stationary covariance in stable directions and a variance of order
kappa in unit-root directions, with the exact cross-coupling. Exact
Durbin-Koopman diffuse recursions are on the roadmap.
A list with elements T, R, Z, d, Qc, H, P1,
vars_all, observables, diffuse, n_unit.
sol <- qpm_solve(qpm_template("bkl")) ss <- state_space(sol, observables = c("pi", "i", "q")) dim(ss$T); ss$dsol <- qpm_solve(qpm_template("bkl")) ss <- state_space(sol, observables = c("pi", "i", "q")) dim(ss$T); ss$d
Steady state of a model or solution
steady_state(x, ...)steady_state(x, ...)
x |
A |
... |
Unused. |
Named vector of steady-state values for the declared variables.
steady_state(qpm_template("bkl"))steady_state(qpm_template("bkl"))
Returns the posterior summary as a data frame — prior, mode, mean, standard deviation, credible interval, R-hat and effective sample size — so it can be used programmatically rather than only read.
## S3 method for class 'qpm_estimate' summary(object, level = 0.9, ...)## S3 method for class 'qpm_estimate' summary(object, level = 0.9, ...)
object |
A |
level |
Credible level for the interval. |
... |
Unused. |
A data frame, one row per estimated parameter.
m <- qpm_model(variables = vars(x = "x"), shocks = shocks(e), equations = eqs(x ~ rho * x[-1] + e), params = list(rho = 0.5)) obs <- simulate(qpm_solve(qpm_calibrate(m, rho = 0.8)), nsim = 120, seed = 1) est <- qpm_estimate(m, obs, priors(rho = beta(0.5, 0.2)), iter = 300, chains = 1, seed = 2, verbose = FALSE) summary(est)m <- qpm_model(variables = vars(x = "x"), shocks = shocks(e), equations = eqs(x ~ rho * x[-1] + e), params = list(rho = 0.5)) obs <- simulate(qpm_solve(qpm_calibrate(m, rho = 0.8)), nsim = 120, seed = 1) est <- qpm_estimate(m, obs, priors(rho = beta(0.5, 0.2)), iter = 300, chains = 1, seed = 2, verbose = FALSE) summary(est)
Summarise a filtration
## S3 method for class 'qpm_filtration' summary(object, ...)## S3 method for class 'qpm_filtration' summary(object, ...)
object |
A |
... |
Unused. |
A data frame with, for each model variable, whether it was observed and the mean, standard deviation and range of its smoothed path, plus the mean estimation standard error for latent states.
sol <- qpm_solve(qpm_template("bkl")) obs <- simulate(sol, nsim = 40, seed = 1, burn = 20) summary(qpm_filter(sol, obs[, c("period", "pi", "i", "q")]))sol <- qpm_solve(qpm_template("bkl")) obs <- simulate(sol, nsim = 40, seed = 1, burn = 20) summary(qpm_filter(sol, obs[, c("period", "pi", "i", "q")]))
Used inside vars() to attach documentation to a variable. Units are
stored and displayed; unit algebra (automatic consistency checking
inside equations) is on the roadmap.
var(label = "", unit = "", ...)var(label = "", unit = "", ...)
label |
Human-readable description, e.g. |
unit |
Unit of measurement, e.g. |
... |
Passed to |
When given numeric data instead of a character label, var()
delegates to stats::var(), so attaching qpmR does not break variance
computations.
An object of class qpm_var, or the result of stats::var().
vars(y_gap = var("Output gap", unit = "pp")) var(rnorm(10)) # still the sample variancevars(y_gap = var("Output gap", unit = "pp")) var(rnorm(10)) # still the sample variance
Declare the endogenous variables of a model
vars(...)vars(...)
... |
Named arguments, one per variable. Each value is a |
A data frame with columns name, label, unit.
vars( y_gap = var("Output gap", unit = "pp"), pi = "CPI inflation, QoQ annualised" )vars( y_gap = var("Output gap", unit = "pp"), pi = "CPI inflation, QoQ annualised" )
Posterior covariance and credible intervals
## S3 method for class 'qpm_estimate' vcov(object, ...) ## S3 method for class 'qpm_estimate' confint(object, parm = NULL, level = 0.9, ...)## S3 method for class 'qpm_estimate' vcov(object, ...) ## S3 method for class 'qpm_estimate' confint(object, parm = NULL, level = 0.9, ...)
object |
A |
... |
Unused. |
parm |
Parameters to report; default all. |
level |
Credible level. |
vcov() returns the posterior covariance matrix of the
estimated parameters; confint() returns equal-tailed posterior
credible intervals (posterior quantiles, not asymptotic intervals).
m <- qpm_model(variables = vars(x = "x"), shocks = shocks(e), equations = eqs(x ~ rho * x[-1] + e), params = list(rho = 0.5)) obs <- simulate(qpm_solve(qpm_calibrate(m, rho = 0.8)), nsim = 120, seed = 1) est <- qpm_estimate(m, obs, priors(rho = beta(0.5, 0.2)), iter = 300, chains = 1, seed = 2, verbose = FALSE) vcov(est) confint(est)m <- qpm_model(variables = vars(x = "x"), shocks = shocks(e), equations = eqs(x ~ rho * x[-1] + e), params = list(rho = 0.5)) obs <- simulate(qpm_solve(qpm_calibrate(m, rho = 0.8)), nsim = 120, seed = 1) est <- qpm_estimate(m, obs, priors(rho = beta(0.5, 0.2)), iter = 300, chains = 1, seed = 2, verbose = FALSE) vcov(est) confint(est)
Re-runs a saved forecast round from its own contents — model, calibration, data vintage, conditions and judgment — and checks that the published numbers come back. This is the audit an institution needs and that a folder of scripts cannot provide: a round that no longer reproduces is a finding, not a mystery.
verify_round(round, store = "rounds", tol = 1e-06)verify_round(round, store = "rounds", tol = 1e-06)
round |
A |
store |
Round store used when |
tol |
Absolute tolerance on the forecast paths. |
When the round is loaded from a store, the human-readable CSV sidecars are checked against the object too, so a hand-edited audit trail is detected.
An object of class qpm_verification: $ok, the largest
deviation, the worst variables, sidecar checks, and the qpmR
versions involved.
m <- qpm_template("bkl") obs <- simulate(qpm_solve(m), nsim = 40, seed = 1, burn = 20) obs$period <- next_quarters("2016-Q1", 40) r <- qpm_round("demo", m, obs[, c("period", "pi", "i", "q")], horizon = 8) verify_round(r)m <- qpm_template("bkl") obs <- simulate(qpm_solve(m), nsim = 40, seed = 1, burn = 20) obs$period <- next_quarters("2016-Q1", 40) r <- qpm_round("demo", m, obs[, c("period", "pi", "i", "q")], horizon = 8) verify_round(r)
Writes the model as Dynare source so that anyone can re-solve it in Dynare and compare. This is deliberate: qpmR's solver is pinned in its own tests to analytic solutions, but "does it agree with Dynare?" is the first question a modeller asks, and the answer should be cheap to check rather than something to take on trust.
write_dynare(model, file = NULL, irf = 40, order = 1, extra = NULL)write_dynare(model, file = NULL, irf = 40, order = 1, extra = NULL)
model |
A |
file |
Output path. The default |
irf |
Horizon for |
order |
Approximation order passed to |
extra |
Character vector of extra Dynare statements appended
after |
The original equations are exported, not qpmR's internal first-order
system: Dynare creates its own auxiliary variables for long leads and
lags, so the two implementations agree only if both handle them
correctly. Expectation wrappers are dropped (E(pi[+1]) becomes
pi(+1)), since leads in Dynare are already model-consistent
expectations.
The file path invisibly, or the Dynare source as a character
vector when file is NULL.
src <- write_dynare(qpm_template("bkl"), file = NULL) cat(head(src, 15), sep = "\n")src <- write_dynare(qpm_template("bkl"), file = NULL) cat(head(src, 15), sep = "\n")