Online Covariance Estimation for Averaged Stochastic Gradient Descent on Riemannian Manifolds

Jingrong Guan; Liwei Jiang; Sen Na

stat.MLOct 24, 2026

Online Covariance Estimation for Averaged Stochastic Gradient Descent on Riemannian Manifolds

Authors: Jingrong Guan, Liwei Jiang, Sen Na

Organizations: Department of Industrial Engineering and Management Sciences, Northwestern University · Edwardson School of Industrial Engineering, Purdue University · School of Industrial and Systems Engineering, Georgia Institute of Technology

Abstract

We study online inference for averaged Riemannian stochastic gradient descent (RSGD) when the unknown parameter lies on a manifold. Although the averaged RSGD iterate is asymptotically normal, practical inference requires estimating the unknown covariance operator in its limiting distribution. Existing online batch-means methods cannot be applied directly because a general Riemannian manifold lacks an addition operation, and the tangent space containing the Riemannian gradient varies across iterations. Through vector transport, we construct overlapping and non-overlapping covariance estimators that adaptively combine information from successive batches while accounting for the changing tangent spaces. The proposed estimators are online, require neither Riemannian Hessian computations nor a prespecified terminal sample size, and attain the same convergence rates as their Euclidean counterparts (Zhu et al., 2023). Together with asymptotic normality, these estimators enable online statistical inference for RSGD on manifolds. We illustrate the proposed framework through a detailed study of streaming principal component analysis on Grassmann manifolds.

Figures & tables

Figure 1: Left: Assumption 2.1 (a). The inverse image

R_{x}^{-1}(\mathcal{X})

lies in the tangent ball

B_{x}(0_{x},\rho_{\mathrm{ret}})

for every

x\in\mathcal{X}

, with a common radius

\rho_{\mathrm{ret}}>\tan(\max_{x,y\in\mathcal{X}}d(x,y))

. Right: Assumption 2.1 (b). For any

x,y\in\mathcal{X}

, the normalization retraction maps the line segment

\{tR_{x}^{-1}(y):t\in[0,1]\}

\mathrm{T}_{x}\mathcal{M}

onto the short great circle arc from

x

y

, which remains in

\mathcal{X}

, illustrating its retraction convexity.

Figure 2: Construction of a confidence region on the manifold. The retraction

R_{\widetilde{x}_{n}}

maps the tangent ellipsoid

\mathcal{E}_{n,\diamond}(1-\delta)

to the shaded region

\mathcal{C}_{n,\diamond}(1-\delta)

, and its center

0_{\widetilde{x}_{n}}

\widetilde{x}_{n}

. The right panel separates the tangent space from the manifold to show corresponding points, connected by the light gray dashed lines. The estimator subscript

\diamond

is omitted in the figure.

Figure 3: Spherical representation of a confidence region on

G_{3,1}

. Left: the black and orange lines represent the estimate

[\widetilde{X}_{n}]

and a subspace

[z]

in the confidence region, respectively. Right: the polar retraction maps each tangent vector

v

to the line represented by

z=(\widetilde{X}_{n}+v)/\sqrt{1+\|v\|^{2}}

. These unit representatives form the shaded spherical cap depicting the confidence region. The subscript

\diamond

is omitted.

Figure 4: Gaussian streaming PCA: log mean covariance error (operator norm) against log sample size. Columns correspond to

G_{3,1}

G_{7,2}

, and

G_{10,3}

, and the top and bottom rows use eigengaps

1

and

0.25

, respectively. The colored curves compare OL and NOL with

c\in\{1,2,5\}

. Each mean is based on

200

independent runs. The black dashed lines indicate the reference slope

-1/8

Figure 5: Gaussian streaming PCA: empirical coverage (top) and mean length (bottom) of

95\%

confidence intervals for

h_{q}(x_{\star})

. Columns correspond to

G_{10,3}

and

G_{20,4}

, with

\Delta_{\lambda}=0.25

and

c=1

. The blue, orange, and purple curves show OL, NOL, and Oracle, respectively, with results based on

200

independent runs. The small plots inside the top panels show the later coverage values, and the red dashed lines mark the nominal coverage level

0.95

Empirical coverage
$G_{10,3}$ ( $d=21$ )	$n=50{,}000$	$n=100{,}000$	$n=200{,}000$	$n=300{,}000$
OL	0.835 (0.02625)	0.905 (0.02073)	0.910 (0.02024)	0.895 (0.02168)
NOL	0.860 (0.02454)	0.905 (0.02073)	0.930 (0.01804)	0.920 (0.01918)
Oracle	0.930 (0.01804)	0.965 (0.01300)	0.955 (0.01466)	0.945 (0.01612)
$G_{20,4}$ ( $d=64$ )	$n=50{,}000$	$n=100{,}000$	$n=200{,}000$	$n=300{,}000$
OL	0.860 (0.02454)	0.860 (0.02454)	0.890 (0.02212)	0.895 (0.02168)

Table 1: Empirical coverage of nominal

95\%

confidence intervals for

h_{q}(x_{\star})

in Gaussian streaming PCA, with

\Delta_{\lambda}=0.25

and

c=1

. Results are based on 200 independent runs, with Monte Carlo standard errors in parentheses. Here

d=k(p-k)

is the dimension of

G_{p,k}

Appendix figures & tables3 assets

Supplementary material from the paper’s appendix.

Appendix

Figure 6: Gaussian streaming PCA: log mean covariance error (operator norm) against log sample size. Columns correspond to

G_{15,3}

G_{20,4}

, and

G_{25,5}

, and the top and bottom rows use eigengaps

1

and

0.25

, respectively. As in Figure 4 , darker and lighter curves represent OL and NOL, with blue, orange, and green indicating

c=1,2,5

, respectively. Each mean is based on

200

independent runs. The black dashed lines indicate the reference slope

-1/8

		Coverage			Mean length $\times 10^{3}$
$G_{p,k}$	$\Delta_{\lambda}$	OL	NOL	Oracle	OL	NOL	Oracle
$G_{3,1}$	1	0.920 (0.019)	0.920 (0.019)	0.980 (0.010)	4.425	4.447	4.953
$G_{3,1}$	0.25	0.920 (0.019)	0.925 (0.019)	0.980 (0.010)	11.305	11.722	13.424
$G_{7,2}$	1	0.915 (0.020)	0.910 (0.020)	0.960 (0.014)	3.267	3.252	3.488
$G_{7,2}$	0.25	0.915 (0.020)	0.895 (0.022)	0.955 (0.015)	4.581	4.628	5.012
$G_{10,3}$	1	0.915 (0.020)	0.925 (0.019)	0.955 (0.015)	2.224	2.273	2.415

Appendix

Table 2: Empirical coverage and mean length of nominal

95\%

confidence intervals for

h_{q}(x_{\star})

in Gaussian streaming PCA, with

n=3\times 10^{5}

and

c=1

. Results are based on

200

independent runs, with Monte Carlo standard errors in parentheses. Mean lengths are multiplied by

10^{3}

Figure 7: Gaussian streaming PCA: empirical coverage and mean length of

95\%

confidence intervals for

h_{q}(x_{\star})

, with

c=1

. Columns correspond to

G_{3,1}

and

G_{25,5}

; the upper and lower pairs of rows use eigengaps

1

and

0.25

, respectively, with coverage above mean length in each pair. The blue, orange, and purple curves show OL, NOL, and Oracle, respectively, with results based on

200

independent runs. The small plots inside the coverage panels show the later coverage values, and the red dashed lines mark the nominal coverage level

0.95

Covariance Estimation

Covariance estimation develops methods for learning dependence and uncertainty structure from noisy, incomplete, or high-dimensional data.

See all 61 papers

Online Covariance Estimation in Averaged SGD: Improved Batch-Mean Rates and Minimax Optimality via Trajectory Regression
Apr 12, 2026Yijin Ni, Xiaoming HuoPolyak-Ruppert Averaging Covariance Estimation
H. Milton Stewart School of Industrial and Systems Engineering, Georgia Institute of Technology
Refining Covariance Matrix Estimation in Stochastic Gradient Descent Through Bias Reduction
Apr 23, 2026Ziyang Wei, Wanrong Zhu, Jingyang Lyu +1Covariance Estimation Stochastic Gradient Descent
Department of Statistics University of Chicago Chicago, IL 60637, USA · Department of Statistics University of California, Irvine Irvine, CA 92617, USA · Department of Statistics University of Wisconsin–Madison Madison, WI 53706, USA
Inversion-Free Natural Gradient Descent on Riemannian Manifolds
Apr 3, 2026Dario Draca, Takuo Matsubara, Minh-Ngoc TranNatural Gradient Methods Riemannian Optimization
School of Mathematics and Statistics, University of Sydney, Australia · University of Sydney Business School, Australia

Apr 12, 2026J/K move · Enter open · S save

Online Covariance Estimation in Averaged SGD: Improved Batch-Mean Rates and Minimax Optimality via Trajectory Regression

Yijin Ni, Xiaoming Huo

H. Milton Stewart School of Industrial and Systems Engineering, Georgia Institute of Technology

Polyak-Ruppert Averaging Covariance Estimation Minimax Estimation Stochastic Gradient Descent

Open paper AskPDF arXiv

We study online covariance matrix estimation for Polyak--Ruppert averaged stochastic gradient descent (SGD). The online batch-means estimator of Zhu, Chen and Wu (2023) achieves an operator-norm convergence rate of

O(n^{-(1-α)/4})

, which yields

O(n^{-1/8})

at the optimal learning-rate exponent

α\rightarrow 1/2^+

. A rigorous per-block bias analysis reveals that re-tuning the block-growth parameter improves the batch-means rate to

O(n^{-(1-α)/3})

, achieving

O(n^{-1/6})

. The modified estimator requires no Hessian access and preserves

O(d^2)

memory. We provide a complete error decomposition into variance, stationarity bias, and nonlinearity bias components. A weighted-averaging variant that avoids hard truncation is also discussed. We establish the minimax rate

Θ(n^{-(1-α)/2})

for Hessian-free covariance estimation from the SGD trajectory: a Le Cam lower bound gives

Ω(n^{-(1-α)/2})

, and a trajectory-regression estimator--which estimates the Hessian by regressing SGD increments on iterates--achieves

O(n^{-(1-α)/2})

, matching the lower bound. The construction reveals that the bottleneck is the sublinear accumulation of information about the Hessian from the SGD drift.

Abstract

Covariance Estimation

Explore similar work

Online Covariance Estimation in Averaged SGD: Improved Batch-Mean Rates and Minimax Optimality via Trajectory Regression

Refining Covariance Matrix Estimation in Stochastic Gradient Descent Through Bias Reduction

Inversion-Free Natural Gradient Descent on Riemannian Manifolds