Isometries preserve geodesics

Your calculation looks like an attempt to prove the naturality of the Levi-Civita connection, the fact that @Zhen Lin implicitly points to. In the settings of the question it can be stated as $$ \nabla^g_X{Y}=f^* \left( \nabla^{(f^{-1})^* g}_{\operatorname{d}f(X)} \operatorname{d}f(Y) \right) $$

Notice also that in fact you are using two different connections: one for vector fields along $\gamma \colon I \rightarrow M$ induced from $\nabla^g$ on $M$, and another one for vector fields along $f \circ \gamma \colon I \rightarrow N$ induced from $\nabla^h$ on $N$. Due to the naturality property they agree, and it may be helpful to distinguish $D^g_t:=\nabla^g_{\frac{d}{dt}{\gamma}}$ and $D^h_t:=\nabla^h_{\frac{d}{dt}(f \circ \gamma)}$ in the present calculation. Indeed, using $$\frac{\operatorname{d}}{\operatorname{d}t}(f\circ \gamma)=\operatorname{d}f(\frac{\operatorname{d}}{\operatorname{d}t}\gamma) $$ we get $$ D^h_t{\frac{\operatorname{d}(f\circ\gamma)}{\operatorname{d}t}} = f_* \left( D^g_t{\frac{\operatorname{d}\gamma}{\operatorname{d}t}} \right) =0 $$

It is actually an exercise in the Lee's book, I try to do it by following the hint. First, you have to understand the naturality of Riemannian connection, then everything will be clear. I like using $\nabla_{\frac{d}{dt}}$ instead of $D_t$ here.

First, Define an operator $\varphi^*\tilde{\nabla}_{\frac{d}{dt}}:\mathcal{J}(\gamma)\rightarrow\mathcal{J}(\gamma)$ by $(\varphi^*\tilde{\nabla}_{\frac{d}{dt}})V=\varphi^*(\tilde{\nabla}_{\frac{d}{dt}}(\varphi_*V))$. Suppose $V$ is a vector field along $\gamma.$ Show this operator satisfies three properties in the Lemma 4.9. in Lee's Riemannian Geometry.

(1) Linearity over $\mathcal{R}$:$$(\varphi^*\tilde{\nabla}_{\frac{d}{dt}})(aV+bW)=a(\varphi^*\tilde{\nabla}_{\frac{d}{dt}})V+b(\varphi^*\tilde{\nabla}_{\frac{d}{dt}})W$$ (2) Product rule: $$(\varphi^*\tilde{\nabla}_{\frac{d}{dt}})(fV)=\frac{df}{dt}V+f(\varphi^*\tilde{\nabla}_{\frac{d}{dt}})V.$$ (3) If $V$ is extendible, then for any extension $\tilde{V}$ of $V$, $$(\varphi^*\tilde{\nabla}_{\frac{d}{dt}})V=(\varphi^*\tilde{\nabla}_{\dot{\gamma}(t)})\tilde{V}.$$ Then by uniqueness, the operator we defined above is just the unique operator $\nabla_{\frac{d}{dt}}$. that is $$\nabla_{\frac{d}{dt}}V=\varphi^*(\tilde{\nabla}_{\frac{d}{dt}}(\varphi_*V)) \ \text{or} \ \ \varphi_*(\nabla_{\frac{d}{dt}}V)=\tilde{\nabla}_{\frac{d}{dt}}(\varphi_*V).$$

Now If $\gamma$ is the geodesic in $M$ with initial $p$ and initial velocity $V$, i.e., $\nabla_{\frac{d}{dt}}\dot{\gamma}(t)=0$. Obviously $\varphi\circ\gamma$ is a curve in $\tilde{M}$ with initial point $\varphi(p)$ and initial velocity $\varphi_*V$, moreover, it is also the geodesic since $$\tilde{\nabla}_{\frac{d}{dt}}\dot{\varphi}(\gamma(t))=\tilde{\nabla}_{\frac{d}{dt}}\varphi_*\dot{\gamma}(t)=\varphi_*(\nabla_{\frac{d}{dt}}\dot{\gamma}(t))=0.$$