For a Gaussian process I run I have chosen to use the Matern52 covariance function, as from research this is practical to use for physical data as it does not overally smooth the function. However is there a numerical way to prove this suits my data?
What would be the benefits of using Matern32 or Matern52, as I am quite sure that the Matern class is the correct thing to use for my data, as it allows for flexibility, whereas the other stationary covariance functions have been described to have no flexibility with regard to local behaviour and act like a prior
For choosing Matern52 or Matern32, you can turn to model selection methods. Perhaps there is a way to do cross-validation or some other information criteria that can help you in choosing. The Rasmussen + Williams GPML book has lots of information about this. You can also look at prior predictive samples too, and visually see which sets of samples best matches the characteristics you see in your data.
Also, you could try this question out somewhere like https://stats.stackexchange.com/, since this isn’t necessarily PyMC related and there you’ll find a much broader community.