Divide-and-Conquer Predictive Coding: a structured Bayesian inference algorithm

Authors: eli  sennesh
hao  wu
tommaso  salvatori
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-14 23:38:47

Abstract: Unexpected stimuli induce "error" or "surprise" signals in the brain. The theory of predictive coding promises to explain these observations in terms of Bayesian inference by suggesting that the cortex implements variational inference in a probabilistic graphical model. However, when applied to machine learning tasks, this family of algorithms has yet to perform on par with other variational approaches in high-dimensional, structured inference problems. To address this, we introduce a novel predictive coding algorithm for structured generative models, that we call divide-and-conquer predictive coding (DCPC); it differs from other formulations of predictive coding, as it respects the correlation structure of the generative model and provably performs maximum-likelihood updates of model parameters, all without sacrificing biological plausibility. Empirically, DCPC achieves better numerical performance than competing algorithms and provides accurate inference in a number of problems not previously addressed with predictive coding. We provide an open implementation of DCPC in Pyro on Github.

Main paper: 
 introduction: 
 In recent decades, the fields of cognitive science, machine learning, and theoretical neuroscience have borne witness to a flowering of successes in modeling intelligent behavior via statistical learning. Each of these fields has taken a different approach: cognitive science has studied probabilistic inverse inference [Chater et al., 2006 , Pouget et al., 2013 , Lake et al., 2017] in models of each task and environment, machine learning has employed the backpropagation of errors [Rumelhart et al., 1986 , Lecun et al., 2015 , Schmidhuber, 2015] , and neuroscience has hypothesized that predictive coding (PC) [Srinivasan et al., 1982 , Rao and Ballard, 1999 , Friston, 2005 , Bastos et al., 2012 , Spratling, 2017 , Hutchinson and Barrett, 2019 , Millidge et al., 2021] may explain neural activity in perceptual tasks. These approaches share in common a commitment to "deep" models, in which task processing emerges from the composition of elementary units.In machine learning, PC-based algorithms have recently gained popularity for their theoretical potential to provide a more biologically plausible alternative to backpropagation for training neural networks [Salvatori et al., 2023 , Song et al., 2024] . However, PC does not perform comparably in these tasks to backpropagation due to limitations in current formulations. First, predictive coding for gradient calculation typically models every node in the computation graph with a Gaussian, and hence fails to express many common generative models. Recent work on PC has addressed this by allowing approximating non-Gaussian energy functions with samples [Pinchetti et al., 2022] . Second, the Laplace approximation to the posterior infers only a maximum-a-posteriori (MAP) estimate and Gaussian covariance for each latent variable, keeping PC from capturing multimodal or correlated distributions. Third, this loose approximation to the posterior distribution results in inaccurate, high-variance updates to the parameters of the generative model. Where classical predictive coding has layers communicate through shared error units, divide-andconquer predictive coding separates recurrent from "bottom-up" error pathways to target complete conditional distributions rather than posterior marginal distributions.In this work we propose a new algorithm, divide-and-conquer predictive coding (DCPC), for approximating structured target distributions with populations of Monte Carlo samples. DCPC goes beyond Gaussian assumptions, and decomposes the problem of sampling from structured targets into local coordinate updates to individual random variables. These local updates are informed by unadjusted Langevin proposals parameterized in terms of biologically plausible prediction errors. Nesting the local updates within divide-and-conquer Sequential Monte Carlo [Lindsten et al., 2017 , Kuntz et al., 2024] ensures that DCPC can target any statically structured graphical model, while Theorem 2 provides a locally factorized way to learn model parameters by maximum marginal likelihood.DCPC also provides a computational perspective on the canonical cortical microcircuit [Bastos et al., 2012 , 2020 , Campagnola et al., 2022] hypothesis in neuroscience. Experiments have suggested that deep laminar layers in the cortical microcircuit represent sensory imagery, while superficial laminar represent raw stimulus information [Bergmann et al., 2024] ; experiments in a predictive coding paradigm specifically suggested that the deep layers represent "predictions" while the shallow layers represent "prediction errors". This circuitry could provide the brain with its fast, scalable, generic Bayesian inference capabilities. Figure 1 compares the computational structure of DCPC with that of previous PC models. The following sections detail the contributions of this work:• Section 3 defines the divide-and-conquer predictive coding algorithm and shows how to use it as a variational inference algorithm;• Section 4 examines under what assumptions the cortex could plausibly implement DCPC, proving two theorems that contribute to biological plausibility;• Section 5 demonstrates DCPC experimentally in head-to-head comparisons against recent generative models and inference algorithms from the predictive coding literature.Section 2 will review the background for Section 3's algorithm: the problem predictive coding aims to solve and a line of recent work adressing that problem from which this paper draws.Typically the marginal and posterior densities have no closed form, so learning and inference algorithms treat the joint distribution as a closed-form unnormalized density over the latent variables; its integral then gives the normalizing constant for approximation γ θ (z) := p θ (x, z), Z θ := z∈Z γ θ (z) dz = p θ (x), π θ (z) := γ θ (z) Z θ .Neal and Hinton [1998] reduced empirical Bayes to minimization of the variational free energy:F(θ, q) := E z∼q(z) -log γ θ (z) q(z) ≥ -log Z(θ).(2)The ratio of densities in Equation 2 is an example of a weight used to approximate a distribution known only up to its normalizing constant. The proposal distribution q(z) admits tractable sampling, while the unnormalized target density γ θ (z) admits tractable, pointwise density evaluation.Predictive Coding Computational neuroscientists now often hypothesize that predictive coding (PC) can optimize the above family of objective functionals in a local, neuronally plausible way [Millidge et al., 2021 [Millidge et al., , 2023]] . More in detail, it is possible to define this class of algorithms as follows:Definition 1 (Predictive Coding Algorithm). Consider approximate inference in a model p θ (x, z) using an algorithm A. Salvatori et al. [2023] calls A a predictive coding algorithm if and only if:1. It maximizes the model evidence log p θ (x) by minimizing a variational free energy;2. The proposal q(z) = z∈z q(z) factorizes via a mean-field approximation; and 3. Each proposal factor is a Laplace approximation (i.e. q µ (z) := N (µ, Σ(µ))).In contrast to predictive coding, particle algorithms approach empirical Bayes problems by setting the proposal to a collection of weighted particles (w k , z k ) drawn from a sampling algorithm meeting certain conditions (see Definition 4 in Appendix B). Any proposal meeting these conditions (see Proposition 1 in Appendix B and Naesseth et al. [2015] , Stites et al. [2021] ) defines a free energy functional, analogous to Equation 2 in upper-bounding the model surprisal:F(θ, q) := E w,z∼q(w,z) [-log w] =⇒ F(θ, q) ≥ -log Z(θ).This paper builds on the particle gradient descent (PGD) algorithm of Kuntz et al. [2023] , that works as follows: At each iteration t, PGD diffuses the particle cloud q K (z) =1 K K k=1 δ z k (z) across the target log-density with a learning rate η and independent Gaussian noise; it then updates the parameters θ by ascending the gradient of the log-likelihood, estimated by averaging over the particles. The update rules are then the following:EQUATIONEQUATION)The above equations target the joint density of an entire graphical model 1 . When the prior p θ (z) factorizes into many separate conditional densities, achieving high inference performance often requires factorizing the inference network or algorithm into conditionals as well [Webb et al., 2018] .Estimating the gradient of the entire log-joint, as in PGD and amortized inference [Dasgupta et al., 2020 , Peters et al., 2024] , also requires nonlocal backpropagation. To provide a generic inference algorithm for high-dimensional, structured models using only local computations, Section 3 will apply Equation 3 to sample individual random variables in a joint density, combine the coordinate updates via sequential Monte Carlo, and locally estimate gradients for model parameters via Equation 4. 
divide-and-conquer predictive coding: 
 The previous section provided a mathematical toolbox for constructing Monte Carlo algorithms based on gradient updates and a working definition of predictive coding. This section will combine those 1 : Comparison of divide-and-conquer predictive coding (DCPC) against other predictive coding algorithms. DCPC provides the greatest flexibility: arbitrary differentiable generative models, an empirical approximation to the posterior, and sampling according to the target's conditional structure. tools to generalize the above notion of predictive coding, yielding the novel divide-and-conquer predictive coding (DCPC) algorithm. Given a causal graphical model, DCPC will approximate the posterior with a population q(z) of K samples, while also learning θ explaining the data. This will require deriving local coordinate updates and then parameterizing them in terms of prediction errors.Let us assume we again have a causal graphical model p θ (x, z) locally parameterized by θ and factorized (as in Equation 1) into conditional densities for each x ∈ x and z ∈ z. DCPC then requires two hyperparameters: a learning rate η ∈ R + , and particle count K ∈ N + , and is initialized (at t = 0) via a population of predictions by ancestor sampling defined as z 0 ∼ z∈z p θ (z 0 | Pa(z 0 )).DCPC aims to minimize the variational free energy (Equation 2). The optimal proposal q * for each random variable would equal, if it had closed form, the complete conditional density for that variable, containing all information from other random variablesq * (z t | z t-1 ) ∝ γ θ (z; z \z ) = p θ (z | Pa(z)) v∈Ch(z) p θ (v | Pa(v)).(5)We observe that the prediction errors ε z in classical predictive coding, usually defined as the precision weighted difference between predicted and actual value of a variable, can be seen as the score function of a Gaussian, where the score is the gradient with respect to the parameter z of the log-likelihood:ε z := ∇ z log N (z, τ ) = τ (x -z) ;When given the ground-truth parameter z, the expected score functionE x∼p(x|z) [∇ z log p(x | z)]under the likelihood becomes zero, making score functions a good candidate for implementing predictive coding. We therefore define ε z in DCPC as the complete conditional's score functionEQUATIONThis gradient consists of a sum of local prediction-error terms: one for the local "prior" on z and one for each local "likelihood" of a child variable. By defining the prediction error as a sum of local score functions, we write Equation 3 in terms of ε z (Equation 6) and the preconditioner of Definition 3:q η (z t | ε t z , z t-1 ) := N z t-1 + η ΣI ε t z , 2η ΣI .The resulting proposal now targets the complete conditional density (Equation 5), simultaneously meeting the informal requirement of Definition 1 for purely local proposal computations while also "dividing and conquering" the sampling problem into lower-dimensional coordinate updates.Since the proposal from which we can sample by predictive coding is not the optimal coordinate update, we importance weight for the true complete conditional distribution that is optimalEQUATIONresampling with respect to these weights corrects for discretization error, yields particles distributed according to the true complete conditional, and estimates the complete conditional's normalizerRESAMPLE z t , u t z ∼ π θ t-1 (z t | z \z ), Ẑθ t-1 (z \z ) t := 1 K K k=1 u t,k z .The recursive step of "Divide and Conquer" Sequential Monte Carlo [Lindsten et al., 2017 , Kuntz et al., 2024] exploits the estimates Ẑθ t-1 (z \z ) t to weigh the samples for the complete target density w t θ t-1 = p θ t-1 (x, z t ) z∈z γ θ (z t ; z \z ) z∈z Ẑθ t-1 (z \z ) t . (8)ε z ← ∇ z log p θ t-1 (z | Pa(z)) ▷ Local prediction error 5: ε z ← ε z + v∈Ch(z) ∇ z log p θ t-1 (v | Pa(v)) ▷ Children's prediction errors 6: ΣI ← ÎK (ε 1:K z ) -1 1 d Tr[ ÎK (ε 1:K z ) -1 ]▷ Estimate precision of prediction errors 7: z t ∼ q η (z t | ε z , z t-F t ← -1 K K k=1 log p θ t-1 (x,z t,k ) z∈z γ θ t-1 (z t,k ;z t,k \z )z∈z Ẑθ t-1 (z \z ) t ▷ Update free energy 12:θ t ← θ t-1 + η 1 K K k=1 ∇ θ t-1 log p θ t-1 (x, z t,k ) ▷ Update parameters13: return z T , θ T , F T ▷ Output: updated particles, parameters, free energy By Proposition 1, log-transforming these weights estimates the free energy (Equation 2):F t (z t-1 , θ t-1 ) := E q * (z t |z t-1 ) -log p θ t-1 (x, z t ) q * (z t | z t-1 ) ≈ E q -log w t θ t-1 .Theorem 3 in Appendix B shows that the gradient∇ θ t-1 F t = E q [-∇ θ t-1 log p θ t-1 (x, z t )]of the above estimator equals the expected gradient of the log-joint distribution. Descending this gradient θ t := θ t-1 -η∇ θ t-1 F t enables DCPC to learn model parameters θ.The above steps describe a single pass of divide-and-conquer predictive coding over a causal graphical model. Algorithm 1 shows the complete algorithm, consisting of nested iterations over latent variables z ∈ z (inner loop) and iterations t ∈ T (outer loop). DCPC satisfies criteria (1) and (2) of Definition 1, and relaxes criterion (3) to allow gradient-based proposals beyond the Laplace assumption. As with Pinchetti et al. [2022] and Oliviers et al. [2024] , relaxing the Laplace assumption enables much greater flexibility in approximating the model's true posterior distribution.Question: For each theoretical result, does the paper provide the full set of assumptions and a complete (and correct) proof? Answer: [Yes] Justification: Proofs can be found in Appendix B. Guidelines:• The answer NA means that the paper does not include theoretical results.• All the theorems, formulas, and proofs in the paper should be numbered and crossreferenced. • All assumptions should be clearly stated or referenced in the statement of any theorems.• The proofs can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide a short proof sketch to provide intuition. • Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material. • Theorems and Lemmas that the proof relies upon should be properly referenced. 
biological plausibility: 
 Different works in the literature consider different criteria for biological plausibility. This paper follows the non-spiking predictive coding literature and considers an algorithm biologically plausible if it performs only spatially local computations in a probabilistic graphical model [Whittington and Bogacz, 2017] , without requiring a global control of computation. However, while in the standard literature locality is either directly defined in the objective function [Rao and Ballard, 1999] , or derived from a mean-field approximation to the joint density [Friston, 2005] , showing that the updates of the parameters of DCPC require only local information is not as trivial. To this end, in this section we first formally show that DCPC achieves decentralized inference of latent variables z (Theorem 1), and then that also the parameters θ are updated via local information (Theorem 2).Gibbs sampling provides the most widely-used algorithm for sampling from a high-dimensional probability distribution by local signaling. It consists of successively sampling coordinate updates to individual nodes in the graphical model by targeting their complete conditional densities π θ (z | x, z \z ). Theorem 1 demonstrates that DCPC's coordinate updates approximate Gibbs sampling.Theorem 1 (DCPC coordinate updates sample from the true complete conditionals). Each DCPC coordinate update (Equation 7) for a latent z ∈ z samples from z's complete conditional (the normalization of Equation 5). Formally, for every measurable h : Z → R, resampled expectations with respect to the DCPC coordinate update equal those with respect to the complete conditionalE z∼qη (z|z t-1 ,ε t z ) E u∼δ(u),z ′ ∼RESAMPLE(z,uz) [h(z)] = z∈Z h(z) π θ (z | z \z ) dz.Proof. See Corollary 4.1 in Appendix B.We follow the canonical cortical microcircuit hypothesis of predictive coding [Bastos et al., 2012 , Gillon et al., 2023] or predictive routing [Bastos et al., 2020] . Consider a cortical column representing z ∈ z; the θ, α/β, and γ frequency bands of neuronal oscillations [Buzsáki and Draguhn, 2004] could synchronize parallelizations (known to exist for simple Gibbs sampling in a causal graphical model [Gonzalez et al., 2011] ) of the loops in Algorithm 1. From the innermost to the outermost and following the neurophysiological findings of Bastos et al. [2015] , Fries [2015] , γ-band oscillations could synchronize the bottom-up conveyance of prediction errors (lines 4-6) from L2/3 of lower cortical columns to L4 of higher columns, β-band oscillations could synchronize the top-down conveyance of fresh predictions (implied in passing from s to s + 1 in the loop of lines 2-9) from L5/6 of higher columns to L1+L6 of lower columns, and θ-band oscillations could synchronize complete attention-directed sampling of stimulus representations (lines 1-11). Figure 5 in Appendix A visualizes these hypotheses for how neuronal areas and connections could implement DCPC.Biological neurons often spike to represent changes in their membrane voltage [Mainen and Sejnowski, 1995 , Lundstrom et al., 2008 , Forkosh, 2022] , and some have even been tested and found to signal the temporal derivative of the logarithm of an underlying signal [Adler and Alon, 2018, Borba et al., 2021] . Theorists have also proposed models [Chavlis and Poirazi, 2021, Moldwin et al., 2021] under which single neurons could calculate gradients internally. In short, if neurons can represent probability densities, as many theoretical proposals and experiments suggest they can, then they can likely also calculate the prediction errors used in DCPC. Theorem 2 will demonstrate that given the "factorization" above, DCPC's model learning requires only local prediction errors.Theorem 2 (DCPC parameter learning requires only local gradients in a factorized generative model).Consider a graphical model factorized according to Equation 1, with the additional assumption that the model parameters θ ∈ Θ = x∈x Θ x × z∈z Θ z factorize disjointly. Then the gradient ∇ θ F(θ, q) of DCPC's free energy similarly factorizes into a sum of local particle averagesEQUATIONProof. See Proposition 5 in Appendix B.Our practical implementation of DCPC, evaluated in the experiments above, takes advantage of Theorem 2 to save memory by detaching samples from the automatic differentiation graph in the forward ancestor-sampling pass through the generative model.Finally, DCPC passes from local coordinate updates to the joint target density via an importance resampling operation, requiring that implementations synchronously transmit numerical densities or log-densities for the freshly proposed particle population. While phase-locking to a cortical oscillation may make this biologically feasible, resampling then requires normalizing the weights. Thankfully, divisive normalization appears ubiquitously throughout the brain [Carandini and Heeger, 2012] , as well as just the type of "winner-take-all" circuit that implements a softmax function (e.g. for normalizing and resampling importance weights) being ubiquitous in crosstalk between superficial and deep layers of the cortical column [Liu, 1999, Douglas and Martin, 2004] .Question: Does the paper fully disclose all the information needed to reproduce the main experimental results of the paper to the extent that it affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)? Answer: [Yes] Justification: The paper includes the learning rates and particle counts, the necessary hyperparameters, for the original experiments, as well as releasing code to run them. Guidelines:• The answer NA means that the paper does not include experiments.• If the paper includes experiments, a No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not. • If the contribution is a dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.• Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is a novel architecture, describing the architecture fully might suffice, or if the contribution is a specific model and empirical evaluation, it may be necessary to either make it possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g., in the case of a large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed. • While NeurIPS does not require releasing code, the conference does require all submissions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. (d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome to describe the particular way they provide for reproducibility.In the case of closed-source models, it may be that access to the model is limited in some way (e.g., to registered users), but it should be possible for other researchers to have some path to reproducing or verifying the results. 
experiments: 
 Divide-and-conquer predictive coding is not the first predictive coding algorithm to incorporate sampling into the inference process, and certainly not the first variational inference algorithm for structured graphical models. This section therefore evaluates DCPC's performance against both models from the predictive coding literature and against a standard deep generative model. Each experiment holds the generative model, dataset, and hyperparameters constant except where noted.We have implemented DCPC as a variational proposal or "guide" program in the deep probabilistic programming language Pyro [Bingham et al., 2019] ; doing so enables us to compute free energy and prediction errors efficiently in graphical models involving neural networks. Since the experiments below involve minibatched subsampling of observations x ∼ B from a dataset D ∼ p(D) of unknown distribution, we replace Equation 9 with a subsampled form (see Welling and Teh [2011] Inference 2 : Negative log-likelihood and mean squared error for MCPC against DCPC on held-out images from the MNISTs. Means and standard deviations are taken across five random seeds.for derivation) of the variational Sequential Monte Carlo gradient estimator [Naesseth et al., 2018 ]EQUATION)We optimized the free energy in all experiments using Adam [Kingma and Ba, 2014] Figure 2 shows that the DLGM structure already requires DCPC to respect hierarchical dependencies. We tested DCPC's performance on elementary reconstruction and generation tasks by using it to train this exact generative model, changing only the likelihood from a discrete Bernoulli to a continuous Bernoulli [Loaiza-Ganem and Cunningham, 2019] . After training we evaluated with a discrete Bernoulli likelihood. Table 2 shows that in terms of both surprise (negative log evidence, with the discrete Bernoulli likelihood) and mean squared reconstruction error, DCPC enjoys better average performance with a lower standard deviation of performance, the latter by an order of magnitude. All experiments used a learning rate η = 0.1 and K = 4 particles.Figure 3 shows an extension of this experiment to EMNIST [Cohen et al., 2017] and Fashion MNIST [Xiao et al., 2017] as well as the original MNIST, with ground-truth images in the top row and their reconstructions from DCPCinferred latent codes below. The ground-truth images come from a 10% validation split of each data-set, on which DCPC only infers particles q K=4 (z).The above datasets do not typically challenge a new inference algorithm. The next experiment will thus attempt to learn representations of color images, as in the widely-used variational autoencoder [Kingma and Welling, 2013] framework, without an encoder network or amortized inference.Image generation with representation learning Zahid et al. [2024] have also recently designed and evaluated Langevin predictive coding (LPC), with differences from both MCPC and DCPC. DN 64 × 64 30 × 150 = 4500 79.0 ± 0.9 Table 3 : FID score comparisons on the CelebA dataset [Liu et al., 2015] . The score for LPC comes from Figure 2 in Zahid et al. [2024] , where they ablated warm-starts and initialized from the prior. We compare DCPC to LPC using the Frechet Inception Distance (FID) [Seitzer, 2020] featured in Zahid et al. [2024] , holding constant the prior, neural network architecture, learning rate on θ, and number of gradient evaluations used to train the parameters θ and latents z. Zahid et al. [2024] evaluated a variety of scenarios and reported that their training could converge quickly when counted in epochs, but they accumulated gradients of θ over inference steps. We compare to the results they report after 15 epochs with 300 inference steps applied to latents initialized from the prior, equivalent to 15×300 = 4500 gradient steps on θ per batch, replicating their batch size of 64. Since Algorithm 1 updates θ only in its outer loop, we set S = 30 and ran DCPC for 150 epochs, with η = 0.001 and K = 16. Table 3 shows that DCPC outperforms LPC in apples-to-apples generative quality, though not to the point of matching other model architectures2 by inference quality alone.Figure 4 shows reconstructed images from the validation set (left) and samples from the posterior predictive generative model (right). There is blurriness in the reconstructions, as often occurs with variational autoencoders, but DCPC training allows the network to capture background color, hair color, direction in which a face is looking, and other visual properties. Figure 4a shows reconstructions over the validation set, while Figure 4b shows samples from the predictive distribution. Kuntz et al. [2023] also reported an experiment on CelebA in terms of FID score, at the lower 32 × 32 resolution. Since they provided both source code and an exact mathematical description, we were able to run an exact, head-to-head comparison with PGD. The line in Table 3 evaluating DCPC with PGD's example neural architecture at the 32 × 32 resolution (with similar particle count and learning rate) demonstrates a significant improvement in FID for DCPC, alongside reduced FID variance.The initial DLGM experiments on the MNIST-alike datasets were performed on a desktop workstation with 128GB of RAM and an NVIDIA Quadro P4000 with 8GB of VRAM. Experiments on CelebA were conducted on an NVIDIA DGX equipped with eight (8) NVIDIA A100's, each with 80GB of VRAM. The latter compute infrastructure was also used for unpublished experiments, on several different datasets, in structured time-series modeling.6 Related Work Pinchetti et al. [2022] expanded predictive coding beyond Gaussian generative models for the first time, applying the resulting algorithm to train variational autoencoders by variational inference and transformer architectures by maximum likelihood. DCPC, in turn, broadens predictive coding to target arbitrary probabilistic graphical models, following the broadening in Salvatori et al. [2022] to arbitrary deterministic computation graphs. DCPC follows on incremental predictive coding [Salvatori et al., 2024] in quickly alternating between updates to random variables and model parameters, giving an incremental EM algorithm [Neal and Hinton, 1998 ]. Finally, Zahid et al. [2024] and Oliviers et al. [2024] also recognized the analogy between predictive coding's prediction errors and the score functions used in Langevin dynamics for continuous random variables.There exists a large body of work on how neurobiologically plausible circuits could implement probabilistic inference. Classic work by Shi and Griffiths [2009] provided a biologically plausible implementation of hierarchical inference via importance sampling; DCPC proceeds from importance sampling as a foundation, while parameterizing the proposal distribution via prediction errors. Recent work by Fang et al. [2022] studied neurally plausible algorithms for sampling-based inference with Langevin dynamics, though only for a Gaussian generative model of sparse coding. Golkar et al. [2022] imposed a whitening constraint on a Gaussian generative model for biological plausibility. Finally, Dong and Wu [2023] and Zahid et al. [2024] both suggest mechanisms for employing momentum to reduce gradient noise in a biologically plausible sampling algorithm; the former intriguingly analogize their momentum term to neuronal adaptation. To conclude, other works have already implemented predictive coding models for image generation tasks, a notable example being the neural generative coding framework Ororbia and Kifer [2022] , Ororbia and Mali [2022] .Question: Does the paper provide open access to the data and code, with sufficient instructions to faithfully reproduce the main experimental results, as described in supplemental material? Answer: [Yes] Justification: The paper provides a link to the code in the abstract and introduction. Guidelines:• The answer NA means that paper does not include experiments requiring code. Justification: The paper notes the hyperparameters and datasets used, gives additional details in the appendices, and releases the code.Guidelines:• The answer NA means that the paper does not include experiments.• The experimental setting should be presented in the core of the paper to a level of detail that is necessary to appreciate the results and make sense of them. • The full details can be provided either with the code, in appendix, or as supplemental material. 
conclusion: 
 This paper proposed divide-and-conquer predictive coding (DCPC), an algorithm that efficiently and scalably approximates Gibbs samplers by importance sampling; DCPC parameterizes efficient proposals for a model's complete conditional densities using local prediction errors. Section 4 showed how Monte Carlo sampling can implement a form of "prospective configuration" [Song et al., 2024] , first inferring a sample from the joint posterior density (Theorem 1) and then updating the generative model without a global backpropagation pass ( Theorem 2). Experiments in Section 5 showed that DCPC outperforms the state of the art Monte Carlo Predictive Coding from computational neuroscience, head-to-head, on the simple generative models typically considered in theoretical neuroscience; DCPC also outperforms the particle gradient descent algorithm of Kuntz et al. [2023] while under the constraint of purely local computation. DCPC's Langevin proposals admit the same extension to constrained sample spaces as applied in Hamiltonian Monte Carlo [Brubaker et al., 2012] ; our Pyro implementation includes this extension via Pyro's preexisting support for HMC.DCPC offers a number of ways forward. Particularly, this paper employed naive Langevin proposals, while Dong and Wu [2023] , Zahid et al. [2024] applied momentum-based preconditioning to take advantage of the target's geometry. Yin and Ao [2006] demonstrated that gradient flows of this general kind can also provide more efficient samplers by breaking the detailed-balance condition necessary for the Metropolis-Hastings algorithm, motivating the choice of SMC over MCMC to correct proposal bias. Appendix C derives a mathematical background for an extension of DCPC to discrete random variables. Future work could follow Marino et al. [2018] , Taniguchi et al. [2022] in using a neural network to iteratively map from particles and prediction errors to proposal parameters.DCPC's main limitations are its longer training time, and greater sensitivity to learning rates, than state-of-the-art amortized variational inference trained end-to-end. Such limitations occur frequently in the literature on neuroscience-inspired learning algorithms, as well as in the literature on particlebased algorithms with no parametric form. This work has no singular ethical concerns specific only to DCPC, rather than the broader implications and responsibilities accompanying advancements in biologically plausible learning and Bayesian inference. Figure 5 : Divide-and-conquer predictive coding provides an algorithmic interpretation for some of the connections mapped in the canonical neocortical microcircuit [Bastos et al., 2012 , 2020 , Campagnola et al., 2022] : prediction errors (red) arrive through ascending pathways into the central laminar layer 4, which transmits them up to layers 2/3 (green). These layers combine the incoming errors with a present posterior estimate (green L5→ L2/3 connection) to generate prediction errors for the next cortical area. Eventually the updated predictions flow back down the cortical hierarchy (blue).Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?Answer: [Yes] Justification: Results are given in terms of sample averages and standard deviations, with the standard deviations serving as error bars. The paper notes the number of replicates over which averages are taken.Guidelines:• The answer NA means that the paper does not include experiments.• The authors should answer "Yes" if the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper. • The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions). • The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.) • The assumptions made should be given (e.g., Normally distributed errors).• It should be clear whether the error bar is the standard deviation or the standard error of the mean. • It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have a 96% CI, if the hypothesis of Normality of errors is not verified. • For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out of range (e.g. negative error rates). • If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text. 
experiments compute resources: 
 Question: For each experiment, does the paper provide sufficient information on the computer resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?Answer: [Yes] Justification: The experiments section now includes a paragraph at the end documenting the specific compute resources used in our experiments. Many of our experiments can be trained with less compute than we used, albeit more slowly.Guidelines:• The answer NA means that the paper does not include experiments.• The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage. • The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute.• The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make it into the paper). 
code of ethics: 
 Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics https://neurips.cc/public/EthicsGuidelines?Answer: [Yes] Justification: The paper does not rely on human subjects, nor release any new dataset. Approximate Bayesian inference is a well-understood problem whose downstream ethical impacts do not differ significantly from those of other machine-learning methods.Guidelines:• The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.• If the authors answer No, they should explain the special circumstances that require a deviation from the Code of Ethics. • The authors should make sure to preserve anonymity (e.g., if there is a special consideration due to laws or regulations in their jurisdiction). 
broader impacts: 
 Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?Answer: [NA] Justification: We propose a novel neuroscience-inspired learning algorithm, which does not have immediate societal impacts, positive or negative.Guidelines:• The answer NA means that there is no societal impact of the work performed.• If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact. • Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations. • The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster. • The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology. • If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML). 
safeguards: 
 Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?Answer: [NA] Justification: We do not use data that is not already available online, nor do we release any models for a task more complex than non-customizable image generation.Answer: [NA] Justification: [NA] Guidelines:• The answer NA means that the paper does not involve crowdsourcing nor research with human subjects. • Including this information in the supplemental material is fine, but if the main contribution of the paper involves human subjects, then as much detail as possible should be included in the main paper. • We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution. • For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review. 
 
 Appendix: Alternate image generation/ representation learning As indicated in Section 2, this paper builds upon the particle gradient descent (PGD) algorithm; Kuntz et al. [2023] demonstrated the algorithm's performance by training a generator network on CelebA. Their network employed a Gaussian likelihood with a fixed standard deviation of 0.01, and evaluated a log-joint objective over 100 epochs on exactly 10,000 subsampled data points. The paper then evaluated mean squared error on an inpainting task and the Frechet Inception Distance over data images.When applied to the resulting target density, DCPC amounts to PGD with a resampling step. Table 4 shows the results of training and evaluating the same model described above with DCPC. Since PGD trained for 100 epochs with a batch size of 128, albeit on a 10,000-image subsample of CelebA, we trained with the entire dataset for 100 epochs with batch-size 128.Log-joint FID ↓ VAE (K = 10) -4.3 × 10 5 171.5 ± 0.1 PGD (K = 10) -3.8 × 10 5 100 ± 2.7 DCPC (ours, K = 10) -1.8 × 10 6 82.7 ± 0.9 Table 4 : Log-joint probabilities and FID metrics show how DCPC performs against the original PGD.We suspect that the supplied code for log-joint calculation averages over either particles or batch items differently from how we have evaluated DCPC (e.g. we call mean() without dividing by any further shape dimensions), accounting for the apparent order-of-magnitude difference between log-joints.At the request of reviewers, we have substituted a simplified Figure 1 in the main text for Figure 5 showing how to map DCPC onto laminar microcircuit structure.B Importance sampling and gradient estimation proofs Titsias [2023] introduced optimal estimators for preconditioning Langevin dynamics to adapt with the Fisher information of the target density. Definition 2 gives the most basic estimator for that Fisher information, defined in terms of the score functions we use as prediction errors.Definition 2 (Bayesian Fisher estimator [Titsias, 2023] ). Denoting by ε z the score function (from Equation 6) and letting λ > 0 be a fixed hyperparameter, the Bayesian Fisher estimatorĨK := E ε 1:K z ∼π θ (z|z\z) ε z ε ⊤ z k + λ K I (11)is an empirical estimator of a Bayesian target density's Fisher information (the target density from which the prediction errors were derived) based on a cloud of K particles. When the particles are not yet distributed around a mode of the target (e.g. the score function does not have an average of zero), substituting the empirical covariance for the expected outer product reduces the estimator's biasEQUATIONThe above matrix does not describe the preconditioner that Titsias [2023] actually recommended applying in a Langevin proposal. Definition 3 provides the fully normalized preconditioner. Definition 3 (Predictive coding Fisher preconditioner [Titsias, 2023] ). Using Definition 2 to parameterize a preconditioner for the Langevin dynamics proposal, the predictive coding Fisher preconditioner is the inverse of Equation 12, normalized to have an average eigenvalue of 1EQUATION)Definition 4 generalizes the definition of importance sampling, suitable for recursively constructing sequential Monte Carlo algorithms with changing target densities. Definition 4 (Strict proper weighting for a density). Given an unnormalized density γ θ (z) with corresponding normalizing constant Z(θ) and normalized density π θ (z)Z(θ) := z∈Z γ θ (z) dz π θ (z) := γ θ (z) Z(θ) ,the random variables w, z ∼ q(w, z) are strictly properly weighted [Naesseth et al., 2015] with respect to γ θ (z) if and only if for any measurable test function h : Z → R, the weighted expectation over the proposal q(w, z) equals the expectation under the target γ θ (z)EQUATION)The following two propositions come from the previous work by Wu et al. [2020] , Stites et al. [2021] and Zimmermann et al. [2021] . The reader looking for foundations can see Naesseth et al. [2015] or Chopin and Papaspiliopoulos [2020] .Proposition 1 (The free energy upper-bounds the surprisal). Given a proposal q ϕ (w, z) strictly properly weighted (Definition 4) for the target γ θ (z), the variational free energy provides an upper bound to the target's surprisalEQUATION)Proof. I begin by writing out the free energy (Equation 2) as an expectation of a negative logarithmF(θ, q) = E z,w∼q(z,w) [-log w] .Jensen's Inequality allows moving the expectation into the negative logarithm by relaxing the definition of the variational free energy from an equality to an upper boundF(θ, q) ≥ -log E z,w∼q(z,w) [w] .Setting h(z) = 1, strict proper weighting for an unnormalized density (Definition 4) says the expected weight will be the normalizing constantE z,w∼q(z,w) [w] = Z(θ)which I substitute back in to obtain the desired inequality F(θ, q) ≥ -log Z(θ).Proposition 2 (Weighted expectations approximate the normalized target up to a constant). Given a proposal q ϕ (w, z) strictly properly weighted (Definition 4) for the target γ θ (z) and a measurable test function h : Z → R, weighted expectations under the proposal equal the target's normalizing constant times the test function's expectation under the normalized targetE (w,z)∼q ϕ (w,z) [wh(z)] = Z(θ)E z∼π θ (•) [h(z)] .Proof. Strict proper weighting (Equation 14) states that weighted expectations under the proposal equal integrals over the unnormalized target, and by definition the normalized target equals the unnormalized density over its normalizing constantE w,z∼q(w,z) [wh(z)] = z∈Z h(z) γ θ (z) dz, π θ (z) := γ θ (z) Z(θ) .The second equation expresses the unnormalized target in terms of the normalized oneZ(θ)π θ (z) = γ θ (z),and substituting this expression into the definition of strict proper weighting leads to the desired resultz∈Z h(z) γ θ (z) dz = z∈Z h(z) Z(θ)π θ (z) dz, = Z(θ) z∈Z h(z) π θ (z) dz E w,z∼q(w,z) [wh(z)] = Z(θ)E π θ (z) [h(z)] .Proposition 3 (DCPC's free energy has a pathwise derivative). The free energyF t+1 = E q -log w t+1 θ tconstructed by the population predictive coding algorithm (Algorithm 1) has a pathwise derivative as the expectation of the negative gradient of the log-joint density∇ θ t F t+1 = E q -∇ θ t log p θ t (x, z t+1 ) .Proof. The free energy has an expression in terms of Equation 8F t+1 = E q -log w t+1 θ t w t+1 θ t = p θ t (x, z) z∈z γ θ (z t+1 b ; z \z ) z∈z Ẑθ t (z \z ) t+1 , Ẑθ t (z \z ) t+1 = 1 K K k=1 u t+1,k b u t+1 z = γ θ (z t+1 ; z \z ) q(z t+1 | ε z (z t )) ,and writing out the free energy itself in full shows that many terms cancelq(z t+1 | z t ) = z t+1 b ∈z t+1 q(z t+1 | z t ), F t+1 = E q(z t+1 |z t )   -log p θ t (x, z) z∈z γ θ (z t+1 b ; z \z ) z∈z ¡ ¡ ¡ 1 K ¡ ¡ ¡ ¡ K k=1 γ θ (z t+1 b ; z \z ) q(z t+1 | ε z (z t ))   = E q(z t+1 |z t ) -log p θ t (x, z)@ @ @ @ @ @ @ @ @ z∈z γ θ (z t+1 b ; z \z ) @ @ @ @ @ @ @ @ @ z∈z γ θ (zt+1 b ; z \z ) z∈z q(z t+1 | ε z (z t )) = E q(z t+1 |z t ) -log p θ t (x, z) q(z t+1 | z t ) .The proposal distribution q is a function of the random variable values themselves through the prediction errors, not of the parameters θ. The above expression therefore admits a pathwise derivative [Schulman et al., 2015] , moving the gradient operator into the expectation∇ θ t F t+1 = ∇ θ t E q(z t+1 |z t ) -log p θ t (x, z t+1 ) q(z t+1 | z t ) = E q(z t+1 |z t ) ∇ θ t -log p θ t (x, z t+1 ) q(z t+1 | z t ) = E q(z t+1 |z t ) ∇ θ t -log p θ t (x, z t+1 ) -log q(z t+1 | z t ) = E q(z t+1 |z t ) -∇ θ t log p θ t (x, z t+1 ) -∇ θ t log q(z t+1 | z t ) ∇ θ t F t+1 = E q(z t+1 |z t ) -∇ θ t log p θ t (x, z t+1 ) .Proposition 4 (DCPC coordinate updates are strictly properly weighted for the complete conditionals).Each DCPC coordinate update (Equation 7) for a latent variable z ∈ z is strictly properly weighted (Definition 4) for z's unnormalized complete conditional. For every measurable h :EQUATION)Proof. Expanding the outer expectation into an integral and replacing the Dirac delta with the expression for the local weights transforms Equation 16into z∈Z γ θ (z; z \z ) @ @ @ @ @ @ @ q η (z| z t-1 , ε t z ) E z ′ ∼RESAMPLE(z,uz) [h(z ′ )]@ @ @ @ @ @ @ q η (z| z t-1 , ε t z ) dz = z∈Z h(z) γ θ (z; z \z ) dz;importance resampling also preserves strict proper weighting (see Naesseth et al. [2015] , Stites et al. [2021] and Chopin and Papaspiliopoulos [2020] for proofs), and so this yieldsz∈Z E z ′ ∼RESAMPLE(z,uz) [h(z ′ )] γ θ (z; z \z ) dz = z∈Z h(z) γ θ (z; z \z ) dz z ′ ∈Z h(z ′ ) γ θ (z ′ ; z \z ) dz ′ = z∈Z h(z) γ θ (z; z \z ) dz.Corollary 4.1 (DCPC coordinate updates sample from the true complete conditionals). Each DCPC coordinate update (Equation 7) for a latent z ∈ z samples from z's complete conditional (the normalization of Equation 5). Formally, for every measurable h : Z → R, resampled expectations with respect to the DCPC coordinate update equal those with respect to the complete conditionalE z∼qη (z|z t-1 ,ε t z ) E u∼δ(u),z ′ ∼RESAMPLE(z,uz) [h(z ′ )] = z∈Z h(z) π θ (z | z \z ) dz.Proof. Proposition 4 in Appendix B provides a lemmaE z∼qη (z t |z t-1 ,ε t z ) E u∼δ(u),z ′ , Ẑ∼RESAMPLE(z,uz) [h(z ′ )] = z∈Z h(z) γ θ (z; z \z ) dz,which we can apply by observing that resampling sums over self-normalized weightsE z∼qη (z|z t-1 ,ε t z ) E u∼δ(u),z ′ ∼RESAMPLE(z,uz) [h(z)] = E z∼qη (z|z t-1 ,ε t z ) E u∼δ(u) E z ′ ∼ uδz (•) u ′ [h(z ′ )] ,which is just a weighted sum that by Definition 4 is itself properly weightedE z∼qη (z|z t-1 ,ε t z ) E u∼δ(u) E z ′ ∼ uδz (•) u ′ [h(z ′ )] = E z∼qη (z|z t-1 ,ε t z ) E u∼δ(u) u u h(z) = E z∼qη (z|z t-1 ,ε t z ) E u∼δ(u) 1 u z∈Z h(z) γ θ (z; x, z \ ) dz = E z∼qη (z|z t-1 ,ε t z ) E u∼δ(u) 1 $ $ $ $ $ Ẑθ (x, z \ ) $ $ $ $ $ Z θ (x, z \ ) z∈Z h(z) π θ (z | x, z \ ) dz = z∈Z h(z) π θ (z | x, z \ ) dz.Proposition 5 (DCPC parameter learning requires only local gradients in a factorized generative model). Consider a graphical model factorized according to Equation 1, with the additional assumption that the model parameters θ ∈ Θ = x∈x Θ x × z∈z Θ z share that factorization. Then the gradient ∇ θ F(θ, q) of DCPC's free energy similarly factorizes into a sum of local particle averages∇ θ F = E q [-∇ θ log p θ (x, z)] = v∈(x,z) E q(v,Pa(v)|εv,ε Pa(v) ) [-∇ θv log p θv (v | Pa(v))] = - v∈(x,z) 1 K K k=1 ∇ θv log p θv (v k | Pa(v) k ).Proof. Proposition 3 provides the lemma that ∇ θ F = E q [-∇ θ log p θ (x, z)], and applying the factorization of the generative model demonstrates that∇ θ F = E q   -∇ θ v∈(x,z) log p θ (v | Pa(v))   .Since the proposal q does not depend on any θ and consists of a particle cloud, we can rewrite it as a mixture over the particles (after sampling is performed)∇ θ F ≈ 1 K K k=1 -∇ θ v∈(x,z) log p θ (v k | Pa(v) k ),and then finally apply the assumption of this theorem that θ ∈ Θ = x∈x Θ x × z∈z Θ z , moving the gradient operation into the sum over individual random variables≈ 1 K K k=1 v∈(x,z) -∇ θ v log p θ v (v k | Pa(v) k ).Contemporaneously to the work of Kuntz et al. [2023] on particle gradient descent, Sun et al.[2023] derived a novel Wasserstein gradient flow and corresponding descent algorithm for discrete distributions. In their setting, each Wasserstein gradient step constructs a D-dimensional, finitely supported distribution over the C-Hamming ball of the starting sample, such that the distribution has DC possible states in total. Let z t+h ∈ N C (z t ) denote the resulting discrete random variable in the C-neighborhood around z t with respect to the Hamming distance. The update rule relies on simulating the gradient flow for time h, sampling from a Markov jump process at time t + hz t+h ∼ d∈[1...D] q(z t+h d | z t d ).A rate matrix Q d (z t ) defined by the entire discrete variable z t parameterizes the proposal distributionEQUATIONthe rate matrix will have nondiagonal entries at indices i ̸= j ∈ [1 . . . C] in the neighborhood N C (z t ), Q d (z t ) i,j = w i,j g π θ (z t \d , z ′ d,j ) π θ (z t \d , z ′ d,i ).The above equation requires that ∀i, j ∈ [1 . . . C], w i,j = w j,i ∈ R and g(a) = ag 1 a . The ratio of normalized target densities π will equal the ratio of unnormalized densities γ πθ (z t \d , z ′ d,j ) π θ (z t \d , z ′ d,i ) = γ θ (z ′ d,j ; z t \d ) $ $ $ $ $ Z z d (z t \d , θ) $ $ $ $ Z z d (z t \d )γ θ (z ′ d,i ; z t \d ) g π θ (z t \d , z ′ d,j ) π θ (z t \d , z ′ d,i ) = g γ θ (z ′ d,j ; z t \d ) γ θ (z ′ d,i ; z t \d ).Based on the experimental recommendations of Sun et al. [2023] , let w i,j = w j,i = 1 and g(a) = √ a. The rate matrix then simplifies to nondiagonal and diagonal termsEQUATIONEquations 17 and 18 give a distribution descending the Wasserstein gradient of the free energy with respect to a particle cloud in a discrete sample space. Applying Equation 18 to γ θ (z; z \z ) yields a factorization in log spaceQ(z t ) i,j = γ θ (z t + i; z t \z ) γ θ (z t + j; z t \z ) log Q(z t ) i,j = 1 2 log γ θ (z t + i; z t \z ) -log γ θ (z t + j; z t \z ) .This difference can be written as a difference of differencesEQUATIONRecent work on efficient sampling for discrete distributions has focused on approximating density ratios, such as the one in Equation 18, with series expansions parameterized by error vectors. When the underlying discrete densities consist of exponentiating a differentiable energy function, as in Grathwohl et al. [2021] , these error vectors have taken the form of gradients and the finite-series expansions have been Taylor series. When they do not, Xiang et al. [2023] showed how they take the form of finite differences and Newton's seriesEQUATIONDiscrete DCPC would therefore use finite differences as discrete prediction errors, breaking each discrete z ∈ z into dimensions and incrementing each dimension separately to construct a vectorEQUATIONwhere ⊖ subtracts the scalar f (z) from the vector elements and f : Z D → R is the target function. This would lead to defining the discrete prediction error as the finite differenceEQUATIONApplying Equation 20to the two terms of Equation 19, we obtain the approximationslog γ θ (z t + i; z t \z ) -log γ θ (z t ; z t \z ) ≈ ∆ 1 log γ θ (z t ; z t \z ) ⊤ • ((z t + i) -z t ) ≈ ε z (z t ) ⊤ • i log γ θ (z t + j; z t \z ) -log γ θ (z t ; z t \z ) ≈ ∆ 1 log γ θ (z t ; z t \z ) ⊤ • ((z t + j) -z t ) ≈ ε z (z t ) ⊤ • j, log Q(z t ) i,j ≈ 1 2 ε z (z t ) ⊤ (i -j) .Discrete DCPC would thus parameterize its discrete proposal (Equation 17) in terms of ε z (Equation 22), so that Equation 18 comes out to the (matrix) exponential of the (elementwise) exponentialq h (z t+h | ε z ) = exp (Q(ε z )h) Q d (ε z ) i,j = exp (ε z ) ⊤ d (i d -j d ) 2 .• The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be. • The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on a few datasets or with a few runs. In general, empirical results often depend on implicit assumptions, which should be articulated. • The authors should reflect on the factors that influence the performance of the approach.For example, a facial recognition algorithm may perform poorly when image resolution is low or images are taken in low lighting. Or a speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon. • The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size. • If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness. • While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an important role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.Kuntz et al. [2023] also interpreted Equation 3 as an update step along the Wasserstein gradient in the space of probability measures. Appendix C extends this perspective to predictive coding of discrete random variables.Such as vision Transformers, denoising diffusion models, etc. 


Decision: The paper presents a novel algorithm, Population Predictive Coding (PPC), which aims to enhance both the performance and biological plausibility of predictive coding by incorporating  Monte Carlo sampling techniques. The reviewers generally agree that the paper makes a meaningful contribution to the field, offering a fresh approach that integrates ideas from machine learning and neuroscience to address the limitations of existing predictive coding methods. While the reviewers raise concerns about the complexity of some explanations, missing results, and the need for clearer documentation and presentation, the authors have effectively responded to these critiques by clarifying their contributions, updating the manuscript (particularly regarding biological plausibility), and providing additional experimental results. Although the overall connection with neural circuitry is intriguing, the paper would benefit from a more direct comparison with neural data. Nevertheless, given the technical soundness, moderate-to-high impact, and the authors' responsiveness to feedback, the paper represents a valuable advancement in understanding predictive coding and structured Bayesian inference.

Reviews: review: paper_summary: paper_summary: The paper introduces a novel algorithm called Population Predictive Coding (PPC)
for structured generative models. The PPC algorithm aims to enhance the 
performance and biological plausibility of predictive coding
approaches by respecting the correlation structure of generative models. The
paper provides theoretical foundations, discusses the biological plausibility of
PPC and performs empirical validation against previous predictive coding and
deep generative modeling algorithms.
strength_weakness: strength_weakness: See text box above.
questions: questions: 1. I might be missing something here, but why is the FID score for PPC missing
   in table 3? On that note, have you considered using `cleanfid` instead of
   `pytorch-fid` to obtain more stable FID score estimates?
2. It is great that you implemented it all in the `Pyro` framework. However, the
   submission seems to be lacking an (anonymous) code attachment. Will you make
   the code publicly available with instructions and all relevant information
   for reproducing your experiments?
limitations: limitations: Yes, the authors have discussed the limitations of their work and potential negative societal impact.
paper_summary: paper_summary: Predictive coding (PC) is a speculative but attractive theory of how certain parts of the brain—especially the so-called 'canonical' cortical microcircuit—might implement Bayesian inference, and hence focus processing effort on the 'surprising' rather than the 'predictable' features of (e.g., sensory) stimuli. PC requires a generative model and an optimization algorithm (i.e., given observations, what is the posterior associated with the model's parameters?). 

While there exists a rich literature related to how PC might be implemented in the brain, previous work has been biologically implausible in a few different senses. In particular, usually (i) a simple (usually Gaussian) generative model is assumed; and/or (ii) inference involves some kind of approximation; and/or (iii) the algorithm is not necessarily fully 'local'. 

The authors propose a novel implementation of PC that they argue is biologically plausible. They call it "population prediction coding" (PPC). They also show that their algorithm can solve some practical machine learning tasks.
strength_weakness: strength_weakness: My main concerns are that the core ideas of the paper are confusingly presented, that the 'biological plausibility' claim may be too strong, and that the model's performance is a bit underwhelming. 

A more minor concern is that figure real-estate is used strangely (why do 80 copies of the same picture take up a huge amount of space in Fig. 4?). It's 'cool', but why is the cortical microcircuit figure (Figure 1) there? It doesn't seem like any of its details are referenced in the main text. Much more helpful would be a diagram of the proposed PPC algorithm, or better yet, two side-by-side diagrams comparing the proposed PPC algorithm to a more 'vanilla' PC algorithm. 

A very minor quibble: "PPC" in neuroscience is already used to refer to "probabilistic population codes" (see Ma et al. 2006) as well as the "posterior parietal cortex", so using PPC here is a bit confusing.

Another relatively minor point. Does vanilla PC really *require* Gaussian densities and the Laplace approximation (Table 1)? My understanding is that the framework itself doesn't *require* this; rather, this is just a useful assumption people use to get something simple and workable. I could be wrong, though. Discussing this point in the paper would be very helpful.

**Confusing presentation.** The meat of the work is in Sections 2 through 4, and I found these sections confusingly written. In Sections 2 and 3, this is probably in part because the authors have to both (a) discuss standard PC in detail, and (b) introduce their modifications to the standard idea. The way things were written made it hard for me to tell what's 'new' and what's not, so it would be extremely helpful if the sections of the text were very obviously separated between 'this is old PC' and 'this is a new modification to PC we are proposing'. As mentioned above, a diagram of vanilla PC vs the proposed modification to it would greatly help. 

In a few places, the motivation for doing things a certain way is also not quite clear. For example, PPC involves using particle-based gradient descent; why? Why this particular approach to empirical Bayes and not some other one? Are there good reasons to believe this to be a good basis for a 'local' solution to empirical Bayes? This does not appear to be discussed in the text. Put differently, is it better to think of this as a *possible* way the brain could implement empirical Bayes, or for some reason a *really good candidate* for that?

**Biological plausibility.** This term is fraught and my view here is one that others may not agree with. But I think the authors are overclaiming here when they discuss biological plausibility. While they acknowledge that it means different things to different people (line 145), I still think this is insufficient. 

The authors really mean *local computations* when they talk about biological plausibility, as they make clear between lines 145 and line 152. But biologically plausible neural computations involve neurons, and hence I think it is reasonable to expect that there is some discussion of whether or not neurons can implement the computations required. Locality is a start, but it is not the end, and in general it can be quite a bit of work to convert the steps of an algorithm into something neurons can plausibly do. There are references suggesting that neurons *might* be able to learn the required log-density gradients (line 175), but this is more of an afterthought than a core part of the contribution of this work, and these thoughts are not well-developed. Without a more specific proposal for how neurons can learn log-density gradients, and some in-silico validation that the proposed strategy works, I really think the term "biological plausibility" should be removed in many places, or else used with extreme caveats included near each use. 

Finally, part of the goal of developing 'biologically plausible' algorithms is to make predictions about biology. What are the specific predictions here? If there are some interesting ones, how do they compare to the predictions of other variants of PC? This should be explicitly discussed.

**Underwhelming performance.** First, I found the tasks used in Sec. 5 confusingly described. It would be helpful if the details of the reconstruction and generation tasks, as well as the training procedure, were better explained in the main text. Second, there are not that many comparisons shown in Table 2; are there other alternatives PPC could be compared to? Third, in Table 3, PPC's FID score is shown as a question mark; I am not sure if this is an error. I assume "undefined" would not be good. 

Finally, there are some issues related to training time and speed. Can plots analyzing these features be included? How do other PC approaches compare on these axes? (Is there some time - performance tradeoff, for example?)
questions: questions: 1. Can the authors briefly summarize the 'new' features of their algorithm relative to other PC ideas?
2. Why particle-based gradient descent as opposed to some other approach to empirical Bayes?
3. Concretely, how might neurons learn log-density gradients?
4. What is the FID score of PPC in Table 3?
5. Is PPC strictly better than other PC algorithms, or is there some kind of tradeoff?
limitations: limitations: The authors discuss some limitations (line 296). I think they could have better discussed limitations related to whether neurons can perform the desired computations, and whether particle-based gradient descent is the only / most appropriate possibility here. I think possible tradeoffs related to (e.g., training) time could also be interesting to discuss.
paper_summary: paper_summary: SETTING: *Biologically plausible* EM with variational inference in directed graphical models.  In particular, the authors aim to implement predictive coding with local parameter updates for learning.


APPROACH: "Variational" EM with sampling (sequential Monte Carlo, SMC) in place of a separate inference/recognition model---a method called "particle gradient descent" (PGD, Kuntz et al., AISTATS 2023).

More precisely, following Neal & Hinton, the authors formulate the objective as minimizing the free energy (i.e. maximizing the ELBo), an upper bound on the marginal cross entropy (the fit of the model to the data).  In place of the recognition model (as in Neal & Hinton, VAEs, etc.), PGD implements Langevin dynamics (LD) in latent space.  More precisely still, a single step of LD is interpreted as providing a (Gaussian) proposal density, for use in importance sampling, with the unnormalized posterior of the generative model as the target distribution.  The importance weights are then used to resample the particles (making this into an SMC algorithm).

To handle structured probabilistic graphical models, this entire procedure is embedded into a Gibbs sampler.  That is, the particle representation for a particular latent variable is computed with all other latent variables fixed; the sampling procedure is repeated for each latent variable in turn; and *this* entire procedure, which constitutes a single step of Gibbs sampling, is repeated S times.  Then the free energy and model parameters are updated once.  *This* entire update procedure is itself repeated T times, with the LD for each particle continuing from where it left off.

The MS aims to produce a biologically plausible algorithm, and identifies its steps with biologically plausible operations:
        
--The Gibbs sampling is proposed to be coordinated with cortical oscillations.
        
--Each cortical column is proposed to encode an individual latent variable, with columnar connectivity reflecting the structure of the directed graph.  Thus columns receive their parents' current (sampled) state via L1 or L6 in order to compute conditional probabilities; locally calculate the energy gradients, and pass it "up" (L2/3 ->) to their parents (-> L4) for accumulation.
        
--Parameter updates can all be performed locally.
        
--Computation of the joint from the local complete conditionals is coordinated by phase-locking to cortical oscillations (at some other frequency?).
strength_weakness: strength_weakness: This MS has two major weaknesses:

(1) The results are very thin/missing.  Indeed, the only results reported in the whole MS are reconstructions (Figs. 2 and 5, Table 2), and these are quantified only for MNIST (and cousins).  They are certainly not close to the performance of modern generative models.  

Several additional results or their descriptions are missing:
        
--Table 3 is missing the FID score!

--There is no Figure 3.

--Figure 4 is not referenced in the text.
  
--Figure 5 refers to "top" and "bottom," but there are four rows of images.  (And it is not clear visually which images are reconstructions of which other images.)

I also expected to see an example with a graphical model that has some actual structure, since this is part of the appeal of the algorithm, but there is none.

(2) It is difficult to determine precisely what aspects of the method are new and not from published work (especially Kuntz 2023, Kuntz 2024, Naesseth 2015).  If the algorithm exists in prior work and the contribution of this MS is supposed to be the mapping to biological circuits, then it does not provide much novelty at all (more on this in the questions below).
questions: questions: --It is difficult to determine precisely what aspects of the method are new and not from published work (especially Kuntz 2023, Kuntz 2024, Naesseth 2015).  Can the authors clarify this?  (I have only read these papers cursorily.)

--Is the mapping to the cortical microcircuit and other neurological constraints supposed to be the (or a) main result of the MS?  These are certainly intriguing but they are not very constraining.  Can the authors make predictions for neuroscience experimentalists based on their mappings?  Can they explain empirical findings in the literature?

--Is it necessary to write the free energy in terms of weights?  Why not just let the recognition model be the product over the normalized "complete conditionals" (with the normalizers computed with Langevin dynamics)?

--The authors propose cortical oscillations as the mechanism for synchronizing a parallel form of Gibbs sampling (that is still guaranteed to have the right stationary distribution).  At what frequency do they hypothesize these oscillations to be?  Given a reasonable number of Gibbs sweeps, does this let the computation happen fast enough?

More generally, can the authors give some more detail on the time complexity of the algorithm and how they would expect it to scale to problems with (say) real graph structure?

--Proposition 5 is repeated as Corollary 5.1, which then references Proposition 5 (should be 3).

--"PPC" is already the name of a (very popular) proposed representation of probability distributions in populations of neurons ("probabilistic population codes," Ma et al., Nature Neuroscience, 2006).  Perhaps the authors can invent another name.
limitations: limitations: N/A
paper_summary: paper_summary: This paper proposes a version of predictive coding that is biologically plausible (in the sense that all computations are local), and can perform both classification and generation (by operating in both directions). This version of predictive coding, dubbed Population Predictive Coding (PPC), overcomes the limitations of traditional PC by enabling a more general class of distributions in the latent layers.
strength_weakness: strength_weakness: Some notation needs to be defined to help the reader. For example:
- eqn (1): the function Pa is not defined.
- eqn (5): Ch(z) is not defined. I infer from Algorithm 1 that it probably means "children of", but it should be made clear when it is first used.


Additionally, there are a few things that could be fixed.
- There is no figure 3.
- figure 4 has 4 rows, but the caption only states what is in the "Top" and "Bottom". Please clarify.
questions: questions: I would have liked to see the classification accuracy of the model on MNIST (as a percent). Predictive coding networks tend not to perform well on classification tasks.
limitations: limitations: Yes

Comments: title: This paper studies generative modeling.
Thank you
comment: Ok, good, we're on the same page now.  To me this formulation is simpler, but I defer to the authors' judgment here.
Dear reviewers,

Thank you again for your time and comments, and for having been amenable to discussion during the rebuttal period. This message is to communicate a nice improvement in our experimental evaluation: in the last days, we have also run an apples-to-apples comparison against the Langevin PC work [1], using the same model, and, more importantly, the discretized Gaussian function as a likelihood function, that they do. The results show that the discretized Gaussian largely improves the performance in our case as well.
In terms of results, the authors report three different numbers for three different methods: prior sampling, and two amortized warm-start techniques, the best one using Jeffrey’s Divergence. However, we decided to not implement warm-start techniques, as they are not biologically plausible in the sense that we mean in our work, and hence only perform prior sampling. Despite that, D&C-PC achieves a FID score of 96.0473 (+/- 0.2667), that is (1) much better than the one of 120 they report for prior sampling, in the apples-to-apples comparison; and (2) better than both of the method that use warm starts, where they report a FID score on the test set of 97.49 (Table 1 in [1]).

We can thus claim that D&C PC is the variation of PC that achieves the best results in the literature when it comes to image generation tasks, even when compared against less bio-plausible PC methods. We have updated the manuscript accordingly.
A second, minor, change that we have applied to the manuscript, is to write the free energy in terms of a product of complete conditionals as the recognition model Q and then consider estimating it by Monte Carlo with the whole Langevin to nested SMC procedure, as suggested by reviewer DHBv. We agree that this would further improve the clarity of our work.

We once again deeply thank the reviewers for their close engagement that has done so much to improve the paper.

[1] Zahid, U., Guo, Q. and Fountas, Z., Sample as you Infer: Predictive Coding with Langevin Dynamics. In Forty-first International Conference on Machine Learning, ICML 2024.
Aha, we seem to have misunderstood what you were saying.  At least this author had always thought of a free-energy as definitionally involving a proposal/recognition model q that admitted closed-form sampling, providing a base-case on which the more elaborate sampling algorithms are written.  From the point of view of how you're putting it, yes, we could write the free energy in terms of a product of complete conditionals as the recognition model Q and then consider *estimating* it by Monte Carlo with the whole Langevin to nested SMC procedure.  We agree with the reviewer that from the point of view of optimization rather than sampling algorithms, that would likely be much clearer, and will adopt that into the manuscript.
The authors have filled in the missing results, and distinguished their contribution from other papers, so I have raised my score.

> It is indeed mathematically necessary to write the free energy in terms of importance weights targeting the generative joint density. By definition, a variational free energy is the expected value of the negative logarithm of an importance weight. Equivalently, the VFE is the cross-entropy of the generative joint distribution, taken with respect to the proposal distribution (recognition model), minus the entropy of the recognition model. The product of normalized complete conditionals (with the normalizers estimated by importance sampling) would only estimate the entropy of the recognition model, the second term.

I was not proposing to drop either term from the free energy, but rather to rewrite it without any reference to "weights" or "strictly properly weighting."  That is, precisely to rewrite it as "the cross-entropy of the generative joint distribution, taken with respect to the proposal distribution (recognition model), minus the entropy of the recognition model," where the recognition model is defined to be the product over the normalized "complete conditionals" (with the normalizers computed with Langevin dynamics).  This would (in this reviewer's opinion) substantially simplify the presentation.  But perhaps there is some obstacle to writing it this way that I am overlooking.
We thank the reviewer for their feedback and the significant improvements it has aided in our manuscript.  We are finely combing the paper for similar small typos.
I thank the authors for their helpful responses, and think they have mostly addressed my concerns. Overall, I think they make a worthy theoretical contribution to how the brain might implement inference over graphical models. I have increased my score.

As a minor point, I encourage the authors to look over the paper for small typos (e.g., "hypotesis", line 25).
Thank you for your quick response.

This paper studies generative modeling, in which predictions flow top-down from a latent (initialized from the prior and then updated towards the posterior by predictive coding) toward the data (images).

In greater detail, this paper studied predictive coding as an inference strategy in Bayesian generative modeling, not PCNs as a neural network architecture.  That is why we have not conducted a classification experiment.
I have read your responses.

In the paper cited above, the PC networks achieved high classification accuracy (approx 98%), but their architectures had the predictive connections projection from the image inputs toward the classification layer, similar to [Whittington, Bogacz; 2017]. Is that the same as the models in this paper? Which directions do the predictions flow? From input toward latent (classification vector), or from latent toward the input (images)?