Joint Neural Posterior Estimation and Model Selection for Chemical Reaction Networks

J. Arturo Esquivel; Jorge A. Campos-Gonzalez-Angulo; Sergio Pablo-García; Alán Aspuru-Guzik; Vianey Leos-Barajas

stat.MLOct 13, 2026

Joint Neural Posterior Estimation and Model Selection for Chemical Reaction Networks

Authors: J. Arturo Esquivel, Jorge A. Campos-Gonzalez-Angulo, Sergio Pablo-García, Alán Aspuru-Guzik, Vianey Leos-Barajas

Organizations: Department of Statistical Sciences, University of Toronto, Toronto, Ontario, Canada · Schools of Mathematics & Physics and Chemistry & Molecular Biosciences, The University of Queensland, St Lucia, Queensland, Australia · Australian Research Council Centre of Excellence in Quantum Biotechnology, St Lucia, Queensland, Australia · Queensland Quantum Decarbonisation Alliance, St Lucia, Queensland, Australia · The Matter Lab, University of Toronto, Toronto, Ontario, Canada · Institute of Micro and Nanotechnology · Departments of Chemistry, Computer Science, Materials Science & Engineering, Chemical Engineering & Applied Chemistry, and Institute of Medical Science, University of Toronto, Toronto, Ontario, Canada · Acceleration Consortium, Toronto, Ontario, Canada · Vector Institute for Artificial Intelligence, Toronto, Ontario, Canada · Canadian Institute for Advanced Research (CIFAR), Toronto, Ontario, Canada · NVIDIA, Toronto, Ontario, Canada · School of the Environment, University of Toronto, Toronto, Ontario, Canada

Abstract

Simulation-based Bayesian inference for Chemical Reaction Networks (CRNs) must simultaneously solve two computationally demanding problems: discriminating between overlapping reaction pathways and estimating reaction rate constants. Classical approaches, such as Sequential Monte Carlo ABC, rely on handcrafted summary statistics and scale poorly across large parameter spaces. To overcome this, we introduce a neural inference architecture. Our pipeline utilizes shared 1D-CNNs and a permutation-invariant Set Attention mechanism to learn low-dimensional summary statistics from variable-length, multi-experiment time-series data. These statistics simultaneously condition a model selection network to compute posterior model probabilities and model-specific conditional Invertible Neural Networks (cINNs) that estimate parameter posterior distributions. We demonstrate the practical utility of our approach on a real protodeboronation example, illustrating how it handles structural ambiguity and extracts kinetic parameters from competing catalyst deactivation pathways.

Figures & tables

Panel of Figure 1 — Figure 1: High-level architecture of the joint neural inference pipeline. Each time-series observation ( $4$ in total) goes through a parallel, shared-weight CNN to extract kinetic features. The architecture bifurcates to evaluate model evidence and condition $\mathcal{M}$ independent model-specific conditional Invertible Neural Networks to estimate parameter posteriors.

Observed Species	Time Points	Global Acc.	$M_{1}$ F1	$M_{2}$ F1	$M_{3}$ F1
Substrate & Product	20	94%	0.94	0.98	0.96
Substrate & Product	10	94%	0.93	0.97	0.93
Product Only	20	93%	0.91	0.97	0.92
Product Only	10	92%	0.89	0.95	0.90

Table 1: Model Selection Performance of Our Deep Learning Pipeline. We evaluate it across varying observational conditions: the impact of hidden chemical species (Product Only) and reduced temporal resolution (10 vs. 20 points). Global accuracy and model-specific F1-scores are reported for the

3,000

independent test trajectories.

Figure 2: Reliability diagram evaluating the calibration of the posterior probabilities across the 3,000 test trajectories. Predictions are grouped into 10 uniformly spaced probability bins. The dashed diagonal represents perfect calibration, where the estimated probability exactly matches empirical accuracy. All three mechanisms closely track the identity line, achieving low Expected Calibration Errors.

Parameter	$M_{1}$	$M_{2}$	$M_{3}$
$k_{1}$	$96.2\%$	$96.2\%$	$96.2\%$
$k_{-1}$	$96.4\%$	$95.2\%$	$96.0\%$
$k_{2}$	$98.1\%$	$97.3\%$	$95.1\%$
$k_{-2}$	$97.0\%$	$94.8\%$	$96.8\%$
$k_{3}$	–	$93.0\%$	$95.2\%$
$k_{-3}$	–	$96.2\%$	$96.0\%$

Table 2: Empirical Coverage of the

95\%

Posterior Credible Intervals. Results for the kinetic parameters and observation noise (

\sigma

) under the

20

-point (Substrate + Product) setting.

Figure 3: Posterior Predictive Check (PPC) for experimental protodeboronation data across three varying initial catalyst loadings. The $300$ predictive trajectories are sampled from the joint posterior, drawing from the candidate mechanisms proportionally to their estimated model probabilities.

Appendix figures & tables20 assets

Supplementary material from the paper’s appendix.

Appendix

Summary Network: CNN-Transformer
Component / Hyperparameter	Specification
1D-CNN filters (3 layers, kernel size 3)	$64\rightarrow 128\rightarrow 256$
Pooling Layer	Global Average Pooling (across time steps, masked for dropout points)
Dense layer ( $d_{model}$ )	$256$
Self-Attention heads	$4$
Feed Forward Network	$512\rightarrow 256$

Appendix

Table 3: Architecture and Training Hyperparameters for our Neural Inference Framework.

Panel of Figure 4 — Figure 4: Detailed architecture of the summary network $T_{\phi}$ . Independent experimental observations are processed through a shared 1D-CNN (indicated by dashed connections) to extract temporal signatures. The sequences are globally average pooled across time to ensure length invariance. The independent features are contextualized using a Set Attention Block (SAB), before being globally pooled across experiments and passed through dense layers to produce the unified embedding $\mathbf{e}$ .

Figure 5: Convergence diagnostics for SMC-ABC inference. The solid lines (left axis) track the shrinkage of the discrepancy threshold (

\epsilon

) across generations for candidate mechanisms M1, M2, and M3. The dashed lines (right axis) denote the corresponding acceptance rates required to maintain the target population of

N=3,000

particles per generation. As the posterior approximation refines and the

\epsilon

tolerance tightens, the acceptance rate exhibits a rapid decline, reflecting the increasing computational cost of proposing valid kinetic parameters from the transition kernel.

Figure 6: PPC for the four experimental observations, where SMC-ABC was implemented. Black lines denote the true observed synthetic data, with solid and dashed lines representing the substrate (S) and product (P) concentrations, respectively. Faint colored lines depict the 100 inferred trajectories generated by drawing from the network’s estimated posterior.

		Predicted Model			Performance Metrics
Setup	Truth	$M_{1}$	$M_{2}$	$M_{3}$	Precision	Recall	$F_{1}$ -score
(1) S & P, $T=20$	$M_{1}$	954	9	37	0.93	0.95	0.94
	$M_{2}$	18	976	6	0.99	0.98	0.98
	$M_{3}$	51	5	944	0.96	0.94	0.95
(2) S & P, $T=10$	$M_{1}$	936	17	47	0.92	0.94	0.93
	$M_{2}$	19	965	16	0.98	0.96	0.97

Appendix

Table 4: Classification Performance. Results for the three-model setting (

\mathcal{M}=\{M_{1},M_{2},M_{3}\}

) across all four observational setups. The left pane displays the confusion matrix across the holdout test set, and the right pane details the corresponding performance metrics.

	(1) S & P, $T=20$			(2) S & P, $T=10$			(3) P only, $T=20$			(4) P only, $T=10$
Param	$M_{1}$	$M_{2}$	$M_{3}$	$M_{1}$	$M_{2}$	$M_{3}$	$M_{1}$	$M_{2}$	$M_{3}$	$M_{1}$	$M_{2}$	$M_{3}$
$k_{1}$	96.2	96.2	96.2	95.0	96.2	96.0	97.7	94.9	98.0	95.9	92.6	97.4
$k_{-1}$	96.4	95.2	96.0	97.0	94.9	96.4	95.8	94.4	94.9	96.5	94.8	94.0
$k_{2}$	98.1	97.3	95.1	97.5	94.3	95.4	97.8	95.2	93.9	96.5	94.6	93.9
$k_{-2}$	97.0	94.8	96.8	97.7	94.8	97.3	97.0	95.4	94.3	96.6	94.6	94.8
$k_{3}$	–	93.0	95.2	–	89.5	93.1	–	97.7	94.0	–	96.5	93.4

Appendix

Table 6: Empirical Coverage for Three Models. Coverage (%) for the

95\%

credible intervals for the kinetic parameters and observation noise across all four observational setups. Values close to

95\%

indicate well-calibrated marginal posteriors.

Figure 8: Estimated posterior distributions, using

2,000

samples, for all parameters from a test case where the true mechanism was

M_{3}

. The true value (in the log-space) is marked with black dashed line.

Figure 9: Inference calibration across all testing cases from

M_{1}

. Left Panel: Posterior mean estimates vs. the true values. Right Panel: Rank histograms for all parameters (normalized to

(0,1)

). The dashed black line indicated the level for uniform density (

1

Figure 10: Inference calibration across all testing cases from

M_{2}

. Left Panel: Posterior mean estimates vs. the true values. Right Panel: Rank histograms for all parameters (normalized to

(0,1)

). The dashed black line indicated the level for uniform density (

1

Figure 11: Inference calibration across all testing cases from

M_{3}

. Left Panel: Posterior mean estimates vs. the true values. Right Panel: Rank histograms for all parameters (normalized to

(0,1)

). The dashed black line indicated the level for uniform density (

1

Figure 12: PPC for a test trajectory generated by mechanism $M_{3}$ under the high-resolution setup (1). Black dots denote the $20$ observations for both species across four distinct initial conditions. The solid lines represent $100$ posterior predictive trajectories sampled from the network’s inferred joint posterior.

Figure 13: PPC for a test trajectory generated by mechanism $M_{1}$ under the low-resolution setup (4). Black dots denote the $10$ observations for the product across four distinct initial conditions. The gray dots indicate the corresponding substrate, which was not observed by the model. The solid lines represent $100$ posterior predictive trajectories sampled from the network’s inferred joint posterior.

	Predicted Model				Performance Metrics
Truth	$M_{1}$	$M_{2}$	$M_{3}$	$M_{5}$	Precision	Recall	$F_{1}$ -score
$M_{1}$	789	6	19	186	0.75	0.79	0.77
$M_{2}$	17	977	6	0	0.98	0.98	0.98
$M_{3}$	54	9	928	9	0.93	0.93	0.93
$M_{5}$	190	4	43	763	0.80	0.76	0.78

Appendix

Table 7: Classification Performance for Four Models. Results for the expanded four-model setting (

\mathcal{M}=\{M_{1},M_{2},M_{3},M_{5}\}

) under the main observational setup (both species observed,

T=20

). The left pane displays the confusion matrix across the holdout test set, and the right pane details the corresponding performance metrics.

Figure 15: PPC for a test trajectory generated by mechanism $M_{1}$ under the high-resolution setup (1) with four candidate models. Black dots denote the $20$ discrete observations for both species across four distinct initial conditions. The solid lines represent $100$ posterior predictive trajectories sampled from the network’s inferred joint posterior.

Figure 16: Inference calibration across all testing cases from $M_{5}$ . Left Panel: Posterior mean estimates vs. the true values. Right Panel: Rank histograms for all parameters (normalized to $(0,1)$ ). The dashed black line indicated the level for uniform density ( $1$ ).

	$M_{14}$ (3.71%)		$M_{15}$ (92.19%)		$M_{16}$ (4.10%)
Parameter	Mean	Median	Mean	Median	Mean	Median
$\log k_{1}$	$3.9572$	$2.2728$	$3.5822$	$1.0382$	$4.2962$	$2.6632$
$\log k_{-1}$	$-5.6564$	$-6.7500$	$-6.1921$	$-7.2830$	$-5.8420$	$-6.8342$
$\log k_{2}$	$-1.8568$	$-3.9131$	$-1.8872$	$-2.5498$	$-2.0608$	$-4.0818$
$\log k_{-2}$	$-0.9392$	$-1.4539$	$-0.2301$	$-1.2468$	$-1.0966$	$-1.8808$
$\log k_{-3}$	$-2.1533$	$-4.2273$	$9.9645$	$10.5861$	$2.2611$	$0.2867$

Appendix

Table 9: Parameter Estimates. Posterior parameter estimates for the protodeboronation structural hypotheses.

Figure 17: Estimated posterior distributions for

M_{15}

, which was estimated to be the most likely mechanism from our candidate models. We use

2,000

samples, for all parameters estimated from real protodeboronation data. The estimated posterior mean (in the log-space) is marked with black dashed line.

Figure 18: Bivariate and marginal posterior distributions (

2,000

samples) for the most probable protodeboronation mechanism. The diagonal panels display the marginal posterior distributions. The lower-triangular panels feature bivariate KDE contour maps, illustrating the pairwise uncertainty between the kinetic constants.

Explore similar work

Reduction of Probabilistic Chemical Reaction Networks
Jun 26, 2026Mauricio Montes, Gregoire Sergeant-PerthuisProbabilistic Graphical Models Reduced-Order Modeling
COSAM, Auburn University, Auburn, AL, USA · CQSB, Sorbonne Universite, Paris, France
Causal Posterior Estimation
May 27, 2025Simon Dirmeier, Antonietta MiraBayesian Inference Simulation-Based Inference
Swiss Data Science Center, ETH Zurich, Switzerland · Università della Svizzera italiana, Switzerland · Insubria University, Italy
Tracking large chemical reaction networks and rare events by neural networks
Dec 11, 2025Jiayu Weng, Xinyi Zhu, Jing Liu +3Dynamical Systems Neural Surrogate Modeling
Institute of Data Science, University of Hong Kong, Hong Kong · Department of Systems Science, Faculty of Arts and Sciences, Beijing Normal University, Zhuhai 519087, China · School of Physics, University of Electronic Science and Technology of China, Chengdu 611731, China +6

Jun 26, 2026J/K move · Enter open · S save

Reduction of Probabilistic Chemical Reaction Networks

Mauricio Montes, Gregoire Sergeant-Perthuis

COSAM, Auburn University, Auburn, AL, USA · CQSB, Sorbonne Universite, Paris, France

Probabilistic Graphical Models Reduced-Order Modeling

Open paper AskPDF arXiv

Programming adaptive behaviors at the cellular level is a long-standing goal that raises the question of how probabilistic computation can be implemented in biochemical systems. Chemical reaction networks (CRNs) provide such a substrate and have been shown to realize probabilistic models, including hidden Markov models and factor graphs, with dynamics reproducing Bayesian inference and belief propagation. However, encoding these algorithms typically requires prohibitively large reaction networks, and classical CRN reduction techniques do not directly apply. By recovering the factor graph structure encoded in Napp--Adams-compiled CRNs, we transport recent factor-graph reduction results to their chemical implementations, obtaining significantly smaller CRNs while preserving the belief-propagation fixed points on surviving variables.

Abstract

Explore similar work

Reduction of Probabilistic Chemical Reaction Networks

Causal Posterior Estimation

Tracking large chemical reaction networks and rare events by neural networks