Summary: Analysis and development of restarted Krylov subspace methods for computing \(f(A){\boldsymbol b}\) have proliferated in recent years. We present an acceleration technique for such methods when applied to Stieltjes functions \(f\) and Hermitian positive definite matrices \(A\). This technique is based on a rank-one modification of the Lanczos matrix derived from a connection between the Lanczos process and Gauss-Radau quadrature. We henceforth refer to the technique paired with the standard Lanczos method for matrix functions as the Radau-Lanczos method for matrix functions. We develop properties of general rank-one updates, leading to a framework through which other such updates could be explored in the future. We also prove error bounds for the Radau-Lanczos method, which are used to prove the convergence of restarted versions. We further present a thorough investigation of the Radau-Lanczos method explaining why it routinely improves over the standard Lanczos method. This is confirmed by several numerical experiments, and we conclude that, in practical situations, the Radau-Lanczos method is superior in terms of iteration counts and timings, when compared to the standard Lanczos method.