It appears to me that one reason why nobody has proved the formula yet is that the formula is still wrong. First, the formula has to depend on $X$ and $Y$. If you rescale $X$ and $Y$, the left side of the formula scales but the right side stays constant. That can't be. Second, the two sides of the equation do not scale the same under a constant scaling of the metric.

I consider the derivation of the correct version to be a reasonable if challenging exercise for a serious graduate student in differential geometry, so I was expecting someone else to provide the details. You can do this using only the basic definitions and properties of a Riemannian metric, its connection, and Riemann curvature with the fundamental theorem of calculus and the product rule for differentiation. Although I learned most of my Riemannian geometry after I was out of graduate school, I spent many, many hours doing calculations and arguments like this over and over again. Almost all of global Riemannian geometry involves working with Jacobi fields using arguments like the one used to prove this local formula.

But I got tired of waiting, so I wrote out all the details. If you're a student, I recommend that you try to read as little of my proof as possible or just scan it quickly and try to finish it yourself.

Warning: I wrote this up very quickly and did not check for typos and errors. It's possible that my final formula is still not right, but I am confident that my argument can be used to obtain a correct formula. I also did not provide every last detail, so, if you're unfamiliar with an argument like this, you need to do a lot of work making sure that everything really works. The key trick is pulling everything back to the unit square, where elementary calculus can be used. I'm sure this trick can be replaced by Stokes' theorem on the manifold itself, but that's too sophisticated for my taste.

[Holonomy calculation][1]

ADDED:

The correct formula, if you assume $|X\wedge Y| = 1$, is

$P_\gamma Z - Z = Area(c) R(X,Y)Z$

This scales properly when you rescale the metric by a constant factor. Notice that the left side is invariant under rescaling of the metric.

I recommend looking at papers written by Hermann Karcher, especially the one with Jost on almost linear functions, the one with Heintze on a generalized comparison theorem, and the one on the Riemannian center of mass. I haven't looked at this or anything else in a long time, but I have the impression that I learned a lot about how to work with Jacobi fields and Riemann curvature from these papers.

Finally, don't worry about citing anything I've said or wrote. Just write up your own proof of whatever you need. If it happens to look very similar to what I wrote, that's OK. I consider all of this "standard stuff" that any good Riemannian geometer knows, even if they would say it differently from me.

  [1]: http://www.math.poly.edu/~yang/holonomy.pdf