1

我想使用 python 传播不确定性。这对于通过不确定性包的简单功能来说相对容易。然而，用用户定义的函数来实现同样的效果并不是那么明显。以下是我正在尝试做的一个例子。

import mcerp as err
import numpy as np


def mult_func(x,xm ,a):
    x[x==0.] = 1e-20
    v = (1.-(xm/x)**a) * (x > xm)
    v[np.isnan(v)] = 0.
    return v


def intg(e,f,cut,s):
    t = mult_func(e,cut,s)
    res = np.trapz(t*f,e)  
    return res  

x=np.linspace(0,1,10000)
y=np.exp(x)

m=0.
mm=0.
N=100000  
for i in range(0,N):
    cut=np.random.normal(0.21,0.02)
    stg=np.random.normal(1.1,0.1)
    v=intg(x,y,cut,stg)
    m=m+v
    mm=mm+v*v



print("avg. %10.5E +/- %10.5E fixed %10.5E"%(m/N,np.sqrt((mm/N-(m/N)**2)),intg(x,y,0.21,1.1)))

上面所做的只是对两个参数进行随机抽样并计算均值和方差。但是，我不确定这种蛮力方法在多大程度上是足够的。我可以使用大数定律并尝试估计需要多少次试验N才能使某个值(P=1-1/(N*k**2))大约k是真实均值的标准差的倍数。

原则上我写的可以工作。然而，我的假设是，作为一种具有许多强大包的灵活语言，python 可以更有效地完成这项任务。我在想uncertainties，mcerp和pymc。由于我使用这些软件包的经验有限，我不确定如何继续。

编辑：我原来的例子没有那么多信息，这就是为什么我决定做一个新的例子来说明我的想法。

4

1 回答 1

-1

Numpy 支持任意数字类型的数组。但是，并非所有函数都支持任意数字类型。

在这种情况下，不支持numpy.exp和。trapz

请注意，不确定性模块包含unumpy包。numpy.exp这里有一个替代品：uncertainties.unumpy.exp

我们将 trapz 定义为ufunc。

在这里查看！

a=un.ufloat(0.3,0.01)
b=un.ufloat(1.2,0.071)



def sample_func(a: un.UFloat, b: un.UFloat) -> np.ndarray:
    x=np.linspace(0,a,100)
    y = un.unumpy.exp(x)
    return utrapz(y, x)

def utrapz(y: np.ndarray, x: np.ndarray) -> np.ndarray:
    Δx = x[1:]-x[:-1]
    avg_y = (y[1:]+y[:-1])/2
    return (Δx*avg_y)


print(sample_func(a, b))

出去：

[0.00026601240063021264+/-nan 0.0005935120815465686+/-6.429403852670308e-06
 0.0006973604419223405+/-3.888235103342809e-06 ...,
 0.002095505706899622+/-6.503985178118233e-05
 0.0021019968633076134+/-6.545802781649068e-05
 0.0021084415802710295+/-6.587387316821736e-05]
于 2017-06-18T18:54:54.990 回答