有没有办法加快这个简单的PyMC模型？在20-40个数据点上,需要约5-11秒才能适应.
import pymc
import time
import numpy as np
from collections import OrderedDict
# prior probability of rain
p_rain = 0.5
variables = OrderedDict()
# rain observations
data = [True, True, True, True, True,
False, False, False, False, False]*4
num_steps = len(data)
p_rain_given_rain = 0.9
p_rain_given_norain = 0.2
p_umbrella_given_rain = 0.8
p_umbrella_given_norain = 0.3
for n in range(num_steps):
if n == 0:
# Rain node at time t = 0
rain = pymc.Bernoulli("rain_%d" %(n), p_rain)
else:
rain_trans = \
pymc.Lambda("rain_trans", …这可能是一个愚蠢的问题.
我正在尝试使用PyMC中的MCMC评估将数据拟合到一个非常奇怪的PDF中.对于这个例子,我只想弄清楚如何适应正常分布,我手动输入普通PDF.我的代码是:
data = [];
for count in range(1000): data.append(random.gauss(-200,15));
mean = mc.Uniform('mean', lower=min(data), upper=max(data))
std_dev = mc.Uniform('std_dev', lower=0, upper=50)
# @mc.potential
# def density(x = data, mu = mean, sigma = std_dev):
# return (1./(sigma*np.sqrt(2*np.pi))*np.exp(-((x-mu)**2/(2*sigma**2))))
mc.Normal('process', mu=mean, tau=1./std_dev**2, value=data, observed=True)
model = mc.MCMC([mean,std_dev])
model.sample(iter=5000)
print "!"
print(model.stats()['mean']['mean'])
print(model.stats()['std_dev']['mean'])
我发现的例子都使用了像mc.Normal,或者mc.Poisson或者诸如此类的东西,但是我想要符合注释掉的密度函数.
任何帮助,将不胜感激.
有人可以提供一些关于如何并行化PyMC MCMC代码的一般性说明.我试图LASSO按照这里给出的例子运行回归.我在某处读到默认情况下并行采样,但是我是否还需要使用类似的功能Parallel Python来使其工作？
这是一些我希望能够在我的机器上并行化的参考代码.
x1 = norm.rvs(0, 1, size=n)
x2 = -x1 + norm.rvs(0, 10**-3, size=n)
x3 = norm.rvs(0, 1, size=n)
X = np.column_stack([x1, x2, x3])
y = 10 * x1 + 10 * x2 + 0.1 * x3
beta1_lasso = pymc.Laplace('beta1', mu=0, tau=1.0 / b)
beta2_lasso = pymc.Laplace('beta2', mu=0, tau=1.0 / b)
beta3_lasso = pymc.Laplace('beta3', mu=0, tau=1.0 / b)
@pymc.deterministic
def y_hat_lasso(beta1=beta1_lasso, beta2=beta2_lasso, beta3=beta3_lasso, x1=x1, x2=x2, x3=x3):
return beta1 * x1 …我试图在pymc3中创建一个三级逻辑回归模型.有一个顶级,中级和一个单独的级别,其中中级系数是从顶级系数估算的.但是,我很难为中级指定正确的数据结构.
这是我的代码:
with pm.Model() as model:
# Hyperpriors
top_level_tau = pm.HalfNormal('top_level_tau', sd=100.)
mid_level_tau = pm.HalfNormal('mid_level_tau', sd=100.)
# Priors
top_level = pm.Normal('top_level', mu=0., tau=top_level_tau, shape=k_top)
mid_level = [pm.Normal('mid_level_{}'.format(j),
mu=top_level[mid_to_top_idx[j]],
tau=mid_level_tau)
for j in range(k_mid)]
intercept = pm.Normal('intercept', mu=0., sd=100.)
# Model prediction
yhat = pm.invlogit(mid_level[mid_to_bot_idx] + intercept)
# Likelihood
yact = pm.Bernoulli('yact', p=yhat, observed=y)
我收到了错误"only integer arrays with one element can be converted to an index"(第16行),我认为这与mid_level变量是一个列表,而不是一个合适的pymc容器有关.(我也没有在pymc3源代码中看到Container类.)
任何帮助,将不胜感激.
编辑:添加一些模拟数据
y = np.array([0, 1, 0, 1, 0, 0, 0, …CrossValidated上有一个关于如何使用PyMC将两个Normal分布拟合到数据的问题.Cam.Davidson.Pilon的答案是使用伯努利分布将数据分配给两个法线之一:
size = 10
p = Uniform( "p", 0 , 1) #this is the fraction that come from mean1 vs mean2
ber = Bernoulli( "ber", p = p, size = size) # produces 1 with proportion p.
precision = Gamma('precision', alpha=0.1, beta=0.1)
mean1 = Normal( "mean1", 0, 0.001 )
mean2 = Normal( "mean2", 0, 0.001 )
@deterministic
def mean( ber = ber, mean1 = mean1, mean2 = mean2):
return ber*mean1 + (1-ber)*mean2
现在我的问题是:如何用三个法线做到这一点？ …
假设我有几个直方图,每个直方图都在不同的 bin位置(在实轴上)进行计数.例如
def generate_random_histogram():
# Random bin locations between 0 and 100
bin_locations = np.random.rand(10,) * 100
bin_locations.sort()
# Random counts between 0 and 50 on those locations
bin_counts = np.random.randint(50, size=len(bin_locations))
return {'loc': bin_locations, 'count':bin_counts}
# We can assume that the bin size is either pre-defined or that
# the bin edges are on the middle-point between consecutive counts.
hists = [generate_random_histogram() for x in xrange(3)]
如何对这些直方图进行标准化,以便获得 PDF,其中每个PDF的积分在给定范围内(例如0和100)加起来为1？
我们可以假设直方图根据预定义的bin大小计算事件(例如10)
我见过的大多数实现都是基于高斯内核(参见scipy …
我有一个多变量蒙特卡罗隐马尔可夫问题要解决:
x[k] = f(x[k-1]) + B u[k]
y[k] = g(x[k])
哪里:
x[k] the hidden states (Markov dynamics)
y[k] the observed data
u[k] the stochastic driving process
PyMC3已经足够成熟以解决这个问题,还是应该继续使用2.3版？其次,非常感谢PyMC框架中对HM模型的任何引用.谢谢.
- 亨克
我有一个贝叶斯网(DAG)模型,我使用pymc 2.3创建.其中的所有变量都是伯努利随机变量.当我在采样之前调用它上面的MAP.fit()方法时,我得到所有随机变量的以下警告:
value is neither numerical nor array with floating-point dtype. Recommend fitting method fmin (default)
从pymc的github repo中,如果随机变量的基础类型不是浮点数,似乎会打印此警告.对于伯努利RV来说,类型是(并且应该是)bool.
这是否意味着MAP步骤会产生不稳定的结果？
要学习PyMC,我正在尝试做一个简单的隐马尔可夫模型,如下所示:
with pymc3.Model() as hmm:
# Transition "matrix"
a_t = np.ones(num_states)
T = [pymc3.Dirichlet('T{0}'.format(i), a = a_t, shape = num_states) for i in xrange(num_states)]
# Emission "matrix"
a_e = np.ones(num_emissions)
E = [pymc3.Dirichlet('E{0}'.format(i), a = a_e, shape = num_emissions) for i in xrange(num_states)]
# State models
p0 = np.ones(num_states) / num_states
# No shape, so each state is a scalar tensor
states = [pymc3.Categorical('s0', p = p0)]
emissions = [pymc3.Categorical('z0',
p = ifelse(eq(states[0], 0), E[0], ifelse(eq(states[0], 1), E[1], E[2])),
observed …当我尝试通过anaconda环境安装模块'pymc'时,它显示错误消息如下:
不满意的错误:发现以下规范存在冲突:
大火 - > pyyaml - > python [version ='> = 2.7,<2.8.0a0'] - > vc = 9
大火 - > pyyaml - > yaml - >*[track_features = vc9]
pymc使用"conda info"查看每个包的依赖项.
我使用的是Python 2.7.14,我在Windows上安装了anaconda 1.6.9.我是Python的新手.我首先尝试使用cmd来安装模块pymc,我遇到了很多问题,比如在windows上安装g77编译器的要求.从MinGW获得编译器并安装了Microsoft Visual C++编译器Python之后,我仍然无法安装模块,因为出现了新的错误.那是当我发现在anaconda环境中列出的pymc模块,我可以手动添加,但它显示此冲突错误.
我不知道冲突是否来自我上面安装的所有其他东西.请帮忙!谢谢!