Expectation propagation, usually shortened to EP, is a method for approximating a difficult Bayesian posterior by repeatedly replacing one awkward likelihood factor with a tractable stand-in. At each update it removes that factor’s old approximation, combines the remainder with the original exact factor, matches selected moments such as mean and variance, and inserts a revised approximation. The remove–recalculate–replace cycle is what distinguishes EP from simply fitting one Gaussian once.
The useful object is a cavity distribution: the current approximate posterior with the chosen approximate site divided away. Multiplying that cavity by the original factor gives a tilted distribution. Matching the tilted distribution’s moments constructs a new tractable approximation; dividing it by the cavity determines the new site. This can support approximate inference when exact integration is costly, but it is not a guarantee of convergence, correct tail probabilities, calibrated uncertainty or an exact posterior.
This guide explains the mechanism from first principles, walks through an exact Gaussian example and an approximate half-normal example with checked numbers, shows how site removal prevents double counting, diagnoses unstable iterations, and gives an independent changed-context calculation. The examples are theoretical calculations with invented observations and explicitly named units; they are not recorded results from a deployed eduKate machine-learning system or experimental benchmark.
Contents and useful reading routes
Begin with the inference problem, factorised sites, the five-step algorithm, moment projection and Gaussian natural parameters.
For arithmetic, use the exact Gaussian control, site removal, the half-normal factor, the meaning of a site, the changed cavity and the independent transfer task.
For engineering, read the multi-site issue, update schedules, damping, invalid cavities, computational cost, implementation tests and failure repair.
For comparisons and boundaries, visit multimodal targets, variational Bayes, Laplace, Monte Carlo, assumed-density filtering, evidence estimates, limits, FAQs, original research and more reading.
Why a Bayesian posterior becomes hard to calculate
Bayesian inference combines prior beliefs about an unknown variable with a likelihood describing observed data. For a latent parameter x and observations y, the posterior is proportional to p(x) multiplied by p(y given x). The normalising constant is the integral of that product across possible x values. In simple conjugate models this is manageable; in nonlinear or nonconjugate models the integral and many requested expectations can become expensive or unavailable in closed form.
Suppose a model contains thousands of observations, each introducing its own likelihood factor. The unnormalised posterior can be written as a prior p0(x) times factors f1(x), f2(x) and so on. A prediction may require integrating this product against a future-outcome model. It is not sufficient to evaluate the product at one attractive value of x when the task needs uncertainty as well as a point estimate.
Several approximation families make different trade-offs. Sampling attempts to represent uncertainty with draws, optimisation may locate a mode, and parametric approximations represent the full posterior using a smaller class of distributions. EP belongs to this last family: it simplifies each factor while trying to preserve selected consequences of the product. That creates an algorithmic loop rather than a one-shot analytical trick.
The goal should be named before choosing a method. If a decision needs a reliable tail probability, a Gaussian approximation with excellent mean and variance may be inadequate. If a model only needs an approximate posterior mean for a smooth prediction, the same approximation may be useful. EP is a procedure with assumptions and diagnostics, not a magical declaration that Bayesian integration has been solved.
Factor the problem before choosing an approximation
Write the unnormalised target as p0(x) times the product over i of f_i(x). Each f_i might represent a likelihood contribution from one observation, a group of records, or a logical constraint. The factorisation tells EP where it can perform a local update. It does not say that the observations are physically independent unless the generative model actually justifies that structure.
The EP approximation replaces each exact factor f_i with a tractable site t_i from a chosen family, often an unnormalised exponential quadratic in x for Gaussian EP. The approximate posterior is q(x) proportional to p0(x) times the product of all these t_i. Each site represents a computational compromise, not a new independent observation. Treating the site as extra evidence on top of its original factor would count the same information twice.
A single global Gaussian q can have a dense covariance matrix, while its individual factors may be local. The expressiveness and cost depend on whether the approximation retains full covariance, diagonal covariance or another structure. A full covariance in d dimensions has on the order of d squared entries; a diagonal one has d. Neither structure wins universally: correlations can matter, and the computational budget is finite.
Before implementing EP, map each factor to an identified portion of the data-generating model and check whether that model captures the decision you care about. An algorithm can perform numerically stable moment matching on a statistically incorrect factorisation. That would be an implementation success but an inference failure.
The exact five-step site update
First select the site index i whose approximation will be revised. Second compute its cavity, q-minus-i(x), proportional to the current q(x) divided by its existing site t_i(x). Third form the tilted distribution, p-tilted-i(x), proportional to f_i(x) times that cavity. At this point the exact original factor is used locally; it is not replaced by its old approximation.
Fourth project the tilted distribution back into the chosen tractable family. If that family is Gaussian, choose a Gaussian with the tilted distribution’s mean and variance, or mean vector and covariance. Fifth solve for the new site by dividing that projected Gaussian by the cavity. The result can be stored as natural parameters, which make products and quotients easier to compute.
Repeat site updates until an appropriate numerical stopping condition is met or an iteration budget is exhausted. Some implementations update sites sequentially; others compute proposals in parallel. Changing the scheduling may change convergence behaviour. Damping can blend an old and proposed site’s natural parameters, but it modifies the iteration path and does not grant a theorem of convergence for arbitrary models.
The cavity operation is indispensable. If you merely multiply the current q by f_i on every visit, the same observation is repeatedly added as though it were new. Posterior precision can then become overconfident for reasons unrelated to fresh data. EP instead tries to replace the previous contribution from that factor, making repeated refinement possible without naive repeated counting.
What moment matching actually preserves
A probability distribution has infinitely many possible properties. Its mean and variance capture location and spread, but not all aspects of shape. When EP projects a tilted density into a Gaussian family by matching mean and variance, the projected density can differ in skewness, tails, number of modes and support. This is an intentional information loss, not a formatting error.
Moment matching is a mathematical property of a particular projection. When the necessary integrals exist, the Gaussian that minimises KL(tilted distribution || candidate Gaussian) has the same mean and covariance as the tilted distribution. That local direction of Kullback–Leibler divergence matters. It does not imply that the full EP fixed point globally minimises KL(true full posterior || q) or that every iteration moves closer to the true target.
Matching must be performed on the normalised tilted distribution, even if computation uses unnormalised factors. If the normalising integral is Z_i, the tilted mean is the integral of x times cavity times f_i divided by Z_i. Omitting Z_i can produce a value on the wrong scale. For variance, the second moment must have its mean squared subtracted correctly.
Analytic moments are available for some factors, and numerical quadrature or other approximations may be required for others. If a purported implementation returns moments of the old site rather than the exact factor multiplied by the cavity, it has stopped implementing the method described here. A good test stores the intended tilted integral, normaliser and resulting moments with clear tolerances.
Why Gaussian natural parameters make EP practical
For a scalar Gaussian with mean mu and variance sigma-squared, the precision is tau = 1 divided by sigma-squared, and the natural linear parameter is nu = mu divided by sigma-squared. A Gaussian density is proportional to exp(nu times x minus one half tau times x squared), up to a normalisation constant. In a product of Gaussian-shaped factors, the tau values add and the nu values add.
Division is equally simple: removing a site’s Gaussian-shaped term subtracts that site’s tau and nu from q’s totals. This is the algebraic advantage that makes a cavity computation manageable. With covariance matrices the precision becomes a matrix and the linear term becomes a vector; positive definiteness and numerical conditioning must then be checked.
Notice that the site is not required to be a normalised Gaussian probability density on its own. It may have a scale constant that affects approximate evidence calculations, and some EP updates produce site parameters that are not themselves positive, normalisable Gaussian densities. The entire approximate posterior must still be meaningful. An implementer who treats every site as a standalone probability distribution may reject otherwise usable representations or silently distort calculations.
Natural parameters also create a useful diagnostic. If q precision is 1.25 and the current site’s precision contribution is 1, the cavity precision is 0.25. If a programmer accidentally adds rather than subtracts, the resulting cavity precision becomes 2.25: a very different distribution that wrongly retains the same evidence. Simple arithmetic unit tests can detect that mistake before a large model is deployed.
Worked example: an exact Gaussian update as a control
Let an unknown one-dimensional measurement x have prior N(1, 4), meaning mean 1 measurement unit and variance 4 square units. Suppose one invented observation reports a noisy value y = 2 measurement units with known Gaussian observation-noise variance 1 square unit. The likelihood as a function of x is proportional to exp(-(x minus 2) squared divided by 2). This conjugate example has an exact posterior, providing a check against an EP implementation.
The prior precision is 1/4 = 0.25 per square unit and its natural linear parameter is 1/4 = 0.25 per unit. The likelihood contributes precision 1 and linear parameter 2. Adding gives posterior precision 1.25 and linear parameter 2.25. The posterior variance is 1/1.25 = 0.8 square units and the posterior mean is 2.25/1.25 = 1.8 units.
If EP begins with a neutral site, the cavity is the prior N(1, 4). Multiplying it by the exact Gaussian likelihood gives a Gaussian tilted distribution N(1.8, 0.8). Moment matching returns that same Gaussian. Dividing by the cavity leaves the likelihood’s original natural contributions, site precision 1 and linear parameter 2. This update is exact because the chosen approximation family already contains the product.
Check the denominator and units: precision is inverse squared measurement units, the posterior mean is measured in ordinary units, and posterior variance is measured in squared units. Comparing the numbers 0.8 and 1.8 without saying which is variance and which is mean would be a category error. The example is not evidence that EP remains exact with nonlinear likelihoods.
Remove a site and recover exactly the previous information
Start from the exact-control posterior N(1.8, 0.8), with precision 1.25 and natural linear parameter 2.25. Remove the likelihood site’s precision 1 and linear parameter 2. The cavity has precision 0.25 and linear parameter 0.25, hence variance 4 and mean 1. It is exactly the original prior, as required because no other factor contributed.
This subtraction provides a sharp implementation test. If a library’s cavity mean is 1.8 after removing the only site, it has not removed that site’s effect. If cavity variance is 0.8, the old likelihood is still present. Either result would cause an update to use evidence that was supposed to be held out from that site’s own refresh.
In a multi-site model the cavity is not generally the prior. It contains the prior and every other site’s approximated information. The correct question is therefore not what the posterior looked like at the beginning of all training, but what remains after removing only the specific site’s current contribution. Confusing these two reference states is a common source of incorrect EP updates.
Treat the cavity as a numerical object subject to validation. In a Gaussian family, a cavity with nonpositive or invalid precision cannot be passed to ordinary Gaussian moment formulas without care. Some EP formulations and damping strategies handle difficult intermediate configurations differently; no robust implementation should silently take the square root of a negative variance.
A non-Gaussian factor: observing that x is positive
Now use a cavity N(0, 1) on a dimensionless scalar x. The exact factor f(x) equals one when x is positive and zero otherwise. It represents a hard inequality observation, not a Gaussian noisy measurement. Multiplying by the standard normal cavity produces the positive half of the normal distribution, renormalised so its total probability is one. That tilted distribution has zero support on the negative half-line.
The half-normal mean is square root of (2 divided by pi), approximately 0.7978845608. Its variance is 1 minus 2 divided by pi, approximately 0.3633802276. These analytic moments use a normalised positive half-normal distribution. A Gaussian family cannot preserve the hard zero probability for negative x, but EP can choose a Gaussian q-new with that mean and variance.
The Gaussian projection’s precision is 1 divided by 0.3633802276, approximately 2.751938394. Its linear natural parameter is 0.7978845608 divided by 0.3633802276, approximately 2.195729157. The cavity N(0, 1) has precision 1 and linear parameter 0. Dividing projected q by the cavity leaves a site with precision approximately 1.751938394 and linear parameter approximately 2.195729157.
This is genuinely approximate inference: the original factor was a step function and the site is a Gaussian-shaped approximation. The projected q assigns nonzero probability to negative x even though the exact tilted distribution forbids it. Matching two moments did not preserve support. If the downstream task is the probability that x is negative, this particular Gaussian replacement answers the wrong qualitative question unless a separate correction is used.
Do not confuse a Gaussian site with an observed datum
The positive-threshold site’s natural linear parameter is approximately 2.19573, but no sensor measured x = 2.19573. It is the coefficient needed so the product of this site and the cavity has the chosen moments. Similarly, the site precision is not a physical sensor precision derived from a noise calibration. It is an algebraic element of an approximation to the indicator factor.
Interpreting site parameters as original observations is a form of reverse-engineering error. The true information is x greater than zero. The moment-matched Gaussian acts like a smooth potential that pushes probability in a positive direction. That does not change the provenance of the observation. A transparent report should record the original likelihood or constraint and mark the approximate site as a computational substitute.
If the cavity later changes, the site parameters should generally be reconsidered. A once-correct local approximation was calculated against a particular distribution of the other factors. Adding another strong constraint can shift the region where the step factor matters most, so EP re-forms its cavity and tilted distribution on later visits. A frozen site may then be less appropriate.
Use this distinction when explaining predictions to colleagues. State what was actually observed, which posterior family was assumed, which moments were preserved and which questions remain unreliable. A smoother density plot may be easier to draw than a hard boundary, but aesthetic convenience is not statistical correctness.
Changed-context calculation: widen the cavity without changing the factor
Keep the same hard threshold x greater than zero but replace the cavity N(0, 1) with N(0, 4). The standard deviation doubles from 1 to 2 units, while the distribution remains symmetric about zero. The tilted distribution is now a half-normal with scale 2. Its mean is 2 times square root of (2/pi), approximately 1.595769122 units.
The tilted variance is 4 times (1 minus 2/pi), approximately 1.453520911 squared units. The projected Gaussian precision is about 0.6879845985 inverse squared units. The cavity precision was only 1/4 = 0.25, so the revised site’s precision contribution is approximately 0.4379845985. The projected linear parameter is mean divided by variance, approximately 1.097864578.
Although the exact factor x greater than zero did not change, the appropriate site approximation changed because the cavity’s spread changed. This is the central iterative idea: the approximation to one factor depends on what all other factors currently say about x. A fixed approximation obtained once in isolation does not capture that adaptation.
Double-check the scale: the mean doubles, variance quadruples, and precision becomes one quarter as large when the entire symmetric cavity’s scale doubles before truncation. Those transformations are consistent with dimensional analysis. This exercise also shows why a comparison of two site precision numbers without the cavity context may be misleading.
One analytic update is not an iterative convergence demonstration
Our first Gaussian and threshold examples each involve one likelihood factor and one cavity. They demonstrate what a single site projection looks like and how natural parameters are recovered. They do not show the challenging feature of EP in a real multi-factor model: a later update to one site changes the cavity used when another site is revisited.
Suppose a second factor rewards values of x near some different region. EP may first approximate the positivity constraint against the prior, then approximate the second factor using a cavity containing the first site’s approximation. When the second site changes q, removing the first site now yields a different cavity, so the first threshold factor’s tilted moments may also change.
The algorithm repeatedly revises approximate sites until a stopping rule is met. Such a rule may check changes in moments, natural parameters, an approximate energy or held-out downstream quantities. Small changes in a site’s parameters are not proof that the true posterior is being represented accurately, and an iteration can exhibit oscillation or fail to approach a valid fixed point.
Do not present one exact conjugate calculation as proof of an algorithm’s performance on nonlinear models. A sound implementation suite should include tractable controls, nonconjugate toy cases, comparisons to dependable high-accuracy reference calculations and diagnostics for cases where the method becomes unstable.
Sequential and parallel updates change the computation
Sequential EP replaces one site, updates q, and then uses that newer q to construct the next cavity. In parallel EP, several site proposals may be based on the same older global approximation before updates are committed. The former carries information forward immediately; the latter can support distributed computation but may introduce stale proposal interactions.
Neither schedule is universally superior. If factors strongly oppose each other, changing many sites at once can produce a large jump in global precision or mean. Sequential updates can also depend on ordering and can oscillate. Useful diagnostics record the schedule, the update order, proposal changes and any damping, not merely the number of passes labelled epochs.
Consider a database split across machines. Each worker sees its factor block and a cavity reflecting the current global approximation minus its own site. Sending back a proposed update is conceptually different from combining independently computed complete posteriors: EP must avoid counting the prior and shared factors multiple times. This distinction matters for distributed Bayesian analysis as well as for local single-machine algorithms.
Latency and bandwidth measurements are separate from inferential accuracy. A parallel implementation can be faster in wall-clock time while converging to a poor approximation, or slower overall if it needs many damping steps. Always report what was measured and the reference standard used to judge statistical adequacy.
What damping changes and what it cannot promise
One practical stabilisation method is damping: instead of accepting a full proposed site update, move partway from the old site’s natural parameters toward the proposal. With damping weight rho between zero and one, a common rule is new-site-parameter = (1 minus rho) times old-site-parameter plus rho times proposed-parameter. Choosing rho = 0.2 takes one fifth of the proposed change.
Take a fictional site’s old precision 0.4, proposed precision 1.4 and rho = 0.2. The damped precision becomes 0.8 times 0.4 plus 0.2 times 1.4 = 0.60. For an old linear parameter 0.5 and proposed 1.5, the damped result is 0.8 times 0.5 plus 0.2 times 1.5 = 0.70. These are changes in natural parameters, not arbitrary blending of raw variances.
Mixing means and variances directly would usually produce a different distribution from mixing their natural parameters. That is why the parameterisation belongs in the method specification. Damping can reduce abrupt changes and improve numerical behaviour in some models, but a tiny update does not establish that the fixed point is unique, desirable or even present.
A robust experiment tries several schedules and damping values on the same data, records whether each reaches stable valid moments and checks predictions against an independent standard. If different initialisations yield very different answers, a practitioner needs to surface that fact rather than silently presenting one run as the posterior.
Negative site precision and invalid cavities are not the same
Because sites are approximation factors rather than standalone normalised Gaussians, an EP site may sometimes have a negative quadratic coefficient relative to the chosen cavity. That fact alone does not necessarily make the global q invalid. What matters for a scalar Gaussian q is that its total precision remains positive; for a multivariate Gaussian it is that the global precision matrix be positive definite.
Removing a site can still produce an invalid cavity in difficult iterations. If the resulting scalar precision is zero or negative, the usual normalised Gaussian cavity formulas do not apply. In higher dimensions, a matrix might have an eigenvalue that prevents a proper covariance. An implementation that silently clips every value to a positive constant changes the algorithm and can hide the source of instability.
Treat non-finite moments, extreme natural parameters and impossible normalisers as first-class errors. Retain the original factors and the last verified valid iterate when a computation fails, and investigate the factor responsible. Methods using damping, fractional EP or alternative projections may help in some cases, but each must be described accurately rather than called a universally safe fix.
For numerical work, prefer stable linear algebra routines to explicit matrix inversion where appropriate. Symmetrise covariance estimates before diagnostic eigenvalue checks only with clearly documented tolerances, and record the condition number when it matters. A posterior that appears normal on a plot can still have numerically unstable precision calculations.
When a single Gaussian cannot represent the answer
Suppose the target posterior has two well-separated modes, one near negative three and one near positive three. Its overall mean might be near zero, where the true density is very low. A moment-matched single Gaussian centred at zero can suggest a plausible value where the model actually assigns little probability. This is a limitation of the approximation family, not an error in algebra.
EP does not automatically repair inadequate expressiveness by iterating longer. Revisiting factors can improve local moment consistency, but if every site and the global q must belong to an overly narrow family, some structural information remains unrepresentable. More expressive mixture or structured approximations may help at extra computational cost.
Hard constraints create another failure. In the half-normal example, exact probability below zero was zero, but the Gaussian projection spread nonzero density into the forbidden region. A prediction requiring that constraint should use the original truncated distribution, a correction, or an inference method that preserves support. A mean estimate may be useful while a safety-critical tail probability is unacceptable.
Examine how the approximation will be used. A posterior for exploratory ranking may tolerate small shape errors; a medical or financial threshold may depend strongly on tails and calibration. No claim that EP is accurate in an average benchmark replaces evaluation of the specific decision-relevant quantity.
EP and variational Bayes minimise different local objectives
Both expectation propagation and variational inference can represent difficult posteriors with tractable families, but they generally use different objectives and update procedures. In common variational Bayes formulations, one minimises KL(q || p) by maximising an evidence lower bound. Local EP moment projection corresponds to matching expectations of a tilted distribution, commonly expressed through KL(tilted || projected q).
KL is asymmetric. In a multimodal example, one direction of divergence can favour concentrating on one mode while another can favour a broader fit spanning several modes, depending on the approximation family and optimisation. This is a useful intuition, not a guarantee that every EP run is mode-covering or every variational run is mode-seeking. Real outcomes depend on structure, parameterisation and fixed points.
EP’s global iteration does not simply perform global gradient ascent on the standard variational evidence lower bound. Its free-energy and consistency interpretations are more nuanced, and convergence can be difficult. The original work by Minka and later expectation-consistent research by Opper and Winther provide more formal accounts.
When comparing methods, fix the same target and test the same posterior quantities. If one algorithm reports only its optimised objective while another reports a predictive log score, the two outputs cannot be compared as though they were the same measure of posterior accuracy. A useful analysis includes computation, stability and uncertainty validity.
EP and the Laplace approximation serve different roles
The Laplace approximation typically expands the log posterior around a mode and uses local curvature to create a Gaussian approximation. It asks what a quadratic approximation near that selected mode looks like. EP instead repeatedly builds factor-wise approximations intended to match moments of locally tilted distributions. Both may return Gaussian summaries, but the route to those summaries differs.
If the posterior is approximately Gaussian and unimodal, both can work well. If it is skewed, a local curvature near the peak need not describe the mean and variance accurately. A moment-matched approximation may capture those two statistics better in some cases, yet neither method can fully represent a hard constraint with a plain full-support Gaussian.
Compute budgets also differ. Laplace may require a mode search and a Hessian or Hessian approximation, while EP may require repeated moment calculations and site management. In high dimensions either can be expensive or challenging. A headline saying that one always beats the other would ignore model class, integration difficulty and hardware constraints.
An honest method comparison asks which quantity the user needs: MAP-like location, posterior mean, a credible interval, a classification probability or an approximate marginal likelihood. It then measures that quantity against a trustworthy reference, rather than treating the common word Gaussian as proof of equivalence.
EP is not Monte Carlo sampling in disguise
Monte Carlo methods approximate integrals using draws; under suitable conditions and sufficient computation, some samplers provide asymptotically exact targets. EP instead builds and repeatedly adjusts a structured parametric density. It may be computationally attractive because it summarises a large amount of information in a limited set of moments, but its approximation error need not vanish simply by performing more iterations.
A sampler can also fail in practice through poor mixing, convergence diagnostics or insufficient exploration. EP can fail through inaccurate moments, unstable sites, a mismatched family or a bad fixed point. The failure modes are different, so a robust study compares both under realistic computation budgets and relevant outcome metrics.
Suppose a posterior is a bounded positive density. An ideal sampler stays within the target’s actual support because its transition and acceptance mechanisms use the target. A single Gaussian EP projection may assign nonzero mass outside the support. That qualitative difference can matter more than an apparently small error in the mean. Choose the evaluation criteria accordingly.
Do not infer superiority from elapsed time alone. If EP returns a mean in one second while sampling takes sixty seconds, the gain is meaningful only if the mean and any required uncertainty are adequate for the job. If the decision turns on a rare-event probability that EP cannot represent reliably, the faster method may be the wrong one.
EP versus assumed-density filtering
Assumed-density filtering processes factors in sequence and projects the updated distribution back into the tractable family after each new factor. Once a factor has been processed, its approximation is typically not revisited. EP extends that intuition by keeping site approximations that can be removed and recalculated when other factors change the cavity.
That difference is why EP’s site bookkeeping matters. In one-pass filtering, the approximation error made when an early factor arrived may persist. In EP, a later change in the global approximation can motivate refining the early site. This may improve consistency, but it also introduces multiple passes, update ordering and potential oscillation.
A streamed data pipeline may prefer an online approximation even when full EP would give a different result after several revisits. If old raw observations cannot be retained or revisited, the exact tilted update for those old factors may be unavailable. A practitioner must then describe the algorithm actually implemented rather than call an irreversible one-pass update repeated EP.
Factor grouping is also a design choice. One site per observation, one site per batch and one site per machine do not create identical approximations or costs. Decisions about grouping should be tested empirically and recorded because they shape which moments are locally preserved.
Posterior moments and marginal likelihood are different outputs
The approximate posterior q answers questions about latent variables conditional on data. The marginal likelihood, or model evidence, integrates the prior times all original likelihood factors over the latent variables. In a Bayesian model comparison that integral can be central. It is not determined merely by writing down q’s mean and covariance.
EP formulations track site normalisation information and have approximate evidence or free-energy calculations, but these require their own precise equations and bookkeeping. Dropping every site’s multiplicative constant may leave some moments unchanged while corrupting the approximate evidence. An implementation that claims a marginal likelihood should document how those constants were handled.
Evidence estimates are not automatically valid lower bounds, unlike the standard variational lower bound interpretation under its assumptions. A numerical value can be accurate for a well-behaved model and misleading for another. Cross-check on tractable cases and do not advertise an EP evidence number as a rigorous bound unless the selected formulation establishes one.
If a product needs only predicted means, evidence tracking might be unnecessary overhead. If the goal is comparing different model structures, omitting normalisers could invalidate the intended decision. State the objective at the outset so that the engineering implementation preserves the quantities the decision actually requires.
Computational cost is more than the number of sites
For a scalar Gaussian EP model, multiplying and dividing site natural parameters is cheap. A high-dimensional full-covariance Gaussian can require costly matrix factorisation, and evaluating tilted moments may dominate runtime. The bottleneck depends on the chosen site family and the structure of the factors, not merely on the label EP.
Diagonal approximations reduce storage and often computation but discard correlations. Low-rank or sparse structures can preserve selected dependencies, at the cost of more complex updates. For many models, the main difficulty is evaluating the moments of a nonlinear tilted distribution accurately and repeatedly. A fast natural-parameter update does not rescue a poorly approximated integral.
A realistic performance record includes the number of sites, latent dimension, matrix representation, factor evaluation cost, integration method, update schedule, damping, convergence criterion and hardware. If one implementation performs a single pass and another performs twenty, comparing raw milliseconds without these details will mislead readers.
Choose a small tractable reference model as an initial correctness test, then scale gradually. Log intermediate cavity precision and tilted-moment quality. When the approximation breaks at a specific factor, inspecting that local problem is often more informative than merely increasing the global iteration cap.
A disciplined implementation and verification plan
Begin with the exact Gaussian control from this guide. Your code should recover posterior mean 1.8 units and variance 0.8 squared units and, after removing its site, the prior mean 1 and variance 4. Add the half-normal threshold case to confirm that moment matching produces the stated positive mean and a projected Gaussian with nonzero negative support.
Next add two factors and deliberately update them in different orders. Compare final approximate moments and stability, then turn on damping with a documented coefficient. Test failure fixtures with zero or invalid cavity precision, non-finite moments and an inconsistent factor that makes numerical integration fail. The program should stop or report a precise failure rather than silently output a plausible-looking posterior.
Keep all input likelihoods, their provenance, approximation choices and versioned code paths separate. If a site update is retried after an error, reconcile the current site state before applying it again. Otherwise a partial update could be accidentally applied twice. For distributed systems, use versioned site ownership or another appropriate consistency control, but do not confuse concurrency engineering with a new theorem about Bayesian inference.
Assess outputs against the user job. For regression, check predictive error and intervals on held-out relevant cases; for classification, examine calibration and subgroup behaviour as appropriate; for scientific inference, compare means, uncertainties and model assumptions. An approximation that passes internal numerical tests can still fail external validity.
Why Gaussian-process models sometimes use EP
A Gaussian process specifies a jointly Gaussian prior over latent function values. With a Gaussian observation model, posterior calculations retain a conjugate Gaussian form. With classification likelihoods such as probit or logistic factors, the posterior over latent function values is generally no longer exactly Gaussian. That is a setting in which approximations such as EP can be considered.
Each classification observation contributes a non-Gaussian likelihood factor. Gaussian EP approximates those factors with site terms so the global latent posterior remains tractable in a Gaussian form. The cavity for a selected observation reflects the other observations and the prior. The tilted distribution then uses the original classification likelihood for that observation, making this an appropriate example of factor-wise replacement.
Numerical implementation of Gaussian-process classification also requires stable matrix calculations, kernel choices and checks for approximated predictive probabilities. The fact that the latent approximation is Gaussian does not mean the final classification probability is a Gaussian random variable. Predictions involve combining the latent distribution with a link function, and calibration still needs empirical evaluation.
For the separate question of kernel-based regression and exact Gaussian observation updates, see eduKate’s already published Gaussian Process Regression. That article owns a different reader job from the EP algorithm’s cavity and moment projection machinery.
A tidy approximation does not certify model fairness
EP makes an inference calculation more manageable under selected assumptions. It does not validate the input data, decide whether a training set represents the deployment population, prove that the likelihood is correct or guarantee that a sensitive decision is fair. A mathematically converged approximate posterior can be confidently centred on the wrong assumptions.
For a classification application, calibration checks whether predicted probabilities align with observed frequencies on a specified evaluation population. Fairness questions may require separate subgroup or conditional analyses and careful decisions about what outcomes and denominators matter. Approximation error can interact with these issues, but computing a moment-matched site does not automatically resolve them.
Protect the distinction between numerical convergence, posterior fidelity and decision quality. A well-behaved EP update can converge yet be too approximate in a tail that controls an operational threshold. A more accurate posterior can still lead to a poor decision if the objective or intervention model is misspecified. Record each layer of the problem independently.
External review is particularly important for high-stakes uses. Compare estimates with alternative inference methods, validate the real-world task, check uncertainty in conditions different from training and avoid claiming certainty beyond the model. This is responsible technical practice, not an argument that EP is inherently unsafe.
Diagnose five characteristic EP failure modes
If results change wildly with update order, investigate oscillation, badly scaled factors and damping. Run reproducible schedules and compare local moment estimates to trustworthy numerical references. If the global covariance becomes invalid, inspect cavity subtraction, sign conventions and the precision matrix before altering tolerances. Returning a variance of negative one is not a slightly inaccurate posterior.
If the projected mean seems plausible but a hard constraint is violated, inspect the approximation family. The half-normal case already demonstrates how a Gaussian can place probability outside support. That problem is not corrected by a better optimisation stopping criterion. Change the representation, use constraint-aware corrections or select a method that retains the necessary support.
If posterior variance becomes suspiciously small after many passes, look for double counting of one site. Verify that the cavity subtracts the old site’s precision and linear term before reinserting the projected replacement. Test with the exact Gaussian control: removing the only site must recover the prior. If it does not, the implementation is algebraically wrong.
If predictions are poor despite apparently stable moments, examine the generative model, data distribution and the moments being matched. Perhaps the target is multimodal, the factor likelihood is wrong or a downstream decision depends on tails. Replace the source of failure rather than simply increase the iteration count.
Finally, if results differ significantly under alternative initialisations, preserve those results and report that the approximation is sensitive. Do not hide disagreement by choosing only a visually attractive trajectory. EP can be a useful method without being a universally convergent black box.
Independent transfer: remove one likelihood after two observations
Try this changed-context calculation before reading the answer. A scalar parameter x has prior N(0, 9), with variance nine squared units. Observation one measures y1 = 3 units with Gaussian noise variance 3 squared units. Observation two measures y2 = negative 1 unit with Gaussian noise variance 1 squared unit. Use Gaussian natural parameters to compute the posterior after both observations, then remove only the first observation’s site to form its cavity.
The prior precision is 1/9 and linear parameter zero. Observation one contributes precision 1/3 and linear parameter 3/3 = 1. Observation two contributes precision 1 and linear parameter negative 1. Thus total precision is 1/9 + 1/3 + 1 = 13/9 per squared unit, while total linear parameter is zero. The posterior mean is 0 units and variance is 9/13, approximately 0.6923077 squared units.
To remove observation one’s site, subtract its precision 1/3 and linear parameter 1. The cavity precision becomes 13/9 minus 1/3 = 10/9 and the cavity linear parameter becomes negative 1. Therefore its variance is 9/10 = 0.9 squared units and its mean is negative 1 divided by (10/9) = negative 0.9 units. This cavity contains the prior and observation two, not the prior alone.
Check by an independent route: combine N(0, 9) directly with observation two N(negative 1, 1). Precision is 1/9 + 1 = 10/9, and the natural linear parameter is negative 1. The result is exactly N(negative 0.9, 0.9). Matching both routes confirms the cavity subtraction, demonstrates why site removal matters and provides a reusable test for an EP implementation. The example is exact Gaussian algebra, not a measured dataset.
What this method can and cannot promise
Expectation propagation is a family of approximation procedures, not one universally specified software implementation. The factors, chosen site family, moment integration, update order, damping and convergence criteria all affect the result. The standard scalar formulas in this guide teach the mechanism but do not determine how every large-scale classifier or structured probabilistic program should be built.
EP may be attractive when factor-wise moments are accessible and the global tractable family supports efficient computation. Its accuracy can be excellent on some problems and poor on others. An accurate mean does not imply an accurate rare-event probability, and a stable fixed point does not necessarily imply the true posterior was recovered.
Do not generalise from a single half-normal or conjugate Gaussian example to a production risk model. Validate under realistic scale, domain shift, data quality and decision costs. The original papers provide formal context and comparisons; a particular deployment still requires its own evidence.
Finally, distinguish the learning algorithm from artificial intelligence as a broad social category. EP is a technical approximate-inference mechanism used in certain probabilistic models; it does not create artificial superintelligence, human-like awareness or guaranteed trustworthy judgments.
Frequently asked questions
Is EP the same as expectation maximisation?
No. EP approximates a posterior distribution by revising local factor approximations and matching moments. Expectation–maximisation is a different optimisation framework usually formulated through alternating inference over latent variables and parameter updates. Similar initials do not make the algorithms interchangeable.
Why must a site be removed before being updated?
The cavity q divided by the old site excludes that factor’s current approximate contribution. Otherwise repeated updates would keep adding the same observation. This is a replacement operation, not another independent data point.
Does moment matching recover the true posterior?
Not generally. Matching a Gaussian’s mean and variance preserves only those moments of the tilted distribution. Skewness, multiple modes, hard support and tail behaviour may remain wrong, as the positive-threshold example shows.
What is the difference between the cavity and the tilted distribution?
The cavity contains the prior and all other approximate factors but removes the selected site’s contribution. The tilted distribution multiplies that cavity by the original selected exact factor, then normalises it for moment calculation.
Are EP sites themselves probability distributions?
Not necessarily. Gaussian EP sites are often unnormalised Gaussian-shaped potentials, and some may have negative precision contributions while the global q remains valid. Validate the cavity and global approximation separately.
Can damping guarantee convergence?
No. Damping can moderate proposed site changes and help some problems, but generic EP convergence is not guaranteed. Vary initialisations and schedules, check validity and compare with reference calculations.
Is EP better than variational inference or sampling?
There is no universal ranking. Compare accuracy, uncertainty quality, computational cost and stability for the particular model and decision. Each method has different assumptions and characteristic failures.
Why can a Gaussian projection violate a hard boundary?
A normal Gaussian has positive density on every real value, while a truncated or constrained target may have zero density on one side. Moment matching does not preserve support, so the approximation can assign probability to impossible values.
How do I test my first EP implementation?
Start with a conjugate Gaussian control for which exact answers are known, verify cavity subtraction, add the half-normal threshold example, then test multi-site updates under different orders and damping. Check both numerical validity and decision-relevant predictions.
Does EP imply a model is intelligent or fair?
No. EP is an approximate-inference algorithm. Its mathematical operation does not certify that the data are representative, the objectives are appropriate or the resulting decisions are fair or safe.
Original primary research sources
Thomas Minka (2001), Expectation Propagation for Approximate Bayesian Inference, UAI 2001 and associated MIT doctoral research. The author’s original project page supplies the paper, thesis, correction notice and related technical notes. The original paper is also available as an author research archive record with the original conference identification. The algorithm’s cavity, tilted distribution and moment-projection structure comes from this work.
Thomas Minka (2001), A Family of Algorithms for Approximate Bayesian Inference, Massachusetts Institute of Technology PhD thesis. This original institutional record offers the broader relationship among assumed-density filtering, approximate message passing and expectation propagation. The examples and editorial diagrams in this article are independent teaching constructions, not excerpts from the thesis.
Manfred Opper and Ole Winther (2005), Expectation Consistent Approximate Inference, Journal of Machine Learning Research, volume 6, pp. 2177–2204; and their Expectation Consistent Free Energies for Approximate Inference conference paper. These original sources illuminate the more general free-energy and moment-consistency perspective.
Aki Vehtari and colleagues, Expectation Propagation as a Way of Life: A Framework for Bayesian Inference on Partitioned Data. This original research investigates divide-and-conquer Bayesian inference and how local partitions are combined. It informs the caution that local posteriors cannot simply be multiplied without regard to shared priors or other dependencies.
The Gaussian and half-normal numerical examples in this article are independently derived from stated distributions. In particular, the half-normal mean sqrt(2/pi), variance 1 minus 2/pi, and transformed site natural parameters were checked algebraically. No performance benchmark, inference run or device test is claimed.
Verified eduKate AI learning routes
Read the existing How AI Works hub for the parent route across machine-learning methods. This article owns the specific EP cavity and moment-matching reader question, not generic Bayesian inference or all AI uncertainty.
For a related but separate model with tractable Gaussian regression posteriors, follow Gaussian Process Regression. For a different practical uncertainty approximation rather than iterative factor-wise inference, see Monte Carlo Dropout. Both destinations are existing specialist owners, not guessed new articles.
A precise next exercise is to calculate a Gaussian posterior in natural parameters, subtract one old site’s contribution, rebuild its tilted distribution from the exact likelihood and explain what moment projection preserves—and what it loses. If you can do that without double-counting an observation, you have the essential mechanism of EP.