Stochastic equations
\[ \begin{aligned} \text{Consumption}\quad C_t &= \beta_1 + \beta_{1,2} C_{t-1} + \gamma_{1,7} Y_t + \varepsilon_{1t},\\ \text{Investment}\quad I_t &= \beta_2 + \beta_{2,3} I_{t-1} + \varepsilon_{2t},\\ \text{Exports}\quad X_t &= \beta_3 + \beta_{3,4} X_{t-1} + \beta_{3,7} WGDP_t + \varepsilon_{3t},\\ \text{Imports}\quad M_t &= \beta_4 + \gamma_{4,8} D_t + \beta_{4,5} I_{t-1} + \varepsilon_{4t},\\ \text{Prices}\quad P_t &= \beta_5 + \beta_{5,9} ER_t + \beta_{5,10} OP_t + \beta_{5,6} P_{t-1} + \varepsilon_{5t},\\ \text{Interest Rate}\quad R_t &= \beta_6 + \beta_{6,6} P_{t-1} + \beta_{6,8} R_{GE,t} + \gamma_{6,5} P_t + \varepsilon_{6t}. \end{aligned} \]
Equilibrium condition
\[ \begin{aligned} \text{Equilibrium Ouput}\quad GDP_t &= 0.6*C_t + 0.6*D_t + 0.5*X_t - 0.4*I_t \end{aligned} \]
Identities
\[ \begin{aligned} \text{Domestic Demand}\quad D_t &= 0.6*C_t + 0.4*I_t \end{aligned} \]
Coefficients
- The \(\beta\)’s and \(\gamma\)’s are structural coefficients to
be estimated.
Endogenous variables
Exogenous variables
Other variables
To estimate and forecast this model, we need to specify the equations, the exogenous variables, and the date ranges for estimation and forecasting.
equations <- "consumption ~ gdp + consumption.L(1),
investment ~ investment.L(1),
exports ~ world_gdp + exports.L(1),
imports ~ domestic_demand + imports.L(1),
prices ~ exchange_rate + oil_price + prices.L(1),
interest_rate ~ prices + interest_rate_germany + prices.L(1),
gdp == 0.6*consumption + 0.6*domestic_demand + 0.5*exports - 0.4*imports,
domestic_demand == 0.6*consumption + 0.4*investment"
exogenous_variables <- c("world_gdp", "interest_rate_germany", "exchange_rate", "oil_price")The equations string contains the model equations, while
exogenous_variables lists the exogenous variables. The
.L(1) notation indicates a lag of one period. For further
details on the syntax, please refer to the equation syntax documentation.
The system_of_equations function is used to create a
system of equations object, which can then be used for estimation and
forecasting.
library(koma)
sys_eq <- system_of_equations(
equations = equations,
exogenous_variables = exogenous_variables
)
print(sys_eq)
##
## ── System of Equations ─────────────────────────────────────────────────────────
## consumption ~ constant + gdp + consumption.L(1)
## investment ~ constant + investment.L(1)
## exports ~ constant + world_gdp + exports.L(1)
## imports ~ constant + domestic_demand + imports.L(1)
## prices ~ constant + exchange_rate + oil_price + prices.L(1)
## interest_rate ~ constant + prices + interest_rate_germany + prices.L(1)
## gdp == 0.6 * consumption + 0.6 * domestic_demand + 0.5 * exports-0.4 * imports
## domestic_demand == 0.6 * consumption + 0.4 * investmentNow, we need to specify the dates for estimation and forecasting.
We will use the small_open_economy dataset, which
contains the required variables for the model. This dataset is a list of
standard time series (ts) objects, where each object
corresponds to a variable in the model. You can explore the dataset by
using the ?small_open_economy command to view its
documentation and structure.
To prepare the dataset for use with the koma package, we
need to convert the ts objects to the ets
format using the as_ets function.
The model is estimated in growth rates. The
small_open_economy dataset is in levels. We use the
ets to tell the model that the data is in levels and how to
convert it to growth rates. ets is generalized to accept
any additional attributes. Below we use this feature to save the series’
value_type (real or nominal).
series <- names(small_open_economy)
series <- series[!series %in% c("interest_rate", "interest_rate_germany")]
ts_data <- lapply(series, function(x) {
as_ets(small_open_economy[[x]],
series_type = "level", method = "diff_log"
)
})
names(ts_data) <- series
ts_data$interest_rate <- as_ets(small_open_economy$interest_rate,
series_type = "rate", method = "none"
)
ts_data$interest_rate_germany <- as_ets(small_open_economy$interest_rate_germany,
series_type = "rate", method = "none"
)With the data prepared, we can now proceed to estimate the model. Before forecasting, however, we need to truncate the endogenous variables to ensure they only include data up to the forecast start date.
The names of the list elements should match the variable names used
in the equations. Now, we can estimate the model using the
estimate function.
estimates <- estimate(
sys_eq,
ts_data = ts_data,
dates = dates
)
##
## ── Gibbs Sampler Settings ──────────────────────────────────────────────────────
## ── System Wide Settings ──
## • Number of draws (`ndraws`): 2000
## • Burn-in ratio (`burnin_ratio`): 0.5
## • Burn-in (`burnin`): 1000
## • Store frequency (`nstore`): 1
## • Number of saved draws (`nsave`): 1000
## • Tau (`tau`): 1.1
##
##
## ── Estimation ──────────────────────────────────────────────────────────────────
##
## ── ⚠ MCMC Acceptance Probability Warnings ──────────────────────────────────────
## • consumption: 61.1%
##
## ℹ Some acceptance probabilities are outside the recommended range (20%-60%).
## Consider revising the equations, tuning each equation's tau, or adjusting your priors.
print(estimates)
##
## ── Estimates ───────────────────────────────────────────────────────────────────
## consumption ~ 0.36 - 0.02 * gdp + 0.11 * consumption.L(1)
## investment ~ 0.45 + 0.22 * investment.L(1)
## exports ~ - 0.29 + 3.18 * world_gdp - 0.29 * exports.L(1)
## imports ~ - 0.04 + 2.59 * domestic_demand - 0.1 * imports.L(1)
## prices ~ 0.03 + 0.03 * exchange_rate + 0.01 * oil_price + 0.56 * prices.L(1)
## interest_rate ~ - 0.51 + 0.52 * prices + 0.57 * interest_rate_germany - 0.22 * prices.L(1)
## gdp == 0.6 * consumption + 0.6 * domestic_demand + 0.5 * exports - 0.4 * imports
## domestic_demand == 0.6 * consumption + 0.4 * investmentsummary(estimates)
##
## ===============================================================================================================
## consumption investment exports imports prices interest_rate
## ---------------------------------------------------------------------------------------------------------------
## constant 0.36 0.45 -0.29 -0.04 0.03 -0.51
## [ 0.27; 0.45] [0.19; 0.72] [-0.88; 0.29] [-0.52; 0.48] [0.00; 0.06] [-0.63; -0.39]
## consumption.L(1) 0.11
## [-0.08; 0.30]
## gdp -0.02
## [-0.11; 0.07]
## investment.L(1) 0.22
## [0.03; 0.40]
## exports.L(1) -0.29
## [-0.46; -0.13]
## world_gdp 3.18
## [ 2.33; 4.05]
## imports.L(1) -0.10
## [-0.27; 0.06]
## domestic_demand 2.59
## [ 1.49; 3.51]
## prices.L(1) 0.56 -0.22
## [0.46; 0.66] [-0.65; 0.20]
## exchange_rate 0.03
## [0.02; 0.05]
## oil_price 0.01
## [0.01; 0.01]
## interest_rate_germany 0.57
## [ 0.52; 0.62]
## prices 0.52
## [ 0.07; 1.03]
## ===============================================================================================================
## Posterior mean (90% credible interval: [5.0%, 95.0%])
## Estimation period: 1996 Q1 - 2019 Q4summary(estimates, variables = "investment")
##
## =============================
## investment
## -----------------------------
## constant 0.45
## [0.19; 0.72]
## investment.L(1) 0.22
## [0.03; 0.40]
## =============================
## Posterior mean (90% credible interval: [5.0%, 95.0%])
## Estimation period: 1996 Q1 - 2019 Q4To see how to run the estimation in parallel, refer to the parallel estimation vignette.
You can visualize MCMC diagnostics for coefficient draws using
trace_plot().
if (requireNamespace("ggplot2", quietly = TRUE)) {
trace_plot(estimates, variables = c("consumption", "investment"))
}You can also inspect running means (cumulative averages) with
running_mean_plot().
if (requireNamespace("ggplot2", quietly = TRUE)) {
running_mean_plot(estimates, variables = c("consumption", "investment"))
}You can inspect autocorrelation diagnostics with
acf_plot().
To forecast the model, we can use the forecast function.
This function will generate forecasts for the specified date range.
forecasts <- forecast(
estimates,
dates = dates
)
##
## ── Forecast ────────────────────────────────────────────────────────────────────
##
## Conditional fill detected after "2019 Q4".
## Missing observations for: "consumption", "investment", "exports", "imports",
## "prices", "interest_rate", "gdp", and "domestic_demand"
## Missing values will be conditionally filled up to "2022 Q4" before forecasting.
## ✔ Forecasting completed.
print(forecasts)
## <koma_ts>
## attributes:
## series_type: list[12]
## method: list[12]
## anker: list[12]
##
## series:
## consumption investment exports imports prices interest_rate gdp
## 2023 Q1 0.3785 0.4633 1.2160 0.8704 -0.0603 0.8378 0.7344
## 2023 Q2 0.3993 0.5732 0.6695 1.1069 -0.0709 1.2625 0.4129
## 2023 Q3 0.3753 0.6284 0.8870 1.0234 0.0126 1.5708 0.5453
## 2023 Q4 0.3873 0.5874 0.3354 1.0933 -0.0132 1.6856 0.2432
## domestic_demand world_gdp interest_rate_germany exchange_rate oil_price
## 2023 Q1 0.4124 0.4827 2.3887 0.9217 -8.7520
## 2023 Q2 0.4688 0.4083 3.1460 -1.3863 -3.3787
## 2023 Q3 0.4766 0.4259 3.6393 -1.7869 9.9001
## 2023 Q4 0.4673 0.2906 3.8850 -0.7502 -3.3373rate(forecasts$mean$gdp)
## <koma_ts>
## attributes:
## series_type: chr "rate"
## method: chr "diff_log"
## anker: num [1:2] 192523 2023
##
## series:
## Qtr1 Qtr2 Qtr3 Qtr4
## 2023 0.7343906 0.4128768 0.5453014 0.2431651
level(forecasts$mean$gdp)
## <koma_ts>
## attributes:
## series_type: chr "level"
## method: chr "diff_log"
##
## series:
## Qtr1 Qtr2 Qtr3 Qtr4
## 2022 192522.6
## 2023 193941.7 194744.1 195808.9 196285.7Conditional forecasting allows you to impose restrictions on certain variables during the forecast period to fix their values. This is useful for scenario analysis and policy evaluation, as the forecasts will reflect these user-defined conditions.