# Fast Controllable Diffusion Models for Undersampled MRI Reconstruction

###### Abstract

Supervised deep learning methods have shown promise in undersampled Magnetic
Resonance Imaging (MRI) reconstruction, but their requirement for paired data
limits their generalizability to the diverse MRI acquisition parameters.
Recently, unsupervised controllable generative diffusion models have been
applied to undersampled MRI reconstruction, without paired data or model
retraining for different MRI acquisitions. However, diffusion models are
generally slow in sampling and state-of-the-art acceleration techniques can
lead to sub-optimal results when directly applied to the controllable
generation process. This study introduces a new algorithm called Predictor-
Projector-Noisor (PPN), which enhances and accelerates controllable generation
of diffusion models for undersampled MRI reconstruction. Our results
demonstrate that PPN produces high-fidelity MR images that conform to
undersampled k-space measurements with significantly shorter reconstruction
time than other controllable sampling methods. In addition, the unsupervised
PPN accelerated diffusion models are adaptable to different MRI acquisition
parameters, making them more practical for clinical use than supervised
learning techniques.

Index Terms—  Undersampled MRI Reconstruction, Deep Learning, Diffusion
Models, Fast Controllable Generation, Inverse Problems.

##  1 Introduction

Deep learning has emerged as a valuable tool for undersampled Magnetic
Resonance Imaging (MRI) reconstruction [1, 2, 3, 4]. However, most of these
approaches are supervised, relying on paired data to learn a mapping between
undersampled k-space measurements and the fully-sampled MR images. Moreover,
supervised methods have poor generalizability when applied to different
undersampling conditions, making them less practical for clinical use.

To address these limitations, Song et al. [5] developed a fully unsupervised
method that can reconstruct MRI and Computed Tomography images from
undersampled signal acquisitions with similar performance but significantly
enhanced generalizability to supervised methods. The idea is to train
unsupervised generative diffusion models [6, 7, 8] using the MR images only,
and then use physical model to guide the sampling process to generate MR
images consistent with the k-space measurements. However, diffusion models in
general come at the cost of slow image generation speed.

Fig. 1: Illustration of a proposed PPN generation step.

Various methods have been proposed to accelerate the diffusion models,
including Denoising Diffusion Implicit Models (DDIM) [9]. Although successful
in speeding up unconditional image generation tasks, they lead to sub-optimal
results when directly applied to accelerate the controllable generation
process. This is because large reverse steps during the controllable sampling
process can cause deviations from the manifolds [10], particularly when
projection algorithms are used to enforce k-space data consistency [5, 11]. It
is well-known that sampling from out-of-distribution regions is challenging
due to the lack of score information in the ambient space [12]. Consequently,
this deviation from the manifolds leads to an inefficient sampling trajectory,
resulting in degraded reconstruction results.

We present a novel method called Predictor-Projector-Noisor (PPN) to
accelerate the controllable generation process of diffusion models in
undersampled MRI reconstruction. Our work follows the fully unsupervised
framework of [5], but contributes a more efficient high-fidelity generative
strategy. Specifically, during each generative step, the PPN method projects
the DDIM-predicted noise-free image to match the k-space measurement and then
noises it back to a higher noise level. This approach ensures the reverse
trajectory stays within a chain of noisy manifolds converging to the data
manifold, enabling fast generation without sacrificing image reconstruction
accuracy. Our experiments demonstrate that PPN outperforms other controllable
generative processes in reconstructing high-fidelity MR images with
significantly improved speed, and shows excellent generalizability across
various under-sampling conditions without model retraining, due to their
unsupervised learning nature.

##  2 methodology

###  2.1 Preliminary: Diffusion Models

Diffusion models (DMs) are a class of deep generative neural networks designed
to sample from an unknown data distribution $p_{\text{data}}$. They consist of
two processes: a forward process that gradually perturbs inputs with noise to
arrive at a tractable distribution (i.e. ${\cal N}({\bf 0},{\bf I})$), and a
reverse process that learns to undo the forward process by denoising. A well-
known DM implementation is Denoising Diffusion Probabilistic Models (DDPM)
[6], while Denoising Diffusion Implicit Models (DDIM) [9] accelerate the
denoising process.

DDPM The forward process of DDPM [6] can be specified as a $T$-step Markov
chain with a transition kernel $q({\bf x}_{t}|{\bf x}_{t-1})={\cal
N}\left({\bf x}_{t}|\sqrt{\alpha_{t}}{\bf x}_{t-1},(1-\alpha_{t}){\bf
I}\right)$, where ${\bf x}_{0}\sim p_{\text{data}}$ and $\alpha_{t}\in(0,1)$
for all $t=1,\ldots,T$ are parameters chosen [6] or trained [7] to ensure that
the distribution of ${\bf x}_{T}$ converges to ${\cal N}({\bf 0},{\bf I})$.
Notably, diffused samples can be generated directly from a marginal
distribution $q({\bf x}_{t}|{\bf x}_{0})={\cal N}\left({\bf
x}_{t};\sqrt{\bar{\alpha}_{t}}{\bf x}_{0},(1-\bar{\alpha}_{t}){\bf I}\right)$,
where $\bar{\alpha}_{t}=\prod_{s=1}^{t}\alpha_{s}$. We can reparameterize it
as:

${\bf x}_{t}=\sqrt{\bar{\alpha}_{t}}{\bf
x}_{0}+\sqrt{1-\bar{\alpha}_{t}}\bm{\epsilon}_{t}$ (1)

where $\bm{\epsilon}_{t}\sim{\cal N}({\bf 0},{\bf I})$. The reverse process of
DDPM utilizes a deep neural network, typically a U-Net, which is trained to
predict the noise added at ${\bf x}_{0}$ given ${\bf x}_{t}$ and time $t$.
This predicted noise is denoted as $\bm{\epsilon}_{\bm{\theta}}({\bf
x}_{t},t)$. One denoising step in DDPM is:

${\bf x}_{t-1}=\frac{1}{\sqrt{\alpha_{t}}}\left({\bf
x}_{t}-\frac{1-\alpha_{t}}{\sqrt{1-\bar{\alpha}_{t}}}\bm{\epsilon}_{\bm{\theta}}({\bf
x}_{t},t)\right)+\sqrt{1-\alpha_{t}}\bm{\epsilon}_{t}$ (2)

DDIM DDIM [9] introduces a non-Markovian process that matching the marginal
distribution $q({\bf x}_{t}|{\bf x}_{0})$ to accelerate the denoising process,
leading to a reformulation of 2 as:

${\bf x}_{t-1}=\sqrt{\bar{\alpha}_{t-1}}{\bf
x}_{0|t}+\sqrt{1-\bar{\alpha}_{t-1}-\sigma_{t}^{2}}\bm{\epsilon}_{\bm{\theta}}({\bf
x}_{t},t)+\sigma_{t}\bm{\epsilon}_{t}$ (3)

where $\bm{\epsilon}_{t}\sim{\cal N}({\bf 0},{\bf I})$ and

${\bf x}_{0|t}=\frac{{\bf
x}_{t}-\sqrt{1-\bar{\alpha}_{t}}\bm{\epsilon}_{\bm{\theta}}({\bf
x}_{t},t)}{\sqrt{\bar{\alpha}_{t}}}$ (4)

obtained by rewriting Eq. 1 with noise $\bm{\epsilon}_{t}$ replaced by
estimated $\bm{\epsilon}_{\bm{\theta}}({\bf x}_{t},t)$. The parameter
$\sigma_{t}$ controls the noise magnitude in $\bm{\epsilon}_{t}$, with DDIM
being deterministic when $\sigma_{t}=0$.

###  2.2 MRI Inverse Problems

In MRI inverse problems, we seek to reconstruct a signal
$\mathbf{x}\in\mathbb{R}^{n}$ from measurements
$\mathbf{y}\in\mathbb{R}^{m}=\mathbf{A}\mathbf{x}+\bm{\epsilon}$, where
$\bm{\epsilon}\in\mathbb{R}^{m}\sim\mathcal{N}(\mathbf{0},\sigma_{e}^{2}\mathbf{I})$
is noise, and $\mathbf{A}\in\mathbb{R}^{m\times n}=\mathbf{M}\mathbf{F}$
represents the forward model with $\mathbf{M}$ as the under-sampling mask and
$\mathbf{F}$ as the Fourier transform.

An inverse problem typically involves optimizing the cost function
$f(\tilde{\mathbf{x}})=\frac{1}{2\sigma_{e}^{2}}\|\mathbf{y}-\mathbf{A}\tilde{\mathbf{x}}\|_{2}^{2}+\bm{s}(\tilde{\mathbf{x}})$.
This combines a fidelity term, ensuring consistency with measurements, with a
prior as a penalty term to ensure $\tilde{\mathbf{x}}$ looks natural. In this
work, we use a pretrained diffusion model as the prior.

Tirer et al. [13] use the Plug-and-Play approach [14] to simplify inverse
problem solving into the Iterative Denoising and Backward Projections (IDBP)
algorithm. IDBP formulates the fidelity term as a projection of
$\tilde{\mathbf{x}}$ onto the affine subspace $\\{{\bf
x}\in\mathbb{R}^{n}|{\bf A}{\bf x}={\bf y}\\}$, and its closed-form solution
is $\tilde{\mathbf{x}}^{\prime}={\bf A}^{\dagger}{\bf y}+({\bf I}-{\bf
A}^{\dagger}{\bf A})\tilde{\mathbf{x}}$ where ${\bf A}^{\dagger}$ is the
pseudoinverse of ${\bf A}$. Using $\mathbf{A}=\mathbf{M}\mathbf{F}$ and
$\mathbf{A}^{\dagger}=\mathbf{F}^{-1}\mathbf{M}$, we can rewrite the closed-
form as:

${\bf P}_{{\bf y}}(\tilde{\mathbf{x}})={\bf F}^{-1}({\bf M}{\bf y}+({\bf
I}-{\bf M}){\bf F}\tilde{\mathbf{x}})$ (5)

This projection-based approach is adopted in Song et al. [5] (We will
henceforth refer to it as ’MedScore’) and DDNM [15]. The other category is
gradient-based approach, including DPS [10]. In this paper, we compare our PPN
with all three methods in section 4 and 5.

Algorithm 1 Predictor-Projector-Noisor (PPN)

1:$S$, ${\bf y}$ $\triangleright$ $S<T$

2:${\bf x}_{\text{zf}}\leftarrow{\bf F}^{-1}{\bf M}{\bf y}$ $\triangleright$
Zero-Filled

3:${\bf x}_{S}\leftarrow\sqrt{\bar{\alpha}_{S}}{\bf
x}_{\text{zf}}+\sqrt{1-\bar{\alpha}_{S}}\bm{\epsilon}_{S}$ $\triangleright$
$\bm{\epsilon}_{S}\sim\mathcal{N}(\bm{0},\mathbf{I})$

4:for $t=S$ to $1$ do

5:
$\mathbf{x}_{0|t}\leftarrow\left(\mathbf{x}_{t}-\sqrt{1-\bar{\alpha}_{t}}\bm{\epsilon}_{\bm{\theta}}(\mathbf{x}_{t},t)\right)/\sqrt{\bar{\alpha}_{t}}$
$\triangleright$ Prediction

6:
$\hat{\mathbf{x}}_{0|t}\leftarrow\mathbf{F}^{-1}(\mathbf{M}\mathbf{y}+(\mathbf{I}-\mathbf{M})\mathbf{F}\mathbf{x}_{0|t})$
$\triangleright$ Projection

7: $\bm{\epsilon}_{t}\sim\mathcal{N}(\bm{0},\mathbf{I})$

8:
$\mathbf{x}_{t-1}\leftarrow\sqrt{\bar{\alpha}_{t-1}}\hat{\mathbf{x}}_{0|t}+\sqrt{1-\bar{\alpha}_{t-1}}\bm{\epsilon}_{t}$
$\triangleright$ Noisor

9:end for

10:return $\mathbf{x}_{0}$

###  2.3 Predictor-Projector-Noisor (PPN) method

Fig. 2: MRI reconstructions for BraTS [16, 17] and fastMRI knee and brain [18]
at 8× acceleration, 50 NFEs.

Iterative generative methods for inverse problems depend heavily on accurate
sampling information. However, the manifold hypothesis suggests that images
exist on a lower-dimensional manifold embedded in a high-dimensional ambient
space, posing a challenge. For example, gradient-based generation can fail due
to the absence of scores in the ambient space [12]. This issue is particularly
prominent in projection-based methods such as DDNM [15] and MedScore [5],
where early iterations often yield ’disharmonized’ images outside the data
manifold, a region that lacks denoising information. Annealed Langevin
Dynamics (ALD) [12] and Diffusion Models (DMs) mitigate this by perturbing
data with Gaussian noise, creating noisy manifolds. The Repaint method [19],
building on DMs, improves controllable generation using precise denoising
information on noisy manifolds, reducing disharmony. However, its iterative
approach increases neural function evaluations (NFEs), slowing the process.

Inspired by Repaint [19], we propose a methood named Predictor-Projector-
Noiser (PPN) for faster sampling in undersampled MRI reconstruction.
Contrasting with Repaint’s method of executing multiple iterations between two
adjacent noise levels, our approach, based on deterministic DDIM [9],
increases the stride of each iteration while reducing the total number. As
depicted in Figure 1 and outlined in Algorithm 1, PPN loops once between ${\bf
x}_{t}$ and ${\bf x}_{0|t}$ per iteration. A single step of conditional
denoising is defined as follows:

$\mathbf{x}_{t-1}=\sqrt{\bar{\alpha}_{t-1}}\cdot{\bf P}_{{\bf y}}({\bf
x}_{0|t})+\sqrt{1-\bar{\alpha}_{t-1}}\cdot\bm{\epsilon}_{t}$ (6)

where ${\bf x}_{0|t}$ is obtained using Eq. 4, and ${\bf P}{\bf y}(\cdot)$ is
defined by Eq. 5. Starting with a noisy Zero-Filled sample ${\bf x}_{S}$, as
supported by Fig. 3, we can further accelerate PPN by iterating only through
the final $S$ steps, such as 50, out of the total 1,000 diffusion steps,
significantly fewer in number.

##  3 experiments

In this experiment section, our goal is to assess the PPN method’s performance
against three baselines: MedScore [5], DDNM [15], and DPS [10]. We utilize the
BraTS dataset [16, 17] for a standard retrospective study (Table 1) and
explore the effectiveness of NFEs on performance (Figure 3). Additionally, we
perform comparative evaluations with two other datasets, fastMRI Knee and
Brain [18] (Table 2).

Table 1: Results of undersampled MRI reconstructions on BraTS [16, 17] for 50 NFEs. Method | 12$\times$ Acceleration | 8$\times$ Acceleration | 4$\times$ Acceleration
---|---|---|---
PSNR$\uparrow$ | SSIM$\uparrow$ | PSNR$\uparrow$ | SSIM$\uparrow$ | PSNR$\uparrow$ | SSIM$\uparrow$
DDNM [15] | ${27.42}\scriptstyle\pm{2.55}$ | ${0.860}\scriptstyle\pm{0.037}$ | ${30.38}\scriptstyle\pm{2.40}$ | ${0.899}\scriptstyle\pm{0.027}$ | ${33.74}\scriptstyle\pm{2.44}$ | ${0.935}\scriptstyle\pm{0.019}$
MedScore [5] | ${27.91}\scriptstyle\pm{2.71}$ | ${0.870}\scriptstyle\pm{0.036}$ | ${31.47}\scriptstyle\pm{2.44}$ | ${0.914}\scriptstyle\pm{0.023}$ | ${34.14}\scriptstyle\pm{2.45}$ | ${0.941}\scriptstyle\pm{0.017}$
DPS [10] | $\mathbf{29.66\scriptstyle\pm 3.06}$ | ${0.892}\scriptstyle\pm{0.034}$ | ${33.23}\scriptstyle\pm{2.41}$ | ${0.928}\scriptstyle\pm{0.021}$ | ${36.93}\scriptstyle\pm{2.29}$ | ${0.957}\scriptstyle\pm{0.013}$
PPN (Ours) | ${29.07}\scriptstyle\pm{3.46}$ | $\mathbf{0.902\scriptstyle\pm 0.033}$ | $\mathbf{37.51\scriptstyle\pm 2.76}$ | $\mathbf{0.964\scriptstyle\pm 0.012}$ | $\mathbf{41.62\scriptstyle\pm 2.83}$ | $\mathbf{0.982\scriptstyle\pm 0.008}$
Table 2: SSIM scores for fastMRI knee and brain [18] reconstructions at 8× acceleration with 50 NFEs. Method | Knee | Brain
---|---|---
DDNM [15] | ${0.675}\scriptstyle\pm{0.087}$ | ${0.733}\scriptstyle\pm{0.060}$
MedScore [5] | ${0.719}\scriptstyle\pm{0.087}$ | ${0.833}\scriptstyle\pm{0.049}$
DPS [10] | ${0.709}\scriptstyle\pm{0.091}$ | ${0.831}\scriptstyle\pm{0.071}$
PPN (Ours) | $\mathbf{0.827\scriptstyle\pm 0.068}$ | $\mathbf{0.918\scriptstyle\pm 0.034}$

Dataset Our study utilized the Brain Tumor Segmentation (BraTS) 2021 dataset
[16, 17], converting 3D MRI volumes to 2D images (240x240 resolution), and
used roughly 257,000 slices for training. We also used the fastMRI brain and
knee datasets [18], selecting ’reconstruction_rss’ 2D images (320x320
resolution) and omitting the last 37.5% of brain slices and the first 5 and
last 3 knee slices, resulting in about 34,698 brain slices and 26,958 knee
slices for training.

Model Architecture We trained our models using an ablated diffusion model
(ADM) [20, 7] with 1,000 steps, enabling learn_sigma and adopting a cosine
noise schedule as ADM suggests. For the BraTS datasets, we employed the U-Net
architecture with 9.6M parameters, following the settings of MedScore [5]. For
the fastMRI knee and brain datasets, we used larger networks with the same
26.8M parameters.

Training All models were trained on NVIDIA GPUs: Tesla V100 was used for BraTS
with a learning rate of $2\times 10^{-4}$ over 75 epochs, while A100 was used
for the fastMRI brain dataset at $2\times 10^{-4}$ for 317 epochs and for the
fastMRI knee dataset at $5\times 10^{-5}$ over 1,400 epochs.

Baselines For fair comparison, we re-implemented MedScore [5], DDNM’s [15],
and DPS’s [10] Algorithm 1 in ADM, testing at the same checkpoint and
adjusting hyperparameters ($\lambda=1.0$ for MedScore, $\zeta=10.0$ for DPS).
The Baselines’ NFEs were configured using the DDIM setting, differing from
PPN’s approach as detailed in Section 2.3.

Unless specified otherwise, we report experiment results on 1,000 randomly
selected test images using the peak signal-to-noise ratio (PSNR) and
structural similarity index (SSIM) metrics, based on 50 NFEs and uniform 1D
Cartesian masks.

Fig. 3: Performance vs. NFEs in 4× acceleration reconstruction on BraTS [16,
17]. Shaded areas is standard deviations.

##  4 Ablation study

In this ablation study, we analyzed PPN’s three main components. First, we
replaced the noisor (Algorithm 1, Line 7) with DDPM’s Eq. 6 [6], resulting in
DDNM’s Algorithm 1 [15]. Second, we applied projection at ${\bf x}_{t}$
instead of ${\bf x}_{0}$, equivalent to MedScore’s Algorithm 2 [5], by setting
$\lambda=1$ in Line 5. Third, we tested the effectiveness of NFEs, as
illustrated in Figure 3. Section 5 shows that PPN achieves optimal results
using the final $S=50$ steps across the majority of experiments on the three
datasets in this study.

##  5 Results

Figure 2 illustrates the reconstructions of various anatomies at 8x
acceleration across all methods, showcasing PPN’s fidelity to the ground
truth. In Table 1, PPN is compared with baseline methods for MRI undersampling
at 4x, 8x, and 12x accelerations using 50 NFEs. Here, PPN outperforms the
baseline methods at all accelerations except for the PSNR metric of the 12x
acceleration, showcasing PPN’s generalizability. This superiority is further
supported by Table 2, where PPN achieves the highest SSIM for fastMRI knee and
brain datasets at 50 NFEs and 8x acceleration. Additionally, Figure 3 shows
PPN’s performance against NFEs in PSNR and SSIM at 4x acceleration. Notably,
PPN achieves better results with 50 NFEs compared to others at 500 NFEs.

##  6 Conclusion

We have introduced a novel approach called PPN for accelerating diffusion
models in the context of undersampled MRI reconstruction. Our proposed method
enables faster sampling while maintaining high fidelity, resulting in
significant time reduction compared to other methods. Furthermore, our
unsupervised model exhibits excellent generalizability and adaptability,
enabling it to reconstruct various undersampling patterns without requiring
retraining. Future work could focus on investigating PPN’s performance in
distribution shift tasks and its applicability to in-vivo scenarios.

## Acknowledgment

WJ acknowledges support from Dr F & Mrs ME Zaccari Scholarship. HS
acknowledges support from the Australian Research Council (DE210101297,
DP230101628).

## Compliance with Ethical Standards

This research study was conducted retrospectively using human subject data
made available in open access by Brain Tumor Segmentation (BraTS) [16, 17] and
fastMRI [18] datasets. Ethical approval was not required as confirmed by the
license attached with the open access data.

## References

  * [1] Kerstin Hammernik, Teresa Klatzer, Erich Kobler, Michael P Recht, Daniel K Sodickson, Thomas Pock, and Florian Knoll,  “Learning a variational network for reconstruction of accelerated mri data,”  Magnetic resonance in medicine, vol. 79, no. 6, pp. 3055–3071, 2018.
  * [2] Yoseo Han, Leonard Sunwoo, and Jong Chul Ye,  “k-space deep learning for accelerated mri,”  IEEE transactions on medical imaging, vol. 39, no. 2, pp. 377–386, 2019.
  * [3] Dongwook Lee, Jaejun Yoo, Sungho Tak, and Jong Chul Ye,  “Deep residual learning for accelerated mri using magnitude and phase networks,”  IEEE Transactions on Biomedical Engineering, vol. 65, no. 9, pp. 1985–1995, 2018.
  * [4] Yang Gao, Martijn Cloos, Feng Liu, Stuart Crozier, G Bruce Pike, and Hongfu Sun,  “Accelerating quantitative susceptibility and r2* mapping using incoherent undersampling and deep neural network reconstruction,”  NeuroImage, vol. 240, pp. 118404, 2021.
  * [5] Yang Song, Liyue Shen, Lei Xing, and Stefano Ermon,  “Solving inverse problems in medical imaging with score-based generative models,”  arXiv preprint arXiv:2111.08005, 2021.
  * [6] Jonathan Ho, Ajay Jain, and Pieter Abbeel,  “Denoising diffusion probabilistic models,”  Advances in neural information processing systems, vol. 33, pp. 6840–6851, 2020.
  * [7] Alexander Quinn Nichol and Prafulla Dhariwal,  “Improved denoising diffusion probabilistic models,”  in International Conference on Machine Learning. PMLR, 2021, pp. 8162–8171.
  * [8] Yang Song, Jascha Sohl-Dickstein, Diederik P Kingma, Abhishek Kumar, Stefano Ermon, and Ben Poole,  “Score-based generative modeling through stochastic differential equations,”  arXiv preprint arXiv:2011.13456, 2020.
  * [9] Jiaming Song, Chenlin Meng, and Stefano Ermon,  “Denoising diffusion implicit models,”  arXiv preprint arXiv:2010.02502, 2020.
  * [10] Hyungjin Chung, Jeongsol Kim, Michael T Mccann, Marc L Klasky, and Jong Chul Ye,  “Diffusion posterior sampling for general noisy inverse problems,”  arXiv preprint arXiv:2209.14687, 2022.
  * [11] Hyungjin Chung and Jong Chul Ye,  “Score-based diffusion models for accelerated mri,”  Medical image analysis, vol. 80, pp. 102479, 2022.
  * [12] Yang Song and Stefano Ermon,  “Generative modeling by estimating gradients of the data distribution,”  Advances in neural information processing systems, vol. 32, 2019.
  * [13] Tom Tirer and Raja Giryes,  “Image restoration by iterative denoising and backward projections,”  IEEE Transactions on Image Processing, vol. 28, no. 3, pp. 1220–1234, 2018.
  * [14] Singanallur V Venkatakrishnan, Charles A Bouman, and Brendt Wohlberg,  “Plug-and-play priors for model based reconstruction,”  in 2013 IEEE global conference on signal and information processing. IEEE, 2013, pp. 945–948.
  * [15] Yinhuai Wang, Jiwen Yu, and Jian Zhang,  “Zero-shot image restoration using denoising diffusion null-space model,”  arXiv preprint arXiv:2212.00490, 2022.
  * [16] Bjoern H Menze, Andras Jakab, Stefan Bauer, Jayashree Kalpathy-Cramer, Keyvan Farahani, Justin Kirby, Yuliya Burren, Nicole Porz, Johannes Slotboom, Roland Wiest, et al.,  “The multimodal brain tumor image segmentation benchmark (brats),”  IEEE transactions on medical imaging, vol. 34, no. 10, pp. 1993–2024, 2014.
  * [17] Spyridon Bakas, Hamed Akbari, Aristeidis Sotiras, Michel Bilello, Martin Rozycki, Justin S Kirby, John B Freymann, Keyvan Farahani, and Christos Davatzikos,  “Advancing the cancer genome atlas glioma mri collections with expert segmentation labels and radiomic features,”  Scientific data, vol. 4, no. 1, pp. 1–13, 2017.
  * [18] Jure Zbontar, Florian Knoll, Anuroop Sriram, Tullie Murrell, Zhengnan Huang, Matthew J. Muckley, Aaron Defazio, Ruben Stern, Patricia Johnson, Mary Bruno, Marc Parente, Krzysztof J. Geras, Joe Katsnelson, Hersh Chandarana, Zizhao Zhang, Michal Drozdzal, Adriana Romero, Michael Rabbat, Pascal Vincent, Nafissa Yakubova, James Pinkerton, Duo Wang, Erich Owens, C. Lawrence Zitnick, Michael P. Recht, Daniel K. Sodickson, and Yvonne W. Lui,  “fastMRI: An open dataset and benchmarks for accelerated MRI,” 2018.
  * [19] Andreas Lugmayr, Martin Danelljan, Andres Romero, Fisher Yu, Radu Timofte, and Luc Van Gool,  “Repaint: Inpainting using denoising diffusion probabilistic models,”  in Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, 2022, pp. 11461–11471.
  * [20] Prafulla Dhariwal and Alexander Nichol,  “Diffusion models beat gans on image synthesis,”  Advances in neural information processing systems, vol. 34, pp. 8780–8794, 2021.


Consider normal coordinates $(u^{i})_{i\in\\{1,\ldots,d\\}}$ centered at
$\theta^{\star}$ with respect to the orthonormal basis
$(\mathbf{e}_{i})_{i\in\\{1,\dots,d\\}}$ of
$\mathrm{T}_{\theta^{\star}}\Theta$. Define for any $i,j\in\\{1,\ldots,d\\}$,
$h^{i}:\Theta\to\mathbb{R}$, $\Sigma_{ij}:\Theta\to\mathbb{R}$ by
$h^{i}=\mathrm{d}u^{i}(h)$ and
$\Sigma_{ij}=[\mathrm{d}u^{i}\otimes\mathrm{d}u^{j}]\\{\Sigma\\}$. Let
$\mathrm{g}:\mathrm{T}_{\theta^{\star}}\Theta\to\mathbb{R}$ be a smooth
function with compact support. Applying Lemma 21 to $\mathrm{g}$ gives (103).
Using MD 3, $\Sigma$ is continuous, which implies that for any
$\overline{u}_{0}\in\mathrm{T}_{\theta^{\star}}\Theta$,

$\Sigma\left(\mathrm{Exp}_{\theta^{\star}}(\upeta^{1/2}\overline{u}_{0})\right)=\sum_{i,j=1}^{d}\left\\{\Sigma_{\star}^{ij}+\mathscr{R}^{ij}_{\Sigma}\left(\upeta^{1/2}\overline{u}_{0}\right)\right\\}\partial
u_{i}\otimes\partial u_{j}\;,$ (135)

where for any $i,j\in\\{1,\dots,d\\}$,
$\Sigma_{\star}^{ij}=\Sigma_{ij}({\theta^{\star}})$,
$\mathscr{R}^{ij}_{\Sigma}$ is continuous over
$\mathrm{T}_{\theta^{\star}}\Theta$ and $\mathscr{R}^{ij}_{\Sigma}(0)=0$.
Using Lemma 22, replacing $\Sigma_{ij}$ and $h^{i}$ in (103) with (116) and
(135) gives for any $\overline{u}_{0}\in\mathrm{T}_{\theta^{\star}}\Theta$,

$\displaystyle\overline{R}_{\upeta}\mathrm{g}(\overline{u}_{0})=\mathrm{g}(\overline{u}_{0})+\upeta\sum_{i=1}^{d}\partial_{i}\mathrm{g}(\overline{u}_{0})\sum_{k=1}^{d}\mathbf{A}_{k}^{i}\overline{u}_{0}^{k}+(\upeta/2)\sum_{i,j=1}^{d}\partial^{2}_{ij}\mathrm{g}(\overline{u}_{0})$
$\displaystyle\Sigma_{\star}^{ij}+\upeta\mathscr{R}_{\mathrm{g},\upeta,\Sigma,h}(\overline{u}_{0})$
(136)
$\displaystyle\,+(\upeta/6)\overline{\mathscr{R}}_{\mathrm{g},\upeta}(\overline{u}_{0})\;,$

where $\overline{u}^{k}_{0}$ are the components of $\overline{u}_{0}$ in
$(\mathbf{e}_{i})_{i\in\\{1,\ldots,d\\}}$,

$\mathscr{R}_{\mathrm{g},\upeta,\Sigma,h}(\overline{u}_{0})=\upeta^{-1/2}\sum_{i=1}^{d}\mathscr{R}_{h}^{i}\left(\upeta^{1/2}\overline{u}_{0}\right)\partial_{i}\mathrm{g}(\overline{u}_{0})\\\
+(1/2)\sum_{i,j=1}^{d}\left\\{\partial^{2}_{ij}\mathrm{g}(\overline{u}_{0})-\upeta^{1/2}\sum_{k=1}^{d}\Gamma_{ij}^{k}\left(\mathrm{Exp}_{\theta^{\star}}(\upeta^{1/2}\overline{u}_{0})\right)\partial_{k}\mathrm{g}(\overline{u}_{0})\right\\}\left[\mathscr{R}_{\Sigma}^{ij}\left(\upeta^{1/2}\overline{u}_{0}\right)\right]\\\
+(1/2)\sum_{i,j=1}^{d}\left\\{\partial^{2}_{ij}\mathrm{g}(\overline{u}_{0})-\upeta^{1/2}\sum_{k=1}^{d}\Gamma_{ij}^{k}\left(\mathrm{Exp}_{\theta^{\star}}(\upeta^{1/2}\overline{u}_{0})\right)\partial_{k}\mathrm{g}(\overline{u}_{0})\right\\}\left[h^{i}h^{j}\left(\mathrm{Exp}_{\theta^{\star}}(\upeta^{1/2}\overline{u}_{0})\right)\right]\\\
-(\upeta^{1/2}/2)\sum_{i,j,k=1}^{d}\Gamma_{ij}^{k}\left(\mathrm{Exp}_{\theta^{\star}}(\upeta^{1/2}\overline{u}_{0})\right)\partial_{k}\mathrm{g}(\overline{u}_{0})\Sigma^{ij}_{\star}\;.$
(137)

By Lemma 23,
$(\overline{\nu}^{\upeta})_{\upeta\in\left(0,\overline{\upeta}\right]}$ is
tight and therefore relatively compact. Therefore, it is enough that for any
limit point $\overline{\nu}^{\star}$,
$\overline{\nu}^{\star}=\mathrm{N}(0,\mathbf{V})$ where
$\mathbf{V}\in\mathbb{R}^{d\times d}$ is the solution of the Lyapunov equation
$\mathbf{A}\mathbf{V}+\mathbf{V}\mathbf{A}^{\top}=\Sigma({\theta^{\star}})$.
Let $(\upeta_{n})_{n\in\mathbb{N}^{*}}$ be a sequence with values in
$(0,\overline{\upeta}]$, such that $\lim_{n\to+\infty}\upeta_{n}=0$, and
$(\overline{\nu}^{\upeta_{n}})_{n\in\mathbb{N}^{*}}$ weakly converges to
$\overline{\nu}^{\star}$.

First by (136), we have

$\displaystyle\int_{\mathrm{T}_{{\theta^{\star}}}\Theta}\overline{\nu}^{\upeta_{n}}(\mathrm{d}\overline{u}_{0})\int_{\mathrm{T}_{{\theta^{\star}}}\Theta}\overline{R}_{\upeta_{n}}(\overline{u}_{0},\mathrm{d}\overline{u}_{1})\mathrm{g}(\overline{u}_{1})$
(138)
$\displaystyle\quad=\int_{\mathrm{T}_{{\theta^{\star}}}\Theta}\overline{\nu}^{\upeta_{n}}(\mathrm{d}\overline{u}_{0})\mathrm{g}(\overline{u}_{0})+\upeta_{n}\int_{\mathrm{T}_{{\theta^{\star}}}\Theta}\overline{\nu}^{\upeta_{n}}(\mathrm{d}\overline{u}_{0})\sum_{i=1}^{d}\partial_{i}\mathrm{g}(\overline{u}_{0})\sum_{k=1}^{d}\mathbf{A}_{k}^{i}\overline{u}_{0}^{k}$
(139)
$\displaystyle\qquad+(\upeta_{n}/2)\int_{\mathrm{T}_{{\theta^{\star}}}\Theta}\overline{\nu}^{\upeta_{n}}(\mathrm{d}\overline{u}_{0})\sum_{i,j=1}^{d}\partial^{2}_{ij}\mathrm{g}(\overline{u}_{0})\Sigma_{\star}^{ij}+\upeta_{n}\int_{\mathrm{T}_{{\theta^{\star}}}\Theta}\overline{\nu}^{\upeta_{n}}(\mathrm{d}\overline{u}_{0})\mathscr{R}_{\mathrm{g},\upeta_{n},\Sigma,h}(\overline{u}_{0})$
(140)
$\displaystyle\qquad+(\upeta_{n}/6)\int_{\mathrm{T}_{{\theta^{\star}}}\Theta}\overline{\nu}^{\upeta_{n}}(\mathrm{d}\overline{u}_{0})\overline{\mathscr{R}}_{\mathrm{g},\upeta_{n}}(\overline{u}_{0})\;.$
(141)

Therefore using that $\overline{\nu}^{\upeta_{n}}$ is stationary with respect
to $\overline{R}_{\upeta_{n}}$, we obtain that

$\limsup_{n\to+\infty}\left|\int_{\mathrm{T}_{{\theta^{\star}}}\Theta}\overline{\nu}^{\upeta_{n}}(\mathrm{d}\overline{u}_{0})\left\\{\sum_{i=1}^{d}\partial_{i}\mathrm{g}(\overline{u}_{0})\sum_{k=1}^{d}\mathbf{A}_{k}^{i}\overline{u}_{0}^{k}+\sum_{i,j=1}^{d}\partial^{2}_{ij}\mathrm{g}(\overline{u}_{0})\Sigma_{\star}^{ij}\right\\}\right|\\\
\leq\limsup_{n\to+\infty}\left|\int_{\mathrm{T}_{{\theta^{\star}}}\Theta}\overline{\nu}^{\upeta_{n}}(\mathrm{d}\overline{u}_{0})\mathscr{R}_{\mathrm{g},\upeta_{n},\Sigma,h}(\overline{u}_{0})\right|+\left|\int_{\mathrm{T}_{{\theta^{\star}}}\Theta}\overline{\nu}^{\upeta_{n}}(\mathrm{d}\overline{u}_{0})\overline{\mathscr{R}}_{\mathrm{g},\upeta_{n}}(\overline{u}_{0})\right|\;.$
(142)

Consider a sequence of independent random variables
$(\mathrm{Y}_{n})_{n\in\mathbb{N}}$ such that for any $n\in\mathbb{N}$, the
law of $\mathrm{Y}_{n}$ is $\overline{\nu}^{\upeta_{n}}$. By Slutsky’s
theorem, since $(\mathrm{Y}_{n})_{n\in\mathbb{N}}$ converges in distribution
and $\lim_{n\to+\infty}\upeta_{n}=0$, we obtain that
$\upeta_{n}^{1/2}\mathrm{Y}_{n}$ converges in distribution towards $0$.
Moreover, using the continuous mapping theorem, we have

$\limsup_{n\to+\infty}\left|\mathbb{E}[\mathscr{R}_{\mathrm{g},\upeta_{n},\Sigma,h}(\mathrm{Y}_{n})]\right|=0\;.$
(143)

Similarly, we use (106) to obtain, for any $n\in\mathbb{N}$ and $K>0$,

$\displaystyle\left|\overline{\mathscr{R}}_{\mathrm{g},\upeta_{n}}(\mathrm{Y}_{n})\right|$
$\displaystyle\leq
8\upeta_{n}^{1/2}\mathbbm{1}_{\mathsf{K}_{K}}(\uptheta_{n})\mathbb{E}\left[\|\mathsf{C}^{3}(\mathrm{g},\upeta_{n})\|_{\upgamma}\mathbbm{1}_{\mathsf{A}_{\uptheta_{n}}^{\complement}}\left\|H_{K}\right\|_{\uptheta_{n}}^{3}\middle|\uptheta_{n}\right]$
(144)
$\displaystyle\quad+16\|\mathsf{C}^{2}(\mathrm{g},\upeta_{n})\|\mathbb{E}\left[\left\|Y_{K}\right\|_{\uptheta_{n}}^{2}\middle|\uptheta_{n}\right]\;,$
(145)

where for any $n\in\mathbb{N}$,
$\uptheta_{n}=\mathrm{Exp}_{\theta^{\star}}(\upeta^{1/2}_{n}\mathrm{Y}_{n})$
are independent random variables and by (94), the distribution of
$\uptheta_{n}$ is $\mu^{\upeta_{n}}$. Thus we obtain for any $K\geq 0$, using
$\mathbbm{1}_{\mathsf{K}_{K}}(\uptheta_{n})\left\|H_{K}\right\|_{\uptheta_{n}}$
is almost surely bounded by
$4[K^{3}+\sup_{\theta\in\mathsf{K}_{K}}\left\|h(\theta)\right\|_{\theta}^{3}]$,
Markov’s inequality and MD 4,

$\displaystyle\limsup_{n\to+\infty}\left|\mathbb{E}[\overline{\mathscr{R}}_{\mathrm{g},{\upeta_{n}}}(\mathrm{Y}_{n})]\right|$
$\displaystyle\leq\limsup_{n\to+\infty}16\|\mathsf{C}^{2}(\mathrm{g},\upeta_{n})\|\mathbb{E}[\left\|e_{\uptheta_{n}}(X_{1})\right\|_{\uptheta_{n}}^{2}\\{1-\chi_{K}(\|e_{\uptheta_{n}}(X_{1})\|_{\uptheta_{n}})\\}]$
(146) $\displaystyle\leq
16\|\mathsf{C}^{2}(\mathrm{g},0)\|K^{-\varepsilon}\\{\tilde{\sigma}_{0}^{2}+\tilde{\sigma}_{1}^{2}\mathbb{E}[V({\theta^{\star}})]\\}\;,$
(147)

using that $(\uptheta_{n})_{n\in\mathbb{N}}$ converges in distribution to
${\theta^{\star}}$. For any smooth function with compact support
$\mathrm{g}:\mathrm{T}_{{\theta^{\star}}}\Theta\to\mathbb{R}$, combining
(142)-(143)-(147), taking $K\to+\infty$ and using the weak convergence of
$(\overline{\nu}^{\upeta_{n}})_{n\in\mathbb{N}}$ to $\overline{\nu}^{\star}$
when $n\to+\infty$ shows that

$\int_{\mathrm{T}_{{\theta^{\star}}}\Theta}\overline{\nu}^{\star}(\mathrm{d}\overline{u}_{0})\left\\{\sum_{i=1}^{d}\partial_{i}\mathrm{g}(\overline{u}_{0})\sum_{k=1}^{d}\mathbf{A}_{k}^{i}\overline{u}_{0}^{k}+\sum_{i,j=1}^{d}\partial^{2}_{ij}\mathrm{g}(\overline{u}_{0})\Sigma_{\star}^{ij}\right\\}=0\;.$
(148)

Finally, by (Horn and Johnson,, 1994, Theorem 2.2.1), there exists a unique
matrix $\mathbf{V}\in\mathbb{R}^{d\times d}$ solution to the Lyapunov equation
$\mathbf{A}\mathbf{V}+\mathbf{V}\mathbf{A}^{\top}=\Sigma({\theta^{\star}})$.
By (Kent,, 1978, Theorem 10.1), $\mathrm{N}(0,\mathbf{V})$ is the unique
probability distribution on $\mathrm{T}_{{\theta^{\star}}}\Theta$ satisfying
(148). This concludes the proof. ∎

##  Appendix D Proofs for Section 4

###  D.1 Proof of Lemma 8

Recall that $f$ is $\lambda_{f}$-strongly geodesically convex, if and only if
for any $\theta_{1},\theta_{2}\in\Theta$,

$f(\theta_{2})\geq
f(\theta_{1})+\left\langle\mathrm{Exp}_{\theta_{1}}^{-1}(\theta_{2}),\mathrm{grad}\,f(\theta_{1})\right\rangle_{\theta_{1}}+\lambda_{f}\rho_{\Theta}^{2}(\theta_{1},\theta_{2})\;.$
(149)

Put $\theta_{1}={\theta^{\star}}$ and $\theta_{2}=\theta$. Since
${\theta^{\star}}$ is a stationary point of $f$, so
$\mathrm{grad}\,f({\theta^{\star}})=0$, it follows from (149) that

$f(\theta)-f({\theta^{\star}})\geq\lambda_{f}\rho_{\Theta}^{2}({\theta^{\star}},\theta)\;,$
(150)

which is the second identity in (15). To obtain the first identity, put
$\theta_{1}=\theta$ and $\theta_{2}={\theta^{\star}}$, in (149), so

$f({\theta^{\star}})-f(\theta)\geq\left\langle\mathrm{Exp}_{\theta}^{-1}({\theta^{\star}}),\mathrm{grad}\,f(\theta)\right\rangle_{\theta}+\lambda_{f}\rho_{\Theta}^{2}({\theta^{\star}},\theta)\;.$
(151)

Since $f({\theta^{\star}})\leq f(\theta)$, this implies

$-\left\langle\mathrm{Exp}_{\theta}^{-1}({\theta^{\star}}),\mathrm{grad}\,f(\theta)\right\rangle_{\theta}\geq\lambda_{f}\rho_{\Theta}^{2}({\theta^{\star}},\theta)=\lambda_{f}\left\|\mathrm{Exp}_{\theta}^{-1}({\theta^{\star}})\right\|_{\theta}^{2}\;.$
(152)

Or, after using the Cauchy-Schwarz inequality,

$\left\|\mathrm{grad}\,f(\theta)\right\|_{\theta}\geq\lambda_{f}\left\|\mathrm{Exp}_{\theta}^{-1}({\theta^{\star}})\right\|_{\theta}\;.$
(153)

Finally, using once more the Cauchy-Schwarz inequality, and (151) and (153),

$f(\theta)-f({\theta^{\star}})\leq-\left\langle\mathrm{Exp}_{\theta}^{-1}({\theta^{\star}}),\mathrm{grad}\,f(\theta)\right\rangle_{\theta}\leq(1/\lambda_{f})\left\|\mathrm{grad}\,f(\theta)\right\|_{\theta}^{2}\;,$
(154)

which is equivalent to the first identity in (15).

###  D.2 Proof of Lemma 10

Without loss of generality, we assume that $f({\theta^{\star}})=0$. First, we
show that for any $\theta\in\Theta$,

$f(\theta)\leq M_{f}\rho_{\Theta}^{2}({\theta^{\star}},\theta)\;.$ (155)

Let $\theta\in\Theta$ and $\upgamma:\left[0,1\right]\to\Theta$ the unique
geodesic such that $\upgamma(0)={\theta^{\star}}$ and $\upgamma(1)=\theta$.
Then since $f$ is continuously differentiable using (Lee,, 2019, Proposition
4.15-(ii) and Theorem 4.24-(iii)), we get that
$f(\theta)=\int_{0}^{1}\left\langle\mathrm{grad}\,f(\upgamma(t)),\dot{\upgamma}(t)\right\rangle_{\upgamma(t)}\mathrm{d}t$.
Therefore, using the Cauchy-Schwarz inequality and for any
$t\in\left[0,1\right]$,
$\left\|\dot{\upgamma}(t)\right\|_{\upgamma(t)}=\rho_{\Theta}({\theta^{\star}},\theta)$
we obtain that
$\left|f(\theta)\right|\leq\rho_{\Theta}({\theta^{\star}},\theta)\left\|\mathrm{grad}\,f(\upgamma(t))\right\|_{\upgamma(t)}$
which shows that (155) holds by assumption.

We now proceed with the proof of the main statement. Since $f$ is twice
continuously differentiable, $\tilde{f}$ has this same property. In addition,
for any $\theta\in\Theta$,

$\mathrm{grad}\,\tilde{f}(\theta)=\mathrm{grad}\,f(\theta)/[2(f(\theta)+1)^{1/2}]\;.$
(156)

Therefore, using the assumption that for any $\theta\in\Theta$,
$\left\|\mathrm{grad}\,f(\theta)\right\|_{\theta}^{2}\leq
M_{f}\rho_{\Theta}^{2}({\theta^{\star}},\theta)$ and the second inequality of
Lemma 8, we get that

$\|\mathrm{grad}\,\tilde{f}(\theta)\|_{\theta}=\left\|\mathrm{grad}\,f(\theta)\right\|_{\theta}/[2(f(\theta)+1)^{1/2}]\leq
M_{f}^{1/2}\rho_{\Theta}({\theta^{\star}},\theta)/[2(\lambda_{f}\rho_{\Theta}^{2}({\theta^{\star}},\theta)+1)^{1/2}]\\\
\leq C_{f}^{1/2}[1\wedge\rho_{\Theta}({\theta^{\star}},\theta)]\;,$ (157)

with $C_{f}^{1/2}\leftarrow(M_{f}^{1/2}/2)[1\wedge\lambda_{f}^{-1/2}]$.

It remains to show that for any $\theta\in\Theta$,
$-\langle\mathrm{Exp}^{-1}_{\theta}({\theta^{\star}}),\mathrm{grad}\,\tilde{f}(\theta)\rangle_{\theta}\geq\tilde{\lambda}_{f}V_{1}(\theta)$,
where $V_{1}$ is defined by (10) with $\delta=1$ and
$\tilde{\lambda}_{f}\leftarrow\lambda_{f}^{1/2}/2$. Using (156) again, F 2 and
(155), we obtain that for any $\theta\in\Theta$,

$-\left\langle\mathrm{Exp}^{-1}_{\theta}({\theta^{\star}}),\mathrm{grad}\,\tilde{f}(\theta)\right\rangle_{\theta}=-\left\langle\mathrm{Exp}^{-1}_{\theta}({\theta^{\star}}),\mathrm{grad}\,f(\theta)\right\rangle_{\theta}/[2(f(\theta)+1)^{1/2}]\\\
\geq\lambda_{f}\rho_{\Theta}^{2}({\theta^{\star}},\theta)/[2(f(\theta)+1)^{1/2}]\geq\lambda_{f}\rho_{\Theta}^{2}({\theta^{\star}},\theta)/[2(M_{f}\rho_{\Theta}^{2}({\theta^{\star}},\theta)+1)^{1/2}]\;.$
(158)

Using that for any $\theta\in\Theta$,
$V_{1}(\theta)=\\{\rho_{\Theta}^{2}({\theta^{\star}},\theta)+1\\}^{1/2}-1\leq\rho_{\Theta}({\theta^{\star}},\theta)$,
we get that

$-\left\langle\mathrm{Exp}^{-1}_{\theta}({\theta^{\star}}),\mathrm{grad}\,\tilde{f}(\theta)\right\rangle_{\theta}\geq\lambda_{f}V_{1}(\theta)\rho_{\Theta}({\theta^{\star}},\theta)/[2(M_{f}\rho_{\Theta}^{2}({\theta^{\star}},\theta)+1)^{1/2}]\geq\lambda_{f}V_{1}(\theta)/(2M_{f}^{1/2})\;.$
(159)

###  D.3 Proof of Proposition 11

The proof consists in an application of Theorem 1-(b). First, by Proposition
5, $V_{1}$ defined by (10) with $\delta=1$, satisfies H 1. In addition, by
(Durmus et al.,, 2020, Lemma 16), $V_{1}$ is continuously differentiable with
gradient given for any $\theta\in\Theta$ by

$\mathrm{grad}\,V_{1}(\theta)=-\mathrm{Exp}^{-1}_{\theta}({\theta^{\star}})/\\{1+\rho_{\Theta}^{2}({\theta^{\star}},\theta)\\}^{1/2}\;.$
(160)

Therefore, for any $\theta\in\Theta$, by F 3 we get

$\left\langle\mathrm{grad}\,V_{1}(\theta),\mathrm{grad}\,f(\theta)\right\rangle_{\theta}=-\left\langle\mathrm{Exp}^{-1}_{\theta}({\theta^{\star}}),\mathrm{grad}\,f(\theta)\right\rangle_{\theta}/\\{1+\rho_{\Theta}^{2}({\theta^{\star}},\theta)\\}^{1/2}\\\
\geq\tilde{\lambda}_{f}V_{1}(\theta)/\\{1+\rho_{\Theta}^{2}({\theta^{\star}},\theta)\\}^{1/2}\;.$
(161)

In addition, $t^{2}\wedge 1-ab\\{(t^{2}+1)^{1/2}-1\\}/(1+t^{2})^{1/2}\leq 0$
for any $t\geq 0$, $b>0$ and $a=4b^{-1}$ using that $(t^{2}+1)^{1/2}-1\geq
t^{2}/[2(1+t^{2})^{1/2}]$. As a result, using F 3 for any $t\geq 0$, $b>0$ and
$a=4b^{-1}$, it follows that H 2 is satisfied with $C_{1}\leftarrow
0,C_{2}\leftarrow 4C_{f}/\tilde{\lambda}_{f}$ for $h=-\mathrm{grad}\,f$ and
$V\leftarrow V_{1}$. Therefore, we obtain using Theorem 1-(b) that for any
$\upeta\in\left(0,\upeta\right]$,

$n^{-1}\sum_{k=0}^{n-1}\mathbb{E}\left[\left\langle\mathrm{grad}\,V_{1}(\theta_{k}),\mathrm{grad}\,f(\theta_{k})\right\rangle_{\theta_{k}}\right]\leq
2V_{1}(\theta_{0})/(n\upeta)+2\upeta(1+\kappa)\sigma_{0}^{2}\;,$ (162)

where
$\overline{\upeta}=[(8C_{f}/\tilde{\lambda}_{f})(1+\kappa)(1+\sigma_{1}^{2})]^{-1}$.
Using (161), we have

$(\tilde{\lambda}_{f}/n)\sum_{k=0}^{n-1}\mathbb{E}\left[V_{1}(\theta_{k})/\\{1+\rho_{\Theta}^{2}({\theta^{\star}},\theta_{k})\\}^{1/2}\right]\leq
2V_{1}(\theta_{0})/(n\upeta)+2\upeta(1+\kappa)\sigma_{0}^{2}\;,$ (163)

which concludes the proof since $(t^{2}+1)^{1/2}-1\geq
t^{2}/[2(1+t^{2})^{1/2}]$ for any $t\geq 0$ implying
$V_{1}(\theta)/\\{1+\rho_{\Theta}^{2}({\theta^{\star}},\theta)\\}^{1/2}\geq
D^{2}_{\Theta}({\theta^{\star}},\theta)/2$ for any $\theta\in\Theta$.

###  D.4 Proof of Proposition 12

Define $\mathsf{X}=\\{\overline{\theta}_{i}\,:\,i\in\\{1,\ldots,M_{\pi}\\}\\}$
and recall that
$\mathrm{D}=\sup\\{\rho_{\Theta}(\theta_{0},\overline{\theta})\,:\,\overline{\theta}\in\mathsf{X}\\}$.
Set $\mathsf{S}=\overline{\operatorname{B}}(\theta_{0},\rm D)$. Note that the
closed ball $\mathsf{S}$, is compact by (Jost,, 2005, Theorem 1.7.1),
geodesically convex, and $\mathsf{X}\subset\mathsf{S}$, as well as
$\theta_{0}\in\mathsf{S}$. We consider in this section, for any
$\theta\in\Theta$ and $x\in\mathsf{X}$,
$H_{\theta}(x)=\mathrm{Exp}^{-1}_{\theta}(x)$.

First note that $\theta_{n}\in\mathsf{S}$, for all $n\in\mathbb{N}$ by a
straightforward induction using that $\mathsf{S}$ is geodesically convex and
$\theta_{0}\in\mathsf{S}$. Indeed, $\theta_{0}\in\mathsf{S}$, and, if
$\theta_{n}\in\mathsf{S}$, then $\theta_{n+1}$ lies on the geodesic segment
connecting $\theta_{n}$ and $X_{n+1}$, two points which belong to
$\mathsf{S}$, and therefore $\theta_{n+1}\in\mathsf{S}$. This means that the
SGD scheme used here is equivalent to

$\theta_{n+1}=\operatorname{proj}_{\mathsf{S}}\left(\mathrm{Exp}_{\theta_{n}}\\!\left(\upeta
H_{\theta_{n}}(X_{n+1})\right)\right)\;.$ (164)

Define $\mathsf{H}$ and $V_{2}$ as in Proposition 4. It is possible to show
that $\mathsf{H}=\mathsf{S}$. Indeed, for $\theta\in\mathsf{S}$, and
$x\in\mathsf{X}$, since $x\in\mathsf{S}$, and $\mathsf{S}$ is convex, the
geodesic segment connecting $\theta$ to $x$ is entirely contained in
$\mathsf{S}$. However, by definition, this geodesic segment is the set of
points $\mathrm{Exp}_{\theta}(tH_{\theta}(x))$, where $t\in[0,1]$. Now, since
$\upeta\leq\overline{\upeta}\leq 1$, Proposition 4 implies that $V_{2}$
verifies H 1-(i)-(ii) where $L\leftarrow CL_{\pi}$,
$L_{\pi}=(\mathrm{D}+1)(1+\kappa\coth(\kappa\mathrm{D}))$ and $C$ is a
universal constant.

The objective function $f$ satisfies F 2 with $\lambda_{f}=1/2$ (that is, $f$
is $1/2$-strongly convex), since by (Jost,, 2005, Theorem 5.6.1)
$f_{i}(\theta)=\rho_{\Theta}^{2}(\theta,\overline{\theta}_{i})/2$ is
$1$-strongly geodesically convex for any $i\in\\{1,\ldots,M_{\pi}\\}$. Thus,
by (149) for all $\theta\in\mathsf{S}$

$\left\langle\mathrm{Exp}^{-1}_{\theta}({\theta^{\star}}),\mathrm{grad}\,f(\theta)\right\rangle_{\theta}\leq-(1/2)\rho_{\Theta}^{2}({\theta^{\star}},\theta)\;.$
(165)

Now, for any $\theta\in\mathsf{S}$, $v\in\mathrm{T}_{\theta}\Theta$, using
(Jost,, 2005, Theorem 5.6.1), we have,

$\displaystyle\left\|\mathrm{Hess}\,f_{\theta}(v,v)\right\|_{\theta}$
$\displaystyle\leq
M_{\pi}^{-1}\sum_{i=1}^{M_{\pi}}\left\|(\mathrm{Hess}\,f_{i})_{\theta}(v,v)\right\|_{\theta}$
(166) $\displaystyle\leq
M_{\pi}^{-1}\sum_{i=1}^{M_{\pi}}\kappa\rho_{\Theta}(\theta,\overline{\theta}_{i})\coth(\kappa\rho_{\Theta}(\theta,\overline{\theta}_{i}))\left\|v\right\|_{\theta}^{2}\leq\tilde{L}_{\pi}\left\|v\right\|_{\theta}^{2}\;,$
(167)

where $\tilde{L}_{\pi}=2\mathrm{D}\kappa\coth(2\kappa\mathrm{D})$, since
$t\mapsto t\coth(t)$ is non-decreasing over $\mathbb{R}_{+}$. Therefore, by
(Durmus et al.,, 2020, Lemma 10), $\mathrm{grad}\,f$ is geodesically
$\tilde{L}_{\pi}$-Lipschitz continuous on $\mathsf{S}$.In particular, for any
$\theta\in\mathsf{S}$,

$\left\|\mathrm{grad}\,f(\theta)\right\|_{\theta}\leq\tilde{L}_{\pi}\rho_{\Theta}({\theta^{\star}},\theta)\;.$
(168)

By (165) and (168), it is straightforward that $V=V_{2}$ and
$h=-\mathrm{grad}\,f$ satisfy H 2, with $C_{1}=0$ and
$C_{2}=2\tilde{L}^{2}_{\pi}\leq 2^{5}L^{2}_{\pi}$. In addition, by Proposition
4, (165) implies $V_{2}$ verifies H 3-$(\emptyset)$, with $\lambda=1/2$.

Finally, MD 1 holds with $\sigma_{0}^{2}={\rm D}^{2}$ and $\sigma_{1}^{2}=0$
since for any $\theta\in\mathsf{S}$ and $x\in\mathsf{X}$,

$\left\|H_{\theta}(x)\right\|_{\theta}=\left\|\mathrm{Exp}^{-1}_{\theta}(x)\right\|_{\theta}\leq
2{\mathrm{D}}\;.$ (169)

Therefore, we can apply Theorem 1-(c) which implies that for any
$\upeta\leq\overline{\upeta}$,

$\mathbb{E}[V_{2}(\theta_{n})]\leq\left\\{1-\upeta/4\right\\}^{n}V_{2}(\theta_{0})+4\upeta\,L_{\pi}{\mathrm{D}}^{2}\;.$
(170)

To conclude, it only remains to note that
$V_{2}(\theta_{n})=\rho_{\Theta}^{2}({\theta^{\star}},\theta_{n})$ and
$V_{2}(\theta_{0})=\rho_{\Theta}^{2}({\theta^{\star}},\theta_{0})$, since
$(\theta_{n})_{n\in\mathbb{N}}$ and ${\theta^{\star}}$ belong to
$\mathsf{H}=\mathsf{S}$.

###  D.5 Proof of Theorem 13

We consider in this section the recursion

$\displaystyle\theta_{n+1}$
$\displaystyle=\mathrm{Exp}_{\theta_{n}}\left[\upeta
H_{\theta_{n}}(X_{n+1})\right]$ (171) $\displaystyle H_{\theta_{n}}(X_{n+1})$
$\displaystyle=\left.\mathrm{Exp}^{-1}_{\theta_{n}}\left(X_{n+1}^{(1)}\right)\middle/\left(2\\{\rho_{\Theta}^{2}(\theta_{n},X_{n+1}^{(2)})/2+1\\}^{1/2}\right)\right.\;,$
(172)

where $X_{n+1}=(X_{n+1}^{(1)},X_{n+1}^{(2)}$ and
$(X_{n}^{(1)},X_{n}^{(2)})_{n\in\mathbb{N}^{*}}$ is an i.i.d. sequence of
pairs of independent random variables with distribution $\pi$. Denote by
$Q_{\upeta}$ the Markov kernel corresponding to (171).

We give first some additional intuition and motivation behind the scheme
(171). It can be interpreted as a stochastic optimization method to minimize

$\tilde{f}_{\pi}=(f_{\pi}+1)^{1/2}\;,$ (173)

in place of $f_{\pi}$. First note that $f_{\pi}$ and $\tilde{f}_{\pi}$ have
the same minimizer, but compared to $f_{\pi}$ it may be shown that
$\mathrm{grad}\,\tilde{f}_{\pi}$, given for any $\theta\in\Theta$ by

$\mathrm{grad}\,\tilde{f}_{\pi}(\theta)=(1/2)\mathrm{grad}\,f_{\pi}(\theta)(f_{\pi}(\theta)+1)^{-1/2}\;,$
(174)

is geodesically Lipschitz. However, note that (171) is not an unbiased
stochastic optimization scheme for the function $\tilde{f}_{\pi}$ since

$\mathbb{E}\left[H_{\theta_{n}}(X_{n+1})\right]=(1/2)\\{\mathrm{grad}\,f_{\pi}(\theta_{n})\\}\mathbb{E}\left[\\{\rho_{\Theta}^{2}(\theta_{n},X_{n+1}^{(2)})/2+1\\}^{-1/2}\right]\;.$
(175)

The proof of Theorem 13 then consists in adapting the proof of Theorem 1 to
deal with this additional difficulty taking for the Lyapunov function $V$,
$V_{1}$ defined by (10) with $\delta=1$. A general theory could be derived but
we believe that this is out the scope of the present document and leave it for
future work. We start by preliminary technical results which are needed to
establish Theorem 13.

######  Lemma 24.

Assume A 2 and MD 5. Let ${\theta^{\star}_{\pi}}$ be the Riemannian barycenter
of the probability measure $\pi$, i.e.
${\theta^{\star}_{\pi}}=\operatorname*{\arg\\!\min}_{\Theta}f_{\pi}$ where
$f_{\pi}$ is defined by (21). Then, for any $\theta\in\Theta$,

$-\int_{\Theta}\left\langle\mathrm{Exp}^{-1}_{\theta}({\theta^{\star}_{\pi}}),\mathrm{Exp}^{-1}_{\theta}(\nu)\right\rangle_{\theta}\pi(\mathrm{d}\nu)\leq-\rho_{\Theta}^{2}(\theta,{\theta^{\star}_{\pi}})/2\;.$
(176)

###### Proof.

Using A 2 and (Jost,, 2005, Theorem 5.6.1), we have that for any
$\nu\in\Theta$, the operator norm of the Riemannian Hessian of
$\theta\mapsto\rho_{\Theta}^{2}(\theta,\nu)/2$ is lower bounded by 1.
Therefore, by (Boumal,, 2020, Theorem 11.19),
$\theta\mapsto\rho_{\Theta}^{2}(\theta,\nu)/2$ is $1/2$-strongly convex.
Applying this to $\theta$ and ${\theta^{\star}_{\pi}}\in\Theta$, we have for
any $\nu\in\Theta$,

$\rho_{\Theta}^{2}({\theta^{\star}_{\pi}},\nu)/2-\rho_{\Theta}^{2}(\theta,\nu)/2\geq-\left\langle\mathrm{Exp}^{-1}_{\theta}({\theta^{\star}_{\pi}}),\mathrm{Exp}^{-1}_{\theta}(\nu)\right\rangle_{\theta}+\rho_{\Theta}^{2}(\theta,{\theta^{\star}_{\pi}})/2\;.$
(177)

Using MD 5, we can integrate this inequality w.r.t. $\pi$, bringing

$f_{\pi}({\theta^{\star}_{\pi}})-f_{\pi}(\theta)\geq-\int_{\Theta}\left\langle\mathrm{Exp}^{-1}_{\theta}({\theta^{\star}_{\pi}}),\mathrm{Exp}^{-1}_{\theta}(\nu)\right\rangle_{\theta}\pi(\mathrm{d}\nu)+\rho_{\Theta}^{2}(\theta,{\theta^{\star}_{\pi}})/2\;.$
(178)

Since by definition of ${\theta^{\star}_{\pi}}$, $0\geq
f_{\pi}({\theta^{\star}_{\pi}})-f_{\pi}(\theta)$, this completes the proof. ∎

######  Lemma 25.

Assume A 2 and MD 5. Let ${\theta^{\star}_{\pi}}$ be the Riemannian barycenter
of the probability measure $\pi$, i.e.
${\theta^{\star}_{\pi}}=\operatorname*{\arg\\!\min}_{\Theta}f_{\pi}$ where
$f_{\pi}$ is defined by (21). Then, for any $\theta\in\Theta$,

$\int_{\Theta}\\{\rho_{\Theta}^{2}(\theta,\nu)/2+1\\}^{-1/2}\pi(\mathrm{d}\nu)\geq\\{\rho_{\Theta}^{2}(\theta,{\theta^{\star}_{\pi}})+2f_{\pi}({\theta^{\star}_{\pi}})+1\\}^{-1/2}\;.$
(179)

###### Proof.

Let $\theta\in\Theta$. Using Jensen’s inequality with the convex function
$t\mapsto(t+1)^{-1/2}$ on $\mathbb{R}_{+}^{*}$, we have

$\int_{\Theta}\\{\rho_{\Theta}^{2}(\theta,\nu)/2+1\\}^{-1/2}\pi(\mathrm{d}\nu)\geq\left\\{f_{\pi}(\theta)+1\right\\}^{-1/2}\;.$
(180)

However, using the triangle and Hölder’s inequalities, we have for any
$\theta$ and $\nu\in\Theta$,
$\rho_{\Theta}^{2}(\theta,\nu)/2\leq\rho_{\Theta}^{2}(\theta,{\theta^{\star}_{\pi}})+\rho_{\Theta}^{2}({\theta^{\star}_{\pi}},\nu)$.
Taking the integral with respect to $\pi$, by MD 5 we get
$f_{\pi}(\theta)\leq\rho_{\Theta}^{2}(\theta,{\theta^{\star}_{\pi}})+2f_{\pi}({\theta^{\star}_{\pi}})$.
Lastly, combining this result with (180) and using that the function
$t\mapsto(t+1)^{-1/2}$ is non-increasing on $\mathbb{R}_{+}^{*}$ completes the
proof. ∎

######  Lemma 26.

Assume A 2 and MD 5. Let ${\theta^{\star}_{\pi}}$ be the Riemannian barycenter
of the probability measure $\pi$, i.e.
${\theta^{\star}_{\pi}}=\operatorname*{\arg\\!\min}_{\Theta}f_{\pi}$ where
$f_{\pi}$ is defined by (21). Then, for any $\theta_{0}\in\Theta$,

$Q_{\upeta}V_{1}(\theta_{0})\leq
V_{1}(\theta_{0})-[\upeta/(4C_{\pi}^{1/2})]D^{2}_{\Theta}(\theta_{0},{\theta^{\star}_{\pi}})+2\upeta^{2}(1+\kappa)\\{1+f_{\pi}({\theta^{\star}_{\pi}})\\}\left(f_{\pi}({\theta^{\star}_{\pi}})+2\right)\;,$
(181)

where $V_{1}$ is defined in (10) with $\delta\leftarrow 1$,
${\theta^{\star}}\leftarrow{\theta^{\star}_{\pi}}$,
$C_{\pi}=1+2f_{\pi}({\theta^{\star}_{\pi}})$ and
$D^{2}_{\Theta}:\Theta^{2}\to[0,1]$ is defined by (17).

###### Proof.

Let $\theta_{0}\in\Theta$, and consider

$H_{\theta_{0}}(X)=\left.(1/2)\mathrm{Exp}^{-1}_{\theta_{0}}\left(X^{(1)}\right)\middle/\left\\{\rho_{\Theta}^{2}\left(\theta_{0},X^{(2)}\right)\middle/2+1\right\\}^{1/2}\right.\;,$
(182)

where $X^{(1)},X^{(2)}$ are independent random variables with distribution
$\pi$.

Let $\upgamma:[0,1]\to\Theta$ be the geodesic curve defined by
$\upgamma:t\mapsto\mathrm{Exp}_{\theta_{0}}[t\upeta H_{\theta_{0}}(X)]$. Using
(Durmus et al.,, 2020, Lemma 1) with $\upgamma$ and $V_{1}$, we get

$\displaystyle V_{1}(\upgamma(1))$ $\displaystyle\leq
V_{1}(\theta_{0})+\left\langle\mathrm{grad}\,V_{1}(\theta_{0}),\dot{\upgamma}(0)\right\rangle_{\theta_{0}}+(L/2)\left\|\dot{\upgamma}(0)\right\|_{\theta_{0}}^{2}$
(183)
$\displaystyle=V_{1}(\theta_{0})+\upeta\left\langle\mathrm{grad}\,V_{1}(\theta_{0}),H_{\theta_{0}}(X)\right\rangle_{\theta_{0}}+((1+\kappa)\upeta^{2}/2)\left\|H_{\theta_{0}}(X)\right\|_{\theta_{0}}^{2}\;,$
(184)

by Proposition 5. We now compute the expectation of the terms in (184). Using
that $(X^{(1)},X^{(2)})$ are independent, we obtain

$\mathbb{E}\left[\left\langle\mathrm{grad}\,V_{1}(\theta_{0}),H_{\theta_{0}}(X)\right\rangle_{\theta_{0}}\right]=(1/2)\left\langle\mathrm{grad}\,V_{1}(\theta_{0}),\mathbb{E}\left[\mathrm{Exp}^{-1}_{\theta_{0}}\left(X^{(1)}\right)\right]\mathbb{E}\left[\left\\{\rho_{\Theta}^{2}\left(\theta_{0},X^{(2)}\right)\middle/2+1\right\\}^{-1/2}\right]\right\rangle_{\theta_{0}}\;.$
(185)

Moreover, using (60) and Lemmas 24 and 25 yields

$\displaystyle\mathbb{E}\left[\left\langle\mathrm{grad}\,V_{1}(\theta_{0}),H_{\theta_{0}}(X)\right\rangle_{\theta_{0}}\right]$
(186)
$\displaystyle\quad=-(1/2)\left\\{\rho_{\Theta}^{2}\left(\theta_{0},{\theta^{\star}_{\pi}}\right)+1\right\\}^{-1/2}\mathbb{E}\left[\left\langle\mathrm{Exp}^{-1}_{\theta_{0}}({\theta^{\star}_{\pi}}),\mathrm{Exp}^{-1}_{\theta_{0}}(X^{(1)})\right\rangle_{\theta_{0}}\right]\mathbb{E}\left[\left\\{\rho_{\Theta}^{2}(\theta_{0},X^{(2)})\middle/2+1\right\\}^{-1/2}\right]$
(187)
$\displaystyle\quad\leq-(1/4)\rho_{\Theta}^{2}\left(\theta_{0},{\theta^{\star}_{\pi}}\right)\left[\left\\{\rho_{\Theta}^{2}(\theta_{0},{\theta^{\star}_{\pi}})+1\right\\}\left\\{\rho_{\Theta}^{2}(\theta_{0},{\theta^{\star}_{\pi}})+2f_{\pi}({\theta^{\star}_{\pi}})+1\right\\}\right]^{-1/2}$
(188)
$\displaystyle\quad\leq-(16C_{\pi})^{-1/2}D^{2}_{\Theta}(\theta_{0},{\theta^{\star}_{\pi}})\;,$
(189)

where $C_{\pi}=1+2f_{\pi}({\theta^{\star}_{\pi}})$ and
$D^{2}_{\Theta}:\Theta^{2}\to[0,1]$ is defined by (17). Looking to bound the
expectation of the last term in (184), we use that
$\|\mathrm{Exp}^{-1}_{\theta_{0}}(X^{(1)})\|_{\theta_{0}}=\rho_{\Theta}(\theta_{0},X^{(1)})$
and that $X^{(1)}$ has distribution $\pi$ to obtain,

$\displaystyle\mathbb{E}\left[\left\|H_{\theta_{0}}(X)\right\|_{\theta_{0}}^{2}\right]$
$\displaystyle=(1/4)\mathbb{E}\left[\rho_{\Theta}^{2}(\theta_{0},X^{(1)})\right]\mathbb{E}\left[\left\\{\rho_{\Theta}^{2}(\theta_{0},X^{(2)})/2+1\right\\}^{-1}\right]$
(190)
$\displaystyle=(f_{\pi}(\theta_{0})/2)\mathbb{E}\left[\left\\{\rho_{\Theta}^{2}(\theta_{0},X^{(2)})/2+1\right\\}^{-1}\right]\;.$
(191)

Denote by $M=\rho_{\Theta}({\theta^{\star}_{\pi}},\theta_{0})/2$. We bound the
expectation in (191) using the event
$\\{\rho_{\Theta}({\theta^{\star}_{\pi}},X^{(2)})\geq M\\}$ and its
complement. On $\\{\rho_{\Theta}({\theta^{\star}_{\pi}},X^{(2)})\geq M\\}$, we
use Markov’s inequality with the increasing map $t\mapsto t^{2}/2+1$,

$\displaystyle\mathbb{E}\left[\mathbbm{1}_{[M,+\infty)}(\rho_{\Theta}({\theta^{\star}_{\pi}},X^{(2)}))/[\rho_{\Theta}^{2}(\theta_{0},X^{(2)})/2+1]\right]$
$\displaystyle\leq\mathbb{P}\left(\rho_{\Theta}({\theta^{\star}_{\pi}},X^{(2)})\geq
M\right)$ (192)
$\displaystyle\leq\left.\left(\mathbb{E}\left[\rho_{\Theta}^{2}({\theta^{\star}_{\pi}},X^{(2)})\right]\middle/2+1\right)\middle/\left(M^{2}/2+1\right)\right.\;.$
(193)

On $\\{\rho_{\Theta}({\theta^{\star}_{\pi}},X^{(2)})<M\\}$, using the triangle
inequality, we have

$\rho_{\Theta}(\theta_{0},X^{(2)})\geq|\rho_{\Theta}(\theta_{0},{\theta^{\star}_{\pi}})-\rho_{\Theta}({\theta^{\star}_{\pi}},X^{(2)})|\geq\rho_{\Theta}(\theta_{0},{\theta^{\star}_{\pi}})-M=M\;.$
(194)

Then, we obtain

$\mathbb{E}\left[\mathbbm{1}_{[0,M)}(\rho_{\Theta}(\theta,X^{(2)}))/\\{\rho_{\Theta}^{2}(\theta_{0},X^{(2)})/2+1\\}\right]\leq\left.1/[M^{2}/2+1]\right.\;.$
(195)

Adding (193) and (195) together and using the definition of $M$ we obtain,

$\mathbb{E}\left[\left\\{\rho_{\Theta}^{2}\left(\theta_{0},X^{(2)}\right)\middle/2+1\right\\}^{-1}\right]\leq\left.\left(f_{\pi}({\theta^{\star}_{\pi}})+2\right)\middle/\left[\rho_{\Theta}^{2}({\theta^{\star}_{\pi}},\theta_{0})/8+1\right]\right.\;.$
(196)

Plugging (196) in (191), we get

$\mathbb{E}\left[\left\|H_{\theta_{0}}(X)\right\|_{\theta_{0}}^{2}\right]\leq(f_{\pi}(\theta_{0})/2)\left(f_{\pi}({\theta^{\star}_{\pi}})+2\right)/\left[\rho_{\Theta}^{2}({\theta^{\star}_{\pi}},\theta_{0})/8+1\right]\;.$
(197)

Using the triangle and Hölder’s inequalities, we have for any $\theta$ and
$\nu\in\Theta$,
$\rho_{\Theta}^{2}(\theta,\nu)/2\leq\rho_{\Theta}^{2}(\theta,{\theta^{\star}_{\pi}})+\rho_{\Theta}^{2}({\theta^{\star}_{\pi}},\nu)$.
Taking the integral with respect to $\pi$, by MD 5 we get
$f_{\pi}(\theta)\leq\rho_{\Theta}^{2}(\theta,{\theta^{\star}_{\pi}})+2f_{\pi}({\theta^{\star}_{\pi}})$.
Combining this result and (197), we obtain

$\mathbb{E}\left[\left\|H_{\theta_{0}}(X)\right\|_{\theta_{0}}^{2}\right]\leq\\{\rho_{\Theta}^{2}({\theta^{\star}_{\pi}},\theta_{0})/2+f_{\pi}({\theta^{\star}_{\pi}})\\}\left(f_{\pi}({\theta^{\star}_{\pi}})+2\right)/\left[\rho_{\Theta}^{2}({\theta^{\star}_{\pi}},\theta_{0})/8+1\right]\leq
4\\{1+f_{\pi}({\theta^{\star}_{\pi}})\\}\left(f_{\pi}({\theta^{\star}_{\pi}})+2\right)\;.$
(198)

Combining this result and (189) in (184) concludes the proof. ∎

###### Proof of Theorem 13.

Let $\theta_{0}\in\Theta,\upeta>0$ and $n\in\mathbb{N}$. Then, for any
$k\in\\{1,\dots,n\\}$, using Markov’s property and Lemma 26 we have,

$\displaystyle[\upeta/(4C_{\pi}^{1/2})]\mathbb{E}\left[D^{2}_{\Theta}(\theta_{k-1},{\theta^{\star}_{\pi}})\right]$
$\displaystyle=[\upeta/(4C_{\pi}^{1/2})]\int_{\Theta}D^{2}_{\Theta}(\theta,{\theta^{\star}_{\pi}})Q_{\upeta}^{k-1}(\theta_{0},\mathrm{d}\theta)$
(199) $\displaystyle\leq
Q_{\upeta}^{k-1}V_{1}(\theta_{0})-Q_{\upeta}^{k}V_{1}(\theta_{0})+2\upeta^{2}(1+\kappa)(1+f_{\pi}({\theta^{\star}_{\pi}}))(f_{\pi}({\theta^{\star}_{\pi}})+2)\;.$
(200)

Summing these inequalities for $k\in\\{1,\dots,n\\}$ implies that

$[\upeta/(4C_{\pi}^{1/2})]\sum_{k=0}^{n-1}\mathbb{E}\left[D^{2}_{\Theta}(\theta_{k},{\theta^{\star}_{\pi}})\right]\leq
V_{1}(\theta_{0})-Q_{\upeta}^{n}V_{1}(\theta_{0})+2n\upeta^{2}(1+\kappa)(1+f_{\pi}({\theta^{\star}_{\pi}}))(f_{\pi}({\theta^{\star}_{\pi}})+2)\;.$
(201)

Finally, dividing both sides by $[n\upeta/(4C_{\pi}^{1/2})]$ and using that
$V_{1}$ is a non-negative function, we obtain

$n^{-1}\sum_{k=0}^{n-1}\mathbb{E}\left[D^{2}_{\Theta}(\theta_{k},{\theta^{\star}_{\pi}})\right]\leq\left.2V_{1}(\theta_{0})C_{\pi}^{1/2}\middle/(\upeta
n)\right.+2\upeta(1+\kappa)(f_{\pi}({\theta^{\star}_{\pi}})+1)(f_{\pi}({\theta^{\star}_{\pi}})+2)(2f_{\pi}({\theta^{\star}_{\pi}})+1)^{-1/2}\;.$
(202)

Which concludes the proof by setting
$B_{\pi}=(1+\kappa)(f_{\pi}({\theta^{\star}_{\pi}})+1)(f_{\pi}({\theta^{\star}_{\pi}})+2)(2f_{\pi}({\theta^{\star}_{\pi}})+1)^{-1/2}$.
∎

Figure 4: Monte Carlo approximations of the mean distance at convergence in
Theorem 13

Similarly to Figure 2, Figure 4 illustrates Theorem 7. To this end, 1000
replications of the experiment derived for Figure 3 are performed, obtaining
$\\{(\theta_{n}^{(i)})\,:\,i\in\\{1,\dots,1000\\}\\}$ for $n=\lceil
50/\upeta\rceil$ and $\upeta\in\\{1,2.8,4.6,6.4,8.2,10\\}\times 10^{-2}$. We
estimate, with these samples, the mean and the variance of
$D_{\Theta}^{2}(\theta,{\theta^{\star}_{\pi}})$, for $\theta$ following the
stationary distribution $\mu^{\upeta}$. We observe that the mean and variance
are both linear w.r.t. the step-size $\upeta$, indicating that the iterates of
the SA scheme remain in a neighborhood of diameter $\mathcal{O}(\upeta^{1/2})$
to the ground truth.

Even though the setting of this experiment goes beyond the assumptions of
Theorem 7, it suggests that such a result may be applicable also in the
setting of Theorem 13. The proof of such a result is left for future work.

##  Appendix E Background on Markov chain theory and Riemannian geometry

We give here some useful definitions and results that are used throughout the
paper.

###  E.1 Markov chain notions

We refer to Meyn and Tweedie, (2009) for a general introduction to Markov
chains in general state space. Let $(\mathsf{Y},\mathcal{Y})$ be a measurable
state space and $P$ be a Markov kernel on $\mathsf{Y}\times\mathcal{Y}$.
Consider for any $y\in\mathsf{Y}$, the distribution $\mathbb{P}_{y}$ of the
canonical Markov chain $(Y_{n})_{n\in\mathbb{N}}$ corresponding to $P$ and
starting from $y$ on the canonical space
$(\mathsf{Y}^{\mathbb{N}},\mathcal{Y}^{\otimes\mathbb{N}})$. Denote by
$\mathbb{E}_{y}$ the corresponding expectation.

Denote for any $\mathsf{A}\in\mathcal{Y}$, $\tau_{\mathsf{A}}=\inf\\{{l}\geq
1\,:\,Y_{l}\in\mathsf{A}\\}$ and
$N_{\mathsf{A}}=\sum_{{l}=1}^{+\infty}\mathbbm{1}_{\\{\mathsf{A}\\}}(Y_{l})$.

We say that $(Y_{n})_{n\in\mathbb{N}}$ is $\psi$-irreducible if there exists a
measure $\psi$ on $\mathcal{Y}$ such that whenever $\psi(\mathsf{A})>0$, we
have $\mathbb{P}_{y}(\tau_{\mathsf{A}}<\infty)>0$ for any $y\in\mathsf{Y}$.
Moreover, a set $\mathsf{A}\in\mathcal{Y}$ is called Harris-recurrent if
$\mathbb{P}_{y}(N_{\mathsf{A}}=\infty)=1$ for any $y\in\mathsf{A}$. Finally, a
chain $(Y_{n})_{n\in\mathbb{N}}$ is called Harris-recurrent if it is
$\psi$-irreducible and every set $\mathsf{A}\in\mathcal{Y}$ such that
$\psi(\mathsf{A})>0$ is Harris-recurrent.

Let $\overline{V}:\mathsf{Y}\to\left[1,+\infty\right)$. We say that $P$ is
$\overline{V}$-uniformly geometrically ergodic if there exist
$\rho\in\left[0,1\right)$ and $C\geq 0$ such that for any $y\in\mathsf{Y}$ and
$k\in\mathbb{N}$, $\left\|\updelta_{y}P^{k}-\mu\right\|_{\overline{V}}\leq
C\rho^{k}\overline{V}(y)$, where $\left\|\cdot\right\|_{\overline{V}}$ is
defined for two probability measures $\nu_{1},\nu_{2}$ on
$(\mathsf{Y},\mathcal{Y})$ by
$\left\|\nu_{1}-\nu_{2}\right\|_{\overline{V}}=\sup\\{\left|\nu_{1}(g)-\nu_{2}(g)\right|\,:\,\sup_{\mathsf{Y}}\\{\left|g\right|/\overline{V}\\}\leq
1\\}$.

###  E.2 Useful results from Riemannian geometry

We now give definitions and auxiliary results related to tensor fields along
curves, their derivatives, and Taylor expansions on Riemannian manifolds.

Let $\mathsf{M}$ be a smooth manifold with or without boundary. Given a smooth
curve $\upgamma:I\to\mathsf{M}$ defined on an interval $I$, and any
$k,l\in\mathbb{N}$, a $(k,l)$-tensor field along $\upgamma$ is a continuous
map $F:I\to\mathrm{T}^{(k,l)}\mathrm{T}\mathsf{M}$, such that
$F(t)\in\mathrm{T}^{(k,l)}(\mathrm{T}_{\upgamma(t)}\mathsf{M})$ for any $t\in
I$, where $\mathrm{T}^{(k,l)}\mathrm{T}\mathsf{M}$ is the bundle of
$(k,l)$-tensors on $\mathsf{M}$, see e.g. (Lee,, 2019, Appendix B). A vector
field $Y$ along $\upgamma$ is a $(1,0)$-tensor field, in which case for any
$t\in I$, $Y(t)$ is just a tangent vector in
$\mathrm{T}_{\upgamma(t)}\mathsf{M}$. We say that a tensor field $F$ along
$\upgamma$ is extendible if there exists a tensor field $\tilde{F}$ defined on
a neighborhood of $\upgamma(I)$ such that $F=\tilde{F}\circ\upgamma$.

We let $\mathfrak{X}^{k,l}(\upgamma)$ denote the set of smooth $(k,l)$-tensor
fields along $\upgamma$, and
$\mathfrak{X}(\upgamma)=\mathfrak{X}^{1,0}(\upgamma)$ denote the set of smooth
vector fields along $\upgamma$. In particular, $\mathfrak{X}^{0,0}(\upgamma)$
is the set of smooth functions $g:I\to\upgamma(I)\times\mathbb{R}$ such that
for any $t\in I$, $g(t)=(\upgamma(t),f(t))$ for some smooth function
$f:I\to\mathbb{R}$ and therefore can be identified with the set of smooth
functions $f:I\to\mathbb{R}$. In the sequel, we adopt if no confusion is
possible this identification. We extend to tensor fields along $\upgamma$ the
following definition of the trace on tensors. For any $(k,l)$-tensor $T$, we
denote by $\operatorname{Tr}_{\square,\triangle}(T)$ the $(k-1,l-1)$-tensor
with component of index $(i_{1},\ldots,i_{k-1},j_{1},\ldots,j_{l-1})$, given
by
$\sum_{m=1}^{d}T_{i_{1},\ldots,i_{\square-1},m,i_{\square},\ldots,i_{k-1}}^{j_{1},\ldots,j_{\triangle-1},m,j_{\triangle},\ldots,j_{l-1}}$.
In particular, for any
$\omega\in\mathfrak{X}^{0,1}(\upgamma),Y\in\mathfrak{X}(\upgamma)$,

$\operatorname{Tr}_{(1,1)}(\omega\otimes Y)=\omega(Y)\;.$ (203)

Also, for any $F\in\mathfrak{X}^{k,l}(\upgamma)$, any
$\omega^{1},\ldots,\omega^{k_{0}}\in\mathfrak{X}^{0,1}(\upgamma)$ and
$Y_{1},\dots,Y_{l_{0}}\in\mathfrak{X}(\upgamma)$, with $k_{0}\leq k,l_{0}\leq
l$, denote by $[F\,:\,\omega^{1}\otimes\cdots\otimes\omega^{k_{0}}\otimes
Y_{1}\otimes\cdots\otimes Y_{l_{0}}]$, the $(k-k_{0},l-l_{0})$ smooth tensor
field along $\upgamma$ defined by the induction:

$\displaystyle[F\,:\,\omega^{\otimes 1:i}]$
$\displaystyle=\operatorname{Tr}_{(1,l+1)}([F\,:\,\omega^{\otimes
1:(i-1)}]\otimes\omega^{i})$ (204) $\displaystyle[F\,:\,\omega^{\otimes
1:{k_{0}}}\otimes Y_{\otimes 1:j}]$
$\displaystyle=\operatorname{Tr}_{(k-{k_{0}}+1,1)}([F\,:\,\omega^{\otimes
1:{k_{0}}}\otimes Y_{\otimes 1:(j-1)}]\otimes Y_{j})\;,$ (205)

setting $\omega^{\otimes 1:i}=\omega^{1}\otimes\cdots\otimes\omega^{i}$,
$Y_{\otimes 1:j}=Y_{1}\otimes\cdots\otimes Y_{j}$. Note that for any
$\omega^{k-{k_{0}}+1},\ldots,\omega^{k}\in\mathfrak{X}^{0,1}(\upgamma)$ and
$Y_{l-{l_{0}}+1},\ldots,Y_{l}\in\mathfrak{X}(\upgamma)$,

$[F\,:\,\omega^{\otimes 1:{k_{0}}}\otimes Y_{\otimes
1:{l_{0}}}](\omega^{k-{k_{0}}+1},\ldots,\omega^{k},Y_{l-{l_{0}}+1},\ldots,Y_{l})=F(\omega^{1},\ldots,\omega^{k},Y_{1},\ldots,Y_{l})\;.$
(206)

######  Proposition 27.

Let $\mathsf{M}$ be a smooth manifold with or without border, $\nabla$ be a
connection on $\mathrm{T}\mathsf{M}$ and $\upgamma:I\to\mathsf{M}$ a smooth
curve defined on an interval $I$. Then, for any $k,l\in\mathbb{N}$, $\nabla$
determines an operator
$\mathrm{D}_{t}:\mathfrak{X}^{k,l}(\upgamma)\to\mathfrak{X}^{k,l}(\upgamma)$,
satisfying the following conditions.

  1. (a)

On $\mathfrak{X}(\upgamma)$, $\mathrm{D}_{t}$ is the usual covariant
derivative along $\upgamma$, see (Lee,, 2019, Theorem 4.24).

  2. (b)

On $\mathfrak{X}^{0,0}(\upgamma)$, $\mathrm{D}_{t}$ is the usual derivative
for real functions, i.e. for any $f\in\mathfrak{X}^{0,0}(\upgamma)$,
$\mathrm{D}_{t}f=\mathrm{d}f/\mathrm{d}t$.

  3. (c)

For any $F\in\mathfrak{X}^{k,l}(\upgamma)$, any
$\omega^{1},\dots,\omega^{k}\in\mathfrak{X}^{0,1}(\upgamma)$ and any
$Y_{1},\dots,Y_{l}\in\mathfrak{X}(\upgamma)$,

$\displaystyle(\mathrm{D}_{t}F)\left(\omega^{1},\dots,\omega^{k},Y_{1},\dots,Y_{l}\right)$
$\displaystyle=\frac{\mathrm{d}}{\mathrm{d}t}\left[F\left(\omega^{1},\dots,\omega^{k},Y_{1},\dots,Y_{l}\right)\right]$
(207)
$\displaystyle\quad-\sum_{i=1}^{k}F\left(\omega^{1},\dots,\omega^{i-1},\mathrm{D}_{t}\omega^{i},\omega^{i+1},\dots,\omega^{k},Y_{1},\dots,Y_{l}\right)$
$\displaystyle\quad-\sum_{j=1}^{l}F\left(\omega^{1},\dots,\omega^{k},Y_{1},\dots,Y_{j-1},\mathrm{D}_{t}Y_{j},Y_{j+1},\dots,Y_{l}\right)\;.$

In particular, $\mathrm{D}_{t}$ satisfies these additional properties.

  1. (i)

$\mathrm{D}_{t}$ satisfies the product rule, i.e. for any
$f\in\mathfrak{X}^{0,0}(\upgamma),F\in\mathfrak{X}^{k,l}(\upgamma),$

$\mathrm{D}_{t}\left(fF\right)=\left(\frac{\mathrm{d}}{\mathrm{d}t}f\right)F+f\mathrm{D}_{t}F\;.$
(208)

  2. (ii)

For any $k_{1},l_{1},k_{2},l_{2}\in\mathbb{N}$, and any
$F\in\mathfrak{X}^{k_{1},l_{1}}(\upgamma),G\in\mathfrak{X}^{k_{2},l_{2}}(\upgamma)$,

$\mathrm{D}_{t}(F\otimes G)=\mathrm{D}_{t}F\otimes
G+F\otimes\mathrm{D}_{t}G\;.$ (209)

  3. (iii)

For any positive integers $k_{0}\leq k,l_{0}\leq l$,
$F\in\mathfrak{X}^{k,l}(\upgamma)$,

$\mathrm{D}_{t}\left\\{\operatorname{Tr}_{(k_{0},l_{0})}(F)\right\\}=\operatorname{Tr}_{(k_{0},l_{0})}\left(\mathrm{D}_{t}F\right)\;.$
(210)

  4. (iv)

Let $F\in\mathfrak{X}^{k,l}$ be an extendible tensor field, i.e., such that
there exists a $(k,l)$-tensor field $\tilde{F}$ defined on a neighborhood of
$\upgamma(I)$ satisfying for any $t\in I$, $F(t)=\tilde{F}(\upgamma(t))$.
Then, for any $t\in I$,

$\mathrm{D}_{t}F(t)=\nabla_{\dot{\upgamma}(t)}\tilde{F}({\upgamma(t)})\;.$
(211)

Finally, if
$\tilde{\mathrm{D}}_{t}:\mathfrak{X}^{k,l}(\upgamma)\to\mathfrak{X}^{k,l}(\upgamma)$
is another operator satisfying (a),(b),(i),(ii) and (iii), then
$\mathrm{D}_{t}=\tilde{\mathrm{D}}_{t}$.

###### Proof.

Let $k,l\in\mathbb{N}$. Note first that (a)-(b) and (207) define
$\mathrm{D}_{t}F$ for any $F\in\mathfrak{X}^{k,l}(\upgamma)$, setting for any
$\omega\in\mathfrak{X}^{0,1}(\upgamma)$ and $Y\in\mathfrak{X}(\upgamma)$,

$\left[\mathrm{D}_{t}\omega\right](Y)={\mathrm{d}}\left[\omega(Y)\right]/{\mathrm{d}t}-\omega(\mathrm{D}_{t}Y)\;.$
(212)

We now show that $\mathrm{D}_{t}F\in\mathfrak{X}^{k,l}$, which will imply that
$\mathrm{D}_{t}:\mathfrak{X}^{k,l}\to\mathfrak{X}^{k,l}$. Second, we establish
that (i)-(ii)-(iii)-(iv) are satisfied. We conclude the proof by proving
uniqueness of $\mathrm{D}_{t}$.

Using (Lee,, 2019, Lemma B.6), to show that
$\mathrm{D}_{t}F\in\mathfrak{X}^{k,l}$ it is enough to prove that
$\mathrm{D}_{t}F$ is multilinear over $\mathfrak{X}^{0,0}(\upgamma)$. For
that, we start proving (i) on $\mathfrak{X}^{0,1}(\upgamma)$. Let
$\omega\in\mathfrak{X}^{0,1}(\upgamma),f\in\mathfrak{X}^{0,0}(\upgamma)$ and
$Y\in\mathfrak{X}(\upgamma)$, then by (212),

$\left[\mathrm{D}_{t}(f\omega)\right](Y)={\mathrm{d}}\left[f\omega(Y)\right]/{\mathrm{d}t}-f\omega\left(\mathrm{D}_{t}Y\right)=\left[{\mathrm{d}}f/{\mathrm{d}t}\right]\omega(Y)+f\left[\mathrm{D}_{t}\omega\right](Y)\;,$
(213)

which proves (i) on $\mathfrak{X}^{0,1}(\upgamma)$. Now, let
$k,l\in\mathbb{N},F\in\mathfrak{X}^{k,l}(\upgamma),\omega^{1},\dots,\omega^{k}\in\mathfrak{X}^{0,1}(\upgamma),Y_{1},\dots,Y_{l}\in\mathfrak{X}(\upgamma)$.
Let $f\in\mathfrak{X}^{0,0}(\upgamma)$ and $k_{0}\in\mathbb{N}^{*},k_{0}\leq
k$. We have, using the multilinearity of $F$ over
$\mathfrak{X}^{0,0}(\upgamma)$, the definition of $\mathrm{D}_{t}$ (207), and
(213)

$\displaystyle\left[\mathrm{D}_{t}F\right]\left(\omega^{1},\dots,\omega^{k_{0}-1},f\omega^{k_{0}},\omega^{k_{0}+1},\dots,\omega^{k},Y_{1},\dots,Y_{l}\right)$
(214)
$\displaystyle\qquad=\frac{\mathrm{d}}{\mathrm{d}t}\left[F\left(\omega^{1},\dots,\omega^{k_{0}-1},f\omega^{k_{0}},\omega^{k_{0}+1},\dots,\omega^{k},Y_{1},\dots,Y_{l}\right)\right]$
(215) $\displaystyle\qquad\qquad\qquad\quad-\sum_{i=1,i\neq
k_{0}}^{k}fF\left(\omega^{1},\dots,\omega^{i-1},\mathrm{D}_{t}\omega^{i},\omega^{i+1},\dots\omega^{k},Y_{1},\dots,Y_{l}\right)$
(216)
$\displaystyle\qquad\qquad\qquad\quad-F\left(\omega^{1},\dots,\omega^{k_{0}-1},\mathrm{D}_{t}(f\omega^{k_{0}}),\omega^{k_{0}+1},\dots,\omega^{k},Y_{1},\dots,Y_{l}\right)$
(217)
$\displaystyle\qquad\qquad\qquad\quad-\sum_{j=1}^{l}fF\left(\omega^{1},\dots,\omega^{k},Y_{1},\dots,Y_{j-1},\mathrm{D}_{t}Y_{j},Y_{j+1},\dots,Y_{l}\right)$
(218)
$\displaystyle\qquad=\left[\frac{\mathrm{d}}{\mathrm{d}t}f\right]\left\\{F\left(\omega^{1},\dots,\omega^{k},Y_{1},\dots,Y_{k}\right)-F\left(\omega^{1},\dots,\omega^{k},Y_{1},\dots,Y_{k}\right)\right\\}$
(219)
$\displaystyle\qquad\qquad\qquad\quad+f\left[\mathrm{D}_{t}F\right]\left(\omega^{1},\dots,\omega^{k},Y_{1},\dots,Y_{l}\right)$
(220)
$\displaystyle\qquad=f\left[\mathrm{D}_{t}F\right]\left(\omega^{1},\dots,\omega^{k},Y_{1},\dots,Y_{l}\right)\;.$
(221)

The same arguments apply if we replace $Y_{l_{0}}$ with $fY_{l_{0}}$, for some
$l_{0}\leq l$. Thus, using (Lee,, 2019, Lemma B.6),
$\mathrm{D}_{t}F\in\mathfrak{X}^{k,l}$.

Next, regarding (i), using the definition of $\mathrm{D}_{t}$,

$\displaystyle\left[\mathrm{D}_{t}fF\right]\left(\omega^{1},\dots,\omega^{k},Y_{1},\dots,Y_{l}\right)$
$\displaystyle=\left[\frac{\mathrm{d}}{\mathrm{d}t}f\right]F\left(\omega^{1},\dots,\omega^{k},Y_{1},\dots,Y_{l}\right)$
(222)
$\displaystyle\quad+f\left[\mathrm{D}_{t}F\right]\left(\omega^{1},\dots,\omega^{k},Y_{1},\dots,Y_{l}\right)\;,$
(223)

thus proving (i). Moreover, we prove (ii). Let
$k_{1},l_{1},k_{2},l_{2}\in\mathbb{N}$ and
$F\in\mathfrak{X}^{k_{1},l_{1}}(\upgamma),G\in\mathfrak{X}^{k_{2},l_{2}}(\upgamma)$,
$\omega^{1},\dots,\omega^{k_{1}+k_{2}}\in\mathfrak{X}^{0,1}(\upgamma),Y_{1},\dots,Y_{l_{1}+l_{2}}\in\mathfrak{X}(\upgamma)$.
Setting

$f=F(\omega^{1},\dots,\omega^{k_{1}},Y_{1},\dots,Y_{l_{1}})\text{ and
}g=G(\omega^{k_{1}+1},\dots,\omega^{k_{1}+k_{2}},Y_{l_{1}+1},\dots,Y_{l_{1}+l_{2}})\;,$

we have

$\displaystyle\left[\mathrm{D}_{t}(F\otimes
G)\right]\left(\omega^{1},\dots,\omega^{k_{1}+k_{2}},Y_{1},\dots,Y_{l_{1}+l_{2}}\right)$
(224)
$\displaystyle\quad=\frac{\mathrm{d}}{\mathrm{d}t}\left[fg\right]-\left[\sum_{i=1}^{k_{1}}F\left(\omega^{1},\dots,\omega^{i-1},\mathrm{D}_{t}\omega^{i},\omega^{i+1},\dots,\omega^{k_{1}},Y_{1},\dots,Y_{l_{1}}\right)\right.$
(225)
$\displaystyle\quad\quad+\left.\sum_{j=1}^{l_{1}}F\left(\omega^{1},\dots,\omega^{k_{1}},Y_{1},\dots,Y_{j-1},\mathrm{D}_{t}Y_{j},Y_{j+1},\dots,Y_{l_{1}}\right)\right]g$
(226)
$\displaystyle\qquad-f\left[\sum_{i=1}^{k_{2}}G\left(\omega^{k_{1}+1},\dots,\omega^{k_{1}+i-1},\mathrm{D}_{t}\omega^{k_{1}+i},\omega^{k_{1}+i+1},\dots,\omega^{k_{1}+k_{2}},Y_{l_{1}+1},\dots,Y_{l_{1}+l_{2}}\right)\right.$
(227)
$\displaystyle\qquad+\left.\sum_{j=1}^{l_{2}}G\left(\omega^{k_{1}+1},\dots,\omega^{k_{1}+k_{2}},Y_{l_{1}+1},\dots,Y_{l_{1}+j-1},\mathrm{D}_{t}Y_{l_{1}+j},Y_{l_{1}+j+1},\dots,Y_{l_{1}+l_{2}}\right)\right]$
(228)
$\displaystyle\quad=\left[\mathrm{D}_{t}F\right]\left(\omega^{1},\dots,\omega^{k_{1}},Y_{1},\dots,Y_{l_{1}}\right)g+f\left[\mathrm{D}_{t}G\right]\left(\omega^{k_{1}+1},\dots,\omega^{k_{1}+k_{2}},Y_{l_{1}+1},\dots,Y_{l_{1}+l_{2}}\right)$
(229) $\displaystyle\quad=\left[\mathrm{D}_{t}F\otimes
G+F\otimes\mathrm{D}_{t}G\right]\left(\omega^{1},\dots,\omega^{k_{1}+k_{2}},Y_{1},\dots,Y_{l_{1}+l_{2}}\right)\;,$
(230)

which proves (ii). Furthermore, to prove (iii), let $t_{0}\in I$ and
$(\mathbf{b}_{i})_{i\in\\{1,\dots,d\\}}$ be a basis of
$\mathrm{T}_{\upgamma(t_{0})}\Theta$. Using (a) and (Lee,, 2019, Theorem
4.32), define for any $i\in\\{1,\dots,d\\}$ and $t\in I$,

$e_{i}(t)=\mathrm{T}_{t_{0},t}^{\upgamma}\mathbf{b}_{i}\;,$ (231)

where $\mathrm{T}_{t_{0},t}^{\upgamma}$ denotes the parallel transport map
along $\upgamma$ from $\mathrm{T}_{\upgamma(t_{0})}\Theta$ to
$\mathrm{T}_{\upgamma(t)}\Theta$. As the parallel transport map is an
isomorphism, $(e_{i}(t))_{i\in\\{1,\dots,d\\}}$ is a basis of
$\mathrm{T}_{\upgamma(t)}\Theta$, for any $t\in I$. Therefore the family of
smooth vector fields $(e_{i})_{i\in\\{1,\dots,d\\}}$ is a parallel frame along
$\upgamma$ (with respect to $\nabla$). Denote
$(\varepsilon^{j})_{j\in\\{1,\dots,d\\}}$ its dual coframe. Using (212) on
$Y=e_{i},\omega=\varepsilon^{j}$, for any $i,j\in\\{1,\dots,d\\}$, shows that
the coframe $(\varepsilon^{j})_{j\in\\{1,\dots,d\\}}$ is parallel along
$\upgamma$. Note that for $(e_{i})_{i\in\\{1,\dots,d\\}}$ and
$(\varepsilon^{j})_{j\in\\{1,\dots,d\\}}$ to be well defined, we have used
$\nabla$, as well as the operator $\mathrm{D}_{t}$ on $\mathfrak{X}(\upgamma)$
and $\mathfrak{X}^{0,1}(\upgamma)$.

Let $k,l\in\mathbb{N}^{*}$ such that $k_{0}\leq k,l_{0}\leq l$, and let
$F\in\mathfrak{X}^{k,l}(\upgamma)$. There exist a family of functions
$\\{F^{j_{1},\dots,j_{l}}_{i_{1},\dots,i_{k}}\in\mathfrak{X}^{0,0}(\upgamma)\,:\,i_{1},\dots,i_{k},j_{1},\dots,j_{l}\in\\{1,\dots,d\\}\\}$
such that

$F=\sum_{i_{1},\dots,i_{k}=1}^{d}\sum_{j_{1},\dots,j_{l}=1}^{d}F^{j_{1},\dots,j_{l}}_{i_{1},\dots,i_{k}}\bigotimes_{\triangle=1}^{k}e_{i_{\triangle}}\bigotimes_{\square=1}^{l}\varepsilon^{j_{\square}}\;.$
(232)

Since the frame and its dual coframe are parallel along $\upgamma$, for any
$i\in\\{1,\dots,d\\}$ $\mathrm{D}_{t}e_{i}=0$ and
$\mathrm{D}_{t}\varepsilon^{i}=0$. Combining this fact with (i) and (ii) gives

$\mathrm{D}_{t}F=\sum_{i_{1},\dots,i_{k}=1}^{d}\sum_{j_{1},\dots,j_{l}=1}^{d}\left[\frac{\mathrm{d}}{\mathrm{d}t}F^{j_{1},\dots,j_{l}}_{i_{1},\dots,i_{k}}\right]\bigotimes_{\triangle=1}^{k}e_{i_{\triangle}}\bigotimes_{\square=1}^{l}\varepsilon^{j_{\square}}\;.$
(233)

Let $k_{0},l_{0}\in\mathbb{N}^{*}$ such that $k_{0}\leq k,l_{0}\leq l$, then
by definition of $\operatorname{Tr}_{(k_{0},l_{0})}$, for any
$i_{1},\dots,i_{k-1},j_{1},\dots,j_{l-1}\in\\{1,\dots,d\\}$,

$\operatorname{Tr}_{(k_{0},l_{0})}(F)^{j_{1},\dots,j_{l-1}}_{i_{1},\dots,i_{k-1}}=\sum_{m=1}^{d}F^{j_{1},\dots,j_{l_{0}-1},m,j_{l_{0}},\dots,j_{l-1}}_{i_{1},\dots,i_{k_{0}-1},m,i_{k_{0}},\dots,i_{k-1}}\;.$
(234)

We remind the reader that $\operatorname{Tr}_{(k_{0},l_{0})}(F)$ does not
depend on the choice of coordinates (Lee,, 2019, Appendix B). Thus, using
(233) and (234), we have

$\displaystyle\mathrm{D}_{t}\left[\operatorname{Tr}_{(k_{0},l_{0})}(F)\right]$
$\displaystyle=\sum_{i_{1},\dots,i_{k-1}=1}^{d}\sum_{j_{1},\dots,j_{l-1}=1}^{d}\frac{\mathrm{d}}{\mathrm{d}t}\left[\operatorname{Tr}_{(k_{0},l_{0})}(F)^{j_{1},\dots,j_{l-1}}_{i_{1},\dots,i_{k-1}}\right]\bigotimes_{\triangle=1}^{k-1}e_{i_{\triangle}}\bigotimes_{\square=1}^{l-1}\varepsilon^{j_{\square}}$
(235)
$\displaystyle=\sum_{i_{1},\dots,i_{k-1}=1}^{d}\sum_{j_{1},\dots,j_{l-1}=1}^{d}\sum_{m=1}^{d}\frac{\mathrm{d}}{\mathrm{d}t}F^{j_{1},\dots,j_{l_{0}-1},m,j_{l_{0}},\dots,j_{l-1}}_{i_{1},\dots,i_{k_{0}-1},m,i_{k_{0}},\dots,i_{k-1}}\bigotimes_{\triangle=1}^{k-1}e_{i_{\triangle}}\bigotimes_{\square=1}^{l-1}\varepsilon^{j_{\square}}$
(236)
$\displaystyle=\operatorname{Tr}_{(k_{0},l_{0})}\left(\mathrm{D}_{t}F\right)\;,$
(237)

thus proving (iii).

To prove (iv), first for any $f\in\mathfrak{X}^{(0,0)}(\upgamma)$, extendible
in $\tilde{f}$, we have by composition and definition of the covariant
derivative, that for any $t\in[0,1]$,

$(\mathrm{d}f/\mathrm{d}t)(t)=\mathrm{d}\tilde{f}_{\upgamma(t)}(\dot{\upgamma}(t))=\nabla_{\dot{\upgamma}(t)}\tilde{f}(\upgamma(t))\;.$
(238)

Also, using (Lee,, 2019, Theorem 4.24-(iii)) gives (iv) for any
$Y\in\mathfrak{X}(\upgamma)$. Combining (238), (212), its counterpart for
tensor fields defined over a manifold (Lee,, 2019, Proposition 4.15-(a)) and
(iv) over $\mathfrak{X}(\upgamma)$, proves (iv) over
$\mathfrak{X}^{(0,1)}(\upgamma)$. Now, for any $k,l\in\mathbb{N}$, using (iv)
over $\mathfrak{X}(\upgamma)$ and $\mathfrak{X}^{(0,1)}(\upgamma)$ combined
with (207) and its counterpart for tensor fields defined over a manifold
(Lee,, 2019, Equation (4.12)) gives (iv) over
$\mathfrak{X}^{(k,l)}(\upgamma)$.

Finally, we address uniqueness. Suppose now that $\tilde{\mathrm{D}}_{t}$ is
an operator on $\mathfrak{X}^{k,l}(\upgamma)$ that satisfies (a),(b),(i),(ii)
and (iii). First, (a) and (b) show that $\mathrm{D}_{t}$ and
$\tilde{\mathrm{D}}_{t}$ coincide on $\mathfrak{X}^{0,0}(\upgamma)$ and
$\mathfrak{X}(\upgamma)$. Second, for any
$Y\in\mathfrak{X}(\upgamma),\omega\in\mathfrak{X}^{0,1}(\upgamma)$, writing
$\omega(Y)=\operatorname{Tr}_{(1,1)}(Y\otimes\omega)$ and using (iii) gives

$\tilde{\mathrm{D}}_{t}\omega=\mathrm{d}[\omega(Y)]/\mathrm{d}t-\omega(\tilde{\mathrm{D}}_{t}Y)=\mathrm{D}_{t}\omega\;,$
(239)

using (212). Thus, $\tilde{\mathrm{D}}_{t}$ and $\mathrm{D}_{t}$ also agree on
$\mathfrak{X}^{0,1}(\upgamma)$. Therefore, the frame
$(e_{i})_{i\in\\{1,\dots,d\\}}$ and its dual coframe
$(\varepsilon^{j})_{j\in\\{1,\dots,d\\}}$ are also parallel with respect to
$\tilde{\mathrm{D}}_{t}$ along $\upgamma$. Let
$F\in\mathfrak{X}^{k,l}(\upgamma)$, then using (i) and (ii) shows that (233)
holds for the operator $\tilde{\mathrm{D}}_{t}$, proving that
$\mathrm{D}_{t}F=\tilde{\mathrm{D}}_{t}F$. This concludes the proof.

∎

######  Lemma 28.

Let $\mathsf{M}$ be a smooth manifold and $\nabla$ be a connection on
$\mathrm{T}\mathsf{M}$. Let $\upgamma:[0,1]\to\mathsf{M}$ be a smooth curve
and denote $\mathrm{D}_{t}$ the covariant derivative operator along $\upgamma$
associated with $\nabla$, defined in Proposition 27. Let
$F\in\mathfrak{X}^{k,l}(\upgamma)$,
$\omega^{1},\ldots,\omega^{k_{0}}\in\mathfrak{X}^{0,1}(\upgamma)$ and
$Y_{1},\dots,Y_{l_{0}}\in\mathfrak{X}(\upgamma)$, with ${k_{0}}\leq
k,{l_{0}}\leq l$. Then, we have

$\displaystyle\mathrm{D}_{t}\left([F:\omega^{\otimes 1:{k_{0}}}\otimes
Y_{\otimes 1:{l_{0}}}]\right)$ $\displaystyle=[\mathrm{D}_{t}F:\omega^{\otimes
1:{k_{0}}}\otimes Y_{\otimes 1:{l_{0}}}]$ (240)
$\displaystyle\qquad+\sum_{i=1}^{k_{0}}[F:\omega^{\otimes
1:(i-1)}\otimes\mathrm{D}_{t}\omega^{i}\otimes\omega^{(i+1):{k_{0}}}\otimes
Y_{\otimes 1:{l_{0}}}]$
$\displaystyle\qquad+\sum_{j=1}^{l_{0}}[F:\omega^{\otimes 1:{k_{0}}}\otimes
Y_{\otimes 1:(j-1)}\otimes\mathrm{D}_{t}Y_{j}\otimes
Y_{\otimes(j+1):{l_{0}}}]\;.$

###### Proof.

Let $F$ be a smooth $(k,l)$-tensor field along $\upgamma$. We show (240) by
induction. Following the recursive definition of the contraction in (204), we
prove it by induction on ${k_{0}}\in\mathbb{N}^{*},{k_{0}}\leq k$, for any
$\omega^{1},\dots,\omega^{k_{0}}\in\mathfrak{X}^{0,1}(\upgamma)$.

The case ${k_{0}}=1$ follows from Proposition 27-(ii) and (iii), combined with
the definition in (204),

$\displaystyle\mathrm{D}_{t}\left[F:\omega^{1}\right]$
$\displaystyle=\mathrm{D}_{t}\operatorname{Tr}_{(1,l+1)}(F\otimes\omega^{1})$
(241)
$\displaystyle=\operatorname{Tr}_{(1,l+1)}(\mathrm{D}_{t}[F\otimes\omega^{1}])$
(242)
$\displaystyle=\operatorname{Tr}_{(1,l+1)}(\mathrm{D}_{t}F\otimes\omega^{1}+F\otimes\mathrm{D}_{t}\omega^{1})$
(243)
$\displaystyle=\left[\mathrm{D}_{t}F:\omega^{1}\right]+\left[F:\mathrm{D}_{t}\omega^{1}\right]\;,$
(244)

where we have used the linearity of $\operatorname{Tr}$. Now assume there
exists ${k_{0}}\in\\{1,\dots,k-1\\}$ such that (240) holds for any smooth $1$
forms $\omega^{1},\dots,\omega^{k_{0}}$ and ${l_{0}}=0$. Moreover, consider
any smooth $1$ forms $\omega^{1},\dots,\omega^{{k_{0}}+1}$. Then, using the
same arguments as for the case ${k_{0}}=1$ and the induction hypothesis, we
obtain

$\displaystyle\mathrm{D}_{t}\left[F:\omega^{\otimes 1:({k_{0}}+1)}\right]$
$\displaystyle=\mathrm{D}_{t}\operatorname{Tr}_{(1,l+1)}\left(\left[F:\omega^{\otimes
1:{k_{0}}}\right]\otimes\omega^{{k_{0}}+1}\right)$ (245)
$\displaystyle=\operatorname{Tr}_{(1,l+1)}\left(\mathrm{D}_{t}\left[F:\omega^{\otimes
1:{k_{0}}}\right]\otimes\omega^{{k_{0}}+1}\right)+\operatorname{Tr}_{(1,l+1)}\left(\left[F:\omega^{\otimes
1:{k_{0}}}\right]\otimes\mathrm{D}_{t}\omega^{{k_{0}}+1}\right)$ (246)
$\displaystyle=\operatorname{Tr}_{(1,l+1)}\left(\left[\mathrm{D}_{t}F:\omega^{\otimes
1:{k_{0}}}\right]\otimes\omega^{{k_{0}}+1}\right)+\left[F:\omega^{\otimes
1:{k_{0}}}\otimes\mathrm{D}_{t}\omega^{{k_{0}}+1}\right]$ (247)
$\displaystyle\quad+\sum_{i=1}^{{k_{0}}}\operatorname{Tr}_{(1,l+1)}\left(\left[F:\omega^{\otimes
1:(i-1)}\otimes\mathrm{D}_{t}\omega^{i}\otimes\omega^{\otimes(i+1):{k_{0}}}\right]\otimes\omega^{{k_{0}}+1}\right)$
(248) $\displaystyle=\left[\mathrm{D}_{t}F:\omega^{\otimes
1:({k_{0}}+1)}\right]+\sum_{i=1}^{{k_{0}}+1}\left[F:\omega^{\otimes
1:(i-1)}\otimes\mathrm{D}_{t}\omega^{i}\otimes\omega^{\otimes(i+1):({k_{0}}+1)}\right]\;.$
(249)

Subsequently, using the recursive definition of the contraction in (205), we
prove (240) by induction on ${l_{0}}\in\mathbb{N}^{*},{l_{0}}\leq l$ for any
${k_{0}}\leq k$ and any
$\omega^{1},\dots,\omega^{k_{0}}\in\mathfrak{X}^{0,1}(\upgamma)$. Let
$Y_{1}\in\mathfrak{X}(\upgamma)$. Then, using once again Proposition 27-(ii)
and (iii), (205), and (240) in the case ${l_{0}}=0$ justified above, the case
${l_{0}}=1$ is proven as follows,

$\displaystyle\mathrm{D}_{t}\left[F:\omega^{\otimes 1:{k_{0}}}\otimes
Y_{1}\right]$
$\displaystyle=\operatorname{Tr}_{(k-{k_{0}}+1,1)}\left(\mathrm{D}_{t}\left\\{\left[F:\omega^{\otimes
1:{k_{0}}}\right]\otimes Y_{1}\right\\}\right)$ (250)
$\displaystyle=\operatorname{Tr}_{(k-{k_{0}}+1,1)}\left(\left[\mathrm{D}_{t}F:\omega^{\otimes
1:{k_{0}}}\right]\otimes Y_{1}\right)+\left[F:\omega^{\otimes
1:{k_{0}}}\otimes\mathrm{D}_{t}Y_{1}\right]$ (251)
$\displaystyle\quad+\sum_{i=1}^{k_{0}}\operatorname{Tr}_{(k-{k_{0}}+1,1)}\left(\left[F:\omega^{\otimes
1:(i-1)}\otimes\mathrm{D}_{t}\omega^{i}\otimes\omega^{\otimes(i+1):{k_{0}}}\right]\otimes
Y_{1}\right)$ (252) $\displaystyle=\left[\mathrm{D}_{t}F:\omega^{\otimes
1:{k_{0}}}\otimes Y_{1}\right]+\sum_{i=1}^{{k_{0}}}\left[F:\omega^{\otimes
1:(i-1)}\otimes\mathrm{D}_{t}\omega^{i}\otimes\omega^{\otimes(i+1):{k_{0}}}\otimes
Y_{1}\right]$ (253) $\displaystyle\quad+\left[F:\omega^{\otimes
1:{k_{0}}}\otimes\mathrm{D}_{t}Y_{1}\right]\;.$ (254)

Furthermore, assume there exists ${l_{0}}\in\\{1,\dots,l-1\\}$ such that (240)
holds for any ${k_{0}}\leq k$, any
$\omega^{1},\dots,\omega^{k_{0}}\in\mathfrak{X}^{0,1}(\upgamma)$ and any
$Y_{1},\dots,Y_{l_{0}}\in\mathfrak{X}(\upgamma)$. Let
$Y_{1},\dots,Y_{{l_{0}}+1}\in\mathfrak{X}(\upgamma)$. Then using the same
arguments as for the case ${l_{0}}=1$ and the induction hypothesis, we obtain

$\displaystyle\mathrm{D}_{t}\left[F:\omega^{\otimes 1:{k_{0}}}\otimes
Y_{\otimes 1:({l_{0}}+1)}\right]$ (255)
$\displaystyle\qquad=\operatorname{Tr}_{(k-{k_{0}}+1,1)}\left(\mathrm{D}_{t}\left\\{\left[F:\omega^{\otimes
1:{k_{0}}}\otimes Y_{\otimes 1:{l_{0}}}\right]\otimes
Y_{{l_{0}}+1}\right\\}\right)$ (256)
$\displaystyle\qquad=\operatorname{Tr}_{(k-{k_{0}}+1,1)}\left(\left[\mathrm{D}_{t}F:\omega^{\otimes
1:{k_{0}}}\otimes Y_{\otimes 1:{l_{0}}}\right]\otimes Y_{{l_{0}}+1}\right)$
(257)
$\displaystyle\qquad\quad+\sum_{i=1}^{k_{0}}\operatorname{Tr}_{(k-{k_{0}}+1,1)}\left(\left[F:\omega^{\otimes
1:(i-1)}\otimes\mathrm{D}_{t}\omega^{i}\otimes\omega^{\otimes(i+1):{k_{0}}}\otimes
Y_{\otimes 1:{l_{0}}}\right]\otimes Y_{{l_{0}}+1}\right)$ (258)
$\displaystyle\qquad\quad+\sum_{j=1}^{l_{0}}\operatorname{Tr}_{(k-{k_{0}}+1,1)}\left(\left[F:\omega^{\otimes
1:{k_{0}}}\otimes Y_{\otimes 1:(j-1)}\otimes\mathrm{D}_{t}Y_{j}\otimes
Y_{\otimes(j+1):{l_{0}}}\right]\otimes Y_{{l_{0}}+1}\right)$ (259)
$\displaystyle\qquad\quad+\operatorname{Tr}_{(k-{k_{0}}+1,1)}\left(\left[F:\omega^{\otimes
1:{k_{0}}}\otimes Y_{\otimes
1:{l_{0}}}\right]\otimes\mathrm{D}_{t}Y_{{l_{0}}+1}\right)$ (260)
$\displaystyle\qquad=\left[\mathrm{D}_{t}F:\omega^{\otimes 1:{k_{0}}}\otimes
Y_{\otimes 1:({l_{0}}+1)}\right]+\sum_{i=1}^{k_{0}}\left[F:\omega^{\otimes
1:(i-1)}\otimes\mathrm{D}_{t}\omega^{i}\otimes\omega^{\otimes(i+1):{k_{0}}}\otimes
Y_{\otimes 1:({l_{0}}+1)}\right]$ (261)
$\displaystyle\qquad\quad+\sum_{j=1}^{{l_{0}}+1}\left[F:\omega^{\otimes
1:{k_{0}}}\otimes Y_{\otimes 1:(j-1)}\otimes\mathrm{D}_{t}Y_{j}\otimes
Y_{\otimes(j+1):({l_{0}}+1)}\right]\;,$ (262)

which concludes the proof. ∎

######  Theorem 29.

Let $\mathsf{M}$ be a smooth manifold and $\nabla$ be a connection on
$\mathrm{T}\mathsf{M}$. Let $\upgamma:[0,1]\to\mathsf{M}$ be a geodesic and
$Y:\mathsf{M}\to\mathrm{T}\mathsf{M}$ a smooth vector field. Then, for any
$t\in[0,1],n\in\mathbb{N}$,

$\displaystyle\mathrm{T}_{t0}^{\upgamma}Y(\upgamma(t))$
$\displaystyle=\sum_{{k_{0}}=0}^{n}(t^{k_{0}}/{k_{0}}!)\left.\nabla^{k_{0}}Y_{\upgamma(0)}\left(\dot{\upgamma}(0),\dots,\dot{\upgamma}(0)\right)\right.$
(263)
$\displaystyle\quad+\int_{0}^{t}\left.[(t-s)^{n}/n!]\mathrm{T}_{s0}^{\upgamma}\nabla^{n+1}Y_{\upgamma(s)}\left(\dot{\upgamma}(s),\dots,\dot{\upgamma}(s)\right)\right.\mathrm{d}s\;,$

where
$\mathrm{T}_{t0}^{\upgamma}:\mathrm{T}_{\upgamma(t)}\mathsf{M}\to\mathrm{T}_{\upgamma(0)}\mathsf{M}$
is the parallel transport map along $\upgamma$, and the $(1,{k_{0}})$-tensor
field $\nabla^{k_{0}}Y$ is the total derivative of order ${k_{0}}$ of the
$(1,0)$-tensor field $Y$.

For a definition of the total covariant derivative, see (Lee,, 2019,
Proposition 4.15). Also, in (263), remark that even though $\dot{\upgamma}$ is
only a vector field along $\upgamma$, and not a vector field, the value of a
vector field $\nabla_{X}Y$ evaluated at $\theta\in\mathsf{M}$ only depends on
$X(\theta)$ and on values of $Y$ along smooth curves
$\mathrm{c}:[0,1]\to\mathsf{M}$ satisfying $\mathrm{c}(0)=\theta$ and
$\dot{\mathrm{c}}(0)=X(\theta)$; by (Lee,, 2019, Proposition 4.26). Therefore
the expression
$\nabla^{k_{0}}Y_{\upgamma(t)}(\dot{\upgamma}(t),\dots,\dot{\upgamma}(t))$ in
Theorem 29 is well defined for any ${k_{0}}\in\mathbb{N},t\in[0,1]$.

###### Proof.

Consider $\mathpzc{V}:[0,1]\to\mathsf{M}$ the smooth vector field along
$\upgamma$ and the function
$\varphi:[0,1]\to\mathrm{T}_{\upgamma(0)}\mathsf{M}$ defined by

$\mathpzc{V}=Y\circ\upgamma\text{ and
}\;\varphi:t\mapsto\mathrm{T}_{t0}^{\upgamma}\mathpzc{V}(t)\;.$ (264)

Then we check by induction on $n\in\mathbb{N}^{*}$ that $\varphi$ is $n$-times
differentiable with derivative of order $n$ given for any $t\in[0,1]$ by
$\varphi^{(n)}(t)=\mathrm{T}_{t0}^{\upgamma}[\mathrm{D}_{t}^{n}\mathpzc{V}(t)]$
and
$\mathrm{D}_{t}^{n}\mathpzc{V}(t)=\nabla^{n}Y_{\upgamma(t)}(\dot{\upgamma}(t),\dots,\dot{\upgamma}(t))$,
where $\mathrm{D}_{t}$ is the covariant derivative operator along $\upgamma$
with respect to the connection $\nabla$, defined in Proposition 27.

First, the case $n=1$ is a direct application of (Lee,, 2019, Theorem 4.34,
Theorem 4.24) since $Y$ is an extension of $\mathpzc{V}$. Assume now that the
property holds for $n\in\mathbb{N}^{*}$. Then, for any $t_{0},t\in[0,1],t\neq
t_{0}$, we have

$\left.\left[\varphi^{(n)}(t)-\varphi^{(n)}(t_{0})\right]\middle/\left(t-t_{0}\right)\right.=\mathrm{T}_{t_{0}0}^{\upgamma}\left.\left[\mathrm{T}_{tt_{0}}^{\upgamma}\mathrm{D}^{n}_{t}\mathpzc{V}(t)-\mathrm{D}^{n}_{t}\mathpzc{V}(t_{0})\right]\middle/\left(t-t_{0}\right)\right.\;.$
(265)

Now (Lee,, 2019, Theorem 4.34) ensures that the limit of the quantity above
exists when $t\to t_{0}$ and in addition this limit is

$\varphi^{(n+1)}(t_{0})=\mathrm{T}_{t_{0}0}^{\upgamma}\mathrm{D}^{n+1}_{t}\mathpzc{V}(t_{0})\;,$
(266)

which shows that $\varphi$ is $n+1$ times differentiable on
$\left[0,1\right]$. We now show that for any $t\in\left[0,1\right]$,
$\mathrm{D}_{t}^{n+1}\mathpzc{V}(t)=\nabla^{n+1}Y_{\upgamma(t)}(\dot{\upgamma}(t),\dots,\dot{\upgamma}(t))$.
Using Lemma 28 on the smooth $(1,n)$-tensor field along $\upgamma$
$F=(\nabla^{n}Y)\circ\upgamma$, taking ${k_{0}}=0$ and $n$ times the vector
field $\dot{\upgamma}$, we have

$\mathrm{D}_{t}\left[F:\dot{\upgamma}\otimes\cdots\otimes\dot{\upgamma}\right]=\left[\mathrm{D}_{t}F:\dot{\upgamma}\otimes\cdots\otimes\dot{\upgamma}\right]\;,$
(267)

since $\mathrm{D}_{t}\dot{\upgamma}=0$ because $\upgamma$ is a geodesic. Also,
by (206),
$[\mathrm{D}_{t}F:\dot{\upgamma}\otimes\cdots\otimes\dot{\upgamma}]=\mathrm{D}_{t}F(\dot{\upgamma},\dots,\dot{\upgamma})$.
Finally, as $\nabla^{n}Y$ is an extension of $F$, using the induction
hypothesis and the definition of the total derivative give for any
$t\in[0,1]$,

$\mathrm{D}^{n+1}_{t}\mathpzc{V}(t)=\mathrm{D}_{t}F\left(\dot{\upgamma},\dots,\dot{\upgamma}\right)(t)=\nabla_{\dot{\upgamma}(t)}(\nabla^{n}Y)_{\upgamma(t)}\left(\dot{\upgamma}(t),\dots,\dot{\upgamma}(t)\right)\\\
=(\nabla^{n+1}Y)_{\upgamma(t)}\left(\dot{\upgamma}(t),\dots,\dot{\upgamma}(t)\right)\;,$
(268)

concluding the induction.

Finally, (263) is simply a consequence of Taylor’s formula with integral
remainder of the vectorial valued function $\varphi$ identifying
$\mathrm{T}_{\upgamma(0)}\mathsf{M}$ with $\mathbb{R}^{d}$. ∎

######  Proposition 30.

Let $\mathsf{M}$ be a smooth manifold, $\nabla$ be a symmetric connection
defined over the smooth vector fields of $\mathsf{M}$. For any smooth function
$f:\mathsf{M}\to\mathbb{R}$ and any local coordinates
$(u_{i})_{i\in\\{1,\dots,d\\}}$, we have

$\nabla\mathrm{Hess}\,f=\sum_{i,j,k=1}^{d}\left\\{\partial^{3}_{kij}f-\sum_{l=1}^{d}\left[\Gamma_{ij}^{l}\partial^{2}_{kl}f+\Gamma_{ki}^{l}\partial^{2}_{jl}f+\Gamma_{kj}^{l}\partial^{2}_{il}f\right]-\sum_{m=1}^{d}\partial_{k}\Gamma_{ij}^{m}\partial_{m}f\right.\\\
\left.+\sum_{l,m=1}^{d}\left[\Gamma_{kj}^{l}\Gamma_{il}^{m}+\Gamma_{ki}^{l}\Gamma_{lj}^{m}\right]\partial_{m}f\right\\}\mathrm{d}u^{i}\otimes\mathrm{d}u^{j}\otimes\mathrm{d}u^{k}\;,$
(269)

where $(\Gamma_{ij}^{k})_{i,j,k\in\\{1,\dots,d\\}}$ are the Christoffel
symbols in these local coordinates, the local frame and its dual coframe are
denoted by $(\partial u_{i})_{i\in\\{1,\dots,d\\}}$ and
$(\mathrm{d}u^{j})_{j\in\\{1,\dots,d\\}}$.

###### Proof.

Let $(u_{i})_{i\in\\{1,\ldots,d\\}}$ be local coordinates. By (Lee,, 2019,
Example 4.22), in this chart, we have

$\mathrm{Hess}\,f=\sum_{i,j=1}^{d}F_{ij}\mathrm{d}u^{i}\otimes\mathrm{d}u^{j}\;,\text{where
for any
}i,j\in\\{1,\dots,d\\}\;,\;F_{ij}=\partial_{ij}^{2}f-\sum_{m=1}^{d}\Gamma_{ij}^{m}\partial_{m}f\;.$
(270)

Applying (Lee,, 2019, Proposition 4.18) on $\mathrm{Hess}\,f$, we obtain that
$\nabla\mathrm{Hess}\,f=\sum_{i,j,k=1}^{d}G_{ijk}\mathrm{d}u^{i}\otimes\mathrm{d}u^{j}\otimes\mathrm{d}u^{k}$,
where for any $i,j,k\in\\{1,\ldots,d\\}$,

$G_{ijk}=\partial_{k}F_{ij}-\sum_{l=1}^{d}\left(\Gamma_{kj}^{l}F_{il}+\Gamma_{ki}^{l}F_{lj}\right)\;.$
(271)

Expanding the expression above using (270) gives for any
$i,j,k\in\\{1,\dots,d\\}$,

$G_{ijk}=\partial^{3}_{ijk}f-\sum_{m=1}^{d}\left(\partial_{k}\Gamma_{ij}^{m}\partial_{m}f+\Gamma_{ij}^{m}\partial^{2}_{km}f\right)-\sum_{l=1}^{d}\Gamma_{kj}^{l}\left(\partial^{2}_{il}f-\sum_{m=1}^{d}\Gamma_{il}^{m}\partial_{m}f\right)\\\
-\sum_{l=1}^{d}\Gamma_{ki}^{l}\left(\partial^{2}_{lj}f-\sum_{m=1}^{d}\Gamma_{lj}^{m}\partial_{m}f\right)\;.$
(272)

The desired result is obtained by reordering this equation, which concludes
the proof. ∎