We consider a certain class of Riemannian submersions π:N→M and study lifted geodesic random walks from the base manifold M to the total manifold N. Under appropriate conditions on the distribution of the speed of the geodesic random walks, we prove an invariance principle, i.e., convergence to horizontal Brownian motion for the lifted walks. This gives us a natural probabilistic proof of the geometric identity relating the horizontal Laplacian ΔH on N and the Laplace–Beltrami operator ΔM on M. In the setting where N is the orthonormal frame bundle O(M), this identity is central in the Malliavin–Eells–Elworthy construction of Riemannian Brownian motion.

We extend Loeper's (Formula presented.) -estimate (Theorem 2.9 in J. Math. Pures Appl. (9) 86 (2006), no. 1, 68–79) relating the force fields to the densities for the Vlasov–Poisson system to (Formula presented.), with (Formula presented.), based on the Helmholtz–Weyl decomposition. This allows us to generalize both the classical Loeper's 2-Wasserstein stability estimate (Theorem 1.2 in J. Math. Pures Appl. (9) 86 (2006), no. 1, 68–79) and the recent stability estimate by the first author relying on the newly introduced kinetic Wasserstein distance (Theorem 3.1 in Arch Rational Mech. Anal. 244 (2022), no. 1, 27–50) to kinetic Wasserstein distances of order (Formula presented.).