diff --git a/NAMESPACE b/NAMESPACE
index 96ea8ba..47064cc 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -1,6 +1,7 @@
 # Generated by roxygen2: do not edit by hand
 
 S3method(estimate,BSVARSIGN)
+S3method(forecast,PosteriorBSVARSIGN)
 export(specify_bsvarSIGN)
 export(specify_identification_bsvarSIGN)
 export(specify_posterior_bsvarSIGN)
diff --git a/R/forecast.R b/R/forecast.R
new file mode 100644
index 0000000..984dd9c
--- /dev/null
+++ b/R/forecast.R
@@ -0,0 +1,132 @@
+
+#' @title Forecasting using Structural Vector Autoregression
+#'
+#' @description Samples from the joint predictive density of all of the dependent 
+#' variables for models from packages \pkg{bsvars}, \pkg{bsvarSIGNs} or 
+#' \pkg{bvarPANELs} at forecast horizons from 1 to \code{horizon} specified as 
+#' an argument of the function.
+#' 
+#' @method forecast PosteriorBSVARSIGN
+#' @param posterior posterior estimation outcome - an object of class 
+#' \code{PosteriorBSVARSIGN} obtained by running the \code{estimate} function.
+#' @param horizon a positive integer, specifying the forecasting horizon.
+#' @param exogenous_forecast a matrix of dimension \code{horizon x d} containing 
+#' forecasted values of the exogenous variables. 
+#' @param conditional_forecast a \code{horizon x N} matrix with forecasted values 
+#' for selected variables. It should only contain \code{numeric} or \code{NA} 
+#' values. The entries with \code{NA} values correspond to the values that are 
+#' forecasted conditionally on the realisations provided as \code{numeric} values.
+#' 
+#' @return A list of class \code{Forecasts} containing the
+#' draws from the predictive density and data. The output list includes element:
+#' 
+#' \describe{
+#'  \item{forecasts}{an \code{NxhorizonxS} array with the draws from predictive density}
+#'  \item{Y}{an \eqn{NxT} matrix with the data on dependent variables}
+#' }
+#' 
+#' @author Tomasz Woźniak \email{wozniak.tom@pm.me} and Xiaolei Wang \email{adamwang15@gmail.com}
+#' 
+#' @examples
+#' # upload data
+#' data(oil)
+#' 
+#' # specify the model and set seed
+#' set.seed(123)
+#' sign_irf       = array(matrix(c(-1, -1, 1, rep(0, 6)), nrow = 3), dim = c(3, 3, 1))
+#' specification  = specify_bsvarSIGN$new(oil, p = 12, sign_irf = sign_irf)
+#' 
+#' # estimate the model
+#' posterior      = estimate(specification, 20)
+#' 
+#' # sample from predictive density 1 year ahead
+#' predictive     = forecast(posterior, 4)
+#' 
+#' # workflow with the pipe |>
+#' ############################################################
+#' set.seed(123)
+#' oil |>
+#'   specify_bsvarSIGN$new(p = 12, sign_irf = sign_irf) |>
+#'   estimate(S = 20) |> 
+#'   forecast(horizon = 4) -> predictive
+#' 
+#' # conditional forecasting 2 quarters ahead conditioning on 
+#' #  provided future values for the Gross Domestic Product 
+#' ############################################################
+#' cf        = matrix(NA , 2, 3)
+#' cf[,3]    = tail(oil, 1)[3]   # conditional forecasts equal to the last gdp observation
+#' predictive    = forecast(posterior, 2, conditional_forecast = cf)
+#' 
+#' # workflow with the pipe |>
+#' ############################################################
+#' set.seed(123)
+#' oil |>
+#'   specify_bsvarSIGN$new(p = 12, sign_irf = sign_irf) |>
+#'   estimate(S = 10) |> 
+#'   forecast(horizon = 2, conditional_forecast = cf) -> predictive
+#' 
+#' @export
+forecast.PosteriorBSVARSIGN = function(
+    posterior, 
+    horizon = 1, 
+    exogenous_forecast = NULL,
+    conditional_forecast = NULL
+) {
+  
+  posterior_Sigma = posterior$posterior$Sigma
+  posterior_A     = posterior$posterior$A
+  T               = nrow(posterior$last_draw$data_matrices$X)
+  X_T             = posterior$last_draw$data_matrices$X[T,]
+  Y               = posterior$last_draw$data_matrices$Y
+  
+  N               = nrow(posterior_Sigma)
+  K               = length(X_T)
+  d               = K - N * posterior$last_draw$p - 1
+  S               = dim(posterior_Sigma)[3]
+  
+  # prepare forecasting with exogenous variables
+  if (d == 0 ) {
+    exogenous_forecast = matrix(NA, horizon, 1)
+  } else {
+    stopifnot("Forecasted values of exogenous variables are missing." = (d > 0) & !is.null(exogenous_forecast))
+    stopifnot("The matrix of exogenous_forecast does not have a correct number of columns." = ncol(exogenous_forecast) == d)
+    stopifnot("Provide exogenous_forecast for all forecast periods specified by argument horizon." = nrow(exogenous_forecast) == horizon)
+    stopifnot("Argument exogenous has to be a matrix." = is.matrix(exogenous_forecast) & is.numeric(exogenous_forecast))
+    stopifnot("Argument exogenous cannot include missing values." = sum(is.na(exogenous_forecast)) == 0 )
+  }
+  
+  # prepare forecasting with conditional forecasts
+  if ( is.null(conditional_forecast) ) {
+    # this will not be used for forecasting, but needs to be provided
+    conditional_forecast = matrix(NA, horizon, N)
+  } else {
+    stopifnot("Argument conditional_forecast must be a matrix with numeric values."
+              = is.matrix(conditional_forecast) & is.numeric(conditional_forecast)
+    )
+    stopifnot("Argument conditional_forecast must have the number of rows equal to 
+              the value of argument horizon."
+              = nrow(conditional_forecast) == horizon
+    )
+    stopifnot("Argument conditional_forecast must have the number of columns 
+              equal to the number of columns in the used data."
+              = ncol(conditional_forecast) == N
+    )
+  }
+  
+  # perform forecasting
+  for_y       = .Call(`_bsvarSIGNs_forecast_bsvarSIGNs`, 
+                      posterior_Sigma,
+                      posterior_A,
+                      X_T,
+                      exogenous_forecast,
+                      conditional_forecast,
+                      horizon
+  ) # END .Call
+  
+  fore            = list()
+  fore$forecasts  = for_y
+  fore$Y          = Y
+  class(fore)     = "Forecasts"
+  
+  return(fore)
+} # END forecast.PosteriorBSVARSIGN
diff --git a/inst/include/bsvarSIGNs_RcppExports.h b/inst/include/bsvarSIGNs_RcppExports.h
index 43061ec..8580342 100644
--- a/inst/include/bsvarSIGNs_RcppExports.h
+++ b/inst/include/bsvarSIGNs_RcppExports.h
@@ -109,6 +109,27 @@ namespace bsvarSIGNs {
         return Rcpp::as<Rcpp::List >(rcpp_result_gen);
     }
 
+    inline arma::cube forecast_bsvarSIGNs(arma::cube& posterior_Sigma, arma::cube& posterior_A, arma::vec& X_T, arma::mat& exogenous_forecast, arma::mat& cond_forecast, const int& horizon) {
+        typedef SEXP(*Ptr_forecast_bsvarSIGNs)(SEXP,SEXP,SEXP,SEXP,SEXP,SEXP);
+        static Ptr_forecast_bsvarSIGNs p_forecast_bsvarSIGNs = NULL;
+        if (p_forecast_bsvarSIGNs == NULL) {
+            validateSignature("arma::cube(*forecast_bsvarSIGNs)(arma::cube&,arma::cube&,arma::vec&,arma::mat&,arma::mat&,const int&)");
+            p_forecast_bsvarSIGNs = (Ptr_forecast_bsvarSIGNs)R_GetCCallable("bsvarSIGNs", "_bsvarSIGNs_forecast_bsvarSIGNs");
+        }
+        RObject rcpp_result_gen;
+        {
+            RNGScope RCPP_rngScope_gen;
+            rcpp_result_gen = p_forecast_bsvarSIGNs(Shield<SEXP>(Rcpp::wrap(posterior_Sigma)), Shield<SEXP>(Rcpp::wrap(posterior_A)), Shield<SEXP>(Rcpp::wrap(X_T)), Shield<SEXP>(Rcpp::wrap(exogenous_forecast)), Shield<SEXP>(Rcpp::wrap(cond_forecast)), Shield<SEXP>(Rcpp::wrap(horizon)));
+        }
+        if (rcpp_result_gen.inherits("interrupted-error"))
+            throw Rcpp::internal::InterruptedException();
+        if (Rcpp::internal::isLongjumpSentinel(rcpp_result_gen))
+            throw Rcpp::LongjumpException(rcpp_result_gen);
+        if (rcpp_result_gen.inherits("try-error"))
+            throw Rcpp::exception(Rcpp::as<std::string>(rcpp_result_gen).c_str());
+        return Rcpp::as<arma::cube >(rcpp_result_gen);
+    }
+
     inline bool match_sign_narrative(const arma::mat& Epsilon, const arma::mat& sign_narrative, const arma::cube& irf) {
         typedef SEXP(*Ptr_match_sign_narrative)(SEXP,SEXP,SEXP);
         static Ptr_match_sign_narrative p_match_sign_narrative = NULL;
diff --git a/inst/tinytest/test_forecast.R b/inst/tinytest/test_forecast.R
new file mode 100644
index 0000000..edb132a
--- /dev/null
+++ b/inst/tinytest/test_forecast.R
@@ -0,0 +1,41 @@
+
+data(oil)
+
+# for bsvar
+set.seed(1)
+suppressMessages(
+  specification_no1 <- specify_bsvarSIGN$new(oil, p = 1)
+)
+run_no1             <- estimate(specification_no1, 3, 1, show_progress = FALSE)
+ff                  <- forecast(run_no1, horizon = 2)
+
+set.seed(1)
+suppressMessages(
+  ff2              <- oil |>
+    specify_bsvarSIGN$new(p = 1) |>
+    estimate(S = 3, thin = 1, show_progress = FALSE) |>
+    forecast(horizon = 2)
+)
+
+
+expect_identical(
+  ff$forecasts[1,1,1], ff2$forecasts[1,1,1],
+  info = "forecast: forecast identical for normal and pipe workflow."
+)
+
+expect_true(
+  is.numeric(ff$forecasts) & is.array(ff$forecasts),
+  info = "forecast: returns numeric array."
+)
+
+
+expect_error(
+  specify_bsvar$new(us_fiscal_lsuw) |> forecast(horizon = 3),
+  info = "forecast: wrong input provided."
+)
+
+expect_error(
+  forecast(run_no1, horizon = 1.5),
+  info = "forecast: specify horizon as integer."
+)
+
diff --git a/man/forecast.PosteriorBSVARSIGN.Rd b/man/forecast.PosteriorBSVARSIGN.Rd
new file mode 100644
index 0000000..3cced11
--- /dev/null
+++ b/man/forecast.PosteriorBSVARSIGN.Rd
@@ -0,0 +1,84 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/forecast.R
+\name{forecast.PosteriorBSVARSIGN}
+\alias{forecast.PosteriorBSVARSIGN}
+\title{Forecasting using Structural Vector Autoregression}
+\usage{
+\method{forecast}{PosteriorBSVARSIGN}(
+  posterior,
+  horizon = 1,
+  exogenous_forecast = NULL,
+  conditional_forecast = NULL
+)
+}
+\arguments{
+\item{posterior}{posterior estimation outcome - an object of class 
+\code{PosteriorBSVARSIGN} obtained by running the \code{estimate} function.}
+
+\item{horizon}{a positive integer, specifying the forecasting horizon.}
+
+\item{exogenous_forecast}{a matrix of dimension \code{horizon x d} containing 
+forecasted values of the exogenous variables.}
+
+\item{conditional_forecast}{a \code{horizon x N} matrix with forecasted values 
+for selected variables. It should only contain \code{numeric} or \code{NA} 
+values. The entries with \code{NA} values correspond to the values that are 
+forecasted conditionally on the realisations provided as \code{numeric} values.}
+}
+\value{
+A list of class \code{Forecasts} containing the
+draws from the predictive density and data. The output list includes element:
+
+\describe{
+ \item{forecasts}{an \code{NxhorizonxS} array with the draws from predictive density}
+ \item{Y}{an \eqn{NxT} matrix with the data on dependent variables}
+}
+}
+\description{
+Samples from the joint predictive density of all of the dependent 
+variables for models from packages \pkg{bsvars}, \pkg{bsvarSIGNs} or 
+\pkg{bvarPANELs} at forecast horizons from 1 to \code{horizon} specified as 
+an argument of the function.
+}
+\examples{
+# upload data
+data(oil)
+
+# specify the model and set seed
+set.seed(123)
+sign_irf       = array(matrix(c(-1, -1, 1, rep(0, 6)), nrow = 3), dim = c(3, 3, 1))
+specification  = specify_bsvarSIGN$new(oil, p = 12, sign_irf = sign_irf)
+
+# estimate the model
+posterior      = estimate(specification, 20)
+
+# sample from predictive density 1 year ahead
+predictive     = forecast(posterior, 4)
+
+# workflow with the pipe |>
+############################################################
+set.seed(123)
+oil |>
+  specify_bsvarSIGN$new(p = 12, sign_irf = sign_irf) |>
+  estimate(S = 20) |> 
+  forecast(horizon = 4) -> predictive
+
+# conditional forecasting 2 quarters ahead conditioning on 
+#  provided future values for the Gross Domestic Product 
+############################################################
+cf        = matrix(NA , 2, 3)
+cf[,3]    = tail(oil, 1)[3]   # conditional forecasts equal to the last gdp observation
+predictive    = forecast(posterior, 2, conditional_forecast = cf)
+
+# workflow with the pipe |>
+############################################################
+set.seed(123)
+oil |>
+  specify_bsvarSIGN$new(p = 12, sign_irf = sign_irf) |>
+  estimate(S = 10) |> 
+  forecast(horizon = 2, conditional_forecast = cf) -> predictive
+
+}
+\author{
+Tomasz Woźniak \email{wozniak.tom@pm.me} and Xiaolei Wang \email{adamwang15@gmail.com}
+}
diff --git a/src/RcppExports.cpp b/src/RcppExports.cpp
index aef846b..85deade 100644
--- a/src/RcppExports.cpp
+++ b/src/RcppExports.cpp
@@ -173,6 +173,45 @@ RcppExport SEXP _bsvarSIGNs_bsvar_sign_cpp(SEXP SSEXP, SEXP lagsSEXP, SEXP YSEXP
     UNPROTECT(1);
     return rcpp_result_gen;
 }
+// forecast_bsvarSIGNs
+arma::cube forecast_bsvarSIGNs(arma::cube& posterior_Sigma, arma::cube& posterior_A, arma::vec& X_T, arma::mat& exogenous_forecast, arma::mat& cond_forecast, const int& horizon);
+static SEXP _bsvarSIGNs_forecast_bsvarSIGNs_try(SEXP posterior_SigmaSEXP, SEXP posterior_ASEXP, SEXP X_TSEXP, SEXP exogenous_forecastSEXP, SEXP cond_forecastSEXP, SEXP horizonSEXP) {
+BEGIN_RCPP
+    Rcpp::RObject rcpp_result_gen;
+    Rcpp::traits::input_parameter< arma::cube& >::type posterior_Sigma(posterior_SigmaSEXP);
+    Rcpp::traits::input_parameter< arma::cube& >::type posterior_A(posterior_ASEXP);
+    Rcpp::traits::input_parameter< arma::vec& >::type X_T(X_TSEXP);
+    Rcpp::traits::input_parameter< arma::mat& >::type exogenous_forecast(exogenous_forecastSEXP);
+    Rcpp::traits::input_parameter< arma::mat& >::type cond_forecast(cond_forecastSEXP);
+    Rcpp::traits::input_parameter< const int& >::type horizon(horizonSEXP);
+    rcpp_result_gen = Rcpp::wrap(forecast_bsvarSIGNs(posterior_Sigma, posterior_A, X_T, exogenous_forecast, cond_forecast, horizon));
+    return rcpp_result_gen;
+END_RCPP_RETURN_ERROR
+}
+RcppExport SEXP _bsvarSIGNs_forecast_bsvarSIGNs(SEXP posterior_SigmaSEXP, SEXP posterior_ASEXP, SEXP X_TSEXP, SEXP exogenous_forecastSEXP, SEXP cond_forecastSEXP, SEXP horizonSEXP) {
+    SEXP rcpp_result_gen;
+    {
+        Rcpp::RNGScope rcpp_rngScope_gen;
+        rcpp_result_gen = PROTECT(_bsvarSIGNs_forecast_bsvarSIGNs_try(posterior_SigmaSEXP, posterior_ASEXP, X_TSEXP, exogenous_forecastSEXP, cond_forecastSEXP, horizonSEXP));
+    }
+    Rboolean rcpp_isInterrupt_gen = Rf_inherits(rcpp_result_gen, "interrupted-error");
+    if (rcpp_isInterrupt_gen) {
+        UNPROTECT(1);
+        Rf_onintr();
+    }
+    bool rcpp_isLongjump_gen = Rcpp::internal::isLongjumpSentinel(rcpp_result_gen);
+    if (rcpp_isLongjump_gen) {
+        Rcpp::internal::resumeJump(rcpp_result_gen);
+    }
+    Rboolean rcpp_isError_gen = Rf_inherits(rcpp_result_gen, "try-error");
+    if (rcpp_isError_gen) {
+        SEXP rcpp_msgSEXP_gen = Rf_asChar(rcpp_result_gen);
+        UNPROTECT(1);
+        Rf_error("%s", CHAR(rcpp_msgSEXP_gen));
+    }
+    UNPROTECT(1);
+    return rcpp_result_gen;
+}
 // match_sign_narrative
 bool match_sign_narrative(const arma::mat& Epsilon, const arma::mat& sign_narrative, const arma::cube& irf);
 static SEXP _bsvarSIGNs_match_sign_narrative_try(SEXP EpsilonSEXP, SEXP sign_narrativeSEXP, SEXP irfSEXP) {
@@ -869,6 +908,7 @@ static int _bsvarSIGNs_RcppExport_validate(const char* sig) {
         signatures.insert("arma::cube(*ir1_cpp)(const arma::mat&,const arma::mat&,int,const int&)");
         signatures.insert("arma::field<arma::cube>(*bsvarSIGNs_ir)(arma::cube&,arma::cube&,const int,const int)");
         signatures.insert("Rcpp::List(*bsvar_sign_cpp)(const int&,const int&,const arma::mat&,const arma::mat&,const arma::field<arma::mat>&,const arma::cube&,const arma::mat&,const arma::mat&,const arma::field<arma::mat>&,const Rcpp::List&,const Rcpp::List&,const bool,const int,const int&)");
+        signatures.insert("arma::cube(*forecast_bsvarSIGNs)(arma::cube&,arma::cube&,arma::vec&,arma::mat&,arma::mat&,const int&)");
         signatures.insert("bool(*match_sign_narrative)(const arma::mat&,const arma::mat&,const arma::cube&)");
         signatures.insert("double(*weight_narrative)(const int&,arma::mat,const arma::cube&)");
         signatures.insert("arma::field<arma::mat>(*ZIRF)(const arma::field<arma::mat>&,const arma::mat&)");
@@ -893,6 +933,7 @@ RcppExport SEXP _bsvarSIGNs_RcppExport_registerCCallable() {
     R_RegisterCCallable("bsvarSIGNs", "_bsvarSIGNs_ir1_cpp", (DL_FUNC)_bsvarSIGNs_ir1_cpp_try);
     R_RegisterCCallable("bsvarSIGNs", "_bsvarSIGNs_bsvarSIGNs_ir", (DL_FUNC)_bsvarSIGNs_bsvarSIGNs_ir_try);
     R_RegisterCCallable("bsvarSIGNs", "_bsvarSIGNs_bsvar_sign_cpp", (DL_FUNC)_bsvarSIGNs_bsvar_sign_cpp_try);
+    R_RegisterCCallable("bsvarSIGNs", "_bsvarSIGNs_forecast_bsvarSIGNs", (DL_FUNC)_bsvarSIGNs_forecast_bsvarSIGNs_try);
     R_RegisterCCallable("bsvarSIGNs", "_bsvarSIGNs_match_sign_narrative", (DL_FUNC)_bsvarSIGNs_match_sign_narrative_try);
     R_RegisterCCallable("bsvarSIGNs", "_bsvarSIGNs_weight_narrative", (DL_FUNC)_bsvarSIGNs_weight_narrative_try);
     R_RegisterCCallable("bsvarSIGNs", "_bsvarSIGNs_ZIRF", (DL_FUNC)_bsvarSIGNs_ZIRF_try);
@@ -916,6 +957,7 @@ static const R_CallMethodDef CallEntries[] = {
     {"_bsvarSIGNs_ir1_cpp", (DL_FUNC) &_bsvarSIGNs_ir1_cpp, 4},
     {"_bsvarSIGNs_bsvarSIGNs_ir", (DL_FUNC) &_bsvarSIGNs_bsvarSIGNs_ir, 4},
     {"_bsvarSIGNs_bsvar_sign_cpp", (DL_FUNC) &_bsvarSIGNs_bsvar_sign_cpp, 14},
+    {"_bsvarSIGNs_forecast_bsvarSIGNs", (DL_FUNC) &_bsvarSIGNs_forecast_bsvarSIGNs, 6},
     {"_bsvarSIGNs_match_sign_narrative", (DL_FUNC) &_bsvarSIGNs_match_sign_narrative, 3},
     {"_bsvarSIGNs_weight_narrative", (DL_FUNC) &_bsvarSIGNs_weight_narrative, 3},
     {"_bsvarSIGNs_ZIRF", (DL_FUNC) &_bsvarSIGNs_ZIRF, 2},
diff --git a/src/forecast_bsvarSIGNs.cpp b/src/forecast_bsvarSIGNs.cpp
new file mode 100644
index 0000000..0966ec9
--- /dev/null
+++ b/src/forecast_bsvarSIGNs.cpp
@@ -0,0 +1,70 @@
+
+#include <RcppArmadillo.h>
+#include "bsvars.h"
+
+using namespace Rcpp;
+using namespace arma;
+
+
+// [[Rcpp::interfaces(cpp)]]
+// [[Rcpp::export]]
+arma::cube forecast_bsvarSIGNs (
+    arma::cube&   posterior_Sigma,    // (N, N, S)
+    arma::cube&   posterior_A,        // (N, K, S)
+    arma::vec&    X_T,                // (K)
+    arma::mat&    exogenous_forecast, // (horizon, d)
+    arma::mat&    cond_forecast,      // (horizon, N)
+    const int&    horizon
+) {
+  
+  const int   N = posterior_Sigma.n_rows;
+  const int   S = posterior_Sigma.n_slices;
+  const int   K = posterior_A.n_cols;
+  const int   d = exogenous_forecast.n_cols;
+  
+  bool        do_exog = exogenous_forecast.is_finite();
+  vec         x_t;
+  if ( do_exog ) {
+    x_t       = X_T.rows(0, K - 1 - d);
+  } else {
+    x_t       = X_T.rows(0, K - 1);
+  } // END if do_exog
+  
+  vec         Xt(K);
+  cube        forecasts(N, horizon, S);
+  
+  for (int s=0; s<S; s++) {
+    
+    if ( do_exog ) {
+      Xt          = join_cols(x_t, trans(exogenous_forecast.row(0)));
+    } else {
+      Xt          = x_t;
+    } // END if do_exog
+    
+    for (int h=0; h<horizon; h++) {
+      
+      mat   Sigma             = posterior_Sigma.slice(s);
+      vec   cond_forecast_h   = trans(cond_forecast.row(h));
+      uvec  nonf_el           = find_nonfinite( cond_forecast_h );
+      int   nonf_no           = nonf_el.n_elem;
+      
+      if ( nonf_no == N ) {
+        forecasts.slice(s).col(h) = mvnrnd( posterior_A.slice(s) * Xt, Sigma );
+      } else {
+        forecasts.slice(s).col(h) = bsvars::mvnrnd_cond( cond_forecast_h, posterior_A.slice(s) * Xt, Sigma );   // does not work if cond_fc_h is all nan
+      } // END if nonf_no
+      
+      if ( h != horizon - 1 ) {
+        if ( do_exog ) {
+          Xt          = join_cols( forecasts.slice(s).col(h), Xt.subvec(N, K - 1 - d), trans(cond_forecast.row(h + 1)) );
+        } else {
+          Xt          = join_cols( forecasts.slice(s).col(h), Xt.subvec(N, K - 1) );
+        }
+      } // END if h
+      
+    } // END h loop
+  } // END s loop
+  
+  return forecasts;
+} // END forecast_bsvarSIGNs
+
diff --git a/src/forecast_bsvarSIGNs.h b/src/forecast_bsvarSIGNs.h
new file mode 100644
index 0000000..55cba7c
--- /dev/null
+++ b/src/forecast_bsvarSIGNs.h
@@ -0,0 +1,18 @@
+
+#ifndef _FORECAST_BSVARSIGNS_H_
+#define _FORECAST_BSVARSIGNS_H_
+
+#include <RcppArmadillo.h>
+
+
+arma::cube forecast_bsvarSIGNs (
+    arma::cube&   posterior_Sigma,    // (N, N, S)
+    arma::cube&   posterior_A,        // (N, K, S)
+    arma::vec&    X_T,                // (K)
+    arma::mat&    exogenous_forecast, // (horizon, d)
+    arma::mat&    cond_forecast,      // (horizon, N)
+    const int&    horizon
+);
+
+
+#endif  // _FORECAST_BSVARSIGNS_H_
\ No newline at end of file