symulacje

##########################################################################################
###                                                                  Zastosowanie algorytmu Gibbsa                                                               ###
##########################################################################################

N=1000
lambda=matrix(0,N,10)
czas=c()
beta=c()
dane=c(5,1,5,14,3,19,1,1,4,22)
czas=c(94.32,15.72,62.88,125.76,5.24,31.44,1.05,1.05,2.10,10.48)
alfa=1.8
beta=c(); beta[1]=0
delta=1
gamma=0.01

for (i in 2:N)
{

for (j in 1:10)
lambda[i,j]=rgamma(1,sh=dane[j]+alfa,ra=czas[j]+beta[i-1])
beta[i]=rgamma(1,sh=gamma+10*alfa,ra=delta+sum(lambda[i,]))

}
par(mfrow=c(1,3))
hist(lambda[,1])
hist(lambda[,2])
hist(beta)