Estimation of single-cell and tissue perturbation effect in spatial transcriptomics via Spatial Causal Disentanglement
Authors: Stathis Megas
Daniel G. Chen
Krzysztof Polanski
Moshe Eliasof
Carola-Bibiane Schönlieb
Sarah A Teichmann
Venue: ICLR 2025 Poster (2025)
Keywords: spatial causal inference
spatially disentangled representations
spatial transcriptomics
mechanistic interpretability
Primary area: applications to physical sciences (physics, chemistry, biology, etc.)
TL;DR: Generative model that identifies and disentangles causal relations in spatial (transcriptomics) data to perform spatial counterfactuals, thus proving insights into experimentally inaccessible states, with potential applications in human health.

Abstract: Models of Virtual Cells and Virtual Tissues at single-cell resolution would allow us to test perturbations in silico and accelerate progress in tissue and cell engineering. 
However, most such models are not rooted in causal inference and as a result, could mistake correlation for causation.
We introduce Celcomen, a novel generative graph neural network grounded in mathematical causality to disentangle intra- and inter-cellular gene regulation in spatial transcriptomics and single-cell data. 
Celcomen can also be prompted by perturbations to generate spatial counterfactuals, thus offering insights into experimentally inaccessible states, with potential applications in human health. 
We validate the model's disentanglement and identifiability through simulations, and demonstrate its counterfactual predictions in clinically relevant settings, including human glioblastoma and fetal spleen, recovering inflammation-related gene programs post immune system perturbation. 
Moreover, it supports mechanistic interpretability, as its parameters can be reverse-engineered from observed behavior, making it an accessible model for understanding both neural networks and complex biological systems.

Published as a conference paper at ICLR 2025
ESTIMATION OF SINGLE-CELL AND TISSUE PERTURBA- TION EFFECT IN SPATIAL TRANSCRIPTOMICS VIA SPA- TIAL CAUSAL DISENTANGLEMENT
Stathis Megas"
University of Cambridge, Cambridge Center for Al in Medicine, Wellcome Sanger Institute em6130cas cam ac. uk
Daniel G. Chen" " & Krzysztof Polanski University of Cambridge, Wellcome Sanger Institute
Moshe Eliasof & Carola-Bibiane Schönlieb University of Cambridge, Department of Applied Mathematics and Theoretical Physics
Sarah A. Teichmann
University of Cambridge, Wellcome Sanger Institute (ec10003ccam.. Bec.mac..uu. ac.uk
ABSTRACT
Models of Virtual Cells and Virtual Tissues at single-cell resolution would al- low to test perturbations in silico and accelerate progress in tissue and cell engineering. However, most such models are not rooted in causal inference and result, could mistake correlation for causation. We introduce Celcomen, novel generative graph neural network grounded in mathematical causality to dis- entangle intra- and inter-cellular gene regulation in spatial transcriptomics and single-cell data. Celcomen can also be prompted by perturbations to generate spa- tial counterfactuals, thus offering insights into experimentally inaccessible states, with potential applications in human health. We validate the model's disentangle- ment and identifiability through simulations, and demonstrate its counterfactual predictions in clinically relevant settings, including human glioblastoma and fetal spleen, recovering ddflamattion-elltter gene programs post immune system per- turbation. Moreover, it supports mechanistic interpretability, reverse-engineered from observed behavior, making an accessible model for understanding both neural networks and complex biological systems.
its parameters can
INTRODUCTION
cell's gene expression profile simultaneously encodes information about intrinsic characteris- tics and extrinsic tissue microenvironment, Disentangling these two effects and understanding the causal links between them are necessary to fully reconstruct the complex interplay of intra- and inter-cellular interactions in human tissues during homeostasis and post-disease therapy-induced perturbation (Rood et al. 2024) Megas et al. [0224b] Bunne et al. 2024). To achieve this, a robust framework for causal disentanglement needed.
Causal inference seeks to uncover the mechanisms that generate the observed data by leveraging the mathematical principle of identifiability (Pearl 2009 Khemakhem et al.) 2020 Xi & Bloem- Reddy 2023 Zhang et al. 2023b). This principle holds when there exists a single unique model that fits the data, and thus model. However, despite the several useful properties that directly follow from identifiability (e.g.
assured that our observations only be explained by the given
"Contributed equally this work.

---

Published as a conference paper at ICLR 2025
robustness, generalizability, and self-consistency), the vast majority of current deep learning models violate this principle (Khemakhem et al., 2020 Somepalli et al. 2022). Moreover, applications of causal inference in spatial context have been limited (Reich et al. 2021 Papadogeorgou & Samanta 2024) and there is a need for causally identifiable methods in spatial causal inference.
One of the goals in representation learning includes disentangling the data into conceptually distinct variables to achieve a representation of it that is more interpretable, generalizable, and suitable for reasoning (Squires al. 2023). The simplest setting of conceptual distinctness of variables is that of statistical independence, which can be approached using independent component analysis (Comon 1994). However, data and our reasoming about it often involves variables that are not statistically ly independent but are linked together in causal acyclic graph (Pearl (009). For example, the expression of ligand gene a in a cell not independent of but causes the expression of a target gene A in a neighboring cell j.
The generalization of learning statistically independent variables to learning sets of variables that are linked together in causal graph is referred to as causal disentanglement (Squires et al.] 2023), and combines the insights of disentangled representation learning and causal inference. In its most ambitious form, causal disentanglement aims to learn both the latent variables and the causality structure of the underlying causality diagram with identifiability guarantees (Zhang et al. (023b)) Although disentangled representations offer more interpretable descriptions of the data and in- sights into the inner workings of neural networks, the emerging field of mechanistic interpretability (Bereska & Gavves 2024 Ferrando & Voita 2024 Rai al. 2024) has the even bigger ambition to completely specify neural network's 's computation. In its most granular form, this entails compre- hensively reverse-engineering the model's weights, which is of particular interest safety alignment of LLMs (Bereska & Gavves 2024).
Concurrently with the theoretical developments in causal inference and machine learning, technol- ogy developments, such as VisiumHD (Nagendran et al. (023), Curio, Stereo-seq (Cheng et al.) 2023) and Xenium (Salas et al. 2023), have allowed profiling the gene expression of cells at single- cell resolution as well as the cell's spatial coordinates in the tissue. These datasets, referred to as spatial transcriptomics, have also enabled us to perform perturbation experiments, such as gene knock-outs, in spatial samples and at large scale (Binan et al. 2023). This advancement has created the need for joint causal modeling of cell and tissue architecture (Rood et al., 2024, Megas et al., 2024b), to better capture the causal links of gene regulation.
In this paper, we introduce Celcomen to address the problem of inferring and disentangling the causal structure diagram of feature interactions from spatial samples, illustrated in Fig. [IA A typical example is gene regulation within cells compared with the behavior across nearest neighbors. That is, we aim to understand the causal links between features within node, but also the causal links across nodes, see Fig. B which might be confounded by spurious correlations see Fig. [C.D. We then use it to perform spatial perturbations, such as spatial gene knock-outs, in selected regions of a tissue to help guide efforts of cell and tissue engineering. In other words, we here aim to introduce a model of Virtual Tissues. Models of Virtual Cells predict the effect that changes in the micro- and macro- eevero-eniirommettt of the cell (such as perturbing the age of the donor, the tissue the cell in, the drug treatment, knock-outs from guide RNAs etc) have on gene expression (Bunne et al. 2024 Roohani et al.) 2024) Megas et al.] 2024b). Conversely, model of Virtual Tissues aims to not only estimate the effect the environment has on the cell but also the effect that the cell has on its environment and overall tissue.
Main Contributions. The contributions of this paper are two-fold, and offer advancements from both (i) computational biology and (ii) machine learning aspects. In terms of (i), the contributions include:
A proof of concept that generalizing models of Virtual Cells to models of Virtual Tissues is possible.
A novel framework of causal structure learning for feature interaction in graph data, such as gene-gene interactions in spatial transcriptomics data.
Gene regulation inference by integrating both dissociated and spatial single-cell data. A framework to perform counterfactuals on spatial transcriptomics graphs, such as answer- ing the question "What would the cells in this tissue have looked like had we performed a gene knock-out in a specific location in the tissue?".
2

---

Published as a conference paper at ICLR 2025
Feature interactions in graph
rudee interazion
cell
call
gate
cell
gena
gene
Ha He
odl
edik
A) Structural causal model (SCM) between feature variables $i, Sg. The directed acyelic graph (DAG) of the SCM has edges of different types, depending on whether a variable causally influences another in the same node, the nearest neighbor nodes (or nodes further away). This causality diagram generates spatial data such as spatial transcriptomics data. B) Our method is able to I) disentangle the causality graph into its two components, and 2) retrieve the undirected Markov equivalence class of the graph of each component. Our model infers and disentangles the multi-type SCM based on samples of those features organized in space. C) Gene 1 causes gene 2 to be expressed in the neighboring cell and gene causes gene in its neighbor. D) Half the time the 1 l-neighbor of a 1-neighbor also 1 1-neighbor leading spurious colocalization of genes 1 and 3 although there is no causal link.
Figure 1:
* To our knowledge, the first causally identifiable model for spatial transcriptomics analysis.
In terms of (ii), this paper offers:
The construction of a k-hop GNN whose network parameters are interpretable, and guar- anteed assume the values of the problem's Lagrange multipliers, which are meaningful, physical quantities.
A toy model for mechanistic interpretability, since Celcomen's weights can be recovered from its input-output computation.
We provide a theoretical understanding of the model, namely, we mathematically prove that the problem of spatial causal disentanglement is identifiable. We also computationally verify our theoretical findings.
2 RELATED WORK
coll-mmmmiiiaiinnn Several previous works on cell communication have relied on prior knowl- edge of protein-protein interactions (PPI) gene regulatory networks (GRN) to distinguish intrinsic and extrinsic circuits of gene regulation; this reliance often excludes key cell-cell interaction part- ners that are unreported (Browaeys et al., 2020, da Rocha al., 2022). Recent deep learning models advance on this limitation by simultaneously modeling intrinsic and extrinsic features; however, these models lack interpretable insight due to their black-box nature (Schaar al.) (2024). Further, virtually all current models lack mathematical (identifiability) guarantees (Khemakhem et al. 2020).
leading to their hyper-sensitivity to input data variability; exceedingly few accept both spatial and single-cell input data (Birk et al. 2024 Bernstein et al. 2022); and many cannot perform in silico perturbation experiments critical to understanding tissue behavior during disease (Cang al.. 2023 Wilk et al. 2024 Jerby-Arnon & Regev] (022). While these works have introduced marked compu- tational leaps in spatial transcriptomics, they often cannot perform causal inferences due their lack of identifiability which mathematically prevents many current models from deriving comprehensive mechanistic insights into cell and tissue biology. Altough applications of causal inference and coun- terfactual predictions in spatial transcriptomic have been limited, there are by now several such methods in the dissociated single-cell field achieving state of the art performance (Tejada-Lapuert
et al.
2023 Aliee et al. 2023 Piran et al. 2024 Zhang et al. 2024).
Graph Neural Networks. K-hop Graph Neural Networks (GNNs) Brossard et al. 2020 Chien
et al. 2020 Wang al. 2020 Abu- al. 2019) are generalization of the message-passing
3

---

Published as a conference paper at ICLR 2025
type of GNNs. Although most paassigg GNNs iteratively aggregate information from the neighbors to update node representations, K-hop GNNs perform message passing from not only the 1st hop but all the neighbors within K. hops of the node, and different weights and activations are used for each hop. K-hop GNNs strictly generalize message-passing GNNs, since an L--laeer messaging passing GNN is L-layer K-hop GNN with K=1. Note that when we say K-hop GNN, we mean that each node receives messages from all the neighbors that have a distance from node In contrast, k-hop GNNs pass messages from only the neighbors that
less than equal exactly distance k from node v. K-hop GNNs have been mathematically shown to be more powerful than K=1 GNNs encoding the topology of the graph and can distinguish almost all regular graphs [eng et al. 2023).
Spatially RNA--88
Correlation between expected and obeerved gene-gume rnterction mark
Vidumm
0.75-
Raw Inpul output p=2.90x10
0.50-
Dissoclated soa te single cal)
Corencoding types allows the dieentanglenme of intra (Hal varsus inter (Ma)
coeftiieet
1.00-
impul CCC
gene-gere
0.75-
inter-callular
intra-celulan
0.50-
p-8 88x10'
(FC)
0.25
0.00
1.00-
Rew Input Intra-cellular p=3.72x10
Gene expression with spatial resolution
Celcomen (DCE)
0.75
collutar gene gene interaction
0.50
0.25
0.00
CCE
CCE
10*
Raw inpul
output
10%
p=1.4tx1011
10*
10%
SCE
Spearmn soneeation (negathre
10
Raw inpul
SCE
SCE
exploration molecular
10%
inter-cellular p-3.38x10*
and celluter per perturbation
104
10%
10
Rew Inpul
SCE
Multile reelietic spatial maancciipooiccss generated from
10-
Intra-cellular
p-1 0000000
10%
Figure 2: Celcomen reproduces its identifiability guarantees in simulations. a) Inference module (CCE) can learn gene-gene relationships from either spatially resolved and (optionally) dissociated RNA-seq data. The highlighted cell-cell pair, in spatial data, and individual cell, in ScRNA-seq data, emphasizes how CCE distinguish gene-gene interactions that intra- (H' (Haa) cellular gene-gens interactions. b) Generative module, called Simcomen (SCE), leverages learned gene-gene relationships from CCE counterfactual tissue behavior after cellular or genetic perturbation. c) Box plots with the x-axis as the dataset being correlated against the raw noise inputted into SCE generate data as magenta and the learned output from SCE or the learned gene- gene relationships from CCE as pink. The y-axis depicts the Spearman correlation coefficient rho (upper three) or correlation p-value (lower three). Mann- Mann-Whitney Whitney U-test p-values are labeled on the center right of each plot and the legend on the upper right of each plot labels each box's dataset.
inter-
METHOD
We describe the method proposed in this paper, called Celcomen. We start by discussing its motivation and inspiration stemming from Lagrangian mechanics in rection[.1] then the biological assumptions behind Celcomen in Section 3.2. and finally the model's derivation in Section 3.3
4

---

Published as a conference paper at ICLR 2025
3.1 MOTIVATION AND INSPIRATION
Causal inference in machine learning aims extract causal structures from data. As such, it stands in between correlation- -bsee methods, and mechanistic models (Peters et al. 2017). To distinguish between correlation and causality, let us assume, for instance, that gene in cell causes gene 2 in its I-neighboring cells and gene in a cell causes gene 3 in its 1-neighboring cells (see C.) Since in some cases (see Fig D). the 1-neighbor of the 1-neighbor of a cell is also 1-neighbor (and not 2 2-neighbor) of that same cell, there will be (spurious) colocalization also of genes and 3 in nearest neighbors. A causal model should be able to de-confound such spurious connections within spatial correlations, even without mechanistic data such as epigenetic information. The inspiration for our work comes from the notion of force in physics. In broad strokes, we aim learn the "least" number of forces (i.e., least number of causal mechanisms that can explain all the observed colocalizations of pairs of genes) that can explain the observed spatial correlations of pairs of genes. "Leaat" here is meant in the sense of entropy, not the absolute number of forces, although, we could additionally impose a & norm penalty on the force matrix.
In the Lagrangian formulation of classical physics, we think of time evolution of physical objects as an optimization problem (optimizing the action) such that certain constraints imposed by Lagrange multipliers obeyed. One can show that Lagrange multipliers to impose the corresponding constraint, which means that they are meaningful. physical quantities Gollsteein T980). At the same time, imposing the constraint via the Lagrange multipliers allows us to remain agnostic about the nature of the force (whether it is electromagnetism, gravity, or nuclear forces) that imposes the constraint. For example, for an forced to walk on a table's surface, the constraint's 's force is electromagnetism, but we do not need to know this in advance to calculate its value. Similarly, in single-cell genomics, measurements are valued in a high-dimensional gene expression space, but they often are hypothesized to lie on lower dimensional surface (Fefferman al. 2013) due to biological mechanisms (already discovered or not) that "force" our measurements to lie on it. Uncovering such causal links is the first step to identifying the underlying molecular mechanisms. We use Lagrange multipliers to impose the observed colocalization of genes.
equal to the force required
Our key idea is to make graph neural network whose parameters are guaranteed to assume the val- ues of the Lagrange multipliers of the problem at hand. This idea has two powerful implications: (i) Because Lagrange multipliers are meaningful and physical quantities, finding them is likely to be well-posed problem, leading to a causally identifiable model. (ii) Since they are meaningful quanti- ties they should be easier recover, thereby providing a new avenue mechanistic interpretability.
3.2
MODEL ASSUMPTIONS
Our model, Celcomen, is the unique model that follows from three assumptions:
The expected gene-gene correlations across Ist neighbors match exactly the observed ones.
The expected gene-gene correlations within spots/cells match exactly the observed ones. Causal sufficiency: there are no unmeasured common causes of any pair of genes that are under consideration, such as different experimental study sites, different tissues which might have more or less signature of long-range hormonal regulation, etc.
These three assumptions can be summarized in the following equation for the entropy
S(P({s))), Iup:
Iugg Sms) = Σ
log(((()
}
Σ
Σ
Σ
8f sf
+
a
i.j mm
(.j
+ Σ
Saye
(1)
5 s )esp)
(1)
a.
where si_ is the spatial gene expression and P(ts is the probability distribution over spatial tran- scriptomics samples, and Sag gay are Lagrange multipliers that enforce our assumptions I and 2. For a summary of our notation, see Appendix A and for a discussion of what the biological limi- tations and implications of these assumptions Appendix Our task is to maximize the entropy
5

---

Published as a conference paper at ICLR 2025
functional in equation I over all possible functions P E (RRN8X) and matrices 9a38 and gasi
max ((P(fs))), 90 Panaa
9a Iaa).
(2)
3.3 MODEL DERIVATION
We should note that the optimization problem in equation 2 particularly hard non-parametric problem because it requires optimizing over the space of normalized functions. On similar note, the entropy is a functional that needs be maximized, not merely a function. To this end, using functional calculus, we perform the maximization of the entropy functional in equation [ over all functions P E T.² (RRN5)) to obtain simpler optimization problem over g.g alone. This relaxed optimization problem is more amenable to neural networks and will reveal the required architecture for Celcomen. The aim of this section to derive the suitable structure Celcomen that directly follows from the assumptions of the previous section.
Theorem (Extremization over P). The following two optimization problems are equivalent
Maximizing the entropy functional in equation 7 over all possible functions P E L (x8)) and matrices 9aß and I
max ((P(fs))) 90 Pag) 9.8,))9))) P.g.a
(3)
where is given by equation 1
Minimizing the experimental/empiriaal log-likelihood over matrices 9oz and Sag
min (log Plopp = =min 9.8"
log ((ga339a + Bap + Ios Cap
(4)
aB
9.8"
where Cog Eis 83 Jiss Cas - E. sq 12. =
H1 e and
ΣΣ H =
ΣΣ ag i
8
Jiggoss
(5)
0,
Theorem whose proof is in Appendix B implies that maximizing equivalent to minimizing
(log P) exp log Z(gaz) + g0,3 Cesp + Iap Cap
where Ca. -i,, 5 Jiis = T Js, and Ia C Trisss JJ) Tr (Jsgs " see Appendix A
for notation reminders. We note that, the term Jsg is standard linear message passing as in GCN (Kipf & Welling, 2017) and other GNN methods, as shown in (Eliasof et al. 2023)
In summary, our non- non-parametria optimization over P. tells that the desired model architecture is
a k-hop graph convolutional network, similar to (Nikolentzos et al. 2020 Feng et al. 2023) with new and simpler loss function, Equation 6
A key challenge in causal structure learning and causal inference is determining whether multiple causal structures could explain the observed data or, more precisely, whether the causal structure identifiable (i.e.. unique explanation exists) from the data. For causal model to be well-defined, must have unique causal structure that fits data; otherwise, perturbation effects cannot be identified, and the model may arbitrarily select one of many equally plausible structures due to minor variations in data noise. We now present a thcoretical result (proved in Appendix F) that guarantees Celcomen's identifiability in causal structure learning and disentanglement.
H((sf()
Theorem (Identifiability). The model defined by P({ P(a Sos3) 7
is identifiable
in the sense that
((s P(fs)))l (a....-..))))))) (((sft)) (hap, bagl))
(7)
Jog8 has and Sag hag-
(8)

---

Published as a conference paper at ICLR 2025
In the following Section, verify Celcomen's identifiability of causal structure learning and dis- entanglement in cenfsssns synthetic and real-world data experiments.
Celcomen overcomes some limitations of existing methods (see section 2 by leveraging causally identifiable framework into a generative graph neural network for learning disentangled representa- tions of intra- and inter- cellular gene regulation in spatial transeriptomies data (Fig. 2a-b).
The inference module of Celcomen, hereby called CCE, finds disentangled representations of gene interactions the cell, 1st neighbor, 2nd neighbor, etc. levels. These representations can then be used by the generative module of Celcomen (Simcomen), hereby called SCE (see Appendix H and Fig. [p). produce single-cell spatially resolved predictions of tissue behavior post perturbation and to derive realistic slides of spatial transcriptomics data from noise.
Finally, thanks carefully deriving powerful approximations computationally heavy tasks in the training process (see Appendices CD) Celcomen easily scales to many graphs, cells, and genes, while the same time having small number of tunable hyperparameters.
4 EXPERIMENTS
We now validate the robustness and insightfulness of CCE and SCE across plethora of simula- tions, input types, and human tissues. summary, we demonstrate Celcomen as a mathematically grounded spatial and single-cell transcriptomics analysis tool that introduces the capability to per- form high-resolution spatially resolved perturbation predictions that are critical for clinically rele- vant disease modeling and tissue engineering efforts.
4.1 MATHEMATICAL IDENTIFIABILITY AND INTERPRETABILITY
We now confirm whether Celcomen's 's identifiability guarantees, discussed in section 3.3 hold in practice by subjecting Celcomen to multitude of self-consistency experiments.
Synthetic Problem Setting. We randomly generate ground truth set of feature-feature (e.g. gene- gene) interactions, which we encode in the networks weights. Next, we utilize Celcomen's gener- ative module, SCE, to generate spatial transcriptomics data representative of these gene-gene inter- actions. Following that, we feed the generated data into Celcomen's 's inference module, CCE, in an attempt to retrieve the originally encoded gene-gene interaction forces, which are also the network weights (see Appendix G for more details).
Results on synthetic data. In agreement with its identifiability guarantee, Celcomen consistently demonstrated strong aligmment between its inferred gene-gene interactions from its simulated data and the ground truth (Fig. (b).. This result suggests that Celcomen possesses strong self colfsconsistency, and thus identifiabiliy. This conclusion emerges because Celcomen can move between encoded gene-gene interactions to simulated spatial transcriptomics and then back to inferred gene-gene in- teractions with minimal, if not no, loss of information. It also demonstrates Celcomen's 's usefulness as a toy model for mechanistic interpretability since its network weights can be reverse-engineered from its observed input-output behavior.
Results on real-world data. Furthermore, we confirm Celcomen's 's identifiability guarantees on real human data, by applying the Celcomen model multiple spatial transcriptomics slides of human fetal spleen (Suo et al. 2022). For each slide, we trained sample-specific model and model trained on the remaining samples. We then correlated the gene-gend cene-gene interaction matrices of these two models. In line with its claimed identifiability, we observed Spearman correlations in the range of 0.5-0.6 between these two gene-gene interaction matrices even though they shared no training samples, as shown in Fig. 6 Moreover, the gene interactions captured in the intra and inter-cellular matrices are biologically sensible because they adhere to known biological intra and inter-cellular processes, shown in Fig. Therefore, through this experiment, demonstrate Celcomen's identifiability, by confirming that its implied stability and robustness extend beyond theory and synthetic data, and also be observed real human samples.

---

Published as a conference paper at ICLR 2025
0.06-
2.57x10
Inter Intoe gane-
gere prsessstt
0.04-
8.70x
0
10x10
-0.02-
Secreted
-0.04
Membrano
Spelial RNAA-eR (Xorum)
Intur- and intrae gemb-- gene imtraction
Interferon kmocked out (KO)
Equiiteium stute idenstte diferential
0.06
Cytoplasmic
40.0.
only
upon
.0.10
sample
IFTMB.
d
Perturbed spol (FTMMB sig. ko)
(
UFIMM sig medium) (F sig medium)
*st neighbore
2nd neighbbrrs
Unperturbed spote
(IFIMM3 sig. hight
Spatal aaxa
TMB sagreture
-1.60
0.9
1.4
a Tall neighters 2nd neighbors
Perturbed
Uncerturbee
12 42201e))
-1700
-1.75-
1.68
0.0
0.5
0.7-
0.8-
0.2
0.6
0.8.
1.2
1.0-
control spots
20
40
10
20
30
20
40
20
40
Spatial X-axis
Epech
Epoch
Epoch
Epoch
Perturted Unperturted
1st neightors Unperturbed
2nd neighbors Unperturbed
Unireected genes up
Hommesstilic prolieration Unstimulated neutrophil
in trymus genes up
genes
Etlccto col opees down
Pastivasoh penas doan
genes up-
Call death genes down-
PD1 Lganc ganes up
Sensscence dygrature
ILA signa ng
tydl
Thymus genas down
Calluier nomensiaals
saasad
Heathy genes up-
Homecstolic proceeces-
Na ve cell sighalure
i nolmmmaory response-
Lymphhoyts acivetion
PIBKAKK signaling
Pial menemmennghh
rojection genes
Pos: Poss-infeciio on
iPka pinduction
immune binding
Effector oel genes up-
Stiukk ressorse-
JNK pescede
Vira
Platelst
PSandduced Infarmations
infect en resporse-
FoR-based phagocytosis
eylokin response-
FOR palnway
3
-1.0 0.0 1.0 Enrichment score
-10 0.0 1.0 Enrichment score
1.0 0.0 1,0 Enrichmentscom soors
KO high low
Figure 3: Celcomen recapitulates known interferon knockout biology in human glioblastoma and disentangles intra- and inter-cellular gene-gene interactions. a) Public spatially resolved RNA-seq data, Xenium, of the human brain during glioblastoma is inputted into CCE to derive intra- and inter- cellular gene-gene relationships. Interferon (IFN) signaling is knocked out (KO) in a previously IFN+ cell and SCE learns the local and global effects of this perturbation. b) Bar plot with the x-axis as the subcellular localization of the gene as acquired from its gene ontology and the y- axis is the difference between the gene's inter- and intra-cellular gene-gene interaction terms. c) Visualization of perturbed spot and its nearest neighbors. d) Scatter plots with the x-axis as the epoch number and the y-axis as the interferon signature score of the given spot(s) at the specified epoch. e) Pathway enrichment analysis in each of the groups of cells from (c). Pathways were derived by first calculating pre- and perttuurrtrbtionn changes in gene expression in each cell, then identifying differentially changed genes between spot(s) of interest and unperturbed controls, this provides a ranking of genes that were differentially upregulated or downregulated in the interferon KO cell, its neighbors, as compared to the unperturbed control cells.
4.2 CAUSAL DISENTANGLEMER AND SPATIAL COUNTERFACTUALE IN REAL DATA
We now seek evaluate two other important properties of Celcomen its ability disentangle intra-cellular from inter-cellular gene regulation programs, as well as perform spatially resolved perturbation modeling. Disentangling such programs is crucial requirement of any Al model of virtual tissue at single-cell resolution, and having such virtual models would enable us to design and test perturbations in silico. thereby accelerating the discovery process biomedicine (as discussed in section D. In this section, we focus on counterfactuals involving the knockout of interferon signaling, crucial pathway the response of the immune system.
Disentangling intra- vs inter-cellular programs. To test the disentanglement abilities, we apply Celcomen in a real, human clinical setting, by analyzing a single-cell resolution spatial transcrip- tomics dataset of human glioblastoma (brain cancer), as illustrated in Fig. 3h. In congruence with our theoretical understanding of Celcomen, discussed in Thm. 2 we find that Celcomen can suc- cessfully disentangle intrinsic from extrinsic sources of transcriptomic variation. This is enabled through Celcomen's 's assignment of gene-gen interactions involving secreted genes as inter-cellular, and those solely involving cytoplasmic genes intra-cellular, as shown in Fig. 3p. It is important to note that the knowledge of which genes are secreted and which are cytoplasmic is not encoded into the model prior information, but rather is learned by the model in an unsupervised manner.
In silico knockout of interferon in a cell prevents its propagation to neighbors. We leverage Cel- comen's perturbation abilities to model interferon signaling, in the context of neurological tumor,
8

---

Published as a conference paper at ICLR 2025
where we investigated the scenario of interferon knockout. We choose model interferon signaling due to its critical role in cancer in inducing antigen presentation, inflammation, and immune acti- vation (Kruse et al. 2023: Haanen 2013 Gocher et al., 2022). First, we quantified the expression gene program by averaging differentially upregulated genes in interferon (IFITM3) high versus low cells. Next, we knocked out interferon expression in ran- domly chosen interferon-hig cell (Fig. Bp.. Utilizing this interferon score, we not only confirmed our in-silico knockout of interferon in the perturbed cell, by observing its marked loss of interferon- associated genes, but we also observed loss of interferon signaling in neighbors of the perturbed cell (Fig. Bi), This behavior is highly consistent with known interferon biology as interferon signaling physically propagates from cell to cell within human tissues; thus, recapitulating this intercellu- lar signaling phenomenon supports the validity of Celcomen's 's perturbation modeling (Green et al.] 2017 Mesev et al. 2019, Lukhele et al. 2019). Moreover, the genes that are active in the intra- and inter-cellular matrices are also associated with processes known to be intra- and inter-cellular, respectively, illusrated in Fig. 7
of our sample's interferon-aasociatee
In silico knock out of interferon in a cell shuts down the immune system in neighbors. To fur- ther confirm the validity of our interferon knockout modeling, we performed pathway enrichment on genes that were differentially changed in perturbed (and perturbed neighboring) compared to un- perturbed cells (see Appendix GG. Indeed, we find that post-interferon knockout, perturbed cells and their neighbors significantly downregulated characteristic interferon response programs compared to unperturbed cells (Fig. (.. For example, we observed the perturbed cells have decreased T cell effector and activation gene programs, as well as greater loss of infection-relate gene sets and marked increases in regulatory programs. The consistency of our model with multiple aspects of known interferon biology strongly affirms Celcomen's 's ability to model perturbations with spatial resolution. Thus, through an in-depth study of Celcomen's 's application on a real human sample, we provide validation to its value in disentangling intra- versus inter-cellular gene regulation programs and in performing high-resolution spatially contextualized perturbation modeling with accuracy.
4.3 VALIDATION OF SPATIAL COUNTERFACTUAL IN-VIVO
To further demonstrate the usefulness of our Celcomen, we now benchmark it on the in-vivo full- transcriptome dataset that has measured gene knockouts in spatial transcriptomics, called Perturb- map (Dhainaut et al.) (2022). and compare it to random baseline. To the best of our knowledge, this the first learning method that attempts to utilize this dataset.
The Perturb-may map dataset consists of mouse model for KP lung cancer where in addition there might be either a Jak2 or a Tgfbr2 knock-out. Their dataset has annotated 5 spatial regions as lesions whose parts either 1) KP wild-type cancer, 2) KP cancer with Jak2 knock-out, 3) KP-cancer with Tgfbr2 knock-out, see Fig. 4h. Notably Perturb-map is based on the Visium technology which is not single-cell resolution because each spot contains cell communities of around seven cells average.
Celcomen correctly captures the perturbation effects in the spatial context. To learn the cancer biology of the KP model, isolated the wild type parts of the lesions by removing all spots with guide RNAS as well as their --neighbor spots. We then applied CCE on this part of the data to learn the gene-gene force matrices of KP cancer biology.
We then chose one of these spots, knocked-out Tgfbr2 (respectively Jak2) in it, and applied SCE to simulate in-silico the effect of this knock-out, see Fig. 44. Spearman correlation was used to compare model Tgfbr2 KO versus WT gene rankings with those directly derived from experimental Tgfbr2 KO spots and WT which was used as ground truth. We derived "random" controls for each lesion by computing correlations on shuffled gene rankings of the observed and predicted differentials between Tgfbr2 KO and WT. Mann-Whitney comparing observed lesion-derived gene rankings with those from random shuffling. Our results show Spearman correlations between our predictions and in-vivo measurements in the range of 0.28-0.47 for all lesions. To assess the significance of this performance, we compare our model to a random baseline where Celcomen is run on randomly shuffled data. We find that Celcomen's performance is significantly higher, with p- valuue of 0.0079, compared to the random baseline, as shown in Fig. 4-f.
test used to derive p-value when

---

Published as a conference paper at ICLR 2025
b
Isolate wild type lesion remove neightors
Rondontly chose Tglto? perturb
Knockout Tifh/22 randomiy chowen spot
Arcctation From Cnaaonan and
that perturbed
2022
g lose
Spata
#KP-2
gfbrr high
Parlurbed
mJlk? KO TTgibi2 KD
Spatal
ePorphory
Spertal
Spatial X-anis
Lesion
Lesion 2
Losion 3 2p
Lesion po4 Op
Lesion
d
All Lesions RRandom
100-p-23 Ba
pos do
de
Ranted owne WVr
(model prodiciio)
.og (p-velule
RD-
9 0
20
p=7 Se
100
100
+02
Spearman corretation
Geres are Iested percentle--oeeed tankinng,, ..s. 0-20%
coofficient rho
Lesion Ae
Ledion 2
Leeion 3
Lesion
Lesion
All Lesions Dbservedd
verpus Wr
(modal predictioo)
100
1-0.58
Am
de
ped
Be
Rardom
Eo-
Fomteed
p-7 Sw
0.0
+0.2
104
Jak2 WT
Spearmnn correlation confficient rho
Genes are Bitted by percent besedd rankinge, 9.9. 0-20%
Figure 4: Counterfactual predictions validated in-vivo in lung cancer model. a) Scatter plot on spatial axes with each Visium spot colored by tumor cell phenotype ("KO" for knockout and "KP" for wild-type tumors). Lesions of interest, large enough for modeling, are labeled with numbers. b) Example workflow with a wild type (WT) lesion trimmed for spots within two- degrees of perturbed clusters, random Tgfbr2+ spot has Tgfbr2 knocked out, Celcomen then predicts the whole transcrip- tome accompanying this perturbation. c,e) Box plots, per lesion, with x-axis as the observed ranked differentially expressed genes (DEGs) between Tgfbr2 (respectively Jak2) KO and WT, and the y- axis as the model predicted gene ranking between our perturbed Tgfbr2 KO spot and wild type spots. Spearman correlation coefficient rhos and p-values are annotated on the plot. d.f) Scatter plot with each dot representing given tumor lesion with Tgfbr2 (respectively Jak2) KO and the x-axis as the Spearman correlation coefficient rho and y-axis as the p-value, the color indicates if the correlation was computed on the lesion's observed gene rankings or a random shuffling of the gene rankings. Mann- Whitney U test p vvluu between observed and randomly shuffled correlations are annotated.
5 CONCLUSIONS AND DISCUSSION
The advent of single-cell resolution spatial transcriptomics has allowed spatial tissue atlases with un- precedented resolution (Farah al., [2024] Megas al., 2024a| To et al.
2024 Yao et al.) 2023) Lin- deboom et al., 2024) Zhang et al.] (2023a). Current computational methods addressing spatial tran-
scriptomics focus on phenotypic characterization but often neglect combined cell- and tissue-level causal perturbation modeling (Rood et al., 2024) which could reveal the mechanisms behind tissue disease states. Here, we present first step towards a model of Virtual Tissues, called Celcomen, which can predict the effect of spatial counterfactuals at the cellular and tissue levels thanks to its strong mathematical foundations. We confirm Celcomen's 's ability to disentangle and recover ground truth gene-gene interactions in real and self-simulated spatial transcriptonics data. Moreover, our model opens new route to mechanistic interpretability via causal inference. As we demonstrate in experiments, thanks to our model's 's causal identifiability, we can recover the values of the parameters of the neural network with high accuracy.
Broader Impact. Celcomen's 's advancements are poised to significantly impact biomedicine by re- vealing how diseases cause tissue failure and enabling testable hypotheses on therapeutic benefits. As technology progresses, Celcomen's 's value will grow, enhancing disease modeling and mecha- nistic understanding. Additionally, its architecture holds promise for advancing mechanistic inter- pretability, contributing to the development of causal foundation models.
10

---

Published as a conference paper at ICLR 2025
REPRODUCIBILITY STATEMENT
Our source code for the method available at https://shoturlattcNNCC and the code for reproducing the results of our experiments synthetic and real data is available at https://shorturlariis44.. The real data can be downloaded from the links provided in Appendix G.1
We also provide clear explanations and step-by-step proofs of all the mathematical claims in the paper. Thm. on rewriting our main optimization problem is proved in Appendix B and Thm. on the idenfiability of our model is proved in Appendix F Several more technical lemmas on deriving the mean-fiel approximation used in Celcomen are proved in Appendices D
ETHICAL STATEMENT
In this work, we do not release any datasets or models that could be misused, and we believe our research carries no direct indirect negative societal implications. We do work with sensitive or privacy-related data, nor do we develop methods that could be applied to harmful purposes. To the best of our knowledge, this study raises ethical concerns risks of negative impact. Additionally, all human samples used in this work were downloaded from appropriately licensed sources (see Reproducibility Statement (5). We also confirm that there are no conflicts of interest external sponsorships influencing the objectivity or results of this study.
REFERENCES
Sami Abu-El-Haile EljHa, Bryan Perozzi, Amol Kapoor, Nazanin Alipourfard, Kristina Lerman, Hrayr Harutyunyan, Greg Ver Steeg, and Aram Galstyan. Mixhop: Higher-order graph convolutional architectures via sparsified neighborhood mixing. In international conference on machine learn- ing. pp. 21-29. PMLR, 2019.
Hananeh Aliee, Ferdinand Kapl, Soroor Hediyeh- Zadeh, and Fabian Theis. Conditionally invariant representation learning for disentangling cellular heterogeneity. arXiv preprint arXiv:2307.0055:: 2023.
Leonard Bereska and Efstratios Gavves. Mechanistic interpretability for ai safety a review, 2024.
URL htttt:////aaxii org/ abs/2404 14082
Matthew N Bernstein, Zijian Ni, Aman Prasad, Jared Brown, Chitrasen Mohanty, Ron Stewart, Michael Newton, and Christina Kendziorski. Spatialcorr identifies gene sets with spatially varying correlation structure. Cell Reports Methods, 2(12), 2022.
Loic Binan, Serwah Danquah, Vera Valakh, Brooke Simonton, Jon Bezney, Ralda Nehme, Brian Cleary, and Samouil L. Farhi. Simultaneous crispr screening and spatial transcriptomics reveals intracellular, intercellular, and functional transcriptional circuits. bioRxiv, 2023. doi: 10.1 1023.11.30.56600 URL https: hetts::/www.bboxxxtt www.::/wioxivvv.o 01/2023. 569494
Sebastian Birk, Irene Bonafonte-Pardass Adib Miraki Feriz, Adam Boxall, Eneritz Agirre, Fani Memi, Anna Maguza, Rong Fan, Gonçalo Castelo-Brande Fabian J. Theis, Omer Ali Bayraktar, Carlos Talevera-López, and Mohammad Lotfollahi. Quantitative characterization of cell niches in spatial atlases. bioRxiv, 2024. doi: 10. 10.11001/22240022..5842222 14288 URL hhttp:: //www. biirxiiv.org/conteett... ott.://eeeettaaael//////////////
Rémy Brossard, Oriel Frigo, and David Dehaene. Graph convolutions that can finally model local structure. arXiv preprint arXiv:2011. 20..200611..,,,,,, 2020.
Robin Browacys, Wouter Saclens, and Yvan Sacys. Nichenet: modeling intercellular communi- cation by linking ligands to target genes. Nature Methods, 17:1-4, 02 2020. doi: 10.1038/ $41192-019-0067-5
Charlotte Bunne, Yusuf Roohani, Yanay Rosen, Ankit Gupta, Xikum Zhang, Marcel Roed, Theo Alexandrov, Mohammed AlQuraishi, Patricia Brennan, Daniel B. Burkhardt. Andrea Califano, Jonah Cool, Abby F. Dernburg, Kirsty Ewing, Emily B. Fox, Matthias Haury, Amy E. Herr,
11

---

Published as a conference paper at ICLR 2025
Eric Horvitz, Patrick D. Hsu, Viren Jain, Gregory R. Johnson, Thomas Kalil, David R. Kel- ley, Shana O. Kelley, Anna Kreshuk, Tim Mitchison, Stephani Otte, Jay Shendure, Nicholas J. Sofroniew, Fabian Theis, Christina V. Theodoris, Srigokul Upadhyayula, Mare Valer, Bo Wang, Eric Xing, Serena Yeung-Levy, Marinka Zitnik, Theofanis Karaletsos, Aviv Regev, Emma Lund- berg, Jure Leskovee, and Stephen R. Quake. How to build the virtual cell with artificial intelli- gence: Priorities and opportunities, 2024. URL https: arxiv. org/ abs / 200.. 11654
Zixuan Cang, Yanxiang Zhao, Axel A. Almet, Adam Stabell, Raul Ramos, Maksim V. Plikus, Scott X. Atwood, and Qing Nie. Screening cell-cell communication in spatial transcrip- tomics via collective optimal transport. Nature Methods, 20(2)2118228,, Feb 2023. ISSN 1548-7105. doi:
URL brtps://eci.org/I00
64 1415522-22201722
Mengnan Cheng. Yujia Jiang, Jiangshan Xu, Alexios-Fotios A Mentis, Shuai Wang, Huiwen Zheng, Sunil Kumar Sahu, Longqi Liu, and Xun Xu. Spatially resolved transcriptomies: comprehensive review of their technological advances, applications, and challenges. Journal of Genetics and Genomics, 50(9):625 640, 2023.
Eli Chien, Jianhao Peng, Pan Li, and Olgica Milenkovie. Adaptive universal generalized pagerank graph neural network. arXiv preprini arXiv:2006.079888, 2020.
Pierre Comon. Independent component analysis, new concept? Signal Processing. 36(3):287-314 1994. ISSN 0165- 1684. doi: http://dc..agg///.ooooooo1111111
URL https://
www. wccsccenin www.eceeiiii...ccc...o..c..cc..c rect. com http:///eeeeeerrrreee//////1111/
Higher
Order Statistics.
Edroaldo da Rocha, Caroline Kubaczka, Wade Sugden, Mohamad Najia, Ran Jing, Arianna Markel, Zachary LeBlanc, Rafael Peixoto, Marcelo Falchetti, James Collins, Trista North, and George Daley, Cellcomm infers cellular crosstalk that drives haematopoietic stem and progenitor cell 105-8/444556-022000355550
development. Nature Cell Biology, 24:1-11, 04 2022. doi:
Maxime Dhainaut, Samuel A Rose, Guray Akturk, Aleksandra Wroblewska, Sebastian R Nielsen, Eun Sook Park, Mark Buckup, Vladimir Roudko, Luisanna Pia, Robert Sweeney, et al. Spatial crispr genomics identifies regulators of the tumor microenvironment. Cell, 185(7):1223-1239 2022.
Moshe Eliasof, Lars Ruthotto, and Eran Treister. Improving graph neural networks with learn- able propagation operators. In International Conference Machine Learning, pp. 9224 9224-9245. PMLR, 2023.
9245.
Elie N Farah, Robert K Hu, Colin Kern, Qingquan Zhang, Ting- Yu Lu, Qixuan Ma, Shaina Tran, Bo Zhang. Daniel Carlin, Alexander Monell, et al. Spatially organized cellular communities form the developing human heart. Nature, 627(8005):554-864, 2024.
Charles Fefferman, Sanjoy Mitter, and Hariharan Narayanan. Testing the manifold hypothesis, 2013.
URL htttp::////tiiii org. abs 1310. 0425
Jiarui Feng, Yixin Chen, Fuhai Li, Anindya Sarkar, and Muhan Zhang. How powerful are k-hop message passing graph neural networks, 2023. URL
13328
Javier Ferrando and Elena Voita. Information flow routes: Automatically interpreting language models at scale, 2024.
John Gardiner and Stathis Megas. 2d tqfts and baby universes. Journal of High Energy Physics, 2021(10):1-43, 2021.
Angela M Gocher, Creg Workman, and Dario AA Vignali. Interferon-y: teammate or opponent in the tumour microenvironment? Nature Reviews Immunology, 22(3):158-172, 2022.
Herbert Goldstein. Classical Mechanics. Addison Addison-Wesley, 1980.
Daniel S Green, Howard A Young, and Julio C Valencia. Current prospects type ii interferon ?
signaling and autoimmunity. Journal of Biological Chemistry, 192(54):13225-1399 2017.
12

---

Published as a conference paper at ICLR 2025
John BAG Haanen. Immunotherapy of melanoma. European Journal of Cancer Supplements, (2):97-105, 2013.
Livnat Jerby-Arno and Aviv Regev. Dialogue maps multicellular programs in tissue from single-cell spatial transcriptomics data. Nature Biotechnology, 40<10):467-1477, Oct 2022. ISSN 1546-1669 doi: 10. 10:1088/44157-022-0188583 URL 6-1187-022-01288-0
llyes Khemakhem, Diederik Kingma, Ricardo Pio Monti, and Aapo Hyvärinen. Variational au- toencoders and nonlinear ica: A unifying framework, 2020. URL https::/axxiv.or
abs/1907. 04809
Thomas N. Kipf and Max Welling. Semi-supervised classification with graph convolutional net- works, 2017. http::tttp/sstteeee wwtp://rgggabbbbbbooooo /ULL:ttt.////ax...oo../////.. 02907
Per Kraus, Stathis Megas, and Allic Sivaramakrishnan. Anomalous dimensions from thermal ads partition functions. Journal of High Energy Physics, 2020(10):1-35, 2020.
Bastian Kruse, Anthony C Buzzai, Naveen Shridhar, Andreas D Braun, Susan Gellert, Kristin Knauth, Joanna Pozniak, Johannes Peters, Paulina Dittmann, Miriam Mengoni, et al. Cd4+ t cell-induced inflammatory cell death controls immune-evasive tumours. Nature, 618(7967): 1033-1040, 2023.
Rik GH Lindeboom, Kaylee B Worlock, Lisa M Dratva, Masahiro Yoshida, David Scobie, Helen Wagstaffe, Laura Richardson, Anna Wilbrey-Clatt Josephine L Barnes, Lorenz Kretschmer, al. Human sars-cov-2 challenge uncovers local and systemic response dynamics. Nature, pp. 1-10, 2024.
Sabelo Lukhele, Giselle M Boukhaled, and David G Brooks. Type interferon signaling, regulation and gene stimulation in chronic virus infection. In Seminars in immunology, volume 43, pp- 101277. Elsevier, 2019.
Stathis Megas, Anna Filbrey-Clark Aidan Maartens, Sarah A Teichmann, and Kerstin B Meyer. Spatial transcriptomics of the respiratory system. Annual Review of Physiology, 87, 2024a.
Stathis Megas, Nadav Yayon, Kerstin B Meyer, and Sarah A Teichmann. Molecular connectomies: Placing cells into morphological tissue context. PLoS biology, 2208):3302203, 2024b.
Emily V Mesev, Robert LeDesma, and Alexander Ploss. Decoding type and ii interferon sig- nalling during viral infection. Nature microbiology, (66)9914-224, 2019.
Mchryar Mohri, Afshin Rostamizadeh, and Ameet Talwalkar. Foundations of Machine Learning. MIT Press, 2018
Monica Nagendran, Jerald Sapida, Joey Arthur, Yifeng Yin, Sultan Doganay Tuncer, Naishita Anaparthy, Anushka Gupta, Marco Serra, David Patterson, and Augusto Tentori. 1457 vi-
sium hd enables spatially resolved, single-cell scale resolution mapping of ffpe human breast cancer tissue. Journal for ImmunoTherapy of Cancer, 11(Suppl 1):A1620-A1620 2023. doi: 20t0:/iie-2222-SS2200000220 Supp1_ 1/21620
URL https: peatte.bmj.. com/ content/1//
Giannis Nikolentzos, George Dasoulas, and Michalis Vazirgiannis. k-hop graph neural networks, 2020. URL htt.httts://ariiv...tt
abs/ 907. 06051
Georgia Papadogeorgou and Srijata Samanta. Spatial causal inference in the presence of unmeasured confounding and interference, 2024. URL httpss::///aatt org/ abs /2003 08218
Judea Pearl. Causality. Cambridge University Press, 2 edition, 2009.
Jonas Peters, Dominik Janzing, and Bernhard Schlkopf. Elements of Causal Inference: Foundations and Learning Algorithms. The MIT Press, 2017. ISBN 0262037319.
Zoe Piran, Niv Cohen, Yedid Hoshen, and Mor Nitzan. Disentanglement single-cell data with biolord. Nature Biotechnology, pp. 1-6, 2024.
13

---

Published as a conference paper at ICLR 2025
Daking Rai, Yilun Zhou, Shi Feng, Abulhair Saparov, and Ziyu Yao. A practical review of mech- anistic interpretability for transformer-bast language models, 2024. URL httpt::/aaxxx org/
abs/ org/abs/ 2407. 02646
Brian J. Reich, Shu Yang, Yawen Guan, Andrew B. Giffin, Matthew J. Miller, and Ana Rappold. A
review of spatial causal inference methods for environmental and epidemiological applications. doi: https://doi.org/00.111 1ins/s 2452. URL http:: //ol i libeayibrrry... wi ey com do /aoo///bs/lie /ttp///b///....//////....
International Statistical Review, 89(3):605-633 634, 2021.
12452 2452
Jennifer Rood, Anna Hupalowska, and Aviv Regev. Toward a foundation model causal cell and
tissue biology with perturbation cell and tissue atlas. Cell, 187(17):4520-444, 2024.
Yusuf Roohani, Kexin Huang, and Jure Leskovec. Predicting transcriptional outcomes of novel multigene perturbations with gears. Nature Biotechnology, 42(6):927-935, 2024.
Sergio Marco Salas, Paulo Czarnewski, Louis Kuemmerle, Saga Helgadottir, Christoffer Matsson- Langseth. Sebastian Tismeyer, Christophe Avenel, Habib Rehman, Katarina Tiklova, Axel Ander- sson, et al. Optimizing xenium in situ data utility by quality assessment and best practice analysis workflows. BioRxiv, pp. 2023-02 2023.
Anna C. Schaar, Alejandro Tejada- Giovanni Palla, Robert Gutgesell, Lennard Halle, Mariia Minaeva, Larsen Vornholz, Leander Dony, Francesca Drummer, Mojtaba Bahrami, and Fabian J. Theis. Nicheformer: a foundation model for single-cell and spatial omics. bioRxiv, 2024. doi: 000110012224.04455555447777 URL https: http://wwwwwom www.:w/wooxvv.ccmm www.booiii...ro/creeee...... cat1k//20..70/////2022...11
589472
Gowthami Somepalli, Liam Fowl, Arpit Bansal, Ping Yeh-Chiang, Yehuda Dar, Richard Baraniuk, Micah Goldblum, and Tom Goldstein. Can neural nets learn the same model twice? investigating reproducibility and double descent from the decision boundary perspective. In Proceedings the JEEE/CVF Conference on Computer Vision and Pattern Recognition, pp. 13699-137000 2022.
Yang Song, Jascha Sohl-Dickstein, Diederik P. Kingma, Abhishek Kumar, Stefano Ermon, and Ben Poole. Score-based generative modeling through stochastic differential equations, 2021. URL
https://arxiv: www.:/iiiirrraeooooooo abs wttp/a//aaa....o/aae/////// 13456
Chandler Squires, Anna Seigal, Salil Bhate, and Caroline Uhler. Linear causal disentanglement via interventions, 2023. http:/tpps///rx..oog//aaaas/e
Chenqu Suo, Emma Dann, Issac Goh, Laura Jardine, Vitalii Kleshchevnikov, Jong-Eun Park, Rachel A Botting, Emily Stephenson, Justin Engelbert, Zewen Kelvin Tuong, et al. Mapping the developing human immune system across organs. Science, 376(659))eab005100, 2022.
Alejandro Tejada-Lapuerta, Paul Bertin, Stefan Bauer, Hananeh Aliee, Yoshua Bengio, and Fabian Theis. Causal machine learning for single-cell genomics. arXiv preprint arXiv:2310. 14935, 2023.
Ken To, Lijiang Fei, J Patrick Pett, Kenny Roberts, Raphael Blain, Krzysztof Polański, Tong Li, Nadav Yayon, Peng He, Chuan Xu, et al. multi-omic atlas of human embryonic skeletal devel- opment. Nature, 635(8039):657-667, 2024.
Martin J. Wainwright. High--Dienssonal Statistics: A Non-Asymplotic Viewpoint. Cambridge Series in Statistical and Probabilistic Mathematics. Cambridge University Press, 2019.
Guangtao Wang, Rex Ying, Jing Huang. and Jure Leskovec. Multi-hop attention graph neural net- work. arXiv preprint 20Xiv:2009. 14332, 2020.
Aaron Wilk, Alex K. Shalek, Susan Holmes, and Catherine A. Blish. Comparative analysis of cell-cell communication at single-cell resolution. Nature Biotechnology, 42(3):470-455 Mar 2024. ISSN 1546-1696. doi: 10. 10.108/4411877.02---777... URL https: http://oooooo///////// or00
1038/ 0538//41587-03-000
Quanhan Xi and Benjamin Bloem- Reddy. Indeterminacy generative models: Characterization and strong identifiability, 2023. U.I.https///axxi.crr//bbb. 00801
14

---

Published as a conference paper at ICLR 2025
Zizhen Yao, Cindy TJ van Velthoven, Michael Kunst, Meng Zhang, Delissa McMillen, Changkyu Lee. Won Jung, Jeff Goldy, Aliya Abdelhak, Matthew Aitken, et al. high-resolttion transcrip- tomic and spatial atlas of cell types in the whole mouse brain. Nature, 624(7991):317-332, 2023.
Bao Zhang, Peng He. John EG Lawrence, Shuaiyu Wang, Elizabeth Tuck, Brian A Williams, Kenny Roberts, Vitalii Kleshchevnikov, Lira Mamanova, Liam Bolt, et al. A human embryonic limb cell atlas resolved in space and time. Nature, pp. 1-11, 2023a.
Jiqqi Zhang, Chandler Squires, Kristjan Greenewald, Akash Srivastava, Karthikeyan Shanmugam, and Caroline Uhler. Identifiability guarantees for causal disentanglement from soft interventions, 2023b.
abs/ 2307. 06250
Ziqi Zhang, Xinye Zhao, Mehak Bindra, Peng Qiu, and Xiuwei Zhang. sedisinfact: disentangled learning for integration and prediction of multi-batch multi-condition single-cell maa-sequeneng data. Nature Communications, 15(1);912, 2024.
15

---

Published as a conference paper at ICLR 2025
A
APPENDIX: NOTATION
For case reference, we provide here a guide to our notation choices:
E RRNSS, count values for spot/cell i and gene where Vi ||s:|| = 1,
H. Hamiltonian of system,
EEee H((()) the partition function,
Σeissnn' sum over pairs of nodes {i,5} that nearest neighbors,
9, the number of nearest neighbors (that assume interacting),
S, the entropy functional,
S, the number of spots/nodes
spatial graph,
N, the number of features/genes in the graph,
9oo E RNXN Lagrange multiplier enforcing gene-gene correlations,
Jij E RRS the spatial adjacency matrix between spots/nodes in the graph,
Pr, the average with respect to the probability distribution P.
(Dep. the empirical/experimental average with respect the observed samples,
P(s))= E NXS)) the probability density of the count matrix of a spatial, transcrip- tomics experiment equals the matrix s."
B APPENDIX: EXTREMIZATION OVER P
Here we prove Theorem]
[ which states that the following two optimization problems are equivalent:
Maximizing the entropy functional in eq over all possible functions P E L Nx and matrices fag and sap
max ((P(ts))).903 Pag)
P.g.9'
where S is given by
Minimizing the experimental/empiriaal log likelihood over matrices Iag and Dag
min (log P)esp min 9:9
log Z (gans Mag) + 9n.9 as + gagG ag
(10)
9.5¹
where
Cao Σ
J
(11)
i
Ca Σ
(12)
Σ H((+)))
(13)
H = a8
+ 3Σ a3 ij
of JigBans
(14)
Proof. Optimizing functional requires taking derivatives with respect to functions. In particular,
8 $ff(x)dz 8(x y), we can maximize S with respect to P:
using
8f(y)
6S
0 8P(s') log p(1s'))) 1
Σ
JaB
Σ
8 5
tn 18
+
9a8
Σ
18
(15)
a.9
i.j nn
H({+?))
(16)
7
16

---

Published as a conference paper at ICLR 2025
where we normalized the probability function and denote
Si Saps
Σ
E
is 9as8)
(17)
<i.jaan
Σ
Σ
+
>>>
82
(18)
0,3
Z = Σ
(19)
Maximizing with respect to the Lagrange multipliers Pag, I'g gives:
5 - Σ
of s cxp
(20)
2.J n
1.j nn
0 =
11 Σ
(21)
Texp
Moreover, by substituting 16 into 1 we get
in
((P((6))), Pas) log 7
gap Σ
Fao Σ
8
0000 exp
(22)
-log P(s))exp
(23)
Therefore maximizing S is equivalent to minimizing
(log P)exp = log ((a + IasC + SasC
(24)
where Ca. Eis 89 Jiof
C APPENDIX: MEAN GENE APPROXIMATION
Even after performing the optimization over functions P, the simpler loss function 4 of our k-hop GCN is still intractable to compute because calculating the partition function (and its derivatives) requires summing over large number of possible spatial transcriptomics datasets.
Several famous algorithms in machine learning circumvent computing the partition function dif- ferent ways. For instance, contrastive learning essentially takes the ratio of probabilities, thereby canceling out the partition functionn score-based diffusion uses score-matching to learn a model of the gradient of the log of the probability density function (Song al. 202I). which again avoids computing the partition function.
In this paper, introduce novel approximation to the partition function, inspired by physics Gardiner & Megas (202I); Kraus et al. (2020), which has never been used before spatial tran- scriptomics. This
Mean Field Theory approximation
ōsf. m 85f
(25)
where we assume that the gene expression does not fluctuate too much around the mean.
Lemma 1 (Mean Field Theory approximation). Using the mean field theory approximation. St 88 + 886 = m + 5sh and retaining terms of order (bss), ), the partition function assumes the form
2 =
(26)
H = Σ
Σ a,
(ghx +
9 (as)) (-m))(mmmm) mm + +m8 + mmm s
(27)
Proof. Our Mean Field Theory approximation is
55 m 85
(28)
17

---

Published as a conference paper at ICLR 2025
where we assume that the gene expression does not fluctuate too much around the mean.
Using this, we can rewrite the exponent as
i gass) = Iao ds + 6s7)
(29)
[B + 0000 + o5
(30)
9 m (s m²) + m (s) (sf- m)
(31)
= gaol mmm m - m 9 of + m m²)
(32)
where in the second line we used the MFT approximation to neglect terms of order higher than 2, and
of Lag = Iapl-m²mo + mmss + m² mn?)))
(33)
This implies that the intercellular term in the exponent can be rewritten as
Easas-mmi m² 8f m (08))==
ΣΣ Jas( m +m m² m
(34)
(1.8)
a.
where g the number of nearest neighbors that we assume are interacting, and therefore
H 2
ΣΣ
Jas( (a m² + m² of + n
(35)
0,
Σ
Hapl (-al-m=" mmmam + m² + m² ),2)
(36)
Σ
Σ
(ohs + 9 (9s8)( (-s) m² 8 mm m s?
(37)
0
since 9a.3 symmetric.
D APPENDIX: COMPUTING THE PARTITION FUNCTION IN THE MFT APPROXIMATION
To sum the exponential of the Hamiltonian over possible count matrices 8f we first prove the following Lemma.
Lemma 2. The following sum can be simplified as follows
(e9H/2 -9-4/22 8
Σ (
exp ΣΣ 2 (m=) of + m =m²;²)) = Vgn-1
(38)
9H/2
where S is the number of spots, Has = Baß + D3ar volume of the dimensional sphere.
VESlE. Hagma)? and Vgn-1 the
18

---

Published as a conference paper at ICLR 2025
Proof.
2
Σ
exp
(m 2 o7 + mm )na)))))
(39)
a.p
Σ exp
9
Σ
(geam* of + 9usm= sf9)
(40)
2
{ap)
a.p
Σ exp {+8)
(gsa + Pag ) Σ Σ
(41)
2
a.p
Σ exp
q
H Hagmas'
(42)
2
(+P)
a.
(J
II
deg
exp
Σ
Σ
Hagma
(43)
BLES
2
a.3
exp
ạ 2
Σ
Hs
dag
(44)
V
exp
Σ H
cos 1 sin 0d0
(45)
/ exp
g
Σ
Hu
du
(46)
2
(e9H/2
8
Vgn
(47)
qH/2
(48)
where the number of spots, Has = 9as + Jāa, v EE Hagma)? and without loss of generality we assumed that the vector Ea Hagm lies only along the first dimension.
Now applying lemma to the results of lemma 2 where we just need to replace Jag + Jag + Jos gives
log 7 -= - (gas ạ )a m
a
+ log Vgn-1
1212 H2
+Slog
(49)
H1/2
where S is the number of spots, H'g 490,3 498a + 29'n 2g H' H' - H' may
we have complete formula for calculating the log-likelihood, and the only optimiza- tion remaining is over the Lagrange multipliers.
Using 4 49
8P({s
0
(50)
8gas
0 8P([s
(51)
8g'n9
In other words, want look for the forces that causing the observed spatial gene expression. Since the Lagrange multipliersfforees are meaningful physical variables, they naturally equip our model with identifiability show in the next section.
19

---

Published as a conference paper at ICLR 2025
A
B
h *
h
h
H
H
h2 *
h **
h
Figure 5: A) An illustration of the decomposition of enror in PAC model. The model explores the fitness of a subset of hypotheses H. B) A non-identifiable model has many hypotheses that have the smallest approximation error. This lead non-robust results small change in the data can cause the model to transition to an alternative hypothesis.
E APPENDIX: REVIEW OF PROBABLY APPROXIMATELY CORRECT LEARNING
Here we offer a brief review of Probably Approximately Correct (PAC) learning, to set up some notation and definitions that will aid our discussion of identifiability. For more in-depth discussion of the relevant theory, readers can consult (Mohri et al. 2018).
In PAC learning, we aim learn the hypothesis h from set of available hypotheses Hi which has the smallest generalization error R(h) to given concept c. Let R = R(h"") be the smallest error achievable any hypothesis in 52, we call R the Bayes error and h² the Bayes optimal
hypothesis. Then the difference between the error of any hypothesis h and the Bayes error can be decomposed as
R(h) R(h**) = (R(A) R(h")) - (R(h*) (R(h")-At) R((h"))
(52)
eiror
approximation crror
where h H. is the hypothesis in H. that minimizes the error, Fig.
There are examples of famous bounds on the overall error which are often referred to as oracle identities but this requires extra assumptions about the distribution of the data expected model. For instance, the Lasso oracle inequality is proved assuming the restricted eigenvalue condition (Wainwright, 2019); and the oracle for graphical Lasso is proved assuming that the precision matrix is a spiky (Wainwright) 2019).
In the context of PAC learning, the identifiability of model is the property of having unique hypothesis As E IDI that minimizes the approximation error,
344" e H. wh H: |[(hh)) (R(h"))((""']1 < [R(h) R[h"")|
(53)
Lack of identifiability can lead to model results that are not robust since tiny differences in the data or random seeds could lead the model to transition to one of the other equally optimal hypotheses.
F APPENDIX: IDENTIFIABILITY RESULT
An important question we want to address is the identifiability our model, ie whether there is unique setting of the forces that leads the observed correlations in the data. If the identifiability property holds then our model will naturally be robust and causal in the sense that can deconfound spurious correlations.
In gauge theory terminology, we want to determine whether there is some gauge symmetry that allows different sets of parameters give the same probability distribution.
20

---

Published as a conference paper at ICLR 2025
We now prove one of the main results of our paper, Theorem 2 which states that the model defined by equation 16 is identifiable since
V(s{): P(Oo))aa faall = (Mas H
(54)
9a has and In h/
(55)
Proof. Let's pick to be a cell/node that has at least one neighbor. If there is not such cell then there wouldn't be cell communication problem to model.
P( M9 faph) = P(ts Hlhap, (aphaap Kaat) HP Pasl) dP Mihas Kagl) dsf
(56)
(57)
dsf
Then we pick to be any of the neighbors of cell i,
d de (Ass
aP Mhap
(58)
dxf dsa
ds; dss
d (P((f Paot)
dP H Nagl)
(59)
ds
0,+0==
ds ds
s 0.0f
9o.t
(60)
Alternatively, taking the second derivative with respect to the same cell i,
d P({s Hoga) dP( [ss])/[aasss (haaa Nasl)
(61)
d8; ds?
ds, ds? (M
P Kgas snab) ds, dsf
(62)
do; dsa
l =0=0,=²=-
No8 N/a8
(63)
Identifiability for model means that unique configuration of our model fits the observed data. This very important for causal models since otherwise, our model could tell us that A either causes B or it doesn't 't cause B. which is unhelpful tautology. Identifiability is point of principle of whether the model could definitively decide that A causes B if provided enough data. For this reason most proofs of identifiability happen in the infinite data limit.
G APPENDIX: METHODS
G.1 SPATIAL TRANSCRIPTOMICS DATASET CURATION AND PREPROCESSING
The fetal spleen datasets were curated from https//deeelommmnaalllllltttttot immune in log-normalized form, which explicitly indicates log-transformation and li- brary size normalization18. The glioblastoma dataset was curated from 10x genomics at https://www. 1ww.100x (0ggenomisscco///atttttttt www.a//wwwweeeiiiiiiioiaaiaannnn wath---hhwwithuumaa -imaano---oo----
and subjected the same library size normalization, counts per million (CPM), and bog-transformation, with base of e; additionally, only genes that were expressed in at least 100 cells were kept. Due to the large size of the Xenium slide, random square portion of the slide was chosen for analysis, this section is defined as cell centroid x-component > 6500 and < 7000 and cell centroid y-component > 8000 and < 8500. The entire fetal spleen slide was kept for each fetal slide sample as they comparatively smaller than the original Xenium slide and post down-sampling, approximately the same size as the analyzed Xenium section. All data normalization was done using Scanpy (v1.9.8 .9..) in Python (v3.9.18)31)
secology-profiliig
21

---

Published as a conference paper at ICLR 2025
G.2 SIMULATIONS TESTING CELCOMEN' IDENTIHIBBILTY GUARANTEES
Simulations were done in Python and completed by first generating random n-genes by n-genes matrix of ground truth gene-gene interactions, for these experiments, four genes were used. We then utilized Celcomen's s generative module to create simulated spatial transcriptomics slide whose val- ues were then fed into Celcomen's 's inference module to decipher back gene-gene interactions. Spear- man correlation was used to compare the original ground-truth gene-gene interaction values and the simulated- seetteed-tt--innfrree gene-gene interaction values to test for model robustness and identifiability. For all exact parameter values utilized during the experiments, see the iaaalyyis.siullationn ippynbbb notebook in the reproducibility GitHub.
G.3 BIOLOGICAL TESTING OF CELCOMEN'S IDENTIFIABILITY GUARANTEES
Biological confirmation of Celcomen's identifiability guarantee was done by training two Celcomen inference module instances at the same time and comparing their derived gene-gene interaction results. The first model instance, which call sample-specific, The second model instance, which we call rest, was training on the remaining samples. Thus, these two model instances are never trained on the same samples. Each model is trained completion utilizing the same model hyperparameters, and their gene-gene interaction matrices are retrieved after the final epoch. We correlate flattened version of their gene-gent gene interaction matrices using Spearman's 's correlation due to the possible non-linear nature of the matrices' values. We repeat this experiment for each of the samples the fetal spleen dataset. The results across each sample's experiments are aggregated together and compared in a bar plot. We derived a "random" control to compare to by shuffling the order the flattened gene-gene interaction matrices and computing a correlation of the shuffled values. are labeled on plot. For the full code utilized, see the "unalysis.boologcca" www.obiioollolyyyy...... ipynb" notebook in the reproducibility GitHub.
trained only on one sample.
U test is used to derive p-values and all p-values
G.4 INTERFERON KNOCKOUT EXPERIMENT ON XENIUM OF HIUMAN GLIOBLASTOMA
Processed Xenium data was subjected to the inference module of Celcomen, CCE, and then these gene-gene interaction values were annotated as containing cytoplasmic, surface membrane (plasma membrane GO ID via GO cellular component), or secreted (extracellular space GO ID also via GO cellular component) genes according their GO IDs from QuickGO32. IFITM3 was knocked out in a randomly selected previously IFTTM3 positive cell. First neighbors were defined as less than 15 jum away and second neighbors were defined less than 30 μm away. Changes in each gene's expression in each cell were calculated and these changes in expression pre- and post- perturbation were compared between different specified cellular subsets. These are the differential genes later used for differential expression analysis and pathway enrichment. Gene set enrichment analysis (GSEA) in (v4.1.2) was utilized to perform pathway enrichment analysis on differentially post- perturbation affected genes. The interferon signature was derived directly from tissue by computing the differentially expressed genes between interferon high and low cells and taking the top 25, ex- cluding the perturbed IFITM3 that would bias analyses. For the full model parameters and code utilized, see the mumalysis.pertubattinn
notebook in the reproducibility GitHub.
H APPENDIX: SIMCOMEN, GENERATION MODULE
Our model is a mathematically robust way of learning the distribution of spatial transeriptomics samples such that there is 1-1 correspondence between a configuration of forces and the learned distribution of spatial transcriptomics samples.
Generating new samples from the learned distribution is classic problem that can be addressed, for instance, by Markov Chain Monte Carlo Methods. However, given the high dimensionality of the space of spatial transcriptomics. MCMC can be very computationally expensive. Therefore, in our generation module, called Simulated Communication Energy, we produce new samples in an adversarial approach by trying find a sample that would trick Celcomen into thinking it is derived from the learned distribution. More specifically, we fix the parameters of our model and optimize the likelihood of the possible datasets. This approach is also followed for generating counterfactual
22

---

Published as a conference paper at ICLR 2025
samples, eg intervene node, and from that starting point find the most likely sample under the learned distribution of spatial transcriptomics.
I COMPUTATIONAL EFFICIENCY
Our derivation of the mean field approximation eq. partition function. In particular, with eq. 26 which is the un-approximated formula that derives from Theorem [ computing would require summing all possible spatial transcriptomies samples, which estimate be roughly equal to 106,000,000 possibilities (for 100 possible gene values, 1000 genes, and 3000 spots), Performing this calculation in every training step is computationally infeasible. On the other hand, using eq. 49] which describes our Celcomen, Z can be calculated by summing over 10 gene- cene-pairs. Using this approach, Celcomen's 's "inference* module runs in 11 seconds for 500 genes and 2500 spots. In the early phases of the design of the model, experimented with sampling from a space of 106.000,000 possibilities using MCMC or Gibbs sampling, which proved to be computationally intractable. That is the reason we implemented the generative module of Celomen to generate the most likely spatial transcriptomicss sample based on its learned distribution of observed samples. Thanks this, running the generative module for 500 genes and 2500 spots takes few hours on CPU or 12 minutes on a GPU.
49 allows to speed the computation of Z, the
J LIMITATIONS AND FUTURE DIRECTIONS
While this work demonstrates the potential of Celcomen in disentangling intra- and inter-cellular regulation, identifying causal links, and performing spatial counterfactuals, several limitations high- light exciting directions for future investigation.
Limitation Celcomen aims identifiably retrieve the undirected form of the DAG underlying intra- and inter-cellular regulation from observational data, i.c.. without interventional data. With- interventional data, impossible to uncover more than the Markov equivalence class of the DAG, which is the undirected form of the graph along with partial information about colliders. Therefore Celcomen's ability is close to being tight against the best possible bound on the abilities of unsupervised models but could be improved by training in supervised fashion on interventional spatial transcriptomics data, which might soon be available high-throughput.
Our modeling of gene regulation relies on DAG, and if there are cycles in gene regulation, then the Bayes optimal hypothesis bound to be outside the set of hypotheses explored by Celcomen. Extending the model to accommodate cyclic regulatory structures perhaps via a causal kinetic model Tejada-Lapuerta et al. (2023) is a promising direction for future work.
Limitation 2:
Limitation 3: A further limitation of our model is that it does not incorporate cell type information which can be incorporated into both the generative and the inference modules of our Celcomen. For instance, in the inference module, gene-CellType--Cllfyye matrices, to allow different gene-gene interactions for each pair of cell types. In the generative module, we could consider cell types when generating the counterfactuals. This amounts to imposing a prior probability for each cell/node in the graph. Incorporating cell type information into Celcomen can also enhance the interpretability of our method, especially in heterogeneous tissues, which have diverse composition in terms of cell types. On the other hand, incorporating cell type information might create challenges when cell type continua are arbitrarily clustered into separate cell types, and requires the generalization of our Theorem [ and Theorem 2 to the case of -ete-geee----eye-Ceeeeee research direction is left for future work.
infer not only gene-gene matrices, but also gene-
matrices. Therefore, this interesting and important
Limitation Our modeling of Batch effects as reflected assumption amounts to the causal sufficiency assumption. Indeed, without it, the identifiability proof would not go through, and the effect of counterfactuals would not be identifiable. We are intrigued by the possibility of weakening the causal sufficiency assumption while preserving the identifiability guarantees of the model and are open to exploring this in future work.
Limitation 5: Our modeling of gene regulation relies only pairwise gene forces and, as such, omits synergistic effects between genes. A straightforward extension Celcomen to include such effects could help understand more complex regulatory dynamics.
23

---

Published as a conference paper at ICLR 2025
Limitation 6: Finally, although technically not limitation, want to explain that for Celcomen to be able to disentangle intra- and inter-cellular gene regulation, it needs data that are of single-cell resolution. When spot and inter-spot gene regulation. Intra-spot gene regulation would involve signatures of intra- cellular, juxta-crine, and short-range paracrine gene regulation, whereas inter-spot gene regulation would involve long-range paracrine gene regulation.
multi-cell resolution data, such as Visium, Celcomen disentangles intra-
K
ENTENDED DATA FIGURES
coefficient rho (intra)) coefficient tho (inter-)
Correlation fetal co-encoded inter- cellular gene-gene Interactions
Spearmnn correaatio Spearman correlation
0.5-
0.5-
0.5-
p=0.03
0.3-
0.3
0.3-
0.3-
-04
0.2
0.2
0.2
0.2
0.2
0.1-
0.1-
Sample
Sample 2
Sample
Sample 4
0.0
0.0
Sample
Sample
Sample
Sample
p=0.03
0.4
0.4-
0.4
0.3-
0.3
0.3-
0.3-
0.2-
0.2
0.2-
0.2
0.1-
0.14
0.0
of epochs the model for Sample versus rest Randomly shuffled
Random
Observed
Figure 6: Extended Data Figure 1: Celcomen recapitulates its identifiability guarantees through strong sample-to-sampl epochs and y-axis the Spearman correlation coefficients between the gene-gene interaction ma- trices of the model trained on the specified sample and the model trained on all other samples. The sample utilized for the sample-specific model is annotated directly the plot. The color of the line, see the lower legend, indicates whether it represents comparisons between the two observed models, pink, between random shuffling of the two gene-gene interactions, black, represent null model. Right: Bar plots with the left black bar representing the average final Spearman correlation coefficient between randomly shuffled gene-gene interaction matrices the sample-specific model and model trained on all other samples, and the right pink bar representing the observed correlation. P-values are derived from the Mann- Mann-Whitney U test and are annotated directly on the plot. Error bars indicate standard error.
correlation on real human samples. Line plots with the x-axis as
24

---

Published as a conference paper at ICLR 2025
400 Pa
0:2
Thes tou Requidio O Cull gnd 101 COMBOO hopal Figa
of hour monrbann Proi healdin 101 IITOROO Happliae hup anfinn To batee 100 E00EE Deg
Phadse Nerria eeter Drea coeeerreeen (G0 20120077500
Rlupdattoo Of Call Nopulaaio doputation Pro unilon LCO CO02 aa
100 LA
of 2y 100000
(O0
from aisaiie brd mil teet
postt Rhep becpaiiie of deacaaae sue menarrttt
Of anl cneae tee
So the Ao aor
the Of th
All coaain
the
Heguis on Exil pod
poster restadig ticinni
Cal
sogganee of
tie an sespin ne slt of
Pusl te FOR coore
postee tos
pech Sool
Pag
Cluder (aaar))
Claar 5 (atraca
Cloter 1 (0M0a-Dekkirr lrr
Cluste 2 (k02)-0m.0)
C Gaste 4 (1889)
0.30
3
6:15-
Depptttoon cal Pro Provutioon Prolleetto (OO 20 2103
Degga Df Call IepAAso Trodlee Troeheetin 100 30 11500
Froll the dalts Per give frome (00 201
100 204100
Ragur d DE Dhes Dy Dymmonm 100 061522
304
Reguld on stenstttt Proftct (00 004
Opokin Irtu aar Peuperde (85 Co0e00
Open doneaiing (B0 CoTTe
5000
liggiie (D0 1500400
to Oofodme pespencee ceneeee 1O0 20544
werooo A0 00000000
Pootie Roge latoon orr Prosee bac Mecansole offire sarie 10 co0 so
Pootiive Ring Artion May of Procussr Oel Pry dilie ood 100 2000
Epen wid Tike pogggmeei (00 100000 Postie Pooy aatten or Pops aatio Pryall cosaenaati one 10 10000
Neggtis Pongy Aarten Of Onl rot vatoe USO 1000000
Nastos flag aad of Call moo mosliteee do condtii
Cyp the tig Mitii hope unteii COD COOOE 16
Hegiis Mon OF sAAr Canra Canaan Co HOROOO
Fleaay beg aaten Call Mealt GC 2000 47
NOW Aapesii Of hour peatto Resooooooo hespooto dO 6130077
Cyrth the seeeii cragaiing Thtt COCOTTIIO seeee
crane loceed To teroch owiie OC 101120
Der recsiaaion Of bood Prodi pod BC (00070
0.00
CCO
Oaster Collddd
HaLtnc Oyem fanniiee
pullee Regu Cfrsiieo megee oI Hesaon
Oerral thh Calk respe Sggai
of Of Wede
Of
huppliic And
supptes Repulll
Clader (ire-Cabsar)
Foother
Poetio
Cluster Collaar
Clster (Mier
Cluder 1 Grae
Figure 7: Extended Data Figure 2: Biological interpretation of Gene-Gene Interaction Matrices and Cluster Analysis for Xenium human glioblastoma. (a) Clustered intra-cellular gene-gene interaction matrix, highlighting clusters 3, 4. and 5 as the most strongly interacting intra-cellular modules. Path- way enrichment analysis identifies these clusters as associated with proliferaaive/apopossi functions (cluster 5), intracellular signaling (cluster 4), and post-cytokine stimulus response (cluster 3). which are predominantly intra-cellular controlled processes. (b) Clustered inter-cellular gene-gene inter- action matrix, emphasizing clusters 2, 4, and 5 as the strongest inter-cellular interactions. Pathway enrichment reveals cytokine response signatures, with inter-cellular cluster 2 capturing markedly more interactions related to cytokine and chemotaxis compared to the analogous ter. This highlights the model's 's ability to distinguish inter-cellular interactions selectively.
clus-
25

---

Published as a conference paper at ICLR 2025
0.50
10000
2
0.00
605
600
cotapor Me (GO
100
E10 Drt roreooo 101101000
Estocoliu Sen Danatogm conattgeeen (SC E
Aed Redicobbl Hunno a Repporse (00
Antiiieatia Hunaral Responon 190
feur Chindton (GD 00020022200
ID (01300
iequlattoo Of Anpogerentt 10D00045 mo)
tigrn ing faem fameeay IDD 0000 1000000
tht Chargen TOp Co0E E0
Altivly indden apd Procedo 1081 suueebo Parp lesttottinn oa Salattate Acteatinn Seperderr
patees seree 100001100
Communicaaaioo on By Coual in Circity kaat 1C0
Neegit on Of meniiiiion pooduc on (C 2013 20102001
seersas Logg 10D 2000 1907
1/3
Entr Wor THort Tide Nerthi pagrrthypoooa pooul 100 0012
Bag Requaaeenn of innelaatoon frokuct (000003100)
Pythudre Buctee baoo pesceel eroee 350 0108000
Dyttndra hie Vie siele trre 500
cnthmeeeee all Tummore inmure helpooe Hel ticltted by Artim crcbal fegios 000020
thes 1010 lecyiig AOOO
Recuuae fottee Thu Pru
Denccuin tan andd see C Teppo
Watter Gren
Rullb Repulation Of Polease Of Secuesters Caliun
by
Cunter: 2 (etrae)
dag sigettve Reguiation of Cyster ne Syee E-000oo
Custer2 2 (itra Caladar)
Cluster 4 Celulad
Cluster Cellaaa
b
Quster3 Unera Celuter)
0.20
1.0
20-
A9000000
0.14
55
6
0.10
Noutughh Ciooretads (CO OONO
Wopoo Frocess (60 (11EO
Encothela Call Difwer tattin UCD CEATE
Collge Wull O picustio 1S0 030
Ardinicosiaa Gumee du Ritratinn (50 0833
granching insteed skood Wosee Moep rei E05770
10000
The Te-Soaple Viari Hesooii Pocces (SO troocceeoongg deter Of Ang nose (CE (54)
tear Werphogents Al taithsla 'ste 00408105
Megalo rocyee De 100 2010000
good Cospration Thin F'R Cormatin CO GO72 Froden
Rysiid Gol Dove Foo 0300000
(o0 515
Choaaay Syttem riem 100 R07 2099 on
cake teetee Orepnd Orppeinn UPD tette
10 no
lapui Pubiie of aepaaton Deviloma CD 601-p
& modduct USD o01e67
Aubiem crrae Hummor tnumua no Resporue Mad by crneeceee peporne contrra
C de Positiie Regustttion Ferese in fesaenee ai de Sppar nonng (CD 1964
(o Wid (5554
hicess UCO C18111
C macc (CD EONOO
0.00
Pulbre Repuitten yoseediinn agga Thundd cunductoon (OD
of Fitretion (O
Custer Cel
Graand ve Dieertioon
barippn
Cel poremmmnn
Oesie Carl cooeertt Eesssi caseoeebty
go mocees Ner fritte to catin
Adteytore Cys Dasdooe Baya led on Of Calcia to foop aponcil moai
Mellogan Nerooadll dode
Res Sipool wel
Propl netuente TAE Signuing
darten Cafoilat
sag
Pox
OF
Iiegudation Of bereec Df
Cuiter 3 (rresCelulari
Custer (IttenGelaarr
Cluster 2 Cehiia)
Claster Cetarr
Cunter 5 LeteeCCllttrr
Figure 8: Extended Data Figure 3: Biological interpretation of Gene-Gene Interaction Matrices and Cluster Analysis for Visium human spleen. a) Clustered intra-cellular gene-gene interaction ma- trix, identifying cluster 2 as having the highest interactions, with cluster 5 as close second. Both clusters contain genes known to be associated with intracellular processes. Specifically, pathway enrichment analysis reveals that cluster 2 is associated with chemotaxis-related intracellular signal- ing, including chemokine response pathways and calcium release. Cluster 5 shows enrichment for nucleotide metabolism, potentially linked to the regulation of proliferation via nucleotide availabil- ity. b) Clustered inter-cellular gene-gene interaction matrix, highlighting cluster 2 and cluster 4 as the most strongly interacting gene clusters in this matrix. Cluster 2 enriched for antimicrobial re- sponse and angiogenesis programs, while cluster 4 is associated with phagocytosis and extracellular matrix-related functions. All of these tend to be associated with inter-cellular rather than intra- cellular functions.
26


Reviews:
[{"review_id": "UmKMCnRf68", "reviewer": "Reviewer_KiqG", "summary": "This paper introduces Celcoman, a physics inspired generative graph neural network to infer intra- and inter-cellular gene regulation in spatial transriptomics and single-cell data. Celcoman can also be used to simulate the perturbation effect (such as gene knock off) in silico, offering insights into experimentally inaccessible states. \n\nThis experiment is validated on both synthetic and real datasets, which includes both glioblastoma and fetal spleen samples.", "questions": "1. In section 3.1, you mentioned that the inspiration of this work comes from the notion of force in physics and the objective is to learn the \"least\" number of forces. Then you said \"'least' here is in the sense of entropy\". However, later on, all we are trying to do is to maximize entropy. After thinking, I think I sort of understand the concept but someone else may get stuck too. Could you add a bit more explanation on how does the objective of minimizing number of forces transform into maximizing the entropy?\n\n2. Could you comment on the runtime of this algorithm? Also, how long does it take (clock time) to finish the experiment you did in your study?", "limitations": "", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 4, "strengths": "1. I agree with the author's claim that this is the first causally identifiable model for spatial transcriptomics data. \n2. The motivation is well explained and the results on the impact of perturbed spot on surrounding cells are very convincing. It shows why we need to consider cell communication on spacial omnics data. \n3. The idea is novel. The author links GRN inference among neighboring cells with entropy maximization under the control of some lagrangian multipliers, which are then modeled by GNN. Theorem 1 is also inspiring.", "weaknesses": "1. Notation issue. The notations in this paper are not consistent. For example, the author used $s^i_{\\alpha}$, $s^\\alpha_{i}$, and $s_{\\alpha i}$ to refer the same concept.\n\n2. While I understand this is one of the first few papers studying, the lack of baseline comparisons in the experiments makes the results section less convincing. For example, while the results does show the impact of perturbed spot on surrounding area, how accurate it is? With GNN, I'm not surprised to see this kind of phenomenon. Are there any ways to provide a quantitative measurement or comparisons with the methods? Could you add some baselines or ground truth comparison?\n\n3. A common acknowledged fact is that for causal models, there is a huge gap between synthetic data and the data we observe in real life. While the idea of in-silico perturbation data generation is very interesting, can you discuss the similarity/difference between the generated data and real data?", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "4b3CyDzt2E", "reviewer": "Reviewer_pnJL", "summary": "The authors present a novel method for learning causal relationships in the context of spatial transcriptomics data, namely Celcomen. The method is devised to distinguish between inter and intra cellular causal interactions between genes. The authors characterize their method both from a theoretical point of view, demonstrating Celcomen identifiability (granted that assumptions hold), as well as from a practical point of view, applying their methods on real-world data.", "questions": "I would invite the authors to add an extensive paragraph regarding their assumptions:\n- Discuss the biological and technological implications of each assumption\n- Compare the assumptions to current limitations of spatial transcriptomics technologies\n- Identify specific conditions or tissues where the assumptions may not hold\n- Relate the assumptions to classical causal learning assumptions like faithfulness, Markov condition, and causal sufficiency\n\nModifying Celcomen so that to include cell type information is an onerous task that might deserve a separate publication. However, I would suggest that the authors:\n- Acknowledge this limitation of their current model\n- Discuss how disregarding cell type information may impact the interpretation of their results\n- Propose specific ways the model could be extended in future work to incorporate cell type or functional cluster information", "limitations": "", "rating": 8, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "The work has very strong theoretical basis. The authors provide proof of the identifiability of their method; moreover, they also honed their implementation so that to tackle large scale problem, for example by introducing a new Mean Field Theory approximation.", "weaknesses": "The main weakness of the manuscript is its poor connection to the underlying biology, as exemplified by the two issues below.\n\n- The authors present a set of three assumptions upon which their whole model relies. However, they never discuss the implications of these assumptions, and without appropriate commentary, I frankly find difficult to understand what these assumptions entails. How do their assumptions compare with respect to the current limits of spatial transcriptomics technologies? Are there specific conditions, tissues, or other situations where these assumptions may not hold? Moreover, how are these assumptions related to more classical assumptions of causal learning theory, namely faithfulness, Markov condition, and causal sufficiency?\n\n- Cells in tissues are usually grouped within functional clusters. Most typically, cells will belong to distinct cell types, which depending on the type of tissue may located together (endothelial cells in the inner layers of blood vessels) or scattered across the connecting matrix (e.g., fibroblast in mussle tissue). Other cases are possible as well: for example, glutamatergic neurons in the cortex will belong to different layers, where each layer will have different roles in cognitive functions. To my understanding, in its current iteration the Celcomen method fully disregards information on cell types or other functional clustering, even if this information can be instrumental in better understanding cell to cell communication mechanisms. For example, causal relationships between cell may be largely dictate on whether these cells are of the same or of different type.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "HwqG3t2bys", "reviewer": "Reviewer_zadU", "summary": "In the manuscript, ‘Estimation of single-cell and tissue perturbation effect in spatial transcriptomics via Spatial Causal Disentanglement’, the authors introduce Celcomen, a generative GNN that incorporates causal inference for the task of predicting perturbation effects and counterfactuals, with applications to spatial transcriptomics data and single-cell data. The authors aim, in addition to predictions, to expose mechanistic interpretability, such as the inference of intra-cellular and inter-cellular gene regulation in such data. The method is demonstrated over multiple synthetic and real biological data.", "questions": "Comments:\n- It’s unclear why the spatial structure of cells in tissues should be captured well by a DAG.\n- The text could be improved at multiple points for readability.\n- In the related work section, it would be better to add a subsection on generative models and disentanglement for single-cell and spatial transcriptomics data.\n - It would be better to give intuition in the text for the meaning of Theorem 2.\n- ‘we observed a strong positive correlation between these two gene-gene interaction matrices’ - there are no results/statistics in the main text, and in the supplement we see a maximum of 0.4/0.5 correlation value, so needs to be rephrased. Also, unclear how these results point to identifiability, as there only seems to be analysis of correlation \n- About the actual gene-gene interaction matrices themselves (corresponding to results in Fig. 5) - do these make sense in terms of the current knowledge of existing gene-gene interactions? There are many ways to evaluate these types of questions, and this is crucially missing from the evaluation and results in the paper.\n- The discussion section is missing a discussion of the limitations of the approach.\n\nMinor comments:\n- some of the references throughout don’t have parentheses around them.\n- when introducing spatial transcriptomics, better to give an additional example (or examples), beyond VisiumHD.\n- The connections between the visualizations in Fig. 1 and both its caption and the reference to it from the text, are not entirely clear.\n- There are some strange phrasing of text throughout. An example (but not the only one): ‘Suppose for instance that nature imposes the colocalization of genes..’\n- There are inexact sentences (again, example, but not the only one): “Since half of the time, the nearest neighbor of a nearest neighbor is also a nearest neighbor”\n- ‘”Least” here is meant in the sense of entropy’ - needs explanation\n- ‘for an ant forced to walk on the surface of a table, this force happens to be electromagnetism’?\n- After theorem 1, there are repetitions in the text relative to the numbered equations.\n- typo: ‘demonstrates that Celcomen’s usefulness as a toy model for mechanistic interpretability since’, ‘resolutionFarah’", "limitations": "", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "- The problem that the authors pose is important and timely\n- The method seems to be well grounded and robust\n- The demonstrations of the method on real data seems to be promising, especially those having to do with disentangling intra vs inter cellular programs in human glioblastoma spatial transcriptomics data, and showing that in silico knockouts affect signal propagation in a biologically expected manner.", "weaknesses": "See the 'Questions' section below for some comments on readability, assumptions made, missing limitations discussion, and presentation of results.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:

Decision: Accept (Poster)