Paper deep dive
Analog-Friendly Predictive Coding without Activation Derivatives
Francesco Innocenti
Intelligence
Status: not_run | Model: - | Prompt: - | Confidence: 0%
Entities (0)
Relation Signals (0)
No relation signals yet.
Cypher Suggestions (0)
No Cypher suggestions yet.
Abstract
Abstract:Predictive coding (PC) is a local, energy-based alternative to backpropagation (BP) whose iterative inference dynamics make it attractive for implementation on analog hardware. However, standard nonlinear PC requires evaluating the derivative of the activation function during both inference and learning, which can be difficult to realise physically. Here, we introduce \textit{activation-matched Bregman PC}, replacing standard squared-error energies with Bregman divergences matched to the activation function. This formulation eliminates activation derivatives and, when combined with inference via mirror descent, yields local inference and learning rules requiring only weighted sums, local prediction errors, state integration, and the activation function. In digital experiments, Bregman PC performs comparably to standard PC and BP on classification and generative tasks, while preserving characteristic learning dynamics of PC and its convergence to BP under stable large-model parameterisations. These results provide a more analog-friendly formulation of nonlinear PC while retaining its key computational properties.
Tags
Links
- Source: https://arxiv.org/abs/2609.32350v1
- Canonical: https://arxiv.org/abs/2609.32350v1
Trouble viewing inline? Open PDF directly β
Full Text
66,542 characters extracted from source content.
Expand or collapse full text
Analog-Friendly Predictive Coding without Activation Derivatives Francesco Innocenti Affiliation: Brain Network Dynamics Unit Affiliation: University of Oxford, UK Email: francesco.innocenti@ndcn.ox.ac.uk Abstract Predictive coding (PC) is a local, energy-based alternative to backpropagation (BP) whose iterative inference dynamics make it attractive for implementation on analog hardware. However, standard nonlinear PC requires evaluating the derivative of the activation function during both inference and learning, which can be difficult to realise physically. Here, we introduce activation-matched Bregman PC, replacing standard squared-error energies with Bregman divergences matched to the activation function. This formulation eliminates activation derivatives and, when combined with inference via mirror descent, yields local inference and learning rules requiring only weighted sums, local prediction errors, state integration, and the activation function. In digital experiments, Bregman PC performs comparably to standard PC and BP on classification and generative tasks, while preserving characteristic learning dynamics of PC and its convergence to BP under stable large-model parameterisations. These results provide a more analog-friendly formulation of nonlinear PC while retaining its key computational properties. 1 Introduction As the energy demands of AI continue to grow (De Vries, 2023; Luccioni et al., 2024), there is increasing interest in alternatives to conventional digital computing and backpropagation-based training (Hooker, 2021). Analog hardware can exploit the intrinsic dynamics of physical systems to perform inference with high energy efficiency (Shainline, 2021; Wright et al., 2022; Aifer et al., 2025; Momeni et al., 2025; Wright et al., 2026; Wang and Achour, 2026; Lockwood et al., 2026), even when model parameters are trained off-chip. Ultimately, however, further gains may come from performing learning directly on the physical substrate itself, reducing analog-to-digital communication and moving closer to the online adaptation of biological systems. This motivates learning algorithms built from local operations that can be realised directly in hardware. Predictive coding (PC) is a promising but underexplored candidate for analog implementation. PC networks (PCNs) define local prediction errors between the activity of each layer and the prediction generated by the preceding layer, minimising an energy given by the sum of these errors (Rao and Ballard, 1999; Friston, 2005). During inference, neural activities evolve to reduce the energy. Once inference has converged, weights are updated using only locally available pre- and postsynaptic information. PC therefore replaces the explicit backward pass of backpropagation with local recurrent dynamics that can, in principle, be realised directly by a physical system. Recent work on scaling PC to larger models and more complex tasks has been encouraging. Several studies have improved the training of deep PCNs (Qi et al., 2025; Goemaere et al., 2025), including networks with more than 100 layers (Innocenti et al., 2025; Seely and Gould, 2026), and stable parameterisations have been derived for PCNs in large-width and depth limits (Ishikawa et al., 2024; Innocenti et al., 2026). More recently, Kerjan et al. (2026) showed that convolutional PCNs trained using equilibrium propagation can achieve BP-competitive performance on ImageNet classification. Despite these promising advances, it remains largely unknown how PC could be implemented on analog hardware. One important obstacle is that standard nonlinear PC formulations require explicit evaluation of activation-function derivatives during both inference and learning (see Figure 1). While an analog device may naturally realise the activation function itself, separately computing its derivative can require additional circuitry or operations that are difficult or inefficient to implement on the same substrate. Motivated by biological plausibility, previous work has explored removing these activation derivatives from PC, but found that doing so significantly degraded performance on simple classification tasks (Millidge et al., 2020). We therefore ask: Can nonlinear PC be formulated without explicit activation derivatives, making it better suited to analog hardware, while preserving its key computational properties? Here, we provide a positive answer to this question. Building on previous activation-matched Bregman constructions (Wang and Benning, 2020; Amid et al., 2022; Wang and Benning, 2023), we show that activation derivatives can be eliminated by replacing the squared Euclidean prediction errors of standard PC with Bregman divergences whose geometry is matched to the activation function. Combined with mirror descent of the neural states, we obtain local inference and learning updates that require only weighted interactions, local prediction errors, state integration, and evaluation of the activation itself (Figure 1). The resulting formulation recovers standard PC for linear networks while providing nonlinear PC dynamics that are more naturally suited to analog implementation. The rest of the paper is organised as follows. We first review standard PC and Bregman divergences (Β§2), highlighting the origin of activation derivatives in nonlinear PC. In Β§3, we introduce activation-matched Bregman-PC, deriving its energy, mirror-descent inference dynamics, and local learning rule. We then evaluate Bregman-PC on standard learning tasks and test whether it preserves key computational properties of PC (Β§4). We conclude by discussing the main limitations and future directions of this work (Β§5). Standard predictive codingBregman predictive coding β±β=12ββββ22=12ββββΟβ‘(β)β22 = 12 \| \|_2^2= 12 \| z -Ο( a ) \|_2^2 β±Bβ=DΞ¨β(β,Οβ‘(β)),βΞ¨β(β)=Οβ1β(β) _B=D_ ( z ,Ο( a ) ), β ( z )=Ο^-1( z ) Error measure dββdβΟ=ββ+Οβ²β(β+1)β(β+1)β€ββ+1 d z dΟ=- + [rgb]0.75,0.375,0Ο ( a +1) ( W +1) +1 dββdβΟ=ββ+β+(β+1)β€ββ+1,β=Οβ‘(β) d u dΟ=- u + a +( W +1) +1, z =Ο( u ) Inference dββdβt=Ξ·βΟβ²β(β)βββ(ββ1)β€ d W dt=Ξ·\, [rgb]0.75,0.375,0Ο ( a ) ( z -1) dββdβt=Ξ·βββ(ββ1)β€ d W dt=Ξ·\, ( z -1) Learning Figure 1: Standard predictive coding versus activation-matched Bregman predictive coding. Both formulations use the same nonlinear prediction ^β=Οβ‘(β) z =Ο( a ) (see Table 1 for notation), but differ in the error measure used at each hidden layer. Standard PC uses a squared Euclidean error (Eq. 3), resulting in explicit activation derivatives Οβ² [rgb]0.75,0.375,0Ο in both inference and learning (Eq. 4 & Eq. 6). Bregman PC instead uses an activation-matched Bregman divergence (Eq. 18) satisfying βΞ¨β(β)=Οβ1β(β)β ( z )=Ο^-1( z ). This removes the explicit activation derivative: inference evolves the dual state β u (Eq. 22), with β=Οβ‘(β) z =Ο( u ), while learning reduces to a local error-activity outer product (Eq. 24). 1.1 Summary of contributions β’ Activation-matched Bregman PC without explicit activation derivatives. We formulate nonlinear PC using activation-matched Bregman divergences as layerwise prediction errors and mirror descent in the corresponding dual coordinates. The resulting local inference and learning rules do not require explicit activation derivatives (Figure 1), while monotonically decreasing the Bregman-PC energy and recovering standard PC exactly in the linear case. β’ Empirical validation on simple learning tasks. Using digital simulations on classification and generative learning tasks, we show that Bregman PC trains nonlinear networks with performance comparable to standard PC and BP (Figure 2), demonstrating that removing explicit activation derivatives does not prevent effective learning in practice. β’ Preservation of established PC learning regimes. We show that Bregman PC retains learning regimes previously established for standard PC (Figure 5), including fast βsaddle-to-saddle dynamicsβ near degenerate saddles and convergence of PC gradients to BP under stable large-model parameterisations. 2 Background Table 1: Key notation. Symbol Name Description β z Activity / primal state State of layer β . The input and output activities are clamped to x and y, respectively, while hidden activities are inferred. β u Dual state Dual coordinate of the activity, β=Οβ1β(β) u =Ο^-1( z ), equivalently β=Οβ‘(β) z =Ο( u ). β a Preactivation Bottom-up linear input to layer β , β=ββββ1 a = W z -1. At feedforward initialisation, β=β a = u . ^β z Predicted activity Bottom-up prediction of the activity at layer β , ^β=Οβ‘(β) z =Ο( a ). ΟΟ Activation function Map Ο:βnβΟβ‘(βn)Ο:R^nβΟ(R^n) applied elementwise. Ο Scalar Bregman potential Strictly convex scalar function matched to the activation through Οβ²β(z)=Οβ1β(z)Ο (z)=Ο^-1(z). Ξ¨ Layerwise Bregman potential Separable extension of Ο to vector-valued activities, Ξ¨β‘()=βiΟβ‘(zi) ( z)= _iΟ(z_i). DΞ¨β(,)D_ ( p, q) Bregman divergence Non-negative discrepancy generated by Ξ¨ between p and q. 2.1 Predictive coding networks While PCNs can be defined on arbitrary graph topologies (Salvatori et al., 2022), here we consider feedforward networks or multi-layer perceptrons (MLPs) with L layers. We denote the state or activity of layer β by β z and the weights connecting layer β to layer β+1 +1 by β W . Without loss of generality, we omit biases and consider a single data sample. For supervised learning, the input and output activities are clamped to some data so that 0β z^0β x and Lβ z^Lβ y, while the hidden states β=1Lβ1\ z \_ =1^L-1 are free to evolve during inference, as explained below. Energy function. We write the bottom-up prediction of the activity of layer β as ^β=Οβ‘(β),β=ββββ1, z =Ο( a ), a = W z -1, (1) where β a is the corresponding preactivation, and Οβ‘(β )Ο(Β·) is an elementwise activation function such as tanh . Standard PC defines the layerwise prediction error as11 1 In some works, the activation function is applied before the weighted sum, ββΟβ(ββ1) W Ο ( z -1 ). This removes the activation derivative from the weight, but not the activity, update. ββββ^β=ββΟβ‘(ββββ1), z - z = z -Ο ( W z -1 ), (2) and minimises the sum of squared prediction errors β±β‘(,)=12βββ=1Lβββ22, ( z, ΞΈ)= 12 _ =1^L \| \|_2^2, (3) where we collect all the states in ββ=0L z \ z \_ =0^L, and the weights in ββ=1L ΞΈ \ W \_ =1^L. Eq. 3 is an energy function in two separate senses: (i) it minimises the so-called variational free energy (equal to the negative evidence lower bound) under certain assumptions (Buckley et al., 2017; Bogacz, 2017); and (i) both the inference and learning dynamics of PC can be derived as gradient flows on the energy, as we review below. Note that prediction errors are measured using the squared Euclidean distance. Inference dynamics. Before updating the weights, the free, hidden activities ββ1,β¦,Lβ1 β\1,β¦,L-1\ are inferred via gradient descent (GD) on the energy (Eq. 3). In continuous time, the inference dynamics are then given by dββdβΟ=ββββ±=ββ+Οβ²β(β+1)β(β+1)β€ββ+1, d z dΟ=- _ z F=- + [rgb]0.75,0.375,0Ο ( a +1 ) ( W +1 ) +1, (4) where Ο represents inference time, Οβ²β(β ) [rgb]0.75,0.375,0Ο (Β·) is the elementwise derivative of the activation function, and β denotes elementwise multiplication. We highlight the activation derivative since it is the object that we aim to avoid explicitly computing. The inference dynamics are local in the sense that updating a hidden layer requires information only from that layer and its immediate neighbours. We typically initialise the activities with a feedforward pass and run the dynamics until an approximate equilibrium ββargβ‘minβ‘β±β‘(,), z β _ zF( z, ΞΈ), (5) which should be numerically close to the unique global minimum for linear networks (Innocenti et al., 2025), with Ο=Ο= I, or any critical point for nonlinear networks. Learning dynamics. At the converged values of the inferred activities (Eq. 5), the weights are updated by stochastic GD on the same energy (Eq. 5). Again in continuous time, the learning dynamics are dββdβt=βΞ·ββββ±β(β)=Ξ·βΟβ²β(β)βββ(ββ1)β€, d W dt=-Ξ·\, _ W F( z )=Ξ·\, [rgb]0.75,0.375,0Ο ( a ) ( z -1 ) , (6) where Ξ· is the learning rate. Note that we do not backpropagate through the inner, inference dynamics (Eq. 5). Β§A explains why this choice is justified when an (approximate) equilibrium is reached. Similar to the states, the weight updates are local in that they require only the presynaptic activity and quantities available at the postsynaptic layer. On analog implementation. The standard formulation of PC reviewed above highlights a potential difficulty for physical implementations of nonlinear PCNs. In particular, although a nonlinear activation function ΟΟ may arise naturally from the input-output response of an analog device, the inference and learning dynamics of Eqs. 4 & 6 always explicitly require its derivative Οβ² [rgb]0.75,0.375,0Ο . In Β§3, we show that this dependence can be removed by replacing the squared Euclidean prediction errors with Bregman divergences matched to the particular activation function. 2.2 Bregman divergences Definition. A Bregman divergence generalises the squared Euclidean distance by using the geometry of a convex function (Bregman, 1967). Using standard notation, let Ξ¨:ββ :D be a differentiable and strictly convex function on a convex domain D. The Bregman divergence generated by Ξ¨ is DΞ¨β(,) D_ ( p, q) βΞ¨β‘()βΞ¨β‘()ββΞ¨β()β€β(β). ( p)- ( q)- _ q ( q) ( p- q). (7) The last two terms in Eq. 7 are the first-order Taylor approximation of Ξ¨β‘() ( p) around q. The Bregman divergence DΞ¨β(,)D_ ( p, q) therefore measures the gap between the value of the convex function at p and the tangent-plane approximation at q (see Figure 2c for an example). Properties. Since the function Ξ¨ is strictly convex by assumption, the divergence cannot be negative DΞ¨β(,)β₯0D_ ( p, q)β₯ 0, and is exactly zero when = p= q. Ξ¨ is also referred to as the Bregman potential, since it is the function from which the divergence is derived (for key notation, see Table 1). Finally, note that Bregman divergences are not, in general, symmetric and therefore need not define a metric. Example: squared Euclidean distance. Consider the scalar quadratic potential Οβ‘(p)=12βp2Ο(p)= 12p^2. Since Οβ²β(q)=qΟ (q)=q, substituting into the Bregman divergenceβs definition (Eq. 7) gives DΟβ(p,q)=12βp2β12βq2βqβ‘(pβq)=12β(pβq)2. D_Ο(p,q)= 12p^2- 12q^2-q(p-q)= 12(p-q)^2. (8) For the squared vector norm Ξ¨β‘()=12ββ22 ( p)= 12\| p\|_2^2, this immediately generalises to DΞ¨β(,)=12ββββ22. D_ ( p, q)= 12\| p- q\|_2^2. (9) Thus, the quadratic prediction errors of the standard PC energy (Eq. 3) can be seen as special Bregman divergences with a squared Euclidean norm as generating or potential function. For other examples of Bregman divergences (e.g. the KullbackβLeibler divergence), see Β§B. 3 Predictive Coding with Bregman Divergences Figure 2: Constructing an activation-matched Bregman divergence for tanh . (a) The activation function Οβ‘(a)=tanhβ‘(a)Ο(a)= (a). (b) Its inverse, Οβ1β(z)=arctanhβ‘(z)Ο^-1(z)=arctanh(z), which defines the derivative of the matching potential through Οβ²β(z)=Οβ1β(z)Ο (z)=Ο^-1(z) (Eq. 10). (c) Integrating the inverse activation gives the matching potential Οβ‘(z)=zβarctanhβ‘(z)+12βlogβ‘(1βz2)Ο(z)=zarctanh(z)+ 12 (1-z^2), up to an additive constant. Since arctanhβ‘(z)arctanh(z) is strictly increasing, Ο is strictly convex. The Bregman divergence DΟβ(z,z^)D_Ο(z, z) is the gap between Οβ‘(z)Ο(z) and the tangent-line approximation of Ο at the prediction z^=Οβ‘(a) z=Ο(a). See Figures 7-8 for divergences with other activations. In this section, we develop a generalisation of PC in which nonlinear prediction errors are measured using Bregman divergences matched to the activation function, instead of squared Euclidean norms. Closely related activation-matched losses have previously been used to avoid explicit activation derivatives in parameter updates (Wang and Benning, 2020; Amid et al., 2022; Wang and Benning, 2023). βLocoPropβ uses these losses to optimise layerwise weights towards targets obtained from backpropagation (Amid et al., 2022), while βlifted Bregman trainingβ introduces hidden activities as auxiliary (PC-like) optimisation variables in a closely related layerwise Bregman objective (Wang and Benning, 2023). Here, we instead formulate these variables as PC states and derive local mirror-descent inference dynamics before parameter learning. Throughout, we use tanh as an example activation (Figure 2a) and provide results for other suitable activations in the Appendix. Activation function assumptions. We assume that ΟΟ is differentiable with Οβ²β(a)>0Ο (a)>0, so that it is strictly increasing and invertible on its range, and the matched potential defined below (Eq. 10) is strictly convex. These assumptions include, for example, tanh , sigmoid, and softplus activations. Non-smooth or non-monotonic activations require a generalised treatment, which we discuss in Β§5. Activation-matched Bregman potential. We define the matched Bregman scalar potential Ο:Οβ‘(β)ββΟ:Ο(R) by requiring its derivative to equal the inverse of the activation function (see Figure 2a-b for an example), Οβ²β(z)βΟβ1β(z). Ο (z) Ο^-1(z). (10) Equivalently, Ο is any antiderivative of Οβ1Ο^-1, since additive constants do not affect the resulting Bregman divergence. Because ΟΟ is strictly increasing by assumption, its inverse Οβ1Ο^-1 is also strictly increasing. Therefore, Οβ²Ο is strictly increasing and Ο is strictly convex (Figure 2c), which is the key property needed to define a Bregman divergence, as reviewed in Β§2.2. Why matching removes the activation derivative. To see the effect of the activation-matching condition of Eq. 10, consider a scalar prediction z^=Οβ‘(a) z=Ο(a) with prediction error Ξ΅βzβΟβ‘(a) z-Ο(a). As clear from the standard PC inference and learning dynamics (Eqs. 4 & 6), differentiating the squared prediction error with respect to the preactivation introduces the activation derivative: βaβ(zβΟβ‘(a))2/2=βΞ΅βΟβ²β(a). β a (z-Ο(a) )^2/2=- \,Ο (a). (11) Thus, the prediction-error signal is modulated by Οβ²β(a)Ο (a), and can be strongly attenuated when the activation saturates (see Figure 3). For the activation-matched Bregman divergence DΟβ(z,Οβ(a))D_Ο(z,Ο(a)), instead, Figure 3: Activation matching prevents attenuation of prediction errors. For tanh , the standard squared error gives βaβ±=βΞ΅βΟβ²β(a) _aF=- Ο (a), which is attenuated near saturation. The matched Bregman divergence contributes 1/Οβ²β(a)1/Ο (a), cancelling this factor and yielding βaDΟ=βΞ΅ _aD_Ο=- . Analogous results for other activations are shown in Figures 9-10. βaβDΟβ(z,Οβ‘(a))=βΟβ²β(Οβ‘(a))βΞ΅βΟβ²β(a). β aD_Ο (z,Ο(a) )=-Ο (Ο(a) ) \,Ο (a). (12) By the matching condition Οβ²=Οβ1Ο =Ο^-1 (Eq. 10) and the inverse-function theorem, Οβ²β(Οβ‘(a))=(Οβ1)β²β(Οβ‘(a))=1Οβ²β(a). Ο (Ο(a) )=(Ο^-1) (Ο(a) )= 1Ο (a). (13) Hence, the factor Οβ²β(a)Ο (a) introduced by differentiating the nonlinear prediction is exactly cancelled by the reciprocal factor arising from the curvature of the matched potential: βaβDΟβ(z,Οβ‘(a))=βΞ΅. β aD_Ο (z,Ο(a) )=- . (14) Activation matching therefore encodes the nonlinearity into the geometry of the error measure itself, so that the preactivation derivative depends directly on the error without requiring explicit evaluation of Οβ²β(a)Ο (a). Vector-valued activities. For multidimensional activities ββn z ^n, we apply the same construction elementwise and define the separable potential Ξ¨:Οβ(β)nββ,Ξ¨β‘()ββiΟβ‘(zi). :Ο(R)^n , ( z) _iΟ(z_i). (15) Its gradient is therefore the elementwise inverse activation, βΞ¨β()=(Οβ1β(z1)Οβ1β(zn))=Οβ1β(). _ z ( z)= pmatrixΟ^-1(z_1)\\ \\ Ο^-1(z_n) pmatrix=Ο^-1( z). (16) Hence, the corresponding activation-matched divergence satisfies βDΞ¨β(,Οβ‘())=β _ aD_ ( z,Ο( a) )=- . Bregman-PC energy. We now use the activation-matched divergence to construct a nonlinear PC energy whose inference and learning dynamics do not explicitly contain Οβ²Ο . In particular, we replace the layerwise squared prediction error of standard PC (Eq. 3) with β±BββDΞ¨β(β,Οβ‘(β)),βΞ¨β()=Οβ1β(). _B D_ ( z ,Ο( a ) ), _ z ( z)=Ο^-1( z). (17) For a network with a conventional linear readout L=LβLβ1 a^L= W^L z^L-1, we allow an arbitrary output loss ββ‘(L,)L( a^L, y) and define the total Bregman-PC energy as β±B=ββ=1Lβ1DΞ¨β(β,Οβ‘(β))+ββ‘(L,). _B= _ =1^L-1D_ ( z ,Ο( a ) )+L ( a^L, y ). (18) Each hidden term is non-negative and exactly zero at the feedforward pass, as in standard PC. This formulation accommodates, for example, a linear readout trained with either mean-squared error (MSE) or softmax cross-entropy.22 2 A Bregman divergence can also be used at the output whenever the output activation and target domain admit an appropriate matching potential. For convenience, we define the output error as the negative output-loss gradient, LβββLβ ^L - _ a^LL, so that for example, L=βL ^L= y- a^L for MSE. The Bregman-PC energy (Eq. 18) itself can be algebraically more involved than the quadratic energy of standard PC (Eq. 3), since the matched potential is obtained by integrating the inverse activation, Οβ‘(z)=β«zΟβ1β(s)βsΟ(z)= ^zΟ^-1(s)ds. For tanh , for example, Οβ‘(z)=zβarctanhβ‘(z)+12βlogβ‘(1βz2)Ο(z)=zarctanh(z)+ 12 (1-z^2), up to an additive constant. Crucially, however, neither this antiderivative nor the inverse activation needs to be explicitly evaluated in the local inference and learning dynamics derived below. Note that our formulation also remains compatible with equilibrium propagation (Scellier and Bengio, 2017), by introducing a nudging parameter Ξ² scaling the output loss βL. Primal and dual neural states. The activation-matched Bregman potential Ξ¨ (Eq. 16) also induces a natural change of coordinates for the neural activities. We refer to β z as the primal state and define its corresponding dual state as β u ββΞ¨β(β)=Οβ1β(β),β=Οβ‘(β) _ z ( z )=Ο^-1 ( z ), z =Ο ( u ) (19) Thus, β u and β z describe the same neural state in two coordinate systems related by the geometry of Ξ¨ . Under forward-pass initialisation, β=Οβ‘(β) z =Ο( a ), and therefore β=β u = a . In other words, the dual state is equal to the feedforward preactivation at initialisation. Mirror-descent inference. We now use these dual coordinates to formulate inference as mirror descent. In particular, rather than taking an ordinary gradient step directly in the primal activity β z , we evaluate the energy gradient with respect to β z but evolves its dual coordinate β=βΞ¨β(β) u = _ z ( z ). In continuous time, dββdβΟ=dβΟββΞ¨β(β)βββββ±B, d u dΟ= ddΟ _ z ( z ) - _ z F_B, (20) where recall that Ο denotes inference time. Since β=Οβ‘(β) z =Ο( u ), the same dynamics can be expressed in primal coordinates as dββdβΟ=β[β2Ξ¨β(β)]β1ββββ±B=βdiagβ‘(Οβ²β(uiβ))ββββ±B, d z dΟ=- [ _ z^2 ( z ) ]^-1 _ z F_B=-diag (Ο (u_i ) ) _ z F_B, (21) where the second equality follows from the activation-matching relation =Οβ‘() z=Ο( u). Mirror inference is therefore equivalent to a preconditioned gradient flow in the primal activities, where the preconditioner is the inverse local curvature of the matched potential Ξ¨ . Thus, mirror flow changes the trajectory and local speed of inference relative to ordinary Euclidean gradient flow, as illustrated in Figure 4. Figure 4: Mirror inference in primal and dual coordinates. Inference trajectories for a two-hidden-state Bregman PCN with Ο=tanhΟ= , initialised from the same state and, in this example, converging to the same equilibrium. (a) In primal coordinates =(z1,z2) z=(z^1,z^2), Euclidean inference follows dβ/dβΟ=βββ±Bd z/dΟ=- _ zF_B, whereas mirror inference follows the preconditioned flow dβ/dβΟ=β[β2Ξ¨β()]β1βββ±Bd z/dΟ=-[β^2 ( z)]^-1 _ zF_B, leading to a different trajectory through activity space. (b) In dual coordinates =βΞ¨β()=Οβ1β() u=β ( z)=Ο^-1( z), the same mirror dynamics take the simpler form dβ/dβΟ=βββ±Bd u/dΟ=- _ zF_B. Grey contours show the Bregman-PC energy in the corresponding coordinates. For other activations, see Figures 11-12. For the Bregman-PC energy (Eq. 18), the dual-coordinate dynamics (Eq. 20) take the particularly simple form: dββdβΟ=ββ+β+(β+1)β€ββ+1. d u dΟ=- u + a + ( W +1 ) +1. (22) See Appendix C.1 for the derivation. The update depends only on the current dual state, its linear bottom-up prediction, and a weighted prediction error from the next layer. Importantly, an analog implementation can evolve β u directly and obtain the primal activity from the activation, β=Οβ‘(β) z =Ο( u ). Inference therefore requires neither an explicit activation derivative Οβ²Ο nor an explicit evaluation of the inverse activation Οβ1Ο^-1. A candidate analog circuit implementing these dynamics is described in Appendix D. Energy descent. Although mirror descent follows a different trajectory from Euclidean gradient flow (Figure 4), it still monotonically decreases the Bregman-PC energy. Since Ξ¨ is convex, its inverse Hessian is positive semidefinite wherever it is defined [β2βΞ¨β(β)]β1βͺ°0 [ _ z^2 ( z ) ]^-1 0, and hence dββ±BdβΟ dF_BdΟ =ββ=1Lβ1(βββ±B)β€dββdβΟ=βββ=1Lβ1(βββ±B)β€[β2Ξ¨(β)]β1(βββ±B)β€0. = _ =1^L-1 ( _ z F_B ) d z dΟ=- _ =1^L-1 ( _ z F_B ) [ _ z^2 ( z ) ]^-1 ( _ z F_B )β€ 0. (23) Thus, while mirror inference changes the path and speed of inference, it preserves descent of the same Bregman-PC energy. Weight updates. The matched PC-Bregman energy (Eq. 18) also removes the activation derivative from the learning dynamics (derivation in Β§C.2): dββdβt=βΞ·ββββ±B=Ξ·βββ(ββ1)β€, d W dt=-Ξ·\, _ W F_B=Ξ·\, ( z -1 ) , (24) which has the same local outer-product form between the postsynaptic prediction error and the presynaptic activity as standard linear PC. Recovery of standard linear PC. Our proposed Bregman-PC formulation recovers standard linear PC exactly. Specifically, if Οβ‘()=Ο( a)= a, then Οβ1β()=Ο^-1( a)= a and the matching potential is the squared norm Ξ¨β‘()=12ββ22 ( z)= 12\| z\|_2^2. As shown in Β§2.2, the layerwise matched Bregman divergence then reduces to the standard squared prediction error: DΞ¨β(β,ββββ1)=12ββββββββ1β22. D_ ( z , W z -1 )= 12 \| z - W z -1 \|_2^2. (25) Moreover, the Hessian is the identity, the primal and dual coordinates coincide β=β u = z , and mirror inference reduces to ordinary Euclidean gradient flow (see Figures 8, 10 & 12). Bregman PC therefore reduces exactly to standard PC in the linear case, while providing an activation-matched nonlinear generalisation in which neither inference nor learning requires explicit derivatives of the activation function. This recovery is useful because most PC theory has been developed for linear networks yet shown to predict and control behaviour of nonlinear ones (Ishikawa et al., 2024; Innocenti et al., 2024; Innocenti et al., 2026; Seely and Gould, 2026), as we further demonstrate experimentally in the next section. 4 Experiments Table 2: Classification and generation performance of Bregman PC vs BP and standard PC. Validation classification accuracy and image-reconstruction MSE on MNIST and Fashion-MNIST (mean Β± SEM over 3 seeds). See Figures 13-14 for loss curves and examples of generated images, and Β§E for experimental details. Classification accuracy (%) β Reconstruction MSE [0,1][0,1] β Algorithm MNIST Fashion-MNIST MNIST Fashion-MNIST BP 98.05Β±0.1098.05Β± 0.10 88.92Β±0.1188.92Β± 0.11 0.05310Β±0.000030.05310Β± 0.00003 0.05263Β±0.000040.05263Β± 0.00004 Standard PC 98.08Β±0.0898.08Β± 0.08 88.87Β±0.1088.87Β± 0.10 0.05320Β±0.000030.05320Β± 0.00003 0.05275Β±0.000060.05275Β± 0.00006 Bregman PC 97.57Β±0.0897.57Β± 0.08 87.30Β±0.2487.30Β± 0.24 0.05334Β±0.000070.05334Β± 0.00007 0.05326Β±0.000060.05326Β± 0.00006 We now show that Bregman PC remains practically competitive with standard PC while still preserving its key computational properties. We focus on two main questions: (i) whether Bregman PC can train nonlinear networks with performance comparable to standard PC; and (i) whether it retains previously established learning regimes for PC. Performance on simple tasks. We first compare Bregman PC against standard PC and BP on standard classification and generative learning tasks using MNIST and Fashion-MNIST, following the experimental setting of Pinchetti et al. (2024). Across these tasks, Bregman PC reaches performance comparable to both standard PC and BP (see Table 2), indicating that replacing squared prediction errors with activation-matched Bregman divergences, and performing inference in the dual coordinates, does not significantly hurt performance. Preservation of established PC learning regimes. We next ask whether Bregman PC preserves some learning regimes previously established for standard PC, ranging from dynamics that differ qualitatively from BP to a regime in which PC approaches BP. First, we reproduced an experiment from Innocenti et al. (2024), training nonlinear MLPs on MNIST from a small, near-origin initialisation. In this regime, BP remains trapped for an extended period near the highly degenerate saddle at the origin, whereas standard PC rapidly escapes, showing fast βsaddle-to-saddleβ dynamics. Bregman PC reproduces this behaviour, escaping on a similar timescale and following a nearly indistinguishable loss trajectory from standard PC (Figure 5a). Figure 5: Bregman PC preserves established PC learning regimes. (a) βSaddle-to-saddle dynamicsβ for 5-layer tanh MLPs trained on MNIST with a small, near-origin initialisation. BP remains on a prolonged high-loss plateau, whereas standard PC and Bregman PC rapidly escape, with nearly overlapping loss trajectories. (b) Cosine similarity between Bregman-PC and BP parameter gradients after 100100 training steps for nonlinear residual networks trained on CIFAR-10 using the width- and depth-stable parameterisation of Innocenti et al. (2026). Similarities are averaged over three random seeds and approach 11 as width becomes large relative to depth. For more details, see Β§E. At the other extreme, we test whether Bregman PC also preserves a recently established regime where the PC gradients converge to those computed by BP under stable large-model parameterisations (Ishikawa et al., 2024; Innocenti et al., 2026). Following Innocenti et al. (2026), we trained nonlinear residual networks using a width- and depth-stable parameterisation on CIFAR-10, varying both network width N and depth L. Consistent with previous results, we find that the Bregman-PC gradients align with the BP gradients as the model width becomes much larger than the depth (Figure 5b). Together, these results show that Bregman PC preserves both the distinctive learning dynamics of standard PC, in the form of fast saddle-to-saddle dynamics, while recovering its asymptotic convergence to BP for large, stably parameterised models. 5 Discussion Summary. We introduced an activation-matched Bregman formulation of PC in which inference and learning avoid explicit evaluation of activation-function derivatives (Figure 1). By choosing the Bregman potential such that βΞ¨=Οβ1β =Ο^-1, the derivative of the nonlinearity is absorbed into the geometry of the prediction error. Combined with mirror-descent inference in the corresponding dual coordinates, this yields local derivative-free dynamics while recovering standard PC exactly in the linear case. Our experiments show that this modification preserves both competitive learning performance and key computational regimes of standard PC. 5.1 Limitations and future directions Hardware implementation and validation. Our main motivation is analog implementation. Our formulation should therefore be understood as making nonlinear PC more hardware-compatible, rather than providing a complete analog implementation. The inference and learning dynamics of Bregman PC require only weighted interactions, local prediction errors, state integration, and evaluation of the activation function itself. Appendix D illustrates how the inference dynamics could be mapped onto a simple analog circuit. Whether the proposed formulation is sufficient for an effective analog implementation of PC remains to be established and should ultimately be validated in hardware, including under device noise and mismatch. We also note that analog friendliness is related to, but distinct from, biological plausibility: while both often favour local and distributed computation, operations that are convenient to implement physically need not be biologically realistic, and vice versa. Activation-function and empirical limitations. A limitation of our construction is the assumption that the activation is differentiable and strictly increasing. This includes common saturating nonlinearities such as tanh and sigmoid, but excludes standard ReLU and globally non-monotonic activations such as GELU and SiLU. More general Bregman constructions based on proximal maps can accommodate non-smooth and non-invertible monotone activations such as ReLU (Wang and Benning, 2023). Extending the mirror-descent PC dynamics developed here to this setting is therefore a natural direction for future work, while globally non-monotonic activations will require a different treatment. Our empirical evaluation is also limited to relatively small feedforward networks and standard benchmark tasks, and testing whether the same behaviour persists in larger and more complex architectures will be important. AI use statement In this work, we used generative AI tools mainly for writing code to implement all the experiments and for editing the paper to improve clarity. Specifically, we used AI tools for formulating certain mathematical claims (the descent result of Eq. 23), designing and providing feedback on experiments, and implementing methods. We have not used generative AI tools for developing theoretical models or conceptual frameworks, providing critical ingredients for proving mathematical claims, proposing or refining hypotheses. Creating synthetic data sets, assisting in the writing of proofs or translation, cleaning or reformatting datasets, and supporting qualitative and thematic data analysis are not applicable to this work. Additionally, we used generative AI tools for suggesting experimental hyperparameters, brainstorming and editing. We have reviewed all AI-assisted work. In particular, LLM-generated code was verified and tested for correctness. We take responsibility for the final content of this work, including text, claims or artifacts produced with the aid of generative AI. Reproducibility statement Experimental details are provided in Β§E, and code to reproduce all the results is available at https://github.com/thebuckleylab/jpc/tree/main/experiments/bregman_pc. Acknowledgments This work was supported by funding from the Wellcome Trust grant 313955/Z/24/Z. I would like to thank Marc Gong Bacvanski, Bryce Primavera, Benjamin Scellier, and Elene Lominadze for discussions that motivated this work. Some of the simulations were supported by Google Cloud Credits awarded by the Google TPU Research Cloud Programme (2026). References Aifer et al. (2025) M. Aifer, Z. Belateche, S. Bramhavar, K. Y. Camsari, P. J. Coles, G. Crooks, D. J. Durian, A. J. Liu, A. Marchenkova, A. J. Martinez, et al. Solving the compute crisis with physics-based asics. arXiv preprint arXiv:2507.10463. Cited by: Β§1. Amid et al. (2022) E. Amid, R. Anil, and M. Warmuth Locoprop: enhancing backprop via local loss optimization. In International conference on artificial intelligence and statistics, p. 9626β9642. Cited by: Β§1, Β§3. Bai et al. (2019) S. Bai, J. Z. Kolter, and V. Koltun Deep equilibrium models. Advances in neural information processing systems 32. Cited by: Appendix A. Bogacz (2017) R. Bogacz A tutorial on the free-energy framework for modelling perception and learning. Journal of mathematical psychology 76, p. 198β211. Cited by: Β§2.1. Bregman (1967) L. M. Bregman The relaxation method of finding the common point of convex sets and its application to the solution of problems in convex programming. USSR computational mathematics and mathematical physics 7 (3), p. 200β217. Cited by: Β§2.2. Buckley et al. (2017) C. L. Buckley, C. S. Kim, S. McGregor, and A. K. Seth The free energy principle for action and perception: a mathematical review. Journal of mathematical psychology 81, p. 55β79. Cited by: Β§2.1. De Vries (2023) A. De Vries The growing energy footprint of artificial intelligence. Joule 7 (10), p. 2191β2194. Cited by: Β§1. Friston (2005) K. Friston A theory of cortical responses. Philosophical transactions of the Royal Society B: Biological sciences 360 (1456), p. 815β836. Cited by: Β§1. Goemaere et al. (2025) C. Goemaere, G. Oliviers, R. Bogacz, and T. Demeester EPC: fast and deep predictive coding in digital simulation. arXiv preprint arXiv:2505.20137. Cited by: Β§1. Hooker (2021) S. Hooker The hardware lottery. Communications of the ACM 64 (12), p. 58β65. Cited by: Β§1. Hoover et al. (2023) B. Hoover, Y. Liang, B. Pham, R. Panda, H. Strobelt, D. H. Chau, M. Zaki, and D. Krotov Energy transformer. Advances in neural information processing systems 36, p. 27532β27559. Cited by: Appendix A. Innocenti et al. (2024) F. Innocenti, E. M. Achour, R. Singh, and C. L. Buckley Only strict saddles in the energy landscape of predictive coding networks?. Advances in Neural Information Processing Systems 37, p. 53649β53683. Cited by: Appendix E, Β§3, Β§4. Innocenti et al. (2026) F. Innocenti, E. M. Achour, and R. Bogacz On the infinite width and depth limits of predictive coding networks. arXiv preprint arXiv:2602.07697. Cited by: Appendix E, Β§1, Β§3, Figure 5, Β§4. Innocenti et al. (2025) F. Innocenti, E. M. Achour, and C. L. Buckley ΞΌ : Scaling predictive coding to 100+ layer networks. arXiv preprint arXiv:2505.13124. Cited by: Β§1, Β§2.1. Ishikawa et al. (2024) S. Ishikawa, R. Yokota, and R. Karakida Local loss optimization in the infinite width: stable parameterization of predictive coding networks and target propagation. arXiv preprint arXiv:2411.02001. Cited by: Β§1, Β§3, Β§4. Kerjan et al. (2026) T. Kerjan, R. HΓΈier, and B. Scellier Training a predictive coding network on imagenet using equilibrium propagation. arXiv preprint arXiv:2606.03584. Cited by: Β§1. Kidger (2022) P. Kidger On neural differential equations. arXiv preprint arXiv:2202.02435. Cited by: Appendix A. Lockwood et al. (2026) O. Lockwood, J. BΓ©janin, J. Bus, C. Chamberland, P. Huembeli, F. SchΓ€fer, and G. Verdon A blueprint for equilibrium-based differentiable continuous-variable thermodynamic computing. arXiv preprint arXiv:2607.16183. Cited by: Β§1. Luccioni et al. (2024) S. Luccioni, Y. Jernite, and E. Strubell Power hungry processing: watts driving the cost of ai deployment?. In Proceedings of the 2024 ACM conference on fairness, accountability, and transparency, p. 85β99. Cited by: Β§1. Millidge et al. (2020) B. Millidge, A. Tschantz, A. Seth, and C. L. Buckley Relaxing the constraints on predictive coding models. arXiv preprint arXiv:2010.01047. Cited by: Β§1. Momeni et al. (2025) A. Momeni, B. Rahmani, B. Scellier, L. G. Wright, P. L. McMahon, C. C. Wanjura, Y. Li, A. Skalli, N. G. Berloff, T. Onodera, et al. Training of physical neural networks. Nature 645 (8079), p. 53β61. Cited by: Β§1. Pinchetti et al. (2024) L. Pinchetti, C. Qi, O. Lokshyn, G. Olivers, C. Emde, M. Tang, A. MβCharrak, S. Frieder, B. Menzat, R. Bogacz, et al. Benchmarking predictive coding networksβmade simple. arXiv preprint arXiv:2407.01163. Cited by: Β§4. Qi et al. (2025) C. Qi, M. Forasassi, T. Lukasiewicz, and T. Salvatori Towards the training of deeper predictive coding neural networks. arXiv preprint arXiv:2506.23800. Cited by: Β§1. Rao and Ballard (1999) R. P. Rao and D. H. Ballard Predictive coding in the visual cortex: a functional interpretation of some extra-classical receptive-field effects. Nature neuroscience 2 (1), p. 79β87. Cited by: Β§1. Salvatori et al. (2022) T. Salvatori, L. Pinchetti, B. Millidge, Y. Song, T. Bao, R. Bogacz, and T. Lukasiewicz Learning on arbitrary graph topologies via predictive coding. Advances in neural information processing systems 35, p. 38232β38244. Cited by: Β§2.1. Scellier and Bengio (2017) B. Scellier and Y. Bengio Equilibrium propagation: bridging the gap between energy-based models and backpropagation. Frontiers in computational neuroscience 11, p. 24. Cited by: Appendix A, Β§3. Seely and Gould (2026) J. Seely and J. Gould Augmented lagrangian predictive coding. arXiv preprint arXiv:2605.31022. Cited by: Β§1, Β§3. Shainline (2021) J. M. Shainline Optoelectronic intelligence. Applied Physics Letters 118 (16). Cited by: Β§1. Wang and Benning (2020) X. Wang and M. Benning Generalised perceptron learning. arXiv preprint arXiv:2012.03642. Cited by: Β§1, Β§3. Wang and Benning (2023) X. Wang and M. Benning Lifted bregman training of neural networks. Journal of Machine Learning Research 24 (232), p. 1β51. Cited by: Β§1, Β§3, Β§5.1. Wang and Achour (2026) Y. Wang and S. Achour Generative models on analog hardware with dynamics. arXiv preprint arXiv:2606.27294. Cited by: Β§1. Wright et al. (2022) L. G. Wright, T. Onodera, M. M. Stein, T. Wang, D. T. Schachter, Z. Hu, and P. L. McMahon Deep physical neural networks trained with backpropagation. Nature 601 (7894), p. 549β555. Cited by: Β§1. Wright et al. (2026) L. G. Wright, T. Wang, T. Onodera, and P. L. McMahon Physical foundation models: fixed hardware implementations of large-scale neural networks. arXiv preprint arXiv:2604.27911. Cited by: Β§1. Appendix Appendix A On predictive coding and βimplicit gradientsβ In our review of PCNs in Β§2.1, we pointed out that the weight update of PC (Eq. 6) does not backpropagate through the inference updates of Eq. 4, effectively treating the activities as constant. This can be justified if (i) the two levels of optimisation minimise the same objective (as in PC), and (i) we converge to a critical point of the inner problem, since ββ±β‘(ββ(),)β ( z ( ΞΈ), ΞΈ)β ΞΈ =ββ±ββexplicit+ββ±βββ()0ββββ()ββimplicit, = 0 β h ( ΞΈ) β ΞΈ_explicit+ 0 β z ( ΞΈ) β z ( ΞΈ)β ΞΈ_implicit, (26) where we define ββ()=argβminβ‘β±β(,) z ( ΞΈ)= _ zF( z, ΞΈ). We emphasise that the activity solution is a function of the weights by writing β‘() z( ΞΈ). Eq. 26 shows the total weight derivative of the energy at a critical point of the activities, where ββ±/ββ= /β z = 0. By the product rule, this decomposes into two terms: (i) the direct or explicit dependence of the energy on the weights (treating the activities as constant), and (i) the indirect or implicit effect of the weights on the energy through the activity equilibrium. Because (i) inference and learning minimise the same objective (i.e. the energy), and because (i) we assume convergence on the inner problem, the implicit gradient vanishes (Eq. 26). We therefore do not need to backpropagate through the PC inference updates and can just take the direct (local) weight gradient. This contrasts with related models such as deep equilibrium models (Bai et al., 2019), neural differential equations (Kidger, 2022) and modern Hopfield networks (Hoover et al., 2023), which use either backpropagation through time or other techniques to estimate the implicit gradient. We also note that this property is incompatible with equilibrium propagation (Scellier and Bengio, 2017), which introduces a cost or loss function in addition to an energy. Appendix B Other examples of Bregman divergences Kullback-Leibler (KL) divergence. Consider two discrete probability distributions p and q, and let the generating function be the negative (Shannon) entropy Ξ¨β‘()=βipiβlogβ‘pi. ( p)= _ip_i p_i. (27) Its gradient has components [βΞ¨β()]i=1+logβ‘qi[β ( q)]_i=1+ q_i. The corresponding Bregman divergence (Eq. 7) then is DΞ¨β(,) D_ ( p, q) =βipiβlogβ‘piββiqiβlogβ‘qi = _ip_i p_i- _iq_i q_i ββi(1+logqi)(piβqi) - _i(1+ q_i)(p_i-q_i) (28) =βipilogpiqi=DKL(β₯), = _ip_i p_iq_i=D_KL( p\| q), (29) where in the second equality we used βipi=βiqi=1 _ip_i= _iq_i=1. Thus, the KL divergence between probability distributions is a Bregman divergence generated by the negative entropy. Bernoulli KL divergence. For a scalar pβ(0,1)pβ(0,1), consider the negative binary entropy Οβ‘(p)=pβlogβ‘p+(1βp)βlogβ‘(1βp). Ο(p)=p p+(1-p) (1-p). (30) Notably, its derivative is the inverse sigmoid or logit function Οβ²β(q)=logβ‘q1βq=Οβ1β(q)=logitβ‘(q), Ο (q)= q1-q=Ο^-1(q)=logit(q), (31) where Οβ‘(p)=1/(1+eβp)Ο(p)=1/(1+e^-p). This is the key property we use in Β§3 for constructing Bregman divergences whose geometry is matched to an activation function. The corresponding Bregman divergence is DΟβ(p,q) D_Ο(p,q) =Οβ‘(p)βΟβ‘(q)βΟβ²β(q)β(pβq) =Ο(p)-Ο(q)-Ο (q)(p-q) (32) =pβlogβ‘pq+(1βp)βlogβ‘1βp1βq =p pq+(1-p) 1-p1-q (33) =DKL(Bern(p)β₯Bern(q)). =D_KL (Bern(p)\,\|\,Bern(q) ). (34) The Bernoulli KL divergence is therefore generated by the negative binary entropy. Appendix C Derivation of Bregman-PC updates C.1 State updates Here we derive the Bregman-PC inference dynamics of Eq. 22. Recall from Eq. 18 that the Bregman-PC energy is β±B=ββ=1Lβ1DΞ¨β(β,Οβ‘(β))+ββ‘(L,). _B= _ =1^L-1D_ ( z ,Ο( a ) )+L( a^L, y). (35) For a hidden layer β , the activity β z appears in two terms of the energy: (i) its local prediction-error term β±BβF_B , and (i) the prediction-error term of the subsequent layer β±Bβ+1F_B +1. Hence, the activity gradient is given by βββ±B=βββ±Bβ+βββ±Bβ+1. _ z F_B= _ z F_B + _ z F_B +1. (36) Local contribution. From the definition of a Bregman divergence (Eq. 7), its gradient with respect to the first argument is βDΞ¨β(,)=βΞ¨β()ββΞ¨β(). _ pD_ ( p, q)=β ( p)-β ( q). (37) Setting =β p= z and =Οβ‘(β) q=Ο( a ) therefore gives βββ±Bβ=βΞ¨β(β)ββΞ¨β(Οβ‘(β))=Οβ1β(β)βΟβ1β(Οβ‘(β))=βββ, _ z F_B =β ( z )-β (Ο( a ) )=Ο^-1( z )-Ο^-1 (Ο( a ) )= u - a , (38) where we used β=βΞ¨β(β)=Οβ1β(β) u =β ( z )=Ο^-1( z ). Feedback contribution. The activity β z enters the subsequent layer only through the next preactivation β+1=β+1ββ a +1= W +1 z . Applying the chain rule therefore gives βββ±Bβ+1 _ z F_B +1 =(ββ+1ββ)β€βββ+1β±Bβ+1=(β+1)β€βββ+1β±Bβ+1. = ( β a +1β z ) _ a +1F_B +1= ( W +1 ) _ a +1F_B +1. (39) Using the activation-matching property of Eq. 14 derived in Β§3, we get βββ±Bβ+1=β(β+1)β€ββ+1. _ z F_B +1=- ( W +1 ) +1. (40) For the final hidden layer, the same expression follows from our definition LβββLβ ^L - _ a^LL. Substituting the two gradient contributions from Eqs. 38 and 40 into the mirror-gradient-flow dynamics dββ/dβΟβββββ±Bd u /dΟ - _ z F_B gives dββdβΟ=ββ+β+(β+1)β€ββ+1, d u dΟ=- u + a + ( W +1 ) +1, (41) which is Eq. 22. C.2 Weight updates We next derive the Bregman-PC weight updates of Eq. 24. For a hidden layer ββ1,β¦,Lβ1 β\1,β¦,L-1\, the weights β W enter the energy only through the corresponding preactivation β=ββββ1 a = W z -1. Therefore, applying the chain rule gives βββ±B=βββ±Bββ(ββ1)β€. _ W F_B= _ a F_B ( z -1 ) . (42) Using the activation-matching property of Eq. 14, βββ±Bβ=ββ _ a F_B =- , we obtain dββdβt=βΞ·ββββ±B=Ξ·βββ(ββ1)β€, d W dt=-Ξ· _ W F_B=Ξ·\, ( z -1 ) , (43) which is Eq. 24. The same update holds for the output-layer weights as LβββLβ ^L - _ a^LL. Appendix D Candidate analog implementation of Bregman-PC inference aiβa_i uiβu_i RRCCΟ βz_i Ο ^iβ z_i β+1 +1(β+1)β€( W +1) weighted sum Ifb,iβI_fb,i no activation derivativerequired Ifb,iβ=1RββjWjβiβ+1βΞ΅jβ+1 I_fb,i = 1R _jW_ji +1 _j +1RβCβdβuiβdβt=βuiβ+aiβ+βjWjβiβ+1βΞ΅jβ+1 RC\, du_i dt=-u_i +a_i + _jW_ji +1 _j +1 Figure 6: Candidate analog implementation of Bregman-PC inference. The capacitor voltage represents the dual state uiβu_i . A resistor connecting it to the bottom-up preactivation aiβa_i generates the local restoring current (aiββuiβ)/R(a_i -u_i )/R, while a weighted-sum circuit supplies a feedback current proportional to βjWjβiβ+1βΞ΅jβ+1 _jW_ji +1 _j +1. The same nonlinear transfer function can be used to obtain the activity ziβ=Οβ‘(uiβ)z_i =Ο(u_i ) and prediction z^iβ=Οβ‘(aiβ) z_i =Ο(a_i ). Kirchhoffβs current law reproduces the Bregman-PC inference dynamics in Eq. 22 up to an overall rescaling of time by RβCRC. Unlike standard nonlinear PC (Figure 1), the feedback pathway requires no explicit evaluation of the activation derivative. The derivative-free dynamics of Bregman PC (Eq. 22) admit a simple analog interpretation. Consider the inference dynamics of a single unit i in layer β , dβuiβdβΟ=βuiβ+aiβ+βjWjβiβ+1βΞ΅jβ+1. du_i dΟ=-u_i +a_i + _jW_ji +1 _j +1. (44) A candidate circuit implementation is shown in Figure 6. The dual state uiβu_i is represented by the voltage across a capacitor C, while the bottom-up preactivation aiβa_i is represented by a second voltage. Connecting these nodes through a resistor R generates the current Ilocal,iβ=aiββuiβR. I_local,i = a_i -u_i R. (45) The weighted feedback term can be supplied as an additional current Ifb,iβ=1RββjWjβiβ+1βΞ΅jβ+1. I_fb,i = 1R _jW_ji +1 _j +1. (46) Kirchhoffβs current law then gives Cβdβuiβdβt C du_i dt =aiββuiβR+Ifb,iβ,RβCβdβuiβdβt=βuiβ+aiβ+βjWjβiβ+1βΞ΅jβ+1, = a_i -u_i R+I_fb,i , RC du_i dt=-u_i +a_i + _jW_ji +1 _j +1, (47) where t denotes physical circuit time. This reproduces Eq. 22 up to an overall rescaling of time by the circuit time constant RβCRC. The same nonlinear transfer function ΟΟ can be applied to the state and preactivation voltages to produce the activity ziβ=Οβ‘(uiβ)z_i =Ο(u_i ) and prediction z^iβ=Οβ‘(aiβ) z_i =Ο(a_i ), from which the local prediction error can be formed. Crucially, the feedback pathway requires only a weighted sum of next-layer prediction errors. Standard nonlinear PC would additionally require each error to be multiplied by the activation derivative Οβ²β(β+1)Ο ( a +1) before this weighted sum is formed. We emphasise that Figure 6 is intended only as a schematic circuit interpretation of the Bregman-PC inference dynamics, not as a complete hardware design including the learning dynamics. The circuitry generating the bottom-up weighted sums, local prediction errors, and programmable synaptic weights is also left unspecified and may be realised using different analog technologies. Appendix E Experimental details Unless otherwise stated, trained networks used no biases, tanh as hidden activation function, and MSE loss. Comparisons between BP, standard PC, and Bregman PC used identical initial weights. Standard PC used Euclidean GD on the activities (Eq. 4), while Bregman PC used the corresponding mirror-descent dynamics (Eq. 22), with both initialised from the feedforward activities. The inference step size Ο and number of steps T are specified separately for each experiment below. Toy Bregman geometry and inference trajectories (Figures 2-4 and Figures 7β12). These figures illustrate the activation-matched construction for Οβtanh,Ο,idΟβ\ ,Ο,id\ and are not trained models. The matched Bregman eneriges (up to an additive constant) are Οβ‘(z)=zβarctanhβ‘(z)+12βlogβ‘(1βz2)Ο(z)=zarctanh(z)+ 12 (1-z^2) for tanh , Οβ‘(z)=zβlogβ‘z+(1βz)βlogβ‘(1βz)Ο(z)=z z+(1-z) (1-z) for sigmoid, and Οβ‘(z)=12βz2Ο(z)= 12z^2 for the identity. Figure 4 compares Euclidean and mirror inference on the same two-hidden-unit Bregman-PC energy β±Bβ(z1,z2)=DΟβ(z1,Οβ‘(w1βx))+DΟβ(z2,Οβ‘(w2βz1))+12β(w3βz2βy)2, _B(z^1,z^2)=D_Ο (z^1,Ο(w_1x) )+D_Ο (z^2,Ο(w_2z^1) )+ 12(w_3z^2-y)^2, (48) Both flows are integrated with classical RK4 (ΞβΟ=0.04 Ο=0.04, 300300 steps) from identical initial activities. Euclidean inference follows dβ/dβΟ=βββ±Bd z/dΟ=- _ zF_B, while mirror inference follows dβ/dβΟ=βββ±Bd u/dΟ=- _ zF_B with =Οβ‘() z=Ο( u). The equilibrium marker is the terminal state of a mirror trajectory, and grey contours are 1212 equally spaced energy levels between the 22nd and 9292nd percentiles of β±BF_B on the plotted grid. Classification and generation tasks (Table 2). We trained fully connected MLPs with two hidden layers of width 256 on MNIST and Fashion-MNIST. Images were flattened and standardised (MNIST: mean 0.13070.1307, std. 0.30810.3081; Fashion-MNIST: mean 0.50.5, std. 0.50.5), while labels were encoded as one-hot vectors. Hyperparameters were selected independently for each algorithm by mean performance on a validation set over three random seeds. All methods used Adam with no weight decay or learning-rate schedule, batch size 6464, and 1515 epochs. We swept the parameter learning rate over Ξ·β5Γ10β4,10β3,5Γ10β3,10β2Ξ·β\5Γ 10^-4,10^-3,5Γ 10^-3,10^-2\. For PC inference, we jointly swept the activity step size Οβ10β3,5Γ10β3,10β2,5Γ10β2,10β1,5Γ10β1,1Οβ\10^-3,5Γ 10^-3,10^-2,5Γ 10^-2,10^-1,5Γ 10^-1,1\ and the number of inference steps Tβ3,5,10,20,50Tβ\3,5,10,20,50\. All methods were evaluated with a single feedforward pass. Table 2 reports the final classification accuracy and the per-pixel MSE, both averaged over the validation set. For generation, the MSE is in the original [0,1][0,1] pixel space so that the two datasets share a common scale. Saddle-to-saddle dynamics (Figure 5a). Following Innocenti et al. (2024), we trained a 55-layer MLP (44 hidden layers of width 500500) on MNIST, comparing BP, standard PC and Bregman PC. Weights were initialised i.i.d. from a zero-mean Gaussian with standard deviation Ο=5Γ10β3Ο=5Γ 10^-3. We used SGD, batch size 6464, learning rate Ξ·=0.05Ξ·=0.05, and 3535 epochs. PC inference used T=50T=50 steps with Ο=0.1Ο=0.1. Figure 5a plots the training loss during training. Bregman-PC convergence to BP (Figure 5b). We trained residual MLPs on a single CIFAR-10 batch of size 6464, using the ΞΌ parameterisation of Innocenti et al. (2026). In particular, weights are initialised from a standard Gaussian, residual branches are scaled by 1/NβL1/ NL, and the network output by 1/N1/N. We compared all algorithms over the grid Nβ2,8,32,128,512,2048Nβ\2,8,32,128,512,2048\ and Lβ2,4,8,16,32Lβ\2,4,8,16,32\. All methods used Adam with base learning rate 10β310^-3, rescaled as Ξ·/NβLΞ·/ NL, for 100100 full-batch updates. PC inference used Ο=0.3Ο=0.3 and T=100β(Lβ1)T=100(L-1) steps. The heatmap shows the cosine similarity between the Bregman-PC and BP parameter gradients after inference, at the last training step, averaged over three random seeds. Appendix F Compute resources The experiments associated with Table 2 and Figure 5 were run on TPUs supported by Google Cloud Credits awarded through the Google TPU Research Cloud Programme (2026). All other simulations were run on a laptop CPU and took no more than a few minutes. Appendix G Supplementary figures Figure 7: Bregman divergence for the sigmoid function Ο. (a) The activation function Οβ‘(a)=Οβ‘(a)Ο(a)=Ο(a). (b) Its inverse, Οβ1β(z)=logitβ‘(z)Ο^-1(z)=logit(z). (c) Integrating the inverse activation gives the matching potential Οβ‘(z)=zβarctanhβ‘(z)+12βlogβ‘(1βz2)Ο(z)=zarctanh(z)+ 12 (1-z^2), up to an additive constant. Figure 8: Bregman divergence for the identity function. (a) The activation function Οβ‘(a)=idβ‘(a)Ο(a)=id(a). (b) Its inverse, Οβ1β(z)=zΟ^-1(z)=z. (c) Integrating the inverse activation gives the matching potential Οβ‘(z)=zβarctanhβ‘(z)+12βlogβ‘(1βz2)Ο(z)=zarctanh(z)+ 12 (1-z^2), up to an additive constant. Figure 9: Activation matching for sigmoid. As in Figure 3, the standard squared error gives βaβ±=βΞ΅βΟβ²β(a) _aF=- Ο (a), which is attenuated as the activation saturates, whereas the matched Bregman divergence yields βaDΟ=βΞ΅ _aD_Ο=- . Figure 10: Activation matching for a linear activation. For Οβ‘(a)=aΟ(a)=a, we have Οβ²β(a)=1Ο (a)=1, so the standard squared error and the matched Bregman divergence coincide, both yielding βaβ±=βaDΟ=βΞ΅ _aF= _aD_Ο=- . Figure 11: Mirror inference with sigmoid activation. Analogous to Figure 4 for tanh , inference trajectories are shown for a two-hidden-state Bregman PCN with a sigmoid activation. Mirror inference follows a curvature-preconditioned gradient flow in primal coordinates and the corresponding simpler dynamics dβ/dβΟ=βββ±Bd u/dΟ=- _ zF_B in dual coordinates. Figure 12: Mirror inference with a linear activation. Analogous to Figure 4, but for Οβ‘(a)=aΟ(a)=a. In this case, the matched Bregman potential is quadratic (Figure 8), primal and dual coordinates coincide, and mirror inference reduces exactly to ordinary Euclidean gradient flow, yielding identical trajectories. Figure 13: Classification learning curves on MNIST and Fashion-MNIST. Validation classification accuracy during training for BP, standard PC, and Bregman PC on MNIST and Fashion-MNIST. Final validation accuracies are reported in Table 2. Figure 14: Generated samples on MNIST and Fashion-MNIST. Representative images generated by models trained with BP, standard PC, and Bregman PC on MNIST and Fashion-MNIST. Corresponding reconstruction losses are reported in Table 2.