我在这里贴了一本python笔记本:http://nbviewer.ipython.org/gist/awellis/9067358
我正在尝试使用pymc3创建一个probit回归模型,使用生成的数据恢复已知参数(参见笔记本)。截距的估计值还可以,但是斜率估计值有点偏差。在
我的模型是这样的:
with pm.Model() as model:
# priors
alpha = pm.Normal('alpha', mu=0, tau=0.001)
beta = pm.Normal('beta', mu=0, tau=0.001)
# linear predictor
theta_p = (alpha + beta * x)
# logic transform (just for comparison - this seems to work ok)
# def invlogit(x):
# import theano.tensor as t
# return t.exp(x) / (1 + t.exp(x))
# theta = invlogit(theta_p)
# Probit transform: this doesn't work
def phi(x):
import theano.tensor as t
return 0.5 * (1 + t.erf(x / t.sqr(2)))
theta = phi(theta_p)
# likelihood
y = pm.Bernoulli('y', p=theta, observed=y)
with model:
# Inference
start = pm.find_MAP() # Find starting value by optimization
print("MAP found:")
print("alpha:", start['alpha'])
print("beta:", start['beta'])
print("Compare with true values:")
print("true_alpha", true_alpha)
print("true_beta", true_beta)
with model:
step = pm.NUTS()
trace = pm.sample(2000,
step,
start=start,
progressbar=True) # draw posterior samples
同样地,从PyX函数中可以看出,用回归函数来定义它是不可能的。在
谁能给我指出正确的方向吗?有更好/更简单的方法吗?在
这可能是很久以前的事了,但我已经尝试了在一个简单的层次二项式模型中实现它,并发现结果可与logit函数相媲美。在
唯一的区别是我使用了tensor sqrt()函数。可能只是你的错别字?在
相关问题 更多 >
编程相关推荐