# Testing for Causal Influence using a Partial Coherence Statistic

Louis Scharf1, , and Yuan Wang2 This work was supported in part by the
National Science Foundation under contract CCF-1712788. Corresponding
author:Y. Wang (email: yuan.wang@wsu.edu). 1Departments of Mathematics and
Statistics, Colorado State University, Fort Collins, CO, USA  2Department of
Mathematics and Statistics, Washington State University, Pullman, WA, USA

###### Abstract

In this paper we explore partial coherence as a tool for evaluating causal
influence of one signal sequence on another. In some cases the signal sequence
is sampled from a time- or space-series. The key idea is to establish a
connection between questions of causality and questions of partial coherence.
Once this connection is established, then a scale-invariant partial coherence
statistic is used to resolve the question of causality. This coherence
statistic is shown to be a likelihood ratio, and its null distribution is
shown to be a Wilks Lambda. It may be computed from a composite covariance
matrix or from its inverse, the information matrix. Numerical experiments
demonstrate the application of partial coherence to the resolution of
causality. Importantly, the method is model-free, depending on no generative
model for causality.

###### Index Terms:

signal, time series, causality, partial coherence, information matrix,
likelihood ratio test, Wilks Lambda, ROC curve

##  I Introduction

Questions of causality are fraught with ambiguities. In fact, we might
paraphrase Norbert Wiener’s discussion of Newtonian physics to say causality
is a concept grounded in a physical picture that can never be fully justified
or fully rejected experimentally [1]. For this reason it is important to state
at the outset that any conclusions are conditioned on available prima facie
evidence and on the method of inference. The emergence of new evidence may
render conclusions moot or false. And, more sophisticated methods of inference
may reveal dependencies not revealed by other methods. For this reason it is
prudent to answer questions of causality in the negative, rather than the
affirmative. In other words, statements like $x$ is not a cause for $y$ are
preferred over statements like $x$ is a cause for $y$. Moreover, a conclusion
based on data is perforce a conclusion based on an experiment and a statistic
whose value is used to resolve the question of causality. This means the
question can only be resolved with respect to prima facie evidence, the
outcome of an experiment, and the value of a statistic. It follows that the
question is only resolved at a prescribed significance level. Change the
significance level, and change the answer to the question of causality. It
might be said that the question of causality is practically a question of null
hypothesis testing, where the null hypothesis states that one signal sequence
has no causal influence on another. One may still speak of the power of a
method to reject the null for a putative causal model, but this does not
validate the causal model, nor does it rule out the possibility that other
methods and evidence would negate the finding of causality.

Granger causality [2, 3] has been widely used for addressing questions of
causality. Loosely speaking, a time series $\\{y_{n}\\}$ is causally dependent
on another time series $\\{x_{n}\\}$ if the past of $\\{x_{n}\\}$ contains
unique information about $\\{y_{n}\\}$. This statement is made precise in
terms of joint probability statements concerning the predictive value of a
time series when composed with prima facie evidence. Geweke [5] has studied
the linear dependence and feedback between two time series and in [6] extended
his results to conditional linear dependence in the presence of a third time
series. Solo [7] analyses the Granger-Geweke causality measure for a variety
of state-space models, and establishes the invariance of this measure to
linear time-invariant filtering. Siggiridou and Kugiumtzis [8] used a
conditional Granger causality index to quantify Granger causality in a vector
autoregression. Kramer [10] and Amblard and Michel [11] define and analyze
directed information in their study of Granger causality graphs. Schamberg and
Coleman [9] recently proposed a novel sample path dependent measure of causal
influence between time series.

Barnett et al. [12, 13, 14] have adapted Granger’s point of view to the study
of causality in state-space models. They derive an estimate of the Granger-
Geweke causality measure for testing causality within state-space models and
demonstrate improved statistical power and reduced bias compared to results
obtained by fitting an autoregressive (AR) model to what should be an
autoregressive-moving average (ARMA) model.

In this paper, we are interested in assessing causal influence between two
signal sequences, using partial correlation or partial coherence. It seems to
us that causal influence is more descriptive than causality for describing the
aims of statistical methods that are applied to signal sequences without any
pretense at identifying or verifying a physical mechanism. As a practical
matter of testing, these signal sequences are often finite-dimensional samples
of two time series, represented as the random vectors ${\bf x}$ and ${\bf y}$.
The partial correlation between ${\bf x}$ and ${\bf y}$, given ${\bf z}$, is
the correlation between the residual for ${\bf x}$ regressed on ${\bf z}$ and
the residual for ${\bf y}$ regressed on ${\bf z}$, where ${\bf z}$ is a third
set of random variables. A zero partial correlation indicates that the random
vector ${\bf y}$ does not participate in a linear minimum mean-squared error
estimator of ${\bf x}$ from $({\bf y},{\bf z})$. In the multivariate normal
case, ${\bf x}$ is then conditionally independent of ${\bf y}$, given ${\bf
z}$. We shall define the prima facie evidence ${\bf z}$ in such a way that the
question of causality can be resolved from an analysis of partial coherence.

In a comprehensive treatment of conditioning and partial correlation for
vector-valued random variables, with no reference to the question of
causality, we identify a quantity we call partial coherence. Partial coherence
is a ratio of determinants that is proportional to the volume of a normalized
error concentration ellipse for linearly estimating a pair of random variables
from a third, or equivalently for linearly estimating one random variable from
two others. It is determined by partial canonical correlations. Partial
coherence is invariant to block diagonal non-singular transformations. In the
case of wide-sense stationary time series, the spectral representation of
partial coherence shows it to be a Riemann sum of narrowband coherences or
information rates, in close analogy with Geweke’s spectral representations.
With appropriate choice of conditioning, partial coherence is a monotone
function of the Granger-Geweke causality measure, which is the same as
Geweke’s transfer entropy.

We then turn to the question of causality, following the language and
philosophy of Granger, and studying problems of the same general class as the
Barnetts and Seth. However, in our study of causality between signal
sequences. we adopt a conservative, model-free, characterization of the
signals to be analyzed, based only on second-order moments. That is, in
contrast to the work of the Barnetts and Seth, where parametric state space
models are assumed for the time series under study, we assume only that
samples of a time series are second-order random variables. We show that the
sampled-data estimate of partial coherence is in fact the model-free ordinary
likelihood ratio for testing the null hypothesis that one signal sequence is
not correlated with another. The statistic of [12] is an approximate
likelihood ratio, based on an approximate ML identification of an underlying
state space model for the signals. That is, our estimate is ordinary
likelihood in a model-free analysis of causality, and theirs is approximate
likelihood in a state-space analysis of causality. Importantly, under the null
hypothesis, ordinary likelihood, or estimated partial coherence, is shown to
be distributed as a Wilks Lambda statistic, allowing for the setting of
thresholds to control significance level in hypothesis testing.

In the examples we analyze, the signal sequences of interest are finite
samples from wide-sense stationary time series, but this is not essential to
our methods. In fact the methodology we advocate extends to the study of
space-time series, as in [15], where the methods of this paper may be used to
test for particularly influential time series in a spatially-distributed set
of time series. Even more generally, partial coherence may be used to
characterize causality between quite arbitrary measurement vectors, as in the
labeling of graph edges by Whittaker [16].

When parametric models are accurate, and the ratio of sample count to
parameter count is large, then one would expect parametric methods to
outperform model-free methods. But when the parametric model is inaccurate,
model bias is introduced. In our re-running of the Barnetts and Seth example,
we find that the statistical power we achieve is generally comparable to
theirs, without using a parametric state-space model. We use a receiver
operating characteristic (ROC) to emphasize the point that resolution of
causality is always a problem in hypothesis testing, where the power to
identify (provisional) causality is dependent upon the significance level of
the test, or equivalently upon the probability of incorrectly rejecting the
hypothesis of non-causality.

Notations: Bold upper case letters denote matrices, boldface lower case
letters denote column vectors, and italics denote scalars. The superscript
$(\cdot)^{T}$ denotes transpose. ${\bf I}_{M}$ is the identity matrix of
dimensions $M\times M$, and ${\bf 0}$ denotes either a zero column vector or
the zero matrix of appropriate dimensions. We use ${\bf A}^{1/2}$ (${\bf
A}^{-1/2}$) to denote the square root matrix of the Hermitian matrix ${\bf A}$
(${\bf A}^{-1}$); ${\rm blkdiag}[{\bf A}_{1},\ldots,{\bf A}_{k}]$ is a block-
diagonal matrix whose diagonal blocks are ${\bf A}_{1},\ldots,{\bf A}_{k}$.
The notation ${\bf x}\sim\mathcal{N}_{M}({\bf 0},{\bf R})$ indicates that
${\bf x}$ is an $M$-dimensional Gaussian random vector of mean ${\bf 0}$ and
covariance ${\bf R}$. We say it is multivariate normal (MVN).

The rest of the paper is organized as follows. In Section II we review a few
motivating notions of predictive causal influence. In Section III, we
introduce partial coherence and demonstrate its usefulness in the estimation
of one random vector from two others. We review how partial coherence may be
extracted from the information matrix, and establish spectral formulas which
bring insight into the way partial canonical correlations additively decompose
partial coherence. We establish, in the case of multivariate normal
experiments, the connection between partial coherence, conditional Kullback-
Leibler divergence, conditional mutual information, conditional information
rate, and a Granger-Geweke causality measure. In Section IV we show that
partial coherence is a scale-invariant likelihood ratio statistic in the case
of multivariate normal experiments. Importantly, we establish that the null
distribution of the coherence statistic is the distribution of a Wilks Lambda
statistic, and show that its stochastic representation is that of a product of
independent beta-distributed random variables. In Section V, we conduct two
experiments to explore the past-causal, future-causal, and mixed-causal
relationship between time series. The first of these experiments demonstrates
the use of pairwise partial coherence as a map of causality as a function of
two time variables. The second reproduces the experiment of Barnett and Seth
[12] in order to compare our model-free estimator of causality to their model-
based estimator. We compare the power of our estimator of partial coherence to
the power of their estimator, and demonstrate the use of the Receiver
Operating Characteristic as a way to view the trade-off between power
(probability of detection, or one minus type II error probability) versus size
(significance level, probability of false alarm, or probability of type I
error). In Section VI we offer concluding remarks on the efficacy of model-
based versus model-free approaches to null hypothesis testing for causal
influence.

##  II Notions of Predictive Causal Influence

Begin with two scalar-valued time series, whose finite pasts may be denoted
$x^{t}=(x_{0},x_{1},\ldots,x_{t})$ and $y^{t}=(y_{0},y_{1},\ldots,y_{t})$.
Assume the existence of a joint probability density function (pdf)
$p(x^{t},y^{t})$, for all $t$, with corresponding marginal pdfs for any subset
of these random variables. A comparison between the conditional pdfs
$p(y_{t+1}|x^{t},y^{t})$ and $p(y_{t+1}|y^{t})$ may be used to test the null
hypothesis ${\makebox{H}}_{0}$ that $x^{t}$ has no causal influence on
$y_{t+1}$. Such a comparison might be based on the error variance
$V_{y_{t+1}|x^{t},y^{t}}$ of a predictor of $y_{t+1}$ from $(x^{t},y^{t})$
versus the error variance $V_{y_{t+1}|y^{t}}$ of a predictor of $y_{t+1}$ from
$y^{t}$ only. There are other possibilities. Conditioned on $y^{t}$, are
$y_{t+1}$ and $x_{s}$ correlated? When this question is asked for $s=t\pm r$,
then hypotheses of past-causal, future-causal, mixed-causal influence may be
answered. These questions are easily generalized to multivariate time series
by comparing the error covariance for predictors of a random vector ${\bf
y}_{t+1}$ from random vectors $({\bf x}^{t},{\bf y}^{t})$.

When prediction error variances are replaced by entropies, then entropy
differences may be used to measure mutual information. Define the following
entropies $H$ and mutual information $I$:

$\displaystyle H_{y_{t+s}}={\mathbb{E}}[-\log p(y_{t+s})],\,\makebox{entropy
of the random variable }y_{t+s},$ $\displaystyle
H_{y_{t+s}|y^{t}}={\mathbb{E}}[-\log p(y_{t+s}|y^{t})],$ $\displaystyle
H_{y_{t+s}|x^{t},y^{t}}={\mathbb{E}}[-\log p(y_{t+s}|x^{t},y^{t}],$
$\displaystyle
I_{y_{t+s},x^{t}|y^{t}}=H_{y_{t+s}|y^{t}}-H_{y_{t+s}|x^{t},y^{t}}.$

In these expressions, the random variable $p(y_{t+s})$ is a pdf random
variable, and ${\mathbb{E}}$ is expectation. The entropy $H_{y_{t+s}}$ is a
measure of uncertainty for the random variable. The other entropies are
entropies conditioned on the past $y^{t}$ or the two pasts $(x^{t},y^{t})$.
The mutual information $I_{y_{t+s},x^{t}|y^{t}}$ measures mutual information
between $y_{t+s}$ and $x^{t}$, after conditioning on the past $y^{t}$. A
special case is $I_{y_{t+s},x_{s}|y^{t}}$, which measures the mutual
information between $y_{t+s}$ and $x_{s}$, after conditioning on the past
$y^{t}$. By construction the conditional entropies and the conditional mutual
information are measures of directed information. For ergodic time series,
they may be averaged over $t$ to determine the minimum codelength for encoding
the time series $(y_{0},y_{1},\ldots)$, with or without the side information
$x^{t}$. Rissanen and Wax [4] are as cautionary as Wiener in their use of the
term causality when they say, “We measure causality by predictability, with
predictability interpreted in terms of the codelength required to encode the
time series with predictive coding. … Obviously, the value of the resulting
measure depends critically on the class of predictive densities selected. A
good selection gives a sharp measure of causal dependence while a bad one
masks a possible causal dependence, which of course is just as it should be.”

In the framing of questions about past-causal, future-causal, or mixed causal-
influence, three random vectors $({\bf x},{\bf y},{\bf z})$ may be extracted
from two signal sequences. For example, ${\bf x}=x^{t},{\bf y}=y_{t+1},{\bf
z}=y^{t}$. Or ${\bf x}=x_{s},{\bf y}=y_{t+1},{\bf z}=y^{t}$. More generally,
${\bf x}$ is a $p$-sample from time series $\\{x_{n}\\}$, ${\bf y}$ is a
$q$-sample from the time series $\\{y_{n}\\}$, and ${\bf z}$ is an $r$-sample
from either $\\{x_{n}\\}$ or $\\{y_{n}\\}$. That is, three finite-dimensional
vectors are extracted from two signal sequences. The construction determines
the nature of the question of causal influence. This suggests a general
framework for analyzing causal influence based on three random vectors ${\bf
x}\in{\mathbb{R}}^{p},{\bf y}\in{\mathbb{R}}^{q},{\bf z}\in{\mathbb{R}}^{r}$.
Then partial coherence may be used as a tool for testing hypotheses of causal
influence. It encodes for conditional prediction error covariance. For zero-
mean Gaussian time series, in which case all finite-dimensional random vectors
are multivariate normal (MVN), partial coherence encodes for conditional
entropies, and mutual information. These points will be further clarified in
due course.

We do not impose a finite-dimensional time series model, such as AR, ARX,
ARMA, or ARMAX for signal sequences. Rather we assume only that the signal
sequences are finite-variance with arbitrary non-negative definite covariances
for every finite windows-worth of the sequences. These are then estimated from
multiple realizations of an experiment, or in the case of ergodic times
series, from consecutive windows of the signal sequences. There is no
parameter estimation for coefficients in a time series model that is based on
a finite-dimensional difference equation for covariances.

##  III  Predictive Causality, Linear Regression, and Partial Coherence

There are two ways to frame questions of causal influence: as a question of
estimating two random vectors from one other, which is to say regressing two
random vectors onto one random vector; or as a question of estimating one
random vector from two others, which is to say regressing one random vector
onto two random vectors. In the first case the random vectors $({\bf x},{\bf
y})$ are to be regressed onto the random vector ${\bf z}$; in the second case,
the random vector ${\bf x}\in{\mathbb{R}}^{p}$ is to be regressed onto the
random vectors ${\bf y}\in{\mathbb{R}}^{q}$ and ${\bf z}\in{\mathbb{R}}^{r}$.
The composite covariance matrix for these three random vectors is

${\bf R}={\mathbb{E}}[\left[\begin{array}[]{ccc}{\bf x}\\\ {\bf y}\\\ {\bf
z}\end{array}\right]\left[\begin{array}[]{ccc}{\bf x}^{T}&{\bf y}^{T}&{\bf
z}^{T}\end{array}\right]]=\left[\begin{array}[]{ccccccccc}{\bf{\bf R}_{{\bf
x}{\bf x}}}&{\bf R}_{{\bf x}{\bf y}}&{\bf R}_{{\bf x}{\bf z}}\\\ {\bf{\bf
R}_{{\bf y}{\bf x}}}&{\bf R}_{{\bf y}{\bf y}}&{\bf R}_{{\bf y}{\bf z}}\\\
{\bf{\bf R}_{{\bf z}{\bf x}}}&{\bf R}_{{\bf z}{\bf y}}&{\bf R}_{{\bf z}{\bf
z}}\end{array}\right].$ (1)

By defining the composite vectors ${\bf u}=({\bf x}^{T},{\bf y}^{T})^{T}$, and
${\bf v}=({\bf y}^{T},{\bf z}^{T})^{T}$ the covariance matrix ${\bf R}$ may be
parsed two ways:

${\bf R}=\left[\begin{array}[]{cccc}{\bf R}_{{\bf u}{\bf u}}&{\bf R}_{{\bf
u}{\bf z}}\\\ {\bf R}_{{\bf u}{\bf z}}^{T}&{\bf R}_{{\bf z}{\bf
z}}\end{array}\right]=\left[\begin{array}[]{cccc}{\bf R}_{{\bf x}{\bf x}}&{\bf
R}_{{\bf x}{\bf v}}\\\ {\bf R}_{{\bf x}{\bf v}}^{T}&{\bf R}_{{\bf v}{\bf
v}}\end{array}\right].$ (2)

The covariance matrix ${\bf R}_{{\bf u}{\bf u}}$ is $(p+q)\times(p+q)$, and
the covariance matrix ${\bf R}_{{\bf x}{\bf x}}$ is $p\times p$. There are two
useful representations for the inverse of the composite covariance matrix
${\bf R}$:

$\displaystyle{\bf R}^{-1}=$ $\displaystyle\left[\begin{array}[]{cccc}{\bf
R}_{{\bf u}{\bf u}|{\bf z}}^{-1}&{\bf R}_{{\bf u}{\bf u}|{\bf z}}^{-1}{\bf
R}_{{\bf u}{\bf z}}\\\ {\bf R}_{{\bf u}{\bf z}}^{T}{\bf R}_{{\bf u}{\bf
u}|{\bf z}}^{-1}&{\bf R}_{{\bf z}{\bf z}}^{-1}+{\bf R}_{{\bf z}{\bf
z}}^{-1}{\bf R}_{{\bf u}{\bf z}}^{T}{\bf R}_{{\bf u}{\bf u}|{\bf z}}^{-1}{\bf
R}_{{\bf u}{\bf z}}{\bf R}_{{\bf z}{\bf z}}^{-1}\end{array}\right]$ (5)
$\displaystyle=$ $\displaystyle\left[\begin{array}[]{cccc}{\bf R}_{{\bf x}{\bf
x}|{\bf v}}^{-1}&{\bf R}_{{\bf x}{\bf x}|{\bf v}}^{-1}{\bf R}_{{\bf x}{\bf
v}}\\\ {\bf R}_{{\bf x}{\bf v}}^{T}{\bf R}_{{\bf x}{\bf x}|{\bf v}}^{-1}&{\bf
R}_{{\bf v}{\bf v}}^{-1}+{\bf R}_{{\bf v}{\bf v}}^{-1}{\bf R}_{{\bf x}{\bf
v}}^{T}{\bf R}_{{\bf x}{\bf x}|{\bf v}}^{-1}{\bf R}_{{\bf x}{\bf v}}{\bf
R}_{{\bf v}{\bf v}}^{-1}\end{array}\right].$ (8)

The matrix ${\bf R}_{{\bf u}{\bf u}|{\bf z}}$ is the error covariance matrix
for estimating the composite vector ${\bf u}$ from ${\bf z}$, and the matrix
${\bf R}_{{\bf x}{\bf x}|{\bf v}}$ is the error covariance matrix for
estimating ${\bf x}$ from ${\bf v}$:

$\displaystyle{\bf R}_{{\bf u}{\bf u}|{\bf z}}={\bf R}_{{\bf u}{\bf u}}-{\bf
R}_{{\bf u}{\bf z}}{\bf R}_{{\bf z}{\bf z}}^{-1}{\bf R}_{{\bf u}{\bf z}}^{T},$
(9) $\displaystyle{\bf R}_{{\bf x}{\bf x}|{\bf v}}={\bf R}_{{\bf x}{\bf
x}}-{\bf R}_{{\bf x}{\bf v}}{\bf R}_{{\bf v}{\bf v}}^{-1}{\bf R}_{{\bf x}{\bf
v}}^{T}.$ (10)

We shall have more to say about these error covariance matrices in due course.
Importantly, the inverses of each may be read out of the inverse for the
composite covariance matrix ${\bf R}^{-1}$ of eqn(8). The dimension of the
error covariance ${\bf R}_{{\bf u}{\bf u}|{\bf z}}$ is $(p+q)\times(p+q)$ and
the dimension of the error covariance ${\bf R}_{{\bf x}{\bf x}|{\bf v}}$ is
$p\times p$.

###  III-A Regressing two random vectors onto one: partial correlation and
partial canonical correlation coefficients

The estimators of ${\bf x}$ and ${\bf y}$ from ${\bf z}$, and their resulting
error covariance matrices are easily read out from the composite covariance
matrix ${\bf R}$:

$\displaystyle\hat{{\bf x}}({\bf z})=$ $\displaystyle{\bf R}_{{\bf x}{\bf
z}}{\bf R}_{{\bf z}{\bf z}}^{-1}{\bf z},$ (11) $\displaystyle{\bf R}_{{\bf
x}{\bf x}|{\bf z}}=$ $\displaystyle{\mathbb{E}}[({\bf x}-\hat{{\bf x}})({\bf
x}-\hat{{\bf x}})^{T}]={\bf R}_{{\bf x}{\bf x}}-{\bf R}_{{\bf x}{\bf z}}{\bf
R}_{{\bf z}{\bf z}}^{-1}{\bf R}_{{\bf x}{\bf z}}^{T},$ (12)
$\displaystyle\hat{{\bf y}}({\bf z})=$ $\displaystyle{\bf R}_{{\bf y}{\bf
z}}{\bf R}_{{\bf z}{\bf z}}^{-1}{\bf z},$ (13) $\displaystyle{\bf R}_{{\bf
y}{\bf y}|{\bf z}}=$ $\displaystyle{\mathbb{E}}[({\bf y}-\hat{{\bf y}})({\bf
y}-\hat{{\bf y}})^{T}]={\bf R}_{{\bf y}{\bf y}}-{\bf R}_{{\bf y}{\bf z}}{\bf
R}_{{\bf z}{\bf z}}^{-1}{\bf R}_{{\bf y}{\bf z}}^{T}.$ (14)

The composite error covariance matrix for the errors ${\bf x}-\hat{{\bf
x}}({\bf z})$ and ${\bf y}-\hat{{\bf y}}({\bf z})$ is the matrix

$\displaystyle{\bf R}_{{\bf u}{\bf u}|{\bf z}}=$
$\displaystyle{\mathbb{E}}\left[\begin{array}[]{cc}{\bf x}-\hat{{\bf x}}({\bf
z})\\\ {\bf y}-\hat{{\bf y}}({\bf
z})\end{array}\right]\left[\begin{array}[]{cc}({\bf x}-\hat{{\bf x}}({\bf
z}))^{T}&({\bf y}-\hat{{\bf y}}({\bf z}))^{T}\end{array}\right]$ (18)
$\displaystyle=$ $\displaystyle\left[\begin{array}[]{cccc}{\bf R}_{{\bf x}{\bf
x}|{\bf z}}&{\bf R}_{{\bf x}{\bf y}|{\bf z}}\\\ {\bf R}_{{\bf x}{\bf y}|{\bf
z}}^{T}&{\bf R}_{{\bf y}{\bf y}|{\bf z}}\end{array}\right],$ (21)
$\displaystyle{\bf R}_{{\bf x}{\bf y}|{\bf z}}=$
$\displaystyle{\mathbb{E}}[({\bf x}-\hat{{\bf x}}({\bf z}))({\bf y}-\hat{{\bf
y}}({\bf z}))^{T}]$ $\displaystyle=$ $\displaystyle{\bf R}_{{\bf x}{\bf
y}}-{\bf R}_{{\bf x}{\bf z}}{\bf R}_{{\bf z}{\bf z}}^{-1}{\bf R}_{{\bf y}{\bf
z}}^{T}.$ (22)

In this matrix the matrix block ${\bf R}_{{\bf x}{\bf y}|{\bf z}}$ is the
$p\times q$ matrix of cross correlations between the random errors $({\bf
x}-\hat{{\bf x}}({\bf z}))$ and $({\bf y}-\hat{{\bf y}}({\bf z}))$. It is
called the partial correlation between the random vectors ${\bf x}$ and ${\bf
y}$, after each has been regressed onto the common vector ${\bf z}$. Under the
null hypothesis $\makebox{H}_{0}$, this partial correlation is zero, and
therefore ${\bf R}_{{\bf u}{\bf u}|{\bf z}}$ is the block-diagonal matrix
$\text{blkdiag}[{\bf R}_{{\bf x}{\bf x}|{\bf z}},{\bf R}_{{\bf y}{\bf y}|{\bf
z}}]$.

Formula (8) shows that ${\bf R}_{{\bf u}{\bf u}|{\bf z}}$ and its determinant
may be read directly out of the $(p+q)\times(p+q)$ Northwest block of the
inverse of the error covariance matrix:

$\displaystyle({\bf R}^{-1})_{NW}$ $\displaystyle={\bf R}_{{\bf u}{\bf u}|{\bf
z}}^{-1}\hskip 7.22743pt\makebox{and}$ $\displaystyle\det[{\bf R}_{{\bf u}{\bf
u}|{\bf z}}]$ $\displaystyle=\frac{1}{\det[({\bf R}^{-1})_{NW}]}.$ (23)

This result was known to Harold Cramer more than 70 years ago [18], and is
featured prominently in the book on graphical models by Whittaker [16]. The
consequence is that the ingredients of a partial coherence statistic we shall
compute from experimental data may be read out of the Northwest blocks of
inverse sample covariance matrices.

The error covariance matrix ${\bf R}_{{\bf u}{\bf u}|{\bf z}}$ may be pre- and
post-multiplied by $\text{blkdiag}[{\bf R}_{{\bf x}{\bf x}|{\bf
z}}^{-1/2},{\bf R}_{{\bf y}{\bf y}|{\bf z}}^{-1/2}]$ to produce the normalized
error covariance matrix and its corresponding determinant:

$\displaystyle{\bf Q}_{{\bf u}{\bf u}|{\bf z}}=$
$\displaystyle\left[\begin{array}[]{cccc}{\bf I}&{\bf R}_{{\bf x}{\bf
x}\mid{\bf z}}^{-1/2}{\bf R}_{{\bf x}{\bf y}\mid{\bf z}}{\bf R}_{{\bf y}{\bf
y}\mid{\bf z}}^{-1/2}\\\ {\bf R}_{{\bf x}{\bf x}\mid{\bf z}}^{-1/2}{\bf
R}_{{\bf x}{\bf y}\mid{\bf z}}^{T}{\bf R}_{{\bf y}{\bf y}\mid{\bf
z}}^{-1/2}&{\bf I}\end{array}\right],$ (26) $\displaystyle\det[{\bf Q}_{{\bf
u}{\bf u}|{\bf z}}]=$ $\displaystyle\frac{\det[{{\bf R}_{{\bf u}{\bf u}|{\bf
z}}}]}{\det[{\bf R}_{{\bf x}{\bf x}|{\bf z}}]\det[{\bf R}_{{\bf y}{\bf y}|{\bf
z}}]}=\det[{\bf I}-{\bf C}_{{\bf x}{\bf y}|{\bf z}}{\bf C}^{T}_{{\bf x}{\bf
y}|{\bf z}}],$ (27) $\displaystyle{\bf C}_{{\bf x}{\bf y}|{\bf z}}=$
$\displaystyle{\bf R}_{{\bf x}{\bf x}\mid{\bf z}}^{-1/2}{\bf R}_{{\bf x}{\bf
y}\mid{\bf z}}{\bf R}_{{\bf y}{\bf y}\mid{\bf z}}^{-1/2}.$ (28)

The NE matrix block of the covariance matrix ${{\bf Q}_{{\bf u}{\bf u}|{\bf
z}}}$ is the partial coherence matrix ${\bf C}_{{\bf x}{\bf y}|{\bf z}}$. The
determinant of ${\bf Q}_{{\bf u}{\bf u}|{\bf z}}$ is proportional to the
volume of the normalized error concentration ellipse for estimating ${\bf u}$
from ${\bf z}$. When specialized for causal influence, this gives an
interpretation of the Granger-Geweke causality measure that is consistent
with, but slightly different than, the discussion given in [12]. Under the
null hypothesis $\makebox{H}_{0}$, the normalized error covariance matrix is
identity, ${\bf Q}_{{\bf u}{\bf u}|{\bf z}}={\bf I}_{p+q}$, its determinant is
one, and the partial coherence matrix is the zero matrix, ${\bf C}_{{\bf
x}{\bf y}|{\bf z}}={\bf 0}_{p\times q}$.

We define partial coherence to be

$\displaystyle\rho^{2}_{{\bf x}{\bf y}|{\bf z}}=$ $\displaystyle 1-\det[{{\bf
Q}_{{\bf u}{\bf u}|{\bf z}}}]=1-\frac{\det[{\bf R}_{{\bf u}{\bf u}|{\bf
z}}]}{\det[{\bf R}_{{\bf x}{\bf x}|{\bf z}}]\det[{\bf R}_{{\bf y}{\bf y}|{\bf
z}}]}$ $\displaystyle=$ $\displaystyle 1-\det[{\bf I}-{\bf C}_{{\bf x}{\bf
y}|{\bf z}}{\bf C}^{T}_{{\bf x}{\bf y}|{\bf z}}].$ (29)

Define the SVD of the partial coherence matrix to be ${\bf C}_{{\bf x}{\bf
y}\mid{\bf z}}={\mathbf{F}}{\mathbf{K}}{\mathbf{G}}^{T}$, where ${\mathbf{F}}$
is a $p\times p$ orthogonal matrix, ${\mathbf{G}}$ is a $q\times q$ orthogonal
matrix, and ${\mathbf{K}}$ is a $p\times q$ diagonal matrix of partial
canonical correlations. The matrix ${\mathbf{K}}$ may be called the partial
canonical correlation matrix. The normalized error covariance matrix of eqn
(26) may be written as

$\displaystyle{\bf Q}_{{\bf u}{\bf u}|{\bf
z}}=\left[\begin{array}[]{cccc}{\mathbf{F}}&{\bf 0}\\\ {\bf
0}&{\mathbf{G}}\end{array}\right]\left[\begin{array}[]{cccc}{\bf
I}&{\mathbf{K}}\\\ {\mathbf{K}}^{T}&{\bf
I}\end{array}\right]\left[\begin{array}[]{cccc}{\mathbf{F}}^{T}&{\bf 0}\\\
{\bf 0}&{\mathbf{G}}^{T}\end{array}\right].$ (36)

As a consequence, partial coherence may be factored as

$\displaystyle\rho^{2}_{{\bf x}{\bf y}|{\bf z}}=1-\makebox{det}[{\bf
I}-{\mathbf{K}}{\mathbf{K}}^{T}]=1-\prod_{i=1}^{\makebox{min}(p,q)}(1-k_{i}^{2}).$
(37)

The partial canonical correlations $k_{i}$ are bounded between $0$ and $1$, as
is partial coherence. When the squared partial canonical correlations
$k_{i}^{2}$ are near to zero, then partial coherence $\rho^{2}_{{\bf x}{\bf
y}|{\bf z}}$ is near to zero, indicating linear independence of ${\bf x}$ and
${\bf y}$, given ${\bf z}$.

These results summarize the error analysis for linearly estimating the random
vectors ${\bf x}$ and ${\bf y}$ from a common random vector ${\bf z}$. The
only assumption is that the random vectors $({\bf x},{\bf y},{\bf z})$ are
second-order random vectors. By constructing the random vectors ${\bf x},{\bf
y},{\bf z}$ appropriately, we shall answer questions of second-order causal
influence.

###  III-B Regressing one random vector onto two: partial correlation and
partial canonical correlation coefficients

Suppose now that the random vector ${\bf x}$ is to be linearly regressed onto
${\bf v}=({\bf y},{\bf z})$:

$\hat{{\bf x}}({\bf v})=\left[\begin{array}[]{cc}{\bf R}_{{\bf x}{\bf y}}&{\bf
R}_{{\bf x}{\bf z}}\end{array}\right]\left[\begin{array}[]{cccc}{\bf R}_{{\bf
y}{\bf y}}&{\bf R}_{{\bf y}{\bf z}}\\\ {\bf R}_{{\bf y}{\bf z}}^{T}&{\bf
R}_{{\bf z}{\bf z}}\end{array}\right]^{-1}\left[\begin{array}[]{cc}{\bf y}\\\
{\bf z}\end{array}\right].$ (38)

Give the matrix inverse in this equation, the following block diagonal LDU
factorization:

$\displaystyle\left[\begin{array}[]{cccc}{\bf R}_{{\bf y}{\bf y}}&{\bf
R}_{{\bf y}{\bf z}}\\\ {\bf R}_{{\bf y}{\bf z}}^{T}&{\bf R}_{{\bf z}{\bf
z}}\end{array}\right]^{-1}=$ (41)
$\displaystyle\left[\begin{array}[]{cccc}{\bf I}&{\bf 0}\\\ -{\bf R}_{{\bf
z}{\bf z}}^{-1}{\bf R}_{{\bf y}{\bf z}}^{T}&{\bf
I}\end{array}\right]\left[\begin{array}[]{cccc}{\bf R}_{{\bf y}{\bf y}\mid{\bf
z}}^{-1}&{\bf 0}\\\ {\bf 0}&{\bf R}_{{\bf z}{\bf
z}}^{-1}\end{array}\right]\left[\begin{array}[]{cccc}{\bf I}&-{\bf R}_{{\bf
y}{\bf z}}{\bf R}_{{\bf z}{\bf z}}^{-1}\\\ {\bf 0}&{\bf I}\end{array}\right].$
(48)

A few lines of algebra produce this result for $\hat{{\bf x}}({\bf v})$, the
linear minimum mean-squared error estimator of ${\bf x}$ from ${\bf v}$:

$\displaystyle\hat{{\bf x}}({\bf v})=\hat{{\bf x}}({\bf z})+{\bf R}_{{\bf
x}{\bf y}\mid{\bf z}}{\bf R}_{{\bf y}{\bf y}\mid{\bf z}}^{-1}\left[{\bf
y}-\hat{{{\bf y}}}({\bf z})\right].$ (49)

It is evident that the vector ${\bf y}$ is not used in a linear minimum mean-
squared error estimator of ${\bf x}$ when the partial covariance ${\bf
R}_{{\bf x}{\bf y}|{\bf z}}$ is zero. That is, the random vector ${\bf y}$
brings no useful second-order information to the problem of linearly
estimating ${\bf x}$.

The error covariance matrix for estimating ${\bf x}$ from ${\bf v}$ is easily
shown to be

$\displaystyle{\bf R}_{{\bf x}{\bf x}\mid{\bf v}}={\mathbb{E}}[({\bf
x}-\hat{{\bf x}}({\bf v})({\bf x}-\hat{{\bf x}}({\bf v}))^{T}]$
$\displaystyle={\bf R}_{{\bf x}{\bf x}\mid{\bf z}}-{\bf R}_{{\bf x}{\bf
y}\mid{\bf z}}{\bf R}_{{\bf y}{\bf y}\mid{\bf z}}^{-1}{\bf R}_{{\bf x}{\bf
y}\mid{\bf z}}^{T}.$ (50)

Thus the error covariance ${{\bf R}_{{\bf x}{\bf x}\mid{\bf z}}}$ is reduced
by a quadratic form depending on the covariance between the errors $({\bf
x}-\hat{{\bf x}}({\bf z}))$ and $({\bf y}-\hat{{\bf y}}({\bf z}))$. If this
error covariance is now normalized by the error covariance matrix achieved by
regressing only on the vector ${\bf z}$, the result is

$\displaystyle{{\bf P}_{{\bf x}{\bf x}|{\bf v}}}=$ $\displaystyle{\bf R}_{{\bf
x}{\bf x}|{\bf z}}^{-1/2}{\bf R}_{{\bf x}{\bf x}|{\bf v}}{\bf R}_{{\bf x}{\bf
x}|{\bf z}}^{-1/2}$ $\displaystyle=$ $\displaystyle{\bf I}-{\bf R}_{{\bf
x}{\bf x}|{\bf z}}^{-1/2}{\bf R}_{{\bf x}{\bf y}\mid{\bf z}}{\bf R}_{{\bf
y}{\bf y}\mid{\bf z}}^{-1}{\bf R}_{{\bf x}{\bf y}\mid{\bf z}}^{T}{\bf R}_{{\bf
x}{\bf x}|{\bf z}}^{-1/2}$ $\displaystyle=$ $\displaystyle{\bf I}-{\bf
C}_{{\bf x}{\bf y}|{\bf z}}{\bf C}_{{\bf x}{\bf y}|{\bf
z}}^{T}={\mathbf{F}}({\bf I}-{\mathbf{K}}{\mathbf{K}}^{T}){\mathbf{F}}^{T}.$
(51)

The determinant of this matrix measures the volume of the normalized error
covariance matrix:

$\displaystyle\makebox{det}[{\bf P}_{{\bf x}{\bf x}|{\bf v}}]$
$\displaystyle=\frac{\makebox{det}[{\bf R}_{{\bf x}{\bf x}|{\bf
v}}]}{\makebox{det}[{\bf R}_{{\bf x}{\bf x}|{\bf z}}]}$
$\displaystyle=\makebox{det}[{\bf
I}-{\mathbf{K}}{\mathbf{K}}^{T}]=\prod_{i=1}^{\makebox{min}(p,q)}(1-k_{i}^{2}).$
(52)

As before, we may define a partial coherence

$\displaystyle\rho^{2}_{{\bf x}|{\bf y}{\bf z}}=1-\det[{\bf P}_{{\bf x}{\bf
x}|{\bf v}}]=1-\prod_{i=1}^{\makebox{min}(p,q)}(1-k_{i}^{2}).$ (53)

When the squared partial canonical correlations $k_{i}^{2}$ are near to zero,
then partial coherence $\rho^{2}_{{\bf x}|{\bf y}{\bf z}}$ is near to zero,
indicating linear independence of ${\bf x}$ on ${\bf y}$, given ${\bf z}$.
Consequently the estimator $\hat{{\bf x}}({\bf v})$ depends only on ${\bf z}$,
and not on ${\bf y}$.

These results summarize the error analysis for estimating the random vector
${\bf x}$ from the composite vector ${\bf v}$. It is notable that, except for
scaling constants dependent only upon the dimensions $p,q,r$, the volume of
the normalized error covariance matrix for estimating ${\bf x}$ from ${\bf v}$
equals the volume of the normalized error covariance matrix for estimating
${\bf u}$ from ${\bf z}$. Both of these volumes are determined by the partial
canonical correlations $k_{i}$. Importantly, for answering questions of causal
influence, based on linear prediction, it makes no difference whether one
considers $\rho^{2}_{{\bf x}{\bf y}|{\bf z}}$ or $\rho^{2}_{{\bf x}|{\bf
y}{\bf z}}$. These two measures of coherence are identical. So, for example,
to analyze $\rho^{2}_{x_{s}y_{t}|{\bf z}}$ is to analyze
$\rho^{2}_{x_{s}|y_{t}{\bf z}}$, for ${\bf x}$ defined to be $x_{s}$ and ${\bf
y}$ defined to be $y_{t}$. In Section III.D we show that in the MVN case, the
partial coherence between ${\bf x}$ and ${\bf y}$, given ${\bf z}$ is a
monotone function of conditional KL divergence, conditional mutual
information, conditional information rate, and a Granger-Geweke causality
measure.

Importantly the partial canonical correlations $k_{i}$ are invariant to
transformation of the random vector $({\bf x}^{T},{\bf y}^{T},{\bf
z}^{T})^{T}$ by a block diagonal, non-singular, matrix ${\bf
T}=\makebox{blkdiag}[{\bf T}_{{\bf x}{\bf x}},{\bf T}_{{\bf y}{\bf y}},{\bf
T}_{{\bf z}{\bf z}}]$. As a consequence partial coherence is invariant to
transformation ${\bf T}$. This result is the finite-dimensional version of
Solo’s finding [7] that the Granger-Geweke causality measure is invariant to
uncoupled linear filtering of wide-sense stationary time series. A slight
variation on Proposition 10.6 in [17] shows partial canonical correlations to
be maximal invariants under group action ${\bf T}$.

###  III-C Spectral formulas

There are spectral formulas. Begin with wide-sense stationary scalar time
series $\\{x_{n}\\},\\{y_{n}\\},\\{z_{n}\\}$, from which error time series are
computed. Assume the error covariance matrices ${\bf R}_{{\bf x}{\bf x}|{\bf
z}}$, ${\bf R}_{{\bf y}{\bf y}|{\bf z}}$, and ${\bf R}_{{\bf x}{\bf y}|{\bf
z}}$ are Toeplitz matrices constructed from the scalar time series of error
variances, $\\{R_{xx|z}[n]\\}$, $\\{R_{yy|z}[n]\\}$, and $\\{R_{xy|z}[n]\\}$.
Then the formulas of Szego [19], [20] may be applied to show that in the limit
$p,q\longrightarrow\infty$,

$\displaystyle\rho^{2}_{xy|z}$ $\displaystyle\longrightarrow
1-\exp\\{\int_{0}^{2\pi}\log\left(1-\hat{k}^{2}(e^{j\theta})\right)\frac{d\theta}{2\pi}\\},$
(54) $\displaystyle\hat{k}^{2}(e^{j\theta})=$
$\displaystyle\frac{|\hat{R}_{xy|z}(e^{j\theta})|^{2}}{\hat{R}_{xx|z}(e^{j\theta})\hat{R}_{yy|z}(e^{j\theta})}.$
(55)

where $\hat{R}_{xx|z}(e^{j\theta})$ is the discrete-time Fourier transform
(DTFT) of $\\{R_{xx|z}[n]\\}$, $\hat{R}_{yy|z}(e^{j\theta})$ is the DTFT of
$\\{R_{yy|z}[n]\\}$, and $\hat{R}_{xy|z}(e^{j\theta})$ is the DTFT of
$\\{R_{xy|z}[n]\\}$. This formula generalizes to vector-valued time series
[15]. We may say the narrowband partial canonical correlations
$\hat{k}(e^{j\theta})$ resolve the broadband partial coherence
$\rho^{2}_{xy|z}$, in direct correspondence to eqn (37). Similar formulas may
be found in [21]. The partial canonical correlation spectrum
$\\{\hat{k}(e^{j\theta}),0<\theta\leq 2\pi\\}$ is invariant to time-invariant
linear filtering of the wide-sense stationary time series
$\\{x_{n}\\},\\{y_{n}\\},\\{z_{n}\\}$.

###  III-D Kullback-Leibler Divergence, Conditional Mutual Information,
Conditional Information Rate, Partial Coherence, and Granger-Geweke Causality
Measure

Begin with the random vectors $({\bf x},{\bf y},{\bf z})$. The Kullback-
Leibler divergence between two probability laws $P_{1}$ and $P_{0}$ for these
is defined to be

$D_{KL}(P_{1}||P_{0})={\mathbb{E}}_{P_{1}}\log\frac{p_{1}({\bf x},{\bf y},{\bf
z})}{p_{0}({\bf x},{\bf y},{\bf z})},$ (56)

where $p_{i}({\bf x},{\bf y},{\bf z})$ is the probabiity density function
under probability law $P_{i}$. We may think of the random variables
$p_{1}({\bf x},{\bf y},{\bf z})$ and $p_{0}({\bf x},{\bf y},{\bf z})$ as pdf
random variables or belief state random variables. Their ratio $p_{1}({\bf
x},{\bf y},{\bf z})/p_{0}({\bf x},{\bf y},{\bf z})$ is a likelihood ratio
random variable. So KL divergence may be considered the expected value of the
log-likelihood ratio for testing the model $P_{1}$ against the model $P_{0}$,
under the hypothesis that measurements are drawn from the probability law
$P_{1}$.

When $({\bf x},{\bf y},{\bf z})$ is composed as $({\bf u},{\bf z})$, then
$p({\bf x},{\bf y},{\bf z})=p({\bf u}|{\bf z})p({\bf z})$ and the ratio is
$p_{1}({\bf u}|{\bf z})/p_{0}({\bf u}|{\bf z})$, provided the marginal
distribution of ${\bf z}$ is the same under probability laws $P_{1}$ and
$P_{0}$. In the case where $({\bf u},{\bf z})$ is multivariate normal (MVN)
with mean ${\bf 0}$ and the covariance matrix ${\bf R}$ of eqn (2), it is not
hard to show that the Kullback-Leibler divergence is

$D_{KL}(P_{1}||P_{0})=-\frac{1}{2}\log\det[{\bf Q}_{{\bf u}{\bf u}|{\bf z}}].$
(57)

That is, after regression of ${\bf x}$ and ${\bf y}$ onto ${\bf z}$, the
problem is reset to the analysis of the errors $({\bf x}-\hat{{\bf x}}({\bf
z}))$ and $({\bf y}-\hat{{\bf y}}({\bf z}))$, and their partial coherence. In
the MVN case, KL divergence for this problem is also the conditional mutual
information, $I_{{\bf x}{\bf y}|{\bf z}}={\mathbb{E}}\left[\log\frac{p({\bf
x},{\bf y}|{\bf z})}{p({\bf x}|{\bf z})p({\bf y}|{\bf z})}\right]$ between the
random vectors ${\bf x}$ and ${\bf y}$, given ${\bf z}$. Thus KL divergence is
also Shannon’s definition of the rate $R$ at which the MVN error ${\bf
y}-\hat{{\bf y}}({\bf z})$ brings information about the MVN error ${\bf
x}-\hat{{\bf x}}({\bf z})$. This is also the rate at which ${\bf y}$ brings
information about ${\bf x}$, given ${\bf z}$, and a monotone function of what
the Granger-Geweke causality measure would say for finite-dimensional random
vectors $({\bf x},{\bf y},{\bf z})$. So partial coherence, conditional KL
divergence, conditional mutual information, conditional information rate, and
a Granger-Geweke causality measure are monotone functions of each other. All
are resolved by the partial canonical correlations $k_{i}$:

$\displaystyle D_{KL}(P_{1}||P_{0})=$ $\displaystyle I_{{\bf x}{\bf y}|{\bf
z}}=-\frac{1}{2}\log(1-\rho^{2}_{{\bf x}{\bf y}|{\bf z}})$ (58)
$\displaystyle=$
$\displaystyle-\frac{1}{2}\sum_{i=1}^{\makebox{min}(p,q)}\log(1-k_{i}^{2})$
(59)

These finite-dimensional connections are consistent with the finding in [13]
that the Granger-Geweke causality measure and transfer entropy are equivalent
for MVN variables. These equivalent measures are invariant to block diagonal
nonsingular transformations of the random vectors ${\bf x},{\bf y},{\bf z}$,
and they are functions of the partial canonical correlations $k_{i}$, which
are maximal invariants.

When the partial canonical correlations $k_{i}$ are near to zero, indicating
the coherence matrix ${\bf C}_{{\bf x}{\bf y}|{\bf z}}$ is near to the zero
matrix, then partial coherence is near to zero, and conditional KL divergence,
conditional mutual information, and conditional rate are near to zero. The
Granger-Geweke causality measure is near to one, and its transfer entropy is
near to zero.

##  IV  Partial Coherence, Likelihood, and a Test for Causality

Partial coherence is inherently connected to likelihood in the case where the
random vectors $({\bf x},{\bf y},{\bf z})$ are jointly normally distributed.
The matrix ${\bf R}_{{\bf u}{\bf u}|{\bf z}}$ is the conditional error
covariance matrix for the error $({\bf u}-\hat{{\bf u}}({\bf z}))$; when the
partial correlation ${\bf R}_{{\bf x}{\bf y}|{\bf z}}$ is zero, then ${\bf
R}_{{\bf u}{\bf u}|{\bf z}}=\makebox{blkdiag}[{\bf R}_{{\bf x}{\bf x}|{\bf
z}},{\bf R}_{{\bf y}{\bf y}|{\bf z}}]$.

The null hypothesis test of interest is $\makebox{H}_{0}:({\bf u}-\hat{{\bf
u}}({\bf z}))\sim\mathcal{N}_{p+q}({\bf 0},\makebox{blkdiag}[{\bf R}_{{\bf
x}{\bf x}|{\bf z}},{\bf R}_{{\bf y}{\bf y}|{\bf z}}])$ versus the alternative
$\makebox{H}_{1}:({\bf u}-\hat{{\bf u}}({\bf z}))\sim\mathcal{N}_{p+q}({\bf
0},{\bf R}_{{\bf u}{\bf u}|{\bf z}})$.

###  IV-A Likelihood

Given the data matrix ${\bf D}^{T}=({\bf X}^{T},{\bf Y}^{T},{\bf Z}^{T})=({\bf
U}^{T},{\bf Z}^{T})$, where each of the matrices ${\bf
X}\in{\mathbb{R}}^{p\times M},{\bf Y}\in{\mathbb{R}}^{q\times M},{\bf
Z}\in{\mathbb{R}}^{r\times M}$ consists of $M$ i.i.d. realizations of the
jointly multivariate normal random vectors ${\bf x},{\bf y},{\bf z}$, compute
the sample covariance matrix ${\bf S}={\bf D}{\bf D}^{T}$, structured as
follows:

$\displaystyle{\bf S}$ $\displaystyle=\left[\begin{array}[]{cccc}{\bf S}_{{\bf
x}{\bf x}}&{\bf S}_{{\bf x}{\bf y}}&{\bf S}_{x{\bf z}}\\\ {\bf S}_{{\bf y}{\bf
x}}&{\bf S}_{{\bf y}{\bf y}}&{\bf S}_{y{\bf z}}\\\ {\bf S}_{{\bf z}x}&{\bf
S}_{{\bf z}y}&{\bf S}_{{\bf z}{\bf
z}}\end{array}\right]=\left[\begin{array}[]{cccc}{\bf S}_{{\bf u}{\bf u}}&{\bf
S}_{{\bf u}{\bf z}}\\\ {\bf S}_{{\bf z}{\bf u}}&{\bf S}_{{\bf z}{\bf
z}}\end{array}\right].$ (65)

For example, ${\bf S}_{{\bf x}{\bf y}}={\bf X}{\bf Y}^{T}$.

The ordinary likelihood for testing $\makebox{H}_{0}$ versus $\makebox{H}_{1}$
is then [22]

$1-\hat{\rho}^{2}_{{\bf x}{\bf y}|{\bf z}}=\frac{\det[{\bf S}_{{\bf u}{\bf
u}|{\bf z}}]}{\det[{\bf S}_{{\bf x}{\bf x}|{\bf z}}]\det[{\bf S}_{{\bf y}{\bf
y}|{\bf z}}]},$ (66)

where

$\displaystyle{\bf S}_{{\bf u}{\bf u}|{\bf z}}={\bf S}_{{\bf u}{\bf u}}-{\bf
S}_{{\bf u}{\bf z}}{\bf S}_{{\bf z}{\bf z}}^{-1}{\bf S}_{{\bf u}{\bf z}}^{T},$
(67) $\displaystyle{\bf S}_{{\bf x}{\bf x}|{\bf z}}={\bf S}_{{\bf x}{\bf
x}}-{\bf S}_{{\bf x}{\bf z}}{\bf S}_{{\bf z}{\bf z}}^{-1}{\bf S}_{{\bf x}{\bf
z}}^{T},$ (68) $\displaystyle{\bf S}_{{\bf y}{\bf y}|{\bf z}}={\bf S}_{{\bf
y}{\bf y}}-{\bf S}_{{\bf y}{\bf z}}{\bf S}_{{\bf z}{\bf z}}^{-1}{\bf S}_{{\bf
y}{\bf z}}^{T}.$ (69)

Of course, as noted previously, the matrix ${\bf S}_{{\bf u}{\bf u}|{\bf
z}}^{-1}$ may be read out of the inverse ${\bf S}^{-1}$ and used to compute
the determinant $\det[{\bf S}_{{\bf u}{\bf u}|{\bf z}}]$.

The following proposition establishes the distribution, and a stochastic
representation, for the likelihood ratio statistic $\hat{\rho}^{2}_{{\bf
x}{\bf y}|{\bf z}}$ under the null hypothesis. This result forms the basis for
setting a detection threshold to achieve a desired type-I error.

######  Proposition 1.

Under the null hypothesis $\makebox{H}_{0}$ that the coherence matrix ${\bf
C}_{{\bf x}{\bf y}|{\bf z}}$ is the zero matrix, the likelihood ratio
statistic $1-\widehat{\rho}^{2}_{{\bf x}{\bf y}|{\bf z}}$ has the Wilks
$\Lambda$ distribution, $\Lambda(p,M-r-q,q)$ with stochastic representation

$\widehat{\rho}^{2}_{{\bf x}{\bf y}|{\bf z}}\sim 1-\prod_{i=1}^{p}b_{i},$

where the $b_{i}$ are independent Beta random variables with
$\makebox{Beta}(\frac{M-r-q-i+1}{2},\frac{q}{2})$ distributions.

###### Proof.

Let ${\bf U}_{1}={\bf X}-{\bf X}{\bf Z}^{T}({\bf Z}{\bf Z}^{T})^{-1}{\bf Z}$,
${\bf U}_{2}={\bf Y}-{\bf Y}{\bf Z}^{T}({\bf Z}{\bf Z}^{T})^{-1}{\bf Z}$, and
${\bf U}=[{\bf U}_{1}^{T},{\bf U}_{2}^{T}]^{T}$. It follows that ${\bf U}{\bf
U}^{T}={\bf S}_{{\bf u}{\bf u}|{\bf z}}$, ${\bf U}_{1}{\bf U}_{1}^{T}={\bf
S}_{{\bf x}{\bf x}|{\bf z}}$, ${\bf U}_{2}{\bf U}_{2}^{T}={\bf S}_{{\bf y}{\bf
y}|{\bf z}}$, and therefore,

$\displaystyle 1-\widehat{\rho}^{2}_{{\bf x}{\bf y}|{\bf z}}=\frac{\det({\bf
U}{\bf U}^{T})}{\det({\bf U}_{1}{\bf U}_{1}^{T})\det({\bf U}_{2}{\bf
U}_{2}^{T})}.$ (70)

Kshirsagar [23] has shown that the error covariance ${\bf U}{\bf U}^{T}$ has a
Wishart distribution $W_{p+q}({\bf R}_{{\bf u}{\bf u}|{\bf z}},M-r)$ where $r$
is the dimension of ${\bf z}$. The statistic $\frac{\det({\bf U}{\bf
U}^{T})}{\det({\bf U}_{1}{\bf U}_{1}^{T})\det({\bf U}_{2}{\bf U}_{2}^{T})}$ is
invariant to left multiplication of $({\bf X},{\bf Y})$ by a block-diagonal
matrix $\left[\begin{array}[]{ccc}{\bf T}_{1}&0\\\ 0&{\bf
T}_{2}\end{array}\right]$ for non-singular $p\times p$ matrix ${\bf T}_{1}$
and $q\times q$ matrix ${\bf T}_{2}$. Choose ${\bf T}_{1}={\bf R}_{{\bf x}{\bf
x}|{\bf z}}^{-1/2}$ and ${\bf T}_{2}={\bf R}_{{\bf y}{\bf y}|{\bf z}}^{-1/2}$,
so that WLOG ${\bf U}\sim\mathcal{N}({\bm{0}},{\bf I}_{p+q}\otimes{\bf
I}_{M-r})$ which yields ${\bf U}_{1}\sim\mathcal{N}({\bm{0}},{\bf
I}_{p}\otimes{\bf I}_{M-r})$ and ${\bf U}_{2}\sim\mathcal{N}({\bm{0}},{\bf
I}_{q}\otimes{\bf I}_{M-r})$. Moreover, ${\bf U}_{1}$ and ${\bf U}_{2}$ are
independent under the null hypothesis

Let $[{\bf V},{\bf V}^{\perp}]\in{\mathbb{R}}^{(M-r)\times(M-r)}$ be the
orthogonal matrix where ${\bf V}\in{\mathbb{R}}^{(M-r)\times q}$ is a unitary
basis for the row space of ${\bf U}_{2}$ and ${\bf
V}^{\perp}\in{\mathbb{R}}^{(M-r)\times(M-r-q)}$ is the unitary basis for its
orthogonal space. We can choose ${\bf V}={\bf U}_{2}^{T}({\bf U}_{2}{\bf
U}_{2}^{T})^{-1}{\bf U}_{2}$ and define ${\bf W}_{1}={\bf U}_{1}{\bf V}$ and
${\bf W}_{2}={\bf U}_{1}{\bf V}^{\perp}$. Conditional on ${\bf V}$ and ${\bf
V}^{\perp}$, it can be shown that ${\bf W}_{1}\sim\mathcal{N}(0,{\bf
I}_{p}\otimes{\bf I}_{q})$, and ${\bf W}_{2}\sim\mathcal{N}(0,{\bf
I}_{p}\otimes{\bf I}_{M-r-q})$ (assuming $M-r>p+q$) which are invariant to
${\bf V}$. Simple linear algebra leads to

$\displaystyle 1-\widehat{\rho}^{2}_{{\bf x}{\bf y}|{\bf z}}=\frac{\det({\bf
W}_{2}{\bf W}_{2}^{T})}{\det({\bf W}_{1}{\bf W}_{1}^{T}+{\bf W}_{2}{\bf
W}_{2}^{T})}.$

This has the Wilk’s $\Lambda$ distribution $\Lambda(p,M-r-q,q)$ [24], which is
identical to the product of $p$ independent beta-distributed random variables:

$\displaystyle 1-\widehat{\rho}^{2}_{{\bf x}{\bf y}|{\bf
z}}\sim\prod_{i=1}^{p}b_{i}\text{ where }b_{i}\sim\makebox{Beta}(\frac{M-r-
q-i+1}{2},\frac{q}{2}).$

∎

It is important to note that this distribution is invariant to the actual
error covariances ${\bf R}_{{\bf x}{\bf x}|{\bf z}}$ and ${\bf R}_{{\bf y}{\bf
y}|{\bf z}}$, which are themselves functions of the underlying composite
covariance matrix ${\bf R}$ of eqn (2). Moreover, this distribution result
specializes to the case $r=0$, in which case partial coherence is ordinary
multivariate-multivariate coherence, and the distribution of the coherence
statistic is $1-\Lambda(p,M-q,q)$. When $p=1$, then partial coherence is
multivariate coherence, with distribution $1-\Lambda(1,M-q,q)$. In [25]
Bartlett derived the asymptotic distribution for large $M$,

$-(M-r-\frac{p+q+1}{2})\log(1-\widehat{\rho}^{2}_{{\bf x}{\bf y}|{\bf
z}})\sim\mathcal{X}^{2}_{pq},$

where $\chi^{2}_{pq}$ denotes a chi-squared sdistribution with $pq$ degrees of
freedom.

##  V Experiments

If partial coherence is to be used to answer questions of causal influence
between signal sequences, then there is required a princpled way to map a
causality question into a question of conditional dependence among random
vectors.

Our approach to causality analysis of time series will be to start with two
channels of signals, call them $\\{x_{n}\\}$ and $\\{y_{n}\\}$, and then to
define a third signal channel $\\{z_{n}\\}$ of prima facie evidence in such a
way that the question of causality may be resolved from an analysis of partial
coherence. For example the random vectors $({\bf x},{\bf y},{\bf z})$ might be
${\bf x}=x_{s}$ and ${\bf y}=y_{t}$, which are readings of the signal
sequences at the two times $s$ and $t$. In some applications, the prima facie
evidence ${\bf z}$ might be defined as the past of the time series
$\\{y_{n}\\}$ at time $t$, i.e., ${\bf z}=\\{y_{n},n<t\\}$. In other
applications the prima facie evidence might be ${\bf z}=\\{x_{n},n<t\\}$. As a
practical matter in testing from experimental data, this past will be finite.
More commonly, ${\bf x}$ is assembled from the finite past of $\\{x_{n}\\}$ up
to time $t-1$, ${\bf y}$ is defined to be $y_{t}$, and ${\bf z}$ is assembled
from the finite past of $\\{y_{n}\\}$ up to time $t-1$. With ${\bf x},{\bf
y},{\bf z}$ so defined, we can compute the partial coherence statistic
$\rho_{{\bf z}}^{2}(s,t)\doteq\rho^{2}_{x_{s}y_{t}|{\bf z}},$ or more
generally $\rho_{y_{t}{\bf x}|{\bf z}}^{2}$. The statistic $\rho_{{\bf
z}}^{2}(s,t)$ vs the pair $(s,t)$ reveals pairwise causality with respect to
prima facie evidence ${\bf z}$. The statistic $\rho_{y_{t}{\bf x}|{\bf
z}}^{2}$ reveals the causal influence of the finite past of $\\{x_{n}\\}$ on
$y_{t}$, conditioned on the finite past of $\\{z_{n}\\}$.

A value of partial coherence below a threshold is taken to be a finding of
non-causality at significance level $\alpha$, whereas a value of coherence
above this threshold is taken only to be an indication of causality with
respect to this specific prima facie evidence and a second-order analysis of
partial coherence. New evidence and/or new methods of analysis may invalidate
an indication of causality. In the multivariate normal case, no new evidence
or method of analysis can invalidate the finding of non-causality, as we have
found prima facie evidence and a method of analysis that indicates non-
causality at significance level $\alpha$. An equivalent view is that the
likelihood ratio test is a null hypothesis test for non-causality.

To demonstrate the use of partial coherence for resolving the question of
causality, we present in this section a demonstration and an experiment. The
demonstration is an analytical calculation of partial coherence for a known
bi-variate time series which may model past-causal, future-causal, or mixed-
causal influence of one time series upon the other, depending on the choice of
coupling parameters between the two time series. The experiment is a re-
running of the Barnett and Seth experiment [12].

###  V-A Demonstration of Partial Coherence for Past-causal, Future-causal,
and Mixed-Causal Models

Consider the following multivariate linear system,

$\displaystyle{\bf x}_{n}$ $\displaystyle=\sum_{k=-\infty}^{+\infty}{\bf
H}_{k}{\bm{\mu}}_{n-k}$ $\displaystyle{\bf y}_{n}$
$\displaystyle=\sum_{k=-\infty}^{+\infty}{\bf
G}_{k}{\bm{\nu}}_{n-k}+\sum_{k=-\infty}^{+\infty}{\bf F}_{k}{\bf x}_{n-k},$

where ${\bm{\mu}},{\bm{\nu}}$ are multivariate unit variance white noises. The
two time series are correlated with auto- and cross-covariance matrices

$\displaystyle{\bf R}_{{\bf x}{\bf x}}[m]$ $\displaystyle={\mathbb{E}}[{\bf
x}_{n}{\bf x}_{n+m}]=\sum_{k}{\bf H}_{k}{\bf H}_{k+m},$ $\displaystyle{\bf
R}_{{\bf y}{\bf y}}[m]$ $\displaystyle={\mathbb{E}}[{\bf y}_{n}{\bf y}_{n+m}]$
$\displaystyle=\sum_{k}{\bf G}_{k}{\bf G}_{k+m}+\sum_{k,l}{\bf F}_{k}{\bf
R}_{{\bf x}{\bf x}}[m+k-l]{\bf F}_{l},$ $\displaystyle{\bf R}_{{\bf x}{\bf
y}}[m]$ $\displaystyle={\mathbb{E}}[{\bf x}_{n}{\bf y}_{n+m}]=\sum_{k}{\bf
F}_{k}{\bf R}_{{\bf x}{\bf x}}[m-k].$

When these infinite-dimensional moving average (MA) models arise as finitely-
parameterized autoregressive-moving average (ARMA) models, then these infinite
sums may be evaluated analytically.

For our demonstration of partial coherence, we assume the time series are
scalar-valued. We set the filter coefficients at $h_{k}=h_{0}a^{k}$ and
$g_{k}=g_{0}b^{k}$ for $k\geq 0$ where $h_{0}=0.8$, $a=0.1$, $g_{0}=0.7$, and
$b=0.7$. The coefficients $\\{f_{k}\\}$, are the filtering coefficients that
determine the causality between $\\{y_{n}\\}$ and $\\{x_{n}\\}$. For any
selected pair of times $s,t$, we evaluate the partial coherence between
$x_{s}$ and $y_{t}$. We consider two sets of prima facie evidence: ${\bf z}$
is the past of $\\{y_{n}\\}$ up to time $t-1$, and ${\bf z}$ is the past of
$\\{x_{n}\\}$ up to time $t$ with $x_{s}$ excluded. Partial coherence is
therefore $\rho^{2}_{{\bf z}}(s,t)\doteq\rho^{2}_{y_{t}x_{s}|{\bf z}}$.

We consider three different scenarios. Case I: conditioned on the finite past
of the input $x$, the output $y$ is linearly independent of any input outside
this finite past which includes the future of $x$; Case II: conditioned on the
finite future of the input $x$ the output $y$ is linearly independent of any
input outside this finite future which includes the past of $x$; Case III;
conditioned the finite past and finite future of the input $x$, the output $y$
is independent of any inputs outside this finite past and future.

Case I: $\\{y_{n}\\}$ is past-causally dependent on $\\{x_{n}\\}$:

$\displaystyle y_{n}=$
$\displaystyle.7\sum_{k=0}^{+\infty}(.7)^{k}\nu_{n-k}+.7x_{n}+.8x_{n-1}+.7x_{n-2}+.6x_{n-3}$
$\displaystyle x_{n}=$ $\displaystyle.8\sum_{k=0}^{+\infty}(.1)^{k}\mu_{n-k}.$

Case II: $\\{y_{n}\\}$ is future-causally dependent on $\\{x_{n}\\}$:

$\displaystyle y_{n}=$
$\displaystyle.7\sum_{k=0}^{+\infty}(.7)^{k}\nu_{n-k}+.7x_{n}+.8x_{n+1}+.7x_{n+2}+.3x_{n+3}$
$\displaystyle x_{n}=$ $\displaystyle.8\sum_{k=0}^{+\infty}(.1)^{k}\mu_{n-k}.$

Case III: the mixed case:

$\displaystyle y_{n}=$
$\displaystyle.7\sum_{k=0}^{+\infty}(.7)^{k}\nu_{n-k}+.4x_{n-2}+.8x_{n-1}$
$\displaystyle\quad+.7x_{n}+.8x_{n+1}+.4x_{n+2}$ $\displaystyle x_{n}=$
$\displaystyle.8\sum_{k=0}^{+\infty}(.1)^{k}\mu_{n-k}.$

Figure 1: Pairwise partial coherence between $x_{s}$ at time $s$ (horizontal
axis) and $y_{t}$ at time $t$ (vertical axis). The prima facie evidence ${\bf
z}$ is the past of $\\{x_{n}\\}$ up to time $t$ in the left column, and the
past of $\\{y_{n}\\}$ up to time $t-1$ in the right column; (top row)
$\\{y_{n}\\}$ is past-causally dependent on $\\{x_{n}\\}$; (middle row)
$\\{y_{n}\\}$ is future-causally dependent on $\\{x_{n}\\}$; (bottom row)
$\\{y_{n}\\}$ has mixed past- and future-causal dependence on $\\{x_{n}\\}$.

The resulting pairwise partial coherence maps are illustrated in Figure 1. The
three rows show the results of case I, II, and III, respectively. Each panel
of Figure 1 illustrates the partial coherence between $y_{t}$ and $x_{s}$,
with the prima facie evidence being the past of $\\{x_{n}\\}$ on the left
column, and the past of $\\{y_{n}\\}$ on the right column. As we already know,
the partial coherence are conditioned on available prima facie evidence. With
${\bf z}=\\{x_{t},x_{t-1},\ldots\\}$, (Notice that $x_{s}$ is excluded from
${\bf z}$ if $s\leq t$.), $\rho^{2}_{{\bf z}}(s,t)=0$ for any $s>t$ or $s<t-3$
in Case I (top row), which suggests that future and long-past values of
$\\{x_{n}\\}$ do not contain any second-order information about $y_{t}$ given
this ${\bf z}$ (In fact, they are conditionally independent in the MVN case).
In Case II the future-causal case (middle row), $\rho^{2}_{{\bf z}}(s,t)=0$
for any $s\leq t$ since in this case the past of $\\{x_{n}\\}$ does not
contain unique information about $y_{t}$ given ${\bf z}$. In Case III (bottom
row) the mixed case where $y_{t}$ depends on both the past and future of $x$,
$\rho^{2}_{{\bf z}}(s,t)>0$ for both $s\leq t$ and $s>t$. In summary, there is
an asymmetry in the past-causal and future-causal cases that bears comment. In
Case I, $\rho^{2}_{{\bf z}}(s,t)$ for $s<t$ is non-zero, but weak, compared
with $\rho^{2}_{{\bf z}}(s,t)$ for $s>t$ in Case II. This effect is explained
by noting that for these demonstrations our conditioning is always on the past
of $\\{x_{n}\\}$ at time $t$. The right column of Figure 1 shows the partial
coherence between $y_{t}$ and $x_{s}$ given the past of $\\{y_{n}\\}$, i.e.,
${\bf z}=\\{y_{t-1},y_{t-2},\ldots\\}$ for the three cases. Here the partial
correlation generally is nonzero for both $s<t$ and $s>t$. The reason is that
$\\{x_{n}\\}$ is an auto-regressive process and $\\{y_{n}\\}$ is a finite
moving-average of $\\{x_{n}\\}$. When the prima facie is the past of
$\\{y_{n}\\}$, $x_{s}$ still contains some unique information about $y_{t}$
due to the infinite memory of $\\{x_{n}\\}$. While all three cases are Granger
causal, the choice of prima facie evidence as the past of input time series
$\\{x_{n}\\}$ leads to clear separation between past-causal, future-causal,
and the mixed-causal influence, whereas the conventional choice of prima facie
evidence as the past of the output time series $\\{y_{n}\\}$ does not resolve
the direction of causal influence.

###  V-B The experiment of Barnett, et al.

In this subsection, we replicate the experiment in Barnett et al. [12] and use
the partial coherence statistic to test for time-domain causality. As in [12],
we define ${\bf y}=y_{t}$, ${\bf z}=[y_{t-1},y_{t-2},\ldots]$, and ${\bf
x}=[x_{t-1},x_{t-2},\ldots]$. The causality from time series $\\{x_{n}\\}$ to
$\\{y_{n}\\}$ is measured by the partial coherence between ${\bf x}$ and ${\bf
y}$ given ${\bf z}$. Straightforward calculation shows that partial coherence
$\rho^{2}_{{\bf x}{\bf y}|{\bf z}}$ is a monotone function of the causality
measure in [12]. However, as a statistic, our estimation of the causality
measure is very different. The test statistic in Barnett et al. is computed by
replacing state-space parameters by their approximate maximum likelihood
estimates, and inserting these estimates into an analytical formula for their
ratio of determinants. The result is an approximate likelihood ratio for an
ARMA model of measurements. Our test statistic is the ordinary likelihood
ratio for testing the null hypothesis that partial correlation is zero. It is
model-free, requiring only the computation of partial coherence from sample
correlations. To compare the two approaches, we replicate the experiment in
[12] by generating bivariate $ARMA(r,1)$ time series data.

Figure 2: Left: The time series model in the simulation study from Barnett et
al. [12] where $\mu$ and $\nu$ are time series of white noise. Right: Pairwise
partial coherence between $x_{s}$ at time $s$ ( horizontal axis) and $y_{t}$
at time $t$ (vertical axis) with the prima facie evidence ${\bf z}$ is the
past of $\\{y_{n}\\}$ up to time $t-1$.

Start with a bivariate auto-regressive (AR) time series
$\\{[\eta_{1t},\eta_{2t}]^{T}\\}$,

$\displaystyle(1-az)\eta_{1t}$ $\displaystyle=cz\eta_{2,t}+\mu_{t},$
$\displaystyle(1-bz)\eta_{2t}$ $\displaystyle=\nu_{t},$

where $z$ is the backward shift or delay operator:
$z\\{\eta_{t}\\}=\\{\eta_{t-1}\\}$. The bivariate auto-regressive moving
average (ARMA) time series are generated as

$\displaystyle\left[\begin{array}[]{cc}y_{t}\\\
x_{t}\end{array}\right]=\left[\begin{array}[]{cc}(1+f_{1}z)^{r}&0\\\
0&(1+f_{2}z)^{r}\end{array}\right]\cdot\left[\begin{array}[]{cc}\eta_{1t}\\\
\eta_{2t}\end{array}\right].$ (77)

The circuit diagram in the left panel of Figure. 2 shows clearly the
dependence of the time series $\\{y_{n}\\}$ on $\\{x_{n}\\}$ through the
coupling parameter $c$. The MA order $r$ varies from $0$ to $10$, but there
are only two parameters, $f_{1}=0.6$ and $f_{2}=0.7$ determining these $r$ MA
parameters. The AR parameters are set to $a=0.9,b=0.8$,
$c=\sqrt{e^{-F}(e^{F}-1)(e^{F}-b^{2})}$ with the transfer entropy $F=0.02$.
The true partial coherence in this case is $\rho_{{\bf x}{\bf y}|{\bf
z}}^{2}=1-e^{-0.02}\approx 0.02$ and $GG=e^{-0.02}$.

To replicate the Barnett, et al. experiment, we choose $1000$ time samples to
estimate $\rho^{2}_{{\bf x}{\bf y}|{\bf z}}$ from second-order statistics, and
this experiment is replicated $10,000$ times to obtain Monte-Carlo estimates
of performance. The vector ${\bf y}$ is defined to be $y_{t}$; ${\bf x}$ and
${\bf z}$ are truncated to finite $T$-dimensional vectors with ${\bf
x}=[x_{t-1},\ldots,x_{t-T}]$ and ${\bf z}=[y_{t-1},\ldots,y_{t-T}]$. So the
sample covariance matrices for these vectors are just $10\times 10$. Using the
null distribution given in Theorem 1, we set the rejection region with
significance level $0.05$. As a practical matter, $(10,000)\times(1000)$
samples are generated and consecutive windows of $1000$ samples are
constructed to obtain the 10,000 replications. Because of the weak dependence
between windows, these are not actually $1000$ independent realizations, so
our analytical formula for the null distribution does not perfectly predict
the experimentally-determined significance level. But it is very close. In
fact, when we design for a significance level of $\alpha=0.05$ we achieve an
experimental level of $0.044.$

Figure 3 displays our results based on the model-free test in Section IV. The
left panel plots the power of our test statistic with respect to the order of
the underlying ARMA model, at significance level $0.05$. This power is
relatively invariant to MA order, whereas the experimental results in [12]
show decreasing power as the MA order increases, as a result of the increased
model complexity. The model-free coherence statistic for memory $T=10$ has
power about $0.90$ regardless of the MA order. This power is slightly lower
than the result of Barnett, et al. at low MA orders, and higher at higher MA
orders. Power is relatively invariant to MA order, and it is not susceptible
to model mismatch, as no state space model is assumed. An appropriate value of
$T$ can be selected by cutting off the empirical auto-correlations at a pre-
determined threshold or comparing AIC or BIC values for an auto-regressive
model fitted to the sampled time series data. In the right panel of Figure 3,
we plot the ROC curve (power vs size) to demonstrate that probability of
correctly rejecting the null hypothesis (the power) depends on the probability
of incorrectly rejecting it (the size). We also plot in the right panel of
Figure 2 the pairwise partial coherence between $y_{t}$ and $x_{s}$ given the
prima facie evidence ${\bf z}=[y_{t-1},\ldots,y_{t-T}]$. These pairwise
partial coherences are not the partial correlation $\rho^{2}_{{\bf x}{\bf
y}|{\bf z}}$ that is used in the hypothesis test. Rather they are pairwise
partial coherences, $\rho^{2}_{x_{s},y_{t}|{\bf z}}$, plotted only for
intuition, for the case $c=1$. It can be seen that the partial coherence given
the past of $\\{y_{n}\\}$ is non-zero for short ranges of $s>t$ and $s<t$.

Figure 3: The simulation study from Barnett et al. [12]. (Left) Power versus
the MA order of the ARMA time series. (Right) Receiver operating
characteristic (ROC) showing power versus size.

##  VI Discussion

In this paper, our objective has been to explore the use of partial coherence
among finite-dimensional vectors, constructed from signal sequences, as a
scale-invariant statistic to assess causal influence of one signal sequence on
another. Of course these signal sequences might well be space series or space-
time series, in which case the idea of past-, future- or mixed- causal
influence is also a physically meaningful question [26]. The partial coherence
statistic we argue for is a monotone function of conditional KL divergence,
conditional mutual information, conditional information rate, and the Granger-
Geweke causality measure in the multivariate normal case. It is a measure of
second-order partial correlation in the Hilbert space of second-order random
variables, and it is resolved by partial canonical correlations. The sampled-
data version of partial coherence is in fact a likelihood ratio statistic
whose null distribution is a Wilks Lambda. It may be computed from the sample
information matrix. Importantly, from this null distribution, a threshold on
partial coherence may be set to control the confidence level $1-\alpha$.

The key idea in the application of partial coherence to questions of causal
influence is to define three channels of signals so that they code for the
question of interest. Then finite-dimensional vectors are constructed from
these channels and their composite covariance matrix is estimated without
appeal to an underlying generative model. Thus the method is model-free.

We have demonstrated the methodology for two particular experiments in causal
influence. One of these experiments demonstrates how partial coherence may be
used to map out a partial coherence surface, parameterized by two time
variables, and the other replicates the experiment of Barnett and Seth [12].
For this latter experiment we demonstrate the power of our method and plot its
experimentally computed Receiver Operating Characteristic (ROC). The
conclusion of the test is that the performance of a model-free coherence test
for causal influence is competitive with a model-based coherence test.
Moreover, the null distribution of the model-free test is invariant to the
underlying generative model for the measurements. Consequently, as a null
hypothesis test, the partial coherence test is not vulnerable to model bias.
Of course in applications where a model might be well-established, and
measurements are plentiful, a method like that of the Barnetts and Seth [12],
[13] might outperform a model-free method such as model-free partial
coherence. In applications where there is no generative model like an ARMA
model, then the method of second-order partial coherence is compelling.

Acknowledgment. We thank Profs. Sonia and Tomas Charleston Villalobos for
piquing our interest in causality, and Prof. Todd Moon for alerting us to the
work of the Barnetts and Seth, which led us to the work of Geweke.

## References

  * [1] N. Wiener, “The Human Use of Human Beings”, Avon Books, New York, NY, 1950.
  * [2] C. W. J. Granger,  “Investigating Causal Relations by Econometric Models and Cross-spectral Methods,”  Econometrica, 37 (3): 424–438, 1969.
  * [3] C. W. J. Granger,  “Testing for Causality: a Personal View,”  Journal of Economic Dynamics and Control, 2: 329–352, 1980.
  * [4] J. Rissanen and M. Wax,  “Measures of mutual and causal dependence between two time series (Corresp.),”  IEEE Trans Information Theory , 33(4), 598–601, 1987.
  * [5] J. Geweke,  “Measurement of linear dependence and feedback between multiple time series,”  J. Am. Stat. Assoc. 77(378): 304–313, 1982.
  * [6] J. Geweke,  “Measures of conditional linear dependence and feedback between time series,”  J. Am. Stat. Assoc.,79(388): 907–915, 1984.
  * [7] V. Solo, “State-Space Analysis of Granger-Geweke Causality Measures with Application to fMRI,” Neural Computation, 28(5), 914–949, 2016.
  * [8] E. Siggiridou and D. Kugiumtzis,  “Granger Causality in Multivariate Time Series Using a Time-Ordered Restricted Vector Autoregressive Model,”  IEEE Transactions on Signal Processing, 64(7):1759–1773, 2016.
  * [9] G. Schamberg and T. P. Coleman,  “Measuring Sample Path Causal Influences With Relative Entropy,”  IEEE Transactions on Information Theory, 66(5):2777–2798, 2020.
  * [10] Kramer, G., “Directed information for channels with feedback,” Ph.D. thesis, Swiss Federal Institute of Technology Zurich, 1998.
  * [11] P-O. Amblard and J.J. Michel, “On directed information theory and Granger causality graphs,” J. Computational Neuroscience,30(1), 7–16, 2011.
  * [12] L. Barnett and A. K. Seth,  “Granger causality for state-space models,”  Physical Review E 91:040101(R), 2015.
  * [13] L. Barnett, A. B. Barrett, and A. K. Seth,  “Granger causality and transfer entropy are equivalent for Gaussian variables,”  Phys. Rev. Lett., 103(23), 238701, 2009.
  * [14] L. Barnett, A. B. Barrett, and A. K. Seth,  “Solved problems for Granger causality in neuroscience: A response to Stokes and Purdon,”  Neuroimage, 178:744–748, 2018.
  * [15] D. Ramirez, J. Via, I. Santamaria, and L.L. Scharf, “Detection of Spatially-Correlated Gaussian Time Series,” IEEE Trans Signal Processing, 58(10), 5006–5015, 2010.
  * [16] J. Whittaker,  “Graphical Models in Applied Multivariate Statistics,”  Wiley, 2009.
  * [17] M. L. Eaton, “Multivariate Statistics: A Vector Space Approach,” Wiley Interscience,New York, 1983.
  * [18] H. Cramer, “Mathematical Methods of Statistics”, Princeton University Press, Princeton, 1946
  * [19] Szego, G, “Ein Grenzwertsatz über die Toeplitzschen Determinanten einer reellen positiven Funktion”. Math. Ann. 76 (4): 490–503, 1915.
  * [20] A. Bottcher, B. Silbermann. “Toeplitz determinants”. Analysis of Toeplitz operators. Berlin: Springer-Verlag. p. 525, 1990
  * [21] L. L. Scharf and C. T. Mullis, “Canonical Coordinates and the Geometry of Inference, Rate, and Capacity,” IEEE Trans Signal Processing, 48(3), 824–831, 2000.
  * [22] K. V. Mardia, J. T. Kent, and J. M. Bibby, Multivariate Analysis, Academic Press, 133–138, 1979.
  * [23] A. M. Kshirsagar  “Multivariate Analysis,”  Marcel Dekker Inc, 1972.
  * [24] S. S. Wilks, “Certain generalizations in the analysis of variance,” Biometrika, vol 24, 471–494, 1932.
  * [25] Bartlett, M. S.  “A Note on the Multiplying Factors for Various $\mathcal{X}^{2}$ Approximations,”  J R Stat Soc Series B. 16(2): 296–298, 1954.
  * [26] I. Ebert-Uphoff, and Y. Deng,  “Causal Discovery for Climate Research Using Graphical Models,”  J. Climate, 25:5648–5665, 2012.