# Validation Diagnostics for SBI algorithms based on Normalizing Flows

Julia Linhart, Alexandre Gramfort, Pedro Luiz Coelho Rodrigues

► To cite this version:

Julia Linhart, Alexandre Gramfort, Pedro Luiz Coelho Rodrigues. Validation Diagnostics for SBI algorithms based on Normalizing Flows. NeurIPS 2022 - the 36th conference on Neural Information Processing Systems - Machine Learning and the Physical Sciences workshop, Nov 2022, New Orleans, United States. pp.1-7. hal-03856444v2

HAL Id: hal-03856444

https://hal.science/hal-03856444v2

Submitted on 24 Nov 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
# Validation Diagnostics for SBI algorithms based on Normalizing Flows

**Julia Linhart**

Université Paris-Saclay, Inria, CEA
Palaiseau 91120, France
julia.linhart@inria.fr

**Alexandre Gramfort**

Université Paris-Saclay, Inria, CEA
Palaiseau 91120, France
alexandre.gramfort@inria.fr

**Pedro L. C. Rodrigues**

Univ. Grenoble Alpes, Inria, CNRS, Grenoble INP, LJK
Grenoble 38000, France
pedro.rodrigues@inria.fr

## Abstract

Building on the recent trend of new deep generative models known as Normalizing Flows (NF), simulation-based inference (SBI) algorithms can now efficiently accommodate arbitrary complex and high-dimensional data distributions. The development of appropriate validation methods however has fallen behind. Indeed, most of the existing metrics either require access to the true posterior distribution, or fail to provide *theoretical guarantees* on the consistency of the inferred approximation beyond the one-dimensional setting. This work proposes easy to interpret validation diagnostics for multi-dimensional conditional (posterior) density estimators based on NF. It also offers theoretical guarantees based on results of local consistency. The proposed workflow can be used to check, analyse and guarantee consistent behavior of the estimator. The method is illustrated with a challenging example that involves tightly coupled parameters in the context of computational neuroscience. This work should help the design of better specified models or drive the development of novel SBI-algorithms, hence allowing to build up trust on their ability to address important questions in experimental science.

## 1 Introduction

Recent advances in computing have led to a new generation of expressive simulators used to study complex systems in many scientific fields [4]. They implicitly encode the intractable likelihood $p(x | \theta)$ of the underlying mechanistic model which relates observed data $x \in \mathbb{R}^d$ to scientifically meaningful internal parameters $\theta \in \mathbb{R}^m$. To perform statistical inference in this setting, one can recur to simulation based inference (SBI) [4] to approximate the posterior distribution $p(\theta | x)$ using samples from the joint pdf $p(x, \theta)$. In this work, we consider SBI methods based on normalizing flows [14], which are invertible neural networks that can be trained via maximum likelihood. Once the flow is trained, and for any new observation $x$, one can directly evaluate the estimated density over the parameter space (i.e. $\theta$-space) [6], draw samples to construct confidence regions [11], etc.

Flow-based SBI has been used in many recent applied works [10, 1, 6], but it lacks an important feature before becoming a technology for experimental science: *validation*. Ideally, one would like to have a method that provides finite-sample guarantees of nominal coverage (or calibration) of the estimated posterior regions, but also ensures that the approximation fits the true underlying posterior of the model when new data is observed. While solutions have been proposed to provide such finite-sample guarantees [5, 11], assessing the convergence and consistency of the underlying
inference remains a challenging task [10]. Simulation based calibration (SBC) [18] is arguably the most popular metric for validating posterior approximations in the applied SBI community [3], but it only provides necessary conditions for consistency and fails to give any insight on the local behavior of the estimator. Moreover, it is a univariate procedure, thus ignoring any information about the coupling between parameters. Zhao et al. [19] recently proposed local coverage tests (LCT). They leverage machine learning to evaluate their test quantities on different locations of the feature space (i.e. x-space). In 1D, the probability integral transform (PIT) provides necessary and sufficient conditions for consistency, which is not the case for the proposed multivariate extensions (e.g. HPD).

In this work, we present a multivariate version of LCT for density models based on normalizing flows. Our method comes with the same theoretical guarantees on local consistency and interpretable diagnostics as provided by [19] in 1D. We also present a workflow that can be used as a practical user guide. Lastly, we provide numerical illustrations on a well known model from computational neuroscience [8], and demonstrate the importance of trustworthy diagnostics for multivariate posterior distributions when correlations between parameter variables can play an important role.

## 2 Methods

Our conditional density estimator $q_\phi$ is a normalizing flow defined for samples $\theta \in \mathbb{R}^m$ and is conditioned on observations $x \in \mathbb{R}^d$. It uses a Gaussian base distribution $p(z) = \mathcal{N}(0, I_m)$ and a bijective transformation defined for every $x$, $T_\phi(\cdot; x) := (T_{\phi,1}(\cdot; x), \dots, T_{\phi,m}(\cdot; x))$, with Jacobian $J_{T_\phi}(\cdot; x)$, such that $q_\phi(\theta | x) = p(T_\phi(z; x)) = p(z) | \det J_{T_\phi}(z; x)|^{-1}$.

Our goal is to evaluate the local consistency [9] of $q_\phi$ with respect to the true posterior density, i.e. whether for a given $x$ the following null hypothesis holds:

$$ \mathcal{H}_0(x) : q_\phi(\theta | x) = p(\theta | x), \quad \forall \theta \in \mathbb{R}^m. \qquad (1) $$

We define the multivariate probability integral transform (PIT) of $\theta$ at $x$ and associated to $q_\phi$ as the vector of $m$ one-dimensional projections:

$$ \text{PIT}_m(\theta, x, q_\phi) = [P_1(\theta, x), \dots, P_m(\theta, x)], \quad P_i(\cdot, x) = F_{\mathcal{N}(0,1)} \circ T_{\phi,i}^{-1}(\cdot; x), \quad \forall i \in [1, m] \quad (2) $$

where $F_{\mathcal{N}(0,1)}$ is the cumulative distribution function (c.d.f.) of the univariate normal distribution.

**Theorem 1: Local Consistency and multivariate PIT.** For any $x \in \mathbb{R}^d$, the null hypothesis $\mathcal{H}_0(x)$ holds if, and only if, the covariates of $\text{PIT}_m(\theta, x, q_\phi)$ conditioned on $x$ are mutually independent and uniformly distributed over $(0, 1)$. We refer to Appendix A for a detailed proof of this result.

$$ \begin{align*} p(\boldsymbol{\theta} | x) = q_{\phi}(\boldsymbol{\theta} | x) &\iff p(T_{\phi}^{-1}(\boldsymbol{\theta}, x) | x) = p(z) \\ &\iff p(T_{\phi,1}^{-1}(\boldsymbol{\theta}, x), \dots, T_{\phi,m}^{-1}(\boldsymbol{\theta}, x) | x) = \mathcal{N}(0, I_m) \\ &\iff \begin{cases} p(P_i(\boldsymbol{\theta}, x) | x) = \mathcal{U}(0, 1) & \forall i \in [1, m] \text{ and} \\ \{P_i(\boldsymbol{\theta}, x) | x\}_{i=1,\dots,m} \text{ are mutually independent} \end{cases} \end{align*} $$

We can now verify the null hypothesis of local consistency $\mathcal{H}_0(x)$ (1) via $m$ statistical tests for the uniformity of the 1D local PIT covariates and an additional test for their mutual independence.

**Multivariate Local Coverage Tests (LCT).** Noting that for every $i = 1, \dots, m$ we have

$$ P_i(\theta, x) | x \sim \mathcal{U}(0, 1) \iff \forall \alpha \in [0, 1] \quad r_{i,\alpha}(x) = \mathbb{P}(P_i(\theta, x) \le \alpha | x) = \alpha \quad (3) $$

and following the same approach as in [19], we propose $m$ test statistics

$$ \mathcal{T}_i(x) := \frac{1}{|G|} \sum_{\alpha \in G} \left( \hat{r}_{i,\alpha}(x) - \alpha \right)^2 \quad \forall i \in [1, m], \qquad (4) $$

where $G$ is a grid of $\alpha$-values and the estimators $\hat{r}_{i,\alpha}$ are obtained by regressing $\mathbb{I}_{\{P_i(\theta,x) \le \alpha\}}$ on $x$, which is optimal for $\mathbb{E}[\mathbb{I}_{\{P_i(\theta,x) \le \alpha\}} | x] = r_{i,\alpha}(x)$ (when using appropriate loss-functions [12]). We refer to Algorithm 1 (resp. 2) in [19] for computing the $p$-values (resp. confidence bands) associated to each test. Since there are $m$ independent tests, we recur to a Bonferroni correction of the $p$-values
(resp. confidence levels) [2]. These tests come with interpretable graphical diagnostics such as PP-plots or histograms that depict distributional deviations in 1D, as shown in Figure 1. If any of the uniformity tests is rejected, we reject $H_0(x)$. If not, we proceed to the mutual independence test.

We now assume that the covariates of PIT$_m$ conditioned on $x$ are uniformly distributed over $(0, 1)$, which also means that every coordinate of the flow-transformation $T_{\phi}^{-1}(\theta, x)$ follows a normal distribution given $x$ (cf. Theorem 1). Their mutual independence is thus characterized by a covariance matrix equal to the identity $I_m$. We are currently working on how to perform this check in practice.

**The workflow in practice.** Let $\mathcal{D} = \{x_n, \theta_n\}_{n=1}^N$ be a calibration dataset with $(x_n, \theta_n) \sim p(x, \theta)$ which were not used to train $q_{\phi}$. We use $\mathcal{D}$ to calculate the PIT values PIT$_m(\theta_n, x_n, q_{\phi})$ (2) and estimate our local test quantities as functions of $x$. We investigate the consistency of $q_{\phi}$ in two parts:

(1) **Global consistency check:** first, we look at the global PIT-distribution, i.e. on average over the entire $x$-space. More specifically, we directly compute the empirical approximation of

$$r_{i,\alpha} = \mathbb{P}(P_i(\theta, x) \le \alpha) \approx \frac{1}{N} \sum_{n=1}^{N} \mathbb{I}_{\{P_i(\theta_n, x_n) \le \alpha\}}, \quad \forall i \in [1, m]$$

with samples from $\mathcal{D}$ to test the global uniformity of each PIT-covariate and check which one(s) might be responsible for making $q_{\phi}$ deviate from the true posterior distribution. Note that passing such global test is only a necessary condition for consistency of $q_{\phi}$, as it is insensitive to covariate transformations in $x$-space [19] (and ignores the condition on mutual independence).

(2) **Local consistency check:** we construct $m \times |G|$ transformed datasets from $\mathcal{D}$

$$\mathcal{D}_{i,\alpha} = \{(x_n, W_n^{i,\alpha})\}_{n=1}^N, \quad \alpha \in G, \quad \forall i \in [1, m]$$

where $W_n^{i,\alpha} = \mathbb{I}_{\{P_i(\theta_n, x_n) \le \alpha\}}$. We can then compute the test statistics defined in (4) for any new observation $x$ and check whether $H_0(x)$ should be rejected or not. In the latter case we proceed to the mutual independence test. Only then can we conclude on the validity of $H_0(x)$ according to Theorem 1. If the check in (1) passes, this allows us to guarantee (or reject) consistency anywhere in $x$-space. Even if the check in (1) does not pass, analyzing local consistency allows to 'open the box' and better understand why and where in $x$-space the estimator fails. In such situations – as the goal is not to guarantee local consistency – it can be enough to test for the covariate-wise uniformity of PIT$_m$, putting the test for mutual independence aside. This approach was adopted for our numerical illustrations in Section 3.

## 3 Numerical illustrations

We apply our method to check the validity of a posterior estimate $q_{\phi}$ of the Jansen and Rit neural mass model (JRNMM) [8]. This well known model from computational neuroscience takes parameters $\theta = (C, \mu, \sigma, g) \in \mathbb{R}^4$ as input and generates time series $x \in \mathbb{R}^{1024}$ with properties similar to brain signals obtained in neurophysiology. Parameter $C$ influences the oscillatory behavior of the signals and $(\mu, \sigma)$ characterize their amplitude. The gain factor $g$ rescales the signals and models the effects of the amplifier used for measuring them in practice. Note that the coupling-effect of $g$ and $(\mu, \sigma)$ on the amplitude leads to intrinsic indeterminacies in the inversion of the model [16]. Our approximation $q_{\phi}$ is a conditioned masked autoregressive flow (MAF) [13] with 10 layers and implemented in the sbi package [20]. Based on the setup in [16], we train $q_{\phi}$ on 50 000 simulations from the JRNMM with a uniform prior defined on physiologically relevant values (cf. pair-plots in Figure 2).

**Results and discussion.** All test quantities are computed with a calibration dataset $\mathcal{D}$ containing 10 000 simulations and a grid of $|G| = 100$ $\alpha$-values in $[0, 1]$. Following the workflow described in Section 2, we first check the global consistency of $q_{\phi}$. In the left part of Figure 1, our graphical diagnostics illustrate how the c.d.f. for every PIT-covariate deviates from the identity function (black dashed line), outside of the 95%-confidence region (in gray), thus rejecting the null hypothesis of global consistency. We compare our method to SBC (right part of Figure 1) as implemented in the sbi package. We observe that SBC is unable to detect any inconsistencies in $q_{\phi}$.

We proceed to investigate the local consistency of $q_{\phi}$ on different locations of $x$-space. Note that since the global test for consistency has not passed, we may focus on the uniformity of the local PIT-covariates (cf. (2) of the workflow in Section 2). We consider a 1D subspace in $\theta$-space, where $(C_0, \mu_0, \sigma_0)$ are fixed and the gain $g_0$ varies in $[-20, 20]$ to generate observations $x_0$. The upper-left part of Figure 2 shows how the test statistics evolve with $g_0$. We observe a strong deviation from
<figure><img src="image_2.png" /><figcaption>Figure 1: Global Consistency of the JRNMM posterior estimator: PP-plots and corresponding histograms obtained for every global PIT covariate $P_i(\theta)$ (right) and for SBC applied to every element $\theta_i$ of $\theta$ (left). The gray zone indicates the 95%-confidence region of acceptance outside of which the uniformity test is rejected. SBC fails to reject the null hypothesis of global consistency.</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>Figure 2: Local Consistency Analysis. (Left) Diagnostics obtained for every PIT-covariate, when evaluated for different $x_0$, simulated via the JRNMM with fixed parameters ($C_0, \mu_0, \sigma_0$) and variable $g_0$ in $[-20, 20]$. The test-statistics are plotted as a function of $g_0$ and indicate inconsistent behavior (values are not constantly close to zero). Below, the PP-plots report the nature (bias/dispersion) of these inconsistencies for $g_0 \in \{-20, 0, 20\}$. The gray zone indicates the 95%-confidence region of acceptance outside of which the uniformity test is rejected. (Right) Pair-plots of $q_\phi(\theta | x_0)$ with ground-truth parameters $\theta_0 = (C_0, \mu_0, \sigma_0, g_0)$ (black dots and dashed lines) for $g_0 \in \{-20, 0, 20\}$.</figcaption></figure>

uniformity for covariates $P_2$ and $P_3$: they vary smoothly in a 'U-shape', with higher values as $g_0$ deviates from zero. We also generate local PP-plots (lower-left part of Figure 2) and observe positive (resp. negative) bias for small (resp. high) values of $g_0$. Our procedure shows that there are certain locations in $x$-space (here $g_0 = 0$) where $q_\phi$ performs well (i.e. test statistics are close to zero and PP-plots show little deviations from the 95%-confidence region in gray), even though global consistency does not hold. These observations can be compared to the pair-plots (right part of Figure 2) which shows the estimated posterior density over the entire $\theta$-space. We observe that gain values at the boundary of the prior support ($g_0 = -20$ and $g_0 = 20$) induce discrepancies (cf. marginals plotted in the diagonals of the pair-plots in Figure 2) that are detected by our PIT-based diagnostics.

The above results are obtained by directly applying the method from [19] using their default regressor, *MLPClassifier* from scikit-learn [15]. Indeed, this model is well suited for the binary target variables $W^{i,\alpha}$ and, based on neural networks, is able to scale to high-dimensional $x$-spaces. Furthermore, it has been shown [12] that for a large enough dataset $\mathcal{D}$ (and number of training iterations), the cross-entropy loss converges to the optimal solution of our regression problem. Note that although such method allows to capture and visualize local discrepancies, it does not cope well with our high-dimensional, complex data-distributions (high variances represented by large confidence regions computed over 100 trials). Also, it requires the training of a large number ($m \times |G|$) of regressors (one for each $\theta$-dimension and $\alpha$-value). Finally, this method is only applicable for the uniformity tests. We are currently investigating algorithms that are computationally more efficient and more accurate in estimating the *full* multivariate local PIT-distribution.
**Conclusions.** Numerical illustrations demonstrate that our diagnostics capture well the inconsistencies in $q_\phi$ with respect to the true JRNMM posterior. To be precise, the proposed validation method is statistically more powerful and computationally more efficient than SBC and, importantly, it allows for a local analysis that reveals where in x-space the estimations should be improved. Our method exploits useful properties of modern normalizing flows. Indeed, contrarily to other SBI-strategies [7], they allow for efficient density evaluation and use bijective transformations involving the mixing between elements of $\theta$. These transformations can be interpreted as a set of 1D-projections of $\theta$, including information about its joint p.d.f. (i.e. interactions between its elements $\theta_i$). The covariates of our multivariate PIT therefore define fast and easy to compute 1D test-quantities on which we can perform univariate LCTs that do not completely ignore correlations in the $\theta$-space. Combined with a test to check their mutual independence, this would provide theoretical guarantees for consistency, which is not the case for HPD, the existing multivariate version of LCT proposed in [19].

Ongoing work will investigate the link between PIT-covariates and their deviations from uniformity with the actual parameters of the model. This could reveal the true nature of their coupling which can be used for the development of better specified models (e.g. hierarchical posterior estimation [16, 17]). Ultimately, the goal is to use our diagnostics as a tool for model selection.

## Broader Impact

This work tackles an important open question in simulation based inference. We introduce theoretically valid and interpretable validation diagnostics that scale to both high-dimensional data and parameter spaces. Our contribution should help to further improve SBI methods and drive the design of better specified models, hence allowing to build up trust on their ability to address important questions in experimental science. Moreover, our work is in line with many other tentatives in the machine learning community of ensuring the quality and calibration of complex models based on neural networks.

## Acknowledgments

Julia Linhart is recipient of the Pierre-Aguilar Scholarship and thankful for the funding of the Capital Fund Management (CFM). Alexandre Gramfort thanks the support of the ANR BrAIN (ANR-20-CHIA0016) grant.

## References

[1] Sean R. Bittner, Agostina Palmigiano, Alex T. Piet, Chunyu A. Duan, Carlos D. Brody, Kenneth D. Miller, and John Cunningham. Interrogating theoretical models of neural computation with emergent property inference. *eLife*, 10, 7 2021. ISSN 2050084X. doi: 10.7554/eLife.56265.

[2] J. Martin Bland and Douglas G. Altman. Multiple significance tests: the bonferroni method. *BMJ (Clinical research ed)*, 310:170, 1 1995. ISSN 0959-8138. doi: 10.1136/BMJ.310.6973.170.

[3] Jan Boelts, Jan-Matthis Lueckmann, Richard Gao, and Jakob H. Macke. Flexible and efficient simulation-based inference for models of decision-making. *bioRxiv*, page 2021.12.22.473472, 6 2022. ISSN 2050084X. doi: 10.1101/2021.12.22.473472.

[4] Kyle Cranmer, Johann Brehmer, and Gilles Louppe. The frontier of simulation-based inference. *Proceedings of the National Academy of Sciences (PNAS)*, 117:30055-30062, 2020. ISSN 0027-8424. doi: 10.1073/pnas.1912789117.

[5] Niccolò Dalmasso, David Zhao, Rafael Izbicki, and Ann B. Lee. Likelihood-free frequentist inference: Confidence sets with correct conditional coverage. 7 2021. doi: 10.48550/arxiv.2107.03920.

[6] Pedro J. Gonçalves, Jan Matthis Lueckmann, Michael Deistler, Marcel Nonnenmacher, Kaan Öcal, Giacomo Bassetto, Chaitanya Chintaluri, William F. Podlaski, Sara A. Haddad, Tim P. Vogels, David S. Greenberg, and Jakob H. Macke. Training deep neural density estimators to identify mechanistic models of neural dynamics. *eLife*, 9:1-46, 9 2020. ISSN 2050084X. doi: 10.7554/ELIFE.56261.
[7] Joeri Hermans, Volodimir Begy, and Gilles Louppe. Likelihood-free MCMC with amortized approximate ratio estimators. In Hal Daumé III and Aarti Singh, editors, *Proceedings of the 37th International Conference on Machine Learning*, volume 119 of *Proceedings of Machine Learning Research*, pages 4239–4248. PMLR, 13–18 Jul 2020. doi: 10.48550/arxiv.1903.04057.

[8] Ben H. Jansen and Vincent G. Rit. Electroencephalogram and visual evoked potential generation in a mathematical model of coupled cortical columns. *Biological Cybernetics* 1995 73:4, 73:357–366, 9 1995. ISSN 1432-0770. doi: 10.1007/BF00199471.

[9] Ilmun Kim, Ann B. Lee, and Jing Lei. Global and local two-sample tests via regression. *Electronic Journal of Statistics*, 13:5253–5305, 12 2018. ISSN 19357524. doi: 10.48550/arxiv.1812.08927.

[10] Jan-Matthis Lueckmann, Jan Boelts, David S Greenberg, Pedro J Gonçalves, and Jakob H Macke. Benchmarking simulation-based inference. *Proceedings of The 24th International Conference on Artificial Intelligence and Statistics (PMLR)*, 130:343–351, 4 2021. doi: 10.48550/arxiv.2101.04653.

[11] Luca Masserano, Tommaso Dorigo, Rafael Izbicki, Mikael Kuusela, and Ann B. Lee. Simulation-based inference with waldo: Perfectly calibrated confidence regions using any prediction or posterior estimation algorithm. 5 2022. doi: 10.48550/arxiv.2205.15680.

[12] John W. Miller, Rod Goodman, and Padhraic Smyth. On loss functions which minimizes conditional expected values and posterior probabilities. *IEEE TRANSACTIONS ON INFORMAITON THEORY*, 39, 1993.

[13] George Papamakarios, Theo Pavlakou, and Iain Murray. Masked autoregressive flow for density estimation. *Advances in Neural Information Processing Systems (NeurIPS)*, pages 2339–2348, 12 2017. ISSN 10495258. doi: 10.48550/arxiv.1705.07057.

[14] George Papamakarios, Eric Nalisnick, Danilo Jimenez Rezende, Shakir Mohamed, and Balaji Lakshminarayanan. Normalizing flows for probabilistic modeling and inference. *Journal of Machine Learning Research*, 22:1–64, 2021. ISSN 15337928. doi: 10.48550/arxiv.1912.02762.

[15] F. Pedregosa, G. Varoquaux, A. Gramfort, V. Michel, B. Thirion, O. Grisel, M. Blondel, P. Prettenhofer, R. Weiss, V. Dubourg, J. Vanderplas, A. Passos, D. Cournapeau, M. Brucher, M. Perrot, and E. Duchesnay. Scikit-learn: Machine learning in Python. *Journal of Machine Learning Research*, 12:2825–2830, 2011.

[16] Pedro L C Rodrigues, Thomas Moreau, and Gilles Louppe. Hnpe: Leveraging global parameters for neural posterior estimation. *Neural Information Proceedings (NeurIPS)*, 12 2021. ISSN 10495258. doi: 10.48550/arxiv.2102.06477.

[17] Louis Rouillard and Demian Wassermann. ADAVI: Automatic Dual Amortized Variational Inference Applied To Pyramidal Bayesian Models. In *ICLR 2022*, Virtual, France, April 2022. doi: 10.48550/arXiv.2106.12248.

[18] Sean Talts, Michael Betancourt, Daniel Simpson, Aki Vehtari, and Andrew Gelman. Validating bayesian inference algorithms with simulation-based calibration. 4 2018. doi: 10.48550/arxiv.1804.06788.

[19] David Zhao, Niccolò Dalmasso, Rafael Izbicki, and Ann B. Lee. Diagnostics for conditional density models and bayesian inference algorithms. *PMLR*, 2 2021. ISSN 2640-3498.

[20] Álvaro Tejero-Cantero, Jan Boelts, Michael Deistler, Jan-Matthis Lueckmann, Conor Durkan, Pedro J Gonçalves, David S Greenberg, Jakob H Macke, Computational Neuroengineering, and TU Munich. Sbi – a toolkit for simulation-based inference. *Journal of Open Source Software*, 5:2505, 7 2020. ISSN 2475-9066. doi: 10.21105/joss.02505.
A Appendix

**Proof of Theorem 1.** Let $x \in \mathbb{R}^d$. We consider the random variables defined on $\mathbb{R}^m$: $\Theta \sim p(\theta | x)$ and $\Theta_{\phi} \sim q_{\phi}(\theta | x)$. Our local null hypothesis (cf. Equation 1), can be rewritten as

$$ \mathcal{H}_0(x) : \mathrm{Prob}(\Theta \in \Omega_{\theta}) = \mathrm{Prob}(\Theta_{\phi} \in \Omega_{\theta}) \quad \forall \Omega_{\theta} \subset \mathbb{R}^m \quad (5) $$

From the definition of our normalizing flow $q_{\phi}$ with bijective transformation $T_{\phi}$ and normal base distribution, there exists $Z \sim p(z) = \mathcal{N}(0, \mathbf{I}_m)$, such that $\Theta_{\phi} = T_{\phi}(Z, x)$.

Replacing $\Theta_{\phi}$ by $T_{\phi}(Z, x)$ in (5) and defining $\Omega_z = T_{\phi}^{-1}(\Omega_{\theta}; x)$, we have

$$ \mathrm{Prob}(\Theta \in \Omega_{\theta}) = \mathrm{Prob}(T_{\phi}(Z; x) \in \Omega_{\theta}) \iff \mathrm{Prob}(T_{\phi}^{-1}(\Theta; x) \in \Omega_{z}) = \mathrm{Prob}(Z \in \Omega_{z}) $$

Considering $T_{\phi,i}^{-1}$, the i-th coordinate of $T_{\phi}^{-1}$, and $Z_1, \dots, Z_m$ the mutually independent and normally distributed covariates of $Z \sim \mathcal{N}(0, \mathbf{I}_m)$, we can rewrite the right part of the above equivalence as:

$$ \mathrm{Prob}(T_{\phi,1}^{-1}(\Theta; x) \in \Omega_{z,1} \cap \dots \cap T_{\phi,m}^{-1}(\Theta; x) \in \Omega_{z,m}) = \prod_{i=1}^{m} \mathrm{Prob}(Z_i \in \Omega_{z,i}) \quad (6) $$

where the $\Omega_{z,i}$ are projections of the set $\Omega_z$ on each of its dimensions.

Note that since the above equivalence is true for any choice of $\Omega_{\theta}$, hence of $\Omega_z$, we can write what happens when we fix a given $\Omega_{z,1}$ and let $\Omega_{z,2} = \dots = \Omega_{z,m} = \mathbb{R}$:

$$ \mathrm{Prob}(T_{\phi,1}^{-1}(\Theta; x) \in \Omega_{z,1}) = \mathrm{Prob}(Z_1 \in \Omega_{z,1}) $$

Doing the same for all other coordinates and using these equalities in (6), we end up with

$$ \forall i \in [1, m], \quad \forall \Omega_{z,i} \subset \mathbb{R}, \quad \mathrm{Prob}(T_{\phi,i}^{-1}(\Theta; x) \in \Omega_{z,i}) = \mathrm{Prob}(Z_i \in \Omega_{z,i}) \quad (7) $$

$$ \text{and} \quad \mathrm{Prob}(T_{\phi,1}^{-1}(\Theta; x) \in \Omega_{z,1} \cap \dots \cap T_{\phi,m}^{-1}(\Theta; x) \in \Omega_{z,m}) = \prod_{i=1}^{m} \mathrm{Prob}(T_{\phi,i}^{-1}(\Theta; x) \in \Omega_{z,i}) \quad (8) $$

where (8) is the definition of mutual independence itself.

We now apply the c.d.f. $F_{\mathcal{N}(0,1)}$ to each random variable in (7) and (8). Mutual independence stays true and the probability integral transform theorem (in 1D) states that $F_{\mathcal{N}(0,1)}(Z_i) \sim \mathcal{U}(0, 1)$. We therefore get that the null hypothesis $\mathcal{H}_0(x)$ holds if, and only if,

$$ \left\{ \begin{array}{l} P_i(\Theta, x) = F_{\mathcal{N}(0,1)}(T_{\phi,i}^{-1}(\Theta; x)) \sim \mathcal{U}(0, 1), \quad \forall i = [1, m] \quad \text{and} \\ \{P_i(\Theta, x)\}_{i=1,\dots,m} \text{ are mutually independent} \end{array} \right. \quad (9) $$

The result in Theorem 1 directly follows from rewriting (5) and (9) with initial notations from Section 2: remember that $\Theta \sim p(\theta | x)$, so $P_i(\Theta, x) \sim p(P_i(\theta, x) | x)$ and we get:

$$ p(\theta | x) = q_{\phi}(\theta | x) \iff \left\{ \begin{array}{l} p(P_i(\theta, x) | x) = \mathcal{U}(0, 1) \quad \forall i \in [1, m] \quad \text{and} \\ \{P_i(\theta, x) | x\}_{i=1,\dots,m} \text{ are mutually independent} \end{array} \right. $$

where $P_i(\theta, x)$ is the $i^{th}$ covariate of PIT$_m(\theta, x, q_{\phi})$, the multivariate PIT of $\theta$ at $x$, associated to $q_{\phi}$.

**Conclusion:** The null hypothesis $\mathcal{H}_0(x)$ holds if, and only if, the covariates of PIT$_m$ conditioned on $x$ are mutually independent and uniformly distributed over $(0, 1)$.