1opgam.iscluster.default<-function(data, idx, idxorder, alpha, ...) 2{ 3 #Indexes of regions in the balls ordered by distance to the centre 4 localidx<-idxorder[idx[idxorder]] 5 6 localO<-sum(data$Observed[localidx]) 7 localE<-sum(data$Expected[localidx]) 8 if(localE ==0) return(c(localO, NA, NA, NA)) 9 10 11 12 localP<-sum(data$Population[localidx]) 13 14 pvalue<-ppois(localO, localE, lower.tail=FALSE) 15 16 return (c(localO, alpha>pvalue, pvalue, sum(idx))) 17} 18 19opgam.iscluster.negbin<-function(data, idx, idxorder, alpha, mle, R=999, ...) 20{ 21 #Indexes of regions in the balls ordered by distance to the centre 22 localidx<-idxorder[idx[idxorder]] 23 24 25 localO<-sum(data$Observed[localidx]) 26 localE<-sum(data$Expected[localidx]) 27 if(localE ==0) return(c(localO, NA, NA, NA)) 28 29# bt<-boot(data[localidx, ], statistic=function(x){sum(x$Observed)}, 30# R=R, sim="parametric",ran.gen=negbin.sim, 31# mle=list(n=sum(idx),size=mle$size,prob=mle$prob[localidx])) 32# pvalue<-sum(bt$t>localO)/(R+1) 33 34 pvalue<-.Call("Ropgam_iscluster_negbin",data$Observed[localidx], data$Expected[localidx], mle$size, mle$prob[localidx], R, PACKAGE="DCluster") 35 36 return (c(localO, alpha>pvalue, pvalue, sum(idx))) 37} 38 39opgam<-function(data, thegrid=NULL, radius=Inf, step=NULL, alpha, iscluster=opgam.iscluster.default, set.idxorder=TRUE, ...) 40{ 41 #If the Grid is null, then create a new grid 42 if(is.null(thegrid)) 43 { 44 if(is.null(step)) 45 step<-.2*radius 46 47 xgrid<-seq(min(data$x), max(data$x), by=step) 48 ygrid<-seq(min(data$y), max(data$y), by=step) 49 50 xlen<-length(xgrid) 51 ylen<-length(ygrid) 52 npoints<-xlen*ylen 53 54 thegrid<-matrix(rep(NA, 2*npoints) , ncol=2) 55 56 thegrid[,1]<-rep(xgrid, times=ylen) 57 thegrid[,2]<-rep(ygrid, each=xlen) 58 } 59 60 61 rr<-radius*radius 62 63 GAM<-apply(thegrid, 1, opgam.intern, data, rr, set.idxorder, iscluster, alpha, ...) 64 65 #Take only those balls which were significant 66 GAM<-GAM[,!is.na(GAM[4,])] 67 GAM<-as.data.frame(t(GAM[, as.logical(GAM[4,])==TRUE ])) 68 #Just the first five names are set because it is possible 69 #to get more than five columns if user-defined functions 70 #are used 71 names(GAM)<-c("x", "y", "statistic", "cluster", "pvalue", "size") 72 73 return(GAM) 74} 75 76 77opgam.intern<-function(point, data, rr, set.idxorder, iscluster, alpha, ...) 78{ 79 xd<-(data$x-point[1]) 80 yd<-(data$y-point[2]) 81 dist<-xd*xd+yd*yd 82 83 idx<-(dist<=rr) 84 if(set.idxorder) idxorder<-order(dist) 85 86 cl<-iscluster(data=data, idx=idx, idxorder=idxorder, alpha=alpha, ...) 87 88 return(c(point, cl)) 89} 90