# Scipy functions in pymc3

**URL:** <https://discourse.pymc.io/t/scipy-functions-in-pymc3/4554>\
**Category:** General\
**Tags:** theano\
**Created:** [February 25, 2020, 12:19am UTC](https://discourse.pymc.io/t/scipy-functions-in-pymc3/4554 "2020-02-25T00:19:10Z")\
**Posts on this page:** 1\
**Page:** 1

<div class="post-metadata">

**Author:** ![balarsen](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/balarsen/32/816_2.png) [@balarsen](https://discourse.pymc.io/u/balarsen)\
**Post date:** [February 25, 2020, 12:19am UTC](https://discourse.pymc.io/t/scipy-functions-in-pymc3/4554/1 "2020-02-25T00:19:10Z")

</div>

Hello all,

I am trying to reproduce some literature results, this time utilizing the special function “[Confluent hypergeometric function 1F1](https://docs.scipy.org/doc/scipy/reference/generated/scipy.special.hyp1f1.html#scipy.special.hyp1f1)”, OK I didn’t know it either.

The paper is: Coath, Christopher D., Robert CJ Steele, and W. Fred Lunnon. “Statistical bias in isotope ratios.” Journal of Analytical Atomic Spectrometry 28.1 (2013): 52-58.

The issue is a non-biased estimator for expected value of the ratio of two Poisson variables.

r'(X, Y)=\frac{X}{Y+1}\cdot M(1, Y+2, 0) where M is that special function above and assuming no noise.

This seems like it should be straightforward enough but I don’t understand something about Theano maybe. My attempt at a simple model is:

```
import pymc3 as pm
import numpy as np
import scipy
import scipy.special as sc
import theano.tensor as tt
from theano.compile.ops import as_op

@as_op(itypes=[tt.lscalar], otypes=[tt.dscalar])
def hyp1f1(value):
    return sc.hyp1f1(1, value+2, 0)

with pm.Model() as model:
    p1 = pm.Poisson('p1', 10)
    p2 = pm.Poisson('p2', 15)
    ratio = pm.Deterministic('ratio', p1/p2)
    quasi = pm.Deterministic('quasi', ((p1)/p2+1)*hyp1f1(p2))
    step = pm.Metropolis()
    trace = pm.sample(10000, step=step)

```

But try as I may I cannot escape from  
`ValueError: expected an ndarray` or `TypeError: For compute_test_value, one input test value does not have the requested type.`

Anyone have any help/fixes for this? Might be simple.

Thanks!
