Complexity and Last-Iterate Convergence of Inexact Riemannian Proximal Momentum Shuffling with Variance Reduction

We develop a pathwise convergence and complexity theory for inexact Riemannian proximal methods that combine full-pass shuffling, variance reduction, projected momentum, and inexact tangent-space proximal solves. The main difficulty is that gradient estimators within a shuffled epoch are dependent, so standard one-step conditional variance recursions do not directly apply. We introduce an epochwise estimator error–memory condition and a projected-momentum lifting principle, and verify them for projection-based SVRG, projection-based SARAH/SPIDER, projected ambient SAGA, and projected ambient SAG. Combined with a prescribed Fenchel-dual residual criterion for the tangent proximal subproblems, this yields a deterministic epochwise descent estimate that controls estimator error, momentum memory, shuffling, and inner inexactness simultaneously. For fixed finite-sum size, we obtain $O({\epsilon }^{-2})$ epoch and retraction complexity and total proximal-evaluation complexities of $O({\epsilon }^{-4})$, $O({\epsilon }^{-3}\log {\epsilon }^{-1})$, and $O({\epsilon }^{-3})$ for FISTA, NFG, and accumulative-regularization inner solvers, respectively. Unlike stochastic analyses that guarantee stationarity in expectation, these bounds hold along every admissible shuffling trajectory. We further establish epochwise nonmonotone descent and relative-error estimates for the natural tangent-bundle lift of the objective. Under its ordinary Kurdyka–\Lojasiewicz property, the complete within-epoch trajectory has finite length and converges to a single critical point, without requiring a common deterministic limiting objective value across trajectories. Under a KL exponent, we obtain explicit last-iterate rates and a pathwise $O({\epsilon }^{-1})$ best-iterate outer complexity.