R/RcppExports.R

In bvhar: Bayesian Vector Heterogeneous Autoregressive Modeling

# Generated by using Rcpp::compileAttributes() -> do not edit by hand
# Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393

#' BVAR(p) Point Estimates based on Minnesota Prior
#' 
#' Point estimates for posterior distribution
#' 
#' @param y Time series data
#' @param lag VAR order
#' @param bayes_spec BVAR Minnesota specification
#' @param include_mean Constant term
#' 
#' @noRd
estimate_bvar_mn <- function(y, lag, bayes_spec, include_mean) {
    .Call(`_bvhar_estimate_bvar_mn`, y, lag, bayes_spec, include_mean)
}

#' BVHAR Point Estimates based on Minnesota Prior
#' 
#' Point estimates for posterior distribution
#' 
#' @param y Time series data
#' @param week VHAR week order
#' @param month VHAR month order
#' @param bayes_spec BVHAR Minnesota specification
#' @param include_mean Constant term
#' 
#' @noRd
estimate_bvhar_mn <- function(y, week, month, bayes_spec, include_mean) {
    .Call(`_bvhar_estimate_bvhar_mn`, y, week, month, bayes_spec, include_mean)
}

#' @noRd
estimate_bvar_mh <- function(num_chains, num_iter, num_burn, thin, x, y, x_dummy, y_dummy, param_prior, param_init, seed_chain, display_progress, nthreads) {
    .Call(`_bvhar_estimate_bvar_mh`, num_chains, num_iter, num_burn, thin, x, y, x_dummy, y_dummy, param_prior, param_init, seed_chain, display_progress, nthreads)
}

#' BVAR(p) Point Estimates based on Nonhierarchical Matrix Normal Prior
#' 
#' Point estimates for Ghosh et al. (2018) nonhierarchical model for BVAR.
#' 
#' @param x Design matrix X0
#' @param y Response matrix Y0
#' @param U Positive definite matrix, covariance matrix corresponding to the column of the model parameter B
#' 
#' @details
#' In Ghosh et al. (2018), there are many models for BVAR such as hierarchical or non-hierarchical.
#' Among these, this function chooses the most simple non-hierarchical matrix normal prior in Section 3.1.
#' 
#' @references
#' Ghosh, S., Khare, K., & Michailidis, G. (2018). *High-Dimensional Posterior Consistency in Bayesian Vector Autoregressive Models*. Journal of the American Statistical Association, 114(526). [https://doi:10.1080/01621459.2018.1437043](https://doi:10.1080/01621459.2018.1437043)
#' 
#' @noRd
estimate_mn_flat <- function(x, y, U) {
    .Call(`_bvhar_estimate_mn_flat`, x, y, U)
}

#' @noRd
estimate_mniw <- function(num_chains, num_iter, num_burn, thin, mn_mean, mn_prec, iw_scale, iw_shape, seed_chain, display_progress, nthreads) {
    .Call(`_bvhar_estimate_mniw`, num_chains, num_iter, num_burn, thin, mn_mean, mn_prec, iw_scale, iw_shape, seed_chain, display_progress, nthreads)
}

#' Forecasting BVAR(p)
#' 
#' @param object A `bvarmn` or `bvarflat` object
#' @param step Integer, Step to forecast
#' @param num_sim Integer, number to simulate parameters from posterior distribution
#' @details
#' n-step ahead forecasting using BVAR(p) recursively.
#' 
#' For given number of simulation (`num_sim`),
#' 
#' 1. Generate \eqn{(A^{(b)}, \Sigma_e^{(b)}) \sim MIW} (posterior)
#' 2. Recursively, \eqn{j = 1, \ldots, h} (`step`)
#'     - Point forecast: Use \eqn{\hat{A}}
#'     - Predictive distribution: Again generate \eqn{\tilde{Y}_{n + j}^{(b)} \sim A^{(b)}, \Sigma_e^{(b)} \sim MN}
#'     - tilde notation indicates simulated ones
#' 
#' @references
#' Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing. [https://doi.org/10.1007/978-3-540-27752-1](https://doi.org/10.1007/978-3-540-27752-1)
#' 
#' Litterman, R. B. (1986). *Forecasting with Bayesian Vector Autoregressions: Five Years of Experience*. Journal of Business & Economic Statistics, 4(1), 25. [https://doi:10.2307/1391384](https://doi:10.2307/1391384)
#' 
#' Bańbura, M., Giannone, D., & Reichlin, L. (2010). *Large Bayesian vector auto regressions*. Journal of Applied Econometrics, 25(1). [https://doi:10.1002/jae.1137](https://doi:10.1002/jae.1137)
#' 
#' Ghosh, S., Khare, K., & Michailidis, G. (2018). *High-Dimensional Posterior Consistency in Bayesian Vector Autoregressive Models*. Journal of the American Statistical Association, 114(526). [https://doi:10.1080/01621459.2018.1437043](https://doi:10.1080/01621459.2018.1437043)
#' 
#' Karlsson, S. (2013). *Chapter 15 Forecasting with Bayesian Vector Autoregression*. Handbook of Economic Forecasting, 2, 791-897. doi:[10.1016/b978-0-444-62731-5.00015-4](https://doi.org/10.1016/B978-0-444-62731-5.00015-4)
#' 
#' @noRd
forecast_bvar <- function(object, step, num_sim, seed) {
    .Call(`_bvhar_forecast_bvar`, object, step, num_sim, seed)
}

#' Forecasting Bayesian VHAR
#' 
#' @param object A `bvharmn` object
#' @param step Integer, Step to forecast
#' @param num_sim Integer, number to simulate parameters from posterior distribution
#' @details
#' n-step ahead forecasting using VHAR recursively.
#' 
#' For given number of simulation (`num_sim`),
#' 
#' 1. Generate \eqn{(\Phi^{(b)}, \Sigma_e^{(b)}) \sim MIW} (posterior)
#' 2. Recursively, \eqn{j = 1, \ldots, h} (`step`)
#'     - Point forecast: Use \eqn{\hat\Phi}
#'     - Predictive distribution: Again generate \eqn{\tilde{Y}_{n + j}^{(b)} \sim \Phi^{(b)}, \Sigma_e^{(b)} \sim MN}
#'     - tilde notation indicates simulated ones
#' 
#' @references Kim, Y. G., and Baek, C. (2024). *Bayesian vector heterogeneous autoregressive modeling*. Journal of Statistical Computation and Simulation, 94(6), 1139-1157.
#' @noRd
forecast_bvharmn <- function(object, step, num_sim, seed) {
    .Call(`_bvhar_forecast_bvharmn`, object, step, num_sim, seed)
}

#' Out-of-Sample Forecasting of BVAR based on Rolling Window
#' 
#' This function conducts an rolling window forecasting of BVAR with Minnesota prior.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param lag BVAR order
#' @param bayes_spec List, BVAR specification
#' @param include_mean Add constant term
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' 
#' @noRd
roll_bvar <- function(y, lag, bayes_spec, include_mean, step, y_test, seed_forecast, nthreads) {
    .Call(`_bvhar_roll_bvar`, y, lag, bayes_spec, include_mean, step, y_test, seed_forecast, nthreads)
}

#' Out-of-Sample Forecasting of BVAR based on Rolling Window
#' 
#' This function conducts an rolling window forecasting of BVAR with Flat prior.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param lag BVAR order
#' @param bayes_spec List, BVAR specification
#' @param include_mean Add constant term
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' 
#' @noRd
roll_bvarflat <- function(y, lag, U, include_mean, step, y_test, seed_forecast, nthreads) {
    .Call(`_bvhar_roll_bvarflat`, y, lag, U, include_mean, step, y_test, seed_forecast, nthreads)
}

#' Out-of-Sample Forecasting of BVHAR based on Rolling Window
#' 
#' This function conducts an rolling window forecasting of BVHAR with Minnesota prior.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param har `r lifecycle::badge("experimental")` Numeric vector for weekly and monthly order.
#' @param bayes_spec List, BVHAR specification
#' @param include_mean Add constant term
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' 
#' @noRd
roll_bvhar <- function(y, week, month, bayes_spec, include_mean, step, y_test, seed_forecast, nthreads) {
    .Call(`_bvhar_roll_bvhar`, y, week, month, bayes_spec, include_mean, step, y_test, seed_forecast, nthreads)
}

#' Out-of-Sample Forecasting of BVAR based on Expanding Window
#' 
#' This function conducts an expanding window forecasting of BVAR with Minnesota prior.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param lag BVAR order
#' @param bayes_spec List, BVAR specification
#' @param include_mean Add constant term
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' 
#' @noRd
expand_bvar <- function(y, lag, bayes_spec, include_mean, step, y_test, seed_forecast, nthreads) {
    .Call(`_bvhar_expand_bvar`, y, lag, bayes_spec, include_mean, step, y_test, seed_forecast, nthreads)
}

#' Out-of-Sample Forecasting of BVAR based on Expanding Window
#' 
#' This function conducts an expanding window forecasting of BVAR with Flat prior.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param lag BVAR order
#' @param bayes_spec List, BVAR specification
#' @param include_mean Add constant term
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' 
#' @noRd
expand_bvarflat <- function(y, lag, U, include_mean, step, y_test, seed_forecast, nthreads) {
    .Call(`_bvhar_expand_bvarflat`, y, lag, U, include_mean, step, y_test, seed_forecast, nthreads)
}

#' Out-of-Sample Forecasting of BVHAR based on Expanding Window
#' 
#' This function conducts an expanding window forecasting of BVHAR with Minnesota prior.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param har `r lifecycle::badge("experimental")` Numeric vector for weekly and monthly order.
#' @param bayes_spec List, BVHAR specification
#' @param include_mean Add constant term
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' 
#' @noRd
expand_bvhar <- function(y, week, month, bayes_spec, include_mean, step, y_test, seed_forecast, nthreads) {
    .Call(`_bvhar_expand_bvhar`, y, week, month, bayes_spec, include_mean, step, y_test, seed_forecast, nthreads)
}

#' Generalized Spillover of Minnesota prior
#' 
#' @param object varlse or vharlse object.
#' @param step Step to forecast.
#' @param num_iter Number to sample MNIW distribution
#' @param num_burn Number of burn-in
#' @param thin Thinning
#' @param seed Random seed for boost library
#' 
#' @noRd
compute_mn_spillover <- function(object, step, num_iter, num_burn, thin, seed) {
    .Call(`_bvhar_compute_mn_spillover`, object, step, num_iter, num_burn, thin, seed)
}

#' Rolling-sample Total Spillover Index of BVAR
#' 
#' @param y Time series data of which columns indicate the variables
#' @param window Rolling window size
#' @param step forecast horizon for FEVD
#' @param num_iter Number to sample MNIW distribution
#' @param num_burn Number of burn-in
#' @param thin Thinning
#' @param lag BVAR order
#' @param bayes_spec BVAR specification
#' @param include_mean Add constant term
#' @param seed_chain Random seed for each window
#' @param nthreads Number of threads for openmp
#' 
#' @noRd
dynamic_bvar_spillover <- function(y, window, step, num_iter, num_burn, thin, lag, bayes_spec, include_mean, seed_chain, nthreads) {
    .Call(`_bvhar_dynamic_bvar_spillover`, y, window, step, num_iter, num_burn, thin, lag, bayes_spec, include_mean, seed_chain, nthreads)
}

#' Rolling-sample Total Spillover Index of BVHAR
#' 
#' @param y Time series data of which columns indicate the variables
#' @param window Rolling window size
#' @param step forecast horizon for FEVD
#' @param num_iter Number to sample MNIW distribution
#' @param num_burn Number of burn-in
#' @param thin Thinning
#' @param week Week order
#' @param month Month order
#' @param bayes_spec BVHAR specification
#' @param include_mean Add constant term
#' @param seed_chain Random seed for each window
#' @param nthreads Number of threads for openmp
#' 
#' @noRd
dynamic_bvhar_spillover <- function(y, window, step, num_iter, num_burn, thin, week, month, bayes_spec, include_mean, seed_chain, nthreads) {
    .Call(`_bvhar_dynamic_bvhar_spillover`, y, window, step, num_iter, num_burn, thin, week, month, bayes_spec, include_mean, seed_chain, nthreads)
}

#' Compute VAR(p) Coefficient Matrices and Fitted Values
#' 
#' This function fits VAR(p) given response and design matrices of multivariate time series.
#' 
#' @param x Design matrix X0
#' @param y Response matrix Y0
#' @param method Method to solve linear equation system. 1: normal equation, 2: cholesky, 3: HouseholderQR.
#' @details
#' Given Y0 and Y0, the function estimate least squares
#' Y0 = X0 A + Z
#' 
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing. doi:[10.1007/978-3-540-27752-1](https://doi.org/10.1007/978-3-540-27752-1)
#' @noRd
estimate_var <- function(y, lag, include_mean, method) {
    .Call(`_bvhar_estimate_var`, y, lag, include_mean, method)
}

#' Compute Vector HAR Coefficient Matrices and Fitted Values
#' 
#' This function fits VHAR given response and design matrices of multivariate time series.
#' 
#' @param x Design matrix X0
#' @param y Response matrix Y0
#' @param week Integer, order for weekly term
#' @param month Integer, order for monthly term
#' @param include_mean bool, Add constant term (Default: `true`) or not (`false`)
#' @param method Method to solve linear equation system. 1: normal equation, 2: cholesky, 3: HouseholderQR.
#' @details
#' Given Y0 and Y0, the function estimate least squares
#' \deqn{Y_0 = X_1 \Phi + Z}
#' 
#' @references
#' Baek, C. and Park, M. (2021). *Sparse vector heterogeneous autoregressive modeling for realized volatility*. J. Korean Stat. Soc. 50, 495-510. doi:[10.1007/s42952-020-00090-5](https://doi.org/10.1007/s42952-020-00090-5)
#' 
#' Corsi, F. (2008). *A Simple Approximate Long-Memory Model of Realized Volatility*. Journal of Financial Econometrics, 7(2), 174-196. doi:[10.1093/jjfinec/nbp001](https://doi.org/10.1093/jjfinec/nbp001)
#' @noRd
estimate_har <- function(y, week, month, include_mean, method) {
    .Call(`_bvhar_estimate_har`, y, week, month, include_mean, method)
}

#' Covariance Estimate for Residual Covariance Matrix
#' 
#' Compute ubiased estimator for residual covariance.
#' 
#' @param z Matrix, residual
#' @param num_design Integer, Number of sample used (s = n - p)
#' @param dim_design Ingeger, Number of parameter for each dimension (k = mp + 1)
#' @details
#' See pp75 Lütkepohl (2007).
#' 
#' * s = n - p: sample used (`num_design`)
#' * k = mp + 1 (m: dimension, p: VAR lag): number of parameter for each dimension (`dim_design`)
#' 
#' Then an unbiased estimator for \eqn{\Sigma_e} is
#' 
#' \deqn{\hat{\Sigma}_e = \frac{1}{s - k} (Y_0 - \hat{A} X_0)^T (Y_0 - \hat{A} X_0)}
#' 
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing.
#' @noRd
compute_cov <- function(z, num_design, dim_design) {
    .Call(`_bvhar_compute_cov`, z, num_design, dim_design)
}

#' Statistic for VAR
#' 
#' Compute partial t-statistics for inference in VAR model.
#' 
#' @param object A `varlse` object
#' @details
#' Partial t-statistic for H0: aij = 0
#' 
#' * For each variable (e.g. 1st variable)
#' * Standard error =  (1st) diagonal element of \eqn{\Sigma_e} estimator x diagonal elements of \eqn{(X_0^T X_0)^(-1)}
#' 
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing. doi:[10.1007/978-3-540-27752-1](https://doi.org/10.1007/978-3-540-27752-1)
#' @noRd
infer_var <- function(object) {
    .Call(`_bvhar_infer_var`, object)
}

#' Statistic for VHAR
#' 
#' Compute partial t-statistics for inference in VHAR model.
#' 
#' @param object A `vharlse` object
#' @details
#' Partial t-statistic for H0: \eqn{\phi_{ij} = 0}
#' 
#' * For each variable (e.g. 1st variable)
#' * Standard error =  (1st) diagonal element of \eqn{\Sigma_e} estimator x diagonal elements of \eqn{(X_1^T X_1)^(-1)}
#' @noRd
infer_vhar <- function(object) {
    .Call(`_bvhar_infer_vhar`, object)
}

#' Forecasting Vector Autoregression
#' 
#' @param object A `varlse` object
#' @param step Integer, Step to forecast
#' @details
#' n-step ahead forecasting using VAR(p) recursively, based on pp35 of Lütkepohl (2007).
#' 
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing. [https://doi.org/10.1007/978-3-540-27752-1](https://doi.org/10.1007/978-3-540-27752-1)
#' @noRd
forecast_var <- function(object, step) {
    .Call(`_bvhar_forecast_var`, object, step)
}

#' Forecasting Vector HAR
#' 
#' @param object A `vharlse` object
#' @param step Integer, Step to forecast
#' @details
#' n-step ahead forecasting using VHAR recursively.
#' 
#' @noRd
forecast_vhar <- function(object, step) {
    .Call(`_bvhar_forecast_vhar`, object, step)
}

#' Out-of-Sample Forecasting of VAR based on Rolling Window
#' 
#' This function conducts an rolling window forecasting of VAR.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param lag VAR order
#' @param include_mean Add constant term
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' @param method Method to solve linear equation system. 1: normal equation, 2: cholesky, 3: HouseholderQR.
#' @param nthreads Number of threads for openmp
#' 
#' @noRd
roll_var <- function(y, lag, include_mean, step, y_test, method, nthreads) {
    .Call(`_bvhar_roll_var`, y, lag, include_mean, step, y_test, method, nthreads)
}

#' Out-of-Sample Forecasting of VHAR based on Rolling Window
#' 
#' This function conducts an rolling window forecasting of VHAR.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param week Integer, order for weekly term
#' @param month Integer, order for monthly term
#' @param include_mean Add constant term
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' @param method Method to solve linear equation system. 1: normal equation, 2: cholesky, 3: HouseholderQR.
#' @param nthreads Number of threads for openmp
#' 
#' @noRd
roll_vhar <- function(y, week, month, include_mean, step, y_test, method, nthreads) {
    .Call(`_bvhar_roll_vhar`, y, week, month, include_mean, step, y_test, method, nthreads)
}

#' Out-of-Sample Forecasting of VAR based on Expanding Window
#' 
#' This function conducts an expanding window forecasting of VAR.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param lag VAR order
#' @param include_mean Add constant term
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' @param method Method to solve linear equation system. 1: normal equation, 2: cholesky, 3: HouseholderQR.
#' @param nthreads Number of threads for openmp
#' 
#' @noRd
expand_var <- function(y, lag, include_mean, step, y_test, method, nthreads) {
    .Call(`_bvhar_expand_var`, y, lag, include_mean, step, y_test, method, nthreads)
}

#' Out-of-Sample Forecasting of VHAR based on Expanding Window
#' 
#' This function conducts an expanding window forecasting of VHAR.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param week Integer, order for weekly term
#' @param month Integer, order for monthly term
#' @param include_mean Add constant term
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' @param method Method to solve linear equation system. 1: normal equation, 2: cholesky, 3: HouseholderQR.
#' @param nthreads Number of threads for openmp
#' 
#' @noRd
expand_vhar <- function(y, week, month, include_mean, step, y_test, method, nthreads) {
    .Call(`_bvhar_expand_vhar`, y, week, month, include_mean, step, y_test, method, nthreads)
}

#' Generalized Spillover of VAR
#' 
#' @param object varlse or vharlse object.
#' @param step Step to forecast.
#' 
#' @noRd
compute_ols_spillover <- function(object, step) {
    .Call(`_bvhar_compute_ols_spillover`, object, step)
}

#' Rolling-sample Total Spillover Index of VAR
#' 
#' @param y Time series data of which columns indicate the variables
#' @param window Rolling window size
#' @param step forecast horizon for FEVD
#' @param lag VAR order
#' @param include_mean Add constant term
#' 
#' @noRd
dynamic_var_spillover <- function(y, window, step, lag, include_mean, method, nthreads) {
    .Call(`_bvhar_dynamic_var_spillover`, y, window, step, lag, include_mean, method, nthreads)
}

#' Rolling-sample Total Spillover Index of VHAR
#' 
#' @param y Time series data of which columns indicate the variables
#' @param window Rolling window size
#' @param step forecast horizon for FEVD
#' @param har VHAR order
#' @param include_mean Add constant term
#' 
#' @noRd
dynamic_vhar_spillover <- function(y, window, step, week, month, include_mean, method, nthreads) {
    .Call(`_bvhar_dynamic_vhar_spillover`, y, window, step, week, month, include_mean, method, nthreads)
}

#' VAR with Shrinkage Priors
#' 
#' This function generates parameters \eqn{\beta, a, \sigma_{h,i}^2, h_{0,i}} and log-volatilities \eqn{h_{i,1}, \ldots, h_{i, n}}.
#' 
#' @param num_chain Number of MCMC chains
#' @param num_iter Number of iteration for MCMC
#' @param num_burn Number of burn-in (warm-up) for MCMC
#' @param thin Thinning
#' @param x Design matrix X0
#' @param y Response matrix Y0
#' @param param_reg Regression specification list
#' @param param_prior Prior specification list
#' @param param_intercept Intercept specification list
#' @param param_init Initialization specification list
#' @param grp_id Unique group id
#' @param grp_mat Group matrix
#' @param include_mean Constant term
#' @param seed_chain Seed for each chain
#' @param display_progress Progress bar
#' @param nthreads Number of threads for openmp
#' 
#' @noRd
estimate_sur <- function(num_chains, num_iter, num_burn, thin, x, y, param_reg, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, seed_chain, display_progress, nthreads) {
    .Call(`_bvhar_estimate_sur`, num_chains, num_iter, num_burn, thin, x, y, param_reg, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, seed_chain, display_progress, nthreads)
}

#' Forecasting predictive density of BVAR
#' 
#' @param num_chains Number of chains
#' @param var_lag VAR order.
#' @param step Integer, Step to forecast.
#' @param response_mat Response matrix.
#' @param sv Use Innovation?
#' @param sparse Use restricted model?
#' @param level CI level to give sparsity. Valid when `prior_type` is 0.
#' @param fit_record MCMC records list
#' @param prior_type Prior type. If 0, use CI. Valid when sparse is true.
#' @param seed_chain Seed for each chain
#' @param stable Filter stable draws
#' @param include_mean Include constant term?
#' @param nthreads OpenMP number of threads
#' 
#' @noRd
forecast_bvarldlt <- function(num_chains, var_lag, step, response_mat, sparse, level, fit_record, prior_type, seed_chain, include_mean, stable, nthreads) {
    .Call(`_bvhar_forecast_bvarldlt`, num_chains, var_lag, step, response_mat, sparse, level, fit_record, prior_type, seed_chain, include_mean, stable, nthreads)
}

#' Forecasting Predictive Density of BVHAR
#' 
#' @param num_chains Number of MCMC chains
#' @param month VHAR month order.
#' @param step Integer, Step to forecast.
#' @param response_mat Response matrix.
#' @param HARtrans VHAR linear transformation matrix
#' @param sv Use Innovation?
#' @param sparse Use restricted model?
#' @param level CI level to give sparsity. Valid when `prior_type` is 0.
#' @param fit_record MCMC records list
#' @param prior_type Prior type. If 0, use CI. Valid when sparse is true.
#' @param seed_chain Seed for each chain
#' @param include_mean Include constant term?
#' @param stable Filter stable draws
#' @param nthreads OpenMP number of threads 
#'
#' @noRd
forecast_bvharldlt <- function(num_chains, month, step, response_mat, HARtrans, sparse, level, fit_record, prior_type, seed_chain, include_mean, stable, nthreads) {
    .Call(`_bvhar_forecast_bvharldlt`, num_chains, month, step, response_mat, HARtrans, sparse, level, fit_record, prior_type, seed_chain, include_mean, stable, nthreads)
}

#' Out-of-Sample Forecasting of VAR-SV based on Rolling Window
#' 
#' This function conducts an rolling window forecasting of BVAR-SV.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param lag VAR order
#' @param num_chains Number of MCMC chains
#' @param num_iter Number of iteration for MCMC
#' @param num_burn Number of burn-in (warm-up) for MCMC
#' @param thinning Thinning
#' @param param_reg SV specification list
#' @param param_prior Prior specification list
#' @param param_intercept Intercept specification list
#' @param param_init Initialization specification list
#' @param get_lpl Compute LPL
#' @param seed_chain Seed for each window and chain in the form of matrix
#' @param seed_forecast Seed for each window forecast
#' @param nthreads Number of threads for openmp
#' @param grp_id Unique group id
#' @param grp_mat Group matrix
#' @param include_mean Constant term
#' @param stable Filter stable draws
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' @param nthreads Number of threads
#' 
#' @noRd
roll_bvarldlt <- function(y, lag, num_chains, num_iter, num_burn, thinning, sparse, level, fit_record, param_reg, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, stable, step, y_test, get_lpl, seed_chain, seed_forecast, display_progress, nthreads) {
    .Call(`_bvhar_roll_bvarldlt`, y, lag, num_chains, num_iter, num_burn, thinning, sparse, level, fit_record, param_reg, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, stable, step, y_test, get_lpl, seed_chain, seed_forecast, display_progress, nthreads)
}

#' Out-of-Sample Forecasting of VAR-SV based on Rolling Window
#' 
#' This function conducts an rolling window forecasting of BVAR-SV.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param lag VAR order
#' @param num_chains Number of MCMC chains
#' @param num_iter Number of iteration for MCMC
#' @param num_burn Number of burn-in (warm-up) for MCMC
#' @param thinning Thinning
#' @param param_sv SV specification list
#' @param param_prior Prior specification list
#' @param param_intercept Intercept specification list
#' @param param_init Initialization specification list
#' @param get_lpl Compute LPL
#' @param seed_chain Seed for each window and chain in the form of matrix
#' @param seed_forecast Seed for each window forecast
#' @param nthreads Number of threads for openmp
#' @param grp_id Unique group id
#' @param grp_mat Group matrix
#' @param include_mean Constant term
#' @param stable Filter stable draws
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' @param nthreads Number of threads
#' 
#' @noRd
roll_bvharldlt <- function(y, week, month, num_chains, num_iter, num_burn, thinning, sparse, level, fit_record, param_reg, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, stable, step, y_test, get_lpl, seed_chain, seed_forecast, display_progress, nthreads) {
    .Call(`_bvhar_roll_bvharldlt`, y, week, month, num_chains, num_iter, num_burn, thinning, sparse, level, fit_record, param_reg, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, stable, step, y_test, get_lpl, seed_chain, seed_forecast, display_progress, nthreads)
}

#' Forecasting predictive density of VAR-SV
#' 
#' @param num_chains Number of chains
#' @param var_lag VAR order.
#' @param step Integer, Step to forecast.
#' @param response_mat Response matrix.
#' @param sv Use Innovation?
#' @param sparse Use restricted model?
#' @param level CI level to give sparsity. Valid when `prior_type` is 0.
#' @param fit_record MCMC records list
#' @param prior_type Prior type. If 0, use CI. Valid when sparse is true.
#' @param seed_chain Seed for each chain
#' @param stable Filter stable draws
#' @param include_mean Include constant term?
#' @param nthreads OpenMP number of threads
#' 
#' @noRd
forecast_bvarsv <- function(num_chains, var_lag, step, response_mat, sv, sparse, level, fit_record, prior_type, seed_chain, include_mean, stable, nthreads) {
    .Call(`_bvhar_forecast_bvarsv`, num_chains, var_lag, step, response_mat, sv, sparse, level, fit_record, prior_type, seed_chain, include_mean, stable, nthreads)
}

#' Forecasting Predictive Density of VHAR-SV
#' 
#' @param num_chains Number of MCMC chains
#' @param month VHAR month order.
#' @param step Integer, Step to forecast.
#' @param response_mat Response matrix.
#' @param HARtrans VHAR linear transformation matrix
#' @param sv Use Innovation?
#' @param sparse Use restricted model?
#' @param level CI level to give sparsity. Valid when `prior_type` is 0.
#' @param fit_record MCMC records list
#' @param prior_type Prior type. If 0, use CI. Valid when sparse is true.
#' @param seed_chain Seed for each chain
#' @param include_mean Include constant term?
#' @param stable Filter stable draws
#' @param nthreads OpenMP number of threads 
#'
#' @noRd
forecast_bvharsv <- function(num_chains, month, step, response_mat, HARtrans, sv, sparse, level, fit_record, prior_type, seed_chain, include_mean, stable, nthreads) {
    .Call(`_bvhar_forecast_bvharsv`, num_chains, month, step, response_mat, HARtrans, sv, sparse, level, fit_record, prior_type, seed_chain, include_mean, stable, nthreads)
}

#' Out-of-Sample Forecasting of VAR-SV based on Rolling Window
#' 
#' This function conducts an rolling window forecasting of BVAR-SV.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param lag VAR order
#' @param num_chains Number of MCMC chains
#' @param num_iter Number of iteration for MCMC
#' @param num_burn Number of burn-in (warm-up) for MCMC
#' @param thinning Thinning
#' @param param_sv SV specification list
#' @param param_prior Prior specification list
#' @param param_intercept Intercept specification list
#' @param param_init Initialization specification list
#' @param get_lpl Compute LPL
#' @param seed_chain Seed for each window and chain in the form of matrix
#' @param seed_forecast Seed for each window forecast
#' @param nthreads Number of threads for openmp
#' @param grp_id Unique group id
#' @param grp_mat Group matrix
#' @param include_mean Constant term
#' @param stable Filter stable draws
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' @param nthreads Number of threads
#' 
#' @noRd
roll_bvarsv <- function(y, lag, num_chains, num_iter, num_burn, thinning, sv, sparse, level, fit_record, param_sv, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, stable, step, y_test, get_lpl, seed_chain, seed_forecast, display_progress, nthreads) {
    .Call(`_bvhar_roll_bvarsv`, y, lag, num_chains, num_iter, num_burn, thinning, sv, sparse, level, fit_record, param_sv, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, stable, step, y_test, get_lpl, seed_chain, seed_forecast, display_progress, nthreads)
}

#' Out-of-Sample Forecasting of VAR-SV based on Rolling Window
#' 
#' This function conducts an rolling window forecasting of BVAR-SV.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param lag VAR order
#' @param num_chains Number of MCMC chains
#' @param num_iter Number of iteration for MCMC
#' @param num_burn Number of burn-in (warm-up) for MCMC
#' @param thinning Thinning
#' @param param_sv SV specification list
#' @param param_prior Prior specification list
#' @param param_intercept Intercept specification list
#' @param param_init Initialization specification list
#' @param get_lpl Compute LPL
#' @param seed_chain Seed for each window and chain in the form of matrix
#' @param seed_forecast Seed for each window forecast
#' @param nthreads Number of threads for openmp
#' @param grp_id Unique group id
#' @param grp_mat Group matrix
#' @param include_mean Constant term
#' @param stable Filter stable draws
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' @param nthreads Number of threads
#' 
#' @noRd
roll_bvharsv <- function(y, week, month, num_chains, num_iter, num_burn, thinning, sv, sparse, level, fit_record, param_sv, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, stable, step, y_test, get_lpl, seed_chain, seed_forecast, display_progress, nthreads) {
    .Call(`_bvhar_roll_bvharsv`, y, week, month, num_chains, num_iter, num_burn, thinning, sv, sparse, level, fit_record, param_sv, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, stable, step, y_test, get_lpl, seed_chain, seed_forecast, display_progress, nthreads)
}

#' Out-of-Sample Forecasting of VAR-SV based on Rolling Window
#' 
#' This function conducts an rolling window forecasting of BVAR-SV.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param lag VAR order
#' @param num_chains Number of MCMC chains
#' @param num_iter Number of iteration for MCMC
#' @param num_burn Number of burn-in (warm-up) for MCMC
#' @param thinning Thinning
#' @param param_sv SV specification list
#' @param param_prior Prior specification list
#' @param param_intercept Intercept specification list
#' @param param_init Initialization specification list
#' @param get_lpl Compute LPL
#' @param seed_chain Seed for each window and chain in the form of matrix
#' @param seed_forecast Seed for each window forecast
#' @param nthreads Number of threads for openmp
#' @param grp_id Unique group id
#' @param grp_mat Group matrix
#' @param include_mean Constant term
#' @param stable Filter stable draws
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' @param nthreads Number of threads
#' 
#' @noRd
expand_bvarldlt <- function(y, lag, num_chains, num_iter, num_burn, thinning, sparse, level, fit_record, param_reg, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, stable, step, y_test, get_lpl, seed_chain, seed_forecast, display_progress, nthreads) {
    .Call(`_bvhar_expand_bvarldlt`, y, lag, num_chains, num_iter, num_burn, thinning, sparse, level, fit_record, param_reg, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, stable, step, y_test, get_lpl, seed_chain, seed_forecast, display_progress, nthreads)
}

#' Out-of-Sample Forecasting of VAR-SV based on Rolling Window
#' 
#' This function conducts an rolling window forecasting of BVAR-SV.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param lag VAR order
#' @param num_chains Number of MCMC chains
#' @param num_iter Number of iteration for MCMC
#' @param num_burn Number of burn-in (warm-up) for MCMC
#' @param thinning Thinning
#' @param param_sv SV specification list
#' @param param_prior Prior specification list
#' @param param_intercept Intercept specification list
#' @param param_init Initialization specification list
#' @param get_lpl Compute LPL
#' @param seed_chain Seed for each window and chain in the form of matrix
#' @param seed_forecast Seed for each window forecast
#' @param nthreads Number of threads for openmp
#' @param grp_id Unique group id
#' @param grp_mat Group matrix
#' @param include_mean Constant term
#' @param stable Filter stable draws
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' @param nthreads Number of threads
#' 
#' @noRd
expand_bvharldlt <- function(y, week, month, num_chains, num_iter, num_burn, thinning, sparse, level, fit_record, param_reg, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, stable, step, y_test, get_lpl, seed_chain, seed_forecast, display_progress, nthreads) {
    .Call(`_bvhar_expand_bvharldlt`, y, week, month, num_chains, num_iter, num_burn, thinning, sparse, level, fit_record, param_reg, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, stable, step, y_test, get_lpl, seed_chain, seed_forecast, display_progress, nthreads)
}

#' Out-of-Sample Forecasting of VAR-SV based on Rolling Window
#' 
#' This function conducts an rolling window forecasting of BVAR-SV.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param lag VAR order
#' @param num_chains Number of MCMC chains
#' @param num_iter Number of iteration for MCMC
#' @param num_burn Number of burn-in (warm-up) for MCMC
#' @param thinning Thinning
#' @param param_sv SV specification list
#' @param param_prior Prior specification list
#' @param param_intercept Intercept specification list
#' @param param_init Initialization specification list
#' @param get_lpl Compute LPL
#' @param seed_chain Seed for each window and chain in the form of matrix
#' @param seed_forecast Seed for each window forecast
#' @param nthreads Number of threads for openmp
#' @param grp_id Unique group id
#' @param grp_mat Group matrix
#' @param include_mean Constant term
#' @param stable Filter stable draws
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' @param nthreads Number of threads
#' 
#' @noRd
expand_bvarsv <- function(y, lag, num_chains, num_iter, num_burn, thinning, sv, sparse, level, fit_record, param_sv, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, stable, step, y_test, get_lpl, seed_chain, seed_forecast, display_progress, nthreads) {
    .Call(`_bvhar_expand_bvarsv`, y, lag, num_chains, num_iter, num_burn, thinning, sv, sparse, level, fit_record, param_sv, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, stable, step, y_test, get_lpl, seed_chain, seed_forecast, display_progress, nthreads)
}

#' Out-of-Sample Forecasting of VAR-SV based on Rolling Window
#' 
#' This function conducts an rolling window forecasting of BVAR-SV.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param lag VAR order
#' @param num_chains Number of MCMC chains
#' @param num_iter Number of iteration for MCMC
#' @param num_burn Number of burn-in (warm-up) for MCMC
#' @param thinning Thinning
#' @param param_sv SV specification list
#' @param param_prior Prior specification list
#' @param param_intercept Intercept specification list
#' @param param_init Initialization specification list
#' @param get_lpl Compute LPL
#' @param seed_chain Seed for each window and chain in the form of matrix
#' @param seed_forecast Seed for each window forecast
#' @param nthreads Number of threads for openmp
#' @param grp_id Unique group id
#' @param grp_mat Group matrix
#' @param include_mean Constant term
#' @param stable Filter stable draws
#' @param step Integer, Step to forecast
#' @param y_test Evaluation time series data period after `y`
#' @param nthreads Number of threads
#' 
#' @noRd
expand_bvharsv <- function(y, week, month, num_chains, num_iter, num_burn, thinning, sv, sparse, level, fit_record, param_sv, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, stable, step, y_test, get_lpl, seed_chain, seed_forecast, display_progress, nthreads) {
    .Call(`_bvhar_expand_bvharsv`, y, week, month, num_chains, num_iter, num_burn, thinning, sv, sparse, level, fit_record, param_sv, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, stable, step, y_test, get_lpl, seed_chain, seed_forecast, display_progress, nthreads)
}

compute_varldlt_spillover <- function(lag, step, fit_record, sparse) {
    .Call(`_bvhar_compute_varldlt_spillover`, lag, step, fit_record, sparse)
}

compute_vharldlt_spillover <- function(week, month, step, fit_record, sparse) {
    .Call(`_bvhar_compute_vharldlt_spillover`, week, month, step, fit_record, sparse)
}

dynamic_bvarldlt_spillover <- function(y, window, step, num_chains, num_iter, num_burn, thin, sparse, lag, param_reg, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, seed_chain, nthreads) {
    .Call(`_bvhar_dynamic_bvarldlt_spillover`, y, window, step, num_chains, num_iter, num_burn, thin, sparse, lag, param_reg, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, seed_chain, nthreads)
}

dynamic_bvharldlt_spillover <- function(y, window, step, num_chains, num_iter, num_burn, thin, sparse, week, month, param_reg, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, seed_chain, nthreads) {
    .Call(`_bvhar_dynamic_bvharldlt_spillover`, y, window, step, num_chains, num_iter, num_burn, thin, sparse, week, month, param_reg, param_prior, param_intercept, param_init, prior_type, ggl, grp_id, own_id, cross_id, grp_mat, include_mean, seed_chain, nthreads)
}

#' Dynamic Total Spillover Index of BVAR-SV
#' 
#' @param lag VAR lag.
#' @param window Rolling window size
#' @param step forecast horizon for FEVD
#' @param response_mat Response matrix.
#' @param phi_record Coefficients MCMC record
#' @param h_record log volatility MCMC record
#' @param a_record Contemporaneous coefficients MCMC record
#' 
#' @noRd
dynamic_bvarsv_spillover <- function(lag, step, num_design, fit_record, sparse, include_mean, nthreads) {
    .Call(`_bvhar_dynamic_bvarsv_spillover`, lag, step, num_design, fit_record, sparse, include_mean, nthreads)
}

#' Dynamic Total Spillover Index of BVHAR-SV
#' 
#' @param month VHAR month order.
#' @param window Rolling window size
#' @param step forecast horizon for FEVD
#' @param response_mat Response matrix.
#' @param HARtrans VHAR linear transformation matrix
#' @param phi_record Coefficients MCMC record
#' @param h_record log volatility MCMC record
#' @param a_record Contemporaneous coefficients MCMC record
#' 
#' @noRd
dynamic_bvharsv_spillover <- function(week, month, step, num_design, fit_record, sparse, include_mean, nthreads) {
    .Call(`_bvhar_dynamic_bvharsv_spillover`, week, month, step, num_design, fit_record, sparse, include_mean, nthreads)
}

get_maxomp <- function() {
    .Call(`_bvhar_get_maxomp`)
}

check_omp <- function() {
    invisible(.Call(`_bvhar_check_omp`))
}

is_omp <- function() {
    .Call(`_bvhar_is_omp`)
}

#' Build Response Matrix of VAR(p)
#' 
#' This function constructs response matrix of multivariate regression model formulation of VAR(p).
#' 
#' @param y Matrix, multivariate time series data.
#' @param var_lag Integer, VAR lag.
#' @param index Integer, Starting index to extract
#' 
#' @details
#' Let s = n - p.
#' \deqn{Y_j = (y_j, y_{j + 1}, \ldots, y_{j + s - 1})^T}
#' is the s x m matrix.
#' 
#' In case of response matrix, t = p + 1 (i.e. \eqn{Y_0 = Y_{p + 1}}).
#' This function is also used when constructing design matrix.
#' 
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing. [https://doi.org/10.1007/978-3-540-27752-1](https://doi.org/10.1007/978-3-540-27752-1)
#' 
#' @noRd
build_response <- function(y, var_lag, index) {
    .Call(`_bvhar_build_response`, y, var_lag, index)
}

#' Build Design Matrix of VAR(p)
#' 
#' This function constructs design matrix of multivariate regression model formulation of VAR(p).
#' 
#' @param y Matrix, time series data
#' @param var_lag VAR lag
#' @param include_mean bool, Add constant term (Default: `true`) or not (`false`)
#' 
#' @details
#' X0 is
#' \deqn{X_0 = [Y_p, \ldots, Y_1, 1]}
#' i.e. (n - p) x (mp + 1) matrix
#' 
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing. [https://doi.org/10.1007/978-3-540-27752-1](https://doi.org/10.1007/978-3-540-27752-1)
#' 
#' @noRd
build_design <- function(y, var_lag, include_mean) {
    .Call(`_bvhar_build_design`, y, var_lag, include_mean)
}

#' Building a Linear Transformation Matrix for Vector HAR
#' 
#' This function produces a linear transformation matrix for VHAR for given dimension.
#' 
#' @param dim Integer, dimension
#' @param week Integer, order for weekly term
#' @param month Integer, order for monthly term
#' @param include_mean bool, Add constant term (Default: `true`) or not (`false`)
#' @details
#' VHAR is linearly restricted VAR(month = 22) in \eqn{Y_0 = X_0 A + Z}.
#' \deqn{Y_0 = X_1 \Phi + Z = (X_0 C_{HAR}^T) \Phi + Z}
#' This function computes above \eqn{C_{HAR}}.
#' 
#' Default VHAR model sets `week` and `month` as `5` and `22`.
#' This function can change these numbers to get linear transformation matrix.
#' 
#' @noRd
scale_har <- function(dim, week, month, include_mean) {
    .Call(`_bvhar_scale_har`, dim, week, month, include_mean)
}

#' Construct Dummy response for Minnesota Prior
#' 
#' Define dummy Y observations to add for Minnesota moments.
#' 
#' @param p Integer, VAR lag. For VHAR, put 3.
#' @param sigma Vector, standard error of each variable
#' @param lambda Double, tightness of the prior around a random walk or white noise
#' @param daily Vector, prior belief about white noise (Litterman sets 1)
#' @param weekly Vector, this was zero in the original Minnesota design
#' @param monthly Vector, this was zero in the original Minnesota design
#' @param include_mean bool, Add constant term (Default: `true`) or not (`false`)
#' 
#' @details
#' Bańbura et al. (2010) defines dummy observation and augment to the original data matrix to construct Litterman (1986) prior.
#' 
#' @references
#' Litterman, R. B. (1986). *Forecasting with Bayesian Vector Autoregressions: Five Years of Experience*. Journal of Business & Economic Statistics, 4(1), 25. [https://doi:10.2307/1391384](https://doi:10.2307/1391384)
#' 
#' Bańbura, M., Giannone, D., & Reichlin, L. (2010). *Large Bayesian vector auto regressions*. Journal of Applied Econometrics, 25(1). [https://doi:10.1002/jae.1137](https://doi:10.1002/jae.1137)
#' 
#' @noRd
build_ydummy_export <- function(p, sigma, lambda, daily, weekly, monthly, include_mean) {
    .Call(`_bvhar_build_ydummy_export`, p, sigma, lambda, daily, weekly, monthly, include_mean)
}

#' Construct Dummy design matrix for Minnesota Prior
#' 
#' Define dummy X observation to add for Minnesota moments.
#' 
#' @param lag_seq Vector, sequence to build Jp = diag(1, ... p) matrix inside Xp.
#' @param sigma Vector, standard error of each variable
#' @param lambda Double, tightness of the prior around a random walk or white noise
#' @param eps Double, very small number
#' 
#' @details
#' Bańbura et al. (2010) defines dummy observation and augment to the original data matrix to construct Litterman (1986) prior.
#' 
#' @references
#' Litterman, R. B. (1986). *Forecasting with Bayesian Vector Autoregressions: Five Years of Experience*. Journal of Business & Economic Statistics, 4(1), 25. [https://doi:10.2307/1391384](https://doi:10.2307/1391384)
#' 
#' Bańbura, M., Giannone, D., & Reichlin, L. (2010). *Large Bayesian vector auto regressions*. Journal of Applied Econometrics, 25(1). [https://doi:10.1002/jae.1137](https://doi:10.1002/jae.1137)
#' 
#' @noRd
build_xdummy_export <- function(lag_seq, lambda, sigma, eps, include_mean) {
    .Call(`_bvhar_build_xdummy_export`, lag_seq, lambda, sigma, eps, include_mean)
}

#' Parameters of Normal Inverted Wishart Prior
#' 
#' Given dummy observations, compute parameters of Normal-IW prior for Minnesota.
#' 
#' @param x_dummy Matrix, dummy observation for X0
#' @param y_dummy Matrix, dummy observation for Y0
#' 
#' @details
#' Minnesota prior give prior to parameters \eqn{B} (VAR matrices) and \eqn{\Sigma_e} (residual covariance) the following distributions
#' 
#' \deqn{B \mid \Sigma_e, Y_0 \sim MN(B_0, \Omega_0, \Sigma_e)}
#' \deqn{\Sigma_e \mid Y_0 \sim IW(S_0, \alpha_0)}
#' (MN: [matrix normal](https://en.wikipedia.org/wiki/Matrix_normal_distribution), IW: [inverse-wishart](https://en.wikipedia.org/wiki/Inverse-Wishart_distribution))
#' 
#' Bańbura et al. (2010) provides the formula how to find each matrix to match Minnesota moments.
#' 
#' @references
#' Litterman, R. B. (1986). *Forecasting with Bayesian Vector Autoregressions: Five Years of Experience*. Journal of Business & Economic Statistics, 4(1), 25. [https://doi:10.2307/1391384](https://doi:10.2307/1391384)
#' 
#' Bańbura, M., Giannone, D., & Reichlin, L. (2010). *Large Bayesian vector auto regressions*. Journal of Applied Econometrics, 25(1). [https://doi:10.1002/jae.1137](https://doi:10.1002/jae.1137)
#' 
#' @noRd
minnesota_prior <- function(x_dummy, y_dummy) {
    .Call(`_bvhar_minnesota_prior`, x_dummy, y_dummy)
}

#' Generate Multivariate Normal Random Vector
#' 
#' This function samples n x muti-dimensional normal random matrix.
#' 
#' @param num_sim Number to generate process
#' @param mu Mean vector
#' @param sig Variance matrix
#' @noRd
sim_mgaussian <- function(num_sim, mu, sig) {
    .Call(`_bvhar_sim_mgaussian`, num_sim, mu, sig)
}

#' Generate Multivariate Normal Random Vector using Cholesky Decomposition
#' 
#' This function samples n x muti-dimensional normal random matrix with using Cholesky decomposition.
#' 
#' @param num_sim Number to generate process
#' @param mu Mean vector
#' @param sig Variance matrix
#' @details
#' This function computes \eqn{\Sigma^{1/2}} by choleksy decomposition.
#' 
#' @noRd
sim_mgaussian_chol_export <- function(num_sim, mu, sig) {
    .Call(`_bvhar_sim_mgaussian_chol_export`, num_sim, mu, sig)
}

#' Generate Multivariate t Random Vector
#' 
#' This function samples n x muti-dimensional normal random matrix.
#' 
#' @param num_sim Number to generate process
#' @param df Degrees of freedom
#' @param mu Location vector
#' @param sig Scale matrix
#' @param method Method to compute \eqn{\Sigma^{1/2}}. 1: spectral decomposition, 2: Cholesky.
#' 
#' @noRd
sim_mstudent <- function(num_sim, df, mu, sig, method) {
    .Call(`_bvhar_sim_mstudent`, num_sim, df, mu, sig, method)
}

#' Generate Matrix Normal Random Matrix
#' 
#' This function samples one matrix gaussian matrix.
#' 
#' @param mat_mean Mean matrix
#' @param mat_scale_u First scale matrix
#' @param mat_scale_v Second scale matrix
#' @param u_prec If `TRUE`, use `mat_scale_u` as its inverse.
#' @details
#' Consider n x k matrix \eqn{Y_1, \ldots, Y_n \sim MN(M, U, V)} where M is n x k, U is n x n, and V is k x k.
#' 
#' 1. Lower triangular Cholesky decomposition: \eqn{U = P P^T} and \eqn{V = L L^T}
#' 2. Standard normal generation: s x m matrix \eqn{Z_i = [z_{ij} \sim N(0, 1)]} in row-wise direction.
#' 3. \eqn{Y_i = M + P Z_i L^T}
#' 
#' This function only generates one matrix, i.e. \eqn{Y_1}.
#' @return One n x k matrix following MN distribution.
#' @export
sim_matgaussian <- function(mat_mean, mat_scale_u, mat_scale_v, u_prec) {
    .Call(`_bvhar_sim_matgaussian`, mat_mean, mat_scale_u, mat_scale_v, u_prec)
}

#' Generate Inverse-Wishart Random Matrix
#' 
#' This function samples one matrix IW matrix.
#' 
#' @param mat_scale Scale matrix
#' @param shape Shape
#' @details
#' Consider \eqn{\Sigma \sim IW(\Psi, \nu)}.
#' 
#' 1. Upper triangular Bartlett decomposition: k x k matrix \eqn{Q = [q_{ij}]} upper triangular with
#'     1. \eqn{q_{ii}^2 \chi_{\nu - i + 1}^2}
#'     2. \eqn{q_{ij} \sim N(0, 1)} with i < j (upper triangular)
#' 2. Lower triangular Cholesky decomposition: \eqn{\Psi = L L^T}
#' 3. \eqn{A = L (Q^{-1})^T}
#' 4. \eqn{\Sigma = A A^T \sim IW(\Psi, \nu)}
#' @return One k x k matrix following IW distribution
#' @export
sim_iw <- function(mat_scale, shape) {
    .Call(`_bvhar_sim_iw`, mat_scale, shape)
}

#' Generate Normal-IW Random Family
#' 
#' This function samples normal inverse-wishart matrices.
#' 
#' @param num_sim Number to generate
#' @param mat_mean Mean matrix of MN
#' @param mat_scale_u First scale matrix of MN
#' @param mat_scale Scale matrix of IW
#' @param shape Shape of IW
#' @param prec If true, use mat_scale_u as its inverse
#' @noRd
sim_mniw_export <- function(num_sim, mat_mean, mat_scale_u, mat_scale, shape, prec) {
    .Call(`_bvhar_sim_mniw_export`, num_sim, mat_mean, mat_scale_u, mat_scale, shape, prec)
}

#' Generate Multivariate Time Series Process Following VAR(p)
#' 
#' This function generates multivariate time series dataset that follows VAR(p).
#' 
#' @param num_sim Number to generated process
#' @param num_burn Number of burn-in
#' @param var_coef VAR coefficient. The format should be the same as the output of [coef()] from [var_lm()]
#' @param var_lag Lag of VAR
#' @param sig_error Variance matrix of the error term. Try `diag(dim)`.
#' @param init Initial y1, ..., yp matrix to simulate VAR model. Try `matrix(0L, nrow = var_lag, ncol = dim)`.
#' @param process Process type. 1: Gaussian. 2: student-t.
#' @param mvt_df DF of MVT
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing. doi:[10.1007/978-3-540-27752-1](https://doi.org/10.1007/978-3-540-27752-1)
#' @noRd
sim_var_eigen <- function(num_sim, num_burn, var_coef, var_lag, sig_error, init, process, mvt_df) {
    .Call(`_bvhar_sim_var_eigen`, num_sim, num_burn, var_coef, var_lag, sig_error, init, process, mvt_df)
}

#' Generate Multivariate Time Series Process Following VAR(p) using Cholesky Decomposition
#' 
#' This function generates VAR(p) using Cholesky Decomposition.
#' 
#' @param num_sim Number to generated process
#' @param num_burn Number of burn-in
#' @param var_coef VAR coefficient. The format should be the same as the output of [coef()] from [var_lm()]
#' @param var_lag Lag of VAR
#' @param sig_error Variance matrix of the error term. Try `diag(dim)`.
#' @param init Initial y1, ..., yp matrix to simulate VAR model. Try `matrix(0L, nrow = var_lag, ncol = dim)`.
#' @param process Process type. 1: Gaussian. 2: student-t.
#' @param mvt_df DF of MVT
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing. doi:[10.1007/978-3-540-27752-1](https://doi.org/10.1007/978-3-540-27752-1)
#' @noRd
sim_var_chol <- function(num_sim, num_burn, var_coef, var_lag, sig_error, init, process, mvt_df) {
    .Call(`_bvhar_sim_var_chol`, num_sim, num_burn, var_coef, var_lag, sig_error, init, process, mvt_df)
}

#' Generate Multivariate Time Series Process Following VHAR
#' 
#' This function generates multivariate time series dataset that follows VHAR.
#' 
#' @param num_sim Number to generated process
#' @param num_burn Number of burn-in
#' @param vhar_coef VHAR coefficient. The format should be the same as the output of [coef.vharlse()] from [vhar_lm()]
#' @param week Order for weekly term. Try `5L` by default.
#' @param month Order for monthly term. Try `22L` by default.
#' @param sig_error Variance matrix of the error term. Try `diag(dim)`.
#' @param init Initial y1, ..., y_month matrix to simulate VHAR model. Try `matrix(0L, nrow = month, ncol = dim)`.
#' @param process Process type. 1: Gaussian. 2: student-t.
#' @param mvt_df DF of MVT
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing. doi:[10.1007/978-3-540-27752-1](https://doi.org/10.1007/978-3-540-27752-1)
#' @noRd
sim_vhar_eigen <- function(num_sim, num_burn, vhar_coef, week, month, sig_error, init, process, mvt_df) {
    .Call(`_bvhar_sim_vhar_eigen`, num_sim, num_burn, vhar_coef, week, month, sig_error, init, process, mvt_df)
}

#' Generate Multivariate Time Series Process Following VHAR using Cholesky Decomposition
#' 
#' This function generates multivariate time series dataset that follows VHAR.
#' 
#' @param num_sim Number to generated process
#' @param num_burn Number of burn-in
#' @param vhar_coef VHAR coefficient. The format should be the same as the output of [coef.vharlse()] from [vhar_lm()]
#' @param week Order for weekly term. Try `5L` by default.
#' @param month Order for monthly term. Try `22L` by default.
#' @param sig_error Variance matrix of the error term. Try `diag(dim)`.
#' @param init Initial y1, ..., y_month matrix to simulate VHAR model. Try `matrix(0L, nrow = month, ncol = dim)`.
#' @param process Process type. 1: Gaussian. 2: student-t.
#' @param mvt_df DF of MVT
#' @details
#' Let \eqn{M} be the month order, e.g. \eqn{M = 22}.
#' 
#' 1. Generate \eqn{\epsilon_1, \epsilon_n \sim N(0, \Sigma)}
#' 2. For i = 1, ... n,
#' \deqn{y_{M + i} = (y_{M + i - 1}^T, \ldots, y_i^T, 1)^T C_{HAR}^T \Phi + \epsilon_i}
#' 3. Then the output is \eqn{(y_{M + 1}, \ldots, y_{n + M})^T}
#' 
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing. doi:[10.1007/978-3-540-27752-1](https://doi.org/10.1007/978-3-540-27752-1)
#' @noRd
sim_vhar_chol <- function(num_sim, num_burn, vhar_coef, week, month, sig_error, init, process, mvt_df) {
    .Call(`_bvhar_sim_vhar_chol`, num_sim, num_burn, vhar_coef, week, month, sig_error, init, process, mvt_df)
}

#' @noRd
VARcoeftoVMA <- function(var_coef, var_lag, lag_max) {
    .Call(`_bvhar_VARcoeftoVMA`, var_coef, var_lag, lag_max)
}

#' Convert VAR to VMA(infinite)
#' 
#' Convert VAR process to infinite vector MA process
#' 
#' @param object A `varlse` object
#' @param lag_max Maximum lag for VMA
#' @details
#' Let VAR(p) be stable.
#' \deqn{Y_t = c + \sum_{j = 0} W_j Z_{t - j}}
#' For VAR coefficient \eqn{B_1, B_2, \ldots, B_p},
#' \deqn{I = (W_0 + W_1 L + W_2 L^2 + \cdots + ) (I - B_1 L - B_2 L^2 - \cdots - B_p L^p)}
#' Recursively,
#' \deqn{W_0 = I}
#' \deqn{W_1 = W_0 B_1 (W_1^T = B_1^T W_0^T)}
#' \deqn{W_2 = W_1 B_1 + W_0 B_2 (W_2^T = B_1^T W_1^T + B_2^T W_0^T)}
#' \deqn{W_j = \sum_{j = 1}^k W_{k - j} B_j (W_j^T = \sum_{j = 1}^k B_j^T W_{k - j}^T)}
#' @return VMA coefficient of k(lag-max + 1) x k dimension
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing.
#' @export
VARtoVMA <- function(object, lag_max) {
    .Call(`_bvhar_VARtoVMA`, object, lag_max)
}

#' @noRd
compute_var_mse <- function(cov_mat, var_coef, var_lag, step) {
    .Call(`_bvhar_compute_var_mse`, cov_mat, var_coef, var_lag, step)
}

#' Compute Forecast MSE Matrices
#' 
#' Compute the forecast MSE matrices using VMA coefficients
#' 
#' @param object A `varlse` object
#' @param step Integer, Step to forecast
#' @details
#' See pp38 of Lütkepohl (2007).
#' Let \eqn{\Sigma} be the covariance matrix of VAR and let \eqn{W_j} be the VMA coefficients.
#' Recursively,
#' \deqn{\Sigma_y(1) = \Sigma}
#' \deqn{\Sigma_y(2) = \Sigma + W_1 \Sigma W_1^T}
#' \deqn{\Sigma_y(3) = \Sigma_y(2) + W_2 \Sigma W_2^T}
#' 
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing. doi:[10.1007/978-3-540-27752-1](https://doi.org/10.1007/978-3-540-27752-1)
#' @noRd
compute_covmse <- function(object, step) {
    .Call(`_bvhar_compute_covmse`, object, step)
}

#' Convert VAR to Orthogonalized VMA(infinite)
#' 
#' Convert VAR process to infinite orthogonalized vector MA process
#' 
#' @param var_coef VAR coefficient matrix
#' @param var_covmat VAR covariance matrix
#' @param var_lag VAR order
#' @param lag_max Maximum lag for VMA
#' @noRd
VARcoeftoVMA_ortho <- function(var_coef, var_covmat, var_lag, lag_max) {
    .Call(`_bvhar_VARcoeftoVMA_ortho`, var_coef, var_covmat, var_lag, lag_max)
}

#' @noRd
VHARcoeftoVMA <- function(vhar_coef, HARtrans_mat, lag_max, month) {
    .Call(`_bvhar_VHARcoeftoVMA`, vhar_coef, HARtrans_mat, lag_max, month)
}

#' Convert VHAR to VMA(infinite)
#' 
#' Convert VHAR process to infinite vector MA process
#' 
#' @param object A `vharlse` object
#' @param lag_max Maximum lag for VMA
#' @details
#' Let VAR(p) be stable
#' and let VAR(p) be
#' \eqn{Y_0 = X_0 B + Z}
#' 
#' VHAR is VAR(22) with
#' \deqn{Y_0 = X_1 B + Z = ((X_0 \tilde{T}^T)) \Phi + Z}
#' 
#' Observe that
#' \deqn{B = \tilde{T}^T \Phi}
#' @return VMA coefficient of k(lag-max + 1) x k dimension
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing.
#' @export
VHARtoVMA <- function(object, lag_max) {
    .Call(`_bvhar_VHARtoVMA`, object, lag_max)
}

#' @noRd
compute_vhar_mse <- function(cov_mat, vhar_coef, har_trans, month, step) {
    .Call(`_bvhar_compute_vhar_mse`, cov_mat, vhar_coef, har_trans, month, step)
}

#' Compute Forecast MSE Matrices for VHAR
#' 
#' Compute the forecast MSE matrices using VMA coefficients
#' 
#' @param object \code{varlse} object by \code{\link{var_lm}}
#' @param step Integer, Step to forecast
#' @details
#' See pp38 of Lütkepohl (2007).
#' Let \eqn{\Sigma} be the covariance matrix of VHAR and let \eqn{W_j} be the VMA coefficients.
#' Recursively,
#' \deqn{\Sigma_y(1) = \Sigma}
#' \deqn{\Sigma_y(2) = \Sigma + W_1 \Sigma W_1^T}
#' \deqn{\Sigma_y(3) = \Sigma_y(2) + W_2 \Sigma W_2^T}
#' 
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing. doi:[10.1007/978-3-540-27752-1](https://doi.org/10.1007/978-3-540-27752-1)
#' @noRd
compute_covmse_har <- function(object, step) {
    .Call(`_bvhar_compute_covmse_har`, object, step)
}

#' Orthogonal Impulse Response Functions of VHAR
#' 
#' Compute orthogonal impulse responses of VHAR
#' 
#' @param vhar_coef VHAR coefficient
#' @param vhar_covmat VHAR covariance matrix
#' @param HARtrans_mat HAR linear transformation matrix
#' @param lag_max Maximum lag for VMA
#' @param month Order for monthly term
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing. [https://doi.org/10.1007/978-3-540-27752-1](https://doi.org/10.1007/978-3-540-27752-1)
#' @noRd
VHARcoeftoVMA_ortho <- function(vhar_coef, vhar_covmat, HARtrans_mat, lag_max, month) {
    .Call(`_bvhar_VHARcoeftoVMA_ortho`, vhar_coef, vhar_covmat, HARtrans_mat, lag_max, month)
}

#' h-step ahead Forecast Error Variance Decomposition
#' 
#' [w_(h = 1, ij)^T, w_(h = 2, ij)^T, ...]
#'
#' @noRd
compute_fevd <- function(vma_coef, cov_mat, normalize) {
    .Call(`_bvhar_compute_fevd`, vma_coef, cov_mat, normalize)
}

#' h-step ahead Normalized Spillover
#'
#' @noRd
compute_spillover <- function(fevd) {
    .Call(`_bvhar_compute_spillover`, fevd)
}

#' To-others Spillovers
#' 
#' @noRd
compute_to_spillover <- function(spillover) {
    .Call(`_bvhar_compute_to_spillover`, spillover)
}

#' From-others Spillovers
#' 
#' @noRd
compute_from_spillover <- function(spillover) {
    .Call(`_bvhar_compute_from_spillover`, spillover)
}

#' Total Spillovers
#' 
#' @noRd
compute_tot_spillover <- function(spillover) {
    .Call(`_bvhar_compute_tot_spillover`, spillover)
}

#' Net Pairwise Spillovers
#' 
#' @noRd
compute_net_spillover <- function(spillover) {
    .Call(`_bvhar_compute_net_spillover`, spillover)
}

#' VAR(1) Representation Given VAR Coefficient Matrix
#' 
#' Compute the VAR(1) coefficient matrix form
#' 
#' @param x VAR without constant coefficient matrix form
#' @details
#' Each VAR(p) process can be represented by mp-dim VAR(1).
#' 
#' \deqn{Y_t = A Y_{t - 1} + C + U_t}
#' 
#' where
#' 
#' \deqn{
#'     A = \begin{bmatrix}
#'       A_1 & A_2 & \cdots A_{p - 1} & A_p \\
#'       I_m & 0 & \cdots & 0 & 0 \\
#'       0 & I_m & \cdots & 0 & 0 \\
#'       \vdots & \vdots & \vdots & \vdots & \vdots \\
#'       0 & 0 & \cdots & I_m & 0
#'     \end{bmatrix}
#' }
#' 
#' \deqn{C = (c, 0, \ldots, 0)^T}
#' 
#' and
#' 
#' \deqn{U_t = (\epsilon_t^T, 0^T, \ldots, 0^T)^T}
#' 
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing. doi:[10.1007/978-3-540-27752-1](https://doi.org/10.1007/978-3-540-27752-1)
#' @noRd
compute_stablemat <- function(x) {
    .Call(`_bvhar_compute_stablemat`, x)
}

#' VAR(1) Representation of VAR(p)
#' 
#' Compute the coefficient matrix of VAR(1) form
#' 
#' @param object Model fit
#' 
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing. doi:[10.1007/978-3-540-27752-1](https://doi.org/10.1007/978-3-540-27752-1)
#' @noRd
compute_var_stablemat <- function(object) {
    .Call(`_bvhar_compute_var_stablemat`, object)
}

#' VAR(1) Representation of VHAR
#'
#' Compute the coefficient matrix of VAR(1) form of VHAR
#'
#' @param object Model fit
#' @details
#' Note that \eqn{A^T = \Phi^T T_{HAR}}.
#' This gives the VAR(1) form of constrained VAR(22).
#'
#' @references Lütkepohl, H. (2007). *New Introduction to Multiple Time Series Analysis*. Springer Publishing. doi:[10.1007/978-3-540-27752-1](https://doi.org/10.1007/978-3-540-27752-1)
#' @noRd
compute_vhar_stablemat <- function(object) {
    .Call(`_bvhar_compute_vhar_stablemat`, object)
}

#' Log of Multivariate Gamma Function
#' 
#' Compute log of multivariate gamma function numerically
#' 
#' @param x Double, non-negative argument
#' @param p Integer, dimension
#' @noRd
log_mgammafn <- function(x, p) {
    .Call(`_bvhar_log_mgammafn`, x, p)
}

#' Numerically Stable Log ML Excluding Constant Term of BVAR and BVHAR
#' 
#' This function computes log of ML stable,
#' in purpose of objective function.
#' 
#' @param object Bayesian Model Fit
#' 
#' @noRd
logml_stable <- function(object) {
    .Call(`_bvhar_logml_stable`, object)
}

#' AIC of VAR(p) using RSS
#' 
#' Compute AIC using RSS
#' 
#' @param object A `varlse` or `vharlse` object
#' 
#' @noRd
compute_aic <- function(object) {
    .Call(`_bvhar_compute_aic`, object)
}

#' BIC of VAR(p) using RSS
#' 
#' Compute BIC using RSS
#' 
#' @param object A `varlse` or `vharlse` object
#' 
#' @noRd
compute_bic <- function(object) {
    .Call(`_bvhar_compute_bic`, object)
}

#' HQ of VAR(p) using RSS
#' 
#' Compute HQ using RSS
#' 
#' @param object A `varlse` or `vharlse` object
#' 
#' @noRd
compute_hq <- function(object) {
    .Call(`_bvhar_compute_hq`, object)
}

#' FPE of VAR(p) using RSS
#' 
#' Compute FPE using RSS
#' 
#' @param object A `varlse` or `vharlse` object
#' 
#' @noRd
compute_fpe <- function(object) {
    .Call(`_bvhar_compute_fpe`, object)
}

#' Choose the Best VAR based on Information Criteria
#' 
#' This function computes AIC, FPE, BIC, and HQ up to p = `lag_max` of VAR model.
#' 
#' @param y Time series data of which columns indicate the variables
#' @param lag_max Maximum Var lag to explore
#' @param include_mean Add constant term
#' 
#' @noRd
tune_var <- function(y, lag_max, include_mean) {
    .Call(`_bvhar_tune_var`, y, lag_max, include_mean)
}