Skip to content
Prev 19151 / 29559 Next

output data from R to SaTScan

Hi Hugh,

Thanks for your answer. Yes I've seen it said that those packages can do most of what SaTScan does.
I'm playing with Kulldorff's Statistic as implemented in section 11.5.6 of ASDAR - but a few things confuse me:

1. How do I decide what the correct value for fractpop is ? I initially had it set to .25 and I was getting cluster of 50% of my cases which made no sense.
2. Is there any correction for multiple testing in the opgam() command ? I have over 3000 areas - do I need to set a very low alpha ?
3. How do I detect more than one cluster ? I have maps where I suspect 2 or 3 separate clusters based on smoothed RR's.


My code currently looks like this (which I based on the code for Fig11.18 in the book and this useful question and answer: http://r-sig-geo.2731867.n2.nabble.com/Strange-results-in-DCluster-package-td7326832.html )

sa<-data.frame(Observed=ED$var1)
sa<-cbind(sa, Expected=mean(ED$var1))
sa<-cbind(sa, x=ED$x, y=ED$y)
sa$Observed<-as.numeric(sa$Observed)

#Kuldorff-Nagarwalla analysis
mle<-calculate.mle(sa, model="poisson") 
knres<-opgam(data=sa, thegrid=sa[,c("x", "y")], alpha=.025, R=99,
             iscluster=kn.iscluster, fractpop=.10, model="poisson", mle=mle, log.v=TRUE)

#Print cluster result
clusters<-get.knclusters(as(ED, "data.frame"), knres)
i<-which.max(knres$statistic)

ED$KNcluster<-""
ED$KNcluster[clusters[[i]]]<-"cluster"
ED$KNcluster[clusters[[i]][1]]<-"centre"
ED$KNcluster<-as.factor(ED$KNcluster)

print(spplot(ED, "KNcluster", main="Kulldorff's method",
              col.regions=c(gray(1), gray(.5), gray(.8))))


Thanks,
James