Volume of small geodesic balls; Uniform bounds vs Landau notation

Question

Let $(M,g)$ be a smooth compact Riemannian manifold without boundary. Then for $\epsilon< inj_{M}$, the injectivity radius of $M$, there is the expansion formula for small geodesic balls, $p \in M$ $$ \, \,\frac{\mbox{Vol}(B_\epsilon(p)\subset M)}{\mbox{Vol}(B_\epsilon(p)\subset\mathbb{R}^n)}=1-\frac{R(p)}{6(n+2)}\epsilon^2+O(\epsilon^4).$$ To get this, a orthonormal basis $e_1,...,e_n $ of $T_pM$ is chosen. Then in the corresponding normal coordinate system in $p$, so that in this coordinate system $$ \int_{B_\epsilon(p)} d\mu= \int_{B_\epsilon(0)} ((\exp_p)^{-1})^{\ast}d\mu= \int_{B_\epsilon(0)} (1 + \frac{\operatorname{Ric}_{ij} x_ix_j}{6} + O(|x|^3)) dx $$ where $d\mu$ is the canonical volume form of $M$ and $$ (*) \,|O(|x|^3)| \leq |x|^3 \sup_{x \in B_\epsilon(0)} |D^3\left(((\exp_p)^{-1})^{\ast}d\mu\right)| $$ How do I argue that the error term $(*)$ is independent of the chosen coordinate system? And, if it possible, is there some way to compute $(*)$ in terms of some geometric invariants of $M$?