Leave One Out Cross-validation and RMSE

I’m trying to compare the results of my model with a model that uses a frequentist approach. This frequentist model uses a LOO Cross-validation and gives the RMSE as result.

Thus, my naive solution is to fit the model with (n-1) training examples, change the shared variables, and use sample_ppc to “predict”. Finally I compute the RMSE.
This is a looong process (suggest me alternatives in case) and I have some errors about shared variabiles.

using this:

X_shared = shared(X[features].values[train_index, :])
mu_phi = CAR2('mu_phi', adjacency=adj_shared, rho=rho, tau=tau, shape=X_shared.get_value().shape[0])

X_shared.set_value(X[features].values[test_group_indexes, :])

I have that mu_phi’s shape doesn’t get updated (n = 76 but it should be n=1). I don’t find how to do this, so I am not sure whether it is a bug.

Update: pm.stats.loo is different from what I need, but maybe I can use it to compute the RMSE between posteriors (in LOO) and observed data. How?


I dont know much about the LOO, @aloctavodia has more experience on that.
However, my comment would be: are you sure you want to use RMSE? It’s not a good measurement for model fit in general except GLMs.

I’m using a GLM (~ Poisson) and I have to compare my predictive results with another paper without Implementing it. I would not use the RMSE otherwise :slight_smile:

Hi @marcodena The LOO implemented in PyMC3 is actually PSIS-LOO (Pareto smoothed importance sampling-LOO) which is an approximation to the exact LOO that you are trying to compute. PSIS-LOO (and WAIC) are both approximation for estimating pointwise out-of-sample prediction accuracy. Under certain circumstances these information criteria should be proportional to the RMSE, I think the proportionality only holds for normal models.

For the comparison you want to perform i think the right approach is what you are doing, I am not sure there is a shortcut here, maybe do you want to check this paper.

I see.
Looking at the LOO paper they compare the RMSE as well, so there might be a way to do it easily.

I read that paper and I saw also this http://www.g3journal.org/content/6/10/3107.full. Despite the fact that it could be proportional, I don’t find easy to get in pymc3 the expected value of the posterior in a cross-validation way (“Importance sampling” section in the paper)


One way to go about of your problem of the shape is to create a new model:

Testdata = X[features].values[test_group_indexes, :]
with Model() as Prediction_Model:
    mu_phi = CAR2('mu_phi', adjacency=adj_shared, rho=rho, tau=tau, shape=Testdata.shape[0])

Then you sample using the training set as before to get trace, and generate sample using the Prediction_Model:

ppc = sample_ppc(trace, model=Prediction_Model, ...)

Yes, although I think there is a bug in pymc3 if it does not get a shared variable updated. This was the reason I posted in github.

I hope there is another way to go (@aloctavodia?) , and a bugfix :))

Hmm, it is indeed a bug if the shared variable is not updated, could you try to reproduce it with a minimal example?

I can’t think of another approach other than the direct comparison (as you were doing), but I will keep thinking about this and I will let you know of any news.

@junpenglao added to github
@aloctavodia thanks!

Saw it, thanks. I am not sure how to approach that, maybe somebody else can give it a stab. I think in pm.minibatch you usually set a total_size but, otherwise I dont think this is a easy to solve issue.