1. Introduction
Time series data from several subjects that can be classified into groups according to a set of features are oftentimes recorded in clinical and non-clinical studies. For instance, a set of brain signals recorded from several individuals may be grouped by the type of stimulus each individual received in a given experimental setting. In such cases, one of the main goals of the time series analysis is to determine if there are differences in the spectral characteristics of the signals across the groups.
More specifically, suppose that a collection of time series \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\{ y_{i_{1},i_{2},h}(t)\}$\end{document}
are recorded during a particular experimental setting in which i 1 indexes the level of some factor, say A, i 2 indexes the level of some other factor, say B, and h indexes the time series within given levels of the A and B factors, where i 1=1:N 1, i 2=1:N 2, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$h=1:H_{i_{1},i_{2}}$\end{document}
, and t=1:T. Such data may arise in designed experiments where factor A is a given treatment or experimental condition—e.g., the type of stimulus in a neuroscience experiment or the disease level—, factor B is another treatment or experimental condition, and h indexes the individuals who were assigned to levels i 1 and i 2 of the A and B factors, respectively. We are interested in fitting models that allow us to characterize the spectral features of the time series in terms of effects associated to the levels of the factors and, possibly also, in terms of effects associated with the interaction of such factors. Assume for instance that each time series can be decomposed as
where the unobserved processes \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{1}}^{(1)}(t)$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{2}}^{(2)}(t)$\end{document}
are assumed to be uncorrelated, i.e., \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\mbox{Cov}(\beta_{i_{1}}^{(1)}(t),\beta_{i_{2}}^{(2)}(t))=0$\end{document}
for all i 1,i 2 and all t=1:T. Furthermore, for a given k with k=1,2 assume that \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\mbox{Cov}(\beta_{i}^{(k)}(t),\beta_{j}^{(k)}(t))=0$\end{document}
when i≠j for all t. The \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\varepsilon_{i_{1},i_{2},h}(t)$\end{document}
s are idiosyncratic—i.e., individual specific—error series assumed to be independent of the other unobserved processes and also mutually independent. This implies that the spectral density of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$y_{i_{1},i_{2},h}(t)$\end{document}
, which is the Fourier transform of the theoretical autocovariance function of the \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$y_{i_{1},i_{2},h}(t)$\end{document}
process and is denoted as \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{y_{i_{1},i_{2},h}}^{*}(\lambda)$\end{document}
, can be written as
where \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{i_{1}}^{*,(1)}(\lambda)$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{i_{2}}^{*,(2)}(\lambda)$\end{document}
are the spectral densities of the unobserved processes associated with levels i 1 and i 2 of factors A and B, respectively; and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f^{*}_{\varepsilon_{i_{1},i_{2},h}}(\lambda)$\end{document}
is the spectral density of the idiosyncratic error component associated to individual h within such levels.
There is a relatively rich literature on Bayesian approaches for estimation of functions that can be used in the context of analyzing multiple signals. DiMatteo, Genovese, and Kass (Reference DiMatteo, Genovese and Kass2001) describe a Bayesian method for fitting curves to data drawn from a distribution in the exponential family. This approach assumes that the curves can be well approximated by splines with an unknown number of knots and unknown knot locations that are then inferred via reversible jump Markov chain Monte Carlo. The method is used to analyze individual (not multiple) time series data obtained from functional magnetic resonance imaging (fMRI) experiments. In such context the method of DiMatteo et al. (Reference DiMatteo, Genovese and Kass2001) is a time-domain approach that assumes that a given time series y(t) can be modeled as y(t)=f(t)+ϵ(t), with \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f(t)=\sum_{j=1}^{k+2}\theta_{j} b_{j}(t)$\end{document}
, where b j(t) is the jth function of a cubic B-spline. In contrast, our proposed frequency-domain method aims to infer the spectral characteristics of multiple, not just a single, time series. A possible way of analyzing multiple times series in the spectral domain would consist on extending the method of DiMatteo et al. (Reference DiMatteo, Genovese and Kass2001) for fitting curves to the periodogram ordinates of the multiple time series recorded at various levels of the factors. Such extension would involve developing and implementing models that can handle multiple series given that the approach of DiMatteo et al. (Reference DiMatteo, Genovese and Kass2001) can only be directly applied to a single data set, e.g., the periodogram ordinates of an individual time series. Furthermore, the method of DiMatteo et al. (Reference DiMatteo, Genovese and Kass2001) applies to independent data (y 1,x 1),…,(y n,x n) that satisfy a model of the form y i|x 1,…,x n∼p(y i|f(x i),σ). However, the periodogram ordinates are only asymptotically independent. Therefore, it would not be advisable to directly apply the free-knots splines method in order to estimate spectral densities. On the other hand, the method proposed here is based on the approaches of Choudhuri, Ghosal, and Roy (Reference Choudhuri, Ghosal and Roy2004) and Macaro (Reference Macaro2010), which provide consistent Bayesian estimates even if the periodogram ordinates are only asymptotically independent (see Choudhuri et al. Reference Choudhuri, Ghosal and Roy2004, pp. 1056–1057, and Macaro Reference Macaro2010, p. 384).
Regarding factorial temporal data, several authors have considered various time-domain and frequency-domain approaches. Some key references include, among others, Shumway (Reference Shumway1970), Brillinger (Reference Brillinger and Krishnaiah1973, Reference Brillinger, Krishnaiah and Brillinger1980), and Stoffer (Reference Stoffer1999). Shumway and Stoffer (Reference Shumway and Stoffer2006) also summarize and illustrate several aspects of some of such approaches. As mentioned above, we follow Choudhuri et al. (Reference Choudhuri, Ghosal and Roy2004) and Macaro (Reference Macaro2010) to obtain Bayesian non-parametric posterior inference of the spectral representation of multiple time series, including the particular case of factorial spectral decompositions, by using Bernstein–Dirichlet prior distributions (Petrone Reference Petrone1999a, Reference Petrone1999b) on the spectral densities in (2). More specifically, Choudhuri et al. (Reference Choudhuri, Ghosal and Roy2004) described a Bayesian approach to estimating the spectral density of a single stationary time series by imposing a non-parametric prior on such density through Bernstein polynomials. Then, Macaro (Reference Macaro2010) proposed a mixture generalization of such approach in which each component in the spectral decomposition is identified using informative prior distributions. We extend the methods and algorithms of Macaro (Reference Macaro2010) to provide a Bayesian non-parametric spectral analysis of multiple time series recorded during designed experiments that may involve several factors, various levels within each factor and several individuals. Bernstein–Dirichlet prior distributions have been successfully applied to analyze other types of data that arise in psychometric applications. For example, Karabatsos and Walker (Reference Karabatsos and Walker2009) use a bivariate Bernstein–Dirichlet prior in the context of modeling test equating. Alternative Bayesian approaches to spectral density estimation of a single time series, including stationary and long-range dependence time series, as well as nonstationary time series with slowly varying dynamics or piecewise stationary time series, can be found in Carter and Kohn (Reference Carter and Kohn1997), Gangopadhyay, Mallick, and Denison (Reference Gangopadhyay, Mallick and Denison1998), Liseo, Marinucci, and Petrella (Reference Liseo, Marinucci and Petrella2001), and Rosen, Stoffer, and Wood (Reference Rosen, Stoffer and Wood2009), among others. A Bayesian approach for estimating the spectral density of multivariate stationary time series is presented in Rosen and Stoffer (Reference Rosen and Stoffer2007).
The article is organized as follows. The details of the modeling approach, including some examples, are provided in Section 2. Section 3 summarizes the Markov chain Monte Carlo algorithmFootnote 1 for posterior inference. In Section 4 we present the results of several simulation studies and discuss aspects of identification through the prior distributions. In Section 5 we analyze multiple time series of fMRI brain responses measured in individuals who participated in an experiment designed to study pain perception in humans. Finally, Section 6 presents the conclusions.
2. Bayesian Spectral Decompositions of Factorial Time Series Data
2.1. General Model Formulation
Let \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\{y_{i_{1},\ldots,i_{D},h}(t)\}$\end{document}
be a set of time series for i d=1:N d, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$h=1:H_{i_{1},\ldots,i_{D}}$\end{document}
and t=1:T. Assume that each time series \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$y_{i_{1},\ldots,i_{D},h}(t)$\end{document}
can be decomposed as a sum of D unobservable components plus an idiosyncratic error term. That is,
where \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{j}}^{(j)}(t)$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{k}}^{(k)}(t)$\end{document}
are assumed to be independent for all j≠k and all t, and the \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\varepsilon_{i_{1},\ldots,i_{D},h}(t)$\end{document}
are assumed to be independent of the \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{j}}^{(d)}(t)$\end{document}
processes for all d and also mutually independent. Furthermore, it is assumed that for a given d, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{j}}^{(d)}(t)$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{k}}^{(d)}(t)$\end{document}
are independent for all i j≠i k and all t, and that \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\varepsilon_{i_{1},\ldots,i_{D},h}(t_{0})$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\varepsilon_{i_{1},\ldots,i_{D},h}(t_{1})$\end{document}
are independent for all t 0≠t 1. This implies that the spectral density of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$y_{i_{1},\ldots,i_{D},h}(t)$\end{document}
, denoted as \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{y_{i_{1},\ldots,i_{D},h}}^{*}(\lambda)$\end{document}
, can be written as
where \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{i_{d}}^{*,(d)}(\lambda)$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{\varepsilon_{i_{1},\ldots,i_{D},h}}^{*}(\lambda)$\end{document}
are, respectively, the spectral densities of the unobserved factors and the idiosyncratic error components. The spectral density is the Fourier transform of the theoretical autocovariance function. Under some regularity conditions it is a real and continuous function over λ∈(−π,π]. A further assumption requires the spectral densities to be bounded and bounded away from zero, ruling out long memory and band limited processes.
2.2. Examples
2.2.1. Baseline Plus Single Factor Model
In this example we consider a one-way model represented as a baseline process plus one factor with two levels. Specifically, let d=1:2, i 1=1, i 2=1:2, h=1:2, and assume the following structure on \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$y_{1,i_{2},h}(t)$\end{document}
:
Therefore, the baseline process \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(1)}(t)$\end{document}
is an autoregressive process of order one, or AR(1), with coefficient 0.9, while the factor has two levels, the first level corresponds to an autoregressive process of order two, or AR(2), with two reciprocal characteristic roots each with modulus 0.949 and frequency 2.006 (period 3.132), and the second level is an AR(2) with two reciprocal characteristic roots each with modulus 0.894 and frequency 2.566 (period 2.449). The spectral representation is
with \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{1}^{*,(1)}(\lambda)=1/[2\pi(1+(\phi^{(1)})^{2}-2\phi^{(1)} \cos(\lambda))]$\end{document}
, and
Therefore, the decomposition in (8) is such that \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{1}^{*,(1)}(\lambda)$\end{document}
captures constant and persistent effects, while the \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{i_{2}}^{*,(2)}(\lambda)$\end{document}
spectra capture quasiperiodic dynamics that shift across the levels, as shown in Figure 1.
2.2.2. Two-Way Models
The probabilistic spectral decomposition in the example below corresponds to that of a two-way temporal model. In particular, assume that a two-way factor model with two levels in each factor is obtained by letting d=1:2, i 1=1:2, i 2=1:2, h=1:H and by assuming the following structure on \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$y_{i_{1},i_{2},h}(t)$\end{document}
:
Therefore, the spectral representation of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$y_{i_{1},i_{2},h}(t)$\end{document}
is given by
with \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{1}^{*,(1)}(\lambda)$\end{document}
the flat spectrum of a white-noise process and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{2}^{*,(1)}(\lambda)$\end{document}
, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{1}^{*,(2)}(\lambda)$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{2}^{*,(2)}(\lambda)$\end{document}
the spectra of quasiperiodic AR(2) processes. The densities in the spectral decomposition of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$y_{i_{1},i_{2},h}(t)$\end{document}
are shown in Figure 2.
2.2.3. Interaction Effects
Oftentimes studies are characterized by the presence of one or more interaction terms. For example, a two-way time-domain model with interactions can be written as
where \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\gamma_{i_{1},i_{2}}^{(1,2)}(t)$\end{document}
are time series processes that represent possible interaction effects between factors 1 and 2 at the i 1 and i 2 levels. In general, if we have D factors, we can extend the model in (3) to allow for second order interactions as follows:
The Bayesian non-parametric approach described in Section 3 is illustrated with one-way and two-way models with no interactions in Sections 4 and 5. However, the methodology can be applied to models that include interaction effects by adding such terms in Whittle’s likelihood approximation and by directly modeling the spectra of the interaction processes using the Bernstein–Dirichlet priors discussed in Section 2.4. Note that higher order interaction processes (e.g., third or higher order) can be added to (13) and modeled similarly.
As it will be shown in Section 4 with simulated data, when the number of factors increases, more restrictions need to be added to guarantee the identifiability of the different processes in the spectral decomposition. Therefore, instead of adding interaction effects and directly modeling the spectral densities of such processes, we propose simple post-processing calculations that can be used to explore if such interaction processes are present. We summarize such calculations in the case of a two-way model. Assume that \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$y_{i_{1},i_{2},h}(t)$\end{document}
is represented as the two-way model process in (12). Such model can be written as \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$y_{i_{1},i_{2},h}(t) = \beta_{i_{1}}^{(1)}(t) + \beta_{i_{2}}^{(2)}(t)+\varepsilon^{\prime}_{i_{1},i_{2},h}(t)$\end{document}
, with \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\varepsilon^{\prime}_{i_{1},i_{2},h}(t)= \gamma_{i_{1},i_{2}}(t)+ \varepsilon_{i_{1},i_{2},h}(t)$\end{document}
. If we were to analyze this scenario without implicitly considering the interaction terms, such terms would appear in the idiosyncratic factors \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\varepsilon^{\prime}_{i_{1},i_{2},h}(t)$\end{document}
. Therefore, we could obtain estimates of the spectral densities of the interaction terms as follows:
-
1. Evaluate the common baseline across all the individuals (i 1,i 2,h). Compute
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$f^{*}_{\bar{\varepsilon}}(\lambda)= \min_{i_1,i_2,h} f_{\varepsilon_{i_1,i_2,h}}^{*\prime}(\lambda).$$\end{document}
-
2. Evaluate the corresponding individual residual. For each λ and each (i 1,i 2,h) compute
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$f_{[\mathrm{residual}]_{i_1,i_2,h}}^{*\prime}(\lambda)=f_{\varepsilon_{i_1,i_2,h}}^{*\prime}(\lambda)-f^{*}_{\bar{\varepsilon}}(\lambda).$$\end{document}
-
3. Evaluate the interaction term across all the individuals within the same group. For each λ and each (i 1,i 2) compute an estimate of the spectrum of the interaction process \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$\gamma_{i_{1},i_{2}}^{(1,2)}(t)$\end{document}
as
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$f^{*}_{\gamma_{i_1,i_2}}(\lambda)=\min_{h}f_{[\mathrm{residual}],i_1,i_2,h}^{*\prime}(\lambda). $$\end{document}
By inspecting the posterior distribution of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{\gamma_{i_{1},i_{2}}}^{*}(\lambda)$\end{document}
for each combination of (i 1,i 2), it can be determined if one or more interaction terms should be included in the analysis.
2.3. Spectral Representation
Following a similar approach to that in Macaro (Reference Macaro2010), we use a discrete version of Whittle’s approximation of the likelihood function (Whittle Reference Whittle1957, Reference Whittle1962) to estimate the factors. That is, the likelihood is approximated as
where \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$I_{y_{i_{1},\ldots,i_{D},h}}(\lambda_{j})$\end{document}
are the periodogram ordinates evaluated at the Fourier frequencies,
These summarize the information obtained from the data, while the \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{\beta_{i_{d}}}^{*,d}(\lambda)$\end{document}
s are spectral ordinates to be estimated.
We now combine Whittle’s approximation to the likelihood in (14) with prior distributions on each \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{\beta_{i_{d}}}^{*,d}(\lambda)$\end{document}
to obtain posterior inference in the spectral domain. As discussed in Macaro (Reference Macaro2010), prior elicitation is difficult since it involves the approximation of continuous spectral densities over the set of Fourier frequencies. Choudhuri et al. (Reference Choudhuri, Ghosal and Roy2004) propose the use of Bernstein–Dirichlet prior distributions (Petrone Reference Petrone1999a, Reference Petrone1999b), while Macaro (Reference Macaro2010) presents a mixture generalization of their work. Specifically, the spectrum of an observable time series is modeled as a weighted sum of several unobserved spectral densities in Macaro (Reference Macaro2010). The resulting procedure is indeed more flexible, although the real novelty of the approach is the possibility of studying unobserved components from a non-parametric viewpoint. Here we extend Macaro (Reference Macaro2010) to consider factorial designs with multiple time series. In other words, the models and methods described in this section lead to a non-parametric spectral representation of the unobserved components underlying multiple time series recorded in factorial experimental settings. Below we describe the prior structure and how posterior inference can be achieved under such prior structure. We also discuss some issues related to prior choice and summarize the steps of a Markov chain Monte Carlo (MCMC) algorithm for posterior inference.
2.4. The Bernstein–Dirichlet Prior Distribution
The Bernstein polynomial prior developed by Petrone (Reference Petrone1999a, Reference Petrone1999b) was used as a non-parametric prior to estimate densities on the [0,1] interval, and was then used by Choudhuri et al. (Reference Choudhuri, Ghosal and Roy2004) for Bayesian estimation of spectral densities. This non-parametric prior hinges on the uniform convergence of
to q(λ), a continuous probability density on [0,1]. Here Q †(ϵ 1,ϵ 2):=Q(ϵ 2)−Q(ϵ 1), with Q(λ) being the probability function associated to q(λ). In particular, the Bernstein–Dirichlet prior induces the weights w 1,…,w K of the mixture through a Dirichlet process, i.e.,
with
where θ 1,…,θ K are the shape parameters, M>0 is the concentration parameter, and K determines how flexible is the prior, with larger values of K leading to more flexible distributions.
We now represent each of the spectral densities, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{\beta_{i_{d}}}^{*,(d)}(\lambda)$\end{document}
in (4), through a Bernstein kernel. In order to implement Bernstein–Dirichlet priors in the spectral domain, such densities must be normalized to the interval [0,1]. Therefore, we can define a pseudo-spectral density function, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$q_{y_{i_{1},\ldots,i_{D},h}}(\lambda)$\end{document}
, and a normalization parameter, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\tau_{y_{i_{1},\ldots,i_{D},h}}$\end{document}
, such that
where
Then, the normalization parameters for the spectral representation of the unobserved factors \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{d}}^{(d)}(t)$\end{document}
and the idiosyncratic errors \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\varepsilon_{i_{1},\ldots,i_{D},h}(t)$\end{document}
can be derived as follows:
such that \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\tau_{y_{i_{1},\ldots,i_{D},h}}=\sum_{d=1}^{D} \tau_{\beta_{i_{d}}^{(d)}} + \tau_{\varepsilon_{i_{1},\ldots,i_{D},h}}$\end{document}
, with
Since \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${\tau}_{\beta_{i_{1}}^{(d)}},\ldots,{\tau}_{\beta_{i_{D}}^{(d)}}$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${\tau}_{\varepsilon_{i_{1},\ldots,i_{D},h}}$\end{document}
are proportional to the variances of the factors, inverse-gamma distributions can be used as priors. That is,
The pseudo-spectral densities \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${q}_{\beta_{i_{d}}}^{(d)}(\lambda)$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${q}_{\varepsilon_{i_{1},\ldots,i_{D},h}}(\lambda)$\end{document}
will be discretized over the set of Fourier frequencies and estimated using Bernstein–Dirichlet priors.
Note that infinitely many choices of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{d}}^{(d)}$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\varepsilon_{i_{1},\ldots,i_{D}}(t)$\end{document}
would lead to the same \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$y_{i_{1},\ldots,i_{D},h}(t)$\end{document}
. The same is true in the spectral domain. Therefore, identifiability conditions need to be imposed. One way of doing so is by adding ANOVA-type restrictions (see Section 4.2). Alternatively, following Macaro (Reference Macaro2010), identifiability can be achieved through proper specification of informative Bernstein–Dirichlet prior as follows. For each \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$q_{\beta_{i_{d}}}^{(d)}(\lambda)$\end{document}
we have the representation in (15) and (16) characterized by weights \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$w_{k_{d}}^{i_{d},d}$\end{document}
for \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$k_{d}=1:K_{i_{d}}^{d}$\end{document}
and concentration parameter \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{i_{d}}^{d}$\end{document}
. Similarly, for each \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$q_{\varepsilon_{i_{1},\ldots,i_{D},h}}(\lambda)$\end{document}
we have a Bernstein–Dirichlet prior characterized by mixture weights \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$w^{\varepsilon_{i_{1},\ldots,i_{D},h}}_{k}$\end{document}
for \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$k=1:K_{i_{1},\ldots,i_{D},h}$\end{document}
and concentration parameter \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{i_{1},\ldots,i_{d},h}^{\varepsilon}$\end{document}
. The concentration parameters control the variance of the Dirichlet processes. In other words, as \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{i_{d}}^{d} \rightarrow 0$\end{document}
the variance of the corresponding Dirichlet process increases; and as \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{i_{d}}^{d} \rightarrow \infty $\end{document}
such variance goes to zero. Therefore, instead of letting \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{i_{d}}^{d} \rightarrow \infty$\end{document}
, we choose informed priors (\documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{i_{d}}^{d}< \infty$\end{document}
) on the \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$f_{i_{d}}^{(d)}(\lambda)$\end{document}
s. In practice, if we use the Dirichlet process characterization proposed by Sethuraman (Reference Sethuraman1994), a natural choice is to fix \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{i_{d}}^{d}=1$\end{document}
(Choudhuri et al. Reference Choudhuri, Ghosal and Roy2004), which is equivalent to saying that the prior for \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$v^{i_{d},d}_{k_{d},l}$\end{document}
for l=1:L is flat (see Section 3 for a detailed description of such parameters).
Theoretically, we could also set a prior on \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$K_{i_{d}}^{d}$\end{document}
as proposed by Choudhuri et al. (Reference Choudhuri, Ghosal and Roy2004), however, we fix these parameters in order to improve the speed and the stability of the MCMC algorithm (see Section 3 for further discussion). This will be illustrated in the analyses of simulated and real data sets in Sections 4 and 5. One of the main reasons for choosing this approach is that in many practical scenarios we have information about the unobserved components in the decomposition of the series that we may want to incorporate into the model through the prior distribution. For instance, electroencephalographic signals typically display activity in four main frequency bands (see, e.g., Prado, Reference Prado, O’Hagan and West2010), which induces a natural way of decomposing each signal into four unobserved components, one per frequency band.
Finally, it is important to notice the difference between \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$K_{i_{d}}^{d}$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{i_{d}}^{d}$\end{document}
. The parameter \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$K_{i_{d}}^{d}$\end{document}
is similar to the bandwidth parameter used when smoothing the periodogram: a smaller value of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$K_{i_{d}}^{d}$\end{document}
leads to a higher degree of smoothing. On the other hand, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{i_{d}}^{d}$\end{document}
indicates the degree of trust that we have in our prior distribution.
3. Posterior Inference
The posterior associated with the likelihood function (14) and with the priors described above is not available in closed form. Nevertheless, samples from this distribution can be obtained with a Metropolis–Hastings Markov chain Monte Carlo (MCMC) algorithm. The original construction of the Dirichlet process is not suitable for MCMC schemes. Sethuraman (Reference Sethuraman1994) proposed an alternative construction which has been implemented by Gelfand and Kottas (Reference Gelfand and Kottas2002) for a more flexible algorithm. Our sampling procedure is a generalization of the algorithm in Macaro (Reference Macaro2010) which, in turn, is an adaptation of the work of Choudhuri et al. (Reference Choudhuri, Ghosal and Roy2004) for the analysis of spectral densities. Particularly, the MCMC algorithm outlined in Section 3.2 is not built for the weight parameters \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$w^{i_{d},d}_{1},\ldots,w^{i_{d},d}_{K_{i_{d}}^{d}}$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$w^{\varepsilon_{i_{1},\ldots,i_{D},h}}_{1},\ldots,w^{\varepsilon_{i_{1},\ldots,i_{D},h}}_{K_{i_{1},\ldots,i_{D},h}}$\end{document}
, but for the set of parameters \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$r^{i_{d},d}_{k_{d},1},\ldots,r^{i_{d},d}_{k_{d},L}$\end{document}
, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$v^{i_{d},d}_{k_{d},1},\ldots,v^{i_{d},d}_{k_{d},L}$\end{document}
, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$r^{\varepsilon_{i_{1},\ldots,i_{D},h}}_{k_{h},1},\ldots,r^{\varepsilon_{i_{1},\ldots,i_{D},h}}_{k_{h},L}$\end{document}
, and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$v^{\varepsilon_{i_{1},\ldots,i_{D},h}}_{k_{h},1},\ldots,v^{\varepsilon_{i_{1},\ldots,i_{D},h}}_{k_{h},L}$\end{document}
such that
with
and
with
where the truncation level of the Dirichlet process is L=max(20,T 1/3) (Choudhuri et al. Reference Choudhuri, Ghosal and Roy2004).
3.1. Prior Choice
The prior probabilities can be chosen by performing the following steps:
-
1. Choose the number and types of factors and idiosyncratic components. This will depend on the data structure and the particular decomposition of interest (e.g., one-way vs. two-way factor structures). For example, in the data analysis presented in Section 5, for each time series we use a model with two factors, where each factor has two levels, and one idiosyncratic component which tries to capture the individual-specific features that are not described by any of the two factors.
-
2. Determine the a priori spectral properties that characterize each factor and transform this information into kernel parameters \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$\theta^{i_{d},d}_{1},\ldots,\theta^{i_{d},d}_{K_{i_{d}}^{d}}$\end{document}
for the corresponding Dirichlet process. For the idiosyncratic components the parameters can be set all equal to one to represent the prior assumption that these components are white noise. The prior for the remaining parameters \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$r^{i_{d},d}_{k_{d},1},\ldots,r^{i_{d},d}_{k_{d},L}$\end{document}, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$v^{i_{d},d}_{k_{d},1},\ldots, v^{i_{d},d}_{k_{d},L}$\end{document}, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$r_{k_{h},1}^{\varepsilon_{i_{1},\ldots,i_{D},h}},\ldots, r_{k_{h},L}^{\varepsilon_{i_{1},\ldots,i_{D},h}}$\end{document}, and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$v_{k_{h},1}^{\varepsilon_{i_{1},\ldots,i_{D},h}},\ldots,v_{k_{h},L}^{\varepsilon_{i_{1},\ldots,i_{D},h}}$\end{document} of the Bernstein representation are chosen according to the guidelines provided by Choudhuri et al. (Reference Choudhuri, Ghosal and Roy2004). -
3. Determine the flexibility of the priors on the factors: the smaller the value of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$K_{i_{d}}^{d}$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$K_{i_{1},\ldots,i_{D},h}$\end{document}, the smoother the corresponding spectral densities. Notice that this is particularly useful to separate the white-noise effects—which are characterized by flat and smooth spectral densities—from other effects. It is worth to remember that although \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$K_{i_{d}}$\end{document} and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$K_{i_{1},\ldots,i_{D},h}$\end{document} control for the smoothness of the spectra, they do not control for the strength of the prior distributions (this is achieved by changing the values of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{i_{d}}^{d}$\end{document} and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{i_{1},\ldots,i_{D},h}$\end{document}). -
4. Choose the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$M_{i_{d}}^{d}$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{i_{1},\ldots,i_{D},h}$\end{document} concentration parameters: large values reduce the variance of the corresponding Dirichlet processes. Note for example that by letting \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{i_{d}}^{d} \rightarrow\infty$\end{document} the underlying components become fixed (for more details see Macaro, Reference Macaro2010). -
5. Determine the prior distributions for the normalization parameters. These parameters are proportional to the variances, therefore inverse gammas seem to be appropriate:
(22)\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \tau_{(.)}\sim \mathrm {InvGamma} \bigl\{ s^\mathrm{prior}_{(.)}, \tilde{\tau}_{(.)}\bigl[s^\mathrm{prior}_{(.)}-1\bigr] \bigr\} $$\end{document}with \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$s^{\mathrm{prior}}_{(.)}$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\tilde{\tau}_{(.)}$\end{document} chosen using the fact that
(23)\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \mathrm {E}[\tau_{(.)} ]=\tilde{\tau}_{(.)}\quad \mbox{and} \quad \mathop {\mathrm {var}}[\tau_{(.)} ]=\frac{\tilde{\tau}_{(.)}^2}{s_{(.)}^\mathrm{prior}-2}. $$\end{document}
3.2. MCMC Algorithm
The posterior distribution associated with the likelihood function (14) and with the prior distributions described above is not available in closed form. The following Metropolis–Hastings MCMC algorithm can be used to obtain samples from the posterior distribution.
-
1. For each factor indexed by d=1:D
-
(a) For each level of by i d=1:N d of the factor d
-
i. Update v. For each l=1:L
-
• Propose a candidate value for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$v_{k_{d},l}^{i_{d},d}$\end{document}
. -
• Evaluate the new set of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$p_{k_{d},1}^{i_{d},d},\ldots,p_{k_{d},L}^{i_{d},d}$\end{document}
and the corresponding new set of mixture weights \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$w_{1}^{i_{d},d},\ldots,w^{i_{d},d}_{K_{i_{d}}^{d}}$\end{document}. -
• Evaluate the new spectral density \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$f^{*,(d)}_{i_{d}}(\lambda)$\end{document}
of the corresponding factor component and the likelihood (14). -
• Decide whether to accept the proposed value of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$v_{k_{d},l}^{i_{d},d}$\end{document}
.
-
-
ii. Update r. For each l=0:L
-
• Propose a candidate value for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$r_{k_{d},l}^{i_{d},d}$\end{document}
. -
• Evaluate the new set of mixture weights \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$w_{1}^{i_{d},d},\ldots,w^{i_{d},d}_{K_{i_{d}}^{d}}$\end{document}
. -
• Evaluate the new spectral density \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$f^{*,(d)}_{i_{d}}(\lambda)$\end{document}
of the corresponding factor component and the likelihood (14). -
• Decide whether to accept the proposed value of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$r_{k_{d},l}^{i_{d},d}$\end{document}
.
-
-
iii. Updateτ. Propose a candidate value of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$\tau_{\beta_{i_{d}}}^{(d)}$\end{document}
and evaluate the corresponding spectral density and the likelihood (14).
-
-
-
2. Given the updates for the levels i d of each process d, update the corresponding idiosyncratic component.
-
(a) For the set of indices (i 1,…,i D,h),
-
i. Update v for the spectral density \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$f_{\varepsilon_{i_{1},\ldots,i_{D},h}}^{*}(\lambda)$\end{document}
as described in 1(a)i. -
ii. Update r for the spectral density \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$f_{\varepsilon_{i_{1},\ldots,i_{D},h}}^{*}(\lambda)$\end{document}
as described in 1(a)ii. -
iii. Updateτ for the spectral density \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$f_{\varepsilon_{i_{1},\ldots,i_{D},h}}^{*}(\lambda)$\end{document}
as described in 1(a)iii.
-
-
The above algorithm was implemented in R (R Development Core Team 2007) using the MCMC package written by Martin and Quinn (2005).Footnote 2 The algorithm is run for I iterations with I as large as needed to obtain MCMC convergence. This algorithm is an extension of the MCMC scheme of Macaro (Reference Macaro2010) that allows us to consider factorial time series data. From the computational point of view this is rather complex; and, so, the computation time required to do the analysis increases with the number of factors, the number of levels within each factor, and the number of time series within each combination of factors. The proposed candidate values in Steps 1(a)i–ii and 2(a)i–ii are sampled from uniform distributions and those in Steps 1(a)iii and 2(a)iii are sampled from inverse-Gamma distributions whose moments are functions of the old parameter values (see Macaro Reference Macaro2010 for further details). In addition, the proposed values are accepted or rejected according to Metropolis–Hastings schemes.
4. Simulation Studies
4.1. Baseline Plus Single Factor Model
A data set with a total of four time series of length T=500 was simulated from the one-way model described in Section 2.2.1. The priors are chosen to be Beta kernels. Specifically, we chose a low frequency band-pass kernel for the baseline, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${q}_{\beta_{1}^{(1)}}(\lambda) \sim \mathrm {Beta} (1,20)$\end{document}
, and an approximately periodic kernel for the two levels of the remaining factor \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${q}_{\beta^{(2)}_{1}}(\lambda) \sim \mathrm {Beta} (50,20)$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${q}_{\beta^{(2)}_{2}}(\lambda) \sim \mathrm {Beta} (50,20)$\end{document}
. Finally, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${q}_{\varepsilon_{i_{1},i_{2},h}}(\lambda) \sim \mathrm {Beta} (1,1)$\end{document}
implies a white-noise assumption for the error term. The parameters controlling the flexibility of the representation are set to \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$K_{i_{1}}^{1}=40$\end{document}
, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$K_{i_{2}}^{2}=40$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$K_{i_{1},i_{2},h}^{3}=1$\end{document}
. This prior choice empowers the belief that the main factors must not capture white-noise effects; and, so, their spectra are not expected to be flat a priori. In addition, the parameters which control for the variances of the Dirichlet distributions are set to \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{i_{d}}^{d}=1$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{\epsilon_{\cdot}}=1$\end{document}
.
Figure 3 shows the prior and posterior distributions (gray and black areas, respectively), as well as estimators based on the smoothed periodograms for the spectral densities of the observed time series and the unobserved processes \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(1)}(t), \beta_{1}^{(2)}(t)$\end{document}
, and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(2)}(t)$\end{document}
. Note that, in spite of the fact that the priors on the spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(2)}(t)$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(2)}(t)$\end{document}
had modes away from the correct periods, the estimated posterior spectral densities adequately capture the frequencies around 2.0 and 2.6.
Spectral analysis of time series data simulated from the one-way model described in Section 2.2.1. Gray areas are central 95 % spectral prior distributions, while black areas are central 95 % spectral posterior distributions. The dotted lines represent estimators of the spectral densities obtained by smoothing the periodograms of the data. First row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(1)}(t)$\end{document}
, y 1,1,1(t), and y 1,1,2(t). Second row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(2)}(t)$\end{document}
, ε 1,1,1(t), and ε 1,1,2(t). Third row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(1)}(t)$\end{document}
, y 1,2,1(t), and y 1,2,2(t). Fourth row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(2)}(t)$\end{document}
, ε 1,2,1(t) and ε 1,2,2(t).
4.2. Two-Way Models
We analyze time series data simulated from the two-way model described in Section 2.2.2 with d=2, i 1=1:2, i 2=1:2, h=1:2, and T=500. We consider two alternative ways of selecting the prior distributions in this example. First we assume that the spectral density of the first level of the first factor, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(1)}(t)$\end{document}
, as well as the spectral density of the ε ⋅,⋅,⋅(t) processes are constant over all the frequencies. This is consistent with the assumption that these processes are white-noise processes. In addition, we assume weak and flat priors on the spectral densities of the remaining time series processes. Specifically, we have \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${q}_{\beta_{2}^{(1)}}(\lambda) \sim \mathrm {Beta} (1.3,1.3)$\end{document}
, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${q}_{\beta^{(2)}_{1}}(\lambda) \sim \mathrm {Beta} (1.3,1.3)$\end{document}
, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${q}_{\beta^{(2)}_{2}}(\lambda) \sim \mathrm {Beta} (1.3,1.3)$\end{document}
. Notice that we purposely chose these kernels to be nearly flat to contrast the boundary bias (see Zhang & Karunamuni, 2010) which seems to affect weak and flat priors. The parameters controlling the flexibility of the representation are set again to \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$K_{i_{1}}^{1}=40$\end{document}
, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$K_{i_{2}}^{2}=40$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$K_{i_{1},i_{2},h}^{3}=1$\end{document}
. The parameters that control the variances of the Dirichlet distributions are set to \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{i_{d}}^{d}=1$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{\epsilon_{\cdot}}=1$\end{document}
.
By fixing the spectral densities of one of the levels of one of the factors and the spectral densities of the error term processes, identification of the remaining components is achieved, as shown in Figures 4 and 5.
Spectral analysis of time series data simulated from the two-way model described in Section 2.2.2. Gray areas are central 95 % spectral prior distributions, while black areas are central 95 % spectral posterior distributions. The dotted lines represent estimators of the spectral densities obtained by smoothing the periodograms of the data. First row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(1)}(t)$\end{document}
, y 1,1,1(t), and y 1,1,2(t). Second row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(2)}(t)$\end{document}
, ε 1,1,1(t), and ε 1,1,2(t). Third row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(1)}(t)$\end{document}
, y 1,2,1(t), and y 1,2,2(t). Fourth row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(2)}(t)$\end{document}
, ε 1,2,1(t) and ε 1,2,2(t).
Spectral analysis of time series data simulated from the two-way model described in Section 2.2.2. Gray areas represent central 95 % spectral prior distributions, while black areas are central 95 % spectral posterior distributions. The dotted lines represent estimators of the spectral densities obtained by smoothing the periodograms of the data. First row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(1)}(t)$\end{document}
, y 2,1,1(t), and y 2,1,2(t). Second row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(2)}(t)$\end{document}
, ε 2,1,1(t), and ε 2,1,2(t). Third row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(1)}(t)$\end{document}
, y 2,2,1(t), and y 2,2,2(t). Fourth row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(2)}(t)$\end{document}
, ε 2,2,1(t) and ε 2,2,2(t).
The second scenario does not fix the spectra of any of the factor processes. It is assumed that the error processes have white-noise spectra a priori, and slightly more informative priors on the factor processes are considered. In particular we assume band-pass filter types of priors that emphasize smaller frequencies for \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(1)}(t)$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(1)}(t)$\end{document}
and larger frequencies for \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(2)}(t)$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(2)}(t)$\end{document}
. The graphs in Figures 6 and 7 show that the implementation of a model with these band-pass filter priors produces results similar to those obtained under the white-noise assumption on \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(1)}(t)$\end{document}
considered above and adequately captures the data structure.
Spectral analysis of time series data simulated from the two-way model described in Section 2.2.2. Gray areas represent central 95 % spectral prior distributions, while black areas are central 95 % spectral posterior distributions. The dotted lines represent estimators of the spectral densities obtained by smoothing the periodograms of the data. First row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(1)}(t)$\end{document}
, y 1,1,1(t), and y 1,1,2(t). Second row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(2)}(t)$\end{document}
, ε 1,1,1(t), and ε 1,1,2(t). Third row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(1)}(t)$\end{document}
, y 1,2,1(t), and y 1,2,2(t). Fourth row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(2)}(t)$\end{document}
, ε 1,2,1(t) and ε 1,2,2(t).
Spectral analysis of time series data simulated from the two-way model described in Section 2.2.2. Gray areas represent central 95 % spectral prior distributions, while black areas are central 95 % spectral posterior distributions. The dotted lines represent estimators of the spectral densities obtained by smoothing the periodograms of the data. First row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(1)}(t)$\end{document}
, y 2,1,1(t), and y 2,1,2(t). Second row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(2)}(t)$\end{document}
, ε 2,1,1(t), and ε 2,1,2(t). Third row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(1)}(t)$\end{document}
, y 2,2,1(t), and y 2,2,2(t). Fourth row: estimated spectral densities of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(2)}(t)$\end{document}
, ε 2,2,1(t) and ε 2,2,2(t).
4.3. Further Identification Issues in Two-Way Models
In order to illustrate some additional features related to prior elicitation and corresponding posterior inference, we consider a larger simulated data set. Specifically, a total of 12 time series were simulated from a two-way model with the following structure:
for t=1:500, i 1=1:2, i 2=1:2, and l=1:3, where
with \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\varepsilon_{i_{1},i_{2},h}(t) \sim N(0,0.5)$\end{document}
for all i 1,i 2 and l, ν 1(t)∼N(0,0.5), and ν 2(t)∼N(0,0.5). In the notation above, AR[2,(0.95,9),0.5] denotes a quasiperiodic autoregressive process of order two, or AR(2), with a pair of complex reciprocal roots with modulus 0.95, wavelength 9—or equivalently, frequency 0.698—and standard deviation 0.5. Similarly, AR[2,(0.95,3),0.5] denotes a quasiperiodic AR(2) process with a pair of complex reciprocal roots with modulus 0.95, wavelength 3—or frequency 2.094—and standard deviation 0.5; AR[1,0.86,0.5] denotes an AR(1) process with a reciprocal root of modulus 0.86 and standard deviation 0.5; and finally, AR[1,0.94,0.5] denotes an AR(1) with a reciprocal root of modulus 0.94 and standard deviation 0.5.
The structure above implies that the spectra of the series y 1,1,h(t) and y 1,2,h(t) will have a peak at the frequency ω 1=0.698, in addition to the peak at zero related to the AR(1) components; while those of the series y 2,1,h(t) and y 2,2,h(t) will have an additional peak at the frequency ω 2=2.094. We proceed to analyze these data with the Bayesian spectral non-parametric approach described in the previous sections. In particular, we assume the following underlying structure on the observed data:
Furthermore, we assume that, a priori, the spectra of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{1}}^{(1)}(t)$\end{document}
for i 1=1:2 have a single peak around the frequency ω 1=0.698, the spectra of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{2}}^{(2)}(t)$\end{document}
for i 2=1:2 have only a peak at zero, and the spectra of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\varepsilon_{i_{1},i_{2},h}(t)$\end{document}
are flat.
Figures 8 and 9 show the prior structure (gray areas). Note that this prior is closely describing the structure underlying the processes y 1,1,h(t) and y 1,2,h(t), but it misspecifies the structure of the processes underlying y 2,1,h(t) and y 2,2,h(t), completely missing the second peak at ω 2=2.094.
Spectral analysis of the time series data simulated from the model described in Section 4.3. Gray areas represent central 95 % of the prior distributions. Black areas represent central 95 % of the posterior distributions from the two-way model. The dotted lines represent periodogram-based estimators of the spectral densities of the individual series without considering the factor structure.
Spectral analysis of the time series data simulated from the model described in Section 4.3. Gray areas represent central 95 % of the prior distributions. Black areas represent central 95 % of the posterior distributions from the two-way model. The dotted lines represent periodogram-based estimators of the spectral densities of the individual series without considering the factor structure.
Note that the error components, assumed to have flat spectra a priori under a white-noise structure, capture, a posteriori, whatever is left unexplained by the factor structure. Therefore, if an informative prior used on the spectra of a given factor misses certain frequency components present in the data, such components will be captured by the posterior distributions of the spectra of the idiosyncratic components, if relatively noninformative priors were given to these components. This is confirmed by the results shown in Figures 8 and 9. The prior distribution for \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(1)}(t)$\end{document}
does not include one of the frequency components in the true \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\alpha_{2}^{(1)}(t)$\end{document}
and the posterior distribution of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(1)}(t)$\end{document}
does not capture the peak at ω 2=2.094. Note, however, that this peak clearly appears in the estimated spectra of the six idiosyncratic error terms ε 2,1,h(t) and ε 2,2,h(t) for h=1:3, indicating that such structure is missing from the factorial decomposition of the process.
The analysis above that assumes informative priors on some of the unobserved components is performed in order to illustrate a couple of features about the spectral decomposition. First, the use of informative priors allow us to choose a particular decomposition of the series where some of the components may have a specific scientifically interpretable meaning a priori. For example, in this simulation study, we use a strong informative prior on the spectrum of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(1)}(t)$\end{document}
that identifies this time series process as one with activity in a particular frequency band. Alternatively, we could have used a more standard time-domain two-way representation of the data by assuming that the spectrum of one of the levels of one of the two factors, say the spectrum of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(1)}(t)$\end{document}
, is flat (consistent with a white-noise prior assumption on this process). Such representation would also imply strong assumptions on the components that would not lead to the true unobserved representation in (24), (25), and (26), but that would allow us to adequately capture the structure of the spectra of the observed time series a posteriori. Researchers can therefore explore several spectral decompositions of the data by making different assumptions on the unobserved components, as long as such assumptions guarantee identifiability. The second model feature that we want to emphasize is that, regardless of the identifiability conditions used on the spectral characterization of the factors \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{d}}^{(d)}(t)$\end{document}
, the model will be able to adequately describe the structure of the data a posteriori if the priors on the spectra of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\varepsilon_{i_{1},i_{2},h}(t)$\end{document}
are assumed to be flat and diffuse. As seen in the simulation study above, none of the priors on the spectra of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{1}}^{(d)}(t)$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\varepsilon_{i_{1},i_{2},h}(t)$\end{document}
had peaks on frequencies above 1.5, however, the estimated posterior spectra of ε 2,1,h(t), ε 2,2,h(t), and, therefore, the estimated posterior spectra of y 2,1,h(t) and y 2,2,h(t), show a peak at the ω 2=2.094 frequency that does not appear in any of the spectra of the \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$y_{i_{1},i_{2},h}(t)$\end{document}
series.
5. Functional Magnetic Resonance Imaging Data
We consider the analysis of a data set from an experiment of Antognini, Buonocore, Disbrow, and Carstens (Reference Antognini, Buonocore, Disbrow and Carstens1998) in which functional magnetic resonance imaging (fMRI) was used to examine pain perception in humans. The data, previously analyzed in Stoffer (Reference Stoffer1999) and Shumway and Stoffer (Reference Shumway and Stoffer2006), consist of consecutive measures of blood oxygenation level dependent (BOLD) signal intensities at nine locations in brain. The signals were measured in 26 awake and mildly anesthetized individuals who were presented with three types of periodic stimuli: brushing, heat, and shock. The stimuli were applied alternately for 32 seconds, with a sampling rate of one point every two seconds, and then stopped for 32 seconds.
Data were recorded by several sensors located on different parts of the brain, namely, Cortex 1, 2, 3 and 4, Caudate, Thalamus 1 and 2, and Cerebellum 1 and 2. Here we analyze only those related to the region Cortex 1. The 26 patients were divided into six groups depending on the type of stimulus they received (Brush, Heat, or Shock) and on whether they were awake or mildly anesthetized. Specifically, we have the following: (1) data from three patients in Low sedation state who were stimulated with a Brush; (2) data from five patients in Low sedation state stimulated by Heat; (3) data from four patients in Low sedation state stimulated by a Shock; (4) data from five Awake patients recorded while they were stimulated with a Brush; (5) data from four Awake patients stimulated by Heat; and finally, (6) data from five Awake patients stimulated by a Shock. Given this structure, a two-way model with one factor associated to the level of consciousness (Low and Awake) and another factor associated to the stimulus type (Brush, Heat, and Shock) was chosen to describe these data. Therefore, we consider a two-way model of the form
where \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{1}}^{(1)}(t)$\end{document}
for i 1=1:2 model the effects related to the consciousness level, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{2}}^{(2)}(t)$\end{document}
for i 2=1:3 model the effects of the type of stimulus, and each \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\varepsilon_{i_{1},i_{2},h}(t)$\end{document}
is the idiosyncratic component modeling those effects specific to each patient that are not captured by the factor structure. Since the number of patients varies within each combination of consciousness level and stimulus type, we have y 1,1,h(t) (Low and Brush) for h=1:3; y 1,2,h(t) (Low and Heat) for h=1:5; y 1,3,h(t) (Low and Shock) for h=1:4; y 2,1,h(t) (Awake and Brush) for h=1:5; y 2,2,h(t), (Awake and Heat) for h=1:4; and y 2,3,h(t) (Awake and Shock) for h=1:5.
The priors are chosen to be Beta kernels. Specifically, we chose \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${q}_{\beta_{1}^{(1)}}(\lambda) \sim \mathrm {Beta} (1.7,8.6)$\end{document}
; \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${q}_{\beta_{2}^{(1)}}(\lambda)\sim \mathrm {Beta} (1.7,10.6)$\end{document}
; \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${q}_{\beta^{(2)}_{1}}(\lambda) \sim \mathrm {Beta} (1.7,1.1)$\end{document}
; \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${q}_{\beta^{(2)}_{2}}(\lambda) \sim \mathrm {Beta} (1.7,1.1)$\end{document}
; \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${q}_{\beta^{(2)}_{3}}(\lambda) \sim \mathrm {Beta} (1.7,1.1)$\end{document}
; and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}${q}_{\varepsilon_{i_{1},i_{2},h}}(\lambda) \sim \mathrm {Beta} (1,1)$\end{document}
for all i 1,i 2 and h. The idea is that, a priori, the \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{1}}^{(1)}(t)$\end{document}
processes aim to describe low frequency components, the \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{i_{2}}^{(2)}(t)$\end{document}
aim to describe high frequency components and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\varepsilon_{i_{1},i_{2},h}(t)$\end{document}
are white-noise processes. For the remaining parameters, we set \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$K_{i_{1}}^{1}=30$\end{document}
(cyclical), \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$K_{i_{2}}^{2}=100$\end{document}
(long term) and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$K_{i_{1},i_{2},h}^{3}=10$\end{document}
(short-term/white noise) for i 1=1:2, i 2=1:3 and h=[(1:3),(1:5),(1:4),(1:5),(1:4),(1:5)]. The parameters that control the variances of the Dirichlet distributions were set to \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{i_{d}}^{d}=1$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$M_{\epsilon_{\cdot}}=1$\end{document}
.
Figure 10 shows the prior and posterior distributions (light and dark gray areas, respectively) of the spectra of the factors. Figure 11 shows the prior and posterior distributions (light and dark gray areas, respectively) of the spectra of the idiosyncratic processes for each of the 26 individuals.
Spectral analysis of the fMRI data. The plots display the spectra of the factor processes. The light gray areas represent central 95 % of the prior distributions and the dark gray areas represent central 95 % of the posterior distributions.
Spectral analysis of the fMRI data. The plots show the spectra of the idiosyncratic components for each of the 26 subjects. The light gray areas represent central 95 % of the prior distributions. The dark gray areas represent central 95 % of the posterior distributions.
From the graphs we can see that the posterior distributions of the spectral density for the process \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(1)}(t)$\end{document}
related to Low sedation level and that for the process \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(2)}(t)$\end{document}
related to Awake consciousness level show a marked peak around 0.19–0.22 which, given the sampling rate of the data, corresponds to the 64-second period band (or 1/64 cycles per second). The posterior spectral density of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(1)}(t)$\end{document}
(Low) shows more uncertainty in this frequency band than the posterior spectral density of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(1)}(t)$\end{document}
(Awake). In addition, the spectral density of \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(1)}(t)$\end{document}
(Low) shows additional activity in higher frequency bands. Similarly, comparing the posterior estimates of the spectral densities of the processes \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{1}^{(2)}(t)$\end{document}
, \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{2}^{(2)}(t)$\end{document}
and \documentclass[12pt]{minimal}
\usepackage{amsmath}
\usepackage{wasysym}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{amsbsy}
\usepackage{mathrsfs}
\usepackage{upgreek}
\setlength{\oddsidemargin}{-69pt}
\begin{document}$\beta_{3}^{(2)}(t)$\end{document}
associated, respectively, with the Brush, Heat and Shock stimuli, we see that the density for Brush shows a marked peaked around the 0.19–0.22 (again, the 64-second period band) that increases the power at this frequency in the estimated spectra of the subjects who received this type stimulus. For the other two types of stimuli, rather small peaks appear at higher frequencies. In summary, the posterior results based on the Bayesian non-parametric analysis for the data collected at the Cortex 1 location indicate that there is a significant difference in the power spectra at the 1/64 frequency band for patients that received the Brush stimulus.
6. Conclusions
This work presents a Bayesian non-parametric framework for the analysis of multiple time series that are collected in factorial experimental designs. This approach is based on representing the prior distributions on the spectral densities with Bernstein–Dirichlet priors. It is an important extension of the univariate methods of Macaro (Reference Macaro2010) in that it allows us to combine data and prior distributions on the unobserved spectral components in the decomposition of multiple time series. In particular, users can explore various decompositions of the observed time series processes and, consequently, various spectral decompositions, which can provide useful representations of the data.
Because of the possibility of extracting unobserved factors from a set of time series, this work is related to the well known literature on dynamic factor models (Geweke Reference Geweke, Aigner and Goldberger1977; Sargent & Sims, Reference Sargent, Sims and Sims1977; Forni, Hallin, Lippi, & Reichlin, Reference Forni, Hallin, Lippi and Reichlin2000; Stock & Watson, Reference Stock and Watson2002, and many others). In our work we emphasize how the priors can be used to guarantee the identifiability of the different model components at the spectral level. Forni et al. (Reference Forni, Hallin, Lippi and Reichlin2000) and Stock and Watson (Reference Stock and Watson2002), on the other hand, derived their identification conditions from the properties of the principal component analysis.
We illustrate the main features and the performance of the Bayesian non-parametric approach in the analyses of simulated data and in functional magnetic resonance (fMRI) brain responses measured in awake and mildly anesthetized individuals who were presented with three types of stimuli. The posterior distributions on the specific spectral decomposition of the data that we chose—i.e., a decomposition in which each observed signal was represented as a sum of two unobserved components, one related to the level of consciousness and one to the stimulus, plus an idiosyncratic term—reveals differences in the power spectra at certain frequency bands for patients who received a particular stimulus. The spectral domain approaches presented here can be used to study data from a broad range of applications in which multiple time series are collected in the context of a designed experiment.
Acknowledgements
This material was based upon work partially supported by the National Science Foundation under grant DMS-0635449 to the Statistical and Applied Mathematical Sciences Institute. Raquel Prado also was partially supported by NSF grant SES-1060911. Any opinions, findings, and conclusions or recommendations expressed in this material are those of the author(s) and do not necessarily reflect the views of the National Science Foundation.
Open access
