I am looking for a reference with a Taylor expansion of the metric tensor in the normal coordinates. The coefficients should be written in terms of $\mathrm{Rm}, \nabla\mathrm{Rm}, \nabla^2\mathrm{Rm},\dots$ at the point.
It is easy to get using Jacobi equation, but I would prefer to have a reference (if it exists).
Comment. I need a tiny bit of this formula, namely the terms with the components of $\nabla^{n-2}\mathrm{Rm}$ in the coefficients for the monomials of degree $n$.
Using the reference https://arxiv.org/pdf/0903.2087.pdf, which agrees with https://arxiv.org/pdf/hep-th/0001078v1.pdf, which agrees with the reference U. Müller, C. Schubert and Anton M. E. van de Ven, J. Gen. Rel. Grav. 31 (1999) 1759-1768 [arXiv], one sees that the expansion in normal coordinates about a point $p\in M$ of a metric $g$ is given by $$ g_{ij} = g_{ij}(0) + \sum_{n\ge 2} \frac1{n!} C_{ijk_1\cdots k_n}\,x^{k_1}\cdots x^{k_n} $$ where $C_{ijk_1\cdots k_n}$ is the value at $p$ of a tensor of the form $$ -\frac{2(n{-}1)}{n{+}1}\,\nabla_{k_3}\cdots\nabla_{k_n} R_{ik_1jk_2} + LOT_n $$ where $LOT_n$ is a polynomial in the curvature tensor and its derivatives of order strictly less than $n{-}2$. (Note that, in particular, $LOT_2 = 0$.)
Note: The first two references only verify the crucial coefficient $-2(n{-}1)/(n{+}1)$ as above for $n=2,3,4,5$, but, in any case, the coefficients for $n=2$ and $n=3$ are well-known.
Meanwhile, if one does not want to rely on these sources, one can note that these coefficients clearly cannot depend on the dimension of the underlying manifold, so it suffices to verify them in the simplest nontrivial case, i.e., when the dimension of the manifold is $2$. In this case, using the Gauss Lemma, a metric $g$ in geodesic normal coordinates $(x,y)$ centered on $p$ takes the form $$ g = \mathrm{d}x^2 + \mathrm{d}y^2 + h(x,y)\bigl(x\,\mathrm{d}y-y\,\mathrm{d}x)^2, $$ where the function $h$ is arbitrary, subject to the condition that $(x^2{+}y^2)h(x,y)+1>0$.
Letting $r^2 = x^2 + y^2$ and letting $T$ be the radial vector field $x\,\partial_x + y\,\partial_y$, one computes the formula for the Gauss curvature of $g$ to be $$ K = -\frac{2(1+r^2h)(TTh) - r^2(Th)^2+2(5+3r^2h)(Th) + 8r^2h^2+12h}{4(1+r^2h)^2}. $$ Of course, one has the formula for the Riemann curvature tensor in the form $$ R = K\,(1{+}r^2h)\,(\mathrm{d}x\wedge\mathrm{d}y)\otimes(\mathrm{d}x\wedge\mathrm{d}y), $$ while the tensor $(1{+}r^2h)\,(\mathrm{d}x\wedge\mathrm{d}y)\otimes(\mathrm{d}x\wedge\mathrm{d}y)$ is parallel with respect to the Levi-Civita connection. Then, by assuming that $h$ vanishes to order $n{-}2\ge0$, one can easily verify that the above coefficient is correct.
