Source-linked AI summary

Causal Discovery Under Hard Selection Bias: A New Robust Score-Matching Approach

Yiwen Qiu, Francesco Montagna, Shimeng Huang, Francesco Locatello

arXiv:2610.18075v1cs.LG

TL;DR

Selection bias can distort observational data, and a general approach for causal discovery from observational data remains missing. The paper proposes a score-based method for hard selection under ANMs, using properties of the score to recover causal structure. It reports robustness across ANM settings except linear Gaussian models, while noting that severe truncation makes score estimation difficult.

  • Problem

    A general method for recovering the original causal graph from observational data under selection bias is missing, despite existing local or interventional approaches.

    Source: §1 · p. 1, Problem Setting: Hard Selection and… · p. 3
  • Method

    The proposed method uses the conditional rank of the Hessian log density in a score-based approach to causal discovery under hard selection.

    Source: §1 · p. 2
  • Results

    CHess is described as robust to hard selection across ANM settings except linear Gaussian, and experiments cover synthetic and real-world datasets.

    Source: §1 · p. 2, §7 · p. 10
  • Takeaways & Limitations

    The method is presented as a hard-selection approach that does not rely on external information about the selection process.

    Source: §1 · p. 2

    The method assumes hard selection, and severe truncation makes the score difficult to estimate in practice.

    Source: §7 · p. 10

Abstract

from arXiv · show

Selection bias is a difficult yet widespread problem in causal discovery, occurring whenever non-random selection processes lead to data that is not representative of the underlying populations. Due to its practical importance, related works have been proposed to address the problem under local [Versteeg et al., 2022] or interventional settings [Dai et al., 2025a], but general algorithms for observational data remain elusive. In this paper, we focus on truncated data (hard selection), and provide a general result for Additive Noise Models (ANM). We show that score-based methods are a surprisingly suitable candidate due to a special property of the score function ($\nabla \log p(x)$): deterministic truncation preserves both the score of the joint distribution as well as the Jacobian of the score. While established score-matching algorithms still fail in most settings, we actively leverage this observation to propose a new method for causal discovery that is robust under truncated data with ANMs. Overall, our approach offers a general, robust solution that is agnostic to external information about the selection process and still achieves comparable performance to the state-of-the-art approaches when selection is not present.

1 Introduction

Selection bias distorts observational data, while existing causal-discovery approaches for it are limited to settings such as interventions or local structure. This paper studies hard selection and introduces a score-based method intended to recover causal structure across ANM settings without external information about selection.

1 Introduction

  • Existing approaches address selection bias through interventional data or local structures, leaving a general observational-data approach missing.FCI and variants recover a graph only up to a Markov equivalence class in the stated setting.Source: §1 · p. 1
  • Hard selection truncates the observed distribution so that some parts of the true support are never observed.The paper distinguishes this setting from selection processes that merely yield nonrepresentative data.Source: §1 · pp. 1–2
  • The paper proposes a score-based causal-discovery method because causal asymmetry persists in the data score under hard selection, regardless of selection location.The method uses the conditional rank of the Hessian log density and is described as robust without external information about selection.Source: §1 · p. 2
  • SCORE, DAS, and NoGAM can be robust to hard selection, but failures persist and motivate a new method.The paper reports that SCORE is robust for nonlinear additive Gaussian noise, whereas more general extensions are not.Source: §1 · p. 2
  • CHess is the only method identified as resilient to hard selection across the listed ANM settings, except the non-identifiable linear Gaussian case.SCORE handles some settings, including latent hard selection, while other methods have more restricted guarantees.Source: §1 · p. 2, Table 1 · p. 2

2 Preliminaries and Problem Setup

The setup asks whether causal structure under the original distribution can be recovered from biased observational data under hard selection. The paper frames this problem for restricted additive noise models and assumes standard DAG conditions.

Problem Setting: Hard Selection and What We Know about it

  • The model is a restricted ANM in which each variable is generated from its parents through a C2 mechanism plus independent noise with strictly positive density.The associated graph is a DAG, and the model satisfies the restricted ANM identifiability condition.Source: Problem Setting: Hard Selection and… · p. 2
  • The paper assumes the DAG and distribution satisfy Markov and faithfulness conditions.These assumptions are used in the remainder of the paper.Source: Problem Setting: Hard Selection and… · p. 2
  • Hard selection is deterministic: observations are retained only when they fall inside a subset of the variables’ support.The observed distribution is p(X | S = 1), which differs from p(X).Source: Problem Setting: Hard Selection and… · pp. 2–3
  • The target is to recover the original-distribution causal graph from observational data drawn under selection, a general setting for which the paper says no method had been developed.Related work covered local, interventional, and latent-variable settings.Source: Problem Setting: Hard Selection and… · p. 3
  • The motivating question is which properties of the original distribution survive hard selection and can support robust causal discovery.The authors revisit score functions and score-matching methods as a first step.Source: Problem Setting: Hard Selection and… · p. 3

3 The Blessing of Score Functions

Hard truncation preserves local score information in the observed region, but score-based causal tests do not all remain valid: variance and Hessian-based results survive while residual testing can fail.

3.1 Score Preservation

  • Within the truncation region's interior, both the score and Hessian of the log density are preserved at observed data points.The score-matching argument uses this pointwise preservation; the boundary is an exception.Source: §3.1 · p. 3
  • At the truncation boundary, an extra indicator-derivative term breaks score preservation, and the selected density violates the vanishing-boundary condition used by Stein-based estimation.The paper proposes diffusion score matching as more robust to boundary effects.Source: §3.1 · p. 3

3.2 Revisiting Score-Matching

  • Score-matching sink tests use a variance criterion for Gaussian noise and a regression-residual criterion for non-Gaussian noise.The residual for a true sink without selection equals its centered noise.Source: §3.2 · p. 4

3.3 When Score-Matching is Surprisingly Resilient and When it is Not

  • The proposed approach therefore uses the score directly, avoiding reliance on residuals or graphical assumptions to identify sinks once the skeleton is known.The motivation is that the score is local and preserved, whereas the residual test depends on density structure that selection can distort.Source: §3.3 · pp. 4–5
  • Under nonlinear ANMs with Gaussian noise, the variance test still identifies sink nodes from truncated data, and mixed-Hessian conditional-independence testing extends to skeleton discovery.The sink result follows from score preservation; the skeleton result is stated as an extension of the earlier Hessian test.Source: §3.3 · pp. 4–5
  • The residual test can fail when selection depends on both a sink and one of its parents, because the selected residual no longer equals the noise residual.When selection depends on a single node, the residual test is a special case that still holds.Source: §3.3 · p. 5

4 Causal Discovery Method Robust to Hard Selection

The method identifies sinks using the conditional rank of a Hessian log-density outer-product matrix. Its rank separates sinks from non-sinks across ANM settings except the linear Gaussian case.

4.1 Intuition and Definitions

  • After conditioning on neighbors, a sink has one independent noise direction in its score, while a non-sink has at least two.The conditional rank of the variance of the score-gradient outer product indicates the number of independently varying factors.Source: §4.1 · p. 5, §4.2 · p. 6

4.2 Main Results: from Linear to Nonlinear ANM

  • For a restricted linear ANM with independent noise terms, a node is a sink if and only if its conditional matrix rank is 1.For generic z, sinks have rank 1 and non-sinks rank at least 2.Source: §4.2 · p. 7
  • In nonlinear ANMs, mechanism curvature generally adds a linearly independent direction, preserving a rank gap of at least one between sinks and non-sinks.The result is stated for generic z.Source: §4.2 · p. 7
  • The conditional rank test separates sink and non-sink nodes across ANM settings except linear Gaussian, where all nodes have rank 0.The nonlinear rank gap is at least one for generic z.Source: §4.2 · p. 7
  • Without knowing the ANM or noise type, the procedure iteratively selects the node with minimum rank when a gap exists.The conditional test avoids the stronger structural assumptions required by the unconditional rank test.Source: §4.2 · pp. 7–8

5 Estimation

The estimation procedure estimates the Hessian log-density, estimates a conditional variance matrix, and uses its rank to identify sinks. It considers Stein and diffusion-based estimators and a residual proxy for the conditional matrix.

Estimation

  • The procedure estimates the Hessian log-density, estimates a conditional variance matrix, and computes its rank to identify sinks.The passage lists four steps, but the fourth is not specified in the supplied text.Source: §5 · p. 8

Score and Hessian Estimation

  • The paper applies Stein and diffusion-based score estimators; the diffusion estimator does not require a strictly positive density, unlike Stein and KDE estimators.Stein and KDE estimators suffer when p(x) is approximately zero.Source: §5 · p. 8

Practical Estimation of M cond

  • In practice, a regression-residual outer product estimates the conditional matrix and follows a similar approach to NoGAM’s residual-based noise estimate.The residual subtracts an estimate based on the node’s neighbors from its Hessian row.Source: §5 · p. 8

6 Experiments

On synthetic data, the method is competitive without selection and generally outperforms baselines under hard selection, with a noted linear-nG exception. It also outperforms baselines on the Sachs dataset, whose aggregated perturbation data warrants careful interpretation.

6.1 Synthetic Experiments

  • Synthetic experiments use Erdős–Rényi DAGs with 5 or 10 nodes across linear-non-Gaussian, nonlinear-Gaussian, and nonlinear-non-Gaussian regimes.The selection mechanism depends on all observed variables, and selected samples are held at n=1000.Source: §6.1 · p. 8

6 Experiments

  • Without selection, the method is competitive across synthetic SEM regimes; with hard selection, it outperforms baselines except that LiNGAM sometimes wins in linear-nG.On the Sachs dataset, it outperforms baselines both with and without selection, but aggregated perturbation data warrants care in interpreting that result.Source: §6.1 · p. 9, §6.2 · p. 9

6.2 Real-world experiments

  • The Sachs protein-signalling dataset contains 7,466 samples over 11 proteins; selection leaves 373 samples, matched by downsampling the unselected setting.The dataset has a verified ground-truth DAG.Source: §6.2 · p. 9

7 Conclusion

The paper concludes that its score-based method identifies causal relationships under hard selection and reports experiments on synthetic and real-world datasets. It assumes hard selection and notes that severe truncation makes score estimation difficult.

7 Conclusion

  • The method leverages score preservation under hard selection to identify causal relationships, with effectiveness demonstrated on synthetic and real-world datasets.The conclusion frames the work as an approach to hard-selection bias in causal discovery.Source: §7 · pp. 9–10
  • The method assumes hard selection, so soft, non-deterministic selection remains outside its stated scope.The authors identify soft-selection extensions as future work.Source: §7 · p. 10
  • Severe truncation makes the score difficult to estimate in practice, despite the theoretical preservation result.The authors call for more robust score-matching estimation techniques.Source: §7 · p. 10

A.1 Definitions

  • The paper assumes the DAG and distribution satisfy Markov and faithfulness conditions.Faithfulness excludes independencies not entailed by the Markov condition.Source: §A.1 · p. 12

A.2 Illustrative Example: Linear Chain X1 →X2 →X3

In the linear chain example, conditional Hessian rank separates the sink from non-sinks. The score-based conditional-independence criterion also preserves skeleton discovery under deterministic truncation.

A.2 Illustrative Example: Linear Chain X1 →X2 →X3

  • For non-Gaussian noise, distinct functions of the noise contribute to the conditional matrix, producing rank 2 for non-sink X2.The passage describes these functions as generally independent.Source: §A.2 · p. 13
  • In the chain, sink X3 is uniquely identified by its conditional rank, while non-sinks X1 and X2 have rank 2.The example assumes non-Gaussian noise.Source: §A.2 · p. 13

B.1 Preservation of Skeleton Discovery

  • Under deterministic truncation, the mixed second derivative of the conditional log density is zero exactly when Xi and Xj are conditionally independent given the remaining variables in G.This result supports preservation of skeleton discovery under hard selection.Source: §B.1 · p. 13

B.2 Proofs of Conditional Hessian Rank Characterization

The conditional Hessian rank distinguishes sinks from non-sinks across the stated ANM settings, except linear Gaussian models. Conditioning on neighbors isolates the contributions that underpin this rank gap.

B.2 Proofs of Conditional Hessian Rank Characterization

  • The rank gap fails in the linear Gaussian case because all nodes have rank 0.The characterization excludes degenerate cases where graph structure is not faithfully reflected in rank.Source: §B.2 · p. 13
  • For linear ANMs, a sink has conditional rank 1, while a non-sink has rank at least 2, giving a sink-identification criterion for generic z.The proposition concerns restricted linear ANMs with independent noise terms.Source: §B.2 · p. 14
  • For nonlinear ANMs, a sink has rank 1 with non-Gaussian noise and 0 with Gaussian noise; non-sinks have rank at least 2 and 1, respectively.The criterion identifies sinks at the corresponding ranks under the stated noise types.Source: §B.2 · p. 15
  • Conditioning on neighbors fixes the relevant parent and child values, leaving the sink’s exogenous noise as the remaining random quantity in its conditional Hessian row.For a sink, the resulting rank is 1 with non-Gaussian noise and 0 with Gaussian noise.Source: §B.2 · p. 15
  • For non-sinks with non-Gaussian noise, own-noise and child-noise contributions are generally functionally distinct, yielding conditional rank at least 2 for almost every z.The child mechanism contributes additional curvature, and the restricted ANM identifiability condition is invoked.Source: §B.2 · p. 16

B.3 Proofs of Unconditional Hessian Rank Test

The unconditional Hessian rank test can identify sinks under a graph-degree condition, but it is less generally applicable than the conditional test. The section also describes estimation properties and how CHess uses conditional ranks.

B.3 Proofs of Unconditional Hessian Rank Test

  • Without conditioning on neighbors, upstream noise contributions make unconditional rank depend on both parent and child counts, so a sink and non-sink can share rank.In the chain example, both X1 and X3 have unconditional rank 2.Source: §B.3 · p. 16
  • For linear mechanisms, sinks have rank 1 and non-sinks rank at least 2; for nonlinear mechanisms, sink rank is at most |PA_l| + 1 and non-sink rank is generically greater than |PA_i| + 1.The nonlinear non-sink rank is also bounded above by |PA_i| + |CH_i| + 1.Source: §B.3 · p. 16
  • An unconditional rank gap identifies a sink only when every non-sink satisfies |PA_i| + |CH_i| − |PA_l| ≥ 1.Tree-like graphs with non-sinks having multiple children and sinks having at most one parent satisfy this condition, but more general graphs may not.Source: §B.3 · pp. 16–17
  • Unlike the unconditional test, the conditional test does not require the graph-degree constraint because conditioning on neighbors blocks parent-noise contributions.It therefore counts noise terms that directly feed into the score of Xi.Source: §B.3 · p. 17

B.4 Proofs of Consistent Estimator of M cond

  • The residual-based estimator is consistent for the conditional matrix when its regression estimate of the conditional mean is consistent.The passage states this result for linear non-Gaussian or nonlinear Gaussian ANMs.Source: §B.4 · p. 17

C Related Works

  • The paper targets full causal-structure learning among observed variables under selection bias without requiring interventional data or assumptions on the selection mechanism.It allows nonlinear and non-Gaussian mechanisms.Source: C Related Works · p. 18

D.1 Main Algorithm

  • CHess estimates Hessians, recovers the skeleton, iteratively identifies and removes a minimum-rank sink, orients its adjacent edges as incoming, then prunes edges.The described procedure accepts data potentially drawn from a truncated distribution.Source: §D.1 · pp. 18–19

D.2 Algorithm for Conditional Hessian Matrix Estimation

  • Conditional rank estimation uses local averaging for nonlinear non-Gaussian ANMs and regression residuals for linear non-Gaussian or nonlinear Gaussian ANMs.These are the two alternative estimators described for sink identification.Source: §D.2 · p. 19

E.1 Datasets

Across the reported selection-bias experiments, the conditional Hessian rank method outperforms as selection becomes severe; on DREAM4 without selection, it performs best on four of five networks. The section also specifies the datasets, score-estimation approach, evaluation setup, and baselines.

E.1 Datasets

  • Synthetic configurations use 1,000 observations averaged over 20 seeds, with linear and nonlinear SEMs generated using the stated noise and mechanism setups.The linear weights and node noise variances are randomly drawn; nonlinear mechanisms use a one-layer sigmoid MLP.Source: §E.1 · pp. 19–20
  • The Sachs dataset contains 7,466 continuous measurements across nine interventional-context files and is log-transformed because its values are highly skewed.The dataset measures concentrations of nine proteins.Source: §E.1 · p. 20

E.2 Estimation Method

  • Diffusion-based score estimation trains a noise predictor whose output approximates the data score, then obtains Hessian entries by differentiating that output.The passage describes denoising score matching and a single backpropagation pass for each Hessian estimate.Source: §E.2 · p. 20

E.3 Evaluation matrics

  • Baseline outputs are converted from CPDAGs or PAGs to the comparison format, and evaluation metrics include precision and recall.The cited table specifies precision = TP/(TP+FP) and recall = TP/(TP+FN).Source: §E.3 · p. 21, Table 6 · p. 21

E.4 Hyperparameters

  • The reported rank-test significance level is 0.05, the KRR ridge parameter is 10^-5, and CAM pruning thresholds are 0.001 for DREAM4 and 0.01 for Sachs.The thresholds are among the method’s stated hyperparameters.Source: §E.4 · p. 21

E.5 Baseline Implementations

  • PC, GES, FCI, GRaSP, and LiNGAM use causallearn defaults, while CAM, SCORE, DAS, and NoGAM use dodiscover defaults.The passage identifies the implementation packages used for the baselines.Source: §E.5 · p. 21

F.1 All Baselines in Linear non-Gaussian ANM; F.2 All Baselines in MLP non-Gaussian ANM

  • Without hard selection bias, LiNGAM-based methods achieve best performance in one result, while another reports similar F1 around 0.7–0.8 across methods; under severe selection, the conditional Hessian rank-based method starts to outperform.Source: §F.1 · p. 21, §F.2 · p. 26

F.3 DREAM4 experiments

  • On DREAM4 without selection, the method achieves best SHD performance on four of five networks; DAS is best on Network 5.The comparison uses full data because each network has only 30 observations.Source: §F.3 · p. 27

F.4 Sachs experiments

On Sachs, the proposed methods achieve best performance under hard selection in aggregated and per-condition evaluations. The experiments vary selection shapes and strengths and report the hardware used.

F.4 Sachs experiments

  • Across aggregated Sachs experiments, the proposed methods achieve best performance on arrow F1 and causal ordering under varied selection shapes and strengths.The results include half-space, annulus, and two-box selection.Source: §F.4 · p. 27, Table 18 · p. 27, Table 21 · p. 29
  • The aggregated dataset combines nine experimental conditions, and random downsampling at ρ = 1 matches the selected sample size for fair comparison.The dataset has n = 7466 before selection bias is applied.Source: §F.4 · p. 27
  • Per-condition experiments score each of nine datasets against its own ground truth; perfect intervention cuts the perturbed protein’s incoming edges.Each subdataset retains 50% of its data under half-space, annulus, or two-box selection.Source: §F.4 · pp. 27, 30
  • On per-condition Sachs data with hard selection, CHess achieves the best performance; CHess (Stein) dominates overall, while CHess (DiffAN) dominates in some conditions.The DiffAN result is reported for psitect, aktinhib, and pma, but is not shown in the averaged table.Source: §F.4 · p. 30, Table 24 · p. 30

G Hardware

  • Experiments ran on Debian GNU/Linux 12 with an Intel Xeon E5-2680 v4 CPU, 503 GB RAM, and ten NVIDIA GeForce GTX 1080 Ti GPUs.The implementation used Python 3.13.Source: G Hardware · p. 30
Loading 2610.18075v1…