DS-UGE: Doubly-Stochastic Unbiased Gradient Estimator

• DS-UGE is a family of gradient estimators that uses two independent stochastic mechanisms to achieve unbiased estimates in complex, high-dimensional settings.
• It leverages randomized multilevel, telescoping, and debiasing techniques to remove bias and control variance, enhancing computational efficiency.
• DS-UGE methodologies are applied in diverse domains such as PINNs, state-space models, and distributionally robust learning, improving scalability and convergence.

Searching arXiv for the cited papers and topic variants to ground the article in the current literature. {"query":"\"Doubly-Stochastic Unbiased Gradient Estimator\" OR DS-UGE arXiv", "max_results": 10} {"query":"(Badolle et al., 2024)", "max_results": 5} The Doubly-Stochastic Unbiased Gradient Estimator (DS-UGE) denotes a class of gradient constructions in which two distinct stochastic mechanisms are used while preserving unbiasedness of the target gradient. In current arXiv usage, the term is explicit in the Stochastic Dimension Implicit Functional Projection (SDIFP) framework for conservative high-dimensional PINNs, where DS-UGE is the backpropagation mechanism that combines spatial mini-batch sampling with differential-operator subsampling (Liang, 31 Mar 2026). Closely related constructions appear under different names in conditional stochastic optimization, continuous-time state-space models, distributionally robust learning, zeroth-order optimization, learning with covariate-dependent missingness, and doubly stochastic SGD more generally; across these settings, the common objective is to remove one or more sources of bias without reverting to full deterministic differentiation (Goda et al., 2022, Ballesio et al., 2021, Ghosh et al., 2020, Ma et al., 22 Oct 2025, Lee et al., 2018, Kim et al., 2024).

1. Terminology and scope

The terminology is not fully standardized. In SDIFP, DS-UGE is the paper’s own name for a composite stochastic gradient estimator used to train a projection-based conservative PDE solver (Liang, 31 Mar 2026). In other literatures, essentially the same design principle appears under different labels: stochastic doubly robust gradient in learning from incomplete observational data, single-term randomized MLMC gradient estimator in conditional stochastic optimization, doubly randomized unbiased estimator of the score function in partially observed diffusions, and doubly stochastic gradients in finite-sum problems whose components are themselves expectations (Lee et al., 2018, Goda et al., 2022, Ballesio et al., 2021, Kim et al., 2024).

A concise way to organize the literature is to view DS-UGE as a family of estimators indexed by the two stochastic sources they exploit and by the quantity they target.

Formulation	Stochastic sources	Target
SDIFP DS-UGE (Liang, 31 Mar 2026)	Spatial mini-batch; differential-operator subsampling	PDE risk gradient
CSO MLMC estimator (Goda et al., 2022)	Outer $\xi$; inner conditional $\eta\mid\xi$, with randomized level	Nested-expectation gradient
State-space score estimator (Ballesio et al., 2021)	Discretization level; CCPF-based increment estimation	Score function
SDRG (Lee et al., 2018)	Minibatch/index; missingness or propensity mechanism	Target gradient under observational bias
ZOO estimator family (Ma et al., 22 Oct 2025)	Random direction $v$; telescoping level $n$	Zeroth-order gradient
Doubly SGD (Kim et al., 2024)	Outer minibatch; inner Monte Carlo samples	Finite-sum expectation gradient

This heterogeneity matters conceptually. Some papers use the phrase to denote a specific estimator; others provide a general stochastic-optimization setting in which unbiased doubly stochastic gradients may or may not be available. A frequent source of confusion is to treat all doubly stochastic gradients as unbiased. That is not correct: the independent-minibatch formulation in doubly SGD is unbiased, whereas random reshuffling is explicitly conditionally biased within an epoch (Kim et al., 2024).

2. Canonical probabilistic structure

Across the literature, a DS-UGE has the defining property that its expectation matches the exact gradient or score of interest. In SDIFP, the unbiasedness theorem is written as

$$\mathbb{E}_{\mathcal{B}_{\mathrm{res}}}\!\left[\mathbb{E}_{I,J}[g_{I,J}(\theta)]\right] = \int_{\mathcal X} \big(\mathcal{L}[\tilde u_\theta]-R\big)\,\nabla_\theta \mathcal{L}[\tilde u_\theta]\,dx \equiv \nabla_\theta \mathcal J(\theta),$$

with independence between the residual batch $\mathcal{B}_{\mathrm{res}}$ and the operator subsets $I,J$ playing a central role (Liang, 31 Mar 2026).

In conditional stochastic optimization, the target objective is

$$F(x)=\mathbb{E}_{\xi}\!\left[f_{\xi}\!\left(\mathbb{E}_{\eta\mid \xi}[g_{\eta}(x,\xi)]\right)\right],$$

and the gradient is itself a nested expectation,

$$\nabla F(x)= \mathbb{E}_{\xi}\!\left[ \left(\mathbb{E}_{\eta\mid\xi}[\nabla g_{\eta}(x,\xi)]\right)^\top \nabla f_{\xi}\!\left(\mathbb{E}_{\eta\mid\xi}[g_{\eta}(x,\xi)]\right) \right].$$

The DS-UGE there is a randomized multilevel estimator of this nested gradient, constructed so that $$\mathbb{E}[\Delta\psi_\ell(x)/\omega_\ell]=\nabla F(x)$$ after summing the telescoping series over levels (Goda et al., 2022).

In zeroth-order optimization, the same principle is expressed through directional derivatives. The unbiased estimator family $\eta\mid\xi$0 is built so that, conditional on a random direction $\eta\mid\xi$1, $\eta\mid\xi$2, and if $\eta\mid\xi$3 satisfies $\eta\mid\xi$4, then

$\eta\mid\xi$5

Here the two stochastic sources are the sampled direction and the sampled telescoping level (Ma et al., 22 Oct 2025).

A broader probabilistic reading is therefore consistent across domains: a DS-UGE is not merely a noisy gradient, but a carefully coupled estimator whose nested expectations, randomized levels, or independent subsamplings reconstruct the exact first-order object after averaging.

3. Multilevel, telescoping, and debiasing constructions

A major branch of DS-UGE methods uses randomized multilevel Monte Carlo. In conditional stochastic optimization, level $\eta\mid\xi$6 uses $\eta\mid\xi$7 inner samples, and the fine estimator is antithetically coupled to two half-sized coarse estimators. The level increment is

$\eta\mid\xi$8

and the single-term randomized estimator is $\eta\mid\xi$9. Unbiasedness follows because the expectations of the antithetic half-level terms equal $v$0, so the series telescopes exactly to $v$1 (Goda et al., 2022).

Distributionally robust learning uses the same logic for the gradient of the inner maximization problem. The randomized level $v$2 selects a subset size $v$3, the subset is split into two halves, and the multilevel correction is

$v$4

The unbiased estimator is then

$v$5

whose expectation equals the full robust gradient by a telescoping argument across subset sizes (Ghosh et al., 2020).

In continuous-time partially observed state-space models, unbiasedness must remove both discretization bias and the bias from estimating each discretized level. The paper therefore uses a doubly randomized scheme: first, sample a discretization level $v$6; second, estimate the level increment with a coupled conditional particle filter (CCPF). The final estimator has the form

$v$7

and is unbiased for the continuous-time score $v$8 even though only discretized dynamics are simulated (Ballesio et al., 2021).

The same architectural pattern appears outside direct gradient estimation. For underdamped Langevin dynamics, a randomized-level multilevel decomposition is combined with meeting-time debiasing at each level to eliminate both finite-time MCMC bias and discretization bias. That construction targets $v$9 rather than a gradient directly, but it is explicitly presented as a doubly randomized unbiased estimation scheme and is structurally adjacent to DS-UGE methodology (Ruzayqat et al., 2022).

4. Structural and domain-specific realizations

In SDIFP, DS-UGE is tied to a specific projection architecture. The raw network output $n$0 is projected by a global affine map

$n$1

where $n$2 are chosen from detached Monte Carlo estimates of raw moments so that the projected field satisfies exact mass and energy constraints. The closed-form coefficients are

$n$3

Because the large Monte Carlo quadrature for $n$4 is detached from the AD graph, DS-UGE reconstructs $n$5 and $n$6 through unbiased mini-batch estimators of $n$7 and $n$8, while independently subsampling the linear PDE operator during reverse-mode differentiation (Liang, 31 Mar 2026).

In incomplete-data learning, the relevant construction is Stochastic Doubly Robust Gradient (SDRG). The estimator

$n$9

combines an inverse-propensity or importance-weight correction with a regression-adjustment or control-variate term. Its defining theorem states that the estimator satisfies the double robustness property: it remains unbiased if either the weighting model or the regression model is correct. The paper does not use the term DS-UGE, but it explicitly fits the same conceptual family because the two stochastic sources are the sampled data index and the missingness or propensity mechanism (Lee et al., 2018).

A conceptually related but structurally different line of work concerns the generator gradient estimator for SDEs. The note on this topic shows that the identity

$$\mathbb{E}_{\mathcal{B}_{\mathrm{res}}}\!\left[\mathbb{E}_{I,J}[g_{I,J}(\theta)]\right] = \int_{\mathcal X} \big(\mathcal{L}[\tilde u_\theta]-R\big)\,\nabla_\theta \mathcal{L}[\tilde u_\theta]\,dx \equiv \nabla_\theta \mathcal J(\theta),$$0

is the stochastic adjoint-state formula for diffusions, and that the estimator is a close analogue of the exact Integral Path Algorithm for CTMCs. The paper does not name this a DS-UGE, but it explicitly characterizes it as an efficient, unbiased, doubly-stochastic-style estimator whose stochasticity comes from SDE path simulation together with auxiliary randomization (Badolle et al., 2024).

5. Complexity, variance, and scaling laws

One of the main attractions of DS-UGE constructions is computational scaling. In SDIFP, the full reverse-mode cost of differentiating through both the global projection quadrature and all differential operator terms is described as $$\mathbb{E}_{\mathcal{B}_{\mathrm{res}}}\!\left[\mathbb{E}_{I,J}[g_{I,J}(\theta)]\right] = \int_{\mathcal X} \big(\mathcal{L}[\tilde u_\theta]-R\big)\,\nabla_\theta \mathcal{L}[\tilde u_\theta]\,dx \equiv \nabla_\theta \mathcal J(\theta),$$1. DS-UGE replaces that by truncated stochastic differentiation, reported in the paper as $$\mathbb{E}_{\mathcal{B}_{\mathrm{res}}}\!\left[\mathbb{E}_{I,J}[g_{I,J}(\theta)]\right] = \int_{\mathcal X} \big(\mathcal{L}[\tilde u_\theta]-R\big)\,\nabla_\theta \mathcal{L}[\tilde u_\theta]\,dx \equiv \nabla_\theta \mathcal J(\theta),$$2, and also as $$\mathbb{E}_{\mathcal{B}_{\mathrm{res}}}\!\left[\mathbb{E}_{I,J}[g_{I,J}(\theta)]\right] = \int_{\mathcal X} \big(\mathcal{L}[\tilde u_\theta]-R\big)\,\nabla_\theta \mathcal{L}[\tilde u_\theta]\,dx \equiv \nabla_\theta \mathcal J(\theta),$$3. The stated mechanism is the decoupling of detached forward quadrature from the backward graph, together with independent subsampling of operator terms (Liang, 31 Mar 2026).

Randomized-MLMC DS-UGEs are typically analyzed through simultaneous variance and cost conditions. In conditional stochastic optimization, the level increments satisfy

$$\mathbb{E}_{\mathcal{B}_{\mathrm{res}}}\!\left[\mathbb{E}_{I,J}[g_{I,J}(\theta)]\right] = \int_{\mathcal X} \big(\mathcal{L}[\tilde u_\theta]-R\big)\,\nabla_\theta \mathcal{L}[\tilde u_\theta]\,dx \equiv \nabla_\theta \mathcal J(\theta),$$4

and choosing

$$\mathbb{E}_{\mathcal{B}_{\mathrm{res}}}\!\left[\mathbb{E}_{I,J}[g_{I,J}(\theta)]\right] = \int_{\mathcal X} \big(\mathcal{L}[\tilde u_\theta]-R\big)\,\nabla_\theta \mathcal{L}[\tilde u_\theta]\,dx \equiv \nabla_\theta \mathcal J(\theta),$$5

gives finite variance and finite expected cost per estimator draw (Goda et al., 2022). In distributionally robust learning, finite expected cost holds when $$\mathbb{E}_{\mathcal{B}_{\mathrm{res}}}\!\left[\mathbb{E}_{I,J}[g_{I,J}(\theta)]\right] = \int_{\mathcal X} \big(\mathcal{L}[\tilde u_\theta]-R\big)\,\nabla_\theta \mathcal{L}[\tilde u_\theta]\,dx \equiv \nabla_\theta \mathcal J(\theta),$$6, bounded variance is proved for $$\mathbb{E}_{\mathcal{B}_{\mathrm{res}}}\!\left[\mathbb{E}_{I,J}[g_{I,J}(\theta)]\right] = \int_{\mathcal X} \big(\mathcal{L}[\tilde u_\theta]-R\big)\,\nabla_\theta \mathcal{L}[\tilde u_\theta]\,dx \equiv \nabla_\theta \mathcal J(\theta),$$7, and the experiments use the Giles-optimal choice $$\mathbb{E}_{\mathcal{B}_{\mathrm{res}}}\!\left[\mathbb{E}_{I,J}[g_{I,J}(\theta)]\right] = \int_{\mathcal X} \big(\mathcal{L}[\tilde u_\theta]-R\big)\,\nabla_\theta \mathcal{L}[\tilde u_\theta]\,dx \equiv \nabla_\theta \mathcal J(\theta),$$8 to balance expected cost and variance (Ghosh et al., 2020).

The zeroth-order literature makes the variance–cost tradeoff unusually explicit. For the estimator family $$\mathbb{E}_{\mathcal{B}_{\mathrm{res}}}\!\left[\mathbb{E}_{I,J}[g_{I,J}(\theta)]\right] = \int_{\mathcal X} \big(\mathcal{L}[\tilde u_\theta]-R\big)\,\nabla_\theta \mathcal{L}[\tilde u_\theta]\,dx \equiv \nabla_\theta \mathcal J(\theta),$$9, the practical $\mathcal{B}_{\mathrm{res}}$0 and $\mathcal{B}_{\mathrm{res}}$1 constructions admit finite-variance bounds, while the one-evaluation $\mathcal{B}_{\mathrm{res}}$2 construction can have infinite variance. The paper derives variance bounds in terms of the perturbation schedule $\mathcal{B}_{\mathrm{res}}$3 and the level distribution $\mathcal{B}_{\mathrm{res}}$4, and then proves that SGD with the unbiased estimators $\mathcal{B}_{\mathrm{res}}$5, $\mathcal{B}_{\mathrm{res}}$6, achieves the optimal zeroth-order complexity for smooth non-convex problems (Ma et al., 22 Oct 2025).

In overparameterized stochastic dynamical models, the principal scaling issue is not sample count but parameter dimension. The generator gradient estimator is attractive because the number of auxiliary pathwise-differentiation estimators scales with the state dimension $\mathcal{B}_{\mathrm{res}}$7, not with the number of parameters $\mathcal{B}_{\mathrm{res}}$8. The note emphasizes that runtime therefore remains stable as the model becomes overparameterized, which is the same scaling advantage associated with adjoint methods for ODEs and PDEs (Badolle et al., 2024).

A complementary theoretical perspective is provided by the analysis of doubly stochastic SGD in the “finite sum with infinite data” regime. Under a fixed per-iteration computational budget $\mathcal{B}_{\mathrm{res}}$9, where $I,J$0 is minibatch size and $I,J$1 is the number of Monte Carlo samples, the paper’s analysis suggests where one should invest most of the budget. Its practical conclusion is that, for a fixed budget $I,J$2, it is generally better to increase $I,J$3 than $I,J$4, and that random reshuffling improves the complexity dependence on subsampling noise from $I,J$5 to $I,J$6, although the within-epoch gradients under reshuffling are conditionally biased (Kim et al., 2024).

6. Applications, assumptions, and recurring misconceptions

DS-UGE methods are now associated with several distinct application areas. In high-dimensional PINNs, they enable exact macroscopic conservation through projection without retaining a full quadrature graph in memory (Liang, 31 Mar 2026). In conditional stochastic optimization, they make standard stochastic-approximation theory available for nested expectations by removing the bias of naive nested Monte Carlo (Goda et al., 2022). In continuous-time state-space models, they support gradient-based parameter estimation, including stochastic-gradient Langevin descent, by delivering an unbiased score even when the hidden diffusion must be discretized (Ballesio et al., 2021). In distributionally robust learning, they replace biased fixed-minibatch robust-gradient estimates by unbiased multilevel corrections suitable for outer SGD (Ghosh et al., 2020). In zeroth-order optimization, they provide truly unbiased gradient surrogates using only function evaluations, with experiments that include language-model fine-tuning (Ma et al., 22 Oct 2025). In missing-data learning, SDRG adapts doubly robust causal-estimation logic to stochastic optimization under covariate-dependent missingness (Lee et al., 2018).

Several misconceptions recur in this literature. First, doubly stochastic does not by itself imply unbiased. SDIFP explicitly distinguishes its unbiased DS-UGE from the faster $I,J$7 “Sampling Only Once” variant, which introduces local variational bias because the forward and backward subsets are no longer independent (Liang, 31 Mar 2026). Likewise, the general theory of doubly stochastic gradients treats independent-minibatch estimators as unbiased but random-reshuffling estimators as conditionally biased within an epoch (Kim et al., 2024). Second, the term DS-UGE is not universal; some papers formulate the same design principle under names tied to their application domain, such as stochastic doubly robust gradient or zero-bias estimator (Lee et al., 2018, Ghosh et al., 2020).

The assumptions required for rigorous unbiasedness are also domain-specific and often nontrivial. SDIFP relies on independence of sampling sources, unbiased operator subsampling, integrability sufficient for Fubini–Tonelli, and unbiased batch-based moment-gradient estimators (Liang, 31 Mar 2026). Conditional stochastic optimization requires Lipschitz and Hölder conditions on $I,J$8, plus moment bounds on $I,J$9 and $$F(x)=\mathbb{E}_{\xi}\!\left[f_{\xi}\!\left(\mathbb{E}_{\eta\mid \xi}[g_{\eta}(x,\xi)]\right)\right],$$0 (Goda et al., 2022). The continuous-time state-space paper expects the required CCPF coupling and summability conditions to follow from related work, but defers the full proof (Ballesio et al., 2021). The adjoint interpretation of the generator gradient estimator is deliberately presented as a conceptual equivalence and does not provide a new variance bound or convergence theorem (Badolle et al., 2024). SDRG is explicit that a full convergence theory is left for future work (Lee et al., 2018).

Taken together, these results establish DS-UGE not as a single algorithm but as a research pattern: identify two stochastic bottlenecks, couple or randomize them so that expectations telescope or cancel correctly, and obtain an estimator that is unbiased for the exact gradient while remaining computationally viable. The main open questions, as reflected in the current literature, concern sharper variance theory, stronger non-asymptotic convergence guarantees, and principled design of the sampling laws or couplings that optimize the cost–variance tradeoff in each application domain.