Probabilistic error analysis of limited-precision stochastic rounding

TL;DR

The paper develops a hardware-faithful model of limited-precision stochastic rounding, $\text{SR}_{p,r}$, showing that its bias depends on the random-bit budget $r$ and vanishes as $r\to\infty$, recovering the classical $\text{SR}_p$ behavior. It then constructs a probabilistic error-analysis framework that yields bounds of the form $\mathcal{O}(\sqrt{n}u_p) + \mathcal{O}(n u_{p+r})$ for fundamental linear-algebraic routines, such as recursive summation and inner products, with both mean-bias and probabilistic high-probability guarantees derived via martingales and Bienaymé–Chebyshev methods. The results are validated by numerical experiments on summation, the Rosenbrock function, and neural-network parameter updates, demonstrating that an appropriately chosen $r$ (roughly $\lceil(\log_2 n)/2\rceil$) yields substantial accuracy improvements with manageable hardware cost. Overall, the work provides a principled bridge between ideal SR theory and practical, bit-budget-constrained implementations, guiding practitioners in selecting $r$ for reliable stochastic rounding in large-scale computations.

Abstract

Classical probabilistic rounding error analysis is particularly well suited to stochastic rounding (SR), and it yields strong results when dealing with floating-point algorithms that rely heavily on summation. For many numerical linear algebra algorithms, one can prove probabilistic error bounds that grow as O(nu), where n is the problem size and u is the unit roundoff. These probabilistic bounds are asymptotically tighter than the worst-case ones, which grow as O(nu). For certain classes of algorithms, SR has been shown to be unbiased. However, all these results were derived under the assumption that SR is implemented exactly, which typically requires a number of random bits that is too large to be suitable for practical implementations. We investigate the effect of the number of random bits on the probabilistic rounding error analysis of SR. To this end, we introduce a new rounding mode, limited-precision SR. By taking into account the number r of random bits used, this new rounding mode matches hardware implementations accurately, unlike the ideal SR operator generally used in the literature. We show that this new rounding mode is biased and that the bias is a function of r. As r approaches infinity, however, the bias disappears, and limited-precision SR converges to the ideal, unbiased SR operator. We develop a novel model for probabilistic error analysis of algorithms employing SR. Several numerical examples corroborate our theoretical findings.

Figures (4)

• Figure 1: Quantities used in the definitions \ref{['eq:sr']} and \ref{['eq:sr-imp']}.
• Figure 2: Left: relative error of RN and $\text{SR}_{11,r}$ in IEEE-754 binary16 arithmetic ($p=11$) for the recursive summation of $n$ floating-point numbers drawn from a uniform distribution between 0 and 1. For each value of $n$, the reported relative error for $\text{SR}_{11,r}$ is the average value over $500$ runs. Right: comparison of deterministic bound and probabilistic bounds (\ref{['thm:bc-proba-sum']}) with the associated random bits.
• Figure 3: Convergence profiles for 5,000 iterations of gradient descent on the Rosenbrock function. The parameter updates are performed using binary64 arithmetic with RN and binary16 arithmetic with RN and $\text{SR}_{11,r}$. The starting value for the iteration is $\mathbf{x}_0 = (0, 0)$ for the profiles on the left and $\mathbf{x}_0 = (0.5, 0.5)$ for those on the right. For both experiments, we average each $\text{SR}_{11,r}$ error over $500$ different runs, and the learning rate is $t_k=0.001$.
• Figure 4: Evolution of training loss (left) and validation accuracy (right) for a ResNet32 image classification network on the CIFAR-10 dataset using various parameter and optimizer update configurations. In the baseline configuration,binary32 arithmetic with RN is used for compute and the same format is used for storage. For the low-precision configurations, parameters are stored and updated using bfloat16 arithmetic with either RN or $\text{SR}_{p,r}$.

Theorems & Definitions (39)

• Lemma 1: Bienaymé--Chebyshev inequality
• Definition 1: Mean independence
• Proposition 1
• Proposition 2
• Definition 2: miup05
• Lemma 2: Azuma--Hoeffding inequality, miup05
• Definition 3: Stochastic rounding
• Remark 1
• Definition 4: limited-precision stochastic rounding
• Remark 2