I have this problem: I have a cohort of individuals grouped in 5 age groups. Initially all of them are susceptible and then they develop disease and finally they have cancers. I have information about the age group distribution of the susceptible and then the cancer carrier. Between the susceptible state and the cancer they pass through 7 stages , with same transition rate. I'm trying to create a model that simulate each transition as a binomial extraction and fit the data I have. I tried something but in the moment of analysing the traces , nothing work You can see the code Where am I getting wrong? Thanks for any help
from pylab import *
from pymc import *
from pymc.Matplot import plot as plt
#susceptible_data = array([647,1814,8838,9949,1920])
susceptible_data = array([130,398,1415,1303,206])
infected_data_100000 = array([0,197,302,776,927])
infected_data = array([0,7,38,90,17])
prior_values=np.zeros(len(infected_data))
for i in range(0,len(infected_data)):
prior_values[i]=infected_data[i]/susceptible_data[i]
# stochastic priors
beta1 = Uniform('beta1', 0., 1.)
lambda_0_temp=susceptible_data[0]
lambda_0_1=pymc.Binomial("lambda_0_1",lambda_0_temp,pow(beta1,1))
lambda_0_2=pymc.Binomial("lambda_0_2",lambda_0_1.value,pow(beta1,1))
lambda_0_3=pymc.Binomial("lambda_0_3",lambda_0_2.value,pow(beta1,1))
lambda_0_4=pymc.Binomial("lambda_0_4",lambda_0_3.value,pow(beta1,1))
lambda_0_5=pymc.Binomial("lambda_0_5",lambda_0_4.value,pow(beta1,1))
lambda_0_6=pymc.Binomial("lambda_0_6",lambda_0_5.value,pow(beta1,1))
lambda_0_7=pymc.Binomial("lambda_0_7",n=lambda_0_6.value,p=pow(beta1,1),value=infected_data[0],observed=True)
lambda_1_temp=susceptible_data[1]
lambda_1_1=pymc.Binomial("lambda_1_1",lambda_1_temp,pow(beta1,1))
lambda_1_2=pymc.Binomial("lambda_1_2",lambda_1_1.value,pow(beta1,1))
lambda_1_3=pymc.Binomial("lambda_1_3",lambda_1_2.value,pow(beta1,1))
lambda_1_4=pymc.Binomial("lambda_1_4",lambda_1_3.value,pow(beta1,1))
lambda_1_5=pymc.Binomial("lambda_1_5",lambda_1_4.value,pow(beta1,1))
lambda_1_6=pymc.Binomial("lambda_1_6",lambda_1_5.value,pow(beta1,1))
lambda_1_7=pymc.Binomial("lambda_1_7",n=lambda_1_6.value,p=pow(beta1,1),value=infected_data[1],observed=True)
lambda_2_temp=susceptible_data[2]
lambda_2_1=pymc.Binomial("lambda_2_1",lambda_2_temp,pow(beta1,1))
lambda_2_2=pymc.Binomial("lambda_2_2",lambda_2_1.value,pow(beta1,1))
lambda_2_3=pymc.Binomial("lambda_2_3",lambda_2_2.value,pow(beta1,1))
lambda_2_4=pymc.Binomial("lambda_2_4",lambda_2_3.value,pow(beta1,1))
lambda_2_5=pymc.Binomial("lambda_2_5",lambda_2_4.value,pow(beta1,1))
lambda_2_6=pymc.Binomial("lambda_2_6",lambda_2_5.value,pow(beta1,1))
lambda_2_7=pymc.Binomial("lambda_2_7",n=lambda_2_6.value,p=pow(beta1,1),value=infected_data[2],observed=True)
lambda_3_temp=susceptible_data[3]
lambda_3_1=pymc.Binomial("lambda_3_1",lambda_3_temp,pow(beta1,1))
lambda_3_2=pymc.Binomial("lambda_3_2",lambda_3_1.value,pow(beta1,1))
lambda_3_3=pymc.Binomial("lambda_3_3",lambda_3_2.value,pow(beta1,1))
lambda_3_4=pymc.Binomial("lambda_3_4",lambda_3_3.value,pow(beta1,1))
lambda_3_5=pymc.Binomial("lambda_3_5",lambda_3_4.value,pow(beta1,1))
lambda_3_6=pymc.Binomial("lambda_3_6",lambda_3_5.value,pow(beta1,1))
lambda_3_7=pymc.Binomial("lambda_3_7",n=lambda_3_6.value,p=pow(beta1,1),value=infected_data[3],observed=True)
lambda_4_temp=susceptible_data[4]
lambda_4_1=pymc.Binomial("lambda_4_1",lambda_4_temp,pow(beta1,1))
lambda_4_2=pymc.Binomial("lambda_4_2",lambda_4_1.value,pow(beta1,1))
lambda_4_3=pymc.Binomial("lambda_4_3",lambda_4_2.value,pow(beta1,1))
lambda_4_4=pymc.Binomial("lambda_4_4",lambda_4_3.value,pow(beta1,1))
lambda_4_5=pymc.Binomial("lambda_4_5",lambda_4_4.value,pow(beta1,1))
lambda_4_6=pymc.Binomial("lambda_4_6",lambda_4_5.value,pow(beta1,1))
lambda_4_7=pymc.Binomial("lambda_4_7",n=lambda_4_6.value,p=pow(beta1,1),value=infected_data[4],observed=True)
model=pymc.Model([lambda_0_7,lambda_1_7,lambda_2_7,lambda_3_7,lambda_4_7,beta1])
mcmc =pymc.MCMC(model)
mcmc.sample(iter=100000, burn=50000, thin=10, verbose=1)
lambda_0_samples=mcmc.trace('lambda_0_7')[:]
lambda_1_samples=mcmc.trace('lambda_1_7')[:]
lambda_2_samples=mcmc.trace('lambda_2_7')[:]
lambda_3_samples=mcmc.trace('lambda_3_7')[:]
lambda_4_samples=mcmc.trace('lambda_4_7')[:]
beta1_samples=mcmc.trace('beta1')[:]