Sampling from Posterior Distribution

[sgbaird]

If I have the mean and covariance functions (via MATLAB's fitrgp() ), how can I sample from the posterior?

[lawrennd]

Hi Sterling, fitrgp didn't exist when we were writing this code (MATLAB didn't have native GP implementations) so I don't know how it works (and since then we moved to python as a group).

But in gpmat the posterior mean and variance are computed through gpPosteriorMeanVar, details here:

https://github.com/SheffieldML/GPmat/blob/master/gp/gpPosteriorMeanVar.m

[sgbaird]

Hi Neil, that makes sense. I'll be following suit since I just joined a group based mostly in python (these questions are related to my previous research). This seems like a good resource. Since I have the mean values and covariance matrix already from MATLAB's fitrgp(), It seems like gsamp.m from NetLab does what I'm looking for in generating samples of the posterior distribution like in the example given in the README file:
>> gpPosteriorSample('rbf', 5, [1 1], [-3 3], 1e5)

I'm now wondering if MATLAB's mvnrnd() function does something similar to gsamp.m. But either way, I can try it out and if it's not what I'm expecting I can probably use gsamp.m. Thanks!

[lawrennd]

Yes, I suspect it does, netlab dates back to around 1996 when MATLAB had very poor support for stats and machine learning.

I suspect they've covered most of the gaps since then. GPmat was developed from about 2003, and again there were big gaps in MATLAB (or silly toolboxes that required extra licenses).

The history also explains why it is written in a particular way ... when I first was implementing it in 2003, the MATLAB object oriented stuff was incredibly slow, so it's all implemented with feval to do late binding of 'methods'. Which can make it quite hard to follow. But the basic structure carries forward to the GPy software, and therefore influences other Python GP packages (like GPflow).

[sgbaird]

Thanks for the comments, I think that's great you had an implementation before it became standard. It also seems that doing some things that are standard in GPmat are undocumented functions in fitrgp(), e.g. obtaining the posterior covariance matrix to sample from the posterior, which has made it a bit difficult 😅, and hence exploring other options like GPmat. Again, thank you for the help.