Block-SSTS: Block-Wise Causal Discovery

• Block-SSTS is a block-wise variant of the Score-Schur Topological Sort that eliminates multiple leaf nodes at once using the score Jacobian information matrix derived from an unconstrained score model.
• It reduces sequential elimination steps by grouping nodes with similar diagonal energies, thereby averaging out noise and mitigating the accumulation of non-linear expectation gaps.
• The method enhances computational efficiency and structural accuracy, offering faster causal discovery in sparse and non-linear additive noise models compared to sequential approaches.

Searching arXiv for Block-SSTS and the associated SSTS paper. Block-SSTS is a block-wise variant of the Score-Schur Topological Sort (SSTS) introduced for causal discovery in additive-noise causal models (ANMs) over directed acyclic graphs (DAGs). It operates on the Score-Jacobian Information Matrix (SJIM), estimated from an unconstrained score model, and replaces sequential single-node Schur-complement eliminations with simultaneous elimination of a “leaf block” whose diagonal energies are within a relative tolerance of the minimum. In the formulation of "Optimization-Free Topological Sort for Causal Discovery via the Schur Complement of Score Jacobians" (Wu et al., 28 Apr 2026), this modification is motivated by non-linear expectation-gap effects and by the presence of parallel leaves in sparse DAGs; it is designed to reduce extraction depth, average out estimation noise, and bound structural error.

1. Position within score-based causal discovery

Block-SSTS is part of a decoupled causal-discovery paradigm that shifts the principal difficulty from constrained structural optimization to statistical score estimation. The underlying SSTS procedure separates density estimation from structure extraction: one first trains an unconstrained score model $s_\theta(x)$ so that $s_\theta \approx \nabla \log p$, then forms an empirical SJIM, and finally extracts a topological order by repeated Schur-complement elimination of leaf variables (Wu et al., 28 Apr 2026).

The central object is the Score-Jacobian Information Matrix. For a continuous random vector $X \in \mathbb{R}^d$ generated by an ANM with structural equations

$$x_i = f_i(x_{pa(i)}) + \epsilon_i, \qquad \epsilon_i \sim N(0,\sigma^2),$$

the score is $s(x)=\nabla_x \log p(x)$, the sample Hessian is $H(x)=-\nabla_x^2 \log p(x)$, and the SJIM is the population Fisher information with respect to $x$,

$$I \coloneqq E_{p(x)}[H(x)] \in \mathbb{R}^{d \times d}.$$

Under the ANM,

$$I_{ii} = (1/\sigma^2) + (1/\sigma^2)\sum_{j \in ch(i)} E\!\left[\left(\partial f_j/\partial x_i\right)^2\right].$$

A topological sink, equivalently a leaf node, satisfies $ch(\ell)=\varnothing$, and therefore $s_\theta \approx \nabla \log p$0, which is minimal among all diagonals (Theorem 1) (Wu et al., 28 Apr 2026).

This diagonal-minimum property is what makes Schur-based topological sorting possible. A plausible implication is that Block-SSTS should be read not as a separate estimator of causal score geometry, but as a structural extraction rule layered on top of the SJIM estimated by the score model.

2. Schur-complement marginalization and the origin of block elimination

In the linear Gaussian ANM, with strictly upper-triangular adjacency matrix $s_\theta \approx \nabla \log p$1 and $s_\theta \approx \nabla \log p$2, the joint precision matrix equals

$s_\theta \approx \nabla \log p$3

Eliminating a leaf node $s_\theta \approx \nabla \log p$4 is exactly equivalent to taking the Schur complement of $s_\theta \approx \nabla \log p$5 with respect to row and column $s_\theta \approx \nabla \log p$6:

$s_\theta \approx \nabla \log p$7

where $s_\theta \approx \nabla \log p$8. In block form, if

$s_\theta \approx \nabla \log p$9

then eliminating block $X \in \mathbb{R}^d$0 yields

$X \in \mathbb{R}^d$1

Under linear ANMs, recursive application of this Schur complement exactly mimics sequential leaf removal without retraining (Theorem 2) (Wu et al., 28 Apr 2026).

Block-SSTS arises when this exact equivalence is perturbed by non-linearity. In non-linear ANMs, $X \in \mathbb{R}^d$2 varies with $X \in \mathbb{R}^d$3, so

$X \in \mathbb{R}^d$4

For elimination of leaf $X \in \mathbb{R}^d$5, the expectation gap is

$X \in \mathbb{R}^d$6

By Remark 3, this gap is localized to the parent-block of $X \in \mathbb{R}^d$7 but cascades through $X \in \mathbb{R}^d$8 sequential eliminations, degrading structural fidelity (Wu et al., 28 Apr 2026).

The motivation for block-wise elimination is therefore twofold. First, sparse DAGs may contain multiple true leaves in parallel, and their theoretical $X \in \mathbb{R}^d$9 values are identical. Second, empirical noise perturbs those diagonals, so eliminating them one-by-one both overreacts to sampling fluctuations and preserves the full depth-$$x_i = f_i(x_{pa(i)}) + \epsilon_i, \qquad \epsilon_i \sim N(0,\sigma^2),$$0 cascade through which $$x_i = f_i(x_{pa(i)}) + \epsilon_i, \qquad \epsilon_i \sim N(0,\sigma^2),$$1 accumulates.

3. Algorithmic specification of Block-SSTS

Block-SSTS groups together variables whose diagonal energies are sufficiently close to the minimum and eliminates them in a single block Schur complement. The required inputs are a score model $$x_i = f_i(x_{pa(i)}) + \epsilon_i, \qquad \epsilon_i \sim N(0,\sigma^2),$$2, data $$x_i = f_i(x_{pa(i)}) + \epsilon_i, \qquad \epsilon_i \sim N(0,\sigma^2),$$3, a ridge parameter $$x_i = f_i(x_{pa(i)}) + \epsilon_i, \qquad \epsilon_i \sim N(0,\sigma^2),$$4, and a block tolerance $$x_i = f_i(x_{pa(i)}) + \epsilon_i, \qquad \epsilon_i \sim N(0,\sigma^2),$$5 (Wu et al., 28 Apr 2026).

The empirical construction begins with

$$x_i = f_i(x_{pa(i)}) + \epsilon_i, \qquad \epsilon_i \sim N(0,\sigma^2),$$6

followed by symmetrization,

$$x_i = f_i(x_{pa(i)}) + \epsilon_i, \qquad \epsilon_i \sim N(0,\sigma^2),$$7

The active set is initialized as $$x_i = f_i(x_{pa(i)}) + \epsilon_i, \qquad \epsilon_i \sim N(0,\sigma^2),$$8 and the order as $$x_i = f_i(x_{pa(i)}) + \epsilon_i, \qquad \epsilon_i \sim N(0,\sigma^2),$$9. At each iteration one computes

$s(x)=\nabla_x \log p(x)$0

forms the block

$s(x)=\nabla_x \log p(x)$1

prepends $s(x)=\nabla_x \log p(x)$2 (sorted internally) to $s(x)=\nabla_x \log p(x)$3, and, if $s(x)=\nabla_x \log p(x)$4, updates the reduced information matrix by

$s(x)=\nabla_x \log p(x)$5

where $s(x)=\nabla_x \log p(x)$6. The process repeats until $s(x)=\nabla_x \log p(x)$7, after which $s(x)=\nabla_x \log p(x)$8 is returned (Wu et al., 28 Apr 2026).

The paper attributes two effects to this block rule. It “averages out estimation noise among parallel leaves,” and it “reduces the number of sequential steps from $s(x)=\nabla_x \log p(x)$9 to roughly $H(x)=-\nabla_x^2 \log p(x)$0, truncating the accumulation of $H(x)=-\nabla_x^2 \log p(x)$1” (Wu et al., 28 Apr 2026). This suggests that the method is best understood as a compression of the elimination trajectory rather than merely a batching heuristic.

A common misunderstanding is to read the term “optimization-free” as applying to the entire learning pipeline. The stated procedure is narrower: it bypasses constrained structure optimization, while retaining unconstrained score-model training as Step A of the overall workflow (Wu et al., 28 Apr 2026).

4. Block-wise expectation gap and structural guarantees

For a block $H(x)=-\nabla_x^2 \log p(x)$2 of simultaneous leaves, the paper defines the block-Schur expectation gap

$H(x)=-\nabla_x^2 \log p(x)$3

and states that it decomposes into a sum of covariances among the gradients of $H(x)=-\nabla_x^2 \log p(x)$4 for $H(x)=-\nabla_x^2 \log p(x)$5 (Wu et al., 28 Apr 2026). By aggregating nodes, the norm $H(x)=-\nabla_x^2 \log p(x)$6 shrinks, and the total number of Schur steps falls from $H(x)=-\nabla_x^2 \log p(x)$7 to approximately $H(x)=-\nabla_x^2 \log p(x)$8.

The formal guarantees are given under an assumption that the entries of $H(x)=-\nabla_x^2 \log p(x)$9 are sub-Gaussian and the noise is homoscedastic. Let $x$0, let $x$1 denote the minimal block size, and let $x$2 denote the dimension. The SJIM concentration lemma states that, with probability at least $x$3,

$x$4

The block-Schur perturbation proposition states that if $x$5 and every block Schur complement is on at most a $x$6 submatrix, then each block elimination incurs a perturbation in the reduced matrix bounded by $x$7, and over at most $x$8 blocks the total deviation from the population sequence is $x$9 (Wu et al., 28 Apr 2026).

The structural error theorem then states that there exists $$I \coloneqq E_{p(x)}[H(x)] \in \mathbb{R}^{d \times d}.$$0 such that, with probability at least $$I \coloneqq E_{p(x)}[H(x)] \in \mathbb{R}^{d \times d}.$$1, the recovered adjacency $$I \coloneqq E_{p(x)}[H(x)] \in \mathbb{R}^{d \times d}.$$2 after Block-SSTS plus downstream pruning satisfies

$$I \coloneqq E_{p(x)}[H(x)] \in \mathbb{R}^{d \times d}.$$3

where

$$I \coloneqq E_{p(x)}[H(x)] \in \mathbb{R}^{d \times d}.$$4

is the per-leaf expectation gap (Wu et al., 28 Apr 2026). In regimes where $$I \coloneqq E_{p(x)}[H(x)] \in \mathbb{R}^{d \times d}.$$5 and nonlinearities are mild, so that $$I \coloneqq E_{p(x)}[H(x)] \in \mathbb{R}^{d \times d}.$$6 is small, exact recovery is ensured.

The decomposition of error into a statistical estimation term and a non-linear expectation-gap term is central. The paper identifies these as two irreducible terms: the statistical estimation error $$I \coloneqq E_{p(x)}[H(x)] \in \mathbb{R}^{d \times d}.$$7 and the non-linear expectation gap $$I \coloneqq E_{p(x)}[H(x)] \in \mathbb{R}^{d \times d}.$$8 (Wu et al., 28 Apr 2026).

5. Complexity, empirical scaling, and parameterization

The computational motivation for Block-SSTS is explicit. Full SSTS, which eliminates one node at a time, requires $$I \coloneqq E_{p(x)}[H(x)] \in \mathbb{R}^{d \times d}.$$9 Schur complements, each costing $$I_{ii} = (1/\sigma^2) + (1/\sigma^2)\sum_{j \in ch(i)} E\!\left[\left(\partial f_j/\partial x_i\right)^2\right].$$0, for total complexity $$I_{ii} = (1/\sigma^2) + (1/\sigma^2)\sum_{j \in ch(i)} E\!\left[\left(\partial f_j/\partial x_i\right)^2\right].$$1. Block-SSTS reduces this to roughly $$I_{ii} = (1/\sigma^2) + (1/\sigma^2)\sum_{j \in ch(i)} E\!\left[\left(\partial f_j/\partial x_i\right)^2\right].$$2 Schur complements times $$I_{ii} = (1/\sigma^2) + (1/\sigma^2)\sum_{j \in ch(i)} E\!\left[\left(\partial f_j/\partial x_i\right)^2\right].$$3, giving $$I_{ii} = (1/\sigma^2) + (1/\sigma^2)\sum_{j \in ch(i)} E\!\left[\left(\partial f_j/\partial x_i\right)^2\right].$$4 (Wu et al., 28 Apr 2026). For block size $$I_{ii} = (1/\sigma^2) + (1/\sigma^2)\sum_{j \in ch(i)} E\!\left[\left(\partial f_j/\partial x_i\right)^2\right].$$5 in sparse graphs, the paper states that one obtains an $$I_{ii} = (1/\sigma^2) + (1/\sigma^2)\sum_{j \in ch(i)} E\!\left[\left(\partial f_j/\partial x_i\right)^2\right].$$6 speedup.

The empirical SJIM construction requires $$I_{ii} = (1/\sigma^2) + (1/\sigma^2)\sum_{j \in ch(i)} E\!\left[\left(\partial f_j/\partial x_i\right)^2\right].$$7 time to accumulate Hessians, while streaming mini-batches bound memory to $$I_{ii} = (1/\sigma^2) + (1/\sigma^2)\sum_{j \in ch(i)} E\!\left[\left(\partial f_j/\partial x_i\right)^2\right].$$8. In experiments up to $$I_{ii} = (1/\sigma^2) + (1/\sigma^2)\sum_{j \in ch(i)} E\!\left[\left(\partial f_j/\partial x_i\right)^2\right].$$9 and sample size $ch(\ell)=\varnothing$0, Block-SSTS required less than $ch(\ell)=\varnothing$1 for the structural extraction phase, with $ch(\ell)=\varnothing$2 for $ch(\ell)=\varnothing$3 on a single GPU; this is contrasted with $ch(\ell)=\varnothing$4 for repeated acyclicity-penalty evaluations in NOTEARS-style methods (Wu et al., 28 Apr 2026).

The two principal tuning quantities are the block tolerance $ch(\ell)=\varnothing$5 and the effective block size $ch(\ell)=\varnothing$6. The paper states that $ch(\ell)=\varnothing$7 should be small enough to avoid grouping non-parallel nodes, which induces structural error, but large enough to merge true parallel leaves and suppress estimation noise; empirically, $ch(\ell)=\varnothing$8–$ch(\ell)=\varnothing$9 works well (Wu et al., 28 Apr 2026). The effective block size is adaptive: in early strata of very disconnected DAGs, blocks can be large, whereas in dense graphs $s_\theta \approx \nabla \log p$00 may shrink to $s_\theta \approx \nabla \log p$01.

The sample-size discussion is similarly concrete. To ensure concentration $s_\theta \approx \nabla \log p$02, one needs $s_\theta \approx \nabla \log p$03, and in practice $s_\theta \approx \nabla \log p$04 on the order of a few times $s_\theta \approx \nabla \log p$05 suffices for moderate $s_\theta \approx \nabla \log p$06 (Wu et al., 28 Apr 2026).

6. Exactness domain, failure modes, and open problems

Block-SSTS is exact for linear ANMs and robust to moderate nonlinearity, specifically “tanh, sigmoidal” nonlinearities (Wu et al., 28 Apr 2026). Its formal derivation depends on the homoscedastic Gaussian ANM framework and on the algebraic relation between graph marginalization and the Schur complement of the SJIM. In this sense, its strongest claims concern settings where the score geometry is well estimated and the expectation gap remains small.

The paper also states explicit limitations. In purely multiplicative noise models or post-nonlinear transforms, the ANM assumptions break, and the method may fail. To reduce $s_\theta \approx \nabla \log p$07, it suggests the exact “covariance patching” of Appendix C.2, but at the cost of $s_\theta \approx \nabla \log p$08 complexity (Wu et al., 28 Apr 2026). This frames a basic trade-off between algebraic correction of the non-linear gap and computational efficiency.

Several open questions are stated. These include whether block-wise covariance patching can be combined with low-rank updates to retain $s_\theta \approx \nabla \log p$09 memory while compensating $s_\theta \approx \nabla \log p$10; what the minimax rates are for score-based SJIM estimation in high dimensions under general ANMs; and how to extend exact algebraic marginalization to post-nonlinear or heteroscedastic noise families beyond the homoscedastic Gaussian ANM (Wu et al., 28 Apr 2026).

A broader interpretive point follows from the paper’s closing summary. By replacing the $s_\theta \approx \nabla \log p$11 sequential Schur complementarities of SSTS with blockwise eliminations, Block-SSTS is presented as both an acceleration mechanism and an error-control device in non-linear settings. This suggests a reframing of continuous causal discovery in which the principal obstacle is no longer non-convex acyclicity enforcement, but the finite-sample estimation variance of the global score geometry (Wu et al., 28 Apr 2026).