Consistent Bayesian causal discovery for structural equation models with equal error variances

1. Problem Statement

The paper addresses the problem of causal discovery: recovering the true underlying Directed Acyclic Graph (DAG) representing causal relationships among a set of variables from purely observational data.

• Model: The data is generated by a linear, recursive Structural Equation Model (SEM).
• Key Assumption: The error terms are independent, have equal variances ($\sigma^2$), but are not necessarily Gaussian.
• Challenge: While it is known that the true DAG is uniquely identifiable under the equal error variance assumption (unlike the general case where only the Markov equivalence class is identifiable), existing Bayesian methods lack theoretical guarantees (specifically, posterior consistency) for this setting. Most Bayesian approaches rely on Gaussian assumptions or require complex non-Gaussian priors that are difficult to analyze asymptotically under model misspecification.

2. Methodology

The authors propose a Bayesian DAG selection method that utilizes a "working model" which may differ from the true data-generating process but yields consistent results.

A. The Working Model

Instead of modeling the true non-Gaussian errors, the authors assume a Gaussian SEM with equal error variances as the working model:
$$X_j = X_{pa_\gamma(j)}^T b_j^\gamma + e_j^\gamma, \quad e_j^\gamma \sim \mathcal{N}(0, \theta^\gamma)$$
where $b_j^\gamma$ are the regression coefficients and $\theta^\gamma$ is the common error variance.

B. Prior Specification

• Coefficients ($b_j^\gamma$): Independent $g$-priors are placed on the SEM coefficients for each node. Specifically, $b_j^\gamma | \theta^\gamma, D_{j,n} \sim \mathcal{N}(0, g\theta^\gamma (D_{j,n}^T D_{j,n})^{-1})$.
• Variance ($\theta^\gamma$): A non-informative Jeffreys prior is used: $\pi(\theta^\gamma) \propto 1/\theta^\gamma$.
• DAG Prior ($\pi(\gamma)$): A generic prior over the space of DAGs ($\Gamma_p$), requiring only that it assigns strictly positive mass to all DAGs (e.g., a uniform prior).

C. Marginal Likelihood

By integrating out the coefficients and variance, the authors derive a closed-form expression for the marginal likelihood (evidence) of a DAG $\gamma$:
$$m(D_n | \gamma) \propto (V_n + g R_n^\gamma)^{-np/2} (1+g)^{(np - |\gamma|)/2}$$
where:

• $V_n$ is a data-dependent term related to the total variance.
• $R_n^\gamma$ is the sum of the residual sum of squares (RSS) from regressing each variable on its parents in $\gamma$.
• $|\gamma|$ is the number of edges in the DAG.

3. Key Theoretical Contributions

Theorem 1: Characterization of Minimum Error

The paper provides an alternative proof for a critical property regarding the sum of minimum expected squared errors ($r_\gamma$).

• Result: For any DAG $\gamma$, let $r_\gamma$ be the sum of the minimum expected squared errors when predicting each variable using its parents in $\gamma$.
• Property: $r_\gamma \geq r_{\gamma^*}$ (where $\gamma^*$ is the true DAG).
• Equality Condition: $r_\gamma = r_{\gamma^*}$ if and only if $\gamma$ is a supergraph of the true DAG ($\gamma \supseteq \gamma^*$).
• Proof Technique: The proof utilizes the Cholesky factorization of the covariance matrix. It shows that the diagonal entries of the Cholesky factor correspond to the minimum prediction errors, and the product of these entries is invariant under permutation, leading to the minimization result via the AM-GM inequality.

Theorem 2: Posterior Consistency

The main theoretical result establishes that the proposed Bayesian method consistently selects the true DAG as the sample size $n \to \infty$.

• Condition: The scaling parameter in the $g$-prior is set to $g = n$.
• Result: The posterior probability of the true DAG $\gamma^*$ converges to 1 in probability:
  $$1 - \pi(\gamma^* | D_n) \leq O_p(n^{-1/2}) + \exp\left(-\frac{np}{2}(\delta^* + o_p(1))\right)$$
  where $\delta^* > 0$ represents the gap in error variance between the true DAG and any non-supergraph.
• Significance: This holds even if the true errors are non-Gaussian. The method is robust to model misspecification regarding the error distribution, relying only on the equal variance assumption.

4. Simulation Study

The authors conducted a simulation study to validate the theoretical findings:

• Setup: A 4-variable DAG was used with non-Gaussian errors (Laplace, Normal, and t-distribution), all with equal variance.
• Procedure: Data was generated for increasing sample sizes ($n = 1600, 3200, 6400, 12800$). The posterior probability of the true DAG was computed for 100 replications at each size.
• Outcome: As $n$ increased, the posterior probability of the true DAG concentrated sharply around 1, confirming the consistency established in Theorem 2.

5. Significance and Impact

1. Bridging Theory and Practice: This work fills a gap in the literature by providing a theoretically principled Bayesian framework for causal discovery under equal error variances. Previous Bayesian methods often lacked consistency guarantees in non-Gaussian settings.
2. Robustness: The method demonstrates that assuming a Gaussian working model is sufficient for consistency, even when the true data is non-Gaussian, provided the equal variance assumption holds. This simplifies implementation significantly.
3. Computational Efficiency: The derivation of a closed-form marginal likelihood allows for efficient computation of posterior probabilities without the need for expensive MCMC sampling, which is often required in general Bayesian structure learning.
4. Foundational Insight: The alternative proof of Theorem 1 using Cholesky factorization offers a new perspective on why supergraphs minimize the sum of squared errors, reinforcing the geometric interpretation of causal discovery.

Conclusion

The paper presents a robust, consistent, and computationally efficient Bayesian approach for discovering causal structures in linear SEMs with equal error variances. By leveraging $g$-priors and a Gaussian working model, the authors prove that the true DAG can be recovered with probability approaching one, regardless of the underlying error distribution, provided the errors are independent and homoscedastic.