I'm still having trouble getting alpha even with the new starting values. I was able to use similar code with proc mcmc, which as I understand it uses sampling rather than optimization. The mcmc code did estimate an alpha, although it is lower than what was in the pdf that I attached to my first email. Maybe you just need to see if the model is fitting the data in a reasonable way. So you fit both Poisson and Negative binomial and choose the distribution that best fits the data. proc mcmc data=koch36 nmc=200000 thin=10 nbi=10000; parms alpha 0.5 intercept 0 b1 0 b2 0 b3 0 b4 0 b5 0 a1 0 ; y=Melanoma; t=population; prior intercept b1 b2 b3 b4 b5 a1 ~normal(mean=0,var=1E6); prior alpha~uniform(0,6); mu=exp(log(t)+Intercept+a1*area+b1*(AgeGroup="35-44")+b2*(AgeGroup="45-54")+b3*(AgeGroup="54-64")+ b4*(AgeGroup="65-74")+b5*(AgeGroup=">74")); llike=lgamma(y+1/alpha)-lgamma(1/alpha)-lgamma(y+1)-log(1+alpha*mu)/alpha-y*log(1+alpha*mu)+y*log ( alpha)+y*log(mu); model y ~ general(llike); run; Parameter N Mean Deviation 95% HPD Interval alpha 20000 0.0563 0.0825 1.218E-7 0.2058 intercept 20000 -10.6468 0.2068 -11.0502 -10.2185 b1 20000 1.8008 0.2615 1.2863 2.3530 b2 20000 1.9040 0.2674 1.3644 2.4317 b3 20000 2.2320 0.2630 1.7101 2.7652 b4 20000 2.3818 0.2652 1.8874 2.9467 b5 20000 2.9064 0.2735 2.3431 3.4402 a1 20000 0.8222 0.1590 0.5171 1.1548
... View more