We are working on a Bayesian model for a space-time process, and are using a No-U-Turn sampler (NUTS) that requires a model for the log-probability and it's gradient with respect to the model parameters. More succinctly, we have fairly complicated log-probability function $f : \mathbb{R}^n \rightarrow \mathbb{R}$, involving statistical distributions, kronecker products, exponentials, ratios, if-else statements etc, and need to provide it and it's gradient to NUTS. Several packages (Stan and Julia's MCMC) use operator-overloading (to the best of my knowledge) to obtain the gradient automatically.
If we were able to create our own gradient function, perhaps using a source-code-transformation auto-diff tool, would we get better performance, or is OO just as good or better?