Let $(M,h)$ be a compact connected smooth Riemannian manifold with $\partial M=\emptyset$ and $h$ smooth metric. Let $p_t(x,y)$ designate its heat kernel $(t>0,x,y\in M)$.
In the Euclidean case, one would have that $$ \nabla_x p_t(x,y) = -\nabla_y p_t(x,y)\,. $$
My question is: does a similar identity hold on $M$, through parallel transport (since the left hand side belongs to $T_xM$ and the right hand side to $T_yM$)?
(I am aware of Bismut formula which enables one to represent $\nabla_x p_t(x,y)$ in a probabilistic setting, but I don't know how to proceed from it.)