Bilinear Gaussian Belief Propagation (BiGaBP)

• BiGaBP is a family of Gaussian message-passing methods enabling inference on factor graphs with bilinear couplings between two groups of uncertain variables.
• It covers explicit bilinear closures, iterative linearization in deep models, and surrogate modeling for rigid-body localization, reflecting diverse algorithmic approaches.
• Key advances include Gaussian moment matching, adaptive damping, and memory reduction techniques that improve convergence and scalability in complex inference tasks.

Searching arXiv for Bilinear Gaussian Belief Propagation and closely related Gaussian BP work on bilinear inference. Bilinear Gaussian Belief Propagation (BiGaBP) is best understood, across the cited literature, as a family of Gaussian message-passing methods for inference problems in which observations depend on bilinear couplings between two groups of uncertain variables. Typical instances include matrix factorization models of the form $x_{ij}=\mathbf u_i^\top \mathbf v_j$, weight-times-activation interactions in deep factor graphs, and rigid-body models whose kinematics induce coupled products after reparameterization or linearization. The literature does not use the label uniformly: one matrix-completion work develops a Gaussian-parameterized loopy BP method for a bilinear latent factor model without naming it BiGaBP, a rigid-body localization paper explicitly invokes a “bilinear Gaussian belief propagation (GaBP) framework” while operationally using linearized and bivariate GaBP stages, and deep Gaussian factor-graph work contains pervasive bilinear structure but handles it through iterative local linearization rather than dedicated bilinear-factor message formulas (Okajima et al., 2021, Führling et al., 2024, Nabarro et al., 2023).

1. Definition and scope

The central object in BiGaBP is a factor graph whose observation factors couple two latent variable families multiplicatively or bilinearly. In the matrix-completion setting, each observed entry depends on $\mathbf u_i^\top \mathbf v_j$; in deep factor graphs, dense and convolutional layers involve products between uncertain parameters and uncertain activations; in 6D rigid-body localization, the underlying physics contains terms such as $\boldsymbol{Q}\boldsymbol{c}_n$ and $$[\boldsymbol{\omega}]^\times \boldsymbol{Q}\boldsymbol{c}_n$$, which are bilinear or coupled before the paper’s linearization steps (Okajima et al., 2021, Nabarro et al., 2023, Führling et al., 2024).

A concise way to organize the literature is to distinguish between explicit bilinear Gaussian closures, generic Gaussianization of bilinear/nonlinear factors, and learned but non-bilinear GBP systems used as engineering contrasts.

Work	Relation to BiGaBP	Main treatment
(Okajima et al., 2021)	Closest explicit bilinear message-passing instance	Gaussian-parameterized loopy BP for bilinear matrix factorization, plus AMP/TAP-like reduction
(Nabarro et al., 2023)	Architecture-level generalization	Iterative linearized GBP for factors containing weight-times-activation couplings
(Führling et al., 2024)	Uses “bilinear GaBP” terminology	Multi-stage linear and bivariate GaBP on linearized rigid-body models
(Tang et al., 29 Jan 2026)	Not BiGaBP	Learned pairwise Gaussian MRF with unrolled GBP

This suggests that BiGaBP is not a single canonical algorithm but a design space. The invariant theme is Gaussian approximate inference on loopy graphs whose semantics are driven by bilinear interactions, while the main axis of variation is how those interactions are made tractable: closed-form Gaussian closure, perturbative reduction, first-order linearization, or linearized surrogate modeling.

2. Bilinear probabilistic models and factor-graph structure

The cleanest bilinear formulation in the cited material is noisy low-rank matrix completion. With partially observed matrix $Y=\{y_{ij}\}$, observed index set $\Omega\subseteq \{1,\dots,N\}\times\{1,\dots,M\}$, and factorization $X=UV^\top$, the optimization surrogate is

$$\min_{U\in \mathbb{R}^{N\times R},\,V\in \mathbb{R}^{M\times R}} \frac12 \sum_{(i,j)\in\Omega} \bigl(y_{ij}-\mathbf u_i^T \mathbf v_j\bigr)^2 +\frac{\lambda}{2}\sum_{i=1}^N \|\mathbf u_i\|^2 +\frac{\lambda}{2}\sum_{j=1}^M \|\mathbf v_j\|^2.$$

The corresponding probabilistic model uses Gaussian likelihood

$$p(y_{ij}\mid \mathbf u_i,\mathbf v_j)\propto \exp\!\left[-\frac{\beta}{2}\bigl(y_{ij}-\mathbf u_i^T\mathbf v_j\bigr)^2\right]$$

and Gaussian priors

$$p(\mathbf u_i)\propto \exp\!\left(-\frac{\beta\lambda}{2}\|\mathbf u_i\|^2\right),\qquad p(\mathbf v_j)\propto \exp\!\left(-\frac{\beta\lambda}{2}\|\mathbf v_j\|^2\right).$$

Its factor graph is a sparse bipartite graph with vector-valued variable nodes $\mathbf u_i^\top \mathbf v_j$0 and $\mathbf u_i^\top \mathbf v_j$1, observation factors for $\mathbf u_i^\top \mathbf v_j$2, and unary prior factors. Every observation factor couples one $\mathbf u_i^\top \mathbf v_j$3 and one $\mathbf u_i^\top \mathbf v_j$4 through the bilinear term $\mathbf u_i^\top \mathbf v_j$5 (Okajima et al., 2021).

The deep-factor-graph formulation generalizes this idea from matrix factorization to layered architectures. Factors are written in energy form

$\mathbf u_i^\top \mathbf v_j$6

and a generic layer factor has the form

$\mathbf u_i^\top \mathbf v_j$7

When $\mathbf u_i^\top \mathbf v_j$8 is dense or convolutional, the preactivation contains sums of products between uncertain weights and uncertain activations. For example, the dense layer factor uses

$\mathbf u_i^\top \mathbf v_j$9

and the convolution factor uses

$\boldsymbol{Q}\boldsymbol{c}_n$0

The bilinear component is the inner product between uncertain patch variables and uncertain filter variables (Nabarro et al., 2023).

The rigid-body localization model illustrates a different route to bilinear structure. Sensor positions and velocities satisfy

$\boldsymbol{Q}\boldsymbol{c}_n$1

so the unknown pose and motion parameters enter through coupled products. The paper’s range and Doppler measurements are then transformed into linear or bivariate Gaussian subproblems by squaring ranges, multiplying range and Doppler, introducing augmented variables, and applying a small-angle approximation for $\boldsymbol{Q}\boldsymbol{c}_n$2 (Führling et al., 2024).

3. Gaussian parameterization and message closure

The primary technical obstacle in BiGaBP is that exact loopy BP on bilinear graphs is generally not closed under Gaussian families. In the matrix-completion graph, the factor-to-variable update

$\boldsymbol{Q}\boldsymbol{c}_n$3

involves the nonlinear interaction

$\boldsymbol{Q}\boldsymbol{c}_n$4

and even a Gaussian incoming cavity message in $\boldsymbol{Q}\boldsymbol{c}_n$5 does not yield a Gaussian outgoing dependence on $\boldsymbol{Q}\boldsymbol{c}_n$6. The cited GPBP approach resolves this by parameterizing cavity messages as multivariate Gaussians in natural form,

$\boldsymbol{Q}\boldsymbol{c}_n$7

$\boldsymbol{Q}\boldsymbol{c}_n$8

with full $\boldsymbol{Q}\boldsymbol{c}_n$9 precision matrices per message rather than scalar variances. Cavity means are $$[\boldsymbol{\omega}]^\times \boldsymbol{Q}\boldsymbol{c}_n$$0 and $$[\boldsymbol{\omega}]^\times \boldsymbol{Q}\boldsymbol{c}_n$$1 (Okajima et al., 2021).

The distinctive step is a moment-matching closure applied to the cavity factor-to-node distribution itself. In the zero-temperature or MAP limit $$[\boldsymbol{\omega}]^\times \boldsymbol{Q}\boldsymbol{c}_n$$2, the algorithm computes first and second moments of the cavity-integrated factor message and fits a Gaussian reproducing them. The resulting edgewise updates involve rank-one outer products such as $$[\boldsymbol{\omega}]^\times \boldsymbol{Q}\boldsymbol{c}_n$$3 and linear terms such as $$[\boldsymbol{\omega}]^\times \boldsymbol{Q}\boldsymbol{c}_n$$4, modulated by uncertainty corrections of the form $$[\boldsymbol{\omega}]^\times \boldsymbol{Q}\boldsymbol{c}_n$$5, where the $$[\boldsymbol{\omega}]^\times \boldsymbol{Q}\boldsymbol{c}_n$$6-terms are derived from cavity means and inverse precisions. These $$[\boldsymbol{\omega}]^\times \boldsymbol{Q}\boldsymbol{c}_n$$7-coefficients are the main distinction from ALS-MP: if they are dropped, GPBP reduces to ALS-MP. The paper interprets $$[\boldsymbol{\omega}]^\times \boldsymbol{Q}\boldsymbol{c}_n$$8 and $$[\boldsymbol{\omega}]^\times \boldsymbol{Q}\boldsymbol{c}_n$$9 as proportional to Fisher information matrices for cavity distributions, so the $Y=\{y_{ij}\}$0-terms attenuate observation influence when cavity uncertainty is large (Okajima et al., 2021).

A different Gaussianization strategy appears in deep factor graphs. Messages are represented in canonical form by precision $Y=\{y_{ij}\}$1 and information vector $Y=\{y_{ij}\}$2, with variable-to-factor updates

$Y=\{y_{ij}\}$3

For nonlinear factors, the factor’s measurement function is linearized every time a factor message is computed,

$Y=\{y_{ij}\}$4

yielding approximate Gaussian factor parameters

$Y=\{y_{ij}\}$5

This is not exact Gaussian BP on the original bilinear or nonlinear model; it is iterative linearized Gaussian BP, closer in spirit to Gauss-Newton or iterated-EKF-style local Gaussian approximation داخل BP (Nabarro et al., 2023).

4. Linearization, bivariate GaBP, and neighboring approximations

An important conceptual boundary in the literature is the distinction between explicit bilinear message derivations and Gaussian BP on linearized surrogates. The rigid-body localization paper uses the phrase “bilinear Gaussian belief propagation (GaBP) framework,” but its operational stages are not exact BP on product factors. Instead, it constructs a multi-stage pipeline: linear GaBP for sensor positions from range, linear GaBP for sensor velocities from Doppler/range, bivariate GaBP for $Y=\{y_{ij}\}$6, and bivariate GaBP for $Y=\{y_{ij}\}$7, each derived after algebraic preprocessing and small-angle linearization (Führling et al., 2024).

For pose estimation, the small-angle approximation

$Y=\{y_{ij}\}$8

yields

$Y=\{y_{ij}\}$9

after which the paper writes a linear observation model

$\Omega\subseteq \{1,\dots,N\}\times\{1,\dots,M\}$0

The associated bivariate GaBP then uses scalar Gaussian messages, soft interference cancellation, effective residual variances, Gaussian denoisers with prior variances $\Omega\subseteq \{1,\dots,N\}\times\{1,\dots,M\}$1 and $\Omega\subseteq \{1,\dots,N\}\times\{1,\dots,M\}$2, and a damping parameter $\Omega\subseteq \{1,\dots,N\}\times\{1,\dots,M\}$3. The motion stage for $\Omega\subseteq \{1,\dots,N\}\times\{1,\dots,M\}$4 is structurally analogous (Führling et al., 2024).

This suggests a useful taxonomy. One branch of BiGaBP keeps the bilinear factor graph explicit and closes messages by Gaussian moment matching, as in matrix completion. Another branch keeps the model architecture explicit but Gaussianizes each local factor through repeated linearization, as in deep factor graphs. A third branch begins from a physically bilinear problem but carries out Gaussian BP only after reformulating the problem into linear or bivariate Gaussian submodels. The three branches share locality, Gaussian messages, and iterative refinement, but they differ sharply in how much bilinear algebra is preserved in the message updates themselves (Okajima et al., 2021, Nabarro et al., 2023, Führling et al., 2024).

A common misconception is therefore that every method involving GaBP and multiplicative unknowns is an explicit BiGaBP algorithm in the narrow product-factor sense. The cited literature shows otherwise: weight-times-activation couplings may be handled by first-order local Gaussianization, and bilinear rigid-body physics may be handled by linear surrogates plus scalar Gaussian interference models, without ever deriving product-specific bilinear factor messages (Nabarro et al., 2023, Führling et al., 2024).

5. Scheduling, damping, memory, and computational structure

Across the literature, damping is treated as a practical necessity rather than an optional refinement. In matrix completion, the GPBP paper emphasizes that damping is crucial for optimal performance. It applies “probabilistic damping” to the factor-to-node contributions entering the natural parameters, with damping coefficient $\Omega\subseteq \{1,\dots,N\}\times\{1,\dots,M\}$5. The authors remark that this is not exactly a weighted average of old and new natural parameters themselves, but rather a damping of the updated contributions, and that it coincides with the formal probabilistic damping of Pretti only for small $\Omega\subseteq \{1,\dots,N\}\times\{1,\dots,M\}$6. Empirically, without damping, both GPBP and ALS-MP can converge very slowly or fail to achieve the performance predicted by population dynamics, and damping is even more important for the approximate low-memory variants (Okajima et al., 2021).

The same paper also derives a memory-friendly approximation by applying the standard BP-to-AMP perturbative idea: when node degree is sufficiently large, a cavity message differs only slightly from the full node marginal. The resulting approxGPBP does not store all edgewise cavity messages, approximates $\Omega\subseteq \{1,\dots,N\}\times\{1,\dots,M\}$7 and $\Omega\subseteq \{1,\dots,N\}\times\{1,\dots,M\}$8 plus first-order cavity corrections, and retains TAP/Onsager-like residual terms. The full GPBP complexity is stated as time $\Omega\subseteq \{1,\dots,N\}\times\{1,\dots,M\}$9 and space $X=UV^\top$0, while the reduced-memory version lowers space from

$X=UV^\top$1

The paper notes that this is beneficial when $X=UV^\top$2, and that Sherman–Morrison updates can become numerically unstable when $X=UV^\top$3 is near zero (Okajima et al., 2021).

Deep Gaussian factor graphs emphasize a different set of engineering tools. Their message schedule alternates forward and backward sweeps through layers, supports layerwise asynchronous updates, and relies on damping and dropout on factor-to-variable messages for stability. For factor updates, the paper exploits low-rank structure and diagonal incoming precisions to reduce complexity using Woodbury identity and shared intermediates from

$X=UV^\top$4

in time and from

$X=UV^\top$5

in space, for a factor connecting $X=UV^\top$6 scalar variables with observation dimension $X=UV^\top$7 (Nabarro et al., 2023).

Although not bilinear, the learned pairwise GBP system for depth completion offers a useful comparison for loopy Gaussian scheduling. It decomposes edges into four directional loop-free local subgraphs and a non-local edge set, performing serial local sweeps over $X=UV^\top$8 and parallel propagation over $X=UV^\top$9. Messages are damped in information form by

$$\min_{U\in \mathbb{R}^{N\times R},\,V\in \mathbb{R}^{M\times R}} \frac12 \sum_{(i,j)\in\Omega} \bigl(y_{ij}-\mathbf u_i^T \mathbf v_j\bigr)^2 +\frac{\lambda}{2}\sum_{i=1}^N \|\mathbf u_i\|^2 +\frac{\lambda}{2}\sum_{j=1}^M \|\mathbf v_j\|^2.$$0

where the damping rate $$\min_{U\in \mathbb{R}^{N\times R},\,V\in \mathbb{R}^{M\times R}} \frac12 \sum_{(i,j)\in\Omega} \bigl(y_{ij}-\mathbf u_i^T \mathbf v_j\bigr)^2 +\frac{\lambda}{2}\sum_{i=1}^N \|\mathbf u_i\|^2 +\frac{\lambda}{2}\sum_{j=1}^M \|\mathbf v_j\|^2.$$1 is learned by the graph-construction network. The appendix reports a final model with a total of 13 iterations and practical compute figures of 62.58 GFLOPs and 11.03 ms for the graph-construction network versus 0.26 GFLOPs and 33.54 ms for GBP at $$\min_{U\in \mathbb{R}^{N\times R},\,V\in \mathbb{R}^{M\times R}} \frac12 \sum_{(i,j)\in\Omega} \bigl(y_{ij}-\mathbf u_i^T \mathbf v_j\bigr)^2 +\frac{\lambda}{2}\sum_{i=1}^N \|\mathbf u_i\|^2 +\frac{\lambda}{2}\sum_{j=1}^M \|\mathbf v_j\|^2.$$2 resolution (Tang et al., 29 Jan 2026).

6. Empirical domains and representative results

The matrix-completion study provides the most direct evidence about explicit bilinear Gaussian closures. On synthetic Gaussian-noise problems with $$\min_{U\in \mathbb{R}^{N\times R},\,V\in \mathbb{R}^{M\times R}} \frac12 \sum_{(i,j)\in\Omega} \bigl(y_{ij}-\mathbf u_i^T \mathbf v_j\bigr)^2 +\frac{\lambda}{2}\sum_{i=1}^N \|\mathbf u_i\|^2 +\frac{\lambda}{2}\sum_{j=1}^M \|\mathbf v_j\|^2.$$3, $$\min_{U\in \mathbb{R}^{N\times R},\,V\in \mathbb{R}^{M\times R}} \frac12 \sum_{(i,j)\in\Omega} \bigl(y_{ij}-\mathbf u_i^T \mathbf v_j\bigr)^2 +\frac{\lambda}{2}\sum_{i=1}^N \|\mathbf u_i\|^2 +\frac{\lambda}{2}\sum_{j=1}^M \|\mathbf v_j\|^2.$$4, and matched regularization $$\min_{U\in \mathbb{R}^{N\times R},\,V\in \mathbb{R}^{M\times R}} \frac12 \sum_{(i,j)\in\Omega} \bigl(y_{ij}-\mathbf u_i^T \mathbf v_j\bigr)^2 +\frac{\lambda}{2}\sum_{i=1}^N \|\mathbf u_i\|^2 +\frac{\lambda}{2}\sum_{j=1}^M \|\mathbf v_j\|^2.$$5, GPBP and ALS-MP have almost the same reconstruction performance. Under Bernoulli–Gaussian sparse corruption,

$$\min_{U\in \mathbb{R}^{N\times R},\,V\in \mathbb{R}^{M\times R}} \frac12 \sum_{(i,j)\in\Omega} \bigl(y_{ij}-\mathbf u_i^T \mathbf v_j\bigr)^2 +\frac{\lambda}{2}\sum_{i=1}^N \|\mathbf u_i\|^2 +\frac{\lambda}{2}\sum_{j=1}^M \|\mathbf v_j\|^2.$$6

GPBP achieves lower nRMSE than ALS-MP at their respective optimal $$\min_{U\in \mathbb{R}^{N\times R},\,V\in \mathbb{R}^{M\times R}} \frac12 \sum_{(i,j)\in\Omega} \bigl(y_{ij}-\mathbf u_i^T \mathbf v_j\bigr)^2 +\frac{\lambda}{2}\sum_{i=1}^N \|\mathbf u_i\|^2 +\frac{\lambda}{2}\sum_{j=1}^M \|\mathbf v_j\|^2.$$7, and ALS-MP prefers much larger $$\min_{U\in \mathbb{R}^{N\times R},\,V\in \mathbb{R}^{M\times R}} \frac12 \sum_{(i,j)\in\Omega} \bigl(y_{ij}-\mathbf u_i^T \mathbf v_j\bigr)^2 +\frac{\lambda}{2}\sum_{i=1}^N \|\mathbf u_i\|^2 +\frac{\lambda}{2}\sum_{j=1}^M \|\mathbf v_j\|^2.$$8, suggesting stronger regularization is needed to prevent overfitting. On MovieLens, GPBP and ALS-MP are similar on denser data, ALS-MP converges faster on 10M but shows stronger overfitting, approxGPBP outperforms approxALS-MP across all tested datasets, and the advantage of GPBP-style uncertainty is especially visible under sparser conditions (Okajima et al., 2021).

The deep factor-graph work shows that Gaussian message passing with uncertain parameters and activations can scale beyond classical bilinear latent-factor models. It reports toy nonlinear tasks including XOR and nonlinear regression, a video denoising setting in which learnable deep factor graphs outperform a classical pairwise smoothing factor graph and deeper models outperform shallow ones, and continual image classification with strong sample efficiency on MNIST. The reported MNIST test accuracy is $$\min_{U\in \mathbb{R}^{N\times R},\,V\in \mathbb{R}^{M\times R}} \frac12 \sum_{(i,j)\in\Omega} \bigl(y_{ij}-\mathbf u_i^T \mathbf v_j\bigr)^2 +\frac{\lambda}{2}\sum_{i=1}^N \|\mathbf u_i\|^2 +\frac{\lambda}{2}\sum_{j=1}^M \|\mathbf v_j\|^2.$$9, comparable to CNN+Adam with a replay buffer of $$p(y_{ij}\mid \mathbf u_i,\mathbf v_j)\propto \exp\!\left[-\frac{\beta}{2}\bigl(y_{ij}-\mathbf u_i^T\mathbf v_j\bigr)^2\right]$$0 examples, and asynchronous training nearly matches synchronous training. These results do not isolate bilinear message derivations, but they support the broader proposition that Gaussian local inference over uncertain weights and activations is operationally viable in large, layered systems (Nabarro et al., 2023).

The rigid-body localization paper evaluates a multi-stage GaBP pipeline in a simulation setting with $$p(y_{ij}\mid \mathbf u_i,\mathbf v_j)\propto \exp\!\left[-\frac{\beta}{2}\bigl(y_{ij}-\mathbf u_i^T\mathbf v_j\bigr)^2\right]$$1 sensors, $$p(y_{ij}\mid \mathbf u_i,\mathbf v_j)\propto \exp\!\left[-\frac{\beta}{2}\bigl(y_{ij}-\mathbf u_i^T\mathbf v_j\bigr)^2\right]$$2 anchors, $$p(y_{ij}\mid \mathbf u_i,\mathbf v_j)\propto \exp\!\left[-\frac{\beta}{2}\bigl(y_{ij}-\mathbf u_i^T\mathbf v_j\bigr)^2\right]$$3 Monte Carlo trials, and Gaussian priors $$p(y_{ij}\mid \mathbf u_i,\mathbf v_j)\propto \exp\!\left[-\frac{\beta}{2}\bigl(y_{ij}-\mathbf u_i^T\mathbf v_j\bigr)^2\right]$$4, $$p(y_{ij}\mid \mathbf u_i,\mathbf v_j)\propto \exp\!\left[-\frac{\beta}{2}\bigl(y_{ij}-\mathbf u_i^T\mathbf v_j\bigr)^2\right]$$5, $$p(y_{ij}\mid \mathbf u_i,\mathbf v_j)\propto \exp\!\left[-\frac{\beta}{2}\bigl(y_{ij}-\mathbf u_i^T\mathbf v_j\bigr)^2\right]$$6, $$p(y_{ij}\mid \mathbf u_i,\mathbf v_j)\propto \exp\!\left[-\frac{\beta}{2}\bigl(y_{ij}-\mathbf u_i^T\mathbf v_j\bigr)^2\right]$$7. It reports improved performance compared to equivalent state-of-the-art techniques for sensor position estimation, translation estimation, rotation estimation, sensor velocity estimation, translational velocity estimation, and angular velocity estimation. Its BiGaBP relevance lies less in a product-factor derivation than in demonstrating that joint bivariate GaBP stages can be effective once a bilinear physics problem has been linearized into Gaussian submodels (Führling et al., 2024).

The depth-completion paper is not a BiGaBP result, but it helps delimit what belongs to the bilinear category and what merely shares GBP infrastructure. It reports SOTA performance on NYUv2 and KITTI, robustness across sparsity levels and sparsity patterns, and strong gains from moving from fixed local GBP-like models to dynamically parameterized and non-local learned Gaussian MRFs. As an engineering reference, it shows how unrolled loopy GBP, learned damping, adaptive graph construction, and uncertainty-aware losses can be embedded end-to-end inside a neural system, even though the underlying factors are pairwise linear-Gaussian rather than bilinear (Tang et al., 29 Jan 2026).

7. Limitations, misconceptions, and conceptual boundaries

A first limitation is theoretical exactness. In matrix factorization, the paper explicitly notes that the matrix factorization model is a singular statistical model, so Gaussian posterior approximations are not asymptotically exact in a strict Bayesian sense. This makes GPBP empirically effective rather than theoretically exact as a posterior approximation (Okajima et al., 2021).

A second limitation is convergence. The deep-factor-graph work states that tree-structured BP converges exactly, loopy BP has no guarantee in general, and in linear-Gaussian models, if GBP converges, posterior means are correct. No comparable guarantee is stated for nonlinear factors handled by iterative linearization. The depth-completion paper likewise places itself in loopy GBP, relying on damping and graph decomposition rather than general convergence theory, and explicitly notes that GBP typically converges to an accurate mean $$p(y_{ij}\mid \mathbf u_i,\mathbf v_j)\propto \exp\!\left[-\frac{\beta}{2}\bigl(y_{ij}-\mathbf u_i^T\mathbf v_j\bigr)^2\right]$$8 but not the exact precision $$p(y_{ij}\mid \mathbf u_i,\mathbf v_j)\propto \exp\!\left[-\frac{\beta}{2}\bigl(y_{ij}-\mathbf u_i^T\mathbf v_j\bigr)^2\right]$$9, motivating a learned residual correction for the precision map (Nabarro et al., 2023, Tang et al., 29 Jan 2026).

A third limitation is model approximation. The rigid-body localization method depends on neglecting second-order range-noise terms, on approximate Gaussianity of residual interference under soft interference cancellation, and on the small-angle approximation

$$p(\mathbf u_i)\propto \exp\!\left(-\frac{\beta\lambda}{2}\|\mathbf u_i\|^2\right),\qquad p(\mathbf v_j)\propto \exp\!\left(-\frac{\beta\lambda}{2}\|\mathbf v_j\|^2\right).$$0

Its “bilinear GaBP” label therefore refers to the original coupled estimation problem more than to the final message equations, which are scalar Gaussian updates on linearized models (Führling et al., 2024).

A final misconception concerns the meaning of “bilinear” in BiGaBP. The cited literature shows that three different situations are often conflated: exact bilinear factor graphs closed by Gaussian moment matching, generic Gaussian BP frameworks that contain bilinear substructures but handle them by Jacobian linearization, and learned loopy GBP systems with no bilinear latent interaction at all. A precise usage of BiGaBP is therefore best reserved for Gaussian belief-propagation methods whose inferential bottleneck is genuinely the coupling of two uncertain variable families through bilinear terms, whether those couplings are handled by moment-matched Gaussian closure, by perturbative AMP/TAP reduction, or by linearized surrogate stages specialized to the application at hand (Okajima et al., 2021, Nabarro et al., 2023, Führling et al., 2024).