python - 如何舍入并与模型中的 PyMC3 先前参数进行比较?
问题描述
如果它看起来比需要的复杂,我提前道歉。试图使它比这更小,但后来有些问题更难以表述。这样就尽可能接近我的情况。我有一些数据,我想用类似于下面描述的函数来拟合。但是在实现 pyMC3 fit 时,如在第二个代码部分中,我遇到了与这两个问题相关的错误。
# some constants and fixed things
rin = 1.
rout = 1000.
tmin = 10.
x = np.arange(0.1, 10, 0.2)
# parameters to be fitted
a = 4.0
b = 0.5
c = 10
# define r_struct, depends on b, but is one value
r_struct = 0.1**(1./b)
# nrad below, the number of radial points,
# depends on r_struct and the constant rout
# QUESTION 1 here
nrad = round(r_struct, -1) + round((rout - r_struct)/10., -1) * 2
# define the radial bins, logarithmically spaced
rbins = np.logspace(np.log10(rin), np.log10(rout), num=nrad+1, endpoint=True)
rad = 10**(( np.log10(rbins[1:]) + np.log10(rbins[:-1]))/2.)
t = rad**(-b)
# QUESTION 2 here
# t cannot go below tmin, so I do this
t[t < tmin] = tmin
# then I want to know where the limit of r_struct is here
iout = rad > r_struct
iin = rad <= r_struct
# depending on if rad is past r_struct or not
# I do two different things.
c1 = c_struct * (rad[iin]/50.)**(-1.*a)
taper = np.exp( -(rad[iout]/rad[iout][0])**(2.-a) )
c2 = c_taper * taper/taper[0]
tau = np.append(c1, c2)
y_true = J0( x[:, np.newaxis] * rad) * (t * (1.0 - np.exp(-1.*tau) ) )
y_true = y_true.sum(axis=1)
y_obs = y_true + np.random.normal(size=len(x))*50
y_error = np.ones_like(y_true)*0.1*y_true.max()
plt.errorbar(x,y_obs, yerr= y_error, color='r', marker='.', ls='None', label='Observed')
plt.plot(x,y_true, 'b', marker='None', ls='-', lw=1, label='True')
plt.xlabel('x'); plt.ylabel('y'); plt.legend()
好的,所以现在我想把它放在 PyMC3 中。我从以下开始。
with pm.Model() as basic_model:
# parameters
a = pm.Uniform("a", 0, 10)
b = pm.Uniform("b", 0, 10)
c = pm.Uniform("c", 0,100)
r_struct = 0.1**(1./b)
nrad = round(r_struct, -1) + round((rout - r_struct)/10., -1) * 2
rbins = np.logspace(np.log10(rin), np.log10(rout), num=nrad+1, endpoint=True)
rad = 10**(( np.log10(rbins[1:]) + np.log10(rbins[:-1]))/2.)
t = rad**(-b)
t[t < tmin] = tmin
iout = rad > r_struct
iin = rad <= r_struct
c1 = c_struct * (rad[iin]/50.)**(-1.*a)
taper = np.exp( -(rad[iout]/rad[iout][0])**(2.-a) )
c2 = c_taper * taper/taper[0]
tau = np.append(c1, c2)
y_true = J0( x[:, np.newaxis] * rad) * (t * (1.0 - np.exp(-1.*tau) ) )
y_true = y_true.sum(axis=1)
func=pm.Deterministic('gauss',y_true)
y = pm.Normal('y', mu=func, tau=1.0/y_error**2, observed=y_obs)
trace = pm.sample(5000)
然而,这已经在第一次“round()”调用中失败了。请参阅下面的问题。
问题 1. “nrad”必须是整数,它取决于 pymc3 之前的“b”(通过 r_struct),运行 PyMC3 时如何将其舍入为整数? 作为参考,我收到此错误:
---------------------------------------------------------------------------
TypeError Traceback (most recent call last)
<ipython-input-11-fcf897a21b89> in <module>()
15 r_struct = 0.1**(1./b)
16
---> 17 nrad = round(r_struct, -1) + round((rout - r_struct)/10., -1) * 2
18 rbins = np.logspace(np.log10(rin), np.log10(rout), num=nrad+1, endpoint=True)
TypeError: type TensorVariable doesn't define __round__ method
问题 2: 如何与使用 PyMC3 先前创建的数组进行比较?代码在比较中失败。
---------------------------------------------------------------------------
AssertionError Traceback (most recent call last)
<ipython-input-12-7f77e70d147a> in <module>()
26 t = rad**(-b)
27
---> 28 t[t < tmin] = tmin
29 iout = rad > r_struct
30 iin = rad <= r_struct
欢迎任何关于如何澄清问题的建议,包括标题!
解决方案
推荐阅读
- javascript - 如何从 C# 绑定到 HTML?
- swift - @Published 变量上的 @Binding 更改 - 但不会更新视图,除非我出去并返回屏幕
- mysql - 加入/合并具有相同列名的多个表
- javascript - 数组中属于同一字段的值的累积和
- python - 如何从python中的嵌套for循环创建数据框?
- c# - 字符串连接元素的顺序
- r - 如何解决 RStudio 中的 TIFFOpen 错误弹出窗口?
- python - 如何为 django rest 框架 DecimalField 序列化器提供选择?
- react-native - ScrollView 不显示反应原生的所有内容
- javascript - 如何在删除前在 toastr 中要求确认