We develop a sketch-based factor reduction and a Nesterov-accelerated projected gradient algorithm (NPGA) with GPU acceleration, yielding a doubly accelerated solver for large-scale constrained mean-variance portfolio optimization. Starting from the sample covariance factor L, the method combines randomized subspace embedding, spectral truncation, and ridge stabilization to construct an effective factor Leff. It then solves the resulting constrained problem with a structured projection computed by scalar dual search and GPU-friendly matrix-vector kernels, yielding one computational pipeline for the baseline, sketched, and Sketch-Truncate-Ridge (STR)-regularized models. We also establish approximation, conditioning, and stability guarantees for the sketching and STR models, including explicit O(ε) bounds for the covariance approximation, the optimal value error, and the solution perturbation under (ε,δ)-subspace embeddings. Experiments on synthetic and real equity-return data show that the method preserves objective accuracy while reducing runtime substantially. On a 5440-asset real-data benchmark with 48374 training periods, NPGA-GPU solves the unreduced full model in 2.80 seconds versus 64.84 seconds for Gurobi, while the optimized compressed GPU variants remain in the low-single-digit-second regime. These results show that the full dense model is already practical on modern GPUs and that, after compression, the remaining bottleneck is projection rather than matrix-vector multiplication.
Figures & tables
Figure 1: Eigenvalue spectrum of the real-world covariance matrix. The spectrum decays rapidly for the first few dozen components before entering a long, slowly decaying bulk spanning several thousand dimensions.
Figure 2: Cumulative explained variance E(r) . Due to the extremely long noisy bulk, 80% of the total variance is reached only after more than 2500 eigenvalues, demonstrating the inadequacy of energy-based rank selection for real-world covariance matrices.
Model
Leff
γ
Cost / iter
Baseline ( 2 )
L∈Rn×T
0
O(nT)
Sketch ( 3 )
L~∈Rn×s
0
O(ns)
STR ( 4 )
L~ℓ∈Rn×ℓ
>0
O(nℓ)
Table 1: Per-iteration complexity of NPGA.
Figure 3: Approximation quality on synthetic instances as the sketch size grows. Both metrics improve as s/ℓ increases, and JL and CountSketch behave similarly once the retained energy level is fixed.
Sketch
Rel. spec. err.
Rel. obj. gap
JL
0.1212
0.0736
CountSketch
0.1192
0.0731
Table 2: Synthetic approximation quality at η=0.98 and s/ℓ=2 .
Figure 4: Synthetic solver behavior. The convex case follows the k−2 benchmark, the ridge-stabilized case is close to the linear-rate benchmark, and GPU acceleration becomes substantial on large unreduced factors.
n
Gurobi (s)
NPGA-GPU (s)
NPGA gap
STR-NPGA-GPU (s)
STR gap
500
0.16
0.47
0.0006
0.59
0.0014
1000
0.56
0.37
0.0012
0.50
0.0025
2000
2.59
0.33
0.0023
0.42
0.0047
4000
19.69
0.64
0.0065
0.42
0.0108
Table 3: Synthetic solver benchmark against Gurobi on the same training problems.
n
Gurobi (s)
NPGA-CPU (s)
NPGA gap
STR-NPGA-CPU (s)
STR gap
STR annualized return (%)
STR annualized volatility (%)
300
0.22
1.28
0.0025
0.03
0.1348
36.15
14.66
600
0.62
2.10
0.0053
0.07
0.0937
39.30
14.26
1000
1.78
4.22
0.0045
0.12
0.0675
39.99
12.79
Table 4: Subset benchmarks on the real 5-minute A-share data.
Method
Factor cols
Build (s)
Solve (s)
Total (s)
Full-model gap
Annualized return (%)
Annualized volatility (%)
Gurobi
48374
0.00
64.84
64.84
0.0000
36.04
11.00
NPGA-CPU
48374
0.00
26.77
26.77
0.0095
36.02
11.12
NPGA-GPU
48374
0.00
2.80
2.80
0.0080
36.04
11.12
Sketch-NPGA-GPU
948
1.15
2.53
3.68
0.2380
35.09
11.32
STR-NPGA-CPU
474
1.15
0.83
1.99
0.3203
35.24
11.17
STR-NPGA-GPU
474
1.15
2.16
3.31
0.3203
35.24
11.17
Table 5: Full-scale benchmark on the 5-minute A-share data with n=5440 assets and T=48374 periods.
Appendix figures & tables3 assets
Supplementary material from the paper’s appendix.
Appendix
Figure 5: The Sketch–Truncate–Ridge (STR) pipeline for constructing a stable covariance approximation.
Figure 6: Condition number and smallest-eigenvalue uplift under ridge regularization on a representative synthetic instance. The Sketch and STR curves coincide on this instance and therefore appear as a single trajectory.
Figure 7: Convergence behavior of the GPU variants on the full 5440-asset real-data problem. The three GPU variants require nearly the same number of outer iterations, whereas their wall-clock times differ substantially.
Randomized-subspace methods reduce the cost of first-order optimization by using only low-dimensional projected-gradient information, a feature that is attractive in forward-mode automatic differentiation and communication-limited settings. While Nesterov acceleration is well understood for full-gradient and coordinate-based methods, obtaining accelerated methods for general subspace sketches that use only projected-gradient information and can improve over full-dimensional Nesterov acceleration in oracle complexity is technically nontrivial. We develop randomized-subspace Nesterov accelerated gradient methods for smooth convex and smooth strongly convex optimization under matrix smoothness and generic sketch moment assumptions. The key technical ingredient is a three-sequence formulation tailored to matrix smoothness, which recovers the corresponding classical Nesterov methods in the full-dimensional case. The resulting theory establishes accelerated oracle-complexity guarantees and makes explicit how matrix smoothness and the sketch distribution enter the complexity. It also provides a unified basis for comparing sketch families and identifying when randomized-subspace acceleration improves over full-dimensional Nesterov acceleration in oracle complexity.
Gaku Omiya, Pierre-Louis Poirion, Akiko Takeda
Department of Mathematical Informatics, The University of Tokyo, Tokyo, Japan · Center for Advanced Intelligence Project, RIKEN, Tokyo, Japan
We study an online sketched Newton method that approximates the Newton direction at each step via a state-of-the-art sketching solver, called the generalized accelerated sketch-and-project solver (GAS), thereby mitigating the computational bottleneck of classical second-order methods. The GAS solver improves upon vanilla, unaccelerated sketch-and-project solvers by achieving accelerated convergence through Nesterov momentum updates, and accommodates a flexible projection metric whose proper choice further reduces computational cost. Building on this design, we establish asymptotic normality of the averaged sketched Newton iterates and characterize their limiting covariance matrix. The resulting covariance recovers that of the unaccelerated sketched Newton method under a specific choice of acceleration parameters, converges more rapidly (in the number of sketching steps) to the minimax-optimal covariance in general, and is smaller than that of the last iterate produced by the accelerated method. Finally, we strengthen these results by establishing a functional central limit theorem for the Newton iterates, which allows us to bypass explicit covariance estimation and develop an online inference procedure based on random scaling. Specifically, we construct a pivotal test statistic by appropriately rescaling the averaged iterates, so that its limiting distribution is free of any unknown parameters, enabling asymptotically valid online inference. Numerical experiments demonstrate superior performance of the proposed inference procedure.
Xinchen Du, Elizaveta Rebrova, Michał Dereziński +1
School of Industrial and Systems Engineering, Georgia Institute of Technology · Department of Operations Research and Financial Engineering, Princeton University
We investigate mean-variance portfolio selection with an ℓ0-penalty to promote sparsity in asset allocations. Uncertainty in the mean return vector is incorporated through an ellipsoidal uncertainty set, yielding a robust sparse optimization framework. We characterize the structure of both local and global minimizers and exploit these properties in the risk minimization and return maximization formulations. Building on this structural insight, we develop a branch-and-bound algorithm tailored to the resulting robust sparse portfolio problems, together with a new pruning rule that can discard exponentially many candidate portfolios in a single step. Extensive computational experiments on real market data, together with comparisons against a mixed-integer second-order cone programming solver, demonstrate the effectiveness and competitiveness of the proposed approach.
Deniz Akkaya, Emre Can Yayla, Buse Şen +1
Department of Industrial Engineering, Bilkent University, Ankara, Türkiye · Risk Analytics and Optimization Chair, EPFL, Lausanne, Switzerland