# Leave One Out Cross-validation and RMSE

**URL:** <https://discourse.pymc.io/t/leave-one-out-cross-validation-and-rmse/88>\
**Category:** Questions\
**Tags:** from\_github\
**Created:** [June 16, 2017, 9:27am UTC](https://discourse.pymc.io/t/leave-one-out-cross-validation-and-rmse/88 "2017-06-16T09:27:53Z")\
**Posts on this page:** 12\
**Page:** 1

<div class="post-metadata">

**Author:** ![marcodena](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/marcodena/32/41_2.png) [@marcodena](https://discourse.pymc.io/u/marcodena)\
**Post date:** [June 16, 2017, 9:27am UTC](https://discourse.pymc.io/t/leave-one-out-cross-validation-and-rmse/88/1 "2017-06-16T09:27:53Z")

</div>

Hello,  
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?

Thanks

---

<div class="post-metadata">

**Author:** ![junpenglao](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/junpenglao/32/8_2.png) [@junpenglao](https://discourse.pymc.io/u/junpenglao)\
**Post date:** [June 16, 2017, 9:39am UTC](https://discourse.pymc.io/t/leave-one-out-cross-validation-and-rmse/88/2 "2017-06-16T09:39:12Z")

</div>

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.

---

<div class="post-metadata">

**Author:** ![marcodena](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/marcodena/32/41_2.png) [@marcodena](https://discourse.pymc.io/u/marcodena)\
**Post date:** [June 16, 2017, 9:41am UTC](https://discourse.pymc.io/t/leave-one-out-cross-validation-and-rmse/88/3 "2017-06-16T09:41:43Z")

</div>

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 🙂

---

<div class="post-metadata">

**Author:** ![aloctavodia](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/aloctavodia/32/8642_2.png) [@aloctavodia](https://discourse.pymc.io/u/aloctavodia)\
**Post date:** [June 16, 2017, 10:29am UTC](https://discourse.pymc.io/t/leave-one-out-cross-validation-and-rmse/88/4 "2017-06-16T10:29:58Z")

</div>

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](https://arxiv.org/abs/1507.04544).

---

<div class="post-metadata">

**Author:** ![junpenglao](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/junpenglao/32/8_2.png) [@junpenglao](https://discourse.pymc.io/u/junpenglao)\
**Post date:** [June 16, 2017, 10:34am UTC](https://discourse.pymc.io/t/leave-one-out-cross-validation-and-rmse/88/5 "2017-06-16T10:34:25Z")

</div>

I see.  
Looking at the [LOO paper](https://arxiv.org/pdf/1507.04544.pdf) they compare the RMSE as well, so there might be a way to do it easily.

---

<div class="post-metadata">

**Author:** ![marcodena](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/marcodena/32/41_2.png) [@marcodena](https://discourse.pymc.io/u/marcodena)\
**Post date:** [June 16, 2017, 10:35am UTC](https://discourse.pymc.io/t/leave-one-out-cross-validation-and-rmse/88/6 "2017-06-16T10:35:36Z")

</div>

I read that paper and I saw also this [http://www.g3journal.org/content/6/10/3107.full](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)

☹

---

<div class="post-metadata">

**Author:** ![junpenglao](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/junpenglao/32/8_2.png) [@junpenglao](https://discourse.pymc.io/u/junpenglao)\
**Post date:** [June 16, 2017, 10:41am UTC](https://discourse.pymc.io/t/leave-one-out-cross-validation-and-rmse/88/7 "2017-06-16T10:41:19Z")

</div>

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

```python
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`:

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

```

---

<div class="post-metadata">

**Author:** ![marcodena](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/marcodena/32/41_2.png) [@marcodena](https://discourse.pymc.io/u/marcodena)\
**Post date:** [June 16, 2017, 10:46am UTC](https://discourse.pymc.io/t/leave-one-out-cross-validation-and-rmse/88/8 "2017-06-16T10:46:21Z")

</div>

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 :))

---

<div class="post-metadata">

**Author:** ![junpenglao](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/junpenglao/32/8_2.png) [@junpenglao](https://discourse.pymc.io/u/junpenglao)\
**Post date:** [June 16, 2017, 10:53am UTC](https://discourse.pymc.io/t/leave-one-out-cross-validation-and-rmse/88/9 "2017-06-16T10:53:11Z")

</div>

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

---

<div class="post-metadata">

**Author:** ![aloctavodia](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/aloctavodia/32/8642_2.png) [@aloctavodia](https://discourse.pymc.io/u/aloctavodia)\
**Post date:** [June 16, 2017, 10:57am UTC](https://discourse.pymc.io/t/leave-one-out-cross-validation-and-rmse/88/10 "2017-06-16T10:57:59Z")

</div>

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.

---

<div class="post-metadata">

**Author:** ![marcodena](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/marcodena/32/41_2.png) [@marcodena](https://discourse.pymc.io/u/marcodena)\
**Post date:** [June 16, 2017, 1:53pm UTC](https://discourse.pymc.io/t/leave-one-out-cross-validation-and-rmse/88/11 "2017-06-16T13:53:41Z")

</div>

@junpenglao added to github  
@aloctavodia thanks!

---

<div class="post-metadata">

**Author:** ![junpenglao](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/junpenglao/32/8_2.png) [@junpenglao](https://discourse.pymc.io/u/junpenglao)\
**Post date:** [June 16, 2017, 2:20pm UTC](https://discourse.pymc.io/t/leave-one-out-cross-validation-and-rmse/88/12 "2017-06-16T14:20:54Z")

</div>

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.