\section{Introduction}

Neural networks have achieved great success in solving diverse problems ranging from image analysis \cite{NIPS2012_4824,NIPS2013_5207}, to natural language processing \cite{cho-etal-2014-learning,NIPS2017_transformer}, to medical applications \cite{annurev-bioeng-071516-044442,MCBEE20181472}.

 An important theoretical as well as computational challenge is calculating the distribution of outputs of a neural network from a given distribution of inputs.
 This task is a crucial ingredient for uncertainty quantification, explainable AI and Bayesian learning.

 More precisely, the task is to describe the distribution $\mu_{\text{out}}$ of
 \begin{equation*}
    f(\mathbf{x})\quad\text{with}\quad\mathbf{x}\sim\mu_{\text{in}},
 \end{equation*}
 where $f$ is a neural network function and $\mathbf{x}$ is a random input vector distributed according to some probability measure $\mu_{\text{in}}$.
 In other words, we want to determine the push-forward
 \begin{equation}\label{eq:prob_prop}
 \mu_{\text{out}} = f_\ast \mu_{\text{in}}
 \end{equation}
 of $\mu_{\text{in}}$ with respect to the transformation $f$.
 Given the fact that neural networks are highly non-linear functions this generally has no closed form solution.
 Obtaining approximations through numerical integration methods is expensive due to the high dimensionality and complexity of the function.

 This lead to the application of other approximation schemes, such as Assumed Density Filtering \cite{gast-roth-8578453, mwhk-2019-rate-dist} (ADF) and in the case of Bayesian neural networks the more general framework of expectation propagation \cite{jylanki2014expectation, soudry2014expectation} (EP) which contains ADF as a special case. Both methods were originally developed for general large-scale Bayesian inference problems \cite{Minka:2001:FAA:935427, Boyen:1998:TIC:2074094.2074099}.

 Assumed Density Filtering makes use of the feed-forward network structure: If we write $f=f_L\circ f_{L-1}\circ\dots \circ f_2\circ f_1$ as a composition of its layers $f_j$, then \eqref{eq:prob_prop} decomposes as
  \begin{alignat}{3}
    \mu_{0} &= \mu_{\text{in}}, \nonumber\\
    \mu_{j} &= (f_j)_\ast\mu_{j-1} &\quad&\text{for}\quad j=1,\dots,L, \nonumber \\
    \mu_{\text{out}} &= \mu_{L}.\label{eq:layerwise}
 \end{alignat}

 The task can then be described as how to propagate distributions through neural network layers.

 If each $\mu_{j}$ was in the same family of probability distributions as the original $\mu_{\text{in}}$ then we could find efficient layer-wise propagation rules for each layer $f_j$ to obtain $\mu_{\text{out}}$. Such a family, if expressable with a possibly large but finite number of parameters, would be extremely useful for applications involving uncertainty quantification or the need for explainable predictions.

 ADF commonly propagates Gaussians, or more generally exponential families, which are not invariant under the action of neural network layers with non-linear activation. Hence, the method relies on a projection step onto the chosen family of probability distributions after each network layer.
 In other words, \eqref{eq:layerwise} is replaced by
  \begin{alignat*}{3}
    \mu_{0} &= \mathrm{proj}(\mu_{\text{in}}), \\
    \mu_{j} &= \mathrm{proj}((f_j)_\ast\mu_{j-1}) &\quad&\text{for}\quad j=1,\dots,L, \nonumber \\
    \mu_{\text{out}} &= \mu_{L},
 \end{alignat*}
  where $\mathrm{proj}(\cdot)$ is a suitable projection.  In this case only $\mu_{\text{out}}\approx f_{\ast}\mu_{\text{in}}$ holds and there is no guarantee how good the approximation is after multiple layers. This heuristic is justified by the implicit assumption that no invariant families exist.

 Due to the flexibility and richness of the class of neural network functions, it seems intuitively clear that such an invariant family of distributions cannot be realised in a meaningful way---a fact that, to the best of our knowledge, remained unproven. In this paper we will fill this gap.

 \subsection{Contribution}
 We give a complete characterisation of families of probability distributions that are invariant under the action of ReLU neural network layers.
 In fact, we prove that all possible invariant families belong to a set of degenerate cases that are either practically irrelevant or amount to sampling.
 Hence, the implicit assumption made by the ADF and EP frameworks is justified. For each of the degenerate cases we give an explicit example of an invariant family of distributions.

 Our novel proof technique is based on very intuitive geometric constructions and ideas from metric dimension theory. This intuitive argument is made rigorous using properties of the Hausdorff dimension of metric spaces.
 As a result, we provide a theoretical underpinning for statistical inference schemes such as ADF for neural networks.

 \subsection{Applications}

  Calculating output distributions of neural networks is a basic and central problem to make neural networks more applicable. Good performance on test data alone is not sufficient to ensure reliability.
  High-stakes applications such as, e.g., medical imaging or autonomous driving, make questions regarding uncertainty quantification (UQ), robustness, and explainability of network decisions inescapable \cite{pmlr-v136-oala20a, begoli2019need,zhang2018opening,holzinger2017we,michelmore2018evaluating,kim2017interpretable}.

 \paragraph{Uncertainty Quantification}
 Studying the propagation of input uncertainties through an already trained (deterministic) network is used to analyse the reliability of their predictions and their sensitivity to changes in the data. ADF is used towards this goal, e.g., for sigmoid networks for automatic speech recognition \cite{Astudillo2011propagation,Astudillo2014AccountingFT} and ReLU networks for image classification and optical flow \cite{gast-roth-8578453}.

 \paragraph{Explainable Artificial Intelligence}
 So far, neural networks are mostly considered opaque ``black-box'' models. Recent approaches to make them human-interpretable have focused on the crucial subtask of relevance attribution, i.e., finding the most relevant input variables for a given prediction.
 A probabilistic formulation of this task was recently proposed for classification networks \cite{fong2017interpretable, wmhk-2021-jair-compl,mwhk-2019-rate-dist}. Macdonald et al. also considered a heuristic solution strategy for deep ReLU networks that is based on ADF with Gaussian distributions \cite{mwhk-2019-rate-dist}.

  \paragraph{Bayesian Networks}
 Instead of point estimates, Bayesian neural networks (BNN) learn a probability distribution over the network weights through Bayesian inference.
 BNNs provide an inbuilt regularisation by choosing appropriate priors. They naturally lead to network compression when using sparsity-promoting priors or encoding weights of high variance with less precision.
 ADF with Gaussian distributions has been employed for the fast calculation of posteriors during the training of BNNs with sigmoidal \cite{opper1998bayesian} and ReLU \cite{hernandez2015probbackprop} activations.
  A similar approach was taken in \cite{pmlr-v28-wang13a} which relies on layer-wise Gaussian approximations to yield a fast alternative to dropout regularized training.

 \subsection{Outline of the Paper}
 In \cref{sec:prelim_and_result} we state our two main results and provide a high-level description of the proof ideas. In \cref{sec:mainproof} we show that no (non-degenerate) network invariant families of distributions exist. Two constructive steps of the proof are deferred to \cref{sec:xi_surjective,sec:xi_loclip} respectively. Finally, in \cref{sec:counterproof} we prove our second main result by constructing several degenerate invariant families. We conclude with a short discussion in \cref{sec:discussion}.

 \section{Characterisation of Invariant Families of Distributions}\label{sec:prelim_and_result}
 Before we rigorously state and prove our main result, we start by introducing some notation and terminology and precisely specify what we mean by \emph{parametrised families} of probability distributions and by their \emph{invariance} with respect to layers of neural networks.

 \subsection{Preliminaries}
 We denote by $\mathcal{B}(\mathbb{R}^d)$ the Borel $\sigma$-Algebra on $\mathbb{R}^d$ and by $\mathcal{D}(\mathbb{R}^d)$ the set of Radon\footnote{In fact $\mathbb{R}^d$ is a separable complete metric space, thus every Borel probability measure is automatically Radon.} Borel probability measures on $\mathbb{R}^d$. We equip $\mathcal{D}(\mathbb{R}^d)$ with the Prokhorov metric
 \begin{equation*}
     d_P(\mu,\nu) = \inf\{\epsilon >0\,:\, \mu(B) \leq \nu(B^\epsilon) + \epsilon
     \text{ and } \nu(B) \leq \mu(B^\epsilon) + \epsilon \text{ for any }B\in\mathcal{B}(\mathbb{R}^d)\},
 \end{equation*}
 where $B^\epsilon = \{\,x\in\mathbb{R}^d\,:\, \exists y\in B \text{ with } \|x-y\|<\epsilon\,\}$ is
 the open $\epsilon$-neighbourhood of a Borel set $B\in\mathcal{B}(\mathbb{R}^n)$.

 \begin{definition}[Family of Distributions]
    Let $\Omega\subseteq\mathbb{R}^n$ and $p\colon\Omega\to\mathcal{D}(\mathbb{R}^d)$. We call $p$ an $n$-parameter family of probability distributions. It is called a continuous family, if $p$ is continuous.
 \end{definition}
 It might seem more intuitive to refer to the set $\{p(\theta)\}_{\theta\in\Omega}$ as the family of distributions. However, a set of distributions can be parametrised in multiple ways and even the number $n$ of parameters describing such a set is not unique. Whenever we speak about a family of distributions, we never think of it as a mere set of distributions but always attach it to a fixed chosen parametrisation $p$. All our results are stated in terms of this parametrisation.
 \begin{example}[$2$-dimensional Gaussian Family]
 Consider the parameter space
 \[\Omega = \mathbb{R}^2\times\skl{(\sigma_1,\sigma_2,\sigma_3) \in \mathbb{R}^3\,:\,\sigma_1,\sigma_2\geq 0\text{ and }\sigma_1\sigma_2-\sigma_3^2 \geq 0},\]
 and
 \[
  p\colon \Omega \to \mathcal{D}(\mathbb{R}^2)\colon (\mu_1,\mu_2, \sigma_1, \sigma_2, \sigma_3) \mapsto \mathcal{N}\kl{\bmat{\mu_1 \\ \mu_2}, \bmat{\sigma_1 & \sigma_3 \\ \sigma_3 & \sigma_2}},
 \]
 where $\mathcal{N}(\bm{\mu},\bm{\Sigma})$ denotes the Gaussian distribution with mean $\bm{\mu}$ and covariance matrix $\bm{\Sigma}$.
 Then $p$ is a continuous 5-parameter family of 2-dimensional Gaussian distributions.
 \end{example}

 \begin{definition}[Invariance]
    Let $p\colon\Omega\to\mathcal{D}(\mathbb{R}^d)$ be a family of probability distributions and $f\colon\mathbb{R}^d\to\mathbb{R}^d$ any measurable function. Then the family is called $f$-invariant, if for any $\theta\in\Omega$ the pushforward of the measure $p(\theta)$ under $f$ is again in the family, i.e.\ there exists $\omega\in\Omega$ such that $p(\omega)=f_\ast p(\theta)$. For a collection $\mathcal{F}\subseteq \{\,f\colon\mathbb{R}^d\to\mathbb{R}^d\,:\,f\text{ measurable}\,\}$ of measurable functions it is called $\mathcal{F}$-invariant, if it is $f$-invariant for all $f\in\mathcal{F}$.
 \end{definition}
 For us the special case of invariance with respect to layers of neural networks will be of interest. A neural network is a composition of (possibly many) building blocks, each consisting of an affine transform and a simple componentwise non-linear transform referred to as activation function. These building blocks are also called the layers of a neural network. One frequently used activation function is the rectified linear unit (ReLU) defined as $\varrho(x)=\max\{0,x\}$. Thus, one layer of a ReLU neural network has the general form $\mathbf{x}\mapsto\varrho(\mathbf{W}\mathbf{x}+\mathbf{b})$,
 with so called weight matrix $\mathbf{W}$ and bias vector $\mathbf{b}$. An $L$-layer ReLU network is the composition of $L$ such functions and has the form
 \[
 \mathbf{x} \mapsto \varrho\kl{\mathbf{W}_L\varrho\kl{\mathbf{W}_{L-1}\varrho\kl{\dots\varrho\kl{\mathbf{W}_1\mathbf{x}+\mathbf{b}_1}\dots} + \mathbf{b}_{L-1}} +\mathbf{b}_L}.
 \]
 \begin{definition}[ReLU-Invariance]
     Let $p\colon\Omega\to\mathcal{D}(\mathbb{R}^d)$ be a family of distributions. The family is called ReLU-invariant, if it is $\mathcal{F}$-invariant for the collection
     \[
     \mathcal{F} = \{\,f\colon\mathbb{R}^d\to\mathbb{R}^d\colon\mathbf{x}\mapsto\varrho(\mathbf{W}\mathbf{x}+\mathbf{b})\,:\,\mathbf{W}\in\mathbb{R}^{d\times d}, \mathbf{b}\in\mathbb{R}^d\,\},
     \]
     where $\varrho(x)=\max\{0,x\}$ is applied componentwise.
 \end{definition}
 We observe that if a family $p$ is $f$-invariant and $g$-invariant for two functions $f$ and $g$ such that the composition $g\circ f$ is well-defined, then it is also $(g\circ f)$-invariant. In particular, a ReLU-invariant family is invariant for all ReLU networks of any depths $L$.

\subsection{The Main Result}
 We can now state the main theorem of this work.
 \begin{restatable}{theorem}{mainresult}\label{thm:main}
 Let $\Omega\subseteq\mathbb{R}^n$ and $p\colon\Omega\to\mathcal{D}(\mathbb{R}^d)$ be a continuous and ReLU-invariant $n$-parameter family of probability distributions. Then at least one of the following restrictions has to hold:
 \begin{description}
     \item[R1. Restricted Dimension:] $d=1$,
     \item[R2. Restricted Support:] $\supp(p(\theta))$ is finite for all $\theta\in\Omega$,
     \item[R3. Restricted Regularity:] $p$ is not locally Lipschitz continuous.
 \end{description}
 \end{restatable}
 This can be interpreted as follows: Besides some rather degenerate cases there can not be any family of probability distributions that is invariant with respect to the layers of a ReLU neural network. The three restrictions characterise which kind of degenerate cases can occur. The first restriction \emph{R1} shows that neural networks with only one neuron per layer are not really powerful function classes, so that the ReLU-invariance is not a strong limitation in this case. However, already two dimensions are enough to unlock the expressive power of neural networks. Restriction \emph{R2} shows that under mild regularity assumptions on the parametrisation the only ReLU-invariant distributions in dimensions $d\geq 2$ are finite mixtures of Dirac distributions. This amounts to sampling distributions of finite size, which of course can work in many scenarios but are often too computationally expensive in high-dimensions. Restriction \emph{R3} shows that continuity alone is not a strong enough assumptions on the parametrisation. This is due to the fact that the class of continuous functions is too flexible and includes non-intuitive examples such as space-filing curves. These can be used to construct examples of invariant distributions. However such a kind of parametrisation has to be be rather wild and would be impractical to work with. A slightly stronger assumption like local Lipschitz continuity is enough to exclude these degenerate examples.

We complement \cref{thm:main} by providing a kind of reverse statement showing that the three restrictions are individually necessary.
 \begin{restatable}{theorem}{maincounter}\label{thm:counter}
    The following statements hold independently of each other.
    \begin{enumerate}[label=(\roman*)]
    \item \label{itm:counter_i} There exists a continuous ReLU-invariant $n$-parameter family of probability distributions that fulfils restriction R1 but not restrictions R2 and R3 from \cref{thm:main}.
    \item \label{itm:counter_ii} There exists a continuous ReLU-invariant $n$-parameter family of probability distributions that fulfils restriction R2 but not restrictions R1 and R3 from \cref{thm:main}.
    \item \label{itm:counter_iii} There exists a continuous ReLU-invariant $n$-parameter family of probability distributions that fulfils restriction R3 but not restrictions R1 and R2 from \cref{thm:main}.
    \end{enumerate}
 \end{restatable}

 \subsection{High-Level Proof Strategy}
  \begin{figure}
    \centering
    \resizebox{\textwidth}{!}{
    \input{./img/highlevel_proof}
    }
    \caption{Main steps of transforming a generic probability distribution to a distribution supported on a polygonal chain. The original distribution (left) is projected onto a single dimension and partitioned into intervals (centre) that each contain a sufficiently large portion of the probability mass (shown as the lightly shaded region). Finally it is ``bent'' into a polygonal chain (right). The number of segments of the chain and its segment lengths can be chosen arbitrarily. }
    \label{fig:highlevel}
 \end{figure}

 The main idea for proving \cref{thm:main} is to use the layers of a neural network to transform simple, more or less arbitrary probability distributions to complicated distributions that need an arbitrarily high number of parameters to be described. But a family of parametrised distributions necessarily has a finite number of parameters leading to a contradiction if the family is assumed to be invariant under neural network layers. More precisely, ReLU neural network layers can be used to transform arbitrary probability distributions to distributions that are supported on certain polygonal chains, which we call arcs. This is visualised in \cref{fig:highlevel}. These arcs can be described by the lengths of their line segments. However, the number of segments can be made larger than the number of parameters describing the family of distributions as long as the support of the initial distribution is large enough. From this a contradiction can be derived. There are only three ways to circumvent this, corresponding to the three restrictions in \cref{thm:main}. Firstly, restricting the allowed neural network layers so that transformations to the arcs are not possible. Secondly, restricting the support of the distributions so that the number of line segments they can be transformed to is limited. Thirdly, using a wild parametrisation that leverages ideas similar to space-filling curves in order to use only few parameters to describe a set that effectively would require more parameters.

 The idea of mapping (a subset of) the parametrisation domain to distributions supported on polygonal arcs, then to the respective arcs, and finally to segment lengths is schematically shown in \cref{fig:overview}. We will give exact definitions of all involved spaces and mappings shortly.
 The full proof for \cref{thm:main} will be given in \cref{sec:mainproof} and the proof for the reverse \cref{thm:counter} in \cref{sec:counterproof}.

 \begin{figure}
    \centering
    \input{./img/overview_improved}
    \caption{Schematic overview of the spaces and functions involved in our proof. }
    \label{fig:overview}
 \end{figure}

\section{Proof of the Main Result (\cref{thm:main})}\label{sec:mainproof}
We will now present the main steps of the proof of \cref{thm:main}. The proofs of the individual steps will be deferred to later sections in order to not distract from the main idea. We will assume towards a contradiction that all three restrictions in \cref{thm:main} are not satisfied.
\begin{assumption}\label{asm:main}
Let $d\geq 2$, $\Omega\subseteq\mathbb{R}^n$, and $p\colon\Omega\to\mathcal{D}(\mathbb{R}^d)$ be a locally Lipschitz continuous and ReLU-invariant $n$-parameter family of probability distributions such that there exists a $\theta\in\Omega$ for which $\supp\kl{p(\theta)}$ is not finite.
\end{assumption}
From this we will derive a contradiction. The proof is based on the idea of transforming arbitrary probability distributions to distributions supported on certain polygonal chains, which we call arcs.
\begin{definition}[Arcs]\label{def:arc}
    An $m$-arc is a subset $A\subseteq \mathbb{R}^2$ defined as a polygonal chain, i.e. a connected piecewise linear curve,
     \[
        A = \bigcup_{i=1}^m \conv\kl{\skl{\mathbf{v}_{i-1}, \mathbf{v}_i}},
    \]
    determined by $m+1$ vertices $\{\mathbf{v}_0, \dots, \mathbf{v}_m\}$, that satisfy
    \[
     \mathbf{v}_i = \mathbf{v}_{i-1} + r_i \begin{bmatrix}\sin(i \phi_m)\\ \cos(i \phi_m)\end{bmatrix}, \quad i=1,\dots,m
    \]
    for some length scales $r_1,\dots, r_m > 0$ and the angle $\phi_m = \frac{\pi}{2m}$.
\end{definition}
We denote the vertex set of an $m$-arc $A$ as $\vertex(A)=\{\mathbf{v}_0,\dots,\mathbf{v}_m\}$, the set of its line segments as $\segment(A) = \{\ell_1,\dots,\ell_m\}$, where $\ell_i=\conv\kl{\skl{\mathbf{v}_{i-1}, \mathbf{v}_i}}$, and the set of its scaling factors as $\scale(A)=\{r_1,\dots,r_m\}$.

It will turn out to be useful to remove some ambiguity from the set of $m$-arcs by standardising their behaviour at the start and end vertices. Some examples can be seen in \cref{fig:arc_examples}.
\begin{definition}[Standard Arcs]
 An $m$-arc $A$ with vertices $\skl{\mathbf{v}_0,\dots,\mathbf{v}_m}$ and length scales $\skl{r_1,\dots,r_m}$ is called standardised (or a standard $m$-arc), if it starts at the origin and has a normalised last line segment, i.e.\ if $\mathbf{v}_0=\mathbf{0}$ and $r_m=1$. The set of all standard $m$-arcs is denoted $\CA_m$ and equipped with the metric
 \[
    d_{\CA}(A_1, A_2) = \max_{i=1,\dots,m} \|\vertex\kl{A_1}_i - \vertex\kl{A_2}_i\|_2.
 \]
\end{definition}

 The metric $d_\CA$ is induced by a mixed $(\ell_2,\ell_\infty)$-norm and thus indeed a proper metric. There is a one-to-one correspondence between the sets $\CA_m$ and $\mathbb{R}^{m-1}_+\times\skl{1}\cong\mathbb{R}^{m-1}_+$ via the scaling factors.

 \begin{figure}
     \centering
     \input{./img/arc_examples}
     \caption{Examples of several standard $m$-arcs for varying $m$ (left). Standardised arcs are always contained within the non-negative orthant and end on a horizontal line segment. The distance between two $m$-arcs is the maximal Euclidean distance of corresponding vertices (right).}
     \label{fig:arc_examples}
 \end{figure}

 \begin{definition}[Arc-Supported Measures]
 A probability measure $\mu\in\mathcal{D}(\mathbb{R}^2)$ is said to be $\delta$-distributed on a standard $m$-arc $A\in\CA_m$ if $\supp(\mu) \subseteq A$ and for each line segment $\ell\in\segment(A)$ we have $\mu(\ell) \geq \delta$.

 For $d>2$, a probability measure $\mu\in\mathcal{D}(\mathbb{R}^d)$ is said to be $\delta$-distributed on a standard $m$-arc $A\in\CA_m$ if $\supp(\mu) \subseteq \linspan\{\mathbf{e}_1,\mathbf{e}_2\}$ and by identifying $\linspan\{\mathbf{e}_1,\mathbf{e}_2\}$ with $\mathbb{R}^2$ it is $\delta$-distributed on a standard $m$-arc in the above sense.

 We denote the set of probability measures that are $\delta$-distributed on standard $m$-arcs by $\mathcal{D}_{\delta,m}(\mathbb{R}^d)\subseteq\mathcal{D}(\mathbb{R}^d)$.
 \end{definition}

\begin{lemma}\label{lem:arc_unique}
    Let $\mu\in\mathcal{D}_{\delta,m}(\mathbb{R}^d)$ be $\delta$-distributed on a standard $m$-arc, then this arc is unique.
\end{lemma}
\begin{proof}
 Towards a contradiction assume that $\mu$ is $\delta$-distributed on $A_1, A_2 \in \CA_m$ and $A_1 \neq A_2$. Let $\skl{\ell_1^j,\dots,\ell_m^j} = \segment(A_j)$ and $\skl{r_1^j,\dots,r_m^j} = \scale(A_j)$ denote the line segments and scaling factors of $A_1$ and $A_2$ respectively.
 Since $A_1\neq A_2$ and $r^1_m=r^2_m=1$ there is a smallest index $1\leq i\leq m-1$ such that $r_i^1\neq r^2_i$. Without loss of generality assume $r_i^1 < r^2_i$. Then $\ell^2_{i+1}$ lies outside of the convex hull of $A_1$ and in particular
 \[
 \ell^2_{i+1} \cap A_1 = \emptyset.
 \]
 But this contradicts the fact that both $\mu(\ell^2_{i+1})\geq \delta$ and $\supp\kl{\mu}\subseteq A_1$ must hold.
 \end{proof}

\Cref{lem:arc_unique} shows that arc-supported measures induce unique standard $m$-arcs. We denote the corresponding function mapping a measure $\mu\in\mathcal{D}_{\delta,m}(\mathbb{R}^d)$ to its induced arc $A_\mu\in\CA_m$ by $\arc\colon\mathcal{D}_{\delta,m}(\mathbb{R}^d)\to\CA_m\colon\mu\mapsto A_\mu$.
For an $n$-parameter family of distributions, $p\colon\Omega\to\mathcal{D}(\mathbb{R}^d)$, we denote by
\[
    \Omega_{\delta, m}=\skl{\,\theta\in\Omega\,:\,p(\theta)\in\mathcal{D}_{\delta,m}(\mathbb{R}^d)\,}\subseteq\Omega
\]
the set of parameters mapping to arc-supported measures. Without further assumptions on $p$ this set might well be empty. We will see however, that this can not happen for ReLU-invariant families as long as they contain at least one measure with a support of cardinality at least $m$.

We are now ready to start discussing the main steps of the proof of \cref{thm:main}. It relies on the following three results.

\begin{restatable}[Surjectivity]{lemma}{Surjectivity}\label{lem:xi_surjective}
Under \cref{asm:main} and for any $m\in\mathbb{N}$ there exists a $\delta>0$ such that the map $\Xi\colon\Omega_{\delta,m}\subseteq\mathbb{R}^n\to\mathbb{R}^{m-1}_+\times\{1\}$ given by $\Xi=\scale\circ\arc\circ\,p\vert_{\Omega_{\delta,m}}$ is surjective. In particular, in this case the domain $\Omega_{\delta,m}$ is non-empty.
\end{restatable}
The proof for this is given in \cref{sec:xi_surjective}. In short, we show that any measure with non finite support can be transformed by ReLU layers into a measure supported on an arbitrary standard arc. Then using the scaling factors to identify arcs with $\mathbb{R}^{m-1}_+\times\{1\}$ yields the claim.

\begin{restatable}[Local Lipschitz continuity]{lemma}{LocLip}\label{lem:xi_loclip}
Under \cref{asm:main} and for any $m\in\mathbb{N}$ let $\delta>0$ and $\Xi\colon\Omega_{\delta,m}\subseteq\mathbb{R}^n\to\mathbb{R}^{m-1}_+\times\{1\}$ be as in \cref{lem:xi_surjective}. Then $\Xi$ is locally Lipschitz continuous.
\end{restatable}
The proof for this is given in \cref{sec:xi_loclip}. In short, we show that all three partial functions $\scale$, $\arc$, and $p$ are locally Lipschitz continuous, hence also their composition $\Xi$ is.

As a final piece before the main theorem we need a rather general result on maps between Euclidean spaces.
\begin{restatable}{lemma}{hdimlemma}\label{lem:loclip_n_to_m}
Let $B\subseteq \mathbb{R}^n$ and $f\colon B\to \mathbb{R}^m$ be locally Lipschitz continuous. If the image $f(B)$ is a Borel set with non-empty interior in $\mathbb{R}^m$, then $m\leq n$. In particular, if $f$ maps surjectively onto $\mathbb{R}^m$ or $\mathbb{R}^m_+$, then $m\leq n$.
\end{restatable}
The full proof of this is given in \ref{apx:hdim}. It essentially uses the fact that locally Lipschitz continuous maps cannot increase the Hausdorff dimension. Since the Hausdorff dimension of $B$ is at most $n$ and the Hausdorff dimension of $f(B)$ is $m$, the claim follows.

We now have all the ingredients to prove the main theorem.
\begin{proof}[Proof of \cref{thm:main}]
Towards a contradiction let \cref{asm:main} hold. Let $m>n+1$ be arbitrary. By \cref{lem:xi_surjective} there exists a $\delta>0$ so that $\Omega_{\delta,m}$ is non-empty and  $\Xi\colon\Omega_{\delta,m}\to\mathbb{R}^{m-1}_+\times\{1\}$ is surjective. By \cref{lem:xi_loclip} it is also locally Lipschitz continuous. Since $\mathbb{R}^{m-1}_+\times\{1\}$ is isometric to $\mathbb{R}^{m-1}_+$ and $m-1>n$ this contradicts \cref{lem:loclip_n_to_m}. Hence, one of the assumptions in \cref{asm:main} cannot hold, which proves \cref{thm:main}.
\end{proof}

 \section{Proof of \cref{lem:xi_surjective}: Surjectivity of \texorpdfstring{$\Xi$}{Xi}}\label{sec:xi_surjective}

 For the remainder of the section let \cref{asm:main} hold. We will prove \cref{lem:xi_surjective} and show that for any $m\in\mathbb{N}$ there exists a $\delta>0$ so that the map $\Xi\colon\Omega_{\delta,m}\to\mathbb{R}^{m-1}_+\times\{1\}$ is surjective. We will do this by transforming (through a series of ReLU layers) any measure with infinite support into a measure that is supported on an arbitrary standard $m$-arc. Such an infinitely supported measure exists within the $n$-parameter family $p$ by assumption. Due to the ReLU-invariance the transformed measure is then also included in the family. As a consequence $\Omega_{\delta,m}$ must be non-empty and in fact we show that $\arc\circ\,p\vert_{\Omega_{\delta,m}}\colon\Omega_{\delta,m}\to\CA_m$ is surjective. The choice of $\delta$ will depend on how ``evenly spread'' the support of the original, untransformed, measure is. Finally, we have already discussed that $\CA_m$ can be identified with $\mathbb{R}^{m-1}_+\times\{1\}$ via the scaling factors.

 As a first step we observe how transforming measures by continuous functions affects their support. The support $\supp(\mu)$ of a Radon measure $\mu$ is defined as the complement of the largest open $\mu$-null set. A point $\mathbf{x}$ belongs to $\supp(\mu)$ if and only if every open neighbourhood of $\mathbf{x}$ has positive measure. A Borel set contained in the complement of $\supp(\mu)$ is a $\mu$-null set. The converse holds for open sets, i.e.\ an open set intersecting $\supp(\mu)$ has positive measure (in fact it suffices if the intersection is relatively open in $\supp(\mu)$).
 \begin{lemma}\label{lem:cont_supp}
 For any $q,r\in\mathbb{N}$ let $\mu\in\mathcal{D}(\mathbb{R}^q)$ and let $f\colon\mathbb{R}^q\to\mathbb{R}^r$ be continuous. Then
 \[
    \supp\kl{f_\ast(\mu)} = \overline{f\kl{\supp(\mu)}},
 \]
 where $f_\ast(\mu)\in\mathcal{D}(\mathbb{R}^r)$ is the pushforward of $\mu$ under $f$.
 \end{lemma}
 \begin{proof}
 We begin by showing $f\kl{\supp(\mu)}\subseteq\supp\kl{f_\ast(\mu)}$. Let $\mathbf{t}\in f\kl{\supp(\mu)}$ and $U\subseteq\mathbb{R}^r$ be any open neighbourhood of $\mathbf{t}$. There exists $\mathbf{x}\in\supp(\mu)$ so that $\mathbf{t}=f(\mathbf{x})$. Since $\mathbf{t}\in U$, we get $\mathbf{x}\in f^{-1}(U)$. Thus $f^{-1}(U)$ is an open neighbourhood of $\mathbf{x}$ and since $\mathbf{x}\in\supp(f)$ we get $f_\ast(\mu)(U) = \mu(f^{-1}(U))>0$. Since $U$ was an arbitrary neighbourhood of $\mathbf{t}$ this shows $\mathbf{t}\in\supp\kl{f_\ast(\mu)}$. This implies $f\kl{\supp(\mu)}\subseteq\supp\kl{f_\ast(\mu)}$. But $\supp\kl{f_\ast(\mu)}$ is closed, so we can even conclude $\overline{f\kl{\supp(\mu)}}\subseteq\supp\kl{f_\ast(\mu)}$. \par
 We will show the converse inclusion by showing $\overline{f\kl{\supp(\mu)}}^c\subseteq\supp\kl{f_\ast(\mu)}^c$. For this let $\mathbf{t}\in \overline{f\kl{\supp(\mu)}}^c$. Then there exists an open neighbourhood $U$ of $\mathbf{t}$ such that $U\cap \overline{f\kl{\supp(\mu)}} = \emptyset$, which implies $f^{-1}(U)\cap\supp(\mu)=\emptyset$. But this is only possible if $f_\ast(\mu)(U)=\mu\kl{f^{-1}(U)}=0$ and therefore $\mathbf{t}\in\supp\kl{f_\ast(\mu)}^c$.
 \end{proof}

 To simplify notations in the following, we begin the transformation by establishing a standardised situation that takes places exclusively in the two-dimensional subspace $\linspan\skl{\mathbf{e}_1,\mathbf{e}_2}$ of $\mathbb{R}^d$. In fact, all $m$-arcs are essentially one-dimensional objects embedded in a two-dimensional space. Therefore, we start by transforming the infinitely supported measure to the non-negative part of the one-dimensional space $\linspan\skl{\mathbf{e}_1}$, which will then subsequently be ``bent'' into the correct arc living in $\linspan\skl{\mathbf{e}_1,\mathbf{e}_2}$.

\subsection{Projection to a  One-Dimensional Subspace}\label{sec:onedim_twodim}

 \begin{lemma}\label{lem:proj_one_dim}
 Let $\mu\in\mathcal{D}(\mathbb{R}^d)$ such that $\supp(\mu)$ is not finite. Then there exists $\mathbf{W}\in\mathbb{R}^{d\times d}$ and $\mathbf{b}\in\mathbb{R}^d$ such that the image of $f\colon\mathbb{R}^d\to\mathbb{R}^d\colon\mathbf{x}\mapsto\varrho(\mathbf{W}\mathbf{x}+\mathbf{b})$ is contained in $\linspan\skl{\mathbf{e}_1}$ and $\supp\kl{f_\ast (\mu)}$ is infinite.
 \end{lemma}
 \begin{proof}
  Consider the coordinate projections $\mathrm{proj}_i\colon\mathbb{R}^d\to\mathbb{R}\colon\mathbf{x}\mapsto\mathbf{e}_i^\top\mathbf{x}$. Clearly we have
  \[
   \supp(\mu) \subseteq \mathrm{proj}_1\kl{\supp(\mu)}\times\dots\times\mathrm{proj}_d\kl{\supp(\mu)}.
  \]
  By assumption $\supp(\mu)$ is infinite and the product can only be infinite if at least of its factors is infinite, for example $\mathrm{proj}_j\kl{\supp(\mu)}$. If $\mathrm{proj}_j\kl{\supp(\mu)}\cap\mathbb{R}_+$ is infinite, we set $\sigma=+1$, otherwise $\mathrm{proj}_j\kl{\supp(\mu)}\cap\mathbb{R}_-$ must be infinite and we set $\sigma=-1$. Now choosing
  \[
    \mathbf{W} = \begin{bmatrix} \sigma\mathbf{e}_j & \mathbf{0} & \dots & \mathbf{0}\end{bmatrix}^\top\quad\text{and}\quad
    \mathbf{b} = \mathbf{0},
  \]
  clearly yields $f(\mathbb{R}^d)\subseteq\linspan\skl{\mathbf{e}_1}$. Further, by \cref{lem:cont_supp} we have \begin{align*}
    \supp\kl{f_\ast(\mu)} &= \overline{f\kl{\supp(\mu)}} \\
    &= \overline{\varrho\kl{\sigma \mathrm{proj}_j\kl{\supp(\mu)}}}\times\skl{0}\times\dots\times\skl{0} \\
    &= \begin{cases}
        \overline{\mathrm{proj}_j\kl{\supp(\mu)}\cap\mathbb{R}_+}\times\skl{0}\times\dots\times\skl{0}, &\quad\text{if }\sigma=+1\\
        \overline{-\mathrm{proj}_j\kl{\supp(\mu)}\cap\mathbb{R}_+}\times\skl{0}\times\dots\times\skl{0}, &\quad\text{if }\sigma=-1
       \end{cases},
  \end{align*}
  which by the choice of $\sigma$ is infinite.
 \end{proof}
 From now on we only need to consider measures with infinite support within the non-negative part of $\linspan\{\mathbf{e}_1\}$. Two dimensions will be enough for our construction to ``bend'' the measures into arc shape. Any two-dimensional ReLU layer $\mathbb{R}^2\to\mathbb{R}^2\colon\mathbf{x}\mapsto\varrho(\mathbf{W}\mathbf{x}+\mathbf{b})$ can easily be extended to a $d$-dimensional ReLU layer
 \[
    \mathbb{R}^d\to\mathbb{R}^d\colon\mathbf{x}\mapsto\varrho\kl{\begin{bmatrix} \mathbf{W} & \mathbf{0}_{2\times(d-2)} \\ \mathbf{0}_{(d-2)\times 2} & \mathbf{0}_{(d-2)\times(d-2)} \end{bmatrix}\mathbf{x} + \begin{bmatrix} \mathbf{b} \\ \mathbf{0}_{(d-2)\times 1}\end{bmatrix}},
 \]
 only operating on $\linspan\skl{\mathbf{e}_1,\mathbf{e}_2}\subseteq\mathbb{R}^d$ by appropriately padding with zeros. To simplify the notation, we leave out these zero paddings and just assume without loss of generality $d=2$ from now on.

 \subsection{Partitioning into Line Segments}

Now we show how the infinite support on the non-negative real line can be partitioned into an arbitrary number of intervals, which in the end will correspond to the arc line segments. We can write $\supp(\mu) = S\times\skl{0}\subseteq\mathbb{R}^2$ for some set $S\subseteq\mathbb{R}_{\geq 0}$.
 \begin{lemma}\label{lem:interval_partition}
 Let $S\subseteq\mathbb{R}_{\geq0}$ be infinite. Then for any $m\in\mathbb{N}$ there exist distinct points $0=b_0<\dots<b_m<\infty$ such that $S\cap \kl {b_{j-1},b_{j}}\neq\emptyset$ for all $j=1,\dots,m$.
 \end{lemma}
 \begin{proof}
 Since $S$ is infinite for any $m\in\mathbb{N}$ there exist distinct points $s_1,\dots,s_m\in S$ with $0<s_1<\dots<s_m$. Setting $b_0=0$,
 \[
    b_j = \frac{s_j+s_{j+1}}{2}\quad\text{for}\quad j=1,\dots,m-1,
 \]
 and $b_m = s_m + 1$ we get $s_j\in\kl{b_{j-1},b_j}$ for $j=1,\dots,m$.
 \end{proof}
 Further, we can restrict the support of the measure to a compact subset by clipping  the set $S$ to the range $\ekl{b_0,b_m}=\ekl{0,b_m}$. This can be achieved using the two layer ReLU network
 \[
 \mathbb{R}^2\rightarrow\mathbb{R}^2 \colon \mathbf{x} \mapsto \varrho\kl{-\varrho\kl{-\mathbf{x}+ \begin{bmatrix} b_m \\ 0 \end{bmatrix}} + \begin{bmatrix} b_m \\ 0\end{bmatrix}}
 \]
 with $b_m$ as in \cref{lem:interval_partition}. Altogether, without loss of generality we can from now on assume that $\mu$ is supported on $[b_0,b_m]\times\skl{0}\subseteq\mathbb{R}^2$ and the support intersects each subset $\kl{b_{j-1},b_j}\times\skl{0}$ for $j=1,\dots,m$.

 \subsection{Resizing the Line Segments}

 After partitioning the support of a distribution into intervals that each contain a certain amount of the probability mass we now want to bend it into arc shape. Each interval in the partition will correspond to a line segment. However, the bending will distort the segment lengths, so in order to obtain an arc with specified segment lengths we first have to resize the intervals. This is also possible with ReLU layers.

  \begin{lemma}\label{lem:nn_relu_b2a}
 For $m\geq 2$ let $0=a_0<\dots<a_m<\infty$ and $0=b_0<\dots<b_m<\infty$. Then there exists a collection of ReLU layer transformations $\skl{F_j\colon\mathbb{R}^2\to\mathbb{R}^2}_{j=0,...,m}$ such that their composition $F = F_{m}\circ F_{m-1}\circ \dots \circ F_{1}\circ F_{0}$ satisfies
 \[
    F\kl{\bmat{b_i \\ 0}} = \bmat{a_i \\ 0}\quad\text{for all}\quad  i=0,\dots,m,
\]
and restricted to $\mathbb{R}_{\geq 0}\times\{0\}$ it is a piecewise linear map with breakpoints (possibly) at $b_0,\dots,b_m$.
 \end{lemma}

\begin{proof}

 The proof consists of two parts. First, we iteratively construct a collection of piecewise linear and strictly increasing functions $\skl{f_j\colon\mathbb{R}\to\mathbb{R}}_{j=0,...,m}$ such that
 \begin{equation}\label{eq:resize_segments}
    \kl{f_j\circ\dots\circ f_0}\kl{b_i} = a_i \quad\text{for all}\quad i=0,\dots,j,
 \end{equation}
 for any $j=0,\dots,m$. In particular for $j=m$ and $f=f_m\circ f_{m-1}\circ\dots\circ f_0$ we obtain $f(b_i)=a_i$ for all $i=0,\dots,m$.
 Afterwards we show how to construct $\skl{F_j}_j$ from $\skl{f_j}_j$.

 Since $b_0=a_0=0$, we can start with
 \[
   f_0(x) = x
 \]
 and see that \eqref{eq:resize_segments} is satisfied for $j=0$. Now assuming we have already constructed strictly increasing and piecewise linear functions $f_0,\dots,f_j$ satisfying \eqref{eq:resize_segments} we now want to construct $f_{j+1}$. The idea is to choose it in such a way that it leaves the points already correctly mapped by $f_j\circ\dots\circ f_0$ unchanged but linearly transforms the remaining points so that $b_{j+1}$ is mapped to $a_{j+1}$. For brevity we denote $b_i^j = \kl{f_j\circ\dots\circ f_0}(b_i)$ for $i=0,\dots,m$. By assumption $b_i^j=a_i$ for $i=0,\dots,j$ and by strict monotonicity $b_{j+1}^j>b_j^j=a_j$. We set
 \[
   f_{j+1}(x)
   = \begin{cases}
   x, &\quad x \leq a_j, \\

   x + \frac{a_{j+1}-b^j_{j+1}}{b^j_{j+1}-a_j}(x-a_j), &\quad x > a_j.
   \end{cases}
 \]
 It is easy to see that this function leaves all points $a_1<a_2\dots<a_j$ unchanged, maps $b^{j}_{j+1}$ to $a_{j+1}$ and is continuous, monotonously increasing and piecewise linear.
 We can also express $f_{j+1}$ using $\varrho$, i.e.
 \begin{equation}\label{eq:fj_resize}
    f_{j+1}(x) = x + \frac{a_{j+1} - b^j_{j+1}}{b^j_{j+1} - a_j} \varrho\kl{x-a_j}.
 \end{equation}
 Then $f_{j+1}$ maps $b^j_{j+1}$ to $a_{j+1}$ and acts as the identity map on the region $x\leq a_j$ where points are already correctly mapped, as illustrated in \ref{fig:mapping}.  Hence clearly \eqref{eq:resize_segments} is also satisfied for $j+1$. To see that $f_{j+1}$ is strictly increasing we observe that its slope is $1$ on the region $x\leq a_j$ and
 \[
    1+\frac{a_{j+1} - b^j_{j+1}}{b^j_{j+1} - a_j} = \frac{a_{j+1} - a_j}{b^j_{j+1} - a_j} >0
 \]
 on the region $x\geq a_j$. This adds at most one new breakpoint at $b_j$ to the composition $f_{j+1}\circ\dots\circ f_0$. Continuing this process until $j=m$ finishes the first part.
 \begin{figure}
     \centering
     \resizebox{\textwidth}{!}{
     \input{./img/mapping}
     }
     \caption{The function $f_{j+1}(x)$ transforms the points $0=a_0=b^j_0, a_1=b^j_1,\dots,a_j=b^j_j$ and $b^j_{j+1},\dots,b^j_m$ (left). The region $x\leq a_j$ is left unchanged and the region $x\geq a_j$ is linearly rescaled so that $f_{j+1}(b^j_{j+1})=a_{j+1}$. The piecewise-linear function $f=f_m\circ\dots\circ f_0$ maps $b_j$ to $a_j$ (right).}
     \label{fig:mapping}
\end{figure}

Next we derive the ReLU layers $\skl{F_j}_j$ from the functions $\skl{f_j}_j$. For $j=0,\dots,m-1$ we denote
\[
 w_{j+1} = \frac{a_{j+1} - b^j_{j+1}}{b^j_{j+1} - a_j}
\]
and rewrite \eqref{eq:fj_resize} as
\[
 f_{j+1}(x) = \bmat{1 & w_{j+1}}\varrho\kl{ \bmat{1 \\ 1} x + \bmat{0 \\ -a_j}}\quad\text{for}\quad x\geq 0.
\]
Two consecutive maps can be combined as
\[
    f_{j+1}\kl{f_{j}(x)} = \bmat{1 & w_{j+1}} \varrho \kl{ \bmat{1 & w_{j} \\ 1 & w_{j}} \varrho \kl{\bmat{1 \\mathbbm{1}}x + \bmat{0 \\ -a_{j-1}}   } + \bmat{0 \\ -a_{j}}  }
\]
for $x\geq 0$, hence we set
\[
 F_{j}\colon\mathbb{R}^2 \to \mathbb{R}^2 \colon\mathbf{x}\mapsto
\varrho \kl{ \bmat{1 & w_{j} \\ 1 & w_{j}} \mathbf{x} + \bmat{0 \\ -a_{j}}  }
\]
for $j=1,...,m-1$. It is now easy to check that with
\[
 F_0(\mathbf{x}) = \varrho\kl{\bmat{ 1 & 0 \\ 1 & 0} \mathbf{x} }
\]
and
\[
  F_m(\mathbf{x}) \mapsto \varrho\kl{ \bmat{1 & w_m \\ 0 & 0} \mathbf{x} }
\]
we finally get
\[
 F\kl{\bmat{x \\ 0}} = \bmat{f\kl{x} \\ 0}\quad\text{for}\quad x\geq 0
\]
for $F=F_m\circ F_{m-1} \circ \dots \circ F_1 \circ F_0$.
\end{proof}

\subsection{Bending the Arc}

 Let us next consider how to ``bend'' a one dimensional subspace into an $m$-arc. This iterative procedure is illustrated in \cref{fig:bending}. It starts with all vertices in a straight line on the $x_2$-axis. First the polygonal chain is rotated by $-\phi_m$ (``to the right'') and translated in negative $x_1$ direction (``to the left'') such that the second last vertex lies on the $x_2$-axis. Then a ReLU is applied, which projects all vertices up until the second last onto the $x_2$-axis. The last last arc segment now has the correct angle and the process is repeated until all segments are correctly ``bent''. Note that each projection shrinks the unbent line segments by a factor of $\cos(\phi_m)$.
 \begin{figure}
     \centering
     \resizebox{\textwidth}{!}{
     \input{./img/bending}
     }
     \caption{Illustration of the 'bending' of a 3-arc in 4 steps. Step (i) shows $G_0$ ($90^{\circ}$ rotation) applied to $[a_m,a_0]\times\{0\}$.  Steps (ii) and (iii) show the image of $G_1\circ G_0$ and $G_2\circ G_1\circ G_0$ respectively (blue) as well as the corresponding transform just before the last ReLU-activation (gray). Step (iv) shows the final arc arising as the image of $G=G_3\circ G_2\circ G_1\circ G_0$.
     }
     \label{fig:bending}
\end{figure}
 To describe this transformation formally, we denote rotation matrices in $\mathbb{R}^2$ by $$R_\alpha=\begin{bmatrix}\cos(\alpha) & -\sin(\alpha) \\ \sin(\alpha) & \cos(\alpha) \end{bmatrix}.$$
 Furthermore, for $m\geq 2$ and $0=a_0<a_1<\dots<a_m<\infty$ we define
 \[
  \mathbf{d}_j = \bmat{\sin\kl{\phi_m}\cos\kl{\phi_m}^j a_{m-j} \\ 0},\quad\text{for}\quad j=1,\dots,m,
 \]
 as well as
 \begin{align*}
     &G_0\colon \mathbb{R}^2 \rightarrow \mathbb{R}^2\colon \mathbf{x} \mapsto \varrho\kl{R_{\frac{\pi}{2}}\mathbf{x}}
     \intertext{and}
     &G_j\colon \mathbb{R}^2 \rightarrow \mathbb{R}^2\colon \mathbf{x} \mapsto \varrho\kl{R_{-\phi_m}\mathbf{x} - \mathbf{d}_j} \quad\text{for}\quad j=1,\dots,m,
 \end{align*}
 which clearly are ReLU-transformations.

  \begin{lemma}\label{lem:nn_relu_a2arc}
 Let $m\geq 2$ and $0=a_0<a_1<\dots<a_m<\infty$ and $G_j$ as above. Then
 \[
    G\colon[a_0,a_m]\times\{0\}\subseteq\mathbb{R}^2\to\mathbb{R}^2\colon \mathbf{x}\mapsto (G_m \circ G _{m-1}\circ \dots\circ G_1\circ G_0)(\mathbf{x})
 \]
 is a parametrisation of a standard $m$-arc with vertices $\mathbf{v}_i=G\left(\begin{bmatrix}a_i \\ 0\end{bmatrix}\right)$, line segments $\ell_i=G([a_{i-1},a_i]\times\{0\})$, and scaling factors $r_i = \cos^{m-i}(\phi_m)(a_i-a_{i-1})$ for $i=1,\dots,m$.
 \end{lemma}

 \begin{proof}
  We need to show that $G$ is (in its first component) a piecewise-linear function with breakpoints $a_i$, and that $\mathbf{v}_0 = \mathbf{0}$ as well as
  \begin{equation}\label{eq:bend_claim}
      \mathbf{v}_i = \mathbf{v}_{i-1} + \cos^{m-i}(\phi_m)(a_i-a_{i-1}) \bmat{\sin(i\phi_m) \\ \cos(i\phi_m)},\quad\text{for}\quad i=1,\dots,m.
  \end{equation}
  For this we denote \[\mathbf{v}_i^j = (G_j\circ\dots\circ G_0)\kl{\bmat{a_i \\ 0}}\quad\text{for}\quad i,j=0,...,m.\] By the choice of $\phi_m$ we know $\sin(\phi_m)>0$ and $\cos(\phi_m)>0$. Hence, the first component of the bias vectors $\mathbf{d}_j$ satisfy $\sin(\phi_m)\cos(\phi_m)^{j}a_{m-j}\geq 0$ for $j=1,\dots,m$.
  Therefore, $\mathbf{v}^j_0=\mathbf{0}$ for all $j$ and in particular $\mathbf{v}_0=\mathbf{v}^m_0=\mathbf{0}$. Further, we will show inductively over $j$ that for any $i,j$ we have
  \begin{equation}\label{eq:bend_induction}
   \mathbf{v}^j_i-\mathbf{v}^j_{i-1} = \begin{cases}
   \cos\kl{\phi_m}^{j} (a_i-a_{i-1}) \bmat{0 \\ 1}, &\quad i\leq m - j, \\
   \cos\kl{\phi_m}^{m-i} (a_i-a_{i-1}) \bmat{\sin\kl{(i+j-m)\phi_m} \\ \cos\kl{(i+j-m)\phi_m}}, & \quad i> m-j.
   \end{cases}
\end{equation}
  Then \eqref{eq:bend_claim} follows from \eqref{eq:bend_induction} with $j=m$.

  To start the inductive proof we observe that for $j=0$ and $i=1,\dots,m$ we have
  \begin{align*}
     \mathbf{v}_i^0 - \mathbf{v}_{i-1}^0 &= G_0\kl{\bmat{a_i \\ 0}} - G_0\kl{\bmat{a_{i-1} \\ 0}} \\
     &= R_{\frac{\pi}{2}}\bmat{a_i \\ 0} - R_{\frac{\pi}{2}}\bmat{a_{i-1} \\ 0} \\
     &= (a_i-a_{i-1})\bmat{0\\mathbbm{1}},
  \end{align*}
  satisfying \eqref{eq:bend_induction}.

  Next, let us assume that \eqref{eq:bend_induction} holds for some $j<m$. For any $i$ we denote the vectors
  $\mathbf{w}^j_i = R_{-\phi_m}\mathbf{v}^j_i - \mathbf{d}_j$ and therefore $\mathbf{v}^{j+1}_i = G_{j+1}(\mathbf{v}^j_i) = \varrho(\mathbf{w}^j_i)$ and
  \begin{equation*}\label{eq:bend_pre_relu}
    \mathbf{w}^j_i - \mathbf{w}^j_{i-1} =
    \begin{cases}
        \cos\kl{\phi_m}^{j} (a_i-a_{i-1}) \bmat{\sin\kl{\phi_m} \\ \cos\kl{\phi_m}}, &\enskip i\leq m - j, \\
        \cos\kl{\phi_m}^{m-i} (a_i-a_{i-1}) \bmat{\sin\kl{(i+j+1-m)\phi_m} \\ \cos\kl{(i+j+1-m)\phi_m}}, &\enskip i> m-j.
   \end{cases}
 \end{equation*}
 We also have $\sin((i+j+1-m)\phi_m)> 0$ and $\cos((i+j+1-m)\phi_m)\geq 0$ for $i>m-j$ by the choice of $\phi_m$ and therefore $\mathbf{w}^j_m\geq \mathbf{w}^j_{m-1}\geq\dots\geq\mathbf{w}^j_1\geq\mathbf{w}^j_0=-\mathbf{d}_j$ component-wise, where even all inequalities are strict in the first component. We immediatley see that all $\mathbf{w}^j_i$ lie in the upper half-plane. Further,
 \begin{align*}
  \mathbf{w}^j_{m-j} &= \mathbf{w}_0^j + \sum_{i=1}^{m-j} (\mathbf{w}^j_i - \mathbf{w}^j_{i-1})\\
   &= -\mathbf{d}_j + \sum_{i=1}^{m-j}\cos\kl{\phi_m}^{j}(a_i-a_{i-1}) \bmat{\sin\kl{\phi_m} \\ \cos\kl{\phi_m}} \\
   &= -\bmat{\sin\kl{\phi_m}\cos\kl{\phi_m}^j a_{m-j} \\ 0}  +
   \cos\kl{\phi_m}^{j} a_{m-j} \bmat{\sin\kl{\phi_m} \\ \cos\kl{\phi_m}} \\
   &= \cos\kl{\phi_m}^{j+1} a_{m-j} \bmat{0 \\ 1},
 \end{align*}
 from which we can conclude that
 \begin{align}
     (\mathbf{w}_i^j)_1 &\geq 0,\quad\text{for}\quad i=1,\dots,m \label{eq:bending1}\\
     (\mathbf{w}_i^j)_2 &< 0,\quad\text{for}\quad i=1,...,m-j-1, \label{eq:bending2}\\
     (\mathbf{w}_i^j)_2 &\geq 0,\quad\text{for}\quad i=m-j,...,m. \label{eq:bending3}
 \end{align}
 Now $\mathbf{v}^{j+1}_i = \varrho(\mathbf{w}^j_i)$ together with \eqref{eq:bending1} to \eqref{eq:bending3} implies
  \begin{alignat*}{2}
     (\mathbf{v}_i^{j+1})_1 &= (\mathbf{w}_i^j)_1,&&\quad\text{for}\quad i=1,...,m, \\
     (\mathbf{v}_i^{j+1})_2 &= 0,&&\quad\text{for}\quad i=1,...,m-j-1, \\
     (\mathbf{v}_i^{j+1})_2 &= (\mathbf{w}_i^j)_2,&&\quad\text{for}\quad i=m-j,...,m,
 \end{alignat*}
 and therefore $\mathbf{v}^{j+1}_i-\mathbf{v}^{j+1}_{i-1}$ is equivalent to
 \[
   \begin{cases}
   \cos\kl{\phi_m}^{j+1} (a_i-a_{i-1}) \bmat{0 \\ 1}, &\enskip i\leq m - (j+1) \\
   \cos\kl{\phi_m}^{m-i} (a_i-a_{i-1}) \bmat{\sin\kl{(i+j+1-m)\phi_m} \\ \cos\kl{(i+j+1-m)\phi_m}}, &\enskip i> m-(j+1),
   \end{cases}
  \]
  which shows that \eqref{eq:bend_induction} holds for $j+1$ and concludes the inductive step.

  Each $G_{j+1}$ adds only one new breakpoint at to the overall function, corresponding to point $\mathbf{w}^j_{m-j}$ where the second component of $(G_j\circ\dots\circ G_0)(\mathbf{x})$ switches sign. This corresponds to the breakpoint $a_{m-j}$ of $G$.
 \end{proof}

 We are now ready to given the main proof of this section.
 \begin{proof}[Proof of \cref{lem:xi_surjective}]
 Let $m\in\mathbb{N}$ be arbitrary. By assumption there exists a measure $\mu\in\{p(\theta)\}_{\theta\in\Omega}$ with infinite support. As discussed in \cref{sec:onedim_twodim} and in particular by using \cref{lem:proj_one_dim} we can without loss of generality assume $d=2$ and $\supp(\mu)$ is compact and contained in the non-negative part of $\linspan\{\mathbf{e}_1\}$. In other words we can rewrite $\supp(\mu)=S\times\{0\}\subseteq[0,b]\times\{0\}$ for an infinite set $S\subseteq\mathbb{R}_{\geq 0}$ and some $b>0$. By \cref{lem:interval_partition} we can choose points $0=b_0<b_1<\dots<b_m=b$ so that $S$ intersects each $\kl{b_{j-1}, b_j}$. By the choice of $S$ and since $\kl{b_{j-1},b_j}\times\{0\}$ is relatively open in $S\times\{0\}$ we have $\mu((b_{j-1},b_j)\times\{0\})>0$ for all $j$. Set $0<\delta\leq\min_j \mu((b_{j-1},b_j)\times\{0\})$.

 Now let $A\in\CA_m$ be an arbitrary standard $m$-arc and denote its set of vertices by $(\mathbf{v}_0,\dots,\mathbf{v}_m)=\vertex(A)$, its line segments by $(\ell_1,\dots,\ell_m)=\segment(A)$, and its scaling factors by $(r_1,\dots,r_m)=\scale(A)$. We set $a_0=0$ and iteratively $a_j=a_{j-1} + r_j \cos^{j-m}(\phi_m)$ for $j=1,\dots,m$. By \cref{lem:nn_relu_b2a} there is a ReLU neural network $F$ satisfying
 \[
 F\kl{\begin{bmatrix}b_j \\ 0\end{bmatrix}} = \begin{bmatrix}a_j \\ 0\end{bmatrix}\quad\text{for}\quad j=1,\dots,m,
 \]
 and then by \cref{lem:nn_relu_a2arc} another ReLU neural network $G$ satisfying
 \[
 G\kl{\begin{bmatrix}a_i \\ 0\end{bmatrix}} = \mathbf{v}_j\quad\text{for}\quad j=1,\dots,m.
 \]
 Altogether $G\circ F$ is a ReLU neural network transforming $\mathbf{b}_j=\begin{bmatrix}b_i & 0\end{bmatrix}^\top$ to the vertex $\mathbf{v}_j$ and $\ekl{b_{j-1},b_j}\times\{0\}$ to the line segment $\ell_j$.
 By \cref{lem:cont_supp} we have
 \[
    \supp\kl{(G\circ F)_\ast\mu} = \overline{(G\circ F)(\supp(\mu))} \subseteq \overline{(G\circ F)(\ekl{b_0,b_m}\times\{0\})} \subseteq A.
 \]
 Also, for each $j$ we get
 \[
 ((G\circ F)_\ast\mu)(\ell_j)=\mu\kl{(G\circ F)^{-1}(\ell_j)} \geq \mu(\kl{b_{j-1},b_j}\times\{0\})\geq\delta
 \]
 and therefore $(G\circ F)_\ast \mu \in\mathcal{D}_{\delta,m}(\mathbb{R}^2)$ and $\arc\kl{(G\circ F)_\ast \mu}=A$.
 Using the ReLU-invariance, we know that $(G\circ F)_\ast\mu\in \{p(\theta)\}_{\theta\in\Omega_{\delta,m}}\subseteq \{p(\theta)\}_{\theta\in\Omega}$. Altogether, since $A\in\CA_m$ was arbitrary we conclude that $\arc\circ\,p\vert_{\Omega_{\delta,m}}$ is surjective onto $\CA_m$. Clearly, the map $\scale\colon\CA_m\to\mathbb{R}^{m-1}_+\times\{1\}$ is a bijection, hence also $\Xi\colon\Omega_{\delta,m}\to\mathbb{R}^{m-1}_+\times\{1\}$ is surjective.
 \end{proof}

 \section{Proof of \cref{lem:xi_loclip}: Local Lipschitz Continuity of \texorpdfstring{$\Xi$}{Xi}}\label{sec:xi_loclip}
 We will now show that the map $\Xi\colon\Omega_{\delta,m}\to\mathbb{R}^m_+$ is locally Lipschitz continuous by showing this for each of its three composite parts $p\vert_{\Omega_{\delta,m}}\colon\Omega_{\delta,m}\to\mathcal{D}_{\delta,m}(\mathbb{R}^d)$, $\arc\colon\mathcal{D}_{\delta,m}(\mathbb{R}^d)\to\CA_m$, and $\scale\colon\CA_m\to\mathbb{R}^m_+$.

 Firstly, $p\vert_{\Omega_{\delta,m}}\colon\Omega_{\delta,m}\to\mathcal{D}_{\delta,m}(\mathbb{R}^d)$ is locally Lipschitz continuous by \cref{asm:main}, so there is nothing to show.

 Secondly, the local Lipschitz continuity of $\arc\colon\mathcal{D}_{\delta,m}(\mathbb{R}^d)\to\CA_m$ can be derived from geometric observations.
 \begin{lemma}\label{lem:loclip1}
    The map $\arc\colon(\mathcal{D}_{\delta,m}(\mathbb{R}^d),d_P)\to(\CA_m,d_\CA)\colon \mu\mapsto A_\mu$, is locally Lipschitz continuous with Lipschitz constant $\frac{2\sqrt{2}}{\sin\kl{\phi_m}}$.
 \end{lemma}
 \begin{proof}
 Let $\mu_0\in\mathcal{D}_{\delta, m}(\mathbb{R}^d)$ be arbitrary. We will show Lipschitz continuity on the open ball of radius $\frac{\delta}{2}$ around $\mu_0$. For this let $\mu_1,\mu_2\in B_{\frac{\delta}{2}}(\mu_0) \cap \mathcal{D}_{\delta, m}(\mathbb{R}^d)$ and therefore $d_P(\mu_1,\mu_2)\leq \delta$. We know that $\mu_1$ and $\mu_2$ are $\delta$-supported on two unique standard $m$-arcs $A_1=\arc\kl{\mu_1}=A_{\mu_1}\in\CA_m$ and $A_2=\arc\kl{\mu_2}=A_{\mu_2}\in\CA_m$ and by definition these lie in the same two-dimensional subspace $\linspan\skl{\mathbf{e}_1,\mathbf{e}_2}$ of $\mathbb{R}^d$. So without loss of generality we carry out the remaining proof in $\mathbb{R}^2$. Let $(\ell_1^1,\dots,\ell_m^1) = \segment(A_1)$ and $(\ell_1^2,\dots,\ell_m^2) = \segment(A_2)$ denote the line segments of $A_1$ and $A_2$ respectively. We denote the affine hulls of individual line segments as $h_i^j = \affspan\kl{\ell_i^j}$ and further denote the Euclidean distance of two corresponding affine hulls by $d_i=\dist(h_i^1, h_i^2)$. An example for this can be seen in \cref{fig:outer-arc}.

 We want to upper bound the Euclidean distances $\|\mathbf{v}_i^1-\mathbf{v}_i^2\|_2$ of all corresponding vertices of the two arcs. For any $i$ one of four cases can occur, as visualised in \cref{fig:vertex-arrangements}. Firstly, if $\mathbf{v}_i^1=\mathbf{v}_i^2$ (cf.\ \Cref{fig:vertex-arrangements}, top left) we trivially get $\|\mathbf{v}_i^1-\mathbf{v}_i^2\|_2=0$.
 Secondly, if $\mathbf{v}_i^1\neq\mathbf{v}_i^2$ but $h_i^1=h_i^2$ (cf.\ \Cref{fig:vertex-arrangements}, top right) it is not hard to see that
 \[
    \|\mathbf{v}_i^1-\mathbf{v}_i^2\|_2 \leq \frac{d_{i+1}}{\sin\kl{\phi_m}}.
 \]
 Thirdly, if $\mathbf{v}_i^1\neq\mathbf{v}_i^2$ but $h_{i+1}^1=h_{i+1}^2$ (cf.\ \Cref{fig:vertex-arrangements}, bottom left) we similarly get
 \[
    \|\mathbf{v}_i^1-\mathbf{v}_i^2\|_2 \leq \frac{d_{i}}{\sin\kl{\phi_m}}.
 \]
 In the fourth case, where $\mathbf{v}_i^1\neq\mathbf{v}_i^2$, $h_i^1\neq h_i^2$, and $h_{i+1}^1\neq h_{i+1}^2$ the four intersections of the affine spaces form a parallelogram (cf.\ \Cref{fig:vertex-arrangements}, bottom right). We obtain
 \[
    \|\mathbf{v}_i^1-\mathbf{v}_i^2\|_2 \leq \frac{\sqrt{2}\kl{d_{i}+d_{i+1}}}{\sin\kl{\phi_m}},
 \]
 using the parallelogram identity. Thus, it remains to bound all the distances $d_i$.

 For this, assume that for some $i$ we have $d_i>0$. We say that $\ell_i^1$ is the \emph{outer} line segment and $\ell_i^2$ is the \emph{inner} line segment if $\dist\kl{\mathbf{0}, h_i^1} > \dist\kl{\mathbf{0}, h_i^2}$ and vice versa if $\dist\kl{\mathbf{0}, h_i^1} < \dist\kl{\mathbf{0}, h_i^2}$. Without loss of generality assume that $\ell_i^1$ is the outer line segment. Then for any point $\mathbf{p}\in A_2$ we have $\dist(\mathbf{p}, \ell_i^1) \geq d_i$. Thus $\mu_2((\ell_i^1)^\epsilon)=0$ for any $\epsilon<d_i$. But $\mu_1(\ell_i^1)\geq \delta$ and therefore by definition of the Prokhorov metric
 \[
  d_i \leq d_P(\mu_1,\mu_2).
 \]
 Altogether, we obtain
 \[
    d_\CA\kl{A_1, A_2} = \max_{i} \|\mathbf{v}_i^1-\mathbf{v}_i^2\|_2 \leq \frac{2 \sqrt{2}}{\sin\kl{\phi_m}} d_P\kl{\mu_1,\mu_2}. \qedhere
 \]
 \end{proof}

\begin{figure}
    \centering
    \input{./img/outer_arc}
    \caption{Example of two standard $m$-arcs $A_1$ (blue) and $A_2$ (orange) for $m=5$. Supporting affine subspaces $h_4^1$ and $h_4^2$ and their distance $d_4$ are shown exemplarily for the fourth line segments. At this segment $A_2$ is the outer arc and $A_1$ is the inner arc.}
    \label{fig:outer-arc}
 \end{figure}

\begin{figure}
     \centering
     \resizebox{0.85\textwidth}{!}{
     \input{./img/vertex_arrangements}
     }
     \caption{Four possible arrangements of corresponding vertices $\mathbf{v}_i^1$ and $\mathbf{v}_i^2$ of two arcs. In case one both vertices are equal (top left). In cases two and three the vertices are different but either the affine subspaces $h_i^1$ and $h_i^2$ or $h_{i+1}^1$ and $h_{i+1}^2$ coincide (top right and bottom left). In the fourth case all affine subspaces differ and their intersections form a parallelogram with $\mathbf{v}_i^1$ and $\mathbf{v}_i^2$ at opposite corners (bottom right).}
     \label{fig:vertex-arrangements}
\end{figure}

 Thirdly, the local Lipschitz continuity of $\scale\colon\CA_m\to\mathbb{R}^m_+$ is a straight-forward calculation if $\mathbb{R}^m_+$ is equipped with the metric induced by the $\ell_\infty$-norm. By norm equivalence the same also holds for $\mathbb{R}^m_+$ viewed as a subspace of Euclidean $\mathbb{R}^m$.
 \begin{lemma}\label{lem:loclip2}
    The map $\scale\colon(\CA_m, d_\CA)\to(\mathbb{R}_+^m,\|\cdot\|_\infty)\colon A\mapsto (r_1,\dots,r_m)$, is Lipschitz continuous with Lipschitz constant 2.
 \end{lemma}

 \begin{proof}
 Let $A_1, A_2\in\CA_m$ be standard arcs with vertices $(\mathbf{v}_1^1,\dots,\mathbf{v}_m^1)=\vertex(A_1)$ and $(\mathbf{v}_1^2,\dots,\mathbf{v}_m^2)=\vertex(A_2)$ and with scaling factors $(r_1^1,\dots,r_m^1)=\scale(A_1)$ as well as $(r_1^2,\dots,r_m^2)=\scale(A_2)$ respectively. Since for all $i=1,\dots,m$ we have
 \[r_i^1 \leq r_i^2 + \dist\kl{\mathbf{v}^1_{i-1}, \mathbf{v}^2_{i-1}} + \dist\kl{\mathbf{v}^1_{i}, \mathbf{v}^2_{i}}
 \]
 and vice versa
 \[r_i^2 \leq r_i^1 + \dist\kl{\mathbf{v}^1_{i-1}, \mathbf{v}^2_{i-1}} + \dist\kl{\mathbf{v}^1_{i}, \mathbf{v}^2_{i}},
\]
 we obtain
 \[
    \nkl{\scale\kl{A_1} - \scale\kl{A_2}}_{\infty} = \max_{{i=1,\dots,m}} \bkl{r_i^1 - r_i^2} \leq 2d_\CA\kl{A_1, A_2}. \qedhere
 \]
 \end{proof}

 Altogether we can conclude this section by combining all pieces and proving the local Lipschitz continuity of $\Xi\colon\Omega_{\delta,m}\to\mathbb{R}^m_+$.
 \begin{proof}[Proof of \cref{lem:xi_loclip}]
  The claim follows directly from \cref{asm:main,lem:loclip1,lem:loclip2} since the composition of locally Lipschitz continuous functions is again locally Lipschitz continuous.
 \end{proof}

\section{Proof of the Reverse Result (\cref{thm:counter})}\label{sec:counterproof}
 We will now come to the proof of \cref{thm:counter} and show that exploiting any of the three restrictions in \cref{thm:main} indeed allows us to find ReLU-invariant families of distributions. We split the proof in three parts and construct examples for each of the three cases.

 \subsection{Families of Distributions in One Dimension}\label{sec:counterproof_i}
 In this section we will prove part \ref{itm:counter_i} of \cref{thm:counter}. We explicitly construct a family of ReLU-invariant probability distributions on $\mathbb{R}$ with a locally Lipschitz continuous parametrisation map.

 The key observation is that ReLU neural networks in one dimension are a rather restricted class of functions unlike networks in higher dimensions. In fact, we will show that increasing the depth of neural networks beyond three layers does not lead to the ability to express more complicated functions. Every one-dimensional ReLU neural network can be rewritten as a three layer network. Thus we only require a constant number of parameters to describe the set of all one-dimensional ReLU networks. It is then straightforward to obtain a ReLU-invariant family of probability distributions by explicitly making these parameters part of the family's parametrisation.

 We will follow a similar idea to obtain ReLU-invariant families in higher dimension when proving part \ref{itm:counter_iii} of \cref{thm:counter} in \cref{sec:counterproof_iii}. However in this case the number of parameters to describe the neural networks grows with their depth leading to the need for an arbitrary large number of parameters to describe the distributions in the family. Although it is possible to generalise our construction to higher dimensions this comes at the cost of losing the local Lipschitz continuity of the parametrisation map.

 We proceed in three steps. We first analyse the set of one-dimensional ReLU networks and show that it can be completely described by only six parameters. We then use these to parametrise all ReLU neural networks as functions in $C(K,\mathbb{R})$ on some compact domain $K\subseteq\mathbb{R}$ via a so called \emph{realisation} map $R\colon\mathbb{R}^6\to C(K,\mathbb{R})$. Finally, we use a pushforward map $Q\colon C(K,\mathbb{R})\to\mathcal{D}(\mathbb{R})$ mapping neural network functions $f\in C(K,\mathbb{R})$ to the pushforward of a fixed \emph{prototype} measure $\mu_0\in\mathcal{D}(\mathbb{R})$ under $f$ to generate the family of probability distribution. The prototype measure is assumed to be supported within $K$. The one-parameter family $p\colon\mathbb{R}^6\to\mathcal{D}(\mathbb{R})$ is then simply given as $p=Q\,\circ R$. We will now discuss each of the steps in more detail.

 \paragraph{One-dimensional ReLU networks:} Let us first show that any one dimensional ReLU network can be rewritten as a three layer network. For this let $L\in\mathbb{N}$ and $f\colon\mathbb{R}\to\mathbb{R}$ be an $L$-layer ReLU network given as
 \[
    f=f_L\circ f_{L-1}\circ\dots\circ f_2\circ f_1
  \]
  with layers
  \[
    f_i\colon\mathbb{R}\to\mathbb{R}\colon x\mapsto\varrho(w_i x + b_i)\quad\text{for}\quad i=1,\dots,L.
  \]
  We observe that $f$ is constant as soon as any weight $w_i$ is zero and a constant function can easily be rewritten as a three layer network. So from now on we assume $w_i\neq 0$ for all layers. We denote by
  \[
    A_i = \skl{\,x\in\mathbb{R}\,:\,w_i x + b_i \geq 0\,}\quad\text{for}\quad i=1,\dots,L
  \]
  the domain on which the $i$-th layer is affine linear with slope $w_i$. It is constant on the complement $A_i^c$. Further we denote
  \[
    B_i = \kl{f_{i-1}\circ\dots\circ f_1}^{-1}(A_i)\quad\text{for}\quad i=2,\dots,L
  \]
  and $B_1=A_1$. Each $A_i$ is a convex and closed subset of $\mathbb{R}$ and $f_{i-1}\circ\dots\circ f_1$ is continuous and monotone as a composition of continuous and monotone functions. Hence, also all $B_i$ are closed and convex. In $\mathbb{R}$ these are exactly the (possibly unbounded) closed intervals.
  \begin{lemma}\label{lem:onedim_relu_const}
  Let $U\subseteq \mathbb{R}\setminus\bigcap_{i=1}^L B_i$ be connected. Then $f\vert_U$ is constant.
  \end{lemma}
  \begin{proof}
   We have $\mathbb{R}\setminus \bigcap_{i=1}^L B_i = \bigcup_{i=1}^L B_i^c$ and therefore $U = \bigcup_{i=1}^L U\cap B_i^c$. Since $B_i$ is closed $B_i^c$ is open and thus $U\cap B_i^c$ is relatively open in $U$. But by the choice of $B_i$ we know $\kl{f_i\circ\dots\circ f_1}\kl{U\cap B_i^c} \equiv 0$, hence $f\vert_{U\cap B_i^c}$ is constant. But since $f$ is continuous and $U$ is connected this already implies that $f\vert_U$ is constant.
  \end{proof}
  The one-dimensional ReLU network $f$ is piecewise constant except for the (possibly empty) region $\bigcap_{i=1}^L B_i$ where all layers are non constant and act affine linearly. On this region $f$ is affine linear with slope $\prod_{i=1}^L w_i$. As discussed above each $B_i$ is convex and closed, hence the same holds for their intersection. Two different cases can occur, which determine the overall shape of the function $f$. This is also visualized in \cref{fig:onedim_relu}.

  \begin{figure}
     \centering
     \resizebox{\textwidth}{!}{
     \input{./img/onedim_relu_compressed}
     }
     \caption{ReLU networks in dimension $d=1$ can only take one of five distinct forms. We show the four non-constant cases. This includes two unbounded cases (top) and two bounded cases (bottom).}
     \label{fig:onedim_relu}
\end{figure}

  Firstly, if $\bigcap_{i=1}^L B_i$ is bounded it must be empty, a single point, or a bounded and closed interval. Then $f$ is bounded. In fact it is constant or a piecewise linear step function with one step and two constant regions. These functions can be represented by three layer networks.

  Secondly, if $\bigcap_{i=1}^L B_i$ is an unbounded closed interval, then $f$ looks like a shifted and scaled ReLU function and can again be represented by a three layer ReLU network.

  \clearpage

  \begin{corollary}\label{cor:onedim_relu}
  The ReLU network $f$ has one of the following forms.
  \begin{itemize}
      \item $f$ is bounded and has two constant and one affine linear region:
      \[
        f(x) = \begin{cases}
        c_1, & \quad x \leq a_1 \\
        c_2, & \quad x \geq a_2 \\
        c_1 + \frac{c_2-c_1}{a_2-a_1}(x-a_1), & \quad a_1<x<a_2
        \end{cases}
      \]
      for some $c_1, c_2\geq 0$ and $a_1<a_2$. This includes the degenerate case where $f$ is constant by choosing $c_1=c_2$.
      \item $f$ is unbounded and has one constant and one affine linear region:
      \[
        f(x) = \begin{cases}
        c, & \quad w(x-a) \leq 0 \\
        c + w(x-a), & \quad w(x-a)>0
        \end{cases}
      \]
      for some $c\geq 0$, $a\in\mathbb{R}$, and $w\neq 0$.
  \end{itemize}
  \end{corollary}
  \begin{proof}
   This directly follows from \cref{lem:onedim_relu_const}, the continuity and monotonicity of the ReLU network $f$ and the fact that $f$ is lower bounded by zero.
  \end{proof}
  \begin{lemma}
  The two types of piecewise linear functions from \cref{cor:onedim_relu} can both be written as one dimensional ReLU networks with three layers.
  \end{lemma}
  \begin{proof}
   We give explicit formulas for choosing the weights $w_1, w_2, w_3$ and biases $b_1, b_2, b_3$ of the three layer network.

   In the first case of a bounded function we choose
   \begin{alignat*}{2}
   w_1 &= 1, \qquad&\qquad b_1 &= -a_1, \\
   w_2 &= -\frac{|c_2-c_1|}{a_2-a_1}, \qquad&\qquad b_2 &= |c_2-c_1|,  \\
   w_3 &= -\sign\kl{c_2-c_1}, \qquad&\qquad b_3 &= c_2,
   \end{alignat*}
   and in the second case of an unbounded function we choose
   \begin{alignat*}{2}
   w_1 &= \sign\kl{w}, \qquad&\qquad b_1 &= -\sign\kl{w}a, \\
   w_2 &= |w|, \qquad&\qquad b_2 &= c, \\
   w_3 &= 1, \qquad&\qquad b_3 &= 0.
   \end{alignat*}
   The claim now follows by straightforward calculations.
  \end{proof}

 \paragraph{Neural Network Parametrisation:}
 Neural networks can be viewed algebraically (as collections of weights and biases) or analytically (as functions in some function space). This distinction was made in \cite{pv2018network-calculus} and used in \cite{prv2018nn-topology-arxiv} to study the dependence of neural network functions on their weights and biases. The mapping from weights and biases to the corresponding neural network function is referred to as the \emph{realisation} map. In \cite{prv2018nn-topology-arxiv} the continuity of the realisation map was shown for neural networks of fixed size with continuous activation function (in particular this includes the ReLU activation), when they are viewed as functions in $C(K,\mathbb{R}^d)$ on a compact domain $K$. If the activation function is Lipschitz continuous (which again is the case for the ReLU activation) also the realisation map is Lipschitz continuous.

 As we have seen all ReLU neural networks in one dimension can be rewritten into a three layer network and can thus be parametrised by three weight and three bias parameters. Let $K\subseteq\mathbb{R}$ be any compact non-empty domain. For any collection of weights and biases $\theta=(w_1,b_1,w_2,b_2,w_3,b_3)\in\mathbb{R}^6$ we denote the function that is the ReLU neural network realisation of these weights and biases as $f_\theta\colon K\to\mathbb{R}$. Then the realisation map is
 \[
    R\colon\mathbb{R}^6\to C(K,\mathbb{R}^d)\colon \theta\mapsto f_\theta.
 \]
 From \cite[Proposition 5.1]{prv2018nn-topology-arxiv} we obtain the following result.
 \begin{lemma}
    The realisation map $R\colon\mathbb{R}^6\to C(K,\mathbb{R}^d)$ is Lipschitz continuous.
 \end{lemma}

 \paragraph{Pushforward Mapping:} We want to connect neural network functions with the effect they have on probability measures and obtain an invariant family of distributions. For this let $\mu_0\in\mathcal{D}(\mathbb{R})$ be any probability measure on $\mathbb{R}$ with $\supp(\mu_0)\subseteq K$. We call this the \emph{prototype} measure and will derive all other measures in the family as pushforwards of $\mu_0$ under ReLU neural networks. We define the pushforward function
 \[
    Q\colon C(K,\mathbb{R})\to\mathcal{D}(\mathbb{R})\colon f\mapsto f_\ast\mu_0.
 \]
 The restriction of $Q$ to the subset of ReLU neural network functions will be used to generate the invariant family of distributions.
 \begin{lemma}
    The pushforward function $Q\colon C(K,\mathbb{R})\to\mathcal{D}(\mathbb{R})$ is Lipschitz continuous.
 \end{lemma}
 \begin{proof}
   Let $f,g\in C(K,\mathbb{R})$ and $A\in\mathcal{B}(\mathbb{R})$. Assuming $f^{-1}(A)\neq\emptyset$, for any $x\in f^{-1}(A)$ we have $f(x)\in A$ and
   $|f(x)-g(x)|\leq\|f-g\|_\infty$. Consequently, for any $\epsilon>\|f-g\|_\infty$ we get $g(x)\in A^\epsilon$ and therefore $x\in g^{-1}\kl{A^\epsilon}$. Thus $f^{-1}(A)\subseteq g^{-1}\kl{A^\epsilon}$. Conversely, if $f^{-1}(A)=\emptyset$ then the same inclusion trivially holds. Analogously $g^{-1}(A)\subseteq f^{-1}\kl{A^\epsilon}$ can be shown for any $\epsilon>\|f-g\|_\infty$. This directly implies
   \begin{align*}
       d_P(Q(f), Q(g)) &= d_P(f_\ast\mu_0, g_\ast\mu_0)\\
       &= \inf\left\{\,\epsilon >0\,:\, f_\ast\mu_0(A) \leq g_\ast\mu_0(A^\epsilon) + \epsilon \text{ and } \right. \\
       &\qquad\qquad\left. g_\ast\mu_0(A) \leq f_\ast\mu_0(A^\epsilon) + \epsilon \text{ for any }A\in\mathcal{B}(\mathbb{R})\,\right\} \\
       &= \inf\left\{\,\epsilon >0\,:\, \mu_0(f^{-1}(A)) \leq \mu_0(g^{-1}(A^\epsilon)) + \epsilon \text{ and } \right. \\
       &\qquad\qquad\left. \mu_0(g^{-1}(A)) \leq \mu_0(f^{-1}(A^\epsilon)) + \epsilon \text{ for any }A\in\mathcal{B}(\mathbb{R})\,\right\} \\
       &\leq \|f-g\|_\infty.\qedhere
   \end{align*}
 \end{proof}

 \paragraph{ReLU-Invariance:}
 We can now define the family of distributions as
 \[
    p\colon\mathbb{R}^6\to\mathcal{D}(\mathbb{R})\colon \theta \mapsto \kl{Q\circ R}\kl{\theta}.
 \]
 The local Lipschitz continuity of $p$ is clear from the previous steps. It remains to establish the ReLU-invariance.
 \begin{lemma}
 The family of distributions $p$ is ReLU-invariant.
 \end{lemma}
 \begin{proof}
   Let $\theta\in\mathbb{R}^6$ be arbitrary and $p(\theta)\in\mathcal{D}(\mathbb{R})$ the corresponding probability measure. Further, let $f\colon\mathbb{R}\to\mathbb{R}\colon x\mapsto\varrho\kl{wx+b}$ for some $w,b\in\mathbb{R}$ be any ReLU-layer. We need to show that there exists a $\omega\in\mathbb{R}^6$ such that $p(\omega)=f_\ast p(\theta)$.

   By construction $p(\theta) = \kl{Q\circ R}(\theta) = \kl{R(\theta)}_\ast\mu_0 = \kl{f_\theta}_\ast\mu_0$, where $f_\theta=R(\theta)$ is a ReLU neural network. Clearly also $f\circ f_\theta$ is a ReLU neural network, and since any ReLU network can be rewritten as a three-layer network, there exists $\omega\in\mathbb{R}^6$ such that $R(\omega)=f\circ f_\theta$. Finally,
   \[
        p(\omega)=(f\circ f_\theta)_\ast \mu_0 = f_\ast \kl{f_\theta}_\ast \mu_0 = f_\ast p(\theta)
    \]
    yields the claim.
 \end{proof}

 \subsection{Families of Distributions with Finite Support}\label{sec:counterproof_ii}
 In this section we will prove part \ref{itm:counter_ii} of \cref{thm:counter} and show that ReLU-invariant families of probability distributions exist in which all distributions are finitely supported.
 A finitely supported Radon probability measure $\mu\in\mathcal{D}(\mathbb{R}^d)$ can be expressed as a mixture of finitely many Dirac measures,
 \[
 \mu = \sum_{i=1}^{N} c_i \delta_{\mathbf{a}_i},
 \]
 for some $N\in\mathbb{N}$, $\mathbf{a}_i\in\mathbb{R}^d$, and $0\leq c_i \leq 1$ with $\sum_{i=0}^N c_i = 1$. For any measurable function $f\colon\mathbb{R}^d\to\mathbb{R}^d$ the pushforward of such a mixture is simply the mixture of the individual pushforwards, that is
 \[
    f_\ast\kl{\mu} =  \sum_{i=1}^{N} c_i f_\ast\kl{\delta_{\mathbf{a}_i}}.
 \]
 The pushforward of a Dirac measure is again a Dirac measure, more precisely we simply have $f_\ast\kl{\delta_{\mathbf{a}_i}} = \delta_{f(\mathbf{a}_i)}$, and therefore
 \[
    f_\ast\kl{\mu} =  \sum_{i=1}^{N} c_i \delta_{f(\mathbf{a}_i)}.
 \]
 In particular, for a ReLU-layer $f\colon\mathbb{R}^d\to\mathbb{R}^d\colon\mathbf{x}\mapsto\varrho\kl{\mathbf{W}\mathbf{x}+\mathbf{b}}$ we get
  \[
    f_\ast\kl{\mu} =  \sum_{i=1}^{N} c_i \delta_{\varrho\kl{\mathbf{W}\mathbf{a}_i+\mathbf{b}}}.
 \]
 Altogether, this shows that for any $N\in\mathbb{N}$ the set
 \[
    \mathcal{D}_N = \skl{\,\mu\in\mathcal{D}(\mathbb{R}^d)\,:\,|\supp(\mu)|\leq N\,}
 \]
 is ReLU-invariant in any dimension $d\in\mathbb{N}$. It can be described by $n=N(d+1)$ parameters by simply specifying all the locations $\mathbf{a}_i$ and mixture coefficients $c_i$. More precisely, let $\Delta_N=\skl{\,\mathbf{c}\in[0,1]^N\,:\,\sum_i c_i=1\,}$ and
 \[
 \Omega_N=\underbrace{\mathbb{R}^d\times\dots\times\mathbb{R}^d}_{N\text{ times}}\times\Delta_N\subseteq \mathbb{R}^n
 \]
 and
 \[
 p_N\colon\Omega_N\to\mathcal{D}(\mathbb{R}^d)\colon\kl{\mathbf{a}_1,\dots,\mathbf{a}_N,\mathbf{c}}\mapsto \sum_{i=1}^{N} c_i \delta_{\mathbf{a}_i}.
 \]
 Then $p_N$ is Lipschitz continuous and parametrises $\mathcal{D}_N$.

 \subsection{Families of Distributions Without Local Lipschitz Continuity}\label{sec:counterproof_iii}
 In this section we will prove part \ref{itm:counter_iii} of \cref{thm:counter}. If we omit the local Lipschitz continuity of the parametrisation we can construct a ReLU-invariant one-parameter family of distributions in $\mathbb{R}^d$ for any dimension $d\in\mathbb{N}$. We proceed similar to the ideas in \cref{sec:counterproof_i}, however we need to take care of the fact that the set of all ReLU networks can not easily be described by a constant number of parameters, unlike in the one-dimensional case.

 We define a family of distributions each of which can be described by a finite but arbitrarily large number of parameters. We use space-filling curves to fuse these arbitrarily many parameters into a single parameter.

 The main idea is quite straightforward and analogous to \cref{sec:counterproof_i}: If we have an arbitrarily large number of parameters available to describe each distribution in the family we can simply use these parameters to specify the weights of the neural network transformations with respect to which we want to obtain invariance. The challenge in the general $d$-dimensional setting is to show that this leads to a continuous parametrisation.

 We start by clarifying what we mean by a finite but arbitrarily large number of parameters. We use the conventional notation $\mathbb{R}^\omega$ for the product of countably many copies of $\mathbb{R}$ equipped with the product topology. In other words, $\mathbb{R}^\omega$ is the space of all real-valued sequences. Further we denote by
 \[
 \mathbb{R}^\infty = \skl{\,\kl{x_1,x_2,\dots}\in\mathbb{R}^\omega\,:\,x_i\neq 0\text{ for only finitely many }i\,}
 \]
 the subset of sequences that are eventually zero. The space $\mathbb{R}^\infty$ will serve as the parameter space. Each element in $\mathbb{R}^\infty$ effectively uses only a finite number of parameters (since all the remaining ones are zero), however this number can be arbitrarily large.

 We proceed in three steps. We first introduce a one-parameter description $\Gamma\colon\mathbb{R}\to\mathbb{R}^\infty$ of the finitely-but-arbitrarily-many parameter set $\mathbb{R}^\infty$. We then use a subset $\Omega_\infty\subseteq\mathbb{R}^\infty$ to parametrise ReLU neural networks as functions in $C(K,\mathbb{R}^d)$ using a realisation map $R\colon\Omega_\infty\to C(K,\mathbb{R}^d)$. Finally, we use a pushforward map $Q\colon C(K,\mathbb{R}^d)\to\mathcal{D}(\mathbb{R}^d)$ mapping neural network functions $f\in C(K,\mathbb{R}^d)$ to the pushforward of a fixed \emph{prototype} measure $\mu_0\in\mathcal{D}(\mathbb{R}^d)$ under $f$ to generate the family of probability distribution. As before the prototype measure is assumed to be supported within the domain $K$. The one-parameter family $p\colon\Omega\to\mathcal{D}(\mathbb{R}^d)$ is then given as $p=Q\,\circ R\,\circ\Gamma$ for some suitable domain $\Omega$. We will now discuss each of the steps in more detail.

 \paragraph{Space-Filling Curves:}
 We want to continuously describe $\mathbb{R}^\infty$ with a single parameter. This can be achieved by gluing together continuous space-filling curves from the unit interval to cubes of arbitrary dimension. The existence of such maps is guaranteed by the Hahn-Mazurkiewicz theorem, see for example \cite[Chapter 6.8]{Sagan1994spacefilling}.
 \begin{lemma}\label{lem:space-filling}
 There exists a continuous and surjective function $\Gamma\colon\mathbb{R}\to\mathbb{R}^\infty$.
 \end{lemma}
 The construction of $\Gamma$ is deferred to \ref{apx:space-filling}.

 \paragraph{Neural Network Parametrisation:}
 We need to extend the realisation map from \cref{sec:counterproof_i} to arbitrary dimensions and arbitrary numbers of layers. Further we will show that the results in \cite{prv2018nn-topology-arxiv} concerning continuous network parametrisations can be extended to networks of arbitrary depth.

 A ReLU neural network in $d$ dimensions is characterised by the number of its layers as well as $d^2+d$ parameters for the weights and biases for each of the layers. Using the first component in $\mathbb{R}^\infty$ to encode the number of layers we can parametrise all such neural networks with the set
 \[
    \Omega_\infty = \skl{\,\kl{L,x_1,x_2,\dots}\in\mathbb{R}^\infty\,:\,L\in\mathbb{N}\text{ and }x_i=0\text{ for all }i>L\kl{d^2+d}\,} \subseteq \mathbb{R}^\infty.
 \]
 It will be convenient to introduce notations for the subsets of parameters with a fixed number of layers. For $L\in\mathbb{N}$ we write $\Omega_L = \skl{L}\times \mathbb{R}^{L\kl{d^2+d}}\times\skl{0}^\infty\subseteq\Omega_\infty$ and observe that
 \[
    \Omega_\infty = \bigcup_{L\in\mathbb{N}} \Omega_L.
 \]
 In fact this is the partition of $\Omega_\infty$ into its connected components.

 For any $\theta=(L,x_1,x_2,\dots)\in\Omega_\infty$ we can regroup its leading non-zero components into alternating blocks of size $d^2$ and $d$ and write \[\theta=(L,\mathbf{W}_1,\mathbf{b}_1,\mathbf{W}_2,\mathbf{b}_2,\dots,\mathbf{W}_L,\mathbf{b}_L,0,\dots)\]
 with $\mathbf{W}_i\in\mathbb{R}^{d\times d}$ and $\mathbf{b}_i\in\mathbb{R}^d$. As before let $K\subseteq\mathbb{R}^d$ be any non-empty compact domain and denote the function that is the ReLU neural network realisation of the weights $\skl{\mathbf{W}_1,\dots,\mathbf{W}_d}$ and biases $\skl{\mathbf{b}_1,\dots,\mathbf{b}_d}$ as $f_\theta\colon K\to\mathbb{R}^d$. Then we can define the extended realisation map
 \[
    R\colon\Omega_\infty\to C(K,\mathbb{R}^d)\colon \theta\mapsto f_\theta.
 \]
 \begin{lemma}
    The extended realisation map $R\colon\Omega_\infty\to C(K,\mathbb{R}^d)$ is continuous.
 \end{lemma}
 \begin{proof}
  To show the continuity of $R\colon\Omega_\infty\to C(K,\mathbb{R}^d)$ it suffices to show the continuity on all connected components $\Omega_L$ of the domain. But for each fixed $L\in\mathbb{N}$ we know from \cite[Proposition 5.1]{prv2018nn-topology-arxiv} that the realisation map is continuous from $\Omega_L$ to $C(K,\mathbb{R}^d)$. This already proves the overall continuity of $R$.
 \end{proof}

 \paragraph{Pushforward Mapping:} We need to extend our pushforward map from \cref{sec:counterproof_i} to arbitrary dimension. As before let $\mu_0\in\mathcal{D}(\mathbb{R}^d)$ be any probability measure on $\mathbb{R}^d$ with $\supp(\mu_0)\subseteq K$. We define the pushforward function analogous to before as
 \[
    Q\colon C(K,\mathbb{R}^d)\to\mathcal{D}(\mathbb{R}^d)\colon f\mapsto f_\ast\mu_0.
 \]
 The restriction of $Q$ to the subset of ReLU neural network functions will be used to generate the invariant family of distributions.
 \begin{lemma}
    The pushforward function $Q\colon C(K,\mathbb{R}^d)\to\mathcal{D}(\mathbb{R}^d)$ is continuous.
 \end{lemma}
 The proof works exactly as in the one-dimensional case in \cref{sec:counterproof_i}.

 \paragraph{ReLU-Invariance:}
 We can now define the one-parameter family of distributions as
 \[
    p\colon\Omega\to\mathcal{D}(\mathbb{R}^d)\colon \theta \mapsto \kl{Q\circ R\circ\Gamma}\kl{\theta}
 \]
 with the domain $\Omega=\Gamma^{-1}(\Omega_\infty)$ chosen so that $\Gamma$ maps it exactly to the feasible neural network parameters $\Omega_\infty$ and thus $R\circ\Gamma$ maps it to all ReLU neural networks in $d$ dimensions. The continuity of $p$ is clear from the three previous steps. It remains to establish the ReLU-invariance.
 \begin{lemma}
 The one-parameter family of distributions $p$ is ReLU-invariant.
 \end{lemma}
 \begin{proof}
   Let $\theta\in\Omega$ be arbitrary and $p(\theta)\in\mathcal{D}(\mathbb{R}^d)$ the corresponding probability measure. Further, let $f\colon\mathbb{R}^d\to\mathbb{R}^d\colon\mathbf{x}\mapsto\varrho\kl{\mathbf{W}\mathbf{x}+\mathbf{b}}$ for some $\mathbf{W}\in\mathbb{R}^{d\times d}$ and $\mathbf{b}\in\mathbb{R}^d$ be any ReLU-layer. We need to show that there exists a $\omega\in\Omega$ such that $p(\omega)=f_\ast p(\theta)$.

   By construction we have that $p(\theta) = \kl{Q\circ R\circ \Gamma}(\theta) = \kl{R(\Gamma(\theta))}_\ast\mu_0 = \widetilde{f}_\ast\mu_0$, where $\widetilde{f}=R(\Gamma(\theta))$ is a ReLU neural network. Clearly also $f\circ \widetilde{f}$ is a ReLU neural network, so there exists $\omega\in\Omega$ such that $R(\Gamma(\omega))=f\circ \widetilde{f}$. Finally,
   \[
        p(\omega)=(f\circ \widetilde{f})_\ast \mu_0 = f_\ast \widetilde{f}_\ast \mu_0 = f_\ast p(\theta)
    \]
    yields the claim.
 \end{proof}

\section{Discussion}\label{sec:discussion}
Our work gives a comprehensive characterisation of distribution families that are invariant under transformations by layers of ReLU neural networks. The only invariant distributions are either sampling distributions or rather degenerate and elaborately constructed distributions that are infeasible for practical applications. This justifies the use of approximation schemes such as Assumed Density Filtering.

We have limited our analysis to the class of functions defined as layers of ReLU neural networks. Similarly, one could ask the same question for other function classes.

Considering different commonly used activation functions, such as $\tanh$, the logistic function, or the Heaviside function, instead of ReLUs is one variation that comes to mind. Since the Heaviside function is not continuous and has a discrete and finite range, it transforms any distribution into a distribution with finite support. So clearly also in this case the only invariant distributions can be sampling distributions. For smooth activation functions, such as $\tanh$ and the logistic function, we are not aware of any characterisation of invariant distributions. Extending our proof strategy to this scenario is not straightforward as it relies on specific properties of the ReLU. We leave this question open for future research.

A second variation on the class of functions would be to put restrictions on the weight matrices or bias vectors. This could be either general constraints, for example non-negativity or positive-definiteness, or more specific restrictions by allowing the weight matrices only to be chosen from a small set of allowed matrices. Both kinds of constraints arise in the
context of iterative reconstruction methods for solving inverse and sparse coding problems \cite{hoyer2002nonneg} or in corresponding unrolled and learnable iterative algorithms \cite{gregor_learning_2010,pock2017varnets}. Here, the restrictions on the weight matrices often stem from physical constraints of the model describing the inverse problem.
We comment on both aspects of this variation in \ref{apx:restricted_weights}.
A brief investigation shows that for any finite or even countable set of weights and biases there exist indeed parametrisations of distributions that circumvent the three restrictions (R1)--(R3) in \cref{thm:main}. These, however, are purely of theoretical interest, since they rely on calculating infinitely many distributions a priori and are thus not useful in practice.

\section*{Acknowledgements}
We thank Peter Jung for providing helpful feedback during the process of writing the first draft of this manuscript, as well as Mones Raslan and John Sullivan for fruitful discussions during the earlier stages of preparing this work.
J. M. and S. W. acknowledge partial support by DFG-GRK-2260 (BIOQIC).

\clearpage

\bibliographystyle{abbrv}