Linear System Challenges of Dynamic Factor Models

Abstract

A survey is provided dealing with the formulation of modelling problems for dynamic factor models, and the various algorithm possibilities for solving these modelling problems. Emphasis is placed on understanding requirements for the handling of errors, noting the relevance of the proposed application of the model, be it for example prediction or business cycle determination. Mixed frequency problems are also considered, in which certain entries of an underlying vector process are only available for measurement at a submultiple frequency of the original process. Certain classes of processes are shown to be generically identifiable, and others not to have this property.

1. Introduction

The purpose of this paper1 is to survey (and provide commentary on) a collection of contributions to do with Dynamic Factor Models developed especially by the authors, principally the first two, over more than a decade, taking the ideas as far as a listing of some problems of current interest. Matters which we are especially concerned to highlight include

• The theoretical underpinning for allowing modelling to focus on AR rather than ARMA modelling, and the simplifications this makes possible for the identification of models.
• Issues arising in practice from numerical problems in attempting to build models; some of these errors can be mitigated through the use of carefully chosen forms for the models.
• The difficulties arising from modelling time series of different but multiply related periodicities, e.g., monthly and quarterly, and existing tools for resolving, at least partly, the difficulties. These build on the ideas developed for single-periodicity modelling.

Feeding in all these matters is the application of ideas of digital signal processing, for which issues of handling numerical problems, and handling data streams of multiple periodicities have been well studied. The relevance of ideas of digital signal processing to dynamic factor modelling has perhaps not been fully recognized.

Dynamic Factor Models were introduced two decades ago, see in particular Forni et al. (2000), Forni and Lippi (2001), Stock and Watson (2002), Bai and Ng (2002), with macroeconometric applications as their main purpose. The basic idea was that a small number of factors, common to all the time-series of a large macroeconomic dataset, could be used both in structural modeling and forecasting of some series of interest like inflation, industrial production, unemployment. For a typical dataset, see the monthly US macroeconomic panel, including 128 series, described in great detail in McCracken and Ng (2016).

We begin by recalling and commenting on a number of ideas in a recent survey we coauthored, Lippi et al. (2022). The survey confined attention to the modelling of multivariate time series with the same sampling frequency for every component. As foreshadowed above, in this paper, besides dealing with such material, we treat modelling where time series data can be available with different periodicities, provided they are integrally related, e.g., monthly and quarterly.

Throughout the paper, we shall confine attention to stationary processes.

We take some time to define what is meant by a dynamic factor model, progressively introducing and commenting on the various assumptions that in aggregate constitute the definition. The comments may be as important as the assumptions. The various assumptions reflect reality to some degree, and they also reflect the requirement to formulate a solvable problem. We also offer commentary on the notion of errors; the comments deal with two issues, first the need to in some way limit the errors involved in determining a process model, given that the model is obtained from sample statistics rather than population data, and second, the need to choose a type of process model that is appropriate to the application: one does not want massive errors in whatever it is desired to compute because of sensitivities intrinsic to the model. An overall conclusion (motivated by the signal processing literature) is that on occasions state-variable models may be preferred to ARMA or AR. Accordingly, some discussion is also offered explaining how state-variable models can be determined without even intermediate use of an ARMA or AR model2. Analogously to the result that AR models can often be determined from covariance data using the Yule-Walker equations is an algorithm involving a finite number of rational calculations for obtaining state-variable models of a canonical spectral factor from lagged covariance data. There may well be numerical advantages in using such a state-variable model.

The last part of the paper traces developments applicable to multiple frequency systems, or systems in which data of different periodicities is collected. More precisely, we postulate that there is an underlying stationary zero mean vector time series $\left(y_t\right)={\left[{\left(y_t^f\right)}^{\top }{\left(y_t^s\right)}^{\top }\right]}^{\top }$ which is defined at time $t=0,1,2,\dots$ and with $y_t^f$ (the subvector of “fast” components) observed at all t while $y_t^s$ (the subvector of “slow” components) is observed at $t=0,J,2J,\dots$ for some integer $J>1$. The task is to build a canonical spectral factor of the process $\left(y_t\right)$ using the measured second-order statistics. The whole set-up can be embedded within a high-dimensional time-series framework of course, but this brief statement encapsulates what is different about the mixed frequency problem.

The conclusions section of the paper sets out a number of issues which have not yet been fully addressed.

2. Dynamic Factor Models

Dynamic factor models can roughly be described as those associated with a high dimension, say N, vector process $\left(y_t^N\right)$ (normally taken with zero mean) defined at $t=0,1,2,\dots$ and possessing an additive decomposition

$$y_t^N={\chi }_t^N+{\xi }_t^N$$

(1)

with $\left({\chi }_t^N\right)$ termed the process of common components or latent variable process, and $\left({\xi }_t^N\right)$ termed the process of idiosyncratic components. Typically, $\left({\chi }_t^N\right)$ is regarded as the true process of interest, with $\left({\xi }_t^N\right)$ some kind of contamination intruding into the measurement process (‘noise’ in electrical engineering parlance, which in general is colored rather than white). However one possible alternative interpretation is that (1) displays the decomposition of stock price observations into a part $\left({\chi }_t^N\right)$ representing the comovements between the variables (for instance representing the “market” effect on N different stock returns) and individual movements $\left({\xi }_t^N\right)$ (representing effects associated with each specific firm). It is assumed that the one-dimensional processes $({\chi }_{it},t\in \mathbb{Z}),i\in \{1,\dots ,N\}$ are also stationary and strongly dependent across the index i, while the processes $({\xi }_{it},t\in \mathbb{Z}),i\in \{1,\dots ,N\}$ are stationary and weakly cross-dependent (the terms strongly dependent and weakly cross-dependent requiring technical definitions). It is standard to allow for variation in N above some $N_0$, with $y_t^N,{\chi }_t^N$ and ${\xi }_t^N$ being nested, for example $y_t^N$ is the subvector of $y_t^{N+1}$ formed from its first N entries. Both $\left({\chi }_t^N\right)$ and $\left({\xi }_t^N\right)$ are zero mean and the two vector processes are mutually independent. The spectral densities ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right),\omega \in [-\pi ,\pi ]$ and ${\mathrm{\Phi }}_{\xi }^N\left(e^{j\omega }\right),\omega \in [-\pi ,\pi ]$ are assumed to exist and $\left(y_t^N\right)$ is then necessarily stationary with

$${\mathrm{\Phi }}_y^N\left(e^{j\omega }\right)={\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)+{\mathrm{\Phi }}_{\xi }^N\left(e^{j\omega }\right).$$

(2)

These spectral matrices obviously inherit the nesting property of the underlying processes.

We sum up the above description as follows.

Assumption 1.

A nested (with N) set of observation processes $\left(y_t^N\right)$ of zero mean arises as an additive combination of a stationary zero mean nested common component process $\left({\chi }_t^N\right)$ and a stationary zero mean nested process of idiosyncratic components $\left({\xi }_t^N\right)$, as in (1). The latter two processes are independent, the first is strongly cross-dependent (a term defined subsequently) and the latter is weakly cross-dependent (also defined subsequently), and the spectra of the three processes are related by (2).

Note that there is no assumption to this point that the spectra are rational, or otherwise finitely parametrized.

To formulate a solvable identification problem, a number of other assumptions are always invoked in considering dynamic factor models, and we will introduce them gradually and attempting a logical progression over the subsequent subsections, albeit with an attempt to justify their relevance.

2.1. The Goal of Modelling

Exactly what should be the goal of modelling in science is as much a philosophical as it is a scientific question, with answers often reflecting on the real-world utility of any model. Much the same question arises in connection with economic modelling; the very highly-cited paper Sims (1980) seeks to peer behind the mathematical formalism to understand exactly what relation the models of economics have to the real world.

For the purposes of this paper, we will simply assert that a model may be desired for various purposes which may not all be applicable for any one model at a particular time. Two of these for example might be forecasting of future values of a time series, and identification of business cycles. For modelling and especially forecasting purposes, often a canonical (stable, miniphase) spectral factor of ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)$ is desired. Such a spectral factor is to be constructed using the data $y_t^N$ with $0\le t\le T$ where T is some finite positive integer. The dimension N may in fact exceed T, indeed theoretical studies focus on letting N rather than T tend to infinity.

To make progress, further assumptions beyond the nesting property are needed. Some of these assumptions bear on the ability to somehow dispose of the contaminating influence of $\left({\xi }_t^N\right)$, and we will deal with them further below, and at that time clarify the strong and weak dependence notions. So assume temporarily that somehow we can work just with $\left({\chi }_t^N\right)$.

The Issue of Errors and Continuity

Observe that there are different ways to describe a stable miniphase spectral factor, see e.g., Hannan and Deistler (2012) and Caines (2018) for the relevant linear system theory background. If it has a rational transfer matrix, state-variable descriptions could be used (and then there is freedom to choose a coordinate basis). There is interest especially in those which are of minimal dimension. Again, if the spectral factor is rational, a left or right matrix fraction description could be used (a left matrix fraction description being an ARMA description, subject to rules about leading coefficient matrices). There is interest especially in those where the fraction is coprime (common factors are generally as unhelpful in a matrix transfer function as they are for the numerator and denominator of a scalar transfer function). If it is known to be autoregressive, an AR description may be used. A very different form again is provided by a series of graphical images of the values of the transfer matrix entries as a function of $\omega$; this might even be the most convenient form of model if the model does not have a rational transfer function matrix.

No matter what means are used to construct a model for the process $\left({\chi }_t^N\right)$, it is clear that one needs to be certain that a small error in (say) the lagged covariance sequence values (due say to use of sample rather than population statistics) will ultimately flow through to a small error in the quantities obtained in whatever particular application is envisaged when using the model. Since the model is the tool underpinning the application, it is reasonable to seek further assurance that the small error in the lagged covariance values (which in almost all cases are the raw data from which any model is constructed) will also give rise to a small error in the actual model of the minimum phase transfer function, be it frequency domain description, AR, ARMA, or general state-variable. Thus, if the model is being used for prediction, one ultimately wants a small error in the estimation variance of the predictions, and an accurate AR model may consequently be sought (assuming existence of an AR model is guaranteed). On the other hand, and even if an AR model should exist, if the model is being used to look for spectral peaks capturing some oscillatory phenomenon, or to understand the different bandwidths of the processes of interest, one is more likely to want a small error in the frequency domain description of the model, or even the pole positions of the model, as opposed to the AR coefficients.

Evidently then, some kind of a continuity is being sought (small errors in data used for generating the model should give small errors in the ultimate answer in the application for which the model is being used), but evidently also the particular parameters relevant to considering continuity, including continuity within the model description, are application dependent.

The most common sort of continuity in a model to which attention is given is the continuity of coefficients in a finite dimensionally parameterised model. One thinks of the class of systems of interest as collection of subsets, where the subsets are obtained by performing (data driven) model selection or by making various assumptions about integer parameters associated with the dimension of the innovations process, or a minimal state variable realization. However, such thinking may exclude certain systems at the start, and on occasions continuity of a frequency domain plot may be of more interest than continuity of a parameter in a model of some assumed order.

Because so much work in econometrics is linked to forecasting, there is a justified preference for using (when they exist) autoregressive models. The crucial part of the modelling to get right is to have accurate AR coefficients, which (in the event that Yule-Walker equations are used to derive the coefficients) probably requires a well-conditioned and thus more easily invertible matrix of lagged covariance coefficients, in addition to some suitable level of accuracy in the lagged covariance coefficients themselves. [More precisely, since on occasions a singular matrix may appear in the Yule-Walker equations, one will want a modified form of well-conditioning as indicated below.] One should note that the conditioning of the matrix in question can be related to an important property of the associated multivariate spectrum. Let ${\nu }_{max}$ denote the maximum over $\omega$ of the maximum eigenvalue of the Hermitian matrix ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)$, and let ${\nu }_{min}$ denote the minimum over $\omega$ of the smallest nonzero eigenvalue of ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)$. (Zero eigenvalues arise if the matrix is singular and these are for this purpose discarded; it is assumed the number of zero eigenvalues is the same for all $\omega$). The condition number of the block Toeplitz matrix inverted in a Yule-Walker calculation is known to approach ${\nu }_{max}/{\nu }_{min}$ as the size of the matrix goes to infinity, Miranda and Tilli (2000); Serra (1999), with modification in the singular case. Thus, ill-conditioning in a Yule-Walker calculation may be an automatic consequence of the nature of the spectrum being modelled.

In a sense, if we use an information criterion such as AIC for estimating the AR order p, such a problem is unlikely to occur Shibata (1980).

In contrast however to the use of a model for prediction, if bandwidth determination and oscillatory phenomena or accuracy of pole positions in the model are issues of concern, it is a well-established, indeed now old, fact in digital signal processing that autoregressive (and ARMA) models (known often in the relevant literature as direct models) can be very dangerous, in the sense that pole positions and frequency peaks can be extremely sensitive to the AR coefficients, e.g., Mitra (2011); Rabiner and Gold (1975). MATLAB (Signal Processing Toolbox) even allows the determination of “lattice-form” state-variable realizations of a prescribed transfer function that will have low or even minimum sensitivity in the frequency response to truncation or round-off errors associated with the coefficients in the realization, or minimum round-off errors in calculations using the coefficients in the realization to compute outputs from inputs.

Using an AR representation for the canonical spectral factor will likely be a bad idea in this situation: the numerical parameters in the model may need to be known to an unreasonable degree of accuracy. At the same time, as the signal processing literature indicates, there may be a state-variable model for which much less demanding accuracy requirements apply to the numerical parameters in the model, in order to learn from the model what is desired.

Note further that if multi-step prediction is in contemplation, rather than one-step-ahead prediction, recursion using AR coefficients will also be risky for those systems with pole positions which are very sensitive functions of the AR coefficients. In the event though that prediction was required for a certain fixed number of steps, h say, ahead, as opposed to an interval extending some distance into the future, one could potentially avoid recursion by directly estimating an h-step prediction model by minimizing the h-step prediction errors by least squares. See for example Bhansali and Ghosh (1999), which distinguishes the two approaches with the names “plug-in method” and “direct method”.

To sum up then, we must be concerned about two things:

• If there is a true underlying process with an associated canonical spectral factor, and a canonical spectral factor obtained using approximate values of lagged covariance coefficients, there will be an error between the two.
• We should use a description of the computed canonical spectral factor that is appropriate to the application, in the sense that the error mechanisms in the relevant computations will not cause significant errors in the quantities computed as part of the application.

In connection with the first issue of continuity of the canonical spectral factor frequency response, reference Anderson (1985) provides a sufficient condition for continuity: the derivative (with respect to $\omega$) of ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)$ must be bounded, to ensure that a small $L^{\infty }$ error in the spectrum gives rise to a small $L^{\infty }$ error in the canonical spectral factor frequency response. To relate this to a covariance sequence condition, observe that

$${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)=\sum _{s=-\infty }^{\infty }\mathbb{E}\left[{\chi }_{t+s}^N{\left({\chi }_t^N\right)}^{\top }\right]e^{-j\omega s}.$$

Formal differentiation yields

$$\frac{d}{d\omega }{\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)=\sum _{s=-\infty }^{\infty }(-j)s\mathbb{E}\left[{\chi }_{t+s}^N{\left({\chi }_t^N\right)}^{\top }\right]e^{-j\omega s}$$

and by the Weierstrass M-test, uniform convergence in $\omega$ will occur provided

$$\sum _{s=-\infty }^{\infty }s\parallel \mathbb{E}\left[{\chi }_{t+s}^N{\left({\chi }_t^N\right)}^{\top }\right]\parallel <\infty .$$

A sufficient condition ensuring this property is incorporated in the following:

Assumption 2.

The derivative (with respect to ω) of ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)$ is bounded, which is guaranteed if the covariance sequence $\mathbb{E}\left[{\chi }_{t+s}^N{\left({\chi }_t^N\right)}^{\top }\right]$ goes to zero faster than ${(1/s)}^{\alpha }$ for any $\alpha >2$.

As explained in Anderson (1985), in the absence of a derivative bound, an arbitrarily small variation of a scalar spectrum can produce an arbitrary phase shift anywhere in ($0,\pi )$ in the associated spectral factor (though the spectral factor magnitude is obviously simply given by taking the square root of the spectrum value at each frequency, and so continuity applies to the magnitude). Evidently, without the Assumption, the spectral factorization problem is almost ill-posed for any real world situation.

As noted above, this assumption is relevant to securing accuracy in frequency domain behavior of the miniphase spectral factor. It is at best indirectly relevant to securing accuracy in say an AR model being used for prediction purposes, if such is known to exist. Nevertheless, if by virtue of some other assumptions or argument, a stable AR model (and surely stability is a requirement) exists, then this imposes an exponential rate of decay (with increasing lag) on the lagged covariances, and the assumption is fulfilled. Thus, the assumption is not offensive to custom. Actually, the ratio ${\nu }_{max}/{\nu }_{min}$ is also partially relevant to the assumption: the greater is this ratio, the more sensitive appears to be the spectral factor frequency response to variations in the power spectrum, Anderson (1985).

2.2. Controlling the Number of Parameters to Be Estimated

It is clear that any canonical spectral factor of practical use has to be parametrised by a finite number of quantities, be they integers or real parameter values. One way to ensure this can occur is simply to lay down an assumption that the spectrum is rational:

Assumption 3.

For each N, ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)$ is rational and therefore bounded in $e^{j\omega }$.

The boundedness assumption rules out the inclusion of a deterministic component.

On the face of it, a rationality assumption is a leap of faith–there seems no reason grounded in physics or economics that suggests this is reasonable, apart perhaps from a belief in the general applicability of an Occam’s Razor principle. What makes it reasonable and less instrumentalist though is to recognize that a very minor and more reasonable variant on it will be satisfactory. Using the crucial fact that any covariance sequence satisfying Assumption 2 (corresponding to a stationary process) can be approximated arbitrarily closely in $L^p$ norm for any $p\in [1,\infty )$ by a covariance with rational spectrum3, it follows that if the true spectrum is not rational, under the decay condition of Assumption 2, the rationality assumption of Assumption 3 becomes in practical terms acceptable, as it can apply to a spectrum arbitrarily close to that of the true process. There is however a caution: while the approximation idea is valid for a fixed N, the assumption says nothing about what happens when N goes to infinity, and the next and less innocent assumption given below is needed to address the issue.

More specifically and given Assumption 3, there will in fact be three critical N-dependent integers associated with the spectral factor of ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)$:

• The dimension n of a minimal state-space realization of the canonical spectral factor (or its McMillan degree)
• The number of columns q in the minimal state-space factor, or the dimension of an innovations sequence (this being a white noise sequence, known in this context as the dynamic factor sequence)
• The rank r of ${\gamma }^N\left(0\right):=\mathbb{E}\left[{\chi }_t^N{\left({\chi }_t^N\right)}^{\top }\right]$, which is the dimension of the space spanned by ${\chi }_{1t},{\chi }_{2t},\dots ,{\chi }_{Nt}$.

While the existence of n (for fixed N) is a direct consequence of Assumption 3, the other two integers deserve some comment. The rationality assumption implies there is a rational canonical spectral factor for $\left({\chi }_t^N\right)$, thus

$$\begin{array}{cc}\hfill x_{t+1}& =F{x}_t+G{w}_{t+1}\hfill \\ \hfill {\chi }_t^N& =H^N{x}_t,\hfill \end{array}$$

(3)

where $F\in {\mathbb{R}}^{n\times n},G\in {\mathbb{R}}^{n\times q}$ with rank q, $H^N\in {\mathbb{R}}^{N\times n}$ with rank r, and the matrices $H^N$ are nested; this means we can write for n-dimensional row vectors $h_1,h_2,\dots ,$

$$H^N=\left[\begin{array}{c}{h}_1\\ h_2\\ ⋮\\ h_N\end{array}\right]$$

(4)

The process $\left(w_t\right)$ is zero mean white noise with $\mathbb{E}\left[w_t{w}_t^{\top }\right]=I_q$. Stability and the miniphase property correspond to

$$\rho \left(F\right)<1,\mathrm{rank}\left[\begin{array}{cc}I-Fz& -G\\ H^N& 0\end{array}\right]=n+q\mathrm{for}\left|z\right|\le 1$$

(5)

We note that one can also write

$${\chi }_t^N=H^N{(I-Fz)}^{-1}G{w}_t=K^N\left(z\right)w_t$$

(6)

with obvious definition of the $N\times q$ transfer matrix $K^N\left(z\right)$. As is well known, such a stable and miniphase transfer function satisfying the spectral factorization equation ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)=K^N\left(e^{j\omega }\right){\left(K^N\left(e^{-j\omega }\right)\right)}^{\top }$ is unique up to right multiplication by a real orthogonal matrix, here of dimension $q\times q$, but through normalization can be made unique.

Since $G{G}^{\top }$ is of size $n\times n$, and the number of columns of G, viz. q, is minimal in a canonical spectral factor, it is immediate that $q\le n$. Further, since

$${\gamma }^N\left(0\right)=\mathbb{E}\left[H^{N}x\left(0\right)x^{\top }\left(0\right){\left(H^N\right)}^{\top }\right],$$

we see that r is bounded by the dimension of the square matrix $\mathbb{E}\left[x\left(0\right)x^{\top }\left(0\right)\right]$, i.e., $r\le n$. However, also, ${\gamma }^N\left(0\right)$ is related to the spectrum by

$${\gamma }^N\left(0\right)=\frac{1}{2\pi }{\int }_{-\pi }^{\pi }{\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)d\omega$$

which shows that $r\ge q$.

All three integers evidently depend on N, and one of our interests is to let N vary. However, using these integers, we introduce by fiat an almost natural strengthening of the nesting property of the underlying processes and their spectral densities and thereby almost ‘define away’ the dependence on N:

Assumption 4.

There exists a positive integer $N_0$ such that for all $N\ge N_0$, the integers $n,q,r$ are independent of N.

By way of heuristic explanation or even justification, we comment that as N increases, the parameter count in the canonical spectral factor increases, linearly with N. However, the assumption means that the number of independent driving signals which are the ultimate cause of the output of the canonical spectral factor (the dimension of what is termed the dynamic factors of the model) remains constant at q, the complexity (minimal state-space dimension) n of the model stays constant, and once the output vector has become sufficiently big (of size r), all further outputs are linear combinations of other ones, and not even delayed versions of others. There is almost no new information in these further outputs4, other than that explaining how they are related to earlier outputs. Implicit then in the assumption is a sort of Occam’s Razor principle, limiting the rate of increase in complexity (the parameter count effectively) of the model to the least amount consistent with increase in N.

A physical analogy arises if one considers ${\chi }_t$ as like a vector of vibrations arising in a piece of machinery, or electrical thermal noise at the ports of a linear circuit. There may be a limited number q of sources of the vibrations (resistors in the circuit analogy), and one could add more and more sensors (corresponding to N) to the set-up. Their presence would not change the underlying dynamical behavior of the machinery or circuit, governed by a set of differential equations with ‘complexity’ n. At some point when adding sensors (but with no associated sensor noise since that is to be bundled in with ${\xi }_t$ in this analogy), say when $N=r$, essentially no new information is provided by an additional sensor: the additional sensor just provides a linear combination of the measurements from other sensors.

The assumption that q is independent of N is shared by all works on dynamic factor models. A large majority in that literature also shares the assumption that r (and n) is independent of N. However, as the following simple example shows, q can be independent of N without r being independent of N. Let

$${\chi }_{it}={(1-{\alpha }_{i}z)}^{-1}{w}_t=w_t+{\alpha }_i{w}_{t-1}+\cdots ,$$

where $w_t$ is a scalar white noise, so that in this case $q=1$. Now, if the ${\alpha }_i\in [0,.9]$ and ${\alpha }_i\ne {\alpha }_j$ for all $i\ne j$, then the dimension of the space spanned by ${\chi }_{it}$, $i=1,\dots ,N$ is N, so that Assumption 4 does not hold, see Forni et al. (2015) for details.

Henceforth, we shall neglect the dependence of $n,q,r$ on N for small values of N.

We summarize the inequalities linking the three key integer parameters and the minimum allowed value $N_0$ of N in the following theorem:

Theorem 1.

With notation as above and under Assumptions 2–4, there holds $q\le r\le n$ and $r\le N_0$, (with $r,q,n$ independent of $N\ge N_0$).

Actually, $r>q$ has always been observed in empirical applications to large macroeconomic data sets, Barigozzi et al. (2021). One can observe also that as N grows, the number of parameters grows linearly with N, in fact as $O\left(nN\right)$.

Given covariance data for $\left({\chi }_t^N\right)$, one way the integer r can be obtained is to use the fact that for large enough N,

$$r=\mathrm{rank}\mathbb{E}\left[{\chi }_t^N{\left({\chi }_t^N\right)}^{\top }\right]$$

Sample averages have to replace the population averages, and methods such as singular value decomposition could be used to execute the calculations.

For future reference, we note that the equation

$${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)=K^N\left(e^{j\omega }\right){\left(K^N\left(e^{-j\omega }\right)\right)}^{\top }$$

implies that the McMillan degree of $K^N$ (the dimension of a minimal realization of $K^N$, i.e., n) is one half the McMillan degree of ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)$, or equal to the McMillan degree of the transfer matrix with associated impulse response the lagged covariance sequence ${\gamma }^N\left(0\right),{\gamma }^N\left(1\right),\dots ,$ with ${\gamma }^N\left(s\right)=\mathbb{E}\left[{\chi }_t^N{\left({\chi }_{t-s}^N\right)}^{\top }\right]$, see, e.g., Anderson (1969).

2.3. Strong and Weak Cross-Dependence

Strong dependence captures the notion that as N tends to infinity, the effect of the input signals, and all the dynamic complexity showing up in the state variable, continue to show up in a nondecreasing way in the progressively added output components. This can be captured formally by

Assumption 5.

For all $\omega \in [-\pi ,\pi ]$, the largest q eigenvalues of ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)$ diverge linearly with N as $N\to \infty$, while the other eigenvalues are zero. In addition, the r largest eigenvalues of $T^{-1}{\sum }_1^T{\chi }_t^N{\left({\chi }_t^N\right)}^{\top }$ diverge to infinity linearly with N, while the other eigenvalues are bounded.

Note that the second part of the assumption is not a consequence of the first part, Lippi et al. (2022).

Weak dependence in a sense is the converse. The different components of ${\xi }_t^N$ do not ‘reinforce’ one another due to some common dependences, but rather have very limited dependence between each other:

Assumption 6.

The largest eigenvalue of ${\mathrm{\Phi }}_{\xi }^N\left(e^{j\omega }\right)$ is bounded for all N and all $\omega \in [-\pi ,\pi ]$.

Evidently, for any fixed T, the larger N becomes, the greater is the ‘signal-to-noise’ ratio between ${\chi }_t^N,t\in \{0,1,2,\dots ,T\}$ and ${\xi }_t^N,t\in \{0,1,2,\dots ,T\}$. Motivated by communication engineering, ${\lambda }_{max}\left({\mathrm{\Phi }}_{\chi }\left(e^{j\omega }\right)\right)/{\lambda }_{max}\left({\mathrm{\Phi }}_{\xi }\left(e^{j\omega }\right)\right)$ can be regarded as a surrogate for signal-to-noise ratio. This might prompt us to contemplate the possibility of ‘denoising’, i.e., separating out the $\left({\chi }_t^N\right)$ process from $\left(y_t^N\right)$ for large enough N. Note that in using terminology with the word ‘noise’, we acknowledge that this interpretation of the $\left({\chi }_t^N\right)$ process will not be appropriate in some applications.

There are similar but not identical ways to introduce assumptions of this type which, taken in conjunction with the earlier assumptions, are just as suitable. Discussion with comparison can be found in Lippi et al. (2022).

This concludes the setting out of assumptions on the processes being examined. We turn now to the procedures used for modelling.

3. Obtaining the Common Component Process from the Measurement Process

The key to achieving an additive decomposition of the measurement process is to appeal to the strong dependence of the processes $\left({\chi }_{it}^N\right)$ and the weak cross-dependence of the processes $\left({\xi }_{it}^N\right)$.

It was shown in Stock and Watson (2002), Forni et al. (2009) that under the earlier assumptions a consistent estimator of ${\chi }_t^N$ can be obtained by using a static PCA of $y_t^N$. At a similar time, Bai and Ng (2002) demonstrated consistent estimation of r; the work relies also on using asymptotic principal components, is nonparametric, and does not require a low value of n. As for the input dimension of the canonical spectral factor of the common component process, Hallin and Liška (2007), Amengual and Watson (2007), Bai and Ng (2007), Onatski (2009) allow estimation of q.

As noted in Forni et al. (2009), it is evident that for all $N\ge N_0$, the dimension of the space, call it ${\mathcal{S}}_t^N$, spanned by ${\chi }_{it},i=1,2,\dots ,N$ is r. Let $\left(f_t\right)$ be an r-dimensional process such that $f_t$ forms a basis for ${\mathcal{S}}_t^N$ and $f_t=S{\chi }_t^N$ where S is an $r\times N$ selector matrix independent of t and picking out a spanning set of entries5 of ${\chi }_{it}$ for i limited to the set $i=1,2,\dots ,N_0$. It follows also that

$${\chi }_t^N=L^N{f}_t$$

for some factor loading matrix $N\times r$ matrix $L^N$ independent of t. The entries of $f_t$ are often termed minimal static factors.

Actually, in place of the selector matrix S, one can use any left multiple $TS$ with T nonsingular and $r\times r$.6 A minimal static factor is then seen to be unique up to premultiplication by a constant nonsingular matrix T, with the factor loading matrix unique up to postmultiplication by $T^{-1}$. Note that the way $L^N$ has been constructed, the successive (with N) matrices are nested. It is possible to eliminate the nonuniqueness using normalizing procedures, see, e.g., Lippi et al. (2022).

The connection with the earlier model of (3) is straightforward to see. Evidently, for some $r\times n$ matrix $\tilde{H}$, there holds

$$f_t=\tilde{H}{x}_t$$

(7)

with

$$H^N=L^N\tilde{H}$$

(8)

The process $\left(f_t\right)$ in a sense is now the real process of interest. Of course, its dimension does not grow with N. Its determination however depends on using large N, in order that $\left({\chi }_t^N\right)$ can be extracted from $\left(y_t^N\right)$, and then $\left(f_t\right)$ is extracted from $\left({\chi }_t^N\right)$.

Importantly, it now becomes computationally practical to consider the determination of the integer n from the lagged covariance sequence determined by $f_t$. By way of background, recall, see, e.g., Hannan and Deistler (2012); Ho and Kalman (1966) that n is available as the rank of all sufficiently large square block Hankel matrices formed from the covariance data

$$\mathbb{E}\left[f_t{f}_t^{\top }\right],\mathbb{E}\left[f_t{f}_{t-1}^{\top }\right],\mathbb{E}\left[f_t{f}_{t-2}^{\top }\right],\dots$$

Singular value decomposition (SVD) is a computationally effective way to determine the rank. Of course, given that the covariances are actually determined using sample rather than population averages, and using finite rather than infinite values of N and T to extract $\left({\chi }_t^N\right)$ from $\left(y_t^N\right)$, judgment will have to be used to deal with the computational errors; SVD really is a tool for estimating rather than determining the rank. Further details appear in a later section.

Two Alternative Approaches to Additive Decomposition of the Spectral Matrix

For completeness, we summarize here two ways, see Lippi et al. (2022), whereby the spectra ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)$ and ${\mathrm{\Phi }}_{\xi }^N\left(e^{j\omega }\right)$ can be derived from ${\mathrm{\Phi }}_y^N\left(e^{j\omega }\right)$.

First, an assumption is made that the idiosyncratic process vector is a vector of independent second order processes. Under such a special assumption, it is in principle possible using partial fraction expansion of the power spectrum ${\mathrm{\Phi }}_y^N\left(e^{j\omega }\right)$ to separate the two spectra ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)$ and ${\mathrm{\Phi }}_{\xi }^N\left(e^{j\omega }\right)$. Second, dropping the rationality assumption but retaining the assumption that the idiosyncratic process vector is a vector of independent second order AR processes, and by considering $(q+1)\times (q+1)$ submatrices of ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)$ which contain precisely one diagonal entry (this requires $N>2q$), and which therefore have all entries known from ${\mathrm{\Phi }}_y^N\left(e^{j\omega }\right)$ except for that diagonal entry, the singularity of the matrix generically enables recovery of the diagonal entry.

The diagonality assumption on ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)$ in both procedures is stronger than the weak cross-dependence and can be seen as a misspecification of the general model. However, Doz et al. (2012) argues that: (i) under diagonal ${\mathrm{\Phi }}_{\xi }^N\left(e^{j\omega }\right)$ the model can be estimated by Maximum Likelihood, (ii) as N tends to infinity the effect of the misspecification vanishes. Whether the effect also vanishes in the procedures outlined above is unknown.

4. Tall Stable Miniphase Spectral Factors and Singular Rational Spectra

An important property is (for $N>q$) the tallness of the transfer function matrix $H^N{(I-Fz)}^{-1}G$, which is linked with the singularity of ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)$. Likewise, the transfer function matrix between the innovations and the minimal static factor process $\left(f_t\right)$, being given by $\tilde{H}{(I-Fz)}^{-1}G$ is $r\times q$, and while always $r\ge q$, as noted earlier $r>q$ is normally encountered in practice.

In this section, we recall results appearing in Anderson and Deistler (2008, 2009); Anderson et al. (2012b); Chen et al. (2012) stemming from this property. The first key result deals with the zeros of the generating model.

Theorem 2

(Anderson and Deistler (2008)). Consider the stable, miniphase model (3) with $N>q$ and assume that the entries of the matrices $F,G,H^N$ take generic values. Then the model has the zero free property7

$$rank\left[\begin{array}{cc}I-Fz& -G\\ -H^N& 0\end{array}\right]=n+q\forall z$$

(9)

In Theorem 2 we make the assumption that the entries of $F,G,H^N$ take generic values, with the requirements of stability and the miniphase property. The inclusion of these requirements is dealt with further below. Genericity of course can only be defined with respect to a well specified parameter space. In Theorem 2, in agreement with the Anderson and Deistler’s papers mentioned above, the parameter space is implicitly defined as a minimum-restriction set. Precisely, the parameter space is an open subset of ${\mathbb{R}}^g$, contained in ${\mathbb{R}}^g-G$, where g is the total number of entries of $F,G,H^N$ and G is the subset of ${\mathbb{R}}^g$ where either the stability or the miniphase condition does not hold. If the system is zero free, that is consistent with the miniphase requirement, and it is clear that if the system is zero free for one set of parameter values, it will be free for almost all, given that it is a rank condition. (Note also that an inessential variant to the count can occur if some of the entries are fixed at zero, or one, as can occur if a square matrix is a companion matrix or a block companion matrix for example. Such fixed entries are not part of the count.) In other words, apart from the stability and miniphase conditions, we assume that each parameter can vary independently of the others. Throughout the paper genericity will always be implicitly defined with respect to such minimum-restriction parameter spaces. See however the Conclusions for further considerations on this issue.

Zero freeness with $N>q$ implies a number of properties:

• Suppose that for polynomial matrices ${\overline{A}}^N\left(z\right)\in {\mathbb{R}}^{N\times N}\left[z\right],{\overline{B}}^N\left(z\right)\in {\mathbb{R}}^{N\times q}\left[z\right]$ with ${\overline{A}}^N$ and ${\overline{B}}^N$ left coprime8 there holds $K^N\left(z\right)={\left({\overline{A}}^N\right)}^{-1}{\overline{B}}^N$. Then zero freeness implies ${\overline{B}}^N\in {\mathbb{R}}^{N\times q}$ has full rank for all z. Accordingly, there exists a polynomial unimodular (constant determinant) common premultiplier $U^N\left(z\right)$ of ${\overline{A}}^N,{\overline{B}}^N$ such that a new polynomial fraction description of $K^N\left(z\right)$ exists with

  $$K^N\left(z\right)={\left(U^N\left(z\right){\overline{A}}^N\left(z\right)\right)}^{-1}\left(U^N\left(z\right){\overline{B}}^N\left(z\right)\right)={\left(A^N\left(z\right)\right)}^{-1}{B}^N$$

  where $B^N$ is a constant matrix. We have stressed earlier that N can increase; note though that $det{A}^N\left(z\right)$ will have degree n as a polynomial in z for all $N\ge N_0$.
• The system (3) is left invertible with unknown initial state, that is, there exists an integer $L\le n$ such that from the sequence ${\chi }_k^N,{\chi }_{k+1}^N,\dots ,{\chi }_{k+L-1}^M$ the state $x_k$ and the sequence $w_{k+1},w_{k+2},\dots ,w_{k+L-1}$ can be determined.

Other characterizations and properties to be found in Anderson and Deistler (2008) deal with the Smith-McMillan form, the impulse response of $K^N\left(z\right)$, and the absence of nontrivial so-called output-nulling subspaces.

There is of course a very minor restatement of the above Theorem and the two following comments based on the canonical spectral factor generating $f_t$ rather than ${\chi }_t^N$. It too is (normally) zero free, due to the tallness inequality $r>q$. However, in working with an AR description of $f_t$ and setting ${\chi }_t^N=L^N{f}_t$, rather than working with an AR description of ${\chi }_t^N$, the reduction in parameter count is substantial, and varies as $O\left(N\right)$ rather than $O\left(N^2\right)$.

To this point, the stability of the model (3) and its connection with a spectrum have not been appealed to. To bring in those ideas, continue to assume the underlying stable miniphase spectral model of ${\mathrm{\Phi }}_{\chi }^N$ is (3) and the matrices in a state-variable description are generic.

Theorem 3

(Anderson and Deistler (2008, 2009); Anderson et al. (2012b); Chen et al. (2012)). Adopt the same hypotheses as Theorem 2. Then

1. 

The spectrum${\mathrm{\Phi }}_{\chi }^N$can be generated by a singular autoregressive process

$$[I+A_1^{N}z+A_2^N{z}^2+\cdots +A_m^N{z}^m]{\chi }_t^N=B^N{w}_t$$

(10)

where$B^N\in {\mathbb{R}}^{N\times q}$has rank q and$w_t$is a unit covariance innovations process; the$A_i^N$and$B^N$are computable in a finite number of rational calculations from the lagged covariances$\mathbb{E}\left[{\chi }_t^N{\left({\chi }_t^N\right)}^{\top }\right],\mathbb{E}\left[{\chi }_t^N{\left({\chi }_{t-1}^N\right)}^{\top }\right],\dots$.

2. 

Suppose the rational spectral matrix${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)$is written in the form

$${\overline{H}}^N{(I-z\overline{F})}^{-1}{\overline{Q}}^N+{\left({\overline{Q}}^N\right)}^{\top }{(I-z^{-1}{\overline{F}}^{\top })}^{-1}{\left({\overline{H}}^N\right)}^{\top }-{\overline{H}}^N{\overline{Q}}^N,z=e^{j\omega }$$

where the triple$\{\overline{F},{\overline{H}}^N,{\overline{Q}}^N\}$is minimal (controllable and observable) in the state-variable sense and all eigenvalues of$\overline{F}$have magnitude less than 1, such matrices being computable in a finite number of rational calculations from the lagged covariances$\mathbb{E}\left[{\chi }_t^N{\left({\chi }_t^N\right)}^{\top }\right],\mathbb{E}\left[{\chi }_t^N{\left({\chi }_{t-1}^N\right)}^{\top }\right],$…. Then

(a) 

there exists${\overline{G}}^N$such that the stable miniphase (canonical) spectral factor can be written in the form${\overline{H}}^N{(I-z\overline{F})}^{-1}{\overline{G}}^N$

(b) 

${\overline{G}}^N$is computable in a finite number of rational calculations from$\overline{F},{\overline{H}}^N$and${\overline{Q}}^N$, and

(c) 

the triple$\{\overline{F},{\overline{H}}^N,{\overline{G}}^N\}$is minimal.

3. 

In the event that${\overline{G}}^N$defined above is independent of N, say$G={\overline{G}}^N$, then${\overline{Q}}^N=P{\left({\overline{H}}^N\right)}^{\top }$for some positive definite P which is independent of N, with P the unique solution of$P-\overline{F}P{\overline{F}}^{\top }=G{G}^{\top }$.

The following statement is a quite obvious consequence of the previous Theorems:

Theorem 4.

Consider the ARMA process

$$A\left(z\right)x_t=B\left(z\right)u_t,$$

where $A\left(z\right)$ is an $m\times m$ stable polynomial matrix, $B\left(z\right)$ is an $m\times q$ polynomial matrix. Suppose that $m>q$, so that $x_t$ is a singular ARMA process. Then for generic values of the coefficients of the entries of $B\left(z\right)$, $B\left(z\right)$ is zeroless, so that $x_t$ has a finite AR representation.

The relevance of AR processes in our context is not a new idea of course, see Stock and Watson (2005), Forni et al. (2009), Stock and Watson (2016), Forni et al. (2020). These papers apply Theorem 4, in which $x_t$ is replaced by $f_t$, while singularity is a consequence of the assumption that $r>q$, so that we are again working with a generically zero-free model.

The first statement of the Theorem 3 above is little more than a restatement of the first statement of Theorem 2. As a side comment, note that $B^N$ can be uniquely specified by multiplying on the right by an orthogonal matrix rendering it triangular. Reference Anderson et al. (2012b) contains further information concerning exploitation of the freedoms that remain after such a normalization, and the associated parameter counts, including the value of m. The row degrees in particular of the matrix $I+A_1^{N}z+\dots +A_m^N{z}^m$ can be minimized in a certain sense. Notwithstanding, the number of nonzero parameters in the singular AR description is essential proportional to $N^2$, which is an unattractive feature, contrasting with the linear dependence on N of the parameter count for a state-variable model. In a different direction, a different aspect of the awkwardness of working with the singular AR model is that the nesting property for $N\ge N_0$ is not simply captured. Of course, in the event one seeks, as noted above, an AR model for $\left(f_t\right)$ rather than $\left({\chi }_t^N\right)$ the difficulty of working with large N is no longer present. Not only does one have ${\chi }_t^N=L^N{f}_t$, with N is now responsible for only linear growth in the number of parameters, but the nesting property holds, in particular being handled by $L^N$.

The second statement is more subtle: the fact that the same pair $\overline{F},{\overline{H}}^N$ can be used for state variable realizations of the covariance sequence associated with ${\mathrm{\Phi }}_{\chi }^N\left(e^{j\omega }\right)$ and the associated stable miniphase spectral factor is an old though not well-known result, see, e.g., Anderson (1969); it is not restricted to the zero-free (autoregressive) case. The multivariate version of the result is not easy to prove, in contrast to the converse (which states that a triple ‘realizing’ the covariance sequence can inherit two of the component matrices of the triple ‘realizing’ the spectral factor). What is novel in the zero-free case is that the determination of the spectral factor requires no more than a finite number of rational calculations. These calculations involve solving a recursive discrete time matrix Riccati equation involving $\overline{F},{\overline{H}}^N$ and ${\overline{Q}}^N$ over a finite interval of length $O\left(n\right)$ or (what is nontrivially equivalent) performing a Cholesky decomposition on a (square) Toeplitz matrix with $O\left(nN\right)$ rows, formed from $\overline{F},{\overline{H}}^N$ and ${\overline{Q}}^N$, see Anderson and Deistler (2009). In the event that N is large and the realization of the covariance sequence is available, this is an attractive idea: with both approaches, the computational complexity is linear in N. Of course, it is the zero-free property of the canonical spectral factor that assures the finite-time convergence of the Riccati equation.

The third statement simply captures, and is also an extension of, the nesting property. Though not contained in Anderson and Deistler (2009), a proof is contained in Appendix A.

5. Obtaining Parameters in the AR and State-Variable Models

We suppose in this section that from the sample data $\left(y_t^N\right)$ collected over a time interval and with a large value of N, using PCA or an equivalent, an estimate of $\left({\chi }_t^N\right)$ is obtained over effectively the same interval. This is not even a sample of a process over a finite interval, but only an estimate of it. Finding a model for the ‘true’ process $\left({\chi }_t^N\right)$ based on the estimates then has to proceed using the assumptions that $n,q,r$ are known, that $N\ge N_0$, and that closeness of the sample-derived lagged covariance data and population lagged covariances assure closeness of the associated canonical spectral factors. Thus, the Assumptions introduced earlier need to all be in force.

Consider first the determination of the coefficient matrices in an AR model for $\left({\chi }_t^N\right)$. What we say will be very little modified if we choose instead to model $f_t$, and so we omit any detailed discussion of this alternative in this section.

The fact that the spectrum of interest is generated by a (singular) autoregression as in (10) suggests that a stable miniphase spectral factor may be obtainable by Yule-Walker equations–the usual tool for obtaining an AR model (when one is known to exist) from covariance data. Indeed, this possibility was examined in Deistler et al. (2010). The Yule-Walker equations in general take the form

$$[A_1^N{A}_2^N\dots A_N^N]{\mathrm{\Gamma }}_N=[{\gamma }_1{\gamma }_2\dots {\gamma }_N]$$

(11)

where ${\gamma }_j$ is the lagged covariance $\mathbb{E}\left[{\chi }_t^N{\left({\chi }_{t-j}^N\right)}^{\top }\right]$ and ${\mathrm{\Gamma }}_N$ is the $N\times N$ block Toeplitz matrix with first row $[{\gamma }_0{\gamma }_1\dots {\gamma }_{N-1}]$. However, a difficulty immediately arises. For a singular covariance, there is no guarantee that ${\mathrm{\Gamma }}_N$ is nonsingular. Nevertheless, it is established in Deistler et al. (2010) that the choice

$$[A_1^N{A}_2^N\dots A_N^N]=[{\gamma }_1{\gamma }_2\dots {\gamma }_N]{\mathrm{\Gamma }}_N^{†}$$

(12)

where use of a pseudoinverse to replace the regular inverse ensures that the resulting AR model defines a stationary process, i.e., the associated matrix polynomial $[I+A_1^{N}z+A_2^N{z}^2+\dots A_N^N{z}^N]$ has all determinantal zeros in $\left|z\right|>1$, which is a rephrasing of the stability condition. The paper Chen et al. (2011), which imposes fewer assumptions than are incorporated here, takes this further: under some circumstances, there may be unit circle determinantal zeros, corresponding to a deterministic component of the process. (The assumptions in this paper rule out this possibility.) Further, if one does not use a pseudoinverse (which delivers a minimum norm solution of (11)) but uses some other procedure to define the matrices $A_i^N$ to satisfy (11), the AR system so defined may not be stable. Note that the quantities ${\gamma }_i$ will actually be replaced by estimates which are consistent as $N,T\to \infty$.

For completeness, we recall the potential numerical hazards arising when the key Toeplitz matrix is ill-conditioned, as discussed in Section 2.1.

The procedure for determining a state-space model is of course different. The starting point is the sample calculations giving estimates of $\mathbb{E}\left[{\chi }_t^N{\left({\chi }_{t-s}^N\right)}^{\top }\right]$ for $s=0,1,2,\dots$. If these were exact population second moments corresponding to a system satisfying the assumptions, particularly the rational spectrum assumption, an algorithm due to Ho and Kalman (1966) will determine a minimal dimension state variable triple realizing the sequence, i.e., a triple $\{\overline{F},{\overline{Q}}^N,{\overline{H}}^N\}$ for which

$${\overline{H}}^N{\overline{F}}^s{\overline{Q}}^N=\mathbb{E}\left[{\chi }_t^N{\left({\chi }_{t-s}^N\right)}^{\top }\right],s=0,1,2,\dots$$

with $\overline{F}$ of minimal dimension. (Such a triple is only unique up to a nonsingular transformation T, replacing $\{\overline{F},{\overline{Q}}^N,{\overline{H}}^N\}$ by $\{T\overline{F}{T}^{-1},T{\overline{Q}}^N,{\overline{H}}^N{T}^{-1}\}$). Subsequently, Zeiger and McEwen (1974) presented a variant that generates approximate finite dimensional realizations of low dimension from approximate data, as needed here. The algorithm uses truncated singular value decomposition, which means that for $O\left(N\right)$-dimensioned matrices of much lower rank than N the complexity appears to be $O\left(N\right)$. Of course, once the triple $\{\overline{F},{\overline{Q}}^N,{\overline{H}}^N\}$ has been determined, the construction of G using the ideas of Theorem 3 can be attempted.

Once again, it may well make more sense to seek a state variable triple realizing the sequence $\mathbb{E}\left[f_t{\left(f_{t-s}\right)}^{\top }\right],s=0,1,2,\dots$, rather than $E\left[{\chi }_t^N{\left({\chi }_{t-s}^N\right)}^{\top }\right],s=0,1,2,\dots$

The bulk of the suggested algorithms in this subsection have not been simulated with real data. However, even if they do not perform particularly well, they offer the opportunity to obtain a first estimate of a model, using which an iterative algorithm to refine the estimate can be used, Shumway and Stoffer (1982). This algorithm, based on expectation maximization, is a maximum likelihood algorithm which at each step always increases the likelihood, and necessarily therefore reaches a maximum, though it may not be a global maximum. The use of a well-chosen initialization can however promote that desired outcome.

The algorithm of Shumway and Stoffer (1982) appears to offer no straightforward technique for incorporating a rank constraint on the input coupling matrix G of the canonical spectral factor, or simply the zero-free property. Consequently, while it uses Kalman filtered and smoothed estimates, it is not immediately clear how one should exploit the notion that such estimates in the zero-free case can be exactly determined using a finite number of measurements as set out in the commentary immediately following Theorem 3.

6. Mixed Frequency Systems

Some econometric time series may be collected monthly, others may be collected quarterly. A methodology is needed to allow graceful combination of the two periodicities, and the modelling of high dimension time series. Starting with this section, there appears a summary of ideas developed by authors Anderson and Deistler with their colleagues in treating such problems. Note importantly that mixed frequency problems are of interest outside the ‘high dimensionality’ context, and some limited conclusions are drawn here that do not appeal to properties such as zero freeness.

The historical roots of the applicable techniques lie in signal processing, which has treated with a technique termed ‘blocking’ or ‘lifting’ problems involving time series with integrally related periodicities, with a view to bringing to bear tools applicable to single frequency systems, such as rational transfer functions, and their descriptions and properties involving for example state-variable realizations, poles and zeros, see, e.g., (Mitra 2011, see Chapters 13 and 14). The book Bittanti and Colaneri (2009) also contains, or contains references to, many of the relevant ideas, and has a more control systems flavour.

6.1. Blocked Linear Systems

We draw on Chen et al. (2012) for initial ideas, which cover change of periodicity as opposed to treatment of simultaneous multiple periodicities. Once change of periodicity is understood, it is straightforward to introduce multiple periodicities.

Consider a time-invariant finite-dimensional linear system with input sequence $u_t$ and output sequence $y_t$ given by $x_{t+1}=A{x}_t+B{u}_t,y_t=C{x}_t+D{u}_t$ with $t=0,1,2,\dots$. Imagine then that blocks of J successive inputs and outputs are collected into N-vectors $U_t={[u_t^{\top }{u}_{t-1}^{\top },\dots ,u_{t-J+1}^{\top }]}^{\top }$ and $Y_t={[y_t^{\top }{y}_{t-1}^{\top },\dots ,y_{t-J+1}^{\top }]}^{\top }$. Then one can write the following equations for a ‘blocked’ or ‘lifted’ system:

$$\begin{array}{cc}\hfill x_{t+1}& =A_b{x}_{t-J+1}+B_b{U}_t\hfill \\ \hfill Y_t& =C_b{x}_{t-J+1}+D_b{U}_t\hfill \end{array}$$

(13)

where

$$\begin{array}{cc}\hfill A_b& =A^J{B}_b=[BAB\dots A^{J-1}B]\hfill \\ \hfill C_b& =\left[\begin{array}{c}C{A}^{J-1}\\ C{A}^{J-2}\\ ⋮\\ CA\\ C\end{array}\right]D_b=\left[\begin{array}{cccc}D& CB& \dots & C^{J-2}B\\ 0& D& \dots & C{A}^{J-3}B\\ ⋮& ⋮& \ddots & ⋮\\ 0& 0& \dots & CB\\ 0& 0& \dots & D\end{array}\right]\hfill \end{array}$$

(14)

These are equations of a time-invariant linear system, of different input and output dimension but the same state dimension as the original. One can regard this second system as updating every J time instants, thus at (for example) $t=0,J,2J,\dots .$ Note that even if the ‘direct feedthrough’ term $D{u}_t$ is absent in the unblocked system as when $D=0$, a nonzero direct feedthrough term does appear in the blocked system in general. The transfer function of the blocked system is $D_b+C_{b}z{(I-z{A}_b)}^{-1}{B}_b$.

The paper Chen et al. (2011) records many connections between the original and ‘blocked’ transfer functions and state-variable realizations (some being novel to the paper), in particular, the carrying over of controllability and observability and minimality, connections between the normal ranks, a formula connecting the first column of the blocked transfer function and the original transfer function, relations between the poles and, separately, the finite and infinite zeros of the two transfer functions, implying that one is stable (or miniphase) if and only if the other is stable (or miniphase), and one is zero free if and only if the other is zero free.

In this connection, note that if the dimension of $y_t$ exceeds that of $u_t$, the same is true of the dimensions of $Y_t$ and $U_t$. Moreover, if the original $A,B,C,D$ matrices are generic, the unblocked and therefore the blocked system are zero free.

The property of zero-free carryover becomes relevant also in considering time series with multiple periodicities.

6.2. Systems with Multiple Periodicities in Their Outputs

Again our starting point is the linear system defined by $x_{t+1}=A{x}_t+B{u}_t,y_t=C{x}_t+D{u}_t$, with $t=0,1,2,\dots$. However we now suppose that the vector $y_t$ is partitioned into two subvectors of dimensions $p_1$ and $p_2$, thus

$$y_t=\left[\begin{array}{c}{y}_t^f\\ y_t^s\end{array}\right]$$

The process $\left(y_t^f\right)$ is observed at all time instants, while the process $\left(y_t^s\right)$ is observed at time instants $0,J,2J,\dots$ even though it is known to exist for all time instants. From the observer’s point of view, $\left(y_t^f\right)$ and $\left(y_t^s\right)$ constitute ‘fast’ and ‘slow’ processes. Evidently the observed output comprises two (in general multivariate) time series with different integrally related periodicities. As set out in, e.g., Anderson et al. (2016a, 2016b), We can establish a single frequency system by blocking the fast process $\left(y_t^f\right)$.

In more detail, consistent with the partition of $y_t$, we also adopt the corresponding partitions

$$C=\left[\begin{array}{c}{C}^f\\ C^s\end{array}\right]D=\left[\begin{array}{c}{D}^f\\ D^s\end{array}\right]$$

To obtain a blocked model, we block $u_t$ as before, but only use the observed values within each $y_t$ for blocking, which means that for $t=0,J,2J,\dots$,

$$U_t=\left[\begin{array}{c}{u}_t\\ u_{t-1}\\ ⋮\\ u_{t-J+1}\end{array}\right]Y_t=\left[\begin{array}{c}{y}_t^f\\ y_{t-1}^f\\ ⋮\\ y_{t-J+1}^f\\ y_t^s\end{array}\right]$$

(15)

It is straightforward to obtain matrices $A_b,C_b,C_b,D_b$ constituting a state-variable realization for the system taking the sequence $U_0,U_J,\dots$ to $Y_0,Y_J,\dots$, see Zamani et al. (2011). In particular,

$$\begin{array}{cc}\hfill A_b& =A^J{B}_b=[BAB\dots A^{J-1}B]\hfill \end{array}$$

(16)

$$\begin{array}{cc}\hfill C_b& =\left[\begin{array}{c}{C}^f{A}^{J-1}\\ C^f{A}^{J-2}\\ ⋮\\ C^{f}A\\ C^f\\ C^s\end{array}\right]D_b=\left[\begin{array}{cccc}{D}^f& C^{f}B& \dots & C^f{A}^{J-2}B\\ 0& D^f& \dots & C^f{A}^{J-3}B\\ ⋮& ⋮& \ddots & ⋮\\ 0& 0& \dots & C^{f}B\\ 0& 0& \dots & D^f\\ D^s& C^{s}B& \dots & C^s{A}^{J-2}B\end{array}\right]\hfill \\ \hfill x_{t+1}& =A_b{x}_{t-J+1}+B_b{U}_t\hfill \\ \hfill Y_t& =C_b{x}_{t-J+1}+D_b{U}_t\hfill \end{array}$$

This system is self-evidently a purely time-invariant system, governed by a single period which happens to be that defined by the observed slow process $\left(y_t^s\right)$, despite the system containing within its outputs all the information from the fast process $\left(y_t^f\right)$. Of course, there is no attempt in this model to interpolate the missing (nonobserved) values of the slow process.

The original system and the blocked system have the same state-space dimension. It is natural then to ask if minimality carries over. This question is examined in Anderson et al. (2016b). If the original system is controllable, then it is straightforward to establish that there can be no nonzero vector w and constant $\lambda$ for which $w^{\top }{A}_b=\lambda w^{\top },w^{\top }{B}_b=0$, i.e., the blocked model is controllable. The argument does not work for observability because there may be modes solely observed by $C^s$ which at the high frequency are distinct, but when sampled every J time instants, merge. This can occur if for example A has two distinct eigenvalues with the same J-th power. Such a situation is nongeneric however, and one can show that if A is generic (and in particular has distinct eigenvalues whose J-th powers are also distinct), observability of the underlying system implies observability of the blocked system. Of course, if the underlying system is observable from $y^f$ alone, the blocked system is observable irrespective of the genericity of A, this following by essentially the same argument as used to study the carrying over of controllability. Furthermore, if the original system is uncontrollable, or unobservable, the same property necessarily holds for the blocked system. Evidently then, for generic A, the original system and the blocked system are either both minimal, or both nonminimal.

Going on from this, it follows that if the underlying system has $\{A,B,C,D\}$ and $\{TA{T}^{-1},TB,C{T}^{-1},D\}$ as two minimal state variable quadruples for some nonsingular T, the corresponding two blocked systems are $\{A_b,B_b,C_b,D_b\}$ and $\{T{A}_b{T}^{-1},T{B}_b,C_b{T}^{-1},D_b\}$, and if any minimal realization of the blocked system is known, it must be expressible in terms of an associated underlying system via the formula (16).

Summing up, we can state:

Theorem 5.

Consider an underlying linear system defined by $x_{t+1}=A{x}_t+B{u}_t,y_t=C{x}_t+D{u}_t$, with $t=0,1,2,\dots$ and where the vector $y_t$ and matrices $C,D$ are partitioned into blocks with $p_1$ and $p_2$ rows, as

$$y_t=\left[\begin{array}{c}{y}_t^f\\ y_t^s\end{array}\right]C=\left[\begin{array}{c}{C}^f\\ C^s\end{array}\right]D=\left[\begin{array}{c}{D}^f\\ D^s\end{array}\right]$$

The process $\left(y_t^f\right)$ is observed at all time instants, while the process $\left(y_t^s\right)$ is observed at time instants $0,J,2J,\dots$. Then the associated observed block model, with input and output given by (15), has a state-variable realization given by (16). Moreover, controllability of either the underlying system or the blocked observed system implies the same property for the other. If A is generic, the same holds true in respect of observability. Coordinate basis change and blocking commute.

We also note by way of a side remark that there is a new sort of prediction problem indirectly presenting itself here: nowcasting, which can be described roughly as prediction of the immediate but unobserved past, present or immediate future. This will apply for those time instants between those at which the slow process is measured.

Following on from the preceding theorem, it is also relevant to consider what happens when the underlying system is a canonical spectral factor generating the process $\left(y_t\right)$. Under some circumstances, the blocked system will not have this property simply because of its input, state and output dimensions, and there is more than one reason why this can happen:

• Depending on the value of J and the dimensions of $u_t,y_t$, the dimension of $U_t$ may end up larger than that of $Y_t$, i.e., the transfer function of the blocked system may be fat, as noted in Anderson et al. (2016b).
• Depending on the value of J and the dimensions of $u_t,x_t$, the dimension of $U_t$ may end up greater than that of $x_t$ (which is the state dimension of the blocked system).

Apart from this, the helpful property of zero freeness, which in many applications will be a property of the underlying system, is not guaranteed to carry over to the blocked system, although, as shown in Zamani et al. (2011), if the dimension of $y_t^f$ alone exceeds that of $u_t$, and the matrices $A,B,C,D$ are generic, then the system defined by $A_b,C_b,C_b,D_b$ is zero free.

The overarching problem is one of identifying the underlying system using the mixed frequency observations. The key is to work via the blocked system. Unsurprisingly, if the blocked system is indeed zero free, the process is much more straightforward.

In econometric applications, the condition for the zero-free property of the blocked system may well be fulfilled. However, this may not be the case, even though the underlying system is zero free. This is the main scenario treated in the next section.

7. Mixed Frequency System Identification

There are two broad thrusts to this section. First, we postulate that the underlying system has an AR description, and exploit that fact. Second, we exploit the structure of the blocked system (but at a certain point, also appeal to the zero-free structure of the underlying system).

7.1. Mixed Frequency System Identificaiton with an AR Underlying System

It is well known that the AR property is not preserved under marginalization. Therefore, even if the underlying process is AR but some of its entries are only observed every J time instants, one cannot expect that the observed process is AR. (This is roughly equivalent to observing that even if the underlying process is zero-free, the blocked process may not have this property). More to the point, even if it is AR, the associated parameters in the defining matrices of the AR equation cannot be straightforwardly related to those of the underlying process. Nevertheless, using an approach based on exploiting Yule-Walker equations, progress can be made on identification. This idea was originally observed in Chen and Zadrozny (1998) and further developed in Anderson et al. (2012a). Subsequently, it was shown in Anderson et al. (2016a) that for systems with generic parameter values, the suggested algorithm could always be executed. The key here was to ensure utilization of all observed second order moments, in contrast to Chen and Zadrozny (1998) which omitted use of some.

Suppose that the underlying process is

$$\begin{array}{cc}\hfill y^t=& \left[\begin{array}{c}{y}_t^f\\ y_t^s\end{array}\right]=\left[\begin{array}{cc}{a}_{ff}\left(1\right)& a_{fs}\left(1\right)\\ a_{sf}\left(1\right)& a_{ss}\left(1\right)\end{array}\right]\left[\begin{array}{c}{y}_{t-1}^f\\ y_{t-1}^s\end{array}\right]+\cdots +\hfill \\ \hfill & +\left[\begin{array}{cc}{a}_{ff}\left(p\right)& a_{fs}\left(p\right)\\ a_{sf}\left(p\right)& a_{ss}\left(p\right)\end{array}\right]\left[\begin{array}{c}{y}_{t-p}^f\\ y_{t-p}^s\end{array}\right]+\left[\begin{array}{c}{w}_t^f\\ w_t^s\end{array}\right]\hfill \end{array}$$

(17)

with innovations covariance

$$\mathbb{E}\left[\left[\begin{array}{c}{w}_t^f\\ w_t^s\end{array}\right]\left[{\left(w_t^f\right)}^{\top }{\left(w_t^s\right)}^{\top }\right]\right]=\mathrm{\Sigma }=\left[\begin{array}{cc}{\mathrm{\Sigma }}_{ff}& {\mathrm{\Sigma }}_{fs}\\ {\mathrm{\Sigma }}_{sf}& {\mathrm{\Sigma }}_{ss}\end{array}\right]$$

The innovations covariance may be singular. The identifiability question is whether the coefficient matrices $a_{ff}\left(i\right),a_{fs}\left(i\right)$, etc and the entries of $\mathrm{\Sigma }$ can be identified from covariance data associated with the time series $\left(y_t^f\right),t=0,1,2,\dots$ and $\left(y_t^s\right),t=0,J,2J,\dots$.

As asserted above, this problem can be tackled by generating ‘extended’ Yule-Walker equations, which are obtained by multiplying (17) on the right by the fast components of $y_{t-k}^{\top },k=1,2,\dots$ and forming expectations. The resulting equations are linear in the AR coefficients and, crucially, it can be proved that they are generically solvable, i.e., for almost all values of the coefficient and entries of $\mathrm{\Sigma }$, a unique solution exists, see Anderson et al. (2016a). (An example is given below). Once these coefficients have been found, the entries of $\mathrm{\Sigma }$ can be computed, except again on a set of parameter values of measure zero. The corresponding estimators are however not “good” and in particular not asymptotically efficient, also because observed lags of the covariances of the slow process are not exploited.

7.2. A Simple AR(1) Example

In this subsection, and following Anderson et al. (2016b), we consider what is almost the simplest possible case of a mixed frequency AR system. The underlying system is AR(1), and has two scalar outputs, one of which is observed at every second time instant. We also assume regularity of the underlying system, so the zero-lag covariance $\mathrm{\Sigma }$ of the driving white noise is nonsingular. We investigate to what extent the parameters of the process before observation loss can be recovered from the observed process statistics.

Before loss of observations, the process is

$$\left[\begin{array}{c}{y}_t^f\\ y_t^s\end{array}\right]=\left[\begin{array}{cc}{a}_{ff}& a_{fs}\\ a_{sf}& a_{ss}\end{array}\right]\left[\begin{array}{c}{y}_{t-1}^f\\ y_{t-1}^s\end{array}\right]+\left[\begin{array}{c}{w}_t^f\\ w_t^s\end{array}\right]\mathbb{E}\left[\left[\begin{array}{c}{w}_t^f\\ w_t^s\end{array}\right]\left[w_t^f{w}_t^s\right]\right]=\mathrm{\Sigma }=\left[\begin{array}{cc}{\sigma }_{ff}& {\sigma }_{fs}\\ {\sigma }_{sf}& {\sigma }_{ss}\end{array}\right]$$

(18)

Set

$$\mathcal{A}=\left[\begin{array}{cc}{a}_{ff}& a_{fs}\\ a_{sf}& a_{ss}\end{array}\right]$$

With every second value of the ‘slow’ process $\left(y_t^s\right)$ being observed, it is easy to see that

$$\left[\begin{array}{c}{y}_t^f\\ y_t^s\\ y_{t-1}^f\end{array}\right]=\left[\begin{array}{cc}{\mathcal{A}}^2& 0\\ a_{ff}{a}_{fs}& 0\end{array}\right]\left[\begin{array}{c}{y}_{t-2}^f\\ y_{t-2}^s\\ y_{t-3}^f\end{array}\right]+\left[\begin{array}{c}\mathcal{A}{w}_{t-1}+w_t\\ w_{t-1}^f\end{array}\right]=\tilde{\mathcal{A}}\left[\begin{array}{c}{y}_{t-2}^f\\ y_{t-2}^s\\ y_{t-3}^f\end{array}\right]+{\tilde{w}}_t$$

(19)

with obvious definitions of $\tilde{\mathcal{A}}$ and ${\tilde{w}}_t$. It is straightforward to check that

$$\tilde{\mathcal{A}}=\left[\begin{array}{ccc}{a}_{ff}^2+a_{fs}{a}_{sf}& a_{ff}{a}_{fs}+a_{fs}{a}_{ss}& 0\\ a_{ss}{a}_{ff}+a_{ss}{a}_{sf}& a_{sf}{a}_{fs}+a_{ss}^2& 0\\ a_{ff}& a_{fs}& 0\end{array}\right]$$

(20)

and

$$\mathbb{E}\left[\tilde{w}{\tilde{w}}^{\top }\right]=\tilde{\mathrm{\Sigma }}=\left[\begin{array}{ccc}{\sigma }_{ff}& {\sigma }_{fs}& 0\\ {\sigma }_{sf}& {\sigma }_{ss}& 0\\ 0& 0& 0\end{array}\right]+\left[\begin{array}{cc}{a}_{ff}& a_{fs}\\ a_{sf}& a_{ss}\\ 1& 0\end{array}\right]\left[\begin{array}{cc}{\sigma }_{ff}& {\sigma }_{fs}\\ {\sigma }_{sf}& {\sigma }_{ss}\end{array}\right]\left[\begin{array}{ccc}{a}_{ff}& a_{sf}& 1\\ a_{fs}& a_{ss}& 0\end{array}\right]$$

(21)

The crucial question at this point is now: given the AR process generated in (19), and an identification using the process measurements of the entries of $\tilde{\mathcal{A}}$ and $\tilde{\mathrm{\Sigma }}$ (presumably using Yule-Walker equations which will yield unique values in the regular case), can the entries of $\mathcal{A}$ and $\mathrm{\Sigma }$ then be inferred?

The high level answer is that for generic values of the parameters, they can indeed be inferred. The more detailed answer is that the values cannot be inferred if and only if all three of the following conditions (defining a semi-algebraic set of parameters of measure zero) are fulfilled (for details, see Anderson et al. (2016b), which builds on a simpler version of the problem assuming diagonal $\mathrm{\Sigma }$ treated in Anderson et al. (2012a)):

$$a_{fs}=0a_{sf}+\frac{{\sigma }_{sf}}{{\sigma }_{ff}}(a_{ss}-a_{ff})=0a_{ss}\ne 0$$

As noted in Anderson et al. (2016b), if the underlying process continues to have a lag of 1, the approach is generalizable to multivariate $y_t^f$ and $y_t^s$, except that the definition of nongeneric values of parameters becomes more involved.

Generalization to lags greater than 1 is however not possible using the framework given in this subsection; this is because the AR structure of the original system no longer yields an AR structure for the blocked system. This is unsurprising, given the nonclosure of AR systems under marginalization. Nevertheless, the general approach based on modified Yule-Walker equations remains available Anderson et al. (2016b).

7.3. Mixed Frequency System Identification Exploiting Blocked System Structure

In this subsection, we record ideas originally presented in Anderson et al. (2016b) for the case where the slow process comprises every second value in the time sequence of the relevant entries of the underlying process. However, here we describe the situation when every J-th value rather than every second value is used. While the underlying system is zero free, the blocked system is not assumed to have this property.

To begin, suppose the underlying system state variable realization $\{A,B,C,D\}$ is minimal. Then generically, the same is true of the blocked system, defined by the quadruple $\{A_b,B_b,C_b,D_b\}$. The blocked system itself may not be a canonical miniphase system, but as argued in Anderson et al. (2016b), there is a minimal state-variable realization of the canonical miniphase system (yielding the same spectrum) with the same $C_b$ and $A_b$. It follows that if we identify a minimal state-variable realization of the canonical miniphase system for the blocked system, and this includes state-update and state-to-output coupling matrices ${\tilde{A}}_b,{\tilde{C}}_b$, then for some nonsingular T, there will hold ${\tilde{A}}_b=T^{-1}{A}_{b}T,{\tilde{C}}_b=T{C}_b$. From these two matrices, one seeks matrices $\tilde{A},{\tilde{C}}^f,{\tilde{C}}^s$ such that

$${\tilde{A}}_b={\tilde{A}}^J,{\tilde{C}}_b=\left[\begin{array}{c}{\tilde{C}}^f{\tilde{A}}^{J-1}\\ {\tilde{C}}^f{\tilde{A}}^{J-2}\\ ⋮\\ {\tilde{C}}^f\tilde{A}\\ {\tilde{C}}^f\\ {\tilde{C}}^s\end{array}\right]$$

(22)

Of course, ${\tilde{C}}^f$ and ${\tilde{C}}^s$ are immediately obtained from the second equation. With A sufficiently generic that no two eigenvalues have the same J-th power, the matrix $\tilde{A}$ (which has the same property) is determined up to a finite number of possibilities corresponding to different possible choices for the J-th root of each eigenvalue of ${\tilde{A}}_b$. It may be that knowledge of ${\tilde{C}}^f\tilde{A},{\tilde{C}}^f{\tilde{A}}^2,\dots$ then suffices to achieve disambiguation. Alternatively the set can be disambiguated to identify just one of the possible choices by invoking a further appeal to genericity (still excluding a semi-algebraic set of measure zero, but one which contains the set ensuring that no two eigenvalues of A have the same J-th power). Exclusion from the larger semi-algebraic set guarantees the observability of $(A,C^f)$ or equivalently of $(\tilde{A},{\tilde{C}}^f)$, which is enough to secure disambiguation.

As noted above, we must be able to identify matrices ${\tilde{A}}_b,{\tilde{C}}_b$ which are part of a quadruple of matrices defining a state-variable realization of the canonical spectral factor of the blocked system. We do not however need to identify all matrices of a realization of the canonical spectral factor. We can in fact identify the relevant matrices by identifying a state-variable realization of the covariance sequence associated with the spectrum of $Y_t$—see the discussion following Theorem 3—and this is a relatively straightforward matter, using the algorithm earlier cited, Zeiger and McEwen (1974).

Now suppose in fact that the underlying system is a VAR process of order p, thus

$$y_t=A_1{y}_{t-1}+\dots +A_p{y}_{t-p}+{\nu }_t$$

(23)

with $\left({\nu }_t\right)$ a white noise process. Of course, the system is assumed to be stable. This system can also be described by the state-variable equations

$$\begin{array}{cc}\hfill \underset{x_{t+1}}{\underbrace{\left[\begin{array}{c}{y}_t\\ ⋮\\ y_{t-p+1}\end{array}\right]}}& =\underset{A}{\underbrace{\left[\begin{array}{cccc}{A}_1& \dots & A_{p-1}& A_p\\ I& \dots & 0& 0\\ ⋮& \ddots & ⋮& ⋮\\ 0& \dots & I& 0\end{array}\right]}}\underset{x_t}{\underbrace{\left[\begin{array}{c}{y}_{t-1}\\ ⋮\\ y_{t-p}\end{array}\right]}}+\underset{B}{\underbrace{\left[\begin{array}{c}b\\ 0\\ ⋮\\ 0\end{array}\right]}}{ϵ}_t\hfill \\ \hfill y_t& =\underset{C}{\underbrace{[A_1\dots A_p]}}{x}_t+\underset{D}{\underbrace{b}}{ϵ}_t\hfill \end{array}$$

(24)

It is explained in Anderson et al. (2016a) how the coefficient matrices $A_i$ can be determined from the observable pair $(\tilde{A},\tilde{C})$ (where $\tilde{C}={\left[{\left({\tilde{C}}^f\right)}^{\top }{\left({\tilde{C}}^s\right)}^{\top }\right]}^{\top }$) under the additional genericity assumptions that $A_p$ is nonsingular, as is the covariance of the white noise process ${\nu }_t$.

In the earlier part of this paper, we stressed the potential numerical hazards when working with AR descriptions of high dimensional processes. The same warning of course applies here if the processes are typical econometric processes. This means that procedures for reducing the problem dimensionality by exploiting the singularity of the zero lag covariance of the output are desirable. How best to do this when measured signals have two periodicities (whether or not at the blocked system level) has not however been addressed in any depth in the literature. Note though that the minimal static factors can be constructed as easily for the mixed frequency case as for the single frequency case, since they are effectively determined using the zero lag covariance of the measured outputs. Another comparatively untreated issue for mixed frequency systems is that of passing to the common components process from a measurement process in which idiosyncratic components are additively embedded.

7.4. VARMA and VMA System Identification

The content of this subsection is not especially relevant to Dynamic Factor Models, but is related to the earlier material of this and the previous section, in that it focusses on further questions of identifiability when not all output values of an underlying process are observed.

An alternative starting point for the underlying process is to assume a different model class, where one postulates a model

$$y_t+A_1{y}_{t-1}+\cdots +A_p{y}_{t-p}=B_0{w}_t+B_1{w}_{t-1}+\cdots +B_s{w}_{t-s}$$

(25)

with $A_i,B_i$ of dimension $N\times N$. For convenience, make the definitions

$$a\left(z\right)=I+A_{1}z+\cdots +A_p{z}^{p}b\left(z\right)=B_0+B_{1}z+\dots B_s{z}^s$$

(26)

It is assumed that $y_t={\left[{\left(y_t^f\right)}^{\top }{\left(y_t^s\right)}^{\top }\right]}^{\top }$ with $y_t^s$ observed at time $t=0,J,2J,\dots$ for some integer $J>1$. The stability and miniphase properties are captured by the assumptions

$$deta\left(z\right)\ne 0,detb\left(z\right)\ne 0\left|z\right|\le 1$$

(27)

With the degrees of $a\left(z\right)$ and $b\left(z\right)$ assumed known, the identification task is to estimate the entries of the coefficient matrices in $a\left(z\right),b\left(z\right)$.

Of course, $\left(w_t\right)$ is a zero mean white noise sequence; we also assume $\mathbb{E}\left[w_t{w}_t^{\top }\right]=\mathrm{\Sigma }$. If rank $\mathrm{\Sigma }=q<N$, the system is singular. Otherwise it is regular.

Note that this problem setting does not have to arise in a context of high-dimensional time series (though it may), which means that non-AR modelling becomes relevant.

Identifiability of such models is linked to the question of uniqueness of the polynomials making up the transfer function $a^{-1}\left(z\right)b\left(z\right)$. In the regular case, one can require left coprimeness of the pair $a\left(z\right),b\left(z\right)$ and then enforce through multiplication of each of $a,b$ by the same unimodular (constant determinant) multiplier a canonical form. In the singular case, the situation is more complicated. For example, even for an AR system, the pair $(a,b)$ with $b\left(z\right)=B_0$ may not be left coprime. Even when they are left coprime, the Yule-Walker equations which are usually used to obtain the $A_i$ from covariances may not have a unique solution, in contrast to the regular case.

In case $s\le p$, the reference Anderson et al. (2016a) establishes that for generic values of the entries of the coefficient matrices, identifiability can be achieved. A two-step process is used. First, analogously to the procedure outlined in Section 7.1, a modified set of Yule-Walker equations can be obtained in which only known covariance data appears, and from these, the denominator matrix $a\left(z\right)$ can be obtained, given generic values of the coefficient matrices in $a\left(z\right)$ and $b\left(z\right)$. Then by using the knowledge of all the $A_i$ and some of the covariances $\mathbb{E}\left[y_t{y}_{t-s}^{\top }\right]$, it turns out that all the missing covariance data can be obtained. Following that, the matrix polynomial $b\left(z\right)$ is straightforward to obtain.

One might wonder whether the condition $s\le p$ is important. It appears so. One can in fact show, see Deistler et al. (2017),

Theorem 6.

Consider the MA model obtained from (25) by setting $A_1=A_2=\cdots =A_p=0$. Suppose that $\mathbb{E}\left[w_t{w}_t^{\top }\right]$ is nonsingular. With $y_t={\left[{\left(y_t^f\right)}^{\top }{\left(y_t^s\right)}^{\top }\right]}^{\top }]$ and $y_t^s$ observed for $t=0,J,2J\dots$ for some integer $J>1$, the parameter matrices $B_i$ are not identifiable.

The underlying reason for this conclusion is that there are more parameters to identify than there is relevant covariance data available to identify them.

A very simple example, drawn from Anderson et al. (2016a), establishes the point. In obvious notation, consider the following system where $y_t$ has dimension 2:

$$\left[\begin{array}{c}{y}_t^f\\ y_t^s\end{array}\right]=\left[\begin{array}{c}{w}_t^1\\ w_t^2\end{array}\right]+B_1\left[\begin{array}{c}{w}_{t-1}^1\\ w_{t-1}^2\end{array}\right]$$

(28)

Suppose first that all covariance data is available, viz. $\gamma \left(0\right)=\mathbb{E}\left[y_t{y}_t^{\top }\right]$ and $\gamma \left(1\right)=\mathbb{E}\left[y_t{y}_{t-1}^{\top }\right]$. It is easy to check that

$$\mathrm{\Sigma }+B_1\mathrm{\Sigma }{B}_1^{\top }=\gamma \left(0\right)B_1\mathrm{\Sigma }=\gamma \left(1\right)$$

(29)

Identification is the task of solving these equations for the seven unknown parameters defining $B_1$ and $\mathrm{\Sigma }$. The number of independent scalar equations is seven; though there are multiple solutions, they can be disambiguated by the miniphase condition requiring $\det b\left(z\right)\ne 0\mathrm{when}\left|z\right|\le 1$.

Now consider the situation where every second entry of $y_t^s$ is observed. Then the 22 entry of $\gamma \left(1\right)$ is simply unavailable, being $\mathbb{E}\left[y_t^s{y}_{t-1}^s\right]$. A continuum of values is possible for this entry, and hence a continuum of values for the entries of $\mathrm{\Sigma }$ and $B_1$.

It is perhaps surprising to find in the mixed frequency case that with an underlying AR process with lag 1 one has identifiability, while with an underlying MA process with lag 1 one does not have identifiability.

8. Conclusions

In the course of this survey of a number of results, we have highlighted several issues that are worth re-emphasising here. Furthermore, some open issues, both in the theory and applications to data are pointed out.

First, we have recorded a number of issues arising in modelling. Though there are exceptions, most problem statements with their assumptions need to guarantee that an estimated model is continuously dependent on the data from which it is obtained. Once an explicit assumption of rationality of the model transfer function has been made, continuity of the real-valued parameters for a given (and possibly estimated) set of integer-valued specification parameters with respect to the data is almost always guaranteed. Before such a rationality assumption is made however, continuity (in the $L_{\infty }$ sense) of a canonical spectral factor is not automatic, and depends on some additional assumption such as boundedness of the derivative of the spectrum. It is also useful to identify which particular form of model of an entity such as a canonical spectral factor is most preferable. One can think of ARMA/AR, state-variable, or even a frequency-based description involving plots of the different entries of the spectral factor against frequency. Moreover, the proposed application for a model (and the calculations using in addressing that application) may be a determinant of the best type of model. For some applications, e.g., those involving business cycles, state-variable models exhibiting pole positions in diagonal blocks of the state update matrix, or frequency domain plots, may be clearly preferable to AR models. For forecasting, AR models are likely to be preferable. Let us recall too that while some of the assumptions used in ensuring solvability of the identification (or modelling) problem are grounded in intuitively reasonable hypotheses, others are of a more instrumentalist character, thus invoked to ensure mathematical solvability rather than to reflect some underlying natural truth.

As a second theme, we have emphasised through much of the paper both the occurrence of and benefit from the zero-free property of canonical spectral factors arising in dynamic factor modelling. It is of great utility, and deserves to be accorded significant recognition for this fact. This observation is also relevant in the mixed frequency case. Another, perhaps less important, idea flowing from linear systems considerations is the fact that state-variable realizations of a covariance sequence and a canonical spectral factor for the associated power spectrum can be assumed to have the state transition matrix and state to output coupling matrices identical. This underpins an ability to compute from a covariance sequence a canonical spectral factor using quite different calculations to those of a Yule-Walker procedure (and thereby possibly avoid numerical problems which can arise).

A third major theme is that of mixed frequency identification problems. Several key conclusions emerge here. First, it is evident that their treatment is not as mature as the treatment of single frequency problems. Second, issues of genericity become more important. Third, certain mixed frequency identification scenarios are provably nonidentifiable (as indicated in our discussion of VARMA and MA modelling).

What now of the future? Among the outstanding issues worthy of investigation, we would certainly include the following:

• The zero-free property has been applied to macroeconomic analysis as a means to solve the fundamentalness problem. Further work is ongoing on this and, we hope, the zero-free property will increasingly become familiar to macroeconomists, see the literature cited below Theorem 4.
• An important observation on the zero-free property in statements like Theorem 4 is that its genericity is usually obtained with respect to a parameterization in which each parameter is free to vary independently of all the others, see the comment below Theorem 2. However, in applications to macroeconomic models, structural restrictions may prevent such free variability and can therefore interfere with the genericity of zero-freeness. For an attempt to solve this issue, see Forni et al. (2020).
• Another observation on Theorem 4 is that an AR representation may exist also for non-generic parameters, for which zerolessness does not hold but the zeros of $B\left(z\right)$ lie outside the unit circle. However, in that case the AR would not be of finite length.
• In some cases, we have outlined algorithms which have not yet been tried on real data. It would be worthwhile to investigate them using real data.
• A full treatment (algorithms through to testing on real data) including the decomposition into common components and idiosyncratic processes is needed for systems with multiple frequencies.