Intrinsic Information-Theoretic Analysis of ReLU Nets

Introduction

Information-theoretic accounts of neural networks span training dynamics [1, 42], loss design [2, 23], and representation analysis [9, 15], but their empirical validation is brittle: classical non-parametric estimators (binning, kernel density, $k$-nearest neighbours [24, 30]) and variational neural estimators [2, 21] all depend on user-chosen hyperparameters — bin widths, bandwidths, neighbourhood radii, or auxiliary architectures — whose biases dominate the signal in deep networks [11, 15, 33]. Can we compute information-theoretic quantities intrinsic to the network, without such choices?

Continuous Piecewise Linear (CPWL) networks [16, 18] — including ReLU networks [13] — partition the input space $\mathcal{X}$ into convex linear regions on which the network acts as an affine map (Figure 1). Each region $\omega$ is uniquely identified by the binary activation pattern $\pi_\omega$ obtained by thresholding pre-activations at zero, so every input carries a discrete region label $\Pi$ alongside its continuous representation $T = \phi_\Pi(X)$, and information-theoretic quantities involving $\Pi$ are well-defined on a finite alphabet without any discretisation choice. We use this to develop a framework for the geometry-based information of CPWL networks. The main contributions are:

Preprint.

• Information decomposition Theorem 3.1 decomposes $I(Y;T)$ into geometrically interpretable parts, where one is a discrete and estimable measure of routing information. Other terms are described in Proposition 3.2 and Corollaries 3.3 and 3.4.

• Routing estimator (§ 3.5). The plug-in estimator is computed from the partition and, in the well-sampled regime, approximates $I(Y;\Pi)$. When it saturates, the functional-equivalence quotient (Proposition 3.5) introduces a coarsening of $\Pi$ that merges regions implementing nearly identical linear operators under a Frobenius norm criterion. This extends the estimator into the fine-resolution regime where the number of regions approaches the number of samples.

We complement the theory with an empirical study (§ 5) on three structurally distinct datasets. Its calibration part is an internal-consistency check of the estimator: does the plug-in estimate on the network’s own routing alphabet agree with three established MI baselines applied to the same networks and data? It does not test whether $I(Y;\Pi)$ approximates $I(Y;T)$; the two can differ in either direction (Corollary 3.4, § 5.2). A saturation-correction experiment then shows that the functional-equivalence quotient recovers meaningful routing estimates at growing input dimensions.

Implications. Replacing aggregate $I(Y;T)$ with its decomposition turns hyperparameter-conditioned IB observations on CPWL networks — whether information rises or falls with depth, whether folding destroys class content, whether MI saturates before generalisation — into measurable claims about specific terms. This is a first step towards an information-theoretic understanding of neural networks in which the base estimate originates from the network geometry.

Related Work

Geometry of CPWL networks. The combinatorial study of linear regions in ReLU networks was established by Montúfar et al. [29] and Pascanu et al. [35], with refined bounds and empirical counts by Hanin and Rolnick [18]; Serra et al. [41] — the latter showing that practical region counts are far below the theoretical maximum — and generalised to arbitrary CPWL activations via convex-partition arrangements by Goujon et al. [16], the framework we adopt in § 3.1. Algorithmic extraction has been pursued through 2D-specialised exact visualisation [22], edge-subdivision [3], and polyhedral marching that prunes unreachable regions [27, 39, 47]. MI estimation in deep networks. The information-bottleneck framework [34, 42, 44, 45] views representation learning as a tradeoff between $I(X;T)$ and $I(Y;T)$, but Saxe et al. [40] and Goldfeld et al. [15] showed that $I(Y;T)$ is ill-defined for deterministic networks on continuous inputs, with estimates dominated by binning and noise-injection choices. Alternative formulations include sliced MI [14, 49] and variational lower bounds [2, 46]. These replace $I(Y;T)$ by a better-behaved surrogate; sliced MI, for instance, averages the MI of random one-dimensional projections and adds a choice of slice distribution and count. We instead exploit the discrete partition $\Pi$ inherent to CPWL networks, on which $I(Y;\Pi)$ is well-defined on finite alphabets and estimable without such choices. It is not a surrogate for $I(Y;T)$ but one term of its exact decomposition (Theorem 3.1): it lower-bounds $I(Y;\Pi,T)$, and its gap to $I(Y;T)$ is the difference between the within-region information and a routing-loss term $I(Y;\Pi\mid T)$ (Corollary 3.4) that has no analogue in the standard IB framework. Sliced MI on the same networks would be a complementary view of $I(Y;T)$, which we leave to future work. Activation-pattern analysis. It is established that the number of activation patterns encountered in practice is far below the theoretical maximum [19], motivating the data-restricted analysis. Lewandowski et al. [26] introduced a space-folding measure based on Hamming-distance non-monotonicity along input-space paths, and Dhayalkar [8] formalised the global combinatorial structure as the ReLU Transition Graph, and later extended it with region-entropy and KL-based analysis [7], further put into a graph framework by Gaines and Bi [10].

Intrinsic information theoretic analysis of CPWL networks

The inherent partition of CPWL networks assigns every input both a continuous representation $T$ and a discrete region label $\Pi$, and the label takes finitely many values regardless of what the per-region affine maps do. This section uses that structure: we set up notation (§ 3.1), propose the exact information decomposition (§ 3.3) with its structural consequences (§ 3.4), and turn to estimation: a plug-in routing estimator (§ 3.5) and a functional-equivalence quotient that extends its reach into the fine-resolution regime (§ 3.6).

Preliminaries

A CPWL neural network $\mathcal{N}:\mathcal{X}\to\mathcal{Y}$ is the composite function $y=\mathcal{N}(x)=(\sigma_L\circ f_{\theta_L}\circ\cdots\circ\sigma_1\circ f_{\theta_1})(x)$, mapping input $x$ to output $y$. It has affine layers $f_{\theta_l}:\mathbb{R}^{n_{l-1}}\to\mathbb{R}^{n_l}$ parametrized by $\theta_l=\{W^l\in\mathbb{R}^{n_l\times n_{l-1}},b^l\in\mathbb{R}^{n_l}\}$ followed by point-wise ReLU activation $\sigma(z)=\max(0,z)$. We denote $\mathcal{N}^l$ the sub-network truncated after layer $l$, write $z^l=W^lt^{l-1}+b^l$ for its pre-activation and $t^l=\sigma_l(z^l)=\mathcal{N}^l(x)$ for its (post-activation) output. When unambiguous, we drop the superscript $l$.

Regions, activation patterns and partition. A region $\omega\subset\mathcal{X}$ is a subset of the input domain on which $\mathcal{N}^l$ is a single affine function. Each region in layer $l$ is uniquely identified by an activation pattern $\pi_\omega=(\pi_\omega^1,\ldots,\pi_\omega^l)$, the concatenation of per-layer binarisations $\pi_\omega^i=\mathbf{1}[z^i>0]\in\{0,1\}^{n_i}$ for layers $1,\ldots,l$, obtained by thresholding the layer-$l$ pre-activation at zero (Figure 1, panel c). We write $\Pi^l$ for the collection of regions induced by the network through layer $l$, so that $\bigcup_{\omega\in\Pi^l}\omega=\mathcal{X}$ partitions the input domain. In practice, only a small fraction of $\Pi^l$ contain data [18], and all our estimators are evaluated on this subset. Given a dataset $\mathcal{D}\subset\mathcal{X}$ of $N:=|\mathcal{D}|$ samples, the data-supported regions are $\Omega^l_{\mathcal{D}}=\{\omega\in\Pi^l:\omega\cap\mathcal{D}\ne\varnothing\}$ with cardinalities $|\Omega^l_{\mathcal{D}}|\leq\min\{|\Pi^l|,|\mathcal{D}|\}$. Within any $\omega\in\Pi^l$ the network behaves as an affine map $\phi^l_\omega:\omega\to\mathbb{R}^{n_l}$ such that $\mathcal{N}^l|_\omega=\phi^l_\omega$, with

where $\prod_{i=1}^l$ denotes the descending product $Q^l_\omega W^l\cdots Q^1_\omega W^1$ (matching feed-forward composition) and $Q^i_\omega=\operatorname{diag}(\pi^i_\omega)$ selects the active neurons in layer $i$. Construction of $\Omega^l_{\mathcal{D}}$ from a trained network and a dataset reduces to evaluating $\pi_\omega(x)$ at every $x\in\mathcal{D}$ and grouping by pattern (App. § 9).

Probabilistic setup for information bottleneck chain

We treat the input $X$ and label $Y$ as random variables with $|\mathcal{Y}|<\infty$. For layer $l$ and a CPWL with given parameters, both $\Pi=\pi(X)\in\Pi^l$ and $T=\phi(X)\in\mathbb{R}^{n_l}$ are deterministic functions of $X$. The standard bottleneck chain $Y\to X\to T$ lifts to the directed acyclic graph (DAG): $Y\to X\to\Pi$, $X\to T$, $\Pi\to T$ in which $\Pi$ and $T$ are parallel descendants of $X$, not sequential. We use two conditional independencies: $Y\perp(\Pi,T)\mid X$ and $H(\Pi\mid X)=H(T\mid X)=0$. Neither $X \to\Pi\to T$ nor $X \to T \to\Pi$ holds in general: within a region $T$ is nontrivial in $X$, and distinct regions can share $T$-values under rank-deficient $A_\omega$. The special case $l=1$, where $\Pi= \mathbf{1}[T > 0]$, is the only setting in which $T$ determines $\Pi$.

We write $\nu_\omega:= p_{T\mid\Pi}(\cdot\mid\omega)$ for the within-region latent-space density – the density of $T$ inside region $\omega$. By the affine map of Eq. (1), $\nu_\omega$ is supported on the affine image $\phi_\omega(\omega) \subset\mathbb{R}^{n_l}$ (possibly singular when $A_\omega$ is rank-deficient), and the marginal density of $T$ is the mixture $p_T = \sum_\omega p_\Pi(\omega)\,\nu_\omega$. A formal pushforward derivation and rank discussion are in App. § 13.

Decomposition of $I(Y; T)$

Overlapping latent supports and routing loss. The routing-loss term has a clean geometric reading: it is positive precisely when distinct regions place mass on the same $t \in\mathbb{R}^{n_l}$ and disagree about the class distribution there.

Structural inequalities

Unlike the standard information-bottleneck chain, neither of $I(Y;\Pi)$ and $I(Y;T)$ dominates unconditionally; the sign of their difference is a measurable property of the network.

The plug-in routing estimator

Of the four terms in Theorem 3.1, only $I(Y;\Pi)$ is a discrete MI on finite alphabets and directly estimable. Given counts $n_{y,\omega},n_y,n_\omega$ over the data-supported regions $\Omega_{\mathcal{D}}^l$, the plug-in estimator is

In the well-sampled regime $\Pr(\Pi\notin\Omega_D) \ll1$, $\hat{I}$ approximates $I(Y;\Pi)$, and the chain

holds up to the truncation gap, bounded by $H(Y)\Pr(\Pi\notin\Omega_D) + h_2(\Pr(\Pi\notin\Omega_D))$ with $h_2$ the binary entropy (App. § 13.§ 6). $\hat{I}$ is the only quantity that needs to be estimated; every inequality has a clean interpretation through Theorem 3.1. The plug-in estimator is biased upwards at finite $N$. We compare it with four standard bias correction, Miller-Madow [28], Grassberger [17], Chao-Shen [5] and Chao-Wang-Jost [6] (App. § 15.§ 2).

Target and occupancy. The chain uses “$\approx$” rather than “$\le$” because $\hat{I}$, whose empirical frequencies are normalised over the observed regions, targets $I(Y;\Pi\mid\Pi\in\Omega_D)$. Conditioning on the observed regions can move MI in either direction, so $\hat{I}$ is not in general a population lower bound for $I(Y;\Pi)$ (App. § 13.§ 6). To judge which regime an estimate is in, we report the occupancy

together with the singleton fraction, the share of regions that hold a single sample. Both are read off the count table of Eq. (4). At small $\rho$ most regions hold several points and $\hat{I}$ is close to $I(Y;\Pi)$. As $\rho\to1$ almost every region is a singleton, the partition becomes an injective map from samples to regions, and $\hat{I}$ approaches the empirical label entropy whether or not the network has learned class structure. Estimates in this regime are not reliable, and standard bias corrections do not repair them (App. § 15.§ 2). The occupancy $\rho$ is distinct from the compression ratio $\rho_{\mathrm{func}}$ of § 3.6.

Computational cost. Constructing $\Omega_D^l$ takes one forward pass per sample, $O(N \cdot\lvert\theta_{\le l}\rvert)$, which is the ordinary inference cost. Hashing the activation patterns costs $O(N\sum_{i\le l} n_i)$, and filling the count table and evaluating Eq. (4) costs $O(N)$ amortised plus $O(\lvert\Omega_D\rvert\cdot\lvert\mathcal{Y}\rvert)$, with $\lvert\Omega_D\rvert\le N$. The estimator therefore adds overhead linear in $N$ to a forward pass that is needed anyway, and is asymptotically cheaper than the pairwise clustering of the functional quotient (App. § 14).

The fine-resolution regime and the functional quotient

Eq. (5) is informative when $N \gg\lvert\Omega_D\rvert$. For typical networks this fails: $\pi_\omega\in\{0,1\}^{\sum_i n_i}$ admits exponentially many patterns and at modest depth $\lvert\Omega_D\rvert$ approaches $N$. In that regime $\hat{I}$ saturates at $H(Y)$ and the bound is vacuous. The issue is structural: the partition $\Omega_D$ becomes too fine to be meaningfully estimated from $D$. We need a coarsening of $\Pi$ that lower-bounds $I(Y;\Pi)$ and exploits the affine structure of Eq. (1). We propose to group regions that implement the same linear operator.

Functional equivalence. Two regions are functionally equivalent if they implement the same linear operator: $\omega\sim_{\mathrm{func}} \omega' \Longleftrightarrow A_\omega^l = A_{\omega'}^l$. The bias $c_\omega$ is a region-dependent offset and is excluded. The quotient is denoted $\Omega_{D,\mathrm{func}} = \Omega_D^l/\sim_{\mathrm{func}}$.

Proposition 3.5 (Tight lower bound via functional quotient). $I(Y;\Omega_{D,\mathrm{func}}) \le I(Y;\Pi)$, with equality iff $p_{Y\mid\Pi=\omega}$ depends on $\omega$ only through $A_\omega^l$.

Proof and remarks (predicted vs. ground-truth label) are in App. § 13.§ 7.

Computing the quotient. Materializing each $A_\omega^l$ is expensive. We use the active sub-network matrix $\tilde{A}_\omega^l = \prod_{i=l}^{1} W^i[S_\omega^i,S_\omega^{i-1}]$, where $S_\omega^i \subseteq\{1,\ldots,n_i\}$ is the set of active neurons in layer $i$, read directly from $\pi_\omega^l$, with $S_\omega^0 := \{1,\ldots,n_0\}$, and $W^i[S_\omega^i,S_\omega^{i-1}]$ is the corresponding submatrix of $W^i$. $A_\omega^l$ is obtained from $\tilde{A}_\omega^l$ by inserting zero rows at inactive output positions. Since $\tilde{A}_\omega^l$ is a product of $l$ matrices, its Frobenius norm scales with depth. We therefore use the relative Frobenius criterion at tolerance $\varepsilon\ge0$,

which is dimensionless, bounded in $[0,2]$, and scale-invariant across layers and architectures. The endpoints have clean meaning: $\varepsilon= 0$ recovers the raw partition; $\varepsilon\ge1$ merges most regions regardless of architecture. § 5.3 identifies $\varepsilon\in[0.2,0.5]$ as the useful range. Details are in App. § 14.

The functional-quotient estimator. (4) on $\Omega_{D,\mathrm{func}}$ yields $\hat{I}_{\mathrm{func}} := \hat{I}(Y;\Omega_{D,\mathrm{func}})$ satisfying

extending the chain Eq. (5) with an additional estimable lower bound.

Based on the functional-quotient estimator, we propose a new diagnostic tool for understanding compression and redundancy in CPWL networks: The compression ratio

diagnoses how much pattern redundancy the network carries: small values indicate many distinct patterns collapsing onto few linear operators, values near 1 indicate $\hat{I}_{\mathrm{func}} \approx\hat{I}(Y;\Omega_D)$. We track $\rho_{\mathrm{func}}^{l}(\varepsilon)$ as a layerwise diagnostic in § 5.3.

Experimental setup

We evaluate the framework along two axes: § 4.1 calibrates the plug-in routing estimator $\hat{I}(Y;\Omega_D)$ against three established MI baselines on three structurally distinct datasets, and § 4.2 stresses the estimator at growing input dimension on PCA-reduced MNIST and tests whether the functional-equivalence quotient at relative tolerance $\varepsilon$ recovers informativeness when the raw partition saturates at $H(Y)$. For the baselines we use uniform binning ($K = 8$ bins per coordinate) and $k$-means clustering ($K = |Y|$), both with the Miller–Madow correction, and the mixed continuous–discrete KSG estimator ($k = 3$) of [38], the discrete-target adaptation of [24]. Further, we use fully-connected ReLU MLPs trained for 150 epochs with SGD with momentum [43], learning rate $\alpha= 10^{-3}$, momentum $\beta= 0.9$, batch size 32, cross-entropy loss; no weight decay, dropout, or augmentation. Hidden weights use He initialisation [20], output weights Glorot initialisation [12], and biases start at zero. These settings were fixed in advance, not tuned. Five seeds per configuration; the data split is fixed across seeds. Networks are trained on 80% of each dataset; every estimate, routing and baseline alike, is computed on the held-out 20% (MNIST: the standard test split; App. § 10.§ 4). Implementation details, licences and compute are in App. § 16.

Study 1 – Calibration against established estimators

Datasets. Composite – a synthetic 7-class 2D mixture of moons, concentric circles, and Gaussian blobs (10 000 points, of which $N = 2{,}000$ are held out for estimation; $H(Y) = 2.639$ bits; App. § 10). WBC – the UCI Wisconsin Breast Cancer tabular dataset [48], binary, 569 points, of which $N = 114$ are held out; $H(Y) = 0.949$ bits on the held-out points. MNIST FC – MNIST [25] flattened to 784 dimensions, projected to 10 PCA components, and rescaled to $[-1,1]^{10}$; $N = 10{,}000$-image test split, $H(Y) = 3.319$ bits. The three datasets cover label entropies from 0.95 to 3.32 bits.

Architectures. Composite and WBC: six fully-connected ReLU architectures spanning depth $\in\{3,5\}$ and width $\in\{5,9,25\}$. MNIST FC: three architectures $[3,3,3]$, $[5,5,5]$, $[7,7,7]$ – wider or deeper networks enter the saturation regime studied in § 4.2. 75 networks in total.

Study 2 – Saturation and the functional-equivalence quotient

We project MNIST onto its top $d$ principal components for $d \in\{2,3,4,5,10,15,20\}$, train three-hidden-layer ReLU MLPs with widths $\{7,15,25,50\}$, and report estimators at the last hidden layer. The threshold Eq. (7) is evaluated on a grid of 15 values, $\{0,0.01,0.05,0.1,0.2,\ldots,1.0,1.5,2.0\}$; Figure 4 shows the subset $\{0,0.1,0.2,0.3,0.4,0.5,1.0,2.0\}$. This gives $7 \times4 \times5 = 140$ networks. The layerwise compression profile $\rho_{\mathrm{func}}^{l}(\varepsilon)$ is computed on the same grid for the five-layer Study-1 networks of Composite and WBC and for $[5,5,5]$ and $[7,7,7]$ on MNIST FC (40 networks). Training metrics are in App. § 17.

Results

Calibration: pointwise agreement across networks

Figure 2 scatters $\hat{I}(Y;\Omega_D)$ against each of the three baselines, with one point per network: the last hidden layer at the last epoch, estimated on the same held-out points as the baselines ($n = 75$ per panel). Within each dataset cluster the points lie along the diagonal: $\hat{I}$ assigns the same bit counts as the baselines without requiring a bin width, bandwidth, or cluster count. The drop in $r$ for KSG has two sources: (i) MNIST FC, where KSG sits systematically above $\hat{I}_{\mathrm{raw}}$ while binning and $k$-means track it closely, and (ii) Composite, where KSG sits at the $H(Y)$ ceiling for 23 of the 30 networks (within 0.02 bits of $H(Y) = 2.639$) while $\hat{I}$ ranges from 1.777 to 2.629 bits.

Estimator	$\rho< 0.05$ ($n = 41$)	$\rho> 0.3$ ($n = 56$)	$\rho> 0.9$ ($n = 19$)
Miller–Madow	0.114	3.739	5.403
Grassberger	0.029	0.067	0.006
Chao–Shen	0.032	0.139	0.054
Chao–Wang–Jost	0.032	0.185	0.065

Finite-sample regime. The estimator is trustworthy when most regions hold several points, i.e. at small occupancy $\rho$ (Eq. (6)). Most calibration networks are in that regime, but the two width-25 WBC architectures are not: with $N = 114$ held-out points they reach $\rho= 0.900$ and 0.961 at the deepest layer, so their points in Figure 2 sit at $H(Y)$ by construction (App. § 15.§ 1). Standard bias corrections do not repair a saturated partition (Table 1): Miller–Madow over-corrects by several bits, into negative values, while Grassberger and the coverage-based Chao–Shen and Chao–Wang–Jost corrections leave the estimate within 0.069 bits of $H(Y)$ at $\rho> 0.9$ (App. § 15.§ 2). We therefore report $\rho$ with every estimate and treat saturation structurally (§ 5.3).

Calibration: layerwise profile at convergence

On all three datasets in Figure 3, $\hat{I}$ rises monotonically with depth, tracks binning and $k$-means at every layer, and stays at or below $H(Y)$ on every layer of every one of the 75 networks (the chain Eq. (5) is empirically satisfied). KSG behaves differently: it starts near or above the routing estimator at layer 1 and then drifts down with depth on Composite ($\sim2.6 \to2.4$ bits), narrowing or inverting the KSG-to-routing gap by the deepest layer. On MNIST FC, KSG stays at $\sim2.7$–2.8 bits across all three layers while $\hat{I}$ rises from $\sim0.95$ to $\sim2.0$, so the gap narrows from $\sim1.8$ bits at layer 1 to $\sim0.65$ bits at layer 3. Per Corollary 3.4, $I(Y;T)-I(Y;\Pi)=I(Y;T\mid\Pi)-I(Y;\Pi\mid T)$, so a positive KSG-above-routing gap signals that the within-region term $I(Y;T\mid\Pi)$ exceeds the routing-loss term $I(Y;\Pi\mid T)$: in narrow MNIST networks the continuous geometry $T$ within each region carries class information that the coarse routing partition $\Pi$ does not yet resolve, and the gap shrinks as deeper layers route more of that information into $\Pi$.

Saturation and the functional-equivalence correction

Saturation. In Figure 4, for $d \le5$, $\hat{I}_{\mathrm{raw}}$ sits well below $H(Y)$ at widths 7 and 15 (by at least 0.345 bits): the partition is informative. Widths 25 and 50 approach $H(Y)$ already at small $d$ (up to 3.317 bits). As $d$ grows ($d \in\{10,15,20\}$) the raw estimator approaches $H(Y)$ for widths 25 and 50: nearly every sample gets its own region, so the empirical conditional entropy $\hat{H}(Y\mid\Omega_D) \to0$ and $\hat{I}_{\mathrm{raw}} \to\hat{H}(Y)$, whatever the network has learned.

Functional equivalence as correction. Merging starts at small $\varepsilon$: at $\varepsilon=0.1$ and 0.2, $\hat{I}_{\mathrm{func}}$ falls up to 0.175 and 0.384 bits below $\hat{I}_{\mathrm{raw}}$. At moderate $\varepsilon\in\{0.3,0.4,0.5\}$ the quotient restores informativeness only where there are redundant linear operators to collapse: at $d \ge10$ and $\varepsilon=0.5$ the width-15 estimate sits 0.266–0.425 bits below $H(Y)$, the width-25 estimate only 0.027–0.074 bits, and width 50 stays saturated (gap $\le0.004$ bits): nearly every region there implements a distinct linear operator. At large $\varepsilon\in\{1,2\}$ the quotient over-compresses; at $d=2$, $\varepsilon=2$ keeps 0.732–0.828 of the raw value.

Compression-ratio diagnostic and layerwise compression. Figure 5 traces $\rho_{\mathrm{func}}(\varepsilon)$ at the last hidden layer. Small-$d$ networks ($d=2,3$) drop sharply at $\varepsilon=0.1$ – many distinct routing paths share nearly identical linear operators – while large-$d$ networks ($d=15,20$) stay near 1 until $\varepsilon\approx0.2$–$0.3$ before dropping. The transition window aligns with the window in Figure 4 where $\hat{I}_{\mathrm{func}}$ provides informative compression: both diagnostics identify $\varepsilon\in[0.2,0.5]$ as the operating range. Across layers of the calibration networks (Figure 7, App. § 14), $\rho_{\mathrm{func}}^{1}=1$ for every network, and for $\varepsilon>0$ the ratio decreases with depth on average: functional equivalence is a phenomenon of compositions.

Discussion

Calibration against established MI estimators

What the calibration tests. $\hat{I}(Y;\Omega_D)$ is not proposed as a proxy for $I(Y;T)$: by Corollary 3.4 the two can differ in either direction, and the KSG–routing gap discussed below shows that they do. The calibration in § 5.1 is an internal-consistency check of the plug-in estimator: applied to the same networks and data, independent discretisation schemes should recover comparable information, and they do. Agreement between finite-sample estimators cannot rule out a bias they share, and the occupancy $\rho$ (Eq. (6)) marks where that matters. The two width-25 WBC architectures sit at $\rho=0.900$ and $0.961$, where every estimate is pushed towards $H(Y)$ and agreement carries little information (App. § 15.§ 1).

Why the routing estimator tracks the discrete baselines. The Pearson coefficients $r=0.911$ (binning) and $r=0.922$ ($k$-means) over 75 networks place $\hat{I}(Y;\Omega_D)$ within the intrinsic plug-in bias of the discrete baselines, across three structurally distinct datasets and a range from 0 bits up to $H(Y)$. The structural distinction is that $\Omega_D$ is read from the network’s piecewise-linear geometry rather than from an external discretisation: no bin width, bandwidth, or cluster count is tuned, and the chain Eq. (5) fixes what the estimator measures (the routing information $I(Y;\Pi)$, lower-bounding $I(Y;\Pi,T)$).

The KSG–routing gap localises within-region information. KSG estimates $I(Y;T)$ from the continuous representation while $\hat{I}(Y;\Omega_D)$ estimates $I(Y;\Pi)$, so by Corollary 3.4 a positive KSG-above-routing gap reads as $I(Y;T\mid\Pi)>I(Y;T)$. The within-region geometry $T\mid\Pi$ carries class information that the partition $\Pi$ has not yet routed. The narrowing of the gap with depth on all three datasets, and the inversion on Composite where the routing estimator overtakes KSG by layer 5, tracks $\Pi$ becoming more class-discriminative. KSG’s known overestimation in low-dimensional continuous representations [15, 24] likely contributes some additional gap on MNIST FC; disentangling the two contributions requires an independent estimator of $I(Y;T\mid\Pi)$.

The $H(Y)$ ceiling holds. Across all 75 networks, every hidden layer, and four estimators, no estimate exceeds $H(Y)$ in Figure 3; the chain Eq. (5) is empirically tight on MNIST FC and WBC and approached on Composite, across an order of magnitude of label entropy.

Functional equivalence as both correction and diagnostic

Saturation correction. When the input dimension is large enough that $|\Omega_D|$ approaches $N$, almost every region is a singleton, $\hat{I}_{\mathrm{raw}} \approx\hat{H}(Y)$, and the chain Eq. (5) carries no information. The functional quotient at $\varepsilon\in[0.3, 0.5]$ collapses routing paths that implement nearly identical linear operators and pulls the estimate back below $H(Y)$ when the network has such redundancy (widths 7–15 in Figure 4); it cannot when nearly every region implements a distinct operator (width 50). Inside the window the estimate still moves with $\varepsilon$, by up to 0.304 bits between $\varepsilon= 0.3$ and 0.5. The criterion Eq. (7) is bounded in $[0, 2]$, dimensionless, and depth- and architecture-invariant, so $\varepsilon$ is a reportable choice rather than a tuned hyperparameter; Figure 5 provides an audit of whether the chosen $\varepsilon$ is in the merging window.

$\rho_{\mathrm{func}}^{l}(\varepsilon)$ as a structural fingerprint. Two structural facts shape Figure 7: $\rho_{\mathrm{func}}^{1} = 1$ identically (a single ReLU layer maps every activation pattern to a distinct linear operator), and on average $\rho_{\mathrm{func}}^{l}$ decreases in $l$ at any fixed $\varepsilon> 0$ (each composition opens new opportunities for products of selection matrices to coincide). The dataset ordering – Composite most compressible, MNIST FC PCA-10 most resistant at small $\varepsilon$ – reads as how much pattern redundancy each network develops when the data supports it: Composite admits a compact functional description by layer 5, while MNIST FC PCA-10, already in the saturation regime, requires non-trivial $\varepsilon$ before any merging. $\rho_{\mathrm{func}} \approx1$ throughout means $\hat{I}_{\mathrm{func}} \approx\hat{I}_{\mathrm{raw}}$ and the quotient is not worth computing; large drops signal that it is. On MNIST FC the routing and functional estimates separate as $\rho_{\mathrm{func}}$ drops, as the diagnostic predicts.

Limitations and scope The MNIST FC calibration was obtained on networks narrow enough that $|\Omega_D|$ remains below $N$. For wider or deeper networks the raw estimator is pushed by $H(Y)$, and the functional quotient is required (Study 2). Moreover, the plug-in estimator Eq. (4) is bounded away from $I(Y;\Pi)$ by $H(Y)\Pr(\Pi\notin\Omega_D) + h_2(\Pr(\Pi\notin\Omega_D))$ (App. § 13.§ 6). The second term dominates near the saturation threshold and surfaces in Figure 4 as the closing gap between $\hat{I}_{\mathrm{raw}}$ and $H(Y)$ at large $d$. A bound that tracks the occupancy $\rho$ explicitly would sharpen this; until then we report $\rho$ with every estimate, and by that measure the two width-25 WBC calibration architectures ($\rho\ge0.900$) lie outside the reliable regime (§ 5.1). High routing information is also not by itself evidence of learned class structure. With randomly permuted training labels (App. § 15.§ 6), the narrow $[5, 5, 5]$ Composite network’s information about the true labels falls from 1.949 to 1.411 bits, and the functional quotient discounts it further (1.828 to 1.096 bits). The wide $[25, 25, 25]$ network keeps 2.433 bits raw and 2.329 bits after the quotient, close to $H(Y) = 2.693$ bits, at 15.9% true-label test accuracy: a fine partition of spatially coherent classes correlates with the class whether or not the network has learned it. The quotient remains a valid lower bound by the data-processing inequality, but it removes this effect only at moderate capacity, so a high estimate should be read together with $\rho$ and, where possible, a label-permutation control. Lastly, the quotient estimate depends on $\varepsilon$ even inside the window $[0.2, 0.5]$: small $\varepsilon$ already merges some regions ($\rho_{\mathrm{func}}(0.1) \ge0.717$) and $\varepsilon\ge1$ over-merges. We report $\varepsilon$ with every estimate and audit it through $\rho_{\mathrm{func}}$ rather than treat it as an opaque hyperparameter.

Conclusion and Outlook

We introduced a framework for the internal information geometry of CPWL networks: an exact decomposition of $I(Y;T)$ into routing, within-region, and routing-loss terms; a hyperparameter-free plug-in routing estimator $\hat{I}(Y;\Omega_D)$ from the piecewise-linear partition; and a functional-equivalence quotient extending the estimator into the fine-resolution regime, at the cost of one bounded parameter $\varepsilon$ that is reported and audited through $\rho_{\mathrm{func}}(\varepsilon)$. As an estimator sanity check, $\hat{I}$ agrees with binning and $k$-means ($r = 0.911, 0.922$) and remains at or below $H(Y)$ across all 75 networks; the KSG-above-routing gap localises within-region information rather than estimator disagreement; and the functional quotient recovers informative estimates on PCA-reduced MNIST where the raw estimator saturates at $H(Y)$. Future work includes estimating the decomposition gap directly and extending the analysis to convolutional and residual architectures, where nearly every data point occupies its own region.

Acknowledgments and Disclosure of Funding

The authors would like to thank associate professor Marta Panizzut at UiT - The Arctic University of Norway for fruitful discussions on algebraic geometry and tropical algebra.

This work was supported by the Research Council of Norway through its Centre of Excellence Integreat – The Norwegian Centre for knowledge-driven machine learning, project number 332645.

Ali Ramezani-Kebrya was supported by the Research Council of Norway through FRIPRO Grant under project number 356103.

Detailed Formulation of Linear Regions for CPWL Networks

Notation. This appendix separates the contribution of a single layer from the partition through that layer. We write $\Pi^l$ for the partition induced by the hyperplanes of layer $l$ alone and $\Pi^{\le l}$ for the common refinement through layer $l$; the latter is the partition written $\Pi^l$ in the main text. As there, $\pi^i_\omega$ is the activation pattern of region $\omega$ at layer $i$, and $\Omega$ is reserved for data-supported subsets $\Omega_{\mathcal{D}}^{\le l} \subset\Pi^{\le l}$. Matrices are written without boldface. The symbol $S^i_\omega$ (§ 3.6) denotes the set of active neurons at layer $i$ in region $\omega$; to avoid collision we use $\sigma^i_\omega$, $\Sigma^i_\omega$ for the sign vector and sign diagonal introduced below.

Let a CPWL network $\mathcal{N}:\mathcal{X}\to\mathbb{R}^k$ with $L$ layers be expressed as

with weights $W^l$, biases $b^l$, and piecewise-linear activation $\sigma$. Each layer induces a partition $\Pi^l$ of the input space into convex linear regions; the arrangement through layer $l$ is the common refinement

Per-region affine form. Within $\omega\in\Pi^{\le l}$ the restriction $\mathcal{N}^l|_\omega$ is affine: $\mathcal{N}^l(x)=A^l_\omega x+c^l_\omega$ with

where the descending product convention $\prod_{i=l}^{1}M_i:=M_lM_{l-1}\cdots M_1$ matches feed-forward composition, and $Q^i_\omega=\operatorname{diag}(\pi^i_\omega)$ is the diagonal matrix that selects active neurons at layer $i$ for region $\omega$. Equivalently, $A^l_\omega$ and $c^l_\omega$ satisfy the recursion

The quantities $W^lA^{l-1}_\omega$ and $W^lc^{l-1}_\omega+b^l$ are the effective weights and biases of layer $l$ restricted to region $\omega$.

Projected inequalities. Each region $\omega\in\Pi^{\le l}$ can be defined by linear inequalities in the input space derived from the hyperplanes of every layer. Construct a sign vector from the activation pattern by mapping $(0,1)\to(1,-1)$

Under this convention an active neuron ($\pi^l_{\omega,i}=1$, with pre-activation $z_i>0$) has $\sigma_i=-1$, so the sign-corrected pre-activation $\sigma_i z_i\le0$. The layer-$l$ constraint defining the closure of $\omega$ (the boundaries $z_i=0$ have measure zero) is therefore

Aggregating the constraints from all layers up to $l$ fully characterizes $\omega$. Writing $D^l_\omega=\Sigma^l_\omega W^lA^{l-1}_\omega$ and $g^l_\omega=-\Sigma^l_\omega(W^lc^{l-1}_\omega+b^l)$, the region takes the half-space form

Data-supported regions in practice. None of the estimators needs the half-space form or an enumeration of $\Pi^{\le l}$. For each $x\in\mathcal{D}$ a single forward pass records the patterns $\pi^1(x),\ldots,\pi^l(x)$; their concatenation is hashed, and points with the same hash lie in the same region of $\Pi^{\le l}$. The distinct hashes are $\Omega_{\mathcal{D}}^{\le l}$, and counting them per label gives the table $n_{y,\omega}$ of Eq. (4). The pattern of the first point seen in a region is stored as its representative and is the input to the functional quotient (App. § 14). The half-space form above describes the geometry of a region but is not used by the estimators. In order to find, characterize and use these regions, we have developed the software parx [31].

The 2D Composite Dataset

The Composite dataset is designed to provide a controlled but heterogeneous testbed for partition-based information analysis. The two canonical 2D classification benchmarks — two-moons and concentric circles — each test a single topological feature (curvature in the moons, nested manifolds in the circles); separately they are too uniform to exercise the routing estimator across regimes of varying density, class entropy, and local linear separability. We combine them, together with isotropic Gaussian blobs, into a single seven-class joint distribution that mixes the three structurally distinct generators within a shared two-dimensional input domain.

Generators

The dataset contains $N = 10\,000$ samples across $\lvert\mathcal{Y}\rvert= 7$ classes drawn from three generative components:

Moons (classes 0–1, 40% of samples). We sample $N_1 = 4\,000$ points (split equally between the two classes) from the standard make_moons generator with noise $\sigma= 0.10$, producing two interlocking half-moons, translated by $(+0.5,+0.5)$. The resulting support is concentrated along two one-dimensional curves with mild isotropic spread; the two classes have non-linearly entangled support that requires at least one ReLU bend to separate.

Concentric circles (classes 2–3, 40% of samples). We sample $N_2 = 4\,000$ points (split equally) from the standard make_circles generator with inner-to-outer scale factor 0.75 and noise $\sigma= 0.05$. The raw output, supported on the unit disk, is then scaled by 3.0 and translated by $(+1,+1)$, so that both circles enclose the moons without overlapping them (Figure 6). The two classes form nested closed loops; the inner class cannot be linearly separated from the outer, and neither circle can be separated from the moons by a single hyperplane.

Gaussian blobs (classes 4–6, 20% of samples). We sample $N_3 = 2\,000$ points (667, 667 and 666 per class) from three isotropic Gaussians with $\sigma= 0.3$ centred at $(-2,4)$, $(3,-3)$, and $(5,4)$, outside the outer circle. The three classes are well separated from each other and from the moons and circles. They serve as an “easy” sub-task.

After concatenation and the train/test split (App. § 10.§ 4), the inputs are min-max scaled component-wise to $[-1,1]^2$, with the scaling fitted on the training split, bounding the input domain for stable neural-network initialisation.

Class statistics and label entropy

The three generators produce non-uniform class proportions: the four moons-and-circles classes each carry $p = 1/5$ of the mass, while the three blob classes each carry $p = 1/15$, summing exactly to 1. The resulting class entropy is

which is the data-side ceiling appearing as a horizontal reference line on the Composite panel of Figure 3 in the main text; the stratified held-out split has the same value, 2.639 bits. Note that $H(Y) < \log_2 7 \approx2.807$: the dataset is deliberately not class-balanced, which gives a non-trivial upper bound that the routing estimator must respect.

Properties relevant to partition-based analysis

Three properties of the Composite dataset matter for our experiments:

Heterogeneous local geometry. The three generators have very different supports in the $[-1,1]^2$ domain: the blobs concentrate mass in three small isotropic regions near the corners, the moons stretch along two curves in the centre, and the circles form two concentric rings around the moons. A trained network must produce class-discriminative routing across all three regimes, so the data-supported regions $\Omega_{\mathcal{D}}$ are not all of the same size or shape. This stresses the routing estimator more than a single-generator dataset would.

Rich class entropy at low ambient dimension. With two input dimensions and $N = 2{,}000$ held-out points, every Composite network stays in the well-sampled regime (occupancy $\rho< 0.411$ at the deepest layer, App. § 15.§ 1), and at the same time the seven classes give a non-trivial $H(Y) = 2.64$ bits to recover. Composite is therefore the dataset on which the routing estimator can be followed through a multi-bit range without approaching saturation.

Mixed difficulty per class. The blob classes are separated from the rest by a few linear cuts, while the moons and circles classes require non-linear bends. The dataset therefore mixes sub-tasks that a narrow network solves early with sub-tasks that need more layers and more training.

Train/test handling

Composite and WBC are split once into a training set (80%) and a held-out test set (20%), stratified by class and fixed across seeds and architectures; MNIST FC uses the standard MNIST training and test splits (Apps. § 11 and § 12). All preprocessing (min-max scaling, and PCA for MNIST) is fitted on the training set and applied unchanged to the test set, whose values are clipped to $[-1, 1]$. Networks are trained on the training set only. Every information-theoretic quantity in the paper — the counts $n_{y,\omega}$, $n_y$, $n_\omega$ of the routing estimator, the data support $\Omega_D$ of the functional quotient, and all three baselines — is computed on the held-out set, so $N = |\mathcal{D}|$ is 2,000 for Composite, 114 for WBC and 10,000 for MNIST FC. The partition is thus evaluated on points it was not fitted to, at the cost of smaller count tables, and hence of larger occupancy $\rho= |\Omega_D|/N$. App. § 15.§ 3 compares this with the in-sample protocol on all points: at the deepest layer the plug-in estimate moves by at most 0.020 bits on Composite and 0.031 bits on WBC.

Figure 6 shows a representative realisation of the dataset.

The WBC Dataset

The Wisconsin Breast Cancer (Diagnostic) dataset [48] contains 569 samples with 30 real-valued features: the mean, standard error and largest value of ten characteristics of cell nuclei measured in digitised images of fine-needle aspirates. The binary label is the diagnosis, with malignant (212 samples) as the positive class and benign (357) as the negative class, so $H(Y) = 0.953$ bits on all points and $0.949$ bits on the 114 held-out points. The features are used without dimensionality reduction and min-max scaled to $[-1, 1]^{30}$ as described in App. § 10.§ 4. WBC is the smallest dataset in the study and the only one with a binary label: it tests the routing estimator at low label entropy and small $N$, where the occupancy $\rho$ grows quickly with the width (App. § 15.§ 1). The architectures and the training protocol are the same as for Composite (§ 4.1).

The MNIST FC Dataset

The MNIST FC dataset complements the synthetic Composite and the small clinical WBC datasets used in our main calibration experiments. It introduces a third, structurally distinct setting: a 10-class image classification task whose label entropy ($H(Y) \approx3.32$ bits) sits well above the other two datasets, and whose input is a genuine image rather than a low-dimensional synthetic vector.

Pipeline

We start from the standard MNIST handwritten-digit dataset [25]: 60 000 training images and 10 000 test images of $28 \times28$ greyscale digits, labelled 0–9. We flatten each image to a vector in $[0, 1]^{784}$ and apply the following two-step preprocessing:

1. PCA to ten dimensions. We fit a principal-component analysis on the training set and project both train and test sets to the top $d_{\mathrm{PCA}} = 10$ principal components. The choice of $d_{\mathrm{PCA}} = 10$ is deliberate: it keeps the input dimension low enough that narrow networks stay in the well-sampled regime ($|\Omega_D| \ll N$) while retaining about half (48.8%) of the training-set pixel variance. Larger $d_{\mathrm{PCA}}$ pushes networks of the same width into the fine-resolution regime ($|\Omega_D| \approx N$) of § 3.6 where $\hat{I}$ becomes uninformative; smaller $d_{\mathrm{PCA}}$ discards class-relevant directions, so the input itself carries less information about the label (test accuracy falls below 50% at $d_{\mathrm{PCA}} = 2$, Figure 11). The study of § 4.2 varies $d_{\mathrm{PCA}}$ from 2 to 20 with the same pipeline.

2. Min-max scaling to $[-1, 1]^{10}$. The PCA-projected coordinates are min-max scaled component-wise so that each of the ten dimensions lies in $[-1, 1]$. This matches the Composite preprocessing and provides comparable input scales across datasets.

The PCA basis and min-max bounds are fitted on the training set and applied unchanged to the test set, to avoid information leakage. All information-theoretic quantities reported in the main paper are computed on the 10 000-sample test set, on which $H(Y) = 3.319$ bits (the test split is approximately balanced across digits but not exactly so).

Architectures

We use three narrow fully-connected ReLU architectures of three hidden layers each: $[3, 3, 3]$, $[5, 5, 5]$, and $[7, 7, 7]$. The narrow widths are dictated by the requirement $|\Omega_D| \ll N$ at the deepest layer: with $d_{\mathrm{PCA}} = 10$ the occupancy at the last hidden layer is $\rho= 0.122$ at width 7, but 0.743 at width 15 and 0.954 at width 25 (Table 3), where most test samples occupy their own region and the routing estimator becomes uninformative. Restricting to the well-sampled regime is therefore not a cherry-pick but a stated scope of the framework; the wider networks are the subject of the saturation study (§ 5.3).

Occupancy $\rho$ at the last hidden layer on MNIST ($N = 10\,000$), per architecture and PCA dimension $d$. The calibration networks of Study 1 are [3, 3, 3], [5, 5, 5] and [7, 7, 7] at $d = 10$.

Architecture	$d = 2$	$d = 3$	$d = 4$	$d = 5$	$d = 10$	$d = 15$	$d = 20$
[3, 3, 3]	–	–	–	–	0.004	–	–
[5, 5, 5]	–	–	–	–	0.030	–	–
[7, 7, 7]	0.012	0.028	0.047	0.082	0.122	0.125	0.103
[15, 15, 15]	0.049	0.195	0.342	0.518	0.743	0.774	0.783
[25, 25, 25]	0.124	0.452	0.734	0.874	0.954	0.965	0.979
[50, 50, 50]	0.317	0.826	0.960	0.987	0.998	0.999	1.000

The $[3, 3, 3]$ architecture is included in the calibration scatter of Figure 2 (where it provides additional data points across the lower bit range), but its partition is very coarse ($\rho= 0.004$, about 40 regions for 10 000 points at the deepest layer) and the network reaches only about 60% test accuracy (Figure 11); we therefore exclude it from the layerwise profile of Figure 3. The saturation study of § 5.3 trains $[7, 7, 7]$ and the wider $[15, 15, 15]$, $[25, 25, 25]$ and $[50, 50, 50]$ at every PCA dimension $d \in\{2, 3, 4, 5, 10, 15, 20\}$ (§ 4.2); its results are reported separately in Figure 4 and Figure 5.

Training follows the same SGD protocol as Composite/WBC (150 epochs, $\alpha= 10^{-3}$, momentum 0.9, batch size 32, cross-entropy, no weight decay or dropout). We use five seeds per architecture, giving 15 networks total.

Probabilistic Details for Section 3

Joint factorization and pushforwards

Under the DAG $Y \to X \to(\Pi,T)$ with $\Pi,T$ deterministic functions of $X$ via $\Pi=\pi(X)$ and $T=\phi_{\Pi}(X)=A_{\Pi}X+c_{\Pi}$, the joint law factorizes as

The induced marginals are

and the within-region activation law is the pushforward of the within-region input law through the region’s affine map,

Concretely, for any Borel $E\subseteq\mathbb{R}^{n_l}$, $\Pr(T\in E\mid\Pi=\omega)=\Pr(X\in\phi_{\omega}^{-1}(E)\mid\Pi=\omega)$. The marginal of $T$ is the mixture

Remark on rank. The pushforward $\nu_{\omega}$ is always well-defined as a measure, regardless of $\operatorname{rank}(A_{\omega})$. If $A_{\omega}$ has rank $r<n_l$, then $\nu_{\omega}$ is supported on an $r$-dimensional affine subspace of $\mathbb{R}^{n_l}$ and is singular with respect to Lebesgue measure; if $A_{\omega}=0$, then $\nu_{\omega}=\delta_{c_{\omega}}$. All MI terms involving $\nu_{\omega}$ remain finite and well-defined in every case because $Y$ takes finitely many values; the pathologies of differential entropy that afflict $I(Y;T)$ in general IB analyses do not apply to the decomposition terms of Theorem 3.1.

Per-region form of the decomposition terms

Because $\Pi$ measurably partitions $\mathcal{X}$, three of the four terms in Theorem 3.1 localize per-region. Let $P_{\omega}$ denote orthogonal projection onto $\operatorname{row}(A_{\omega})$ and $P_{\omega}^{\perp}$ onto $\ker A_{\omega}$. Within region $\omega$, $T=A_{\omega}X+c_{\omega}$ depends on $X$ only through $P_{\omega}X$, so conditioning on $T$ inside region $\omega$ is equivalent (up to the affine bijection $A_{\omega}|_{\operatorname{row}(A_{\omega})}$) to conditioning on $P_{\omega}X$.

with $I_{\omega}$ denoting MI under the within-region law $p_{\cdot\mid\Pi=\omega}$. The conditioning on $P_{\omega}X$ in the last line is essential: a region $\omega$ is a bounded convex polytope, so the row-space and kernel components $P_{\omega}X$ and $P_{\omega}^{\perp}X$ are in general statistically dependent under $p_{X\mid\Pi=\omega}$, and $I_{\omega}(Y;P_{\omega}^{\perp}X\mid P_{\omega}X)\ne I_{\omega}(Y;P_{\omega}^{\perp}X)$. The representation insufficiency $I(Y;X\mid\Pi,T)$ vanishes iff $Y\perp P_{\omega}^{\perp}X\mid P_{\omega}X$ under $p_{X\mid\Pi=\omega}$ for every region $\omega$.

The routing loss does not localize per-region because it couples regions via pushforward overlap:

where $S(t)=\{\omega:t\in\operatorname{supp}(\nu_{\omega})\}$.

Proof of Theorem 3.1

Applying the chain rule to $I(Y;\Pi,T)$ in two orders,

Equating and solving for $I(Y;T)$ gives (2).

For (3): since $H(\Pi,T\mid X)=0$, chain rule gives $I(Y;X,\Pi,T)=I(Y;X)$ and $I(Y;X,\Pi,T)=I(Y;\Pi,T)+I(Y;X\mid\Pi,T)$, whence $I(Y;\Pi,T)=I(Y;X)-I(Y;X\mid\Pi,T)$.

Proofs of Corollaries 3.3 and 3.4

Corollary 3.3. By the chain rule, $I(Y;\Pi,T)=I(Y;T)+I(Y;\Pi\mid T)\ge I(Y;T)$, with equality iff $I(Y;\Pi\mid T)=0$. By (3), $I(Y;\Pi,T)=I(Y;X)-I(Y;X\mid\Pi,T)\le I(Y;X)$, with equality iff $I(Y;X\mid\Pi,T)=0$, i.e. iff $Y\perp X\mid(\Pi,T)$, which is sufficiency of $(\Pi,T)$ for $Y$.

Corollary 3.4. Subtracting $I(Y;\Pi)$ from both sides of (2) gives $I(Y;T)-I(Y;\Pi)=I(Y;T\mid\Pi)-I(Y;\Pi\mid T)$, so $I(Y;T)>I(Y;\Pi)$ iff $I(Y;T\mid\Pi)>I(Y;\Pi\mid T)$.

Proof of Proposition 3.2

Setup. From (20), the posterior over regions given $T=t$ is

where densities $\nu_{\omega}(t)$ are taken with respect to $p_T$ as the dominating measure. Define $S(t)=\{\omega:t\in\operatorname{supp}(\nu_{\omega})\}$ and $O=\{t:|S(t)|\ge2\}$; the posterior is non-degenerate exactly on $O$.

KL form. The standard expansion of conditional mutual information gives

where $p_{Y\mid\Pi=\omega,T}(\cdot\mid t)$ is the joint region-and-activation conditional law of $Y$, and $p_{Y\mid T}(\cdot\mid t)=\sum_{\omega''\in S(t)}p_{\Pi\mid T}(\omega''\mid t)p_{Y\mid\Pi=\omega'',T}(\cdot\mid t)$ is its mixture over regions in the active set.

$(\Leftarrow)$. Suppose distinct $\omega,\omega'$ have $\operatorname{supp}(\nu_{\omega})\cap\operatorname{supp}(\nu_{\omega'})$ of positive $p_T$-mass and $p_{Y\mid\Pi=\omega,T=t}\ne p_{Y\mid\Pi=\omega',T=t}$ on a positive-mass subset $B\subseteq O$. On $B$ the mixture $p_{Y\mid T}(\cdot\mid t)$ cannot equal both $p_{Y\mid\Pi=\omega,T=t}$ and $p_{Y\mid\Pi=\omega',T=t}$ since these differ; hence at least one KL term in (28) is positive on $B$, and $I(Y;\Pi\mid T)>0$.

$(\Rightarrow)$. If $I(Y;\Pi\mid T)>0$, then by (28) there is a set $B$ with $p_T(B)>0$ on which the inner sum is strictly positive. For each $t\in B$: (i) positivity requires $|S(t)|\ge2$, so $B\subseteq O$; (ii) positivity requires at least two $\omega,\omega'\in S(t)$ with $p_{Y\mid\Pi=\omega,T=t}\ne p_{Y\mid\Pi=\omega',T=t}$ (else every conditional equals the mixture and all KLs vanish). Both conditions hold on a positive-mass set.

Marginal sufficient condition. The proposition is stated in terms of the conditional laws $p_{Y\mid\Pi=\omega,T=t}$ at the overlap. A simpler sufficient condition is that the overlapping regions have differing marginal region-conditional label laws $p_{Y\mid\Pi=\omega}\ne p_{Y\mid\Pi=\omega'}$, provided the within-region label distribution is independent of $T$, i.e. $Y\perp T\mid\Pi$. Under that independence the two formulations coincide, and the marginal version recovers the older statement that “two overlapping regions with differing class content contribute to routing loss”. In the general case the conditional version is required: two regions can have identical marginal label laws but differ at a specific overlapping coordinate $t$, which is enough for the routing loss to be positive there.

Geometric reading. The proposition factors label-relevant pushforward overlap into two independent ingredients: image overlap ($p_T(O)>0$) — two regions with distinct affine maps produce overlapping activation-space images — and label-content disagreement at the overlap — the conditional label distributions at overlapping points $t$ differ between regions. Routing loss is positive iff both hold simultaneously. Overlap between regions whose conditional class content at the overlap agrees (e.g. regions carrying the same class from disconnected parts of the input distribution) contributes zero to $I(Y;\Pi\mid T)$: such overlap preserves the conditional label distribution at $t$ and is invisible to the routing loss, even though it is visible to geometric measures of partition complexity [26].

Truncation gap of the plug-in routing estimator

The plug-in estimator (4) is the standard empirical MI, obtained by substituting maximum-likelihood counts into the Shannon formula. Its target depends on whether one views the count tables as a sample from $p_{Y,\Pi}$ (giving $I(Y;\Pi)$ in the limit) or as the distribution of $(Y,\Pi)$ restricted to the observed support $\Omega_D$ (giving $I(Y;\Pi\mid\Pi\in\Omega_D)$). The two limits differ by the truncation gap analysed below.

Truncation identity. Let $U=\{\Pi\notin\Omega_D\}$ be the event that a fresh sample hits an unobserved region, and let $\mathbf{1}_U\in\{0,1\}$ be its indicator. The classical truncation identity gives

The plug-in estimator $\hat I(Y;\Omega_D)$, computed from empirical frequencies that sum to one over $\Omega_D$, targets the unscaled term $I(Y;\Pi\mid U^c)$ (the population MI conditional on $\Pi$ being observed) rather than the scaled term $\Pr(U^c)I(Y;\Pi\mid U^c)$ that appears in (29). Hence $\hat I$ is in general not a population lower bound for $I(Y;\Pi)$: because conditioning on a subset can change MI in either direction, $I(Y;\Pi\mid U^c)$ can lie above or below $I(Y;\Pi)$ in different regimes.

Gap bound. For the third term in (29), $0\le I(Y;\Pi\mid U)\le H(Y\mid U)\le H(Y)$, so $\Pr(U)I(Y;\Pi\mid U)\le H(Y)\Pr(U)$. For the first term, $I(Y;\mathbf{1}_U)\le H(\mathbf{1}_U)=h_2(\Pr(U))$ where $h_2$ is the binary entropy, which scales as $-\Pr(U)\log_2\Pr(U)$ as $\Pr(U)\to0$. Together,

In the well-sampled regime $\Pr(U)\to0$ both terms vanish and $\hat I(Y;\Omega_D)\to I(Y;\Pi\mid U^c)\approx I(Y;\Pi)$. The chain (5) therefore holds up to (30); we use “$\approx$” in (5) to mark this approximate (rather than provable) lower-bound relationship. The right-hand inequalities $I(Y;\Pi,T)\le I(Y;\Pi)\le H(Y)$ are exact: the first is Corollary 3.3, the second holds because $Y$ is discrete.

Proof of Proposition 3.5

Let $q:\Pi^l\to\Pi_{\mathrm{func}}$ be the quotient map. Since $q$ is deterministic, $Y\to\Pi\to\Pi_{\mathrm{func}}$ is Markov and the data-processing inequality gives $I(Y;\Pi_{\mathrm{func}})\le I(Y;\Pi)$.

For equality: $I(Y;\Pi)=I(Y;\Pi_{\mathrm{func}})+I(Y;\Pi\mid\Pi_{\mathrm{func}})$, so the two coincide iff $I(Y;\Pi\mid\Pi_{\mathrm{func}})=0$, i.e., iff $Y\perp\Pi\mid\Pi_{\mathrm{func}}$. This holds iff $p_{Y\mid\Pi=\omega}$ is constant on each equivalence class of $\sim_{\mathrm{func}}$, i.e., iff $p_{Y\mid\Pi=\omega}$ depends on $\omega$ only through $A^l_\omega$.

Remark. The condition asks $p_{Y\mid X=\omega}$ to be constant on each functional-equivalence class. Two cases are worth distinguishing.

(i) Predicted label. If $Y$ in the proposition is taken to be the network’s own predicted label $\widehat Y=\arg\max(Y)$, the condition is still not automatic. Two functionally-equivalent regions share the linear part $A^l_\omega$ but in general not the offset $c^l_\omega$, and they receive different inputs, so their images $\phi^l_\omega(\omega)$ differ and the head can assign them different predicted labels. The condition holds for the predicted label when, in addition, equivalent regions have equal offsets and the head is constant on the union of their images.

(ii) Ground-truth label. If $Y$ is the data label, two functionally-equivalent regions can still differ: the regions can lie in disjoint parts of the input space and carry different ground-truth class proportions in $p_{X,Y\mid\Pi=\omega}$ even though the linear operator $A^l_\omega$ is the same. The equality condition then additionally requires $p_{Y\mid X,\Pi=\omega}$ to depend on $X$ only through $A^l_\omega X$ — i.e. no $Y$-relevant variation along $\ker A^l_\omega$. – a property of the data distribution, not of the network alone. All reported estimates use the ground-truth label, for which the quotient is a lower bound but equality is not guaranteed. The gap $I(Y;\Pi) - I(Y;\Pi_{\mathrm{func}}) = I(Y;\Pi\mid\Pi_{\mathrm{func}})$ is what the difference $\hat{I} - \hat{I}_{\mathrm{func}}$ estimates, and it is reported with every quotient estimate (Figure 4 and Table 8). □

Dataset	Architecture	Labels	Test acc.	Train-label acc.	$\rho$	$\hat{I}$	$\hat{I}_{\mathrm{func}}$	$\rho_{\mathrm{func}}$
Composite	$\left[5,5,5\right]$	clean	97.9%	–	0.030	1.949	1.829	0.765
Composite	$\left[5,5,5\right]$	permuted	19.6%	20.8%	0.022	1.411	1.096	0.656
Composite	$\left[25,25,25\right]$	clean	99.5%	–	0.304	2.535	2.362	0.377
Composite	$\left[25,25,25\right]$	permuted	15.9%	22.5%	0.245	2.433	2.329	0.645
WBC	$\left[5,5,5\right]$	clean	97.4%	–	0.184	0.767	0.700	0.680
WBC	$\left[5,5,5\right]$	permuted	62.5%	62.8%	0.211	0.469	0.295	0.709
WBC	$\left[25,25,25\right]$	clean	97.2%	–	0.900	0.949	0.885	0.745
WBC	$\left[25,25,25\right]$	permuted	63.3%	65.6%	0.933	0.946	0.932	0.910

Relation to prior IB literature

The pathology of $I(Y;T)$ for deterministic continuous networks has been discussed by Saxe et al. [40] and Goldfeld et al. [15], who show that estimates of it are dominated by binning choices; Goldfeld et al. [15] inject noise into the network to make the quantity well-posed. Our approach is orthogonal: we analyze the joint $(\Pi,T)$ whose combinatorial component is entropy-regular by construction, giving finite exact identities without altering the model. The price is tracking the routing-loss term $I(Y;\Pi\mid T)$, whose geometric interpretation via Proposition 3.2 is of independent interest.

Clustering Regions by Functional Equivalence

Given $\Omega_D$, network weights, and tolerance $\varepsilon\in[0,2]$, we group regions by the relative Frobenius criterion of (7): $\omega\sim_{\varepsilon} \omega'$ iff $S_{\omega}^{l}=S_{\omega'}^{l}$ and $\lVert\widetilde{A}_{\omega}^{l}-\widetilde{A}_{\omega'}^{l}\rVert_F/\frac{1}{2}(\lVert\widetilde{A}_{\omega}^{l}\rVert_F+\lVert\widetilde{A}_{\omega'}^{l}\rVert_F)\leq\varepsilon$. The active-set equality $S_{\omega}^{l}=S_{\omega'}^{l}$ is required for the matrix difference to be well-defined ($\widetilde{A}_{\omega}^{l},\widetilde{A}_{\omega'}^{l}$ otherwise have incompatible row counts), so regions are first grouped by $S_{\omega}^{l}$ and compared only within a group. Within a group, each quotient class is represented by its first member, and a new region joins the first class whose representative is within $\varepsilon$ of it.

For $\varepsilon=0$ this is exact functional equivalence (raw partition); larger $\varepsilon$ merges progressively more pairs. Because the relation $\sim_{\varepsilon}$ is not transitive when $\varepsilon>0$ (regions $\omega_1,\omega_3$ may be more than $\varepsilon$ apart even if $\omega_1\sim_{\varepsilon}\omega_2\sim_{\varepsilon}\omega_3$), the greedy clustering in Algorithm 1 can produce slightly different partitions depending on the order in which regions are visited. We iterate over $\Omega_D$ in first-encounter order, i.e. the order in which the regions are first reached by the data points, which is deterministic and reproducible. App. G.7 measures the effect of the order over 15 random visiting orders: at $\varepsilon=0.3$ the quotient size varies by a median of 2 classes and at most 9.4% of $\lvert\Omega_D\rvert$, and the Miller–Madow-corrected quotient MI by at most 0.148 bits.

Algorithm 1	Functional-quotient clustering
Require:	$\Omega_D$, weights $\{W^i,\mathbf{b}^i\}_{i=1}^{l}$, tolerance $\varepsilon\in[0,2]$
Ensure:	Quotient $\Omega_{D,\mathrm{func}}$
1:	$\Omega_{D,\mathrm{func}}\leftarrow\varnothing$
2:	for each $\omega\in\Omega_D$ in first-encounter order do
3:	Read $S_{\omega}^{i}$ from $\pi_{\omega}^{i}$ for $i=1,\ldots,l$
4:	$\widetilde{A}_{\omega}^{l}\leftarrow\prod_{i=l}^{1}W^i[S_{\omega}^{i},S_{\omega}^{i-1}]$
5:	Search the classes of $\Omega_{D,\mathrm{func}}$ with active set $S_{\omega}^{l}$, in creation order, for the first class $[\omega']$ whose representative $\omega'$ satisfies $\omega\sim_{\varepsilon}\omega'$
6:	if found then
7:	Add $\omega$ to $[\omega']$
8:	else
9:	Create new class $\{\omega\}$ in $\Omega_{D,\mathrm{func}}$ with representative $\omega$
10:	end if
11:	end for
12:	return $\Omega_{D,\mathrm{func}}$

Complexity. Each comparison $\omega\sim_{\varepsilon}\omega'$ costs $O(\lvert S_{\omega}^{l}\rvert\cdot n_0)$ for the Frobenius distance. Each region is compared only with the representatives of existing classes in its active-set group, so the number of comparisons is at most $\lvert\Omega_D\rvert\cdot\lvert\Omega_{D,\mathrm{func}}\rvert$, with an $O(\lvert\Omega_D\rvert^2)$ worst case at small $\varepsilon$. In our experiments $\lvert\Omega_D\rvert$ reaches 9,997 (MNIST, width 50); for large supports the comparison can be accelerated by locality-sensitive hashing on rounded $\widetilde{A}_{\omega}^{l}$, reducing the cost to near-linear in $\lvert\Omega_D\rvert$.

Operating range. The relative criterion is bounded in $[0,2]$ and admits a clean interpretation: $\varepsilon=0$ recovers the raw partition; $\varepsilon\in[0.2,0.5]$ is the window in which $\rho_{\mathrm{func}}$ drops and the quotient pulls a saturated estimate back below $H(Y)$, provided the network has redundant linear operators to collapse (Sec. 5.3 and Figure 4 and Figure 5); and $\varepsilon>1$ over-merges. Merging already starts below the window, and the estimate still moves with $\varepsilon$ inside it (Sec. 6.2), so $\varepsilon$ is reported with every quotient estimate. Where a single value is reported (Table 8), we use $\varepsilon=0.3$ as a representative operating point.

Layerwise compression ratio. Figure 7 reports $\rho_{\mathrm{func}}^{l}(\varepsilon)$ as a function of layer depth on the calibration networks. $\rho_{\mathrm{func}}^{1}=1$ for every dataset, architecture and $\varepsilon$: the first hidden layer realises every distinct activation pattern as a distinct linear operator, so functional equivalence is a phenomenon of compositions. For $\varepsilon>0$, $\rho_{\mathrm{func}}^{l}$ then decreases with depth on average — deeper compositions admit more activation-pattern combinations whose products of selection matrices coincide — at a rate set by $\varepsilon$. The decrease is not monotone for every network: of the 560 per-network curves (40 networks, 14 values of $\varepsilon>0$), 232 rise somewhere with depth, by at most 0.216.

Finite-Sample Validation of the Routing Estimator

This appendix collects the evidence on when the plug-in routing estimator $\hat{I}(Y;\Omega_D)$ of Eq. (4) can be trusted. All numbers are computed on the held-out points (App. § 10), at the last epoch, and averaged over the five seeds; the tables are generated from the stored results of the released code. The organising quantity is the occupancy

the number of data-supported regions per sample. At small $\rho$ most regions hold several points and the counts $n_{y,\omega}$ are informative; as $\rho\to1$ almost every region holds a single point, the partition becomes an injective map from samples to regions, and $\hat{I}$ is pushed towards $H(Y)$ whatever the network has learned (App. § 13.§ 6).

Occupancy of the evaluated networks

Table 2 and Table 3 report $\rho$ at the deepest layer for every configuration, together with the fraction of regions that hold a single point and the fraction of points that sit alone in their region. The Composite networks stay below $\rho= 0.411$. On WBC, with only $N = 114$ held-out points, the two width-25 architectures are close to saturation: $\rho= 0.900$ for $[25,25,25]$ and $\rho= 0.961$ for $[25,25,25,25,25]$ (0.741 and 0.893 on all 569 points). Their calibration points in Figure 2 and Figure 3 lie in the regime where $\hat{I}$ is not trustworthy and should be read with that caveat. On MNIST, $\rho$ grows quickly with the PCA dimension and the width (Table 3); this is the saturation that the functional quotient addresses in § 5.3. Figure 8 shows the full distribution of region sizes behind these numbers. At low $\rho$ most points sit in regions holding tens of points or more; in the saturated networks (WBC width 25, MNIST widths 25 and 50) almost every point sits in a region of size one or two.

Occupancy at the deepest layer, Composite ($N = 2,000$) and WBC ($N = 114$), held-out points, mean over seeds. Regions: $|\Omega_D|$; singleton regions: fraction of regions holding one point; singleton samples: fraction of points alone in their region.

Dataset	Architecture	Regions	$\rho$	Singleton regions	Singleton samples
Composite	[5, 5, 5]	59.2	0.030	0.123	0.004
	[5, 5, 5, 5, 5]	108.2	0.054	0.138	0.009
	[9, 9, 9]	128.4	0.064	0.133	0.009
	[9, 9, 9, 9, 9]	294.2	0.147	0.291	0.044
	[25, 25, 25]	607.4	0.304	0.438	0.133
	[25, 25, 25, 25, 25]	822.2	0.411	0.540	0.222
WBC	[5, 5, 5]	21.0	0.184	0.438	0.081
	[5, 5, 5, 5, 5]	23.0	0.202	0.530	0.111
	[9, 9, 9]	41.0	0.360	0.627	0.230
	[9, 9, 9, 9, 9]	64.0	0.561	0.705	0.398
	[25, 25, 25]	102.6	0.900	0.931	0.840
	[25, 25, 25, 25, 25]	109.6	0.961	0.970	0.933

Bias corrections

The plug-in estimator is biased upwards at finite $N$. We compare it with four standard corrections, each applied to the same contingency table: Miller–Madow [33], which subtracts $(|\Omega_D| - 1)(|\mathcal{Y}| - 1)/(2N \ln2)$; Grassberger [17]; Chao–Shen [5] and Chao–Wang–Jost [6], which are designed for undersampled categorical data. The last three are computed through the decomposition $\hat{H}(\Omega) + \hat{H}(Y) - \hat{H}(\Omega, Y)$ with each entropy corrected separately. All four are implemented in the released code and tested against the infomeasure package [4]. We do not report the asymptotic NSB estimator [32]. It relies on coincidences among the samples, and their number vanishes as $\rho\to1$, so it is unstable in the regime where a correction would matter. It is computed in the released code. Table 4 averages over the 138 (dataset, architecture, PCA dimension, layer) cells of both studies, grouped by $\rho$.

$\rho$ range	Estimator	Value	Correction	Seed std	Gap to $H(Y)$
$\rho< 0.05$ ($n = 41$)	Plug-in	1.439	0.000	0.148	1.698
	Miller–Madow	1.325	0.114	0.131	1.812
	Grassberger	1.410	0.029	0.144	1.727
	Chao–Shen	1.407	0.032	0.143	1.730
	Chao–Wang–Jost	1.407	0.032	0.143	1.730
$\rho> 0.3$ ($n = 56$)	Plug-in	2.586	0.000	0.025	0.146
	Miller–Madow	−1.153	3.739	0.092	3.886
	Grassberger	2.519	0.067	0.031	0.214
	Chao–Shen	2.447	0.139	0.037	0.285
	Chao–Wang–Jost	2.402	0.185	0.055	0.331
$\rho> 0.9$ ($n = 19$)	Plug-in	2.941	0.000	0.002	0.004
	Miller–Madow	−2.462	5.403	0.038	5.407
	Grassberger	2.935	0.006	0.004	0.010
	Chao–Shen	2.888	0.054	0.027	0.058
	Chao–Wang–Jost	2.876	0.065	0.037	0.069

Three regimes appear. At small occupancy ($\rho< 0.05$) every correction is small (Miller–Madow: 0.114 bits on average) compared with the seed-to-seed variability of the plug-in value (0.148 bits), so the plug-in estimate is reliable. At $\rho> 0.3$ the Miller–Madow correction grows to 3.739 bits on average, up to 6.488 bits, and the corrected value becomes negative (as low as −3.169 bits): a first-order correction applied far outside its regime. The seed-to-seed spread of the plug-in value moves the other way, from 0.148 bits at $\rho< 0.05$ to 0.025 bits at $\rho> 0.3$: the estimate is most stable where it is most biased, so its variance gives no warning of its bias. At $\rho> 0.9$ (19 cells) the other corrections are small, but for the wrong reason: their corrected values stay within 0.010 (Grassberger), 0.058 (Chao–Shen) and 0.069 (Chao–Wang–Jost) bits of $H(Y)$, the same ceiling as the plug-in value. When almost every region holds one point, the region marginal and the joint table have nearly the same support, so the corrections of $\hat{H}(\Omega)$ and $\hat{H}(\Omega, Y)$ are nearly equal and cancel in the difference. Miller–Madow fails visibly, Grassberger and the coverage-based corrections fail silently. No statistical correction we tried turns a saturated partition into an informative estimate; this is why we treat $\rho$ as a diagnostic to report and use the functional quotient as a structural remedy, which works where the network has redundant linear operators to collapse (§ 5.3). Table 1 in the main text summarises the corrections.

Held-out versus in-sample estimation

All reported estimates use held-out points. To check how much this choice matters, we recompute the routing estimator in-sample, on all points of Composite ($N = 10{,}000$) and WBC ($N = 569$), for the same trained networks (Table 5). At the deepest layer the plug-in estimates differ by at most 0.020 bits on Composite and 0.031 bits on WBC (over all layers: 0.023 and 0.075 bits), even though $\rho$ is lower at the larger $N$. The estimate is therefore not an artefact of the particular points used to build the partition. For the two saturated WBC architectures this agreement is weaker evidence: once an estimator sits at its ceiling, its value no longer moves with $N$.

Dataset	Architecture	$\rho$ held-out	$\rho$ in-sample	$\hat{I}$ held-out	$\hat{I}$ in-sample	$\Delta\hat{I}$
Composite	[5, 5, 5]	0.030	0.007	1.949	1.939	−0.010
	[5, 5, 5, 5, 5]	0.054	0.013	2.261	2.245	−0.015
	[9, 9, 9]	0.064	0.016	2.188	2.174	−0.014
	[9, 9, 9, 9, 9]	0.147	0.043	2.550	2.538	−0.012
	[25, 25, 25]	0.304	0.104	2.535	2.515	−0.020
	[25, 25, 25, 25, 25]	0.411	0.163	2.620	2.609	−0.011
WBC	[5, 5, 5]	0.184	0.068	0.767	0.765	−0.002
	[5, 5, 5, 5, 5]	0.202	0.101	0.707	0.718	0.011
	[9, 9, 9]	0.360	0.187	0.740	0.727	−0.012
	[9, 9, 9, 9, 9]	0.561	0.316	0.904	0.885	−0.020
	[25, 25, 25]	0.900	0.741	0.949	0.919	−0.031
	[25, 25, 25, 25, 25]	0.961	0.893	0.949	0.943	−0.006

Table 5: Held-out versus in-sample routing estimates at the deepest layer (mean over seeds, bits). Held-out / in-sample $N$: Composite 2,000 / 10,000, WBC 114 / 569. $\Delta\hat{I}$: in-sample minus held-out.

Sensitivity of the baselines to their hyperparameters

The calibration in § 5.1 uses one setting per baseline. Table 6 repeats the correlation with the routing estimator for every setting we computed, both for the plug-in and the Miller–Madow-corrected binning and $k$-means estimates. On Composite and WBC the correlation is high for every setting except the finer Miller–Madow binnings (16 and 30 bins; lowest $r = 0.655$), and among the Miller–Madow binnings $K = 8$ is not the best-correlated bin count ($K = 4$ is). On MNIST FC the conclusion depends on the baseline’s hyperparameter: Miller–Madow binning with 16 and 30 bins per coordinate correlates negatively with the routing estimator ($r = -0.570$ and $-0.888$), because the correction for the many occupied bins dominates the estimate, while the plug-in variants of the same binnings stay positively correlated. This dependence on a free parameter is the problem the routing estimator avoids.

Baseline	Setting	Composite + WBC ($n = 60$)		MNIST FC ($n = 15$)
		plug-in	MM	plug-in	MM
Binning ($K$)	2	0.923	0.946	0.877	0.875
	4	0.966	0.973	0.758	0.728
	8*	0.981	0.945	0.807	0.690
	16	0.982	0.852	0.835	−0.570
	30	0.975	0.655	0.824	−0.888
$k$-means ($K$)	$K Y$*	0.970	0.970	0.720	0.720
	2$K Y$	0.974	0.974	0.866	0.866
	4$K Y$	0.974	0.974	0.868	0.868
	16	0.976	0.976	0.834	0.834
	64	0.972	0.971	0.861	0.861
	256	0.971	0.964	0.829	0.829
KSG ($k$)	3*	0.944		0.759
	5	0.948		0.762
	10	0.953		0.772

Table 6: Pearson $r$ between the routing estimator and each baseline setting, one point per network at the last hidden layer. MM: Miller–Madow-corrected. KSG is not a plug-in estimator and has a single column. *The setting used in Figure 2.

KSG on networks with collapsed representations

On five WBC networks the KSG baseline in Figure 2 reads zero or nearly zero while the routing estimator and the binning baseline report 0.2–0.9 bits (Table 7). In each of them a large group of held-out points, 22.8–78.1% of the 114, shares one identical representation at the last hidden layer. Every point in such a group passes through an earlier layer in which all units are inactive, so everything downstream of that layer is the same constant vector. In four of the five networks the group is almost pure in one class.

Architecture	Seed	Test acc.	$\hat{I}$	Binning	KSG raw	KSG	Tied points
[5, 5, 5, 5, 5]	101	98.2%	0.791	0.736	−0.469	0.000	41 (36%)
[5, 5, 5, 5, 5]	103	97.4%	0.758	0.795	0.148	0.148	26 (23%)
[5, 5, 5, 5, 5]	104	63.2%	0.392	0.217	−3.070	0.000	89 (78%)
[5, 5, 5, 5, 5]	105	97.4%	0.911	0.576	−1.743	0.000	64 (56%)
[5, 5, 5]	104	96.5%	0.779	0.720	−1.161	0.000	53 (46%)

The tie breaks the KSG estimator rather than the network. For a tied point the $k$-th nearest same-class neighbour lies at distance zero, so the neighbourhood the estimator counts contains the whole group instead of $k$ points. This violates the assumption behind the estimator and biases it downwards: the raw estimate is negative on four of the networks, and the reported value is clipped at zero, as for every KSG value in the paper. The routing estimator is unaffected, because a tied group simply contributes its label counts to the regions it occupies.

The effect on the calibration is small. Removing the five networks leaves the pooled KSG correlation of Figure 2 unchanged (0.892 with them, 0.893 without). In the WBC panel of Figure 3, four of them ([5, 5, 5, 5, 5], seeds 101, 103, 104, 105) are tied from layer 3 onwards and pull the mean KSG curve down from 0.837 to 0.619 bits; without them, KSG stays between 0.823 and 0.840 bits at every layer while $\hat{I}$ rises from 0.675 to 0.905 bits. The next-largest tied share among the WBC networks is 15.8% ([5, 5, 5], seed 105, KSG 0.436 bits); on Composite and MNIST FC it is at most 7.7% and 5.3%, with no visible effect on KSG.

Routing information without learned structure: label permutation

High routing information need not reflect learned class structure. As $\rho\to1$ the partition approaches an injective map of the samples, and $\hat{I}$ approaches $H(Y)$ for any labelling. To separate the two, we train the [5, 5, 5] and [25, 25, 25] architectures on Composite and WBC with the training labels randomly permuted, and estimate the routing information about the true labels on the held-out points (Table 8). The permuted networks learn neither the true labels nor the permuted ones: on Composite their accuracy is near chance (test accuracy on true labels 19.6% and 15.9%, accuracy on their own permuted training labels 20.8% and 22.5%), and on WBC they predict close to the majority class.

For the narrow Composite network the routing information drops from 1.949 to 1.411 bits, and the functional quotient discounts it further (1.828 to 1.096 bits). For the wide network it barely moves: 2.535 to 2.433 bits raw and 2.362 to 2.329 bits after the quotient, close to $H(Y) = 2.639$ bits, although the network predicts the true label of only 15.9% of the test points. A fine partition of Composite’s spatially coherent classes correlates with the true class whether or not the network has learned it. On WBC the same holds for [25, 25, 25], whose permuted network sits at $\rho> 0.9$. High routing information is therefore not by itself evidence of learned structure, and the functional quotient, while a valid lower bound by the data-processing inequality, removes this effect only at moderate capacity.

Sensitivity of the functional quotient to the visiting order

The greedy clustering of Algorithm 1 depends on the order in which regions are visited (App. § 14). For every clean network and layer we rerun it with 15 random orders in addition to the first-encounter order and record the spread (maximum minus minimum) of the quotient size and of the Miller–Madow quotient MI (Table 9). At $\epsilon= 0.3$ the median spread is 2 classes; the largest is 84 classes, or 9.4% of $\lvert\Omega_D\rvert$, and the quotient MI varies by at most 0.148 bits (0.100 bits at the deepest layer). The spread grows with $\epsilon$: over all four values, $\epsilon\in\{0.1, 0.3, 0.5, 1\}$, it reaches 0.474 bits. The order is thus a small but not negligible source of variation at large $\epsilon$, which is one more reason to stay in the moderate range.

Layers	$\epsilon$	Cells	Median	95th pct.	Max	Max / $R$	Max $\Delta\hat{I}_{\mathrm{func}}^{\mathrm{MM}}$
all layers	0.1	690	0	9.0	17	6.7%	0.053
all layers	0.3	690	2	38.5	84	9.4%	0.148
all layers	0.5	690	4	43.0	83	18.8%	0.348
all layers	1	690	1	19.0	47	14.3%	0.474
deepest layer	0.1	210	0	11.6	17	6.2%	0.053
deepest layer	0.3	210	11	45.5	60	6.5%	0.100
deepest layer	0.5	210	14	53.1	83	8.3%	0.348
deepest layer	1	210	5	28.0	47	10.5%	0.474

Implementation, Assets and Compute

Code. All code, the stored results and instructions to reproduce every figure and table are released under the MIT licence at https://github.com/Johanmkr/Intrinsic-IT-analysis-of-ReLU-nets, with a locked environment (uv.lock) and a Dockerfile. Networks are trained with PyTorch 2.11.0 (CPU build) [36]; the Composite generators (make_moons, make_circles), PCA, min-max scaling and $k$-means come from scikit-learn 1.8.0 [37].

Baselines. Binning quantises every unit of $T$ into $K$ uniform bins on $[-c,c]$, with $c$ the largest absolute activation in the layer, and uses the joint bin index as a discrete variable. $k$-means uses scikit-learn with 10 initialisations, seeded by the network seed. KSG is our own implementation of the mixed continuous–discrete estimator [38] with $k=3$; negative values are clipped at zero (App. § 15.§ 5). Binning and $k$-means are scored with the plug-in and the Miller–Madow estimators (App. § 15.§ 4).

Licences of existing assets. MNIST [25]: CC BY-SA 3.0. WBC [48]: CC BY 4.0 (UCI Machine Learning Repository). PyTorch, scikit-learn, NumPy and SciPy: BSD-style licences. The infomeasure package 0.6.2 [4] (AGPL-3.0) is used only in tests of our bias-correction implementations and is not part of the default installation.

Compute. Everything runs on CPU. On a 22-core x86-64 workstation with 20 parallel single-threaded workers, the full pipeline takes about 2 h of wall time, about 42 CPU-hours: training the 230 networks of both studies and the label-permutation control takes 71 min (about 24 CPU-hours, mostly the MNIST networks at about 3.5 min each), the routing estimates 39 min, the baselines 8 min and the diagnostics 6 min. Memory use is small. The full research project required additional CPU compute for preliminary experiments and earlier versions of the pipeline that are not reported.

Training Metrics

Test accuracy and test cross-entropy loss over training, on the held-out points, for every network in the paper: the six Composite and six WBC architectures of Study 1 (Figure 9 and Figure 10) and all MNIST networks of both studies (Figure 11). Curves are means over the five seeds, shaded bands $\pm1\sigma$.

Composite

WBC One network does not converge: seed 104 of [5, 5, 5, 5, 5] ends at 63% test accuracy, close to the majority-class rate of 62.7%, which produces the wide band in that panel.

MNIST