#7.7 Radioactive Counts Example
x<-rbind(
c(0,  57),
c(1, 203),
c(2, 383),
c(3, 525),
c(4, 532),
c(5, 408),
c(6, 273),
c(7, 139),
c(8,  45),
c(9,  27),
c(10, 10),
c(11,  4),
c(12,  0),
c(13,  1),
c(14,  1))

#Find mean to work out expected counts from Poisson distribution
v<-sum(x[,1]*x[,2])/sum(x[,2])
v
probs<-v^(x[,1])*exp(-v)/factorial(x[,1])
probs[length(probs)]<-1-sum(probs[1:(length(probs)-1)])
tab<-cbind(x[,2],sum(x[,2])*probs)
colnames(tab)<-c("Observed","Expected")
rownames(tab)<-0:(length(probs)-1)
tab

#totals & plot
apply(tab,2,sum)
barplot(t(tab),beside=T,legend.text=c("Observed","Expected"))

#DDR and plot
ddr<-sqrt(4*tab[,1]+2)-sqrt(4*tab[,2]+2)
ddr[tab[,1]==0]<-1-sqrt(4*tab[,2]+1)[tab[,1]==0]
plot(x[,1],ddr,xlab="Category",ylab="DDR",type="n")
points(x[,1],ddr,cex=2)

#group so expected>1
pooledtab<-tab
while (any(pooledtab[,2]<=1))
{
n<-dim(pooledtab)[1]
pooledtab<-rbind(pooledtab[1:(n-2),],apply(pooledtab[(n-1):n,],2,sum))
}
pooledtab

#gof stuff - note have to specify how many parameters were estimated
GOF(pooledtab,q=1)

#gof stuff - attribute the estimated parameter to solving for the location component to be zero
GOF(pooledtab,q=1,qcomps=c(1))