Hello,
I am newbie in Julia, trying to run a MCMC algorithm, but I getting some complaints about the parameters of the gamma distribution (which are supposed to be > 0), but they are, so I don’t understand this complaints. My code is below:
using Random, Distributions
# Data
M=52
NN=30
n=26
m=34
x = [49., 5., 17., 2., 39., 84., 7., 0., 35., 36., 1400., 5., 34., 15., 11.,
2., 1., 39., 8., 101., 2., 148., 1., 68., 31., 1., 20., 118., 91.,
427.]
dc = [1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 0., 1., 1., 1., 1., 1., 1., 1.,
1., 1., 1., 1., 1., 1., 1., 1., 1., 0., 0., 0.]
y = [0., 35., 50., 11., 25., 16., 36., 27., 19., 17., 7., 11., 2., 82., 24.,
70., 15., 16., 50., 22., 45., 18., 4., 1., 40., 57., 0., 1., 20.,
35., 82., 31., 40., 9., 66., 20., 77., 2., 26., 32., 13., 56., 2.,
9., 4., 30., 3., 26., 4., 45., 25., 5.]
z = [15., 3., 624., 46., 127., 61., 1350., 312., 24., 10., 1024., 39., 730.,
136., 1379., 1., 836., 60., 1140., 1153., 54., 47., 0., 43., 971.,
868., 44., 780., 51., 710., 663., 253., 147., 51., 479., 322., 442.,
65., 419., 362., 64., 228., 65., 264., 25., 193., 196., 63., 12.,
103., 60., 43.]
dt = [1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 0., 1., 1., 1.,
0., 0., 1., 1., 1., 1., 0., 0., 1., 0., 1., 0., 0., 1., 1., 1., 0.,
1., 0., 1., 0., 0., 1., 1., 1., 0., 1., 0., 0., 1., 1., 0., 0., 0.]
#defining vectors and matrices:
lambda=zeros(Float64,1000)
tau=zeros(Float64,1000)
p=zeros(Float64,1000)
theta=zeros(Float64,(1000, NN))
eta=zeros(Float64,(1000, M))
# setting initial values:
lambda[1]=0.1
tau[1]=0.1
p[1]=3
theta[1,:].=0.1
eta[1,:].=0.1
alpha=0.01
beta=0.01
gama=0.01
delta=0.01
# simulation
for k in 2:1000
lambda[k]=rand(Gamma(2+NN*p[k-1]+M*p[k-1]+alpha,2+beta+sum(theta[k-1,:])+sum(eta[k-1,:])))
tau[k]=rand(Gamma(m+gama, delta+sum(eta[k-1,:].*z)))
p[k]=rand(Gamma(2+lambda[k]*((prod(theta[k-1,:])*prod(eta[k-1,:]))^(1/(M+NN)))^(0.999*(M+NN)),1/gamma(0.999)))
for i in 1:NN
theta[k,i]=rand(Gamma(p[k]+dc[i],lambda[k]+x[i]))
end
for j in 1:M
eta[k,j]=rand(Gamma(p[k]+dt[j],lambda[k]+y[j]+tau[k]*z[j]))
end
end
It complains that some of the gamma distributions within the for have parameters not satisfying the definition, but we compute the initial parameters manually they are greater than zero. I can’t see why it is not working. I would appreciate any help. thanks.