# 2d Gaussian process IndexError

**URL:** <https://discourse.pymc.io/t/2d-gaussian-process-indexerror/3244>\
**Category:** Questions\
**Created:** [May 13, 2019, 9:18am UTC](https://discourse.pymc.io/t/2d-gaussian-process-indexerror/3244 "2019-05-13T09:18:06Z")\
**Posts on this page:** 14\
**Page:** 1

<div class="post-metadata">

**Author:** ![madarshahian](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/madarshahian/32/722_2.png) [@madarshahian](https://discourse.pymc.io/u/madarshahian)\
**Post date:** [May 13, 2019, 9:18am UTC](https://discourse.pymc.io/t/2d-gaussian-process-indexerror/3244/1 "2019-05-13T09:18:06Z")

</div>

This is an example similar to an example for applying 2d gaussian process by @fonnesbeck. X\_.shape=(10970, 2) and Y\_.shape = (10970,). After running it I get this error at the line related to posterior predictive:

> if samples.shape[0] == 1 and size == 1:
> 
> IndexError: tuple index out of range

Wondering if you could help me to resolve it?thanks!

```
#%%
X_shared=th.shared(X_)
Y_shared=th.shared(Y_)
#%%
nd = 30
z1, z2 = np.meshgrid(np.linspace(X_[:,0].min(), X_[:,0].max(), nd), np.linspace(X_[:,1].min(), X_[:,1].max(), nd))
Z = np.concatenate([z1.reshape(nd*nd, 1), z2.reshape(nd*nd, 1)], 1)

nd = 15
xu1, xu2 = np.meshgrid(np.linspace(X_[:,0].min(), X_[:,0].max(), nd), np.linspace(X_[:,1].min(), X_[:,1].max(), nd))
Xu = np.concatenate([xu1.reshape(nd*nd, 1), xu2.reshape(nd*nd, 1)], 1)
#%%
with pm.Model() as spatial_model:
    
    l = pm.HalfCauchy("l", beta=3, shape=(2,))
    sf2 = pm.HalfCauchy("sf2", beta=3)
    sn2 = pm.HalfCauchy("sn2", beta=3)

# K = pm.gp.cov.ExpQuad(2, l) * sf2**2
    K = pm.gp.cov.Matern32(2,l) * sf2**2
    
    gp_spatial = pm.gp.MarginalSparse(cov_func=K, approx="FITC")
    obs = gp_spatial.marginal_likelihood("obs", X=X_shared, Xu=Xu, y=Y_shared, noise=sn2)

    trace = pm.sample(cores=1)
#%%

with spatial_model:
    f_pred = gp_spatial.conditional('f_pred2', Z)
    samples = pm.sample_posterior_predictive([trace], vars=[f_pred], samples=100)
#%%

```

---

<div class="post-metadata">

**Author:** ![madarshahian](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/madarshahian/32/722_2.png) [@madarshahian](https://discourse.pymc.io/u/madarshahian)\
**Post date:** [May 15, 2019, 2:12pm UTC](https://discourse.pymc.io/t/2d-gaussian-process-indexerror/3244/2 "2019-05-15T14:12:23Z")

</div>

To add more to this question here is what I get after running the code:  
 ![1_pic](https://canada1.discourse-cdn.com/flex036/uploads/pymc3/original/2X/1/184c29092d4d5f1cf47324d2e2a4e407fd569c93.png)

 ![2_pic](https://canada1.discourse-cdn.com/flex036/uploads/pymc3/original/2X/6/6400755be39c33e8e9340dac23242aa50be70401.png)

---

<div class="post-metadata">

**Author:** ![Gon\_F](https://avatars.discourse-cdn.com/v4/letter/g/f14d63/32.png) [@Gon\_F](https://discourse.pymc.io/u/Gon_F)\
**Post date:** [May 16, 2019, 8:25pm UTC](https://discourse.pymc.io/t/2d-gaussian-process-indexerror/3244/3 "2019-05-16T20:25:55Z")

</div>

You would really have to post where your data comes from so anyone that wants to help could actually run your code.

---

<div class="post-metadata">

**Author:** ![madarshahian](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/madarshahian/32/722_2.png) [@madarshahian](https://discourse.pymc.io/u/madarshahian)\
**Post date:** [May 16, 2019, 10:06pm UTC](https://discourse.pymc.io/t/2d-gaussian-process-indexerror/3244/4 "2019-05-16T22:06:24Z")

</div>

@Gon\_FThank you for the message, but the problem should not be the data, I explained the data shape in my explanation. To explain more, X\_ , input is longitude and altitude, and output is a measurement for each of them. You can assume a grid that y is 2_X\_[:,0]+3_X\_[:,1]+noise. If it really matters I would be happy to share the actual data as well.

---

<div class="post-metadata">

**Author:** ![Gon\_F](https://avatars.discourse-cdn.com/v4/letter/g/f14d63/32.png) [@Gon\_F](https://discourse.pymc.io/u/Gon_F)\
**Post date:** [May 17, 2019, 5:20pm UTC](https://discourse.pymc.io/t/2d-gaussian-process-indexerror/3244/5 "2019-05-17T17:20:21Z")

</div>

I just mentioned it because it gives me peace of mind to know the error is not coming from the data : P , and it allows me to work along with the code and run it, instead of just abstractly analyzing it.

---

<div class="post-metadata">

**Author:** ![madarshahian](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/madarshahian/32/722_2.png) [@madarshahian](https://discourse.pymc.io/u/madarshahian)\
**Post date:** [May 17, 2019, 6:23pm UTC](https://discourse.pymc.io/t/2d-gaussian-process-indexerror/3244/6 "2019-05-17T18:23:45Z")

</div>

[data.csv](https://canada1.discourse-cdn.com/flex036/uploads/pymc3/original/2X/b/b5f7cab871e96c68b5f6898c71273e2047143e7c.csv) (1.5 MB)  
@Gon_F here is the data, there are some missing values and I used (pandas dropna) to discard them.  
Thanks!

---

<div class="post-metadata">

**Author:** ![madarshahian](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/madarshahian/32/722_2.png) [@madarshahian](https://discourse.pymc.io/u/madarshahian)\
**Post date:** [May 20, 2019, 3:19pm UTC](https://discourse.pymc.io/t/2d-gaussian-process-indexerror/3244/7 "2019-05-20T15:19:37Z")

</div>

I am wondering if the code looks good? or there is an error in it that I couldn’t fine? What is the reason for the error?

---

<div class="post-metadata">

**Author:** ![fonnesbeck](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/fonnesbeck/32/13_2.png) [@fonnesbeck](https://discourse.pymc.io/u/fonnesbeck)\
**Post date:** [May 23, 2019, 1:18am UTC](https://discourse.pymc.io/t/2d-gaussian-process-indexerror/3244/8 "2019-05-23T01:18:31Z")

</div>

The model does run on my system, using current master (macOS and Python 3.7.3). Can you try that (or the release candidate that was just pushed)?

---

<div class="post-metadata">

**Author:** ![madarshahian](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/madarshahian/32/722_2.png) [@madarshahian](https://discourse.pymc.io/u/madarshahian)\
**Post date:** [June 3, 2019, 1:36pm UTC](https://discourse.pymc.io/t/2d-gaussian-process-indexerror/3244/9 "2019-06-03T13:36:22Z")

</div>

Hey,  
Sorry I missed your message, and thanks so much for your attention; after seeing it I updated my pymc3. For some reason I couldn’t use the developer version using git clone. So I have just updated it to PyMC3 version v.3.7 from v.3.6.  
The code now works, but I still have two problems after this update:  
1- this warning still is there, and I do not know how to solve it:  
`

> C:\Users\xxxxx\AppData\Local\Continuum\anaconda3\envs\pymc3env\lib\site-packages\theano\tensor\basic.py:6611: FutureWarning: Using a non-tuple sequence for multidimensional indexing is deprecated; use `arr[tuple(seq)]` instead of `arr[seq]`. In the future this will be interpreted as an array index, `arr[np.array(seq)]`, which will result either in an error or a different result.

`  
2-

> File “C:\Users\xxxxx\AppData\Local\Continuum\anaconda3\envs\pymc3env\lib\site-packages\pymc3\plots\__init_\_.py”, line 22, in **call**  
> “ArviZ is not installed. In order to use `{0.attr}`:\npip install arviz”.format(self)
> 
> ImportError: ArviZ is not installed. In order to use `plot_forest`:

This Arviz error is really annoying! Last time when this happened I had to re install whole Anaconda. Even when I install it, and installation is done successfully, it keep giving me this error.

---

<div class="post-metadata">

**Author:** ![Gon\_F](https://avatars.discourse-cdn.com/v4/letter/g/f14d63/32.png) [@Gon\_F](https://discourse.pymc.io/u/Gon_F)\
**Post date:** [June 3, 2019, 6:57pm UTC](https://discourse.pymc.io/t/2d-gaussian-process-indexerror/3244/10 "2019-06-03T18:57:01Z")

</div>

So, using the following code on Pymc3 3.7, everything works out without any errors if I run the code on a subset of the first 400 points in the dataset.

Apologies for taking so long to reply: I kept trying to run the whole model but even on my fast computer, it kept taking 20+ hrs., and my python session would quit before that based on some settings I could not change.

If the code below also works out for you, I imagine the sheer size of your dataset stretches the GP methods to their limits, and whatever breaks first is the error you end up getting.

Data import + pre-preprocessing:

```
data_discourse = pd.read_csv('data_discourse.csv')

np.shape(data_discourse)

data_discourse.head(10)

data_discourse.dropna(inplace=True)

data_discourse.reset_index(inplace=True, drop=True)

np.shape(data_discourse)

data_discourse.head(10)

X_ = data_discourse.iloc[:, :2].values
Y_ = data_discourse.iloc[:, 2].values

X_[:5]
Y_[:5]

```

Data preprocessing:

```
import theano as th
import numpy as np

#%%
X_shared=th.shared(X_)
Y_shared=th.shared(Y_)

nd = 30
z1, z2 = np.meshgrid(np.linspace(X_[:,0].min(), X_[:,0].max(), nd), np.linspace(X_[:,1].min(), X_[:,1].max(), nd))
Z = np.concatenate([z1.reshape(nd*nd, 1), z2.reshape(nd*nd, 1)], 1)

nd = 15
xu1, xu2 = np.meshgrid(np.linspace(X_[:,0].min(), X_[:,0].max(), nd), np.linspace(X_[:,1].min(), X_[:,1].max(), nd))
Xu = np.concatenate([xu1.reshape(nd*nd, 1), xu2.reshape(nd*nd, 1)], 1)

```

Model:

```
with pm.Model() as spatial_model:
    
    l = pm.HalfCauchy("l", beta=3, shape=(2,))
    sf2 = pm.HalfCauchy("sf2", beta=3)
    sn2 = pm.HalfCauchy("sn2", beta=3)

# K = pm.gp.cov.ExpQuad(2, l) * sf2**2
    K = pm.gp.cov.Matern32(2,l) * sf2**2
    
    gp_spatial = pm.gp.MarginalSparse(cov_func=K, approx="FITC")
    obs = gp_spatial.marginal_likelihood("obs", X=X_shared[:200], Xu=Xu, y=Y_shared[:200], noise=sn2)

    trace = pm.sample(cores=2, chains=2)

with spatial_model:
    f_pred = gp_spatial.conditional('f_pred2', Z)
    samples = pm.sample_posterior_predictive([trace], vars=[f_pred], samples=100)
```

---

<div class="post-metadata">

**Author:** ![madarshahian](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/madarshahian/32/722_2.png) [@madarshahian](https://discourse.pymc.io/u/madarshahian)\
**Post date:** [June 4, 2019, 3:42pm UTC](https://discourse.pymc.io/t/2d-gaussian-process-indexerror/3244/11 "2019-06-04T15:42:52Z")

</div>

Thanks so much. Actually it worked for me after updating pymc3. One of the things I need to do is to manually write a code to generate posterior predictive for this model. Since at the end of the day I will pass this to people who don’t have pymc3. For simple models like Multiple regression, this is easy for me to construct, but for this model is more confusing. Specially I see instead of changing x\_shared, we used “`gp_spatial.conditional('f_pred2' , Z)`”. I am wondering if you can advice me how to write posterior predictive using my posterior for this model? Is there any Layman reference? ~Thanks!

---

<div class="post-metadata">

**Author:** ![madarshahian](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/madarshahian/32/722_2.png) [@madarshahian](https://discourse.pymc.io/u/madarshahian)\
**Post date:** [June 7, 2019, 1:50pm UTC](https://discourse.pymc.io/t/2d-gaussian-process-indexerror/3244/12 "2019-06-07T13:50:05Z")

</div>

To continue the discussion above, the parameters are `sn2, sf2, l`. Can I claim “`sn2`” is something like standard deviation of residuals? and how can I make posterior predictive (without using pymc3) if I have posterior samples for `sn2, sf2, l` from pymc3? Is there any easy read reference?

---

<div class="post-metadata">

**Author:** ![Gon\_F](https://avatars.discourse-cdn.com/v4/letter/g/f14d63/32.png) [@Gon\_F](https://discourse.pymc.io/u/Gon_F)\
**Post date:** [June 7, 2019, 4:55pm UTC](https://discourse.pymc.io/t/2d-gaussian-process-indexerror/3244/13 "2019-06-07T16:55:14Z")

</div>

I am confused by what you mean just now, but if you want to write A) an analytical expression for your posterior, that would be difficult, and I would hopefully direct you to Rasmussen’s [GP book](http://www.gaussianprocess.org/gpml/chapters/RW.pdf), or B) create predictions without Pymc3, after gaining sampled chains from Pymc3, that would be recreating a complex component manually, which imho would not be easy, and I have no good reference for that.

The best I would hope is maybe one of the developers can chime in on what posterior predicting based on new data precisely does in Pymc3.

---

<div class="post-metadata">

**Author:** ![madarshahian](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/madarshahian/32/722_2.png) [@madarshahian](https://discourse.pymc.io/u/madarshahian)\
**Post date:** [June 14, 2019, 8:56pm UTC](https://discourse.pymc.io/t/2d-gaussian-process-indexerror/3244/14 "2019-06-14T20:56:52Z")

</div>

> [@Gon\_F](#):
>
> imho

Thanks! You got it. So I want to use chains I have for `sn2` and `sf2` and` l` to generate a posterior predictive if someone ask me what is the function value for new x1 and x2. However, I do not want to use pymc3 anymore, I just want to get chains for `sn2` and `sf2` and` l` from pymc3 and then write a code to use them to generate posterior predictive. This is easy for regression models, but I am not sure how it gonna be exactly for Gaussian process.