我是 pyMC 的新手，我仍然无法使用 pyMC 构建我的 MCMC 的结构。我想建立一个链，我很困惑如何一起定义我的参数和对数似然函数。我的卡方函数由下式给出：
其中
和
分别是观测数据和对应误差，
是具有四个自由参数的模型，参数是非线性的。
X和的先验Y是统一的，例如：
import pymc as pm
import numpy as np
import math
import random
@pm.stochastic(dtype=np.float, observed=False, trace=True)
def Xpos(value=1900,x_l=1851,x_h=1962):
"""The probable region of the position of halo centre"""
def logp(value,x_l,x_h):
if ((value>x_h) or (value<x_l)):
return -np.inf
else:
return -np.log(x_h-x_l+1)
def random(x_l,x_h):
return np.round((x_h-x_l)*random.random())+x_l
@pm.stochastic(dtype=np.float, observed=False, trace=True)
def Ypos(value=1900,y_l=1851,y_h=1962):
"""The probable region of the position of halo centre"""
def logp(value,y_l,y_h):
if ((value>y_h) or (value<y_l)):
return -np.inf
else:
return -np.log(y_h-y_l+1)
def random(y_l,y_h):
return np.round((y_h-y_l)*random.random())+y_l
但对于M和C给出如下：
其中的平均值C是通过计算的
对于M和C，先验应该如下所示：
M=math.pow(10,15)*pm.Exponential('mass', beta=math.pow(10,15))
@pm.stochastic(dtype=np.float, observed=False, trace=True)
def concentration(value=4, zh, M200):
"""logp for concentration parameter"""
def logp(value=4.,zh, M200):
if (value>0):
x = np.linspace(math.pow(10,13),math.pow(10,16),200 )
prob=expon.pdf(x,loc=0,scale=math.pow(10,15))
conc = [5.26/(1.+zh)*math.pow(x[i]/math.pow(10,14),-0.1) for i in range(len(x))]
mu_c=0
for i in range(len(x)):
mu_c+=prob[i]*conc[i]/sum(prob)
if (M200 < pow(10,15)):
tau=1./(0.09*0.09)
else:
tau=1./(0.06*0.06)
return pm.lognormal_like(value, mu_c, tau)
else
return -np.inf
def random(mu_c,tau):
return np.random.lognormal(mu_c, tau, 1)
该参数z也是C先验常数。我想知道如何定义我的可能性
，并且应该将其称为@Deterministic variable？我是否以正确的方式将M和定义C为先验信息？
如果有人给我一些提示，告诉我如何将这些参数与给定的先验结合起来，我将不胜感激。