Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
117,889 characters · 23 sections · 89 citation commands
\setcounter{page}{0}
\thispagestyle{empty}
{ \qquad Daniel Dzikowski\footnote{TU Dortmund University, Department of Statistics, D-44221 Dortmund, Germany; [email removed]; corresponding author} \qquad \qquad \qquad Carsten Jentsch\footnote{TU Dortmund University, Department of Statistics, D-44221 Dortmund, Germany; [email removed]}}
{ \quad TU Dortmund University \qquad \quad \; \; \; TU Dortmund University}
{\bf Keywords:} Periodic Vector Autoregressions, Impulse Response Analysis, Linear Restrictions, Seasonal Adjustment, Residual-Based Seasonal Bootstrap
{\bf JEL Codes:} C30, C32
\setcounter{footnote}{0}
Macroeconomic time series often exhibit seasonal structure defined as recurrent intra-year movement. The reasons for seasonal fluctuations in economic time series are manifold. Classical literature as, for instance, hylleberg1992modelling, hylleberg1993seasonality, canova1995seasonal and franses1996recent argue that seasonality of time series may be caused by calendar and weather effects or by seasonally varying behavior of economic agents, among others. To illustrate, we consider a seasonally unadjusted monthly data set consisting of industrial production (IP), CPI inflation (INF) and the federal funds rate (FFR). In Figure (ref), we present estimated spectral densities (SDs) and autocovariance functions (ACFs) for all three (univariate) time series. We find that IP and INF exhibit strong periodic annual patterns, while FFR shows no periodic structure at all. The standard approach to deal with seasonality is to remove the seasonal structure of macroeconomic time series by using seasonal adjustment methods prior to a further analysis. The most popular approaches are the classical X-11 seasonal adjustment program established by shishkin1967x, the TRAMO/SEATS program by gomez1996programs and model-based upgrades X-12-ARIMA findley1998new, X-13-ARIMA-SEATS monsell2007x and Demetra$+$ eurostat2009guidelines. However, by removing seasonality in raw macroeconomic data before it is used for structural inference, this may distort valuable information in the data. Hence, in the literature, seasonal adjustment methods are often controversial due to their high complexity, lack of linkage to economic theory and the resulting lack of transparency; see, for instance, gersovitz1978seasonality, bell1984issues, osborn1989performance, franses1996recent and mcelroy2022review. ghysels1996seasonal and ghysels2001econometric show that seasonal adjustment methods may distort the structure of the data leading to misspecified standard errors and inference. Recently, doppelt2021should argue that the results of SVAR analyses are highly dependent on the used seasonal adjustment method and that structural analysis of seasonally unadjusted data should be preferred. For instance, when relying on the common technique of seasonal demeaning, this may not capture the entire periodic structure in the data. This becomes visible in Figure (ref), where the estimated SDs and ACFs of seasonally demeaned IP, INF and FFR are displayed. Apparently, while seasonal demeaning removes almost all periodicity of INF, it is clearly not able to adequately remove the entire periodic structure of IP. That is, even after seasonal demeaning, IP still exhibits pronounced periodic structure. While this leftover periodicity is not surprising as only a simple seasonal demeaning is used to get Figure (ref), even after applying more sophisticated seasonal adjustment methods, the seasonally adjusted data may still display seasonal structure; see del2004consequences, bell2012unit and findley2017detecting. As an example, rudebusch2015puzzle and lunsford2017lingering found residual or leftover seasonality in seasonally adjusted gross domestic product estimates by the Bureau of Economic Analysis (BEA). Even after BEA introduced a new estimation strategy to remove residual seasonality of macro variables in 2018, consolvo2019residual and lunsford2025residual detected remaining residual seasonality in several adjusted macro variables after BEAs improvements.
In this paper, to address the issue of often highly complex and not transparent seasonal adjustment methods that lack a reasonable economic interpretation, we propose an alternative to the commonly used vector autoregressions (VARs) fitted after seasonal adjustment. More precisely, we propose the direct modeling of seasonally unadjusted macroeconomic data based on periodic autoregressions. This approach enables a general framework for structural impulse response analyses that allows to take seasonal patterns directly into account. Periodic autoregressions are useful to capture not only the seasonal structure present in the mean of time series data, but also in their covariances in the sense that both the mean and the (auto)covariances become periodic functions of time. Processes with periodically varying covariances, often called periodically correlated stochastic processes, have a long history. First introduced by gladyshev1961periodically, pioneering work on (univariate) periodic autoregressive (PAR) processes was done by jones1967time, pagano1978periodic and troutman1979some, who derived inference techniques on periodic autoregressive parameter estimates, predictions using PAR models and the connection between PAR and the related stationary multivariate autoregressive process. PARs for modeling seasonally unadjusted macroeconomic data have been used by osborn1988seasonality, osborn1989performance, novales1997forecasting, franses2005forecasting and ghysels2006forecasting, who mainly focused on forecasting studies. However, although PAR models have received some attention in economics and certain methodological developments have been made franses2004periodic, they have not received much popularity for structural analysis. This is mainly because, for analyzing seasonally unadjusted macroeconomic data, a vector-valued extension of the PAR model class is required that naturally comes with a lot of parameters.
Asymptotic theory for such periodic vector autoregressive (PVAR) models have been developed by ula1990periodic,ula1993forecasting, franses2004periodic and aknouche2007causality, who derived periodic stationarity conditions and parameter estimation methods. Testing for stationarity against periodic stationary alternatives is considered in jentsch2012periodic. Further, P(V)AR models are often considered in stochastic volatility modeling of financial data, see e.g. bollerslev1996periodic, aknouche2009quasi and aknouche2017periodic. Thanks to their large number of parameters, PVAR models are extremely flexible. However, unconstrained PVAR estimation does often not work well in applications even with moderate sample sizes as typically encountered in macroeconometrics due to pronounced overfitting issues. To address these, ursu2009modelling and boubacar2023estimating derived least squares estimators under linear parameter restrictions. However, they only allow to impose linear restrictions on parameters within a given season to control the number of parameters to be estimated. In this paper, we extend their approach by explicitly allowing also for practically relevant linear constraints across seasons. The importance of constraints across seasons in a PVAR setup is also illustrated in Figures (ref) and (ref). Both figures show that for INF and FFR the autoregressive coefficients may be restricted to be the same across seasons, given evidence that there is hardly any seasonality left after allowing for seasonal mean shifts. This renders the model more parsimonious, facilitating the estimation of the PVAR model.
This paper is organized as follows. Section (ref) introduces PVAR processes and their properties. In Section (ref), we derive (contraint) least-squares estimation for the PVAR parameters that extend the approaches of ursu2009modelling and boubacar2023estimating by allowing also for practically relevant constraints across seasons. We provide asymptotic theory for the (constrained) estimators of the PVAR parameters under (possibly weak) periodic white noise assumptions. Section (ref) introduces the structural PVAR model, which enables a direct structural analysis of periodic time series without prior seasonal adjustment. We discuss two identification schemes for SPVAR models and derive asymptotic theory for corresponding (structural) impulse response analysis. Section (ref) proposes suitable residual-based seasonal bootstrap methods for the construction of confidence intervals and provides bootstrap consistency results. Further, a test for seasonality in impulse responses is introduced. A real data application on seasonally unadjusted IP, INF and FFR is provided in Section (ref), while Section (ref) summarizes our findings. All proofs and additional figures are deferred to the appendix.
In this section, we discuss periodic vector autoregressive (PVAR) processes and their properties. Let $\{ y_t \}_{t\in \mathbb{Z}} = \{ y_{Sn+s} \}_{n\in \mathbb{Z}, s = 1,\dots, S}$ be an $m$-dimensional PVAR process ursu2009modelling described by the model equation
where $S$ is the number of seasons within one cycle and $y_{Sn+s}$ denotes the observation in season $s$ and cycle $n$. The model order $p(s)$, the $m$-dimensional seasonal intercept $\nu(s)$, and the ($m \times m$) seasonal autoregressive coefficient matrices $A_{i}(s)$, $i=1,\dots, p(s)$ are allowed to vary across the seasons $s=1,\dots,S$. The innovation process $\{ \epsilon_{Sn+s} \}_{n\in \mathbb{Z}, s = 1,\dots, S}$ is assumed to be an $m$-dimensional periodic white noise process with period $S$ and non-singular covariance matrices $E(\epsilon_{Sn+s} \epsilon^{\prime}_{Sn+s})=\Sigma_{\epsilon}(s)$, $s=1\dots, S$ meaning that
where $\{w_{t} \}_{t\in \mathbb{Z}}$ is a white noise process with $E( w_{t} w_{t}^{\prime} ) = \boldsymbol{I}_{m}$ for all $t\in\mathbb{Z}$ and the non-singular matrix $H_0(s)$ satisfies $ \Sigma_{\epsilon}(s) = H_0(s)H_0(s)^{\prime}$. Here, $\boldsymbol{I}_{m}$ denotes the $m$-dimensional identity matrix and $ \Sigma_{\epsilon}(s)$ is allowed to vary across the seasons $s=1,\dots,S$. Since the model order $p(s)$ of a PVAR process can vary across seasons as well, we introduce the notation PVAR($\boldsymbol{p}$) with $\boldsymbol{p} = (p(1), \dots, p(S))$ and PVAR($p$), if $p(s) = p$ for all $s=1,\dots,S$. For $S=1$, the PVAR($\boldsymbol{p}$)=PVAR($p$) process collapses to a common VAR($p$) process. However, the PVAR($\boldsymbol{p}$) in (ref) can also be represented with constant order $p$, where $p =\max\{p(1), \dots, p(S)\}$ by imposing zero restrictions on $A_{j}(s)$ for $p(s) < j \leq p$ lund2000recursive.
Throughout the paper, we assume that $\{ y_{Sn+s} \}_{n\in \mathbb{Z}, s = 1,\dots, S}$ is periodically stationary and a causal infinite-order moving-average representation of (ref) exists. franses2004periodic use the higher-dimensional VAR representation of PVAR($\boldsymbol{p}$) processes to obtain conditions for periodic stationarity and for causality of PVAR processes. These conditions can be found in Appendix C. The moving-average representation of a periodically stationary (and causal) process $\{y_{Sn+s}\}_{n\in \mathbb{Z}, s = 1,\dots, S}$ is given by
where $\Phi_{k}(s)$ are $(m \times m)$-dimensional matrices consisting of the moving average coefficients and $\mu(s) = E(y_{Sn+s})$, $s=1,\dots, S$. Due to periodic stationarity, the moving average coefficient matrices $\Phi_{k}(s)$ are absolutely summable for all $s=1,\dots,S$. Following ursu2009modelling, they can be obtained recursively by the relations
The moving average coefficient matrices are periodic in the sense that $\Phi_{k}(Sn+s) = \Phi_{k}(s), \; n \in \mathbb{Z},\; k \in \mathbb{N}_{0}.$ Reduced-form impulse responses are defined as responses of the system $k$ periods after a one-time unit shock at time $Sn+s$ kilian2017structural, while moving average coefficient matrices can be interpreted as the effect of a one-time unit shock at time $Sn + s -k$ on $y_{Sn+s}$. For a standard VAR, the reduced-form impulse responses are given by its moving average coefficient matrices due to the time-invariant structure of VARs. However, in PVARs, this equality does not hold since the shocks' effect on the system of variables depends on the season the shock occurs. Hence, in order to get PVAR reduced-form impulse responses, denoted as $\Phi_{k}^{IR}(s)$, we need to shift the seasonal index of periodic moving average coefficient matrices $k$-steps further, that is \[ \Phi_{k}^{IR}(s) = \Phi_{k}(s+k)= \sum\limits_{j=1}^{k} A_{j}(s+k) \Phi_{k-j}(s+k-j) , \; \; \;k \in \mathbb{N}_{0} . \] For $i,j = 1,\dots, m$, the element in the $i$-th row and $j$-th column of $\Phi_{k}^{IR}(s)$ quantifies the effect of a (one-time) unit shock in the $j$-th component of $\epsilon_{Sn+s}$ on the $i$-th component of $y_{Sn+k}$, while all other components of $\epsilon_{Sn+s}$ are not shocked at all. Please note that the reduced-form errors are typically contemporaneously correlated with each other and that this contemporaneous correlation is not reflected in the reduced-form impulse responses. Section (ref) provides an overview of how this contemporaneous structure can be taken into account.
In this section, we provide asymptotic theory for a least squares estimation approach for the PVAR model parameters, which allows for flexible linear restrictions. In this regard, we extend the setup of ursu2009modelling and cover more general forms of linear restrictions, which also allow for constraints across seasons\footnote{In comparison to ursu2009modelling, we find a representation of the PVAR($\boldsymbol{p}$) process that allows for imposing general linear restrictions on the parameter vector. ursu2009modelling consider each season of the PVAR($\boldsymbol{p}$) process individually and are therefore only able to restrict the seasonal parameters separately such that no constraints across seasons can be imposed. For example, in our approach, the PVAR can be constrained in such a way that all autoregressive parameters are constant across seasons, and thus a VAR is estimated.}. In addition, it should be noted that unconstrained PVAR estimation does not work well in most applications due to pronounced overfitting issues. However, in view of Figures (ref) and (ref), it may be often not necessary to fit a full PVAR model, i.e., without any constraints across seasons, to explain the dependence structure of the data. Consequently, general linear restrictions play a crucial role in our setup. To address these issues, we derived a PVAR setup which is also able to collapse to a non-periodic VAR under suitable linear restrictions.
In order to estimate the PVAR parameters suitably under general linear restrictions, a multivariate least squares approach is used. We assume that we observe $N$ complete cycles consisting of $S$ seasons, that is, we have $y_1, \dots ,y_{SN}$. Additionally, for notational convenience, suppose that pre-sample values $y_{s^{\prime}},\dots,y_0$, where $s^{\prime} = \min \{1-p(1), \dots, S-p(S)\}$ are available. Then, the PVAR($\boldsymbol{p}$) process in (ref) can be rewritten as
where $Z, E, B$ and $X$ are defined as
which are of dimension $(m \times SN)$, $(m \times SN)$, $(m \times \sum_{s=1}^S (mp(s)+1))$ and $(\sum_{s=1}^S (mp(s)+1) \times SN)$, respectively. The periodic autoregressive parameters are included in $B$, while the entries $X_{n}(s)$ of $X$ are $(mp(s)+1)$-dimensional random vectors given by $X_{n}(s) = (1,y_{Sn+s-1}^{\prime}, \dots, y_{Sn+s-p(s)}^{\prime})^{\prime}$. Vectorizing the model equation (ref) yields
where $z= (y_1^{\prime} , \dots, y_{SN}^{\prime})^{\prime}$, $\beta = vec\{B\}$ and $e = vec\{E\} = (\epsilon_1^{\prime} , \dots, \epsilon_{SN}^{\prime})^{\prime}$, respectively. Note that $\beta$ is given by $\beta = (\beta^{\prime}(1), \dots, \beta^{\prime}(S))^{\prime}$ with $\beta(s) = vec\{(\nu(s),A_{1}(s), \dots, A_{p(s)}(s)) \}$ and that the covariance matrix $\Sigma_E $ of $e$ is block-diagonal with $\Sigma_E = I_N \otimes \Sigma_{\epsilon} $, where $\Sigma_{\epsilon} = \text{diag}\bigl( \Sigma_{\epsilon}(1),\Sigma_{\epsilon}(2) , \dots, \Sigma_{\epsilon}(S) \bigr) \in \mathbb{R}^{Sm \times Sm}$ is also block-diagonal.
ursu2009modelling and boubacar2023estimating use linear restrictions that can be imposed on parameters within one season and do not allow for restrictions across seasons. However, restrictions across seasons can play a crucial role, especially in macroeconomic applications, since they allow parameters to be constrained as constant across the seasons. In general, any set of linear constraints imposed on $\beta$ can be represented by
where $R$ is a known $(m \sum_{s=1}^S (mp(s)+1) \times M)$-dimensional matrix with rank $M$, $r$ is a $m \sum_{s=1}^S (mp(s)+1)$-dimensional known vector and $\gamma$ is an $M$-dimensional unconstrained parameter vector of interest. Substituting $ \beta = R\gamma + r$ into (ref) yields
where $\bold{z}_r =z -\{X^{\prime} \otimes \boldsymbol{I}_m\} r$. Hence, the least squares estimator $\widehat{\gamma}_{LS}$ of $\gamma$ is defined as the minimizer of the residual sum of squares $RSS(\gamma)$, where
Hence, according to lutkepohl2005new, the corresponding least squares estimator is given by $\widehat{\gamma} = \left[R^{\prime} \{X X^{\prime} \otimes \boldsymbol{I}_m \} R \right]^{-1} R^{\prime} \{ X \otimes \boldsymbol{I}_m \} \bold{z}_r$ and the restricted least squares estimator $\widehat{\beta}_{res}$ for $\beta$ is obtained by replacing $\gamma$ by $\widehat{\gamma}$ in $\beta = R \gamma + r$, that is
The unrestricted least-squares estimator $\widehat{\beta}_{LS}$ is obtained by setting $R = \boldsymbol{I}_{m \sum_{s=1}^S (mp(s)+1)}$ and $r=\boldsymbol{0}$ in (ref) leading to $\widehat{\beta}_{LS} = \left[ \{X X^{\prime} \otimes \boldsymbol{I}_m \} \right]^{-1} \{ X \otimes \boldsymbol{I}_m \} z$. On the other hand, if $p(s) = p$ for all $s=1,\dots,S$ and $R$ and $r$ are defined as $R =\boldsymbol{1}_S\otimes \boldsymbol{I}_{m (mp+1)} $ and $r = \boldsymbol{0}$, where $\boldsymbol{1}_l$ denotes the $l$-dimensional vector of ones, the seasonal PVAR($p$) estimators $\beta(s)$ are restricted to be the same across seasons. Once the periodic autoregressive parameters are estimated for all $s=1,\dots,S$, a natural candidate to estimate the periodic covariance matrices $\Sigma_{\epsilon}(s)$ is \[ \widehat{\Sigma}_{\epsilon}(s) = \frac{1}{N- k(s)} \sum\limits_{n=0}^{N-1} \widehat{\epsilon}_{Sn+s} \widehat{\epsilon}_{Sn+s}^{\prime}, \] where $k(s)= \lfloor M(s)/m \rfloor$, with $M(s)$ denoting the number of freely varying parameters in season $s$. In the following, let $\sigma(s) = vech\{ \Sigma_{\epsilon}(s)\}$ and denote by $\widehat{\sigma}(s) = vech\{\widehat{\Sigma}_{\epsilon} (s)\}$ its empirical counterpart, where the $vech$ operator stacks the entries on the lower-triangular part of a square matrix columnwise below each other. Further, let $\sigma = (\sigma(1)^{\prime}, \dots, \sigma(S)^{\prime})^{\prime}$ and $\widehat{\sigma} = (\widehat{\sigma}(1)^{\prime}, \dots, \widehat{\sigma}(S)^{\prime})^{\prime}$, which are both of dimension $S\widetilde{m}$ with $\widetilde{m} = m(m+1)/2$.
In order to derive the joint asymptotic properties of $\widehat{\beta}_{res}$ under general linear restrictions and periodic covariance matrix estimators $\widehat{\sigma}$, some notation is introduced. At first, for $s=1\dots, S$ and $k\in\mathbb{N}$, we define $mp(s) \times m$ matrices $C_k(s)=(\Phi_{k-1}^\prime(s-1), \ldots,\Phi_{k-p(s)}^\prime(s-p(s))^{\prime}$, where $\Phi_{k}(s)$ denote the moving average coefficient matrices in (ref). In addition, let $L_{mS}$ be the duplication matrix that satisfies $ (vech\{ B_1 \}^{\prime}, \dots, vech\{ B_S \}^{\prime} )^{\prime} = L_{mS} (vec\{ B_1 \}^{\prime}, \dots, vec\{ B_S \}^{\prime})^{\prime} $ for symmetric ($m\times m$) matrices $B_1, \dots,B_S$. Further, for $a,b,c \in \mathbb{Z}$ and $s_1,s_2=1,\ldots,S$, we define the ($m^2 \times m)$ and ($m^2 \times m^2$) dimensional matrices
and the $(Sm^2 \times Sm^2)$ matrices $\tau_{a,b,c}=(\tau_{a,b,c}(s_1,s_2))_{s_1, s_2 = 1, \dots, S}$. The symbols $\overset{d}{\to}$ and $\overset{p}{\to}$ denote convergence in distribution and in probability, respectively, and $ \mathcal{N}(\boldsymbol{\mu}, \boldsymbol{\Sigma})$ denotes a normal distribution with mean vector $\boldsymbol{\mu}$ and covariance matrix $\boldsymbol{\Sigma}$.
Note that some of the matrices in (ref) and (ref) simplify under an independence assumption imposed on the error process. In contrast to this case, also called strong PVAR, we aim to cover so-called weak PVARs as well that do explicitly not make any independence assumption, but allow for an uncorrelated, but possibly dependent periodic white noise process. To enable the derivation of asymptotic theory for weak PVARs, we impose the following assumptions on the process $\{y_{Sn+s}\}_{n \in \mathbb{Z}, s = 1,\dots,S}$.
Instead of the common i.i.d. assumption for the white noise process $\{w_{t}\}_{t\in\mathbb{Z}}$, we impose a less restrictive $\alpha$-mixing assumption. As demonstrated in our real-data application in Figure (ref) in Section (ref) below, where the (P)VAR residuals are likely to be non-linearly dependent, the derivation of asymptotic theory under this more general assumption of uncorrelated, but possibly dependent innovations, is important in practice. Following jentsch2022asymptotically, we impose a typical $\alpha$-mixing assumption, which implies absolute summability of (joint) cumulants up to order $4$. The resulting periodic weak white noise process $\{\epsilon_{Sn+s} \}_{n\in\mathbb{Z}, s=1,\dots,S}$ cover a large class of uncorrelated, but possibly dependent periodically (strictly) stationary processes and allows, e.g., for periodic conditional heteroscedasticity, see francq2011asymptotic, bibi2016periodic. While boubacar2023estimating already showed asymptotic normality in the weak PVAR case for least squares estimators with linear restrictions within a given season, we generalize their result and enable linearly constrained least squares estimators that also allow for restrictions across seasons. Additionally, we provide asymptotic results for covariance estimators. Under Assumption (ref), we obtain a joint CLT of $\widehat{\beta}_{res}$ and $\widehat{\sigma}$.
Due to the seasonal structure of the estimators $\widehat{\beta}_{res}$ and $\widehat{\sigma}$, the submatrices of $V$ are composed of further season-wise submatrices leading to a tedious and lengthy notation. In particular, the submatrices of $V$ depend on infinite sums of $\tau_{a,b,c}(s_1,s_2)$ and $\kappa_{a,b}(s_1,s_2)$, $s_1,s_2 = 1, \dots, S$ and $a,b,c \in \mathbb{Z}$ given in (ref) and (ref), which are cumbersome to estimate. The submatrix $V^{(2,1)} = V^{(1,2)\prime}$ is generally not the zero-matrix, meaning that the estimators $\widehat{\beta}_{res}$ and $\widehat{\sigma}$ become asymptotically dependent. We recommend bootstrap methods to be discussed in Section (ref) to approximate the limiting variances when conducting inference. If we impose stronger assumptions such as e.g. martingale difference sequences or i.i.d. assumptions, the submatrices of $V$ simplify a lot. Nevertheless, as demonstrated in bruggemann2016inference and jentsch2019dynamic, this will not allow for valid wild bootstrap inference. The simplifications in the i.i.d. case are discussed in Remark (ref) below. For the special case of $S=1$, where the weak (restricted) PVAR collapses to a weak (restricted) VAR, we obtain the same results as in francq2007multivariate and bruggemann2016inference.
In this section, we discuss structural analyses in PVAR setups and propose two novel identification methods for structural-form PVAR models (SPVARs). Further, we derive asymptotic properties of the periodic (structural) impulse responses.
As in the VAR setup, PVAR models do not capture contemporaneous effects between economic variables. This structure remains in the error term as (possibly) periodic contemporaneous correlation. However, in contrast to VAR setups, this periodic contemporaneous correlation in the PVAR error term consists not only of the periodic contemporaneous effects between the variables, but also of the periodic heteroscedastic effects of the variables. This raises the question of whether the periodic heteroscedasticity in the error component should be explained by the SPVAR model or whether it should remain in the structural shocks in the form that structural shocks occur more/less pronounced in some seasons than in others. To compare the periodic dynamics within the data, we generally propose to standardize the structural shocks across seasons, so that the periodic heteroscedastic structures of the variables are also captured by the SPVAR model. There is also the question of whether and to what extent the periodic contemporaneous correlations between the variables are explainable. However, as discussed in the following, a complete explanation of these correlations requires a very large number of restrictions, which may be too restrictive in practice.
In view of commonly used identification approaches used for structural VARs, a natural way to identify the PVAR system structurally is to simply generalize the SVAR identification by imposing restrictions for each period. That is, we impose the following assumptions on the variance-covariance matrices of the periodic white noise error term
where $H_{0}(s) \in \mathbb{R}^{m \times m}$ is invertible due to the non-singularity of $\Sigma_{\epsilon}(s)$. For each season $s = 1,\dots, S$, the system of equations consists of $m(m+1)/2$ independent equations and $m^2$ unknowns. Thus, (ref) is uniquely solvable, if at least $m(m-1)/2$ restrictions are imposed on $H_{0}(s)$ for each season $s=1,\dots, S$.\footnote{A trivial way to identify structural shocks is defining $H_{0}(s)$ as the lower-triangular Cholesky factor of $\Sigma_{\epsilon}(s)$, $s=1,\dots, S$. However, this recursive short-run identification method does not necessarily lead to economically interpretable structural shocks.} Transforming the reduced-form PVAR model (ref) by left-multiplying with the inverse SPVAR impact matrix $H^{-1}_{0}(s)$ yields
where $w_{Sn+s} = H_{0}^{-1}(s)\epsilon_{Sn+s}$ denotes the structural error satisfying $\Sigma_w(s) =Var(w_{Sn+s}) = Var(H_{0}^{-1}(s) \epsilon_{Sn+s}) = \boldsymbol{I}_m$ for all $s=1,\dots,S$. This means that, in line with Assumption (ref) (ii), the structural errors become a strictly stationary white noise process with the identity matrix as variance-covariance matrix. Further, in this setup, the SPVAR model in (ref) is indeed able to explain all the periodic contemporaneous effects and the periodic heteroscedasticity of the variables, since the endogenous variables are driven entirely by strictly stationary white noise shocks.
However, in comparison to the SPVAR setup, the usual SVAR identification techniques such as, e.g., short-run, long-run, proxy, sign or inequality restrictions (see kilian2017structural for an overview) are always based on $\widebar{\Sigma}_{\epsilon} = 1/S \sum_{s=1}^S \Sigma_{\epsilon}(s)$, which can be interpreted as the pooled (average) version of the variance-covariance matrices $\Sigma_{\epsilon}(s)$ in our periodic setup. Note that, to achieve full identification of (ref), identifying restrictions have to be imposed on each season $s=1,\dots ,S$. Hence, structural identification techniques on the PVAR based on (ref) can be more restrictive than the same identification schemes used in an SVAR setup. This is also made clear by the fact that the implicit assumption $Var(w_{Sn+s}) =\Sigma_w(s) = \boldsymbol{I}_m$ for all $s$ in (ref) implies the corresponding assumption in SVARs $ \widebar{\Sigma}_w = 1/S \sum_{s=1}^S \Sigma_{w}(s) = \boldsymbol{I}_m$, while the reverse direction is generally not true. Accordingly, it is questionable whether these harsher restrictions are appropriate in practical issues. However, it must be made clear that they have to be imposed to capture the entire periodic contemporaneous structure of the data by the SPVAR model.
The general idea behind our proposed structural identification methods for PVARs is based on the assumption that periodic contemporaneous structures cancel out on average over all seasons. That is, instead of imposing $Var(w_{Sn+s}) =\Sigma_w(s) = \boldsymbol{I}_m$ for all periods $s=1,\ldots,S$, we restrict (only) the pooled (average) variance-covariance matrix of the structural shocks to a unit matrix. That is, we impose that
holds. By construction $H_0 \in \mathbb{R}^{m \times m}$ is invertible and can be interpreted as the non-seasonal impact matrix of the structural shocks. Note that (ref) is the natural counterpart to the SVAR identification problem.
For structural analyses in PVAR setups, we propose two different identification techniques called full identification and approximate in-average identification, which are discussed in detail in the following. The full identification approach for SPVARs addresses both sources of periodic contemporaneous correlation, that is, the periodic contemporaneous effect as well as the heteroscedasticity of the variables. However, many restrictions have to be made to achieve full identification. The (approximate) in-average identification approach allows the structural shocks to have possibly fully occupied seasonal variance-covariance matrices, while their pooled variance-covariance matrix is (approximately) a unit matrix. Approximate in-average identification has the advantage over full identification that identification can be achieved with much fewer restrictions. However, it cannot explain periodic contemporaneous effects of the variables and these remain in the structural error term.
The most general SPVAR identification problem based on (ref) that is able to explain the complete periodic contemporaneous correlation of the variables can be formulated as follows
where $H_0 \in \mathbb{R}^{m \times m}$ and $W(s) \in \mathbb{R}^{m \times m}, s = 1,\dots, S$ are invertible. In (ref), we naturally impose $m(m+1)/2$ restrictions by $1/S \sum_{s=1}^S W(s) W(s)^{\prime} = \boldsymbol{I}_m$ to achieve (ref). Therefore, $W(s)$ can be interpreted as a seasonal matrix containing all the periodic heteroscedastic and contemporaneous effects around the mean level, while $H_0$ can be interpreted as the non-seasonal impact matrix. As we have $m^2+Sm^2$ parameters from $H_0$ and $W(s), s = 1,\dots, S$, additional $m(m-1)/2 + Sm(m-1)/2$ restrictions on $H_0$ and $W(s), s = 1,\dots,S$ are needed, in order to solve (ref) uniquely. Using $H_0(s) = H_0 W(s)$ as SPVAR impact matrix yields
where $\Sigma_w(s) = \boldsymbol{I}_m$ for all $s=1,\dots,S$. This means that both the periodic heteroscedasticity and periodic contemporaneous effects are fully captured by the SPVAR model. However, the issue of this approach is that due to the large number of restrictions that are needed to fully identify (ref), there is no way around restricting the periodic contemporaneous structure between structural shocks and variables, i.e., restricting $W(s),s=1,\dots,S$. Since the periodic contemporaneous structure between structural shocks and variables is unclear, it can be questionable to constrain this structure in any way. Trivial identification strategies such as setting $W(s)$ as the Cholesky factor of $W(s)W(s)^{\prime}$ or restricting $W(s)$ to be a symmetric matrix for all seasons $s=1,\dots,S$ combined with imposing $m(m-1)/2$ restrictions on $H_0$ solve the identification problem, but rely on highly questionable assumptions on the periodic contemporaneous structure of the data. For example, a Cholesky identification would assume that there are no periodic contemporaneous effects between some structural shocks and variables, while e.g. the assumption of a symmetric matrix $W(s)$ would mean that the periodic contemporaneous effects between the $j$-th shock and the $i$-th variable and between the $i$-th shock and the $j$-th variable are equal for $i \neq j$, $i,j = 1,\dots, m $. \\
Since restrictions on the periodic contemporaneous dependence structure between structural shocks and variables can be highly questionable, we propose an alternative identification method for SPVARs that is based on (ref). The idea behind this identification method is to identify seasonal structural shocks whose variances are normalized for all seasons, but whose contemporaneous effects may vary seasonally, but cancel out on average across all seasons. This alternative identification method is not as flexible as (ref), but can be identified without restricting the periodic contemporaneous dependence structure. The system of equations can be formulated by
where $\Lambda(s) \in \mathbb{R}^{m \times m}$ is a symmetric seasonal matrix. In total, we have $Sm(m+1)/2$ independent equations and $m^2+Sm(m+1)/2$ parameters from $H_0$ and $\Lambda(s), s = 1,\dots, S$. As $m(m+1)/2$ restrictions are naturally imposed by $1/S \sum_{s=1}^S \Lambda(s) = \boldsymbol{I}_m$ to achieve (ref), additional $m(m-1)/2$ restrictions are needed (as in SVAR identification), in order to solve (ref) uniquely\footnote{From a mathematical perspective, restrictions can also be imposed on the seasonal covariance matrix $\Lambda(s)$, but this is not really practicable as such restrictions are highly questionable. Accordingly, in most applications, we restrict the non-seasonal impact matrix $H_0$ exclusively to obtain identification.}. Using $\widetilde{H}_{0}(s) = H_0 \; diag[\Lambda(s)]^{1/2}$ as SPVAR transformations matrix, where $diag[\Lambda(s)]$ is a diagonal matrix containing the diagonal elements of $\Lambda(s)$, yields
The covariance matrix of the structural shocks is then given by $\Sigma_{\widetilde{w}}(s) = \Lambda_0(s)$ for all $s=1,\dots,S$ with $\Lambda_0(s) = diag[\Lambda(s)]^{-1/2} \; \Lambda(s) \; diag[\Lambda(s)]^{-1/2}$. We introduce the notation of $\widetilde{H}_0(s)$ and $\widetilde{w}_{Sn+s}$, because the structural shocks $\widetilde{w}_{Sn+s}$ do not necessarily correspond to those $w_{Sn+s}$ from Assumption (ref), due to $\Sigma_{\widetilde{w}}(s) = \Lambda_0(s) \neq \boldsymbol{I}_m$. The construction of $\widetilde{H}_0(s)$ ensures that the variances of $\widetilde{w}_{Sn+s}$ are normalized to one for all seasons. Accordingly $\Lambda_0(s)$ can be described as the seasonal correlation matrix of the structural shocks meaning that periodic heteroscedasticity of the PVAR shocks is completely captured by the SPVAR model (ref). The main advantage of (ref) compared to (ref) is that significantly fewer restrictions need to be imposed on $\widetilde{H}_{0}(s)$ in order to identify the SPVAR system. However, using this approach, periodic contemporaneous effects are explicitly not captured by the SPVAR in (ref) and in general, due to the variance normalization of $\widetilde{w}_{Sn+s}$, it is not true that the periodic correlations across the seasons cancel each other out completely. That is, because $1/S \sum_{s=1}^S \Lambda(s) = \boldsymbol{I}_m$ does not necessarily imply $\widebar{\Sigma}_{\widetilde{w}} = 1/S \sum_{s=1}^S \Lambda_0(s) = \boldsymbol{I}_m$. That is, with greater periodic heteroscedasticity and periodic contemporaneous effects, the off-diagonal elements of the pooled variance matrix $\bar{\Sigma}_{\widetilde{w}}$ of the structural shocks in (ref) deviate more from zero. Hence, if the reduced-form PVAR shocks either show weak periodic heteroscedasticity or weak periodic (contemporaneous) correlation, the identification method from (ref) is very suitable for structural identification since $\widebar{\Sigma}_{\widetilde{w}} \approx \boldsymbol{I}_m$ holds in this case. This means that the contemporaneous structure in the structural shocks is typically negligibly small such that the structural shocks can be considered orthogonal.
Now, we suppose that the structural PVAR impact matrix $H_0(s),s=1,\dots,S$ is fully identified. Structural impulse responses at lag $k$ are defined as a reaction of $y_{Sn+s+k}$ in response to a one-time impulse in $w_{Sn+s}$. Hence, following kilian2017structural, structural PVAR impulse responses for each season $s=1,\dots, S$ are defined as \[ \Theta^{SIR}_{k}(s) = \frac{\partial y_{Sn+s+k}}{\partial w_{Sn+s}} , \quad k\in\mathbb{N}_0.\] For $i,j = 1,\dots,m$, the element in the $i$-th row and $j$-th column of $\Theta^{SIR}_{k}(s)$ denotes the structural impulse response of the $i$-th component of $y_{Sn+s+k}$ to a unit shock in the $j$-th component of $w_{Sn+s}$ at season $s$. The structural impulse responses can also be deduced by transforming the moving average representation of $y_{Sn+s}$ as follows
where $\Theta_{k}(s) = \Phi_{k}(s)H_{0}(s-k)$ are the structural moving average coefficient matrices indicating the effect of $w_{Sn+s-k}$ on $y_{Sn + s}$. Note that $H_{0}(s) = H_{0}(Sn+s)$ for all $n \in \mathbb{Z}$ and $s=1,\dots,S$. The structural impulse responses $\Theta^{SIR}_{k}(s)$ can be expressed analytically as \[ \Theta^{SIR}_{k}(s) = \Theta_{k}(s+k) = \Phi_{k}(s+k)H_{0}(s) = \Phi^{IR}_{k}(s) H_{0}(s). \] In writing down the asymptotic properties of the impulse responses $\Phi_{k}^{IR}(s)$ and the structural impulse responses $\Theta_{k}^{SIR}(s)$, $s =1,\dots, S$, we use the higher-dimensional VAR representation of the PVAR which is given in Appendix C. It can be easily shown that the reduced-form impulse responses of the higher-dimensional VAR form, denoted by $\Pi_{h}^{IR}$, $h \in \mathbb{N}_0$, in Appendix C, consist of the PVAR impulse responses $\Phi_{k}^{IR}(s)$ for all $k,s$ satisfying $hS < k+s \leq (h+1)S$. The structural impulse responses of the higher-dimensional VAR form, denoted by $\Psi_{h}^{SIR}$, $h \in \mathbb{N}_0$, consist of $\Theta_{k}^{SIR}(s)$ for all $k,s$ that satisfy $hS < k+s \leq (h+1)S$. Precisely, we have \[ \Psi_{h}^{SIR} = \Pi_h^{IR}\boldsymbol{H}_0=
, \; \; \; h\in\mathbb{N}_0, \] where $\boldsymbol{H}_0 = \text{diag}\bigl( H_0(1),H_0(2), \dots, H_0(S) \bigr) \in \mathbb{R}^{Sm \times Sm}$. Let $\widehat{\Pi}_{h}^{IR}$ and $\widehat{\Psi}_{h}^{SIR}$ be the empirical counterparts of $\Pi_{h}^{IR}$ and $\Psi_{h}^{SIR}$ for $h\in\mathbb{N}_0$, respectively. Note that $\Pi_{h}^{IR}$, $h \in \mathbb{N}_0$ can be represented as a continuously differentiable function of the PVAR coefficients $\beta$, while the structural impulse responses $\Psi_{h}^{SIR}$ can be represented as a continuously differentiable function depending on $\beta$ and $\sigma$. Consequently, from Theorem (ref) and following lutkepohl2005new, Proposition 3.6, we immediately get the following result.
The closed-form solutions of the derivatives $G_{h}, F_{h}$ and $D_{h}$ are stated in Appendix C. In the strong PVAR case, the same statements as in Theorem (ref) hold with the matrices $V^{(1,1)}, V^{(2,1)}$ and $V^{(2,2)}$ replaced by $V^{(1,1)}_{iid}, V^{(2,1)}_{iid}$ and $V^{(2,2)}_{iid}$, respectively. The i.i.d. limiting variance matrices can be found in Appendix B.
Analytical expressions of the limiting variances of the PVAR estimators and the (structural) impulse responses are very complex and cumbersome to estimate. Therefore, we propose residual-based bootstrap methods to approximate the (limiting) distributions of the SPVAR estimators and prove their bootstrap consistency. We also construct a bootstrap-based test for seasonality in SPVAR impulse responses.
In this section, we discuss different residual-based bootstrap methods that are suitable for conducting inference in different PVAR setups. In Section (ref), we describe residual-based bootstrap schemes that are tailored for the assumption of a periodic weak white noise, while Section (ref) contains simplified versions that are sufficient for periodic strong white noise.
We describe residual-based seasonal block bootstrap methods for PVAR processes that work under Assumption (ref), where the PVAR error process $\{\epsilon_{Sn+s} \}_{n\in\mathbb{Z}, s = 1,\dots,S}$ is assumed to be periodic, possibly weak white noise. In this context, block bootstrapping is required to mimic the potential non-linear serial dependence, while a seasonal version is required because $\{\epsilon_{Sn+s} \}_{n\in\mathbb{Z}, s = 1,\dots,S}$ is allowed to be periodically correlated. Based on sample values $y_1,\dots, y_{SN}$ and pre-sample values $y_{s^{\prime}},\dots,y_0$, suppose that estimates of the PVAR($\boldsymbol{p})$ parameters $\beta$ and $\sigma$ are given under general linear constraints $\beta= R\gamma + r$ as stated in Section (ref).
The algorithm is as follows. We initialize the algorithm by choosing a positive integer block size $b < SN$ and define $l = \lceil SN/b \rceil$, where $\lceil\cdot\rceil$ denotes the ceiling function. Next, we define $(m \times b)$ blocks $\widehat{\mathcal{E}}_i = (\widehat{\epsilon}_i, \dots, \widehat{\epsilon}_{i+b-1})$ for $i= 1,\dots, SN-b+1$. Then, we proceed as follows:
In the bootstrap algorithm, a seasonal block bootstrap proposed by dudek2014generalized,dudek2016generalized is applied directly to the PVAR residuals. Alternatively, under Assumption (ref), a non-seasonal resampling variant can also be used. To initialize this bootstrap, we choose a positive integer block size $b < SN$ and define $l = \lceil SN/b \rceil$. Next, we define $(m \times b)$ blocks $\widehat{\mathcal{U}}_i = (\widehat{u}_i, \dots, \widehat{u}_{i+b-1})$, $i= 1,\dots, SN-b+1$, where $\widehat{u}_{Sn+s} = \widehat{\Sigma}_{\epsilon}^{-1/2}(s) \widehat{\epsilon}_{Sn+s}$ denote the standardized residuals for $s= 1,\dots,S$ and $n = 0,\dots,N-1$. Then, the non-seasonal variant is obtained by replacing steps $(1)-(2)$ by $(1^{\prime})-(2^{\prime})$, where
In this alternative bootstrap method, the periodic residuals are first standardized and then a classical moving-block bootstrap approach, as described, e.g., in bruggemann2016inference and jentsch2022asymptotically, is applied to the adequately standardized residuals. As both bootstrap algorithms need a suitably chosen block length $b$ as initialization, we refer to nordman2009note and bertail2024optimal, who developed optimal choices for the block length of non-seasonal and seasonal block bootstrap methods, respectively. With $b=1$, both block bootstrap methods collapse to the residual-based independent bootstraps to be discussed in Section (ref).
By repeating the bootstrap schemes $L$ times, where $L$ is large, we can get bootstrap quantiles for all estimators $\widehat \beta_{res}$, $\widehat \sigma$, $\widehat \Pi_h^{IR}$ and $\widehat \Psi_h^{SIR}$, which enables in particular the construction of (standard percentile) confidence intervals for the structural impulse responses.
If Assumption (ref) with $(ii)$ and $(iii)$ replaced by $(ii)^\prime$ and $(iii)^\prime$ holds, where the PVAR innovation process $\{\epsilon_{Sn+s} \}_{n\in\mathbb{Z}, s = 1,\dots,S}$ is assumed to be periodic strong white noise, simplified residual-based independent bootstrap techniques can be used. In Section (ref), by setting $b=1$ in steps $(1)-(2)$, the seasonal version of the residual-based independent bootstrap is obtained, while the non-seasonal version is obtained by setting $b=1$ in steps $(1^{\prime})-(2^{\prime})$.
The residual-based independent bootstrap methods provide valid confidence intervals only when Assumption (ref) with $(ii)$ and $(iii)$ replaced by $(ii)^\prime$ and $(iii)^\prime$ hold. However, as illustrated in Figure (ref), the periodic strong white noise assumptions may be violated in applications. Hence, we recommend to use residual-based block bootstrap methods for inference of periodically correlated time series.
We show bootstrap consistency for the residual-based seasonal block bootstrap methods discussed in Sections (ref) and (ref) under periodic weak and strong white noise assumptions, respectively. Addressing first the weak PVAR case, we need the following assumption.
For weak white noise processes, such an assumption is made, e.g., in gonccalves2007asymptotic to guarantee summability of cumulants up to order eight to prove consistency of wild and pairwise bootstrap methods for AR($\infty$) processes, in bruggemann2016inference to prove consistency of a residual-based moving block bootstrap for VAR($p$) models, and in jentsch2019dynamic,jentsch2022asymptotically to prove consistency of a residual-based moving-block bootstrap for Proxy SVARs. We now show that the residual-based seasonal block bootstrap can approximate the limiting distributions of $\sqrt{N}\bigl( \bigl(\widehat{\beta}_{res}- \beta\bigr)^{\prime}, \bigl(\widehat{\sigma}- \sigma \bigr)^{\prime}\bigr)^{\prime}$, $\sqrt{N}\bigl( vec \{ \widehat{\Pi}_{h}^{IR} - \Pi_{h}^{IR} \}\bigr),$ and $\sqrt{N} \bigl(vec \{ \widehat{\Psi}_{h}^{SIR} - \Psi_{h}^{SIR} \}\bigr),$ $h\in \mathbb{N}_0$ derived in Theorems (ref) and (ref), respectively.
The proof of Theorem (ref) can be found in Appendix D. It is worth noting that in the case, where the more restrictive periodic strong white noise assumption applies, the independent bootstrap, i.e., the block bootstrap with block length $b=1$, is sufficient.
In this section, we derive a bootstrap-based test for (seasonal) differences in the SPVAR and SVAR impulse response functions. This approach allows to test for the null hypothesis that there is no seasonal difference between SPVAR and SVAR impulse response functions, while the alternative finds a difference in at least one season. Formally, for fixed $i,j = 1,\dots, m$, we consider the null hypothesis
where $\Theta^{SIR}_{ij}(s)$ is the $(n_{IR}+1)$-dimensional vector consisting of seasonal impulse responses up to lag $n_{IR}$ of the $i$-th variable to a unit shock in the $j$-th structural shock, that is, $\Theta^{SIR}_{ij}(s) = (\Theta^{SIR}_{0,ij}(s), \dots, \Theta^{SIR}_{n_{IR},ij}(s))^{\prime}$. $\Theta^{VAR}_{ij}$ is the non-seasonal counterpart from the SVAR model on seasonally adjusted data. Then, given seasonally adjusted data $y_1^{SA},\ldots,y_{SN}^{SA}$ and seasonally unadjusted data $y_1^{SU},\ldots,y_{SN}^{SU}$, the test statistic $\widehat Q_{N}(i,j)$ is defined by the sum of squared differences in the estimated SPVAR and SVAR impulse response functions, that is,
We use a bootstrap-based method to mimic the distribution of $\widehat{Q}_N(i,j)$ under the null hypothesis $K_0$ in (ref). For this purpose, the bootstrap algorithm has to address the following points. Firstly, $\widehat{\Theta}^{SIR}_{ij}(s)$ and $\widehat{\Theta}^{VAR}_{ij}$ are generally not independent to each other and secondly, $\widehat{\Theta}^{SIR}_{ij}(s)$ is determined from seasonally unadjusted data, while $\widehat{\Theta}^{VAR}_{ij}$ is determined from seasonally adjusted data. Hence, the impulse response functions $\widehat{\Theta}^{SIR}_{ij}(s)$ and $\widehat{\Theta}^{VAR}_{ij}$ are therefore not only determined by two different models, but also originate from different data sets that are obviously highly interdependent. Accordingly, the bootstrap can only mimic the distribution of $\widehat{Q}_N(i,j)$ properly under the null if this dependence between the seasonally adjusted and unadjusted data is adequately captured.
In order to take these dependencies into account, we propose to generate bootstrap pseudo observations of the seasonally adjusted and unadjusted data by the following steps. First, we initialize the algorithm by choosing a VAR model order $p_{var}$ as the reduced-form VAR model and a (possibly) restricted PVAR of order $p = \max\{p(1), \dots, p(S)\}$ as the reduced-form PVAR model. If the PVAR model order varies across the seasons, note that we impose zero restrictions on $A_{j}(s)$ for $p(s) < j \leq p$. Further, we suppose that their structural model variants are identified using the same identification strategies. We also suppose that $p_{var}$ starting values for seasonally adjusted data and $p$ starting values for seasonally unadjusted data are available. Then, the algorithm works as follows:
In Step (1), we model the seasonally adjusted and unadjusted data independently, which means that the dependence structure between the data is preserved in the residuals. In Step (2), we calculate the pooled autoregressive estimates and the pooled impact matrix, which ensure that the generated pseudo observations in Step (4) come from a model that satisfies the null hypothesis $K_0$. In order to mimic the dependence between the SPVAR and SVAR model estimates, we merge the residuals from both models in Step (3) and apply a seasonal block bootstrap that is able to capture this dependence structure as well as the potential seasonal structure in the SPVAR residuals. Step (4) produces the bootstrap pseudo samples of seasonally adjusted and unadjusted data. For both samples, the same pooled autoregressive coefficients and impact matrix are used to ensure that the true seasonal and non-seasonal impulse responses are the same. Note that the pooled autoregressive coefficients are computed by the mean of the VAR estimates and the pooled PVAR estimates. To calculate the bootstrap test statistic $\widehat{Q}_N^*(i,j)$ under the null, in Step (5), we first fit (separately) the SVAR($p_{var}$) model to $y^{SA*}$ and the SPVAR($p$) model to $y^{SU*}$. Then, we obtain the bootstrap versions of the corresponding impulse responses by identifying both structural models by the same restrictions as used in Step (1). We calculate the bootstrap test statistics $\widehat{Q}_N^{*(l)}(i,j)$, $l=1,\ldots,L$ by inserting the bootstrap versions of the impulse responses and we reject $K_0$ if the test statistic $\widehat{Q}_N(i,j)$ exceeds the critical value $\widehat{q}^*_{1-\alpha}$.
We also derive a bootstrap-based test for the statistical significance of differences in SPVAR and SVAR impulse response functions, both based on seasonally unadjusted data. The null of this test can be formulated by $\widetilde K_0: \Theta_{ij}^{SIR}(s) = \widebar{\Theta}^{SIR}_{ij}, \; \forall s=1,\dots,S$, where $\widebar{\Theta}^{SIR}_{ij}$ is the impulse response function from the SVAR($p$) with seasonal means. The test algorithm can be found in Appendix E.
In this section, we compare the results of applying a structural PVAR to seasonally unadjusted macroeconomic data with those of applying a structural VAR to seasonally adjusted data. The variables included in the model are monthly industrial production (IP) as a measure of real economic activity, CPI inflation (INF) and federal funds rate (FFR) from Jan 1968 - Dec 2019 taken from the FRED database. The ADF test rejects unit roots in INF and the FFR, while it does not reject a unit root in industrial production. Hence, we address this by taking log differences of industrial production, whereas INF and FFR remain untouched. Furthermore, HEGY tests hylleberg1990seasonal reject the null of seasonal unit roots for different frequencies in the variables. This result is also confirmed by the OCSB test osborn1988seasontest, which also rejects the presence of seasonal unit roots. For structural PVAR analysis, IP and INF are considered on a seasonally unadjusted basis, while their seasonally adjusted counterparts are considered for SVAR analysis. As FFR appears to be non-periodic, we use it for both analyses.
Due to the monthly frequency, an obvious choice for the number of seasons is $S =12$, resulting in $N=52$ cycles (years). The choice $S =12$ is strongly supported by the estimates of the SDs and ACFs of IP and INF given in Figure (ref). As demonstrated in Figure (ref) and (ref) in the introductory section, we find that IP and INF exhibit strong periodic patterns, while FFR has no periodic structure at all. The periodic structure of IP cannot be eliminated by seasonal demeaning alone, while the periodic structure of INF can be explained mainly by seasonal variation in the mean. Consequently, seasonal demeaning appears to be not sufficient to remove the entire seasonal component of IP, while in the case of INF the periodic structure is almost eliminated completely.
On the seasonally adjusted data basis, we use a VAR model of order $p= 9$ as the benchmark reduced-form model. The model order is chosen based on BIC. Based on Figures (ref) and (ref), we allow the seasonal intercepts of IP and INF to be periodically varying over the seasons and set the intercepts of FFR to be non-periodic within our PVAR setup. For choosing a suitable model order for the reduced-form PVAR, we set the model order of the PVAR for each season to be the same. For the sake of simplicity, we then choose the model order that minimizes the seasonal BIC for all PVARs with seasonal means and non-periodic autoregressive effects\footnote{Since we do not know which autoregressive coefficients are periodic and which are not without knowing the model order, we proceed in this approach by determining the model order from VAR models with seasonal intercept. Alternatively, the bottom-up approach can be run for different PVAR orders, with the selected model determined by the seasonal BIC.}. Based on BIC, the PVAR model order is selected as $p(s)=12$ for $s=1,\dots,S = 12$.
In order to capture the periodic structure, which goes beyond seasonality in means, and to choose the linear restrictions of the reduced-form PVAR model, we use a slightly more sophisticated bottom-up strategy. We define $\widehat{\beta}^{(PVAR)}_i = \bigl(vec \{ \bigl(\widehat{A}_{1}(s), \dots, \widehat{A}_{12}(s) \bigr) \}_i, s = 1,2, \dots, S = 12 \bigr)^{\prime} $ as the periodic $S$-dimensional vector of the $i$-th autoregressive parameter for $i = 1, \dots, pm^2 = 108$, where $vec \{ \bigl(\widehat{A}_{1}(s), \dots, \widehat{A}_{12}(s) \bigr) \}_i$ is the $i$-th component of $vec \{ \bigl(\widehat{A}_{1}(s), \dots, \widehat{A}_{12}(s) \bigr) \}$. $\widebar{\widehat{\beta}}^{(PVAR)}_i$ is defined as the mean of $\widehat{\beta}^{(PVAR)}_i$ for $i=1,\dots,pm^2$. Then, the algorithm for the bottom-up approach is as follows:
In order to get $\widehat{\Sigma}_{\widehat{\beta}^{(PVAR)}_i}^{-1}$ for each PVAR model in Step (1) above, we propose to use residual-based seasonal block bootstrap procedures from Section (ref). Note that, in the first iteration of the bottom-up strategy, we begin by using a VAR with seasonal intercepts, before gradually allowing the model to incorporate periodically varying autoregressive parameters at the most periodic points. We proceed this way until the most periodic component is no longer significant. Note that in each step, assuming that the initial PVAR model is the true model, the Wald-type statistic divided by $S$ is asymptotically $\chi^2_{1}$-distributed.
In general, we propose $\alpha \in [0.05,0.10]$. Note that too small $\alpha$ can lead to periodically inflexible models, while too large $\alpha$ can lead to models that overfit the periodic structure. In our real data analysis, we use the level $\alpha = 0.075$ and stop after four iterations of the bottom-up strategy. The resulting seasonally varying autoregressive parameters are the one-month delayed effects of IP on INF and on IP, the five-month delayed effects of IP on IP and the twelve-month delayed effects of IP on IP. All the remaining autoregressive coefficients are restricted to be constant. This finding is in line with the findings from Figures (ref) and (ref) that IP is the variable in our system that contains most of the periodic structure, which goes beyond seasonality in means, while INF is slightly periodic and FFR is non-periodic.
Consequently, the reduced-form PVAR(12) with period $S=12$, where $y^{ip}, \pi$ and $i$ denote IP, INF and FFR, for $s = 1,\dots, 12$ and $n = 0,\dots ,51$, is as follows
where $y_{12n+s} = (y^{ip}_{12n+s}, \pi_{12n+s}, i_{12n+s})^{\prime}$, $\nu(s) = (\nu^{y^{ip}}(s), \nu^{\pi}(s), \nu^{i})^{\prime}$, $\mathcal{A} = \{ 1,5,12 \}$, $\mathcal{A}^c = \{1,\dots,12 \} \setminus \mathcal{A}$ and
To investigate the validity of the imposed restrictions, we compare the residuals resulting from fitting the restricted PVAR(12) from (ref) with those of a VAR(12) with seasonal intercept. From Figures (ref) and (ref), it can be clearly seen that in particular the IP residuals from the PVAR have much less periodic structure than those from the VAR. Thus, the PVAR model (ref) is able to capture the periodicity of the data better than the VAR model with seasonal intercept, especially for IP, but also somewhat for INF. Note that the FFR is modeled in the same way in both reduced-form models, namely not periodically. Further, based on the seasonal BIC, the PVAR is preferred over the VAR with seasonal mean.
In our analysis, we use full and approximate in-average identification to identify unit variance aggregate demand $(ad)$, aggregate supply $(as)$ and monetary policy $(mp)$ shocks. In both identification methods, a mixture of short- and long-run restrictions following shapiro1988sources, gali1992well and rubio2010structural are used to identify structural $ad, as$ and $mp$ shocks.
To identify the shocks by the full identification approach, we assume that $ad$ and $mp$ shocks do not have an effect on IP in the long-run for all seasons. Additionally, we impose that $mp$ shocks do not have contemporaneous effects on IP in all seasons. These restrictions result in a total of $Sm(m-1)/2 = 36$ restrictions in identification scheme (ref), with two long-run and one short-run restriction applying to all $S$ seasons. As (ref) is a special case of (ref), these restrictions can also be formulated as restrictions on $H_0$ and $W(s),s=1,\dots, S$ in (ref). The short-run restriction is imposed as described in Remark (ref) by setting the corresponding static effect in $H_0$ to zero and by setting the cumulative contemporaneous effects of all other structural shocks to zero for all $S$ seasons. Each long-run restriction is imposed by setting the corresponding element of the pooled (average) long-run impulse response \[ \frac{1}{S} \sum_{s=1}^S \Biggl[\Bigl(\sum_{k=0}^{\infty}\Phi^{IR}_{k}(s) \Bigr) H_{0} \Biggl]\] to zero and by setting the corresponding element of the periodic long-run impulse response \[ \Biggl[\Bigl(\sum_{k=0}^{\infty}\Phi^{IR}_{k}(s) \Bigr) H_{0} W(s) \Biggl]\] to zero for each season $s=1,\dots,S$. Consequently, each short- and long-run restriction results in $S+1$ restrictions on (ref), totalling $(S+1)m(m-1)/2 = 39$ restrictions, which are needed to fully identify the system.
To achieve approximate in-average identification, we only impose the above restrictions on average, rather than on each season individually. Based on (ref), we identify the non-seasonal impact matrix $H_0$ and thus also $\widetilde{H}_0(s) = H_0 \; diag[\Lambda(s)]$ for $s = 1,\dots, S$ by imposing the short-run restriction on $H_0$ and the long-run restrictions on the pooled long-run impulse response $\widebar{\Theta}^{LR}$ given in Remark (ref). Consequently, $m(m-1)/2 = 3$ restrictions are imposed in (ref) which are needed to fully identify the system.
In the following, seasonal structural impulse responses using the restricted SPVAR(12) on seasonally unadjusted data and structural impulse responses using an SVAR(9) on seasonally adjusted data are analyzed, in which the aggregate demand $(ad)$, aggregate supply $(as)$ and monetary policy shocks $(mp)$ are identified. To identify the structural $ad$, $as$ and $mp$ shocks, we use a mixture of short- and long-run restrictions discussed in Section (ref). In this Section we only compare SPVAR impulse responses identified by the full identification approach with SVAR impulse responses. Impulse responses identified by approximate in-average identification method can be found in Appendix F. SVAR impulse responses are identified using the same short- and long-run restrictions.
The black bold line indicates the structural impulse response of the SVAR on seasonally adjusted data which can also be used as the benchmark, while each of the red lines indicates a structural impulse response generated by the SPVAR for the $S=12$ seasons, i.e., from January to December. The dashed lines give the (pointwise) 68% confidence intervals of the structural impulse responses of the SPVAR (in red) and the SVAR (in black).
For the construction of confidence intervals of the seasonal impulse responses, we use the non-seasonal variant of the residual-based seasonal block bootstrap described in Section (ref), because the bottom row in Figure (ref) shows significant autoregressive structure in the squared structural shocks such that the strong periodic white noise assumption appears to not hold. The block length of the residual-based seasonal block bootstrap is set to $b=5$, which is in line with the choice in jentsch2022asymptotically, who considered a similar data generating process, and with the results in bertail2024optimal. We use standard percentile intervals for the construction of the confidence intervals. For the SVAR analysis, we use a non-seasonal block bootstrap to construct the confidence intervals. $L = 1000$ bootstrap repetitions are used for each bootstrap method.
Figures (ref), (ref), and (ref) present the seasonal structural impulse responses generated by a monetary policy shock. This shock is normalized as an unexpected increase in FFR of 25 basis points in each season. The figures show the effects on IP, INF and FFR, respectively. The impulse response functions generated by positive demand and negative supply shocks are deferred to the Appendix G. It appears that monetary policy shocks do not have significant periodic effects on the macro variables, suggesting that the timing of unexpected rate hikes or cuts by the Fed does not affect the macro variables. However, this may also be because the short- and long-run restrictions of the $mp$ shock to IP naturally limit the periodic flexibility of the impulse response. It is noticeable that the confidence bands for the SPVAR responses, especially for the INF and FFR impulse responses, are narrower compared to the non-seasonal SVAR impulse responses. Accordingly, it appears that we can make more precise statements on the basis of the SPVAR impulse responses in the individual seasons than at the non-seasonal SVAR level.
In terms of periodic patterns, the structural impulse responses of macro variables triggered by positive demand or negative supply shocks exhibit significantly more periodicity, particularly in May and September. After a positive demand shock, we usually observe a strong upward shift in September, while observing an downward shift in May. The opposite pattern is observed after negative supply shocks.
All in all, we find that the SVAR impulse responses are approximately given by the mean of the seasonal structural impulse responses. Interestingly, we do not see pronounced efficiency losses in the seasonal structural impulse responses that would materialize in terms of considerably increased confidence intervals compared to the SVAR impulse responses. Instead, we even notice some efficiency gains in certain seasons resulting in more accurate statements; see, e.g., the December panels in Figures G.1 - G.6 in Appendix G.
We also perform the bootstrap-based test for the statistical significance of the seasonal differences in the SPVAR impulse response functions described in Section (ref). Using the test statistic $\widehat Q_N(i,j)$, we test the null hypothesis $K_0: \Theta_{ij}^{SIR}(s) = \Theta^{VAR}_{ij}, \text{ for all } s = 1,\dots, S$, i.e., whether there are significant differences in the SPVAR impulse response function $\Theta_{ij}^{SIR}(s)$ and the SVAR impulse response function $\Theta^{VAR}_{ij}$ from seasonally adjusted data for $i = y^{ip}, \pi, i$ and $ j = ad, as, mp$. In order to calculate the critical value under the null, we use the reduced-form PVAR model from Section (ref), the VAR model order $p_{var} = 9$ , the block length $b = 5$ and $L = 1000$ bootstrap repetitions. Both, SPVAR and SVAR models, are identified using a mixture of short- and long-run restrictions discussed in Section (ref). We use the full identification approach for SPVAR identification.
Table (ref) shows the $p$-values of the test for seasonality of the impulse responses for different maximum number of lags $n_{IR}$. Overall, the test shows signficant periodicities in the impulse responses following demand and supply shocks, while monetary policy shocks do not appear to have any significant periodic effects on the macro variables. It should be noted that the contemporaneous effects of a monetary policy shock on the FFR are standardized for all seasons. However, as shown in Table (ref), the impulse responses of the variables following monetary policy shocks are not seasonally significant even without standardization.
In Monte Carlo simulations, we examine the seasonal impulse response intervals constructed by the residual-based seasonal bootstrap methods discussed in Section (ref). For this purpose, we generate seasonal data from the restricted structural PVAR that was obtained in the real data application when modeling the dynamics of IP, INF and FFR. Hence, in the simulation study, we generate a 3-dimensional sample $y_1, \dots, y_{12N}$ by
for $s = 1,\dots, 12$ and $n = 0,\dots ,N-1$, where we use the estimated reduced-form PVAR(12) model described in (ref) and the estimated SPVAR impact matrix identified by the full identification approach as the true DGP. We generate the structural shocks $\{ w_{t} \}_{t \in \mathbb{Z}}$ using different DGPs, where one variant relies on i.i.d. errors and three variants make use of differently parametrized GARCH error processes to cover both, the periodic strong and weak white noise cases. Using $v_t = (v_{1,t}, v_{2,t}, v_{3,t})^{\prime} \sim i.i.d. \; \mathcal{N}(\boldsymbol{0}, \boldsymbol{I}_3)$ for all $t \in \mathbb{Z}$, we define $w_{i,t} = \sigma_{i,t} v_{i,t}$ with $\sigma^2_{i,t} = a_0 + a_1 w_{i,t}^2 + b_1 \sigma^2_{i,t-1}, i = 1,2,3$ and $a_0 = 1 - a_1 - b_1$ implying that the components of the structural shocks follow univariate independent GARCH(1,1) processes with $E(w_{i,t}^2) = 1,i = 1,2,3$. We use the following GARCH(1,1) specifications of $\{w_t\}_{t\in\mathbb{Z}}$ for simulations:
where the case G0 represents the periodic strong white noise case, and (G1)-(G3) represent variants of the periodic weak white noise case. All these cases are also considered in bruggemann2016inference.
We simulate a total of $M=1000$ time series of length $N = 20,50,100$. Further, we bootstrap the structural impulse responses by using the non-seasonal variant of the residual-based seasonal block bootstrap for PVAR processes discussed in Section (ref) with block sizes $b \in \{1,3\}$ for $N=20$, $b\in\{1,5\}$ for $N=50$ and $b\in\{1,7\}$ for $N=100$ to obtain approximations of the confidence intervals. The block sizes are adopted from jentsch2022asymptotically, who considered similar block and sample size ratios. Note that the block size has to increase with increasing sample size and that for $b=1$, the residual-based seasonal block bootstrap collapses to a seasonal independent variant. The confidence intervals are constructed by standard percentile intervals. The nominal coverage rate is 90% and we use $L = 1000$ bootstrap repetitions.
Tables (ref) and (ref) show the simulated coverage rates for lag $k = 12$ and $N \in \{50,100 \}$ for GARCH specifications G0 and G3. The tables for lag $k = 0$, the tables for sample size $N=20$ and the tables for the other GARCH specifications G1 and G2 can be found in Appendix H, respectively. Due to a simplified presentation, we only provide coverage rates for lags $k \in\{0, 12\}$. For other lags, the coverage rate is not particularly different. Generally, it can be clearly seen that the simulated coverage rates approach the nominal coverage rate of 90% for increasing $N$. Table (ref) shows small efficiency losses when using the seasonal block bootstrap compared to its independent variant ($b=1$). This is to be expected, as the i.i.d. case G0 contains the assumption of strong periodic white noise. When we consider Table (ref), the simulated coverage rates approach the nominal coverage rate somewhat more closely when using the seasonal block bootstrap. We also find this result for cases G1 and G2. This is because the seasonal independent bootstrap (for $b=1$) is not able to capture the non-linear dependencies of the GARCH errors, while its block variant is able to do so. All in all, the standard percentile intervals hit the nominal coverage very well. In the case of a rather small data basis consisting of $N = 20$ years, the coverage rates look moderate and deviate quite a bit from the nominal coverage for some entries. For $N = 50, 100$, the coverage rates improve significantly and come very close to nominal coverage.
When dealing with seasonal data, deseasonalization using standard seasonal adjustment procedures is the standard methodology nowadays. This is the case even though seasonal adjustment procedures may distort the structure of the data and remove useful information. In this paper, we offer an alternative direct approach to analyze seasonally unadjusted macroeconomic data by structural models. Instead of first seasonally adjusting the data and then fitting a structural VAR model, we propose to fit a structural periodic VAR directly on raw macroeconomic data that is seasonally unadjusted. We provide a PVAR representation that enables linearly restricted estimation of PVAR models under more general linear constraints than ursu2009modelling, which allows to impose also constraints across seasons. In addition, different identification methods for structural PVAR models are proposed, which, compared to the identification in SVARs, also take the (contemporaneous) seasonality of the data into account. Further, we show consistency and asymptotic normality of (constrained) least squares estimators for PVAR coefficients and structural impulse responses under periodic weak white noise assumption. For standard error estimation and confidence interval construction, we introduce different residual-based bootstrap methods for PVARs and prove their bootstrap consistency. We also discuss bootstrap-based statistical testing for seasonality in impulse responses.
Our empirical results illustrate the analysis of seasonally unadjusted industrial production, CPI inflation and federal funds rate using a linearly constrained structural PVAR. For the structural analysis of the data, we identify three structural shocks, i.e., monetary policy shock, aggregate supply and demand shock based on a mixture of short- and long-run restrictions. We determine seasonal structural impulse responses based on the linearly restricted structural PVAR and compare the results with those from a standard SVAR analysis on the same US macro variables in seasonally adjusted form. We find that monetary policy shocks do not seem to have significant periodic effects on the macroeconomic variables. However, we observe significant periodic effects caused by aggregate supply and demand shocks.
Moreover, we do not find pronounced efficiency losses of suitably constrained structural PVARs compared to seasonally adjusted SVARs, but even efficiency gains in some individual seasons. In total, we clearly see that useful insights into the dynamics of the macro variables are lost, if we seasonally adjust the data first. In the simulation study, we show that the simulated coverage rates constructed by the proposed residual-based bootstrap methods for PVARs closely match the nominal coverage rate, even for realistic sample sizes.
Financial support by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation; Project-ID 520388526; TRR 391: Spatio-temporal Statistics for the Transition of Energy and Transport) and by the Mercator Research Center Ruhr (MERCUR) with project number Pe-2019-0044 is gratefully acknowledged.