  We develop new algorithms for Riemannian bilevel optimization. We focus in
particular on batch and stochastic gradient-based methods, with the explicit
goal of avoiding second-order information such as Riemannian hyper-gradients.
We propose and analyze $\mathrm{RF^2SA}$, a method that leverages first-order
gradient information to navigate the complex geometry of Riemannian manifolds
efficiently. Notably, $\mathrm{RF^2SA}$ is a single-loop algorithm, and thus
easier to implement and use. Under various setups, including stochastic
optimization, we provide explicit convergence rates for reaching
$\epsilon$-stationary points. We also address the challenge of optimizing over
Riemannian manifolds with constraints by adjusting the multiplier in the
Lagrangian, ensuring convergence to the desired solution without requiring
access to second-order derivatives.
