如何使用 pymc3 拟合属于实例的方法?

how to fit a method belonging to an instance with pymc3?

我未能将属于 class 实例的方法作为确定性函数与 PyMc3 相匹配。你能告诉我怎么做吗?

为了简单起见,下面用一个简单的例子总结了我的案例。实际上,我的限制是一切都是通过 GUI 进行的,并且像“find_MAP”这样的操作应该在链接到 pyqt 按钮的方法中。

我想在数据点上拟合函数“FunctionIWantToFit”。问题,代码如下:

import numpy as np
import pymc3 as pm3
from scipy.interpolate import interp1d
import theano.tensor as tt
import theano.compile

class cprofile:
    def __init__(self):
        self.observed_x = np.array([0.3,1.4,3.1,5,6.8,9,13.4,17.1])
        self.observations = np.array([6.25,2.75,1.25,1.25,1.5,1.75,1.5,1])
        self.x = np.arange(0,18,0.5)

    @theano.compile.ops.as_op(itypes=[tt.dscalar,tt.dscalar,tt.dscalar],
                              otypes=[tt.dvector])
    def FunctionIWantToFit(self,t,y,z):
        # can be complicated but simple in this example
        # among other things, this FunctionIWantToFit depends on a bunch of 
        # variables and methods that belong to this instance of the class cprofile,
        # so it cannot simply be put outside the class ! (like in the following example)
        val=t+y*self.x+z*self.x**2
        interp_values = interp1d(self.x,val)
        return interp_values(self.observed_x)

    def doMAP(self):
        model = pm3.Model()
        with model:
            t = pm3.Uniform("t",0,5)
            y = pm3.Uniform("y",0,5)
            z = pm3.Uniform("z",0,5)
            MyModel = pm3.Deterministic('MyModel',self.FunctionIWantToFit(t,y,z))
            obs = pm3.Normal('obs',mu=MyModel,sd=0.1,observed=self.observations)
            start = pm3.find_MAP()
            print('start: ',start)

test=cprofile()
test.doMAP()

给出以下错误:

Traceback (most recent call last):

  File "<ipython-input-15-3dfb7aa09f84>", line 1, in <module>
    runfile('/Users/steph/work/profiles/GUI/pymc3/so.py', wdir='/Users/steph/work/profiles/GUI/pymc3')

  File "/Users/steph/anaconda/lib/python3.5/site-packages/spyder/utils/site/sitecustomize.py", line 866, in runfile
    execfile(filename, namespace)

  File "/Users/steph/anaconda/lib/python3.5/site-packages/spyder/utils/site/sitecustomize.py", line 102, in execfile
    exec(compile(f.read(), filename, 'exec'), namespace)

  File "/Users/steph/work/profiles/GUI/pymc3/so.py", line 44, in <module>
    test.doMAP()

  File "/Users/steph/work/profiles/GUI/pymc3/so.py", line 38, in doMAP
    MyModel = pm3.Deterministic('MyModel',self.FunctionIWantToFit(x,y,z))

  File "/Users/steph/anaconda/lib/python3.5/site-packages/theano/gof/op.py", line 668, in __call__
    required = thunk()

  File "/Users/steph/anaconda/lib/python3.5/site-packages/theano/gof/op.py", line 912, in rval
    r = p(n, [x[0] for x in i], o)

  File "/Users/steph/anaconda/lib/python3.5/site-packages/theano/compile/ops.py", line 522, in perform
    outs = self.__fn(*inputs)

TypeError: FunctionIWantToFit() missing 1 required positional argument: 'z'

怎么了?

备注 1:我系统地收到有关“FunctionIWantToFit”最后一个参数的错误消息。这里是“z”,但如果我从签名中删除 z,则错误消息涉及“y”(除变量名称外相同)。如果我在签名中添加第四个变量“w”,则错误消息涉及“w”(除变量名称外相同)。

rk2:看起来我错过了“theano”或“pymc3”中的一些非常基本的东西,因为当我将“FunctionIWantToFit”放在 class 之外时,它起作用了。请参阅以下示例。

class cprofile:
    def __init__(self):
        self.observations = np.array([6.25,2.75,1.25,1.25,1.5,1.75,1.5,1])

    def doMAP(self):
        model = pm3.Model()
        with model:
            t = pm3.Uniform("t",0,5)
            y = pm3.Uniform("y",0,5)
            z = pm3.Uniform("z",0,5)
            MyModel = pm3.Deterministic('MyModel',FunctionIWantToFit(t,y,z))
            obs = pm3.Normal('obs',mu=MyModel,sd=0.1,observed=self.observations)
            start = pm3.find_MAP()
            print('start: ',start)

@theano.compile.ops.as_op(itypes=[tt.dscalar,tt.dscalar,tt.dscalar],
                              otypes=[tt.dvector])
def FunctionIWantToFit(t,y,z):
        observed_x = np.array([0.3,1.4,3.1,5,6.8,9,13.4,17.1])
        x = np.arange(0,18,0.5)
        val=t+y*x+z*x**2
        interp_values = interp1d(x,val)
        return interp_values(observed_x)

test=cprofile()
test.doMAP()

给出:

Warning: gradient not available.(E.g. vars contains discrete variables). MAP estimates may not be accurate for the default parameters. Defaulting to non-gradient minimization fmin_powell.
WARNING:pymc3:Warning: gradient not available.(E.g. vars contains discrete variables). MAP estimates may not be accurate for the default parameters. Defaulting to non-gradient minimization fmin_powell.
Optimization terminated successfully.
         Current function value: 1070.673818
         Iterations: 4
         Function evaluations: 179
start:  {'t_interval_': array(-0.27924150484602733), 'y_interval_': array(-9.940000425802811), 'z_interval_': array(-12.524909223913992)}

除了我不知道如何在不对几个模块进行大的修改的情况下做到这一点,因为真正的“FunctionIWantToFit”依赖于属于 class 配置文件的这个实例的一堆变量和方法.

事实上,我什至不确定我知道该怎么做,因为“FunctionIWantToFit”应该在参数中包含对象(我目前通过 self 使用)而且我不确定该怎么做使用 theano 装饰器。

所以我宁愿避免这种解决方案...除非有必要。然后我需要解释如何实现它...


2017 年 4 月 9 日添加:

即使没有插值问题,它也不起作用,因为我一定是错过了 theano and/or pymc3 中一些明显的东西。请你能解释一下这个问题吗?我只想比较模型和数据。首先,坚持使用 pymc2 真是太可惜了。 ;其次,我确定我不是唯一遇到这种基本问题的人。

例如,让我们考虑围绕这个非常基本的代码的变体:

import numpy as np
import theano
import pymc3
theano.config.compute_test_value = 'ignore'
theano.config.on_unused_input = 'ignore'

class testclass:
    x = np.arange(0,18,0.5)
    observed_x = np.array([0.3,1.4,3.1,5,6.8,9,13.4,17.1])
    observations = np.array([6.25,2.75,1.25,1.25,1.5,1.75,1.5,1])

    def testfunc(self,t,y,z):
        t2 = theano.tensor.dscalar('t2')
        y2 = theano.tensor.dscalar('y2')
        z2 = theano.tensor.dscalar('z2')
        val = t2 + y2 * self.observed_x + z2 * self.observed_x**2
        f = theano.function([t2,y2,z2],val)
        return f

test=testclass()
model = pymc3.Model()
with model:
    t = pymc3.Uniform("t",0,5)
    y = pymc3.Uniform("y",0,5)
    z = pymc3.Uniform("z",0,5)

with model:
   MyModel = pymc3.Deterministic('MyModel',test.testfunc(t,y,z))

with model:
   obs = pymc3.Normal('obs',mu=MyModel,sd=0.1,observed=test.observations)

此代码在最后一行失败并显示错误消息:TypeError: unsupported operand type(s) for -: 'TensorConstant' and 'Function'

如果我把'testfunc'改成:

def testfunc(self,t,y,z):
    t2 = theano.tensor.dscalar('t2')
    y2 = theano.tensor.dscalar('y2')
    z2 = theano.tensor.dscalar('z2')
    val = t2 + y2 * self.observed_x + z2 * self.observed_x**2
    f = theano.function([t2,y2,z2],val)
    fval = f(t,y,z,self.observed_x)
    return fval

代码在 'MyModel =' 行失败,错误为 TypeError: ('Bad input argument to theano function with name "/Users/steph/work/profiles/GUI/pymc3/theanotest170409.py:32" at index 0(0-based)', 'Expected an array-like object, but found a Variable: maybe you are trying to call a function on a (possibly shared) variable instead of a numeric array?')

如果我返回原来的 'testfunc' 但将最后 'with model' 行更改为:

with model:
   fval = test.testfunc(t,y,z)
   obs = pymc3.Normal('obs',mu=fval,sd=0.1,observed=test.observations)

错误同第一个错误

我在这里只展示了 3 次尝试,但我想强调的是,我在几个小时内尝试了很多组合,越来越简单,直到这些组合。与 pymc2 相比,我感觉 pymc3 显示出巨大的精神变化,我没有得到,而且记录很少......

theano.compile.ops.as_op 只是定义简单 Theano Ops 的缩写。如果要编写更多涉及的代码,最好单独定义它class。如果确实有必要,此 class 的对象当然可以引用您的 cprofile 实例。

http://deeplearning.net/software/theano/extending/extending_theano.html

好的,让我们分部分来做。首先,我将解释您收到的错误消息,然后我将告诉您我将如何进行。

关于第一个问题,您抱怨缺少参数的直接原因是因为您在 class 中定义的函数将 (self, t, y, z) 作为输入,而您在 op 装饰器中将其声明为只有三个输入(t、y、z)。您必须在装饰器中将输入声明为四个,以说明 class 实例本身。

在 "added on april 9, 2017:" 上,第一个代码将不起作用,因为 test.testfunc(t,y,z) 的输出本身就是一个 theano 函数。 pymc3.Deterministic 期望它输出 theano 变量(或 python 变量)。相反,直接使 test.testfun 输出 val = t2 + y2 * self.observed_x + z2 * self.observed_x**2。

然后,在 "if I change 'testfunc' into:",由于 pymc3 尝试使用 theano 函数的方式,您会收到该错误。长话短说,问题是当 pymc3 使用这个函数时,它会向它发送 theano 变量,而 fval 期望数值变量(numpy 数组或其他)。和上一段一样,直接输出val即可:不需要为此编译任何theano函数。

至于我将如何进行,我会尝试将 class 实例声明为 theano 装饰器的输入。不幸的是,我找不到任何关于如何执行此操作的文档,而且它实际上可能是不可能的(例如,参见 this old post)。

然后我会尝试将函数需要的所有内容作为输入传递,并在 class 之外定义它。这可能非常麻烦,如果它需要方法作为输入,那么您 运行 会遇到其他问题。

另一种方法是创建 theano.gof.Op 的子 class,其 init 方法采用您的 class(或者更确切地说它的实例)作为输入,然后定义您的 perform() 方法。这看起来像这样:

class myOp(theano.gof.Op):
    """ These are the inputs/outputs you used in your as_op
    decorator.
    """
    itypes=[tt.dscalar,tt.dscalar,tt.dscalar]
    otypes=[tt.dvector]
    def __init__(self, myclass):
        """ myclass would be the class you had from before, which
        you called cprofile in your first block of code."""
        self.myclass = myclass
    def perform(self,node, inputs, outputs):
        """ Here you define your operations, but instead of
        calling everyting from that class with self.methods(), you
        just do self.myclass.methods().

        Here, 'inputs' is a list with the three inputs you declared
        so you need to unpack them. 'outputs' is something similar, so
        the function doesn't actually return anything, but saves all
        to outputs. 'node' is magic juice that keeps the world
        spinning around; you need not do anything with it, but always
        include it.
        """
        t, y, z = inputs[0][0], inputs[0][1], inputs[0][2]
        outputs[0][0] = t+y*self.myclass.x+z*self.myclass.x**2
myop = myOp(myclass)

完成此操作后,您可以将 myop 用作其余代码的操作。请注意,某些部分丢失了。您可以查看 了解更多详情。

至于,不需要定义grad()方法。因此,如果有帮助,您可以在纯 python 中执行 perform() 中的所有操作。

或者,如果您可以访问您正在使用的 class 的定义,那么我说这话的时候,您也可以让它继承自 theano.gof.Op,创建 perform() 方法(,您在其中留言)并尝试像那样使用它。它可能会与您正在做的任何其他事情发生冲突 class 并且可能很难做到正确,但尝试可能会很有趣。

我终于收敛到下面成功的代码:

import numpy as np
import theano
from scipy.interpolate import interp1d
import pymc3 as pm3
theano.config.compute_test_value = 'ignore'
theano.config.on_unused_input = 'ignore'

class cprofile:
    observations = np.array([6.25,2.75,1.25,1.25,1.5,1.75,1.5,1])
    x = np.arange(0,18,0.5)
    observed_x = np.array([0.3,1.4,3.1,5,6.8,9,13.4,17.1])    

    def doMAP(self):
        model = pm3.Model()
        with model:
            t = pm3.Uniform("t",0,5)
            y = pm3.Uniform("y",0,5)
            z = pm3.Uniform("z",0,5)
            obs=pm3.Normal('obs',
              mu=FunctionIWantToFit(self)(t,y,z),
              sd=0.1,observed=self.observations)
            start = pm3.find_MAP()
            print('start: ',start)

class FunctionIWantToFit(theano.gof.Op):
    itypes=[theano.tensor.dscalar,
            theano.tensor.dscalar,
            theano.tensor.dscalar]
    otypes=[theano.tensor.dvector]

    def __init__(self, cp):
        self.cp = cp # note cp is an instance of the 'cprofile' class

    def perform(self,node, inputs, outputs):
        t, y, z = inputs[0], inputs[1], inputs[2]

        xxx = self.cp.x
        temp = t+y*xxx+z*xxx**2
        interpolated_concentration = interp1d(xxx,temp)   
        outputs[0][0] = interpolated_concentration(self.cp.observed_x)

testcp=cprofile()
testcp.doMAP()

感谢达里奥的回答,因为我太慢了,无法自己理解第一个答案。我回顾了一下,但我强烈认为 pymc3 文档非常不清楚。它应该包含非常简单和说明性的示例。

然而,在克里斯的评论之后,我没有成功地做任何有用的事情。谁能解释一下 and/or 举个例子?

还有一点:我不知道我上面的例子是否有效或者是否可以简化。特别是它给我的印象是实例“testcp”在内存中被复制了两次。更多comments/answers欢迎更进一步