Star graphs induce tetrad correlations:
for Gaussian as well as for binary variables

Since the seminal work by Bartlett (1935) and Birch (1963), viewed in a larger perspective by Cox (1972), correlation coefficients were barely used as measures of dependence for categorical data. Instead, functions of the odds-ratio emerged as the relevant parameters in log-linear probability models for joint distributions and in logistic models, that is for regressions with a binary response. Compared with other possible measures of dependence, the outstanding advantage of functions of the odds-ratio is their variation independence of the marginal counts: odds-ratios are unchanged under different types of sampling schemes that result by fixing either the total number, $n$, of observed individuals or the counts at the levels of one of the variables; see Edwards (1963).

As one consequence of the importance of odds-ratios for discrete random variables, it is no longer widely known that Pearson’s simple, observed correlation coefficient, $r_{12}$ say, coincides in $2\times 2$ contingency tables with the so-called Phi-coefficient, so that $\sqrt{n}r_{12}$ is, asymptotically and under independence of the two binary variables, a realization of a standard Gaussian distribution. As we shall see, some properties of correlation coefficients for binary variables, make them important for data generating processes that incorporate many conditional independences.

In particular, we look here at directed star graphs such as the one shown in Figure 1(a). Such graphs have one inner node, $L$, from which $Q$ arrows start and point to the uncoupled, outer nodes $1,\dots,Q$. To simplify notation, the inner node $L$ denotes also the corresponding random variable and both, the node and the variable, are called root. The random variables $X_{i}$ corresponding to the outer nodes, also called the leaves of the graph, are identified just by their index $i$ taken from $\{1,\dots,Q\}$.

The independence structure of the directed star graph is mutual independence of the leaves given the root. In the condensed node notation, this is written as

In general, the types of variables can be of any kind. Densities, $f_{N}$, that are said to be generated over the directed star graph with node set $N=\{1,\dots,Q,L\}$, may also be of any kind, provided that they have the above independence structure.

Each generated density is defined by $Q$ conditional densities, $f_{i|L}$, that are called regressions, and a marginal density, $f_{L}$, of the root. In the condensed node notation, a joint density with the independence structure of a directed star graph, factorizes as

Directed star graphs belong to the class of regression graphs and to their subclass of directed acyclic graphs. Distributions generated over regression graphs have been named and studied as sequences of regressions and one knows when two regression graphs capture the same independence structure, that is when they are Markov equivalent, even though defined differently; see Theorem 1 in Wermuth and Sadeghi (2012). In particular, each directed star graph is Markov equivalent to an undirected star graph with the same node and edge set, such as in Figure 1(b), since it does not contain a collision V: $\circ\mbox{$\hskip 0.50003pt\frac{\hskip 5.7361pt\hskip 5.7361pt}{\hskip 5.7361pt}\!\!\!\!\!\succ\!\hskip 1.07639pt$}\circ\mbox{$\hskip 0.50003pt\prec\!\!\!\!\!\frac{\hskip 5.7361pt\hskip 5.7361pt}{\hskip 5.7361pt}$}\circ$, that is two uncoupled arrows meeting head-on.

Sequences of regressions are traceable if one can use the graph alone to trace pathways of dependence, that is to decide when a non-vanishing dependence is induced for an uncoupled node pair. For this, each edge present in the regression graph corresponds to a non-vanishing, conditional dependence that is considered to be strong in a given substantive context. In evolutionary biology, the required changes in dependences that are strong enough to lead to mutations, were for instance called ‘drastic’ by Neyman (1971). In general, special properties of the generated distributions are required in addition; see Wermuth (2012).

In the last century, regression graphs were not yet defined and properties of traceable regression unknown. Then, distributions generated over directed star graphs, in which the root $L$ corresponds to a variable that is never directly observed, have been studied separately under different distributional assumptions and with changing main objectives. Many of them were named item response models. In these contexts, the unobserved root $L$ is also called latent or hidden, the items are the observed variables.

For instance, Spearman (1904) suggested to measure general intelligence with $Q$ similar quantitative tasks, a method now called factor analysis with a single latent variable. For confirmatory factor analyses, typically, a Gaussian distribution is assumed, after the observed variables are standardized to have mean zero and unit variance. Heywood (1931) and Anderson and Rubin (1956) derived necessary and sufficient conditions for the existence of one or more latent variables, without the constraint of non-vanishing, positive dependences of each item on the latent variables. The latter assumption arises however naturally when the items are to measure a specific latent ability or attitude or when they are the symptoms of a given disease.

For instance, when psychologists try to measure for children in a given age range, what is called the working memory capacity, then each item is the successful repetition – in reverse order – of a sequence of numbers. Typical tasks of the same difficulty are the sequences (3, 5, 7) and (2, 4, 6). The item difficulty increases, for instance, by presenting more numbers or numbers of two digits. For children with a higher capacity of the working memory, one expects more successes for tasks of the same difficulty as well as for tasks of increasing difficulty.

When instead, the leaves and the root of the star graph are both categorical, the resulting model is a latent class model, as proposed by Lazarsfeld (1950), again with an extensive follow-up literature. An important warning was given by Holland and Rosenbaum (1986): such a model can never be falsified when the number of levels of the latent variable is large compared to the number of cells in the observed contingency table. Expressed in other words, merely requiring conditional independence of categorical items given the latent variable imposes then no constraints on the observed distributions.

A general, testable constraint suggested by Holland and Rosenbaum (1986), that is now widely adopted in nonparametric item response theory, is to have a monotonically non-decreasing association of each item on $L$; see van der Ark (2012) or Mair and Hatzinger (2007). However, the underlying notion of ‘conditional association’, that had been proposed in probability theory, includes conditionally independent variables as being ‘conditionally associated’; see Esary, Proschan and Walkup (1967). By contrast, when one uses traceable regressions, each edge present in a graph excludes explicitly a corresponding conditional independence statement.

In this paper, we study similarities of joint Gaussian and of binary distributions generated over directed star graphs where all dependences of leaves, $1,\dots,Q$, on the root, $L$, are positive and the root is unobserved but does not coincide with any leaf or with any combination of the leaves. In Section 2, we summarize results for star graphs and for their traceable regressions. In Section 3, we describe joint Gaussian distributions, so generated, and in Section 4, we study joint binary distributions, especially correlation constraints on the distribution of the leaves. Section 5, gives applications to binary distributions, Section 6 a general discussion. The Appendix contains a technical proof.

A regression graph is said to be edge-minimal when each of its edges, that is present in the graph, indicates a non-vanishing dependence.

Thus, for traceable regressions with exclusively strong dependences, $i\pitchfork L$, in the generating star graph, each $ij$-edge in the induced complete covariance graph of the leaves, drawn as in Figure 2(a), indicates a non-vanishing dependence, $i\pitchfork j$. Probability distributions that do not induce such a dependence violate singleton-transitivity; see Table 1 in Wermuth (2012) for an example of such a family of distributions.

More consequences can be deduced, using strong dependences $i\pitchfork L$ and the Markov equivalence of the directed to the undirected star graph with the same node and edge set, as in Figures 1(a) and 1(b). For a traceable distribution with an undirected star graph, such as in Figure 1(b), each edge in the induced complete concentration graph, obtained by marginalising over the root $L$ and drawn as in Figure 2(b), indicates a non-vanishing conditional dependence of each pair of leaves given the remaining leaves, denoted by $i\pitchfork j|\{1,\dots,Q\}\setminus\{i,j\}$.

In exponential family distributions, a covariance graph corresponds to bivariate central moments, a concentration graph to joint canonical parameters and a directed star graph to regression parameters. If Markov equivalent models are also parameter equivalent, then their parameter sets are in a one-to-one relation. This property does not hold in general, but for instance in joint Gaussian distribution, and in joint binary distributions that are quadratic exponential. More generally, Markov equivalence and parameter equivalence coincide when only a single parameter determines for each variable pair whether the pair is conditionally dependent or independent.

For traceability of a given sequence of regressions, the generated family of distributions needs to have three properties of joint Gaussian distributions, that are not common to all probability distributions. In addition to the dependence-inducing property (singleton-transitivity) stated in Definition 1, these are the intersection (downward combination) and the composition property (upward combination of independences); see Lněnička and Matúš (2007) equations (9) to (10) and Wermuth (2012) section 2.4.

Pairwise independences are said to combine downwards and upwards if

respectively, for $c$ any subset of the remaining nodes, of $\{1,\dots Q\}\setminus\{i,j,k\}$. Both of these properties are already needed if one is using graphs of mixed edges just to decide on implied independences; see Sadeghi and Lauritzen (2014).

In the context of directed star graphs, these two properties are a consequence of the special type of generating process. After removing any two arrows for $Q\geq 5$, a directed star graph in $Q-2$ arrows remains. Thus, by introducing two additional pairwise independences in a directed star graph, these independence statement combine downwards. The Markov equivalence of the directed to the undirected star graph implies that for any two nodes, the statement $i\mbox{{ $\perp\hskip-9.90276pt\perp$ }}j|L$ can be modified to $i\mbox{{ $\perp\hskip-9.90276pt\perp$ }}j|Lc$ with $c\subseteq\{1,\dots Q\}\setminus\{i,j\}$ so that pairwise independences combine upwards.

In Gaussian distributions, each dependence is by definition linear and proportional to some (partial) correlation coefficient. Furthermore, there are no higher than two-factor interactions. In a traceable Gaussian distribution generated over a directed star graph, each directed edge indicates a simple, strong correlation, $\rho_{iL}$, called the loading of item $\bm{i}$ on $\bm{L}$. Then for each item pair $(i,j)$ a simple correlation, $\rho_{ij}$ is induced via

where $\rho_{ij|L}$ denotes the partial correlation coefficient

In factor analysis, the latent root $L$ is assumed, without loss of generality, to have mean zero and unit variance. Typically, when the observed variances of the items are nearly equal, the items are transformed to have mean zero and unit variance, so that their correlation matrix is analyzed.

In general, the model parameters are known to be identifiable for $Q>2$; see for instance Stanghellini (1997). For $Q=3$ items, the positive loadings $\rho_{1L},\rho_{2L},\rho_{3L}$ can be completely recovered using equations (3) for the positive, partial correlations $\rho_{ij|k}$ of each leaf pair. The maximum-likelihood equations (Lawley, 1967) reduce to the same type of three equations so that a unique, closed form solution can be obtained, provided it exists, and be written as $\hat{\bm{\lambda}}^{\mathrm{T}}=(\hat{\rho}_{1L}\hskip 3.50006pt\hat{\rho}_{2L}\hskip 3.50006pt\hat{\rho}_{3L})$ with

Clearly, these equations require for permissible estimates: $0<\hat{\rho}_{iL}<1$. In that case, there can be no zero and no negative marginal or partial item correlation that cannot be removed by recoding some items. We give in Table 1 three examples of $3\times 3$ invertible correlation matrices, showing marginal correlations in the lower half and partial correlations in the upper half. The first two examples have no permissible solution for $\hat{\bm{\lambda}}$ and illustrate so-called Heywood cases.

In the example on the left, $\rho_{23|1}=0$. By equation (4), the estimated loading of item 1 is then equal to one, i.e. $\hat{\rho}_{1L}=1$, so that item 1 cannot be distinguished from the latent variable itself. For the example in the middle, $\hat{\rho}_{1L}=\sqrt{0.6\cdot 0.5/0.1}=\sqrt{3}$ is larger than one, hence leads to an infeasible solution for a correlation coefficient. Thus even for a positive, invertible item correlation matrix, there may not exist a generating process via positive loadings on a single latent variable. The example on the right, has a perfect fit for the vector of estimated positive loadings in equation (4) as $\hat{\bm{\lambda}}^{\mathrm{T}}=[0.9\hskip 6.99997pt0.8\hskip 6.99997pt0.7]$.

For $Q>3$ items, proper positive loadings in equation (3), that is $0<\rho_{iL}<1$ for all items $i$, lead directly to a positive correlation matrix of the items, denoted in this paper by $\bm{P}>0$, that is to exclusively positive correlations, $\rho_{ij}=\rho_{iL}\rho_{jL}$ for $i\neq j$ in $\{1,\dots,Q\}$, and to a tetrad structure. Vanishing tetrads were defined by Spearman and coauthors and nowadays a popular search algorithm for models with possibly several latent variables is named Tetrad. This algorithm is based on and extends work by Spirtes, Glymour and Scheines (1993, 2001).

Thus for any pair $i,j$, the ratio of its correlations to variables in row $h$ of $\bm{S}$ is the same as to variables in row $k$, or, equivalently, there are vanishing tetrads: $s_{ih}s_{jk}-s_{jh}s_{ik}=0$. Let $\bm{\lambda}^{\mathrm{T}}$ denote a row vector of proper positive loadings, $\bm{\Delta}$ a diagonal matrix of elements $1-\rho_{iL}^{2}$, and $\bm{\delta}^{\mathrm{T}}$ a row vector of elements $-\rho_{iL}/\{\sqrt{s}(1-\rho_{iL}^{2})\}$, where $s=1+\sum_{i}\rho_{iL}^{2}/(1-\rho_{iL}^{2})$, then, as proven in the Appendix, the correlation matrix of the leaves $\bm{P}$ and its inverse, the concentration matrix $P^{-1}$, are

Some important direct consequences of (6) are given next. Lemma 1 uses the notion of M(inkowski)-matrices that were named and studied by Ostrowski (1937, 1956) and discussed in connection to totally positive Gaussian distributions much later using the name MTP₂ distributions; see e.g. Karlin and Rinott (1983). General MTP₂ distributions are characterized by having no variable pair negatively associated given all remaining variables. We will use a more strict form of MTP₂ that also excludes any variable pair being conditionally independent given all remaining variables.

When we denote the elements of $\bm{P}^{-1}$ by $\rho^{ij}$, then by Lemma 1 and Definition 3, all the precisions, $\rho^{ii},$ are positive and all concentrations are negative, that is $\rho^{ij}<0$, for all $i\neq j$, if $\bm{P}^{-1}$ is a complete M-matrix. The following important properties of complete M-matrices result from Ostrowski (1956), Section 1.

Thus, if the concentration matrix $P^{-1}$ of observed Gaussian items has exclusively negative off-diagonal elements, then all concentrations in every marginal distribution of the items are negative as well and ${\bm{P}}>0$.

Thus, for a complete concentration matrix of Gaussian items, negative concentrations mean positive dependence and every item pair is positively dependent, no matter which conditioning set is chosen. Now, a known rank-one condition for the existence of a Gaussian item correlation matrix $\bm{P}^{*}\!$ of a single factor simplifies as follows.

Conditions $(i)$ to $(iii)$ in Proposition 1 involve dependences of the items on the latent variable $L$, that is the unobserved loadings. In some applications, prior knowledge may be so strong that a positive loading can be safely predicted for each selected item, but otherwise, these characterizations involve unknown parameters.

By contrast, conditions $(iv)$ and $(v)$ in Proposition 1 concern directly the distribution of only the observed items. Equivalence to the former conditions for existence become possible by the added special properties of the concentration matrix $({\bm{P}}^{*})^{-1}$ being a complete M-matrix with vanishing tetrads.

The following example shows that a positive, tetrad correlation matrix alone, does not assure that its inverse is a M-matrix. The example is another Heywood case, conditions $(iv)$ and $(v)$ of Proposition 1 are violated:

where the dot-notation indicates symmetric entries.

This correlation matrix cannot have been generated over a directed star graph with proper positive loadings. If these correlations were observed, one would get with equation $\eqref{estload}$ that $\hat{\rho}_{1L}>1$, that is not a permissible solution. By contrast, every ${\bm{P}}^{-1}$ that is a complete M-matrix with vanishing tetrads implies a positive tetrad correlation matrix.

When $Q$ binary items are mutually independent given a binary variable $L$ and the factorization of the joint probability in equation (2) cannot be further simplified, since each item has a strong dependence on $L$, then it is a traceable regression generated over a directed star graph. The reason is that binary distributions are, just like Gaussian distributions via equation (3), dependence-inducing, that is

The property assures for star graphs with the independences of equation (1), $i\mbox{{ $\perp\hskip-9.90276pt\perp$ }}j|L$, that $i\pitchfork j$ is implied for each item pair. In applications, strong dependences of each item $i$ on $L$ are needed to obtain relevant dependences for each leaf pair $i,j$.

The joint binary distribution generated over a directed star graph is quadratic exponential since the $Q$ largest cliques in this type of a decomposable graph contain just two nodes, an item and $L$. Expressed equivalently, in the log-linear model of an undirected star graph, as in Figure 1(b), the largest, non-vanishing log-linear interactions are positive 2-factor terms, $\alpha_{iL}$. These are canonical parameters in the generated binary quadratic exponential family. Marginalizing over $L$ in such distributions with $\alpha_{iL}>0$ gives a tetrad form for the canonical parameters in the observed item distribution; see Cox and Wermuth (1994), Section 3:

With $\alpha_{iL}>0$ for all $i$, the joint binary distributions of the items have exclusively positive dependences. In general, the parameters $\alpha_{iL}$ are identifiable for $Q\geq 3$, see Stanghellini and Vantaggi (2013). But equation (9) does not lead to a an explicit form of the induced bivariate dependences for the item pairs.

For this, we write for instance for any two items $A$, $B$, and $L$, each with levels $0$ or $1$,

as well as

There are three equivalent expressions of Pearson’s correlation coefficient for the binary items 1 and 2. To present these, we use $\kappa_{12}=\sqrt{\pi_{0++}\pi_{1++}\pi_{+0+}\pi_{+1+}}$ and abbreviate the operation of taking a determinant by $\det(.)$ :

Equation (10) shows that the correlation is given by the cross-product difference, that it is zero if and only if the odds-ratio, defined as the cross-product ratio, equals one and that a positive correlation is equivalent to a positive log-odds ratio. Equation (11) gives as numerator the covariance and as denominator the product of two standard deviations. This is the usual definition of Pearson’s correlation coefficient, here for binary variables, and it implies that $\rho_{12}$ together with the one-dimensional frequencies of items 1 and 2 give the counts in their $2\times 2$ table. Equation (12) leads best to our main new result.

From equations (2) and (12), we have for each source V of the star graph, that is for each configuration $A\mbox{$\hskip 0.50003pt\prec\!\!\!\!\!\frac{\hskip 5.7361pt\hskip 5.7361pt}{\hskip 5.7361pt}$}L\mbox{$\hskip 0.50003pt\frac{\hskip 5.7361pt\hskip 5.7361pt}{\hskip 5.7361pt}\!\!\!\!\!\succ\!\hskip 1.07639pt$}B$, a trivariate binary distribution with $A\mbox{{ $\perp\hskip-9.90276pt\perp$ }}B|L$ and

Hence, the existence of a correlation matrix $\bm{P}^{*}$ becomes relevant for binary variables.

This result corrects a claim in Cox and Wermuth (2002) that a tetrad condition does not show in correlations of binary items but only in their canonical parameters in the induced distribution for the $Q$ items: just as in Gaussian distributions, a positive, invertible tetrad correlation matrix, $\bm{P}$, is induced for binary items if the dependence of each item on the latent binary $L$ is positive, that is if $\rho_{iL}>0$ for all $i$.

Nevertheless, the directly relevant dependences for the induced, complete concentration graph model are the canonical parameters, the log-linear interaction terms that are functions of the odds-ratios. For instance, for a binary distribution generated over a star graph to have a general MTP₂ distribution, the condition is $\alpha_{iL}\geq 0$, while for a strictly positive subclass, the condition is $\alpha_{iL}>0$ for all items $i$.

Necessary and sufficient conditions, that a density of equation (2) has generated the observed distribution of $Q=3$ items, have been derived as nine inequality constraints on the probabilities of the item by Allman et al. (2013) without noting their relation to MTP₂ distributions. For $\pi_{ijk}>0$, the first three reduce to

so that for each leaf, the odds for level 1 to 0 when the levels of the other two leaves match at level 0 do not exceed those with matches at level 1. Their last six constraints require nonnegative odds-ratios for each leaf pair, that is a binary MTP₂ distribution.

It can be shown that the above three inequalities are satisfied, whenever each leaf pair $i,j$ has a positive conditional dependence given the $Q-2$ remaining leaves, that is if the leaves have a strictly positive distribution. Therefore, this strict form of a MTP₂ binary distribution of the leaves is also sufficient for just $Q=3$ leaves to have been generated over a star graph with positive dependences of the leaves on the root. This implies but is not equivalent to ${\bm{P}}^{-1}$ having exclusively negative off-diagonal elements.

More complex characterizing inequality constraints on probabilities of the leaves for $Q>3$ categorical variables, when leaves and root have the same number of levels, are due to Zwiernik and Smith (2011). It remains to be seen how they simplify for binary MTP₂ distributions of the leaves or with a complete, tetrad ${\bm{P}}^{-1}$.

We use here three sets of binary items. The first is a medical data set, the last two are psychometric data sets where the questions were chosen to expect strong positive dependences of each item on a latent variable. As discussed above, to check conditions for the existence of a single latent variable, that might have generated the observed item dependences, we use here mainly the observed item correlation matrices and the observed marginal tables of all item triples. In the first two cases, no violations are detected. In the third case, item correlation matrices alone provide already enough evidence against the hypothesized generating process. Algorithms to compute maximum-likelihood estimates for latent class models, are widely available; see for instance Linzer and Lewis (2011).

It was known already to Bartlett (1951), that an invertible tetrad correlation matrix implies a tetrad concentration matrix. His proof is in terms of the general form of the inverse of sums of matrices. It is more direct to give the overall concentration matrix in explicit form.

For this, we use the partial inversion operator of Wermuth, Wiedenbeck and Cox (2006), described in the context of Gaussian parameter matrices in Marchetti and Wermuth (2009). It can be viewed as a Gaussian elimination technique (for some history see Grcar, 2011) and as a minor extension of the sweep operator for symmetric matrices discussed by Dempster (1972). This extension is to invertible, square matrices so that an operation is undone by just reapplying the operator to the same set.

Let $\bm{M}$ be a square matrix of dimension $d$ for which all principal submatrices are invertible. To describe partial inversions on $d$, we partition $\bm{M}$ into a matrix $\bm{m}$ of dimension $d-1$, column vector $\bm{v}$, row vector $\bm{w}^{\mathrm{T}}$ and scalar $s$

The transformation of $\bm{m}$ is also known as the vector form of a Schur complement; see Schur (1917). For $\bm{M}$ a covariance matrix, ${\bm{v}}/s$ contains $d-1$ linear, least-squares regression coefficients and $\bm{m}-\bm{vw}^{\mathrm{T}}/s$ is a residual covariance matrix.

For partial inversion on a set $a\subseteq\{1,\ldots,d\}$, one may conceptually apply equation (16) repeatedly for each index $k$ of $a$: one first reorders the matrix $\bm{M}$ so that $k$ is the last row and column, applies equation (16) and returns to the original ordering. Several useful and nice properties of this operator have been derived, such as commutativity and symmetric difference.

For $Q=3$ leaves and the correlation matrices of a star graph models, more detail is

where the $.\hskip 3.50006pt$-notation indicates an entry that is symmetric, the $\scriptstyle\sim\hskip 2.8681pt$-notation an entry that is symmetric up to the sign.

For $Q$ items, a single root $L$ and $\bm{\Psi}$ denoting their joint correlation matrix, we have for instance

For $Q>3$, the structure of these matrices is preserved in the sense that there is a diagonal matrix $\bm{\Delta}$ containing $1-\rho_{iL}^{2}$ as elements, a row vector $\bm{\lambda}^{\mathrm{T}}$ with loadings $\rho_{iL}$, the precision of $L$ as $s=1+\sum_{i}\rho_{iL}^{2}/(1-\rho_{iL}^{2})$, and a row vector $\bm{\delta}^{\mathrm{T}}$ with elements $-\rho_{iL}/\{\sqrt{s}(1-\rho_{iL}^{2})\}$.

For $\bm{P}$ the correlation matrix of the $Q$ items that are uncorrelated given $L$, one gets $\bm{P}$ as the submatrix of rows and columns $\{1,\ldots,Q\}$ of

and $\bm{P}^{-1}$ as the submatrix of rows and columns $\{1,\ldots,Q\}$ of $\mathrm{inv}_{L}\,\bm{\Psi}^{-1}$ so that

Thus, for $\bm{P}^{-1}$ a complete M-matrix with vanishing tetrads, $\bm{P}>0$ has tetrad form.