# Use of scipy functions in pyMC3

**URL:** <https://discourse.pymc.io/t/use-of-scipy-functions-in-pymc3/1621>\
**Category:** Questions\
**Tags:** theano\
**Created:** [July 31, 2018, 2:57am UTC](https://discourse.pymc.io/t/use-of-scipy-functions-in-pymc3/1621 "2018-07-31T02:57:11Z")\
**Posts on this page:** 4\
**Page:** 1

<div class="post-metadata">

**Author:** ![lesthesia](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/lesthesia/32/1015_2.png) [@lesthesia](https://discourse.pymc.io/u/lesthesia)\
**Post date:** [July 31, 2018, 2:57am UTC](https://discourse.pymc.io/t/use-of-scipy-functions-in-pymc3/1621/1 "2018-07-31T02:57:11Z")

</div>

Hi everyone.  
I’ve struggle to use PyMC3 but got stuck in a problem.  
A likelihood function I want to define involves the exponential integral function (expi) imported from scipy library.  
But it doesn’t work when I defined the model in PyMC3 and the I see an error message as below

TypeError: ufunc ‘expi’ not supported for the input types, and the inputs could not be safely coerced to any supported types according to the casting rule ‘‘safe’’

How can I use some special functions such as exponential integral and Bessel functions in PyMC3?

Please find the attached python script I am tackling with.  
Thank you in advance  
[trialExpInt.py](https://canada1.discourse-cdn.com/flex036/uploads/pymc3/original/2X/1/1022702c7f94c25759e4943f1bae7fca7d55eb60.py) (751 Bytes)

---

<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:** [July 31, 2018, 4:13am UTC](https://discourse.pymc.io/t/use-of-scipy-functions-in-pymc3/1621/2 "2018-07-31T04:13:28Z")

</div>

You need to either find a implementation in theano, or wrap the scipy function in [theano using @as\_op](http://docs.pymc.io/notebooks/getting_started.html?highlight=as_op#Arbitrary-deterministics).

---

<div class="post-metadata">

**Author:** ![lesthesia](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/lesthesia/32/1015_2.png) [@lesthesia](https://discourse.pymc.io/u/lesthesia)\
**Post date:** [July 31, 2018, 5:58am UTC](https://discourse.pymc.io/t/use-of-scipy-functions-in-pymc3/1621/3 "2018-07-31T05:58:47Z")

</div>

Thank you for your answer.  
Actually, I’ve already tried the decorator @as\_op  
However, It showed me an error message as below.

AttributeError: ‘int’ object has no attribute ‘type’

Is here something I am missing?

```
# In[1]:

import os
import numpy as np
import matplotlib.pyplot as plt
import pymc3 as pm
from scipy.special import expi
get_ipython().magic('matplotlib inline')

# In[2]:

Cs = 2.5E6
rb = 0.0825
Leff = 1.9
Rb = 0.162
q = 50
t = np.array([i*3600 for i in range(51)])
temp_mea = -(q/(4.0*Leff))*expi(-Cs*(rb**2)/(4.0*t*Leff)) + 6

# In[3]:

import theano
import theano.tensor as tt

from theano.compile.ops import as_op

@as_op(itypes=[tt.dscalar, tt.dvector, tt.dscalar, tt.dscalar, tt.dscalar],
       otypes=[tt.dvector])

def ils(q, t, Cs, rb, Leff):
    temp_bore = -(q/(4.0*Leff))*expi(-Cs*(rb**2)/(4.0*t*Leff))
    return temp_bore

# In[4]:

basic_model = pm.Model()

with basic_model:

    Leff = pm.Normal('Leff', mu=1.8, sd=0.5)
    Rb = pm.Normal('Rb', mu=0.16, sd=0.1) 

    temp_cal = ils(q, t, Cs, rb, Leff) + Rb*q

    Y_obs = pm.Normal('Y_obs', mu=temp_cal, sd=0.4, observed=temp_mea)
```

---

<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:** [July 31, 2018, 1:34pm UTC](https://discourse.pymc.io/t/use-of-scipy-functions-in-pymc3/1621/4 "2018-07-31T13:34:41Z")

</div>

The input to the function need to have the same type as you specified in the decorator, something like this should work:

```python
# Make sure the input being shared are float
q = 50.
t = np.array([i*3600 for i in range(51)], dtype=np.float64)
...
q_shared, t_shared, Cs_shared, rb_shared = theano.shared(q), theano.shared(t), theano.shared(Cs), theano.shared(rb)
with basic_model:
    ...
    temp_cal = ils(q_shared, t_shared, Cs_shared, rb_shared, Leff) + Rb*q
    ...

```
