R : Copyright 2005, The R Foundation for Statistical Computing
Version 2.1.1  (2005-06-20), ISBN 3-900051-07-0

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for a HTML browser interface to help.
Type 'q()' to quit R.

> ### * <HEADER>
> ###
> attach(NULL, name = "CheckExEnv")
> assign(".CheckExEnv", as.environment(2), pos = length(search())) # base
> ## add some hooks to label plot pages for base and grid graphics
> setHook("plot.new", ".newplot.hook")
> setHook("persp", ".newplot.hook")
> setHook("grid.newpage", ".gridplot.hook")
> 
> assign("cleanEx",
+        function(env = .GlobalEnv) {
+ 	   rm(list = ls(envir = env, all.names = TRUE), envir = env)
+            RNGkind("default", "default")
+ 	   set.seed(1)
+    	   options(warn = 1)
+ 	   delayedAssign("T", stop("T used instead of TRUE"),
+ 		  assign.env = .CheckExEnv)
+ 	   delayedAssign("F", stop("F used instead of FALSE"),
+ 		  assign.env = .CheckExEnv)
+ 	   sch <- search()
+ 	   newitems <- sch[! sch %in% .oldSearch]
+ 	   for(item in rev(newitems))
+                eval(substitute(detach(item), list(item=item)))
+ 	   missitems <- .oldSearch[! .oldSearch %in% sch]
+ 	   if(length(missitems))
+ 	       warning("items ", paste(missitems, collapse=", "),
+ 		       " have been removed from the search path")
+        },
+        env = .CheckExEnv)
> assign("..nameEx", "__{must remake R-ex/*.R}__", env = .CheckExEnv) # for now
> assign("ptime", proc.time(), env = .CheckExEnv)
> grDevices::postscript("DCluster-Examples.ps")
> assign("par.postscript", graphics::par(no.readonly = TRUE), env = .CheckExEnv)
> options(contrasts = c(unordered = "contr.treatment", ordered = "contr.poly"))
> options(warn = 1)    
> library('DCluster')
Loading required package: boot
Loading required package: spdep
Loading required package: tripack
Loading required package: maptools
Loading required package: foreign
Loading required package: SparseM
[1] "SparseM library loaded"
> 
> assign(".oldSearch", search(), env = .CheckExEnv)
> assign(".oldNS", loadedNamespaces(), env = .CheckExEnv)
> cleanEx(); ..nameEx <- "achisq.boot"
> 
> ### * achisq.boot
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: achisq.boot
> ### Title: Bootstrap replicates of Pearson's Chi-square statistic
> ### Aliases: achisq.boot achisq.pboot
> ### Keywords: htest
> 
> ### ** Examples
> 
> library(boot)
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> 
> niter<-100
> 
> #Permutation  model
> chq.perboot<-boot(sids, statistic=achisq.boot, R=niter)
> plot(chq.perboot)#Display results
> 
> #Multinomial model
> chq.mboot<-boot(sids, statistic=achisq.pboot, sim="parametric", ran.gen=multinom.sim,  R=niter)
> plot(chq.mboot)#Display results
> 
> #Poisson model
> chq.pboot<-boot(sids, statistic=achisq.pboot, sim="parametric", ran.gen=poisson.sim,  R=niter)
> plot(chq.pboot)#Display results
> 
> #Poisson-Gamma model
> chq.pgboot<-boot(sids, statistic=achisq.pboot, sim="parametric", ran.gen=negbin.sim, R=niter)
> plot(chq.pgboot)#Display results
> 
> 
> 
> cleanEx(); ..nameEx <- "achisq.stat"
> 
> ### * achisq.stat
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: achisq.stat
> ### Title: Another implementation of Pearson's Chi-square statistic
> ### Aliases: achisq.stat
> ### Keywords: htest
> 
> ### ** Examples
> 
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> 
> #Compute the statistic under the assumption that lambda = 1.
> achisq.stat(sids, lambda=1)
$T
[1] 225.5723

$df
[1] 100

$pvalue
[1] 1.087968e-11

> 
> 
> 
> cleanEx(); ..nameEx <- "besagnewell"
> 
> ### * besagnewell
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: besagnewell
> ### Title: Besag and Newell's statistic for spatial clustering
> ### Aliases: besagnewell
> ### Keywords: spatial
> 
> ### ** Examples
> 
> #B&N must use the centroids as grid.
> #The size of teh cluster is 20.
> #100  bootstrap simulations are performed
> #Poisson is the model used in the bootstrap simulations to generate the
> #observations.
> #Signifiance level is 0'05, even though multiple tests are made.
> 
> library(boot)
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> sids<-cbind(sids, x=nc.sids$x, y=nc.sids$y)
> 
> bnresults<-opgam(sids, thegrid=sids[,c("x","y")], alpha=.05, 
+         iscluster=bn.iscluster, set.idxorder=TRUE, k=20, model="poisson", 
+         R=100, mle=calculate.mle(sids) )
> 
> #Plot all the centroids
> plot(sids$x, sids$y)
> 
> #Plot signifiant centroids in red
> points(bnresults$x, bnresults$y, col="red", pch=19)
> 
> 
> 
> cleanEx(); ..nameEx <- "besagnewell.boot"
> 
> ### * besagnewell.boot
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: besagnewell.boot
> ### Title: Generate boostrap replicates of Besag and Newell's statistic
> ### Aliases: besagnewell.boot besagnewell.pboot
> ### Keywords: spatial
> 
> ### ** Examples
> 
> library(boot)
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> sids<-cbind(sids, x=nc.sids$x, y=nc.sids$y)
> 
> niter<-100
> 
> #Permutation  model
> besn.perboot<-boot(sids, statistic=besagnewell.boot, R=niter, k=20)
> plot(besn.perboot)#Display results
> 
> 
> 
> cleanEx(); ..nameEx <- "besagnewell.stat"
> 
> ### * besagnewell.stat
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: besagnewell.stat
> ### Title: Besag and Newell's statistic for spatial clustering
> ### Aliases: besagnewell.stat
> ### Keywords: spatial
> 
> ### ** Examples
> 
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> sids<-cbind(sids, x=nc.sids$x, y=nc.sids$y)
> 
> besagnewell.stat(sids, k=20)
value  size 
    6     6 
> 
> 
> 
> cleanEx(); ..nameEx <- "bn.iscluster"
> 
> ### * bn.iscluster
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: bn.iscluster
> ### Title: Clustering function for Besag and Newell's method
> ### Aliases: bn.iscluster
> ### Keywords: spatial
> 
> ### ** Examples
> 
> library(boot)
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> sids<-cbind(sids, x=nc.sids$x, y=nc.sids$y)
> 
> #B&N's method
> bnresults<-opgam(data=sids, thegrid=sids[,c("x","y")], alpha=.05, 
+         iscluster=bn.iscluster, k=20, R=100, model="poisson", 
+         mle=calculate.mle(sids))
> 
> #Plot all centroids and significant ones in red
> plot(sids$x, sids$y, main="Besag & Newell's method")
> points(bnresults$x, bnresults$y, col="red", pch=19)
> 
> 
> 
> cleanEx(); ..nameEx <- "calculate.mle"
> 
> ### * calculate.mle
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: calculate.mle
> ### Title: Calculate parameters involved in smapling procedures
> ### Aliases: calculate.mle
> ### Keywords: distribution
> 
> ### ** Examples
> 
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> sids<-cbind(sids, x=nc.sids$x, y=nc.sids$y)
> 
> #Carry out simulations
> datasim<-multinom.sim(sids, mle=calculate.mle(sids, model="multinomal") )
> 
> #Estimators for Poisson distribution
> datasim<-poisson.sim(sids, mle=calculate.mle(sids, model="poisson") )
> 
> #Estimators for Negative Binomial distribution
> datasim<-negbin.sim(sids, mle=calculate.mle(sids, model="negbin") )
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "empbaysmooth"
> 
> ### * empbaysmooth
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: empbaysmooth
> ### Title: Empirical Bayes Smoothing
> ### Aliases: empbaysmooth
> ### Keywords: models
> 
> ### ** Examples
> 
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> 
> smth<-empbaysmooth(sids$Observed, sids$Expected)
> 
> 
> 
> cleanEx(); ..nameEx <- "gearyc.boot"
> 
> ### * gearyc.boot
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: gearyc.boot
> ### Title: Generate bootstrap replicates of Moran's I autocorrelation
> ###   statistic
> ### Aliases: gearyc.boot gearyc.pboot
> ### Keywords: spatial
> 
> ### ** Examples
> 
> library(boot)
> library(spdep)
> 
> data(nc.sids)
> col.W <- nb2listw(ncCR85.nb, zero.policy=TRUE)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> 
> niter<-100
> 
> #Permutation model
> gc.perboot<-boot(sids, statistic=gearyc.boot, R=niter, listw=col.W,
+         n=length(ncCR85.nb), n1=length(ncCR85.nb)-1, S0=Szero(col.W) )
> plot(gc.perboot)#Display results
> 
> #Multinomial model
> gc.mboot<-boot(sids, statistic=gearyc.pboot, sim="parametric", 
+         ran.gen=multinom.sim, R=niter, listw=col.W,
+         n=length(ncCR85.nb), n1=length(ncCR85.nb)-1, S0=Szero(col.W) )
> plot(gc.mboot)#Display results
> 
> #Poisson model
> gc.pboot<-boot(sids, statistic=gearyc.pboot, sim="parametric", 
+         ran.gen=poisson.sim, R=niter, listw=col.W,
+         n=length(ncCR85.nb), n1=length(ncCR85.nb)-1, S0=Szero(col.W) )
> plot(gc.pboot)#Display results
> 
> #Poisson-Gamma model
> gc.pgboot<-boot(sids, statistic=gearyc.pboot, sim="parametric", 
+         ran.gen=negbin.sim, R=niter, listw=col.W,
+         n=length(ncCR85.nb), n1=length(ncCR85.nb)-1, S0=Szero(col.W) )
> plot(gc.pgboot)#Display results
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "gearyc.stat"
> 
> ### * gearyc.stat
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: gearyc.stat
> ### Title: Compute Moran's I autocorrelation statistic
> ### Aliases: gearyc.stat
> ### Keywords: spatial
> 
> ### ** Examples
> 
> library(spdep)
> data(nc.sids)
> col.W <- nb2listw(ncCR85.nb, zero.policy=TRUE)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> 
> gearyc.stat(data=sids, listw=col.W, n=length(ncCR85.nb), n1=length(ncCR85.nb)-1,
+         S0=Szero(col.W) )
[1] 0.8282572
> 
> gearyc.stat(data=sids, applyto="SMR", listw=col.W, n=length(ncCR85.nb), 
+         n1=length(ncCR85.nb)-1,S0=Szero(col.W) )
[1] 0.7110422
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "kn.iscluster"
> 
> ### * kn.iscluster
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: kn.iscluster
> ### Title: Clustering function for Kulldorff and Nagarwalla's statistic
> ### Aliases: kn.iscluster
> ### Keywords: spatial
> 
> ### ** Examples
> 
> library(boot)
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> sids<-cbind(sids, Population=nc.sids$BIR74, x=nc.sids$x, y=nc.sids$y)
> 
> #K&N's method over the centroids
> mle<-calculate.mle(sids, model="poisson")
> knresults<-opgam(data=sids, thegrid=sids[,c("x","y")], alpha=.05, 
+         iscluster=kn.iscluster, fractpop=.5, R=100, model="poisson", mle=mle)
> 
> #Plot all centroids and significant ones in red
> plot(sids$x, sids$y, main="Kulldorff and Nagarwalla's method")
> points(knresults$x, knresults$y, col="red", pch=19)
> 
> 
> 
> cleanEx(); ..nameEx <- "kullnagar.boot"
> 
> ### * kullnagar.boot
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: kullnagar.boot
> ### Title: Generate bootstrap replicates of Kulldorff and Nagarwalla's
> ###   statistic
> ### Aliases: kullnagar.boot kullnagar.pboot
> ### Keywords: spatial
> 
> ### ** Examples
> 
> library(boot)
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> sids<-cbind(sids, Population=nc.sids$BIR74, x=nc.sids$x, y=nc.sids$y)
> 
> niter<-100
> 
> #Permutation  model
> kn.perboot<-boot(sids, statistic=kullnagar.boot, R=niter, fractpop=.2)
> plot(kn.perboot)#Display results
> 
> #Multinomial model
> kn.mboot<-boot(sids, statistic=kullnagar.pboot, sim="parametric", 
+         ran.gen=multinom.sim,  R=niter, fractpop=.2)
> plot(kn.mboot)#Display results
> 
> #Poisson model
> kn.pboot<-boot(sids, statistic=kullnagar.pboot, sim="parametric", 
+         ran.gen=poisson.sim,  R=niter, fractpop=.2)
> plot(kn.pboot)#Display results
> 
> #Poisson-Gamma model
> kn.pgboot<-boot(sids, statistic=kullnagar.pboot, sim="parametric", 
+         ran.gen=negbin.sim, R=niter, fractpop=.2)
> plot(kn.pgboot)#Display results
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "kullnagar.stat"
> 
> ### * kullnagar.stat
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: kullnagar.stat
> ### Title: Kulldorff and Nagarwalla's statistic for spatial clustering.
> ### Aliases: kullnagar.stat kullnagar.stat.poisson kullnagar.stat.bern
> ### Keywords: spatial
> 
> ### ** Examples
> 
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> sids<-cbind(sids, Population=nc.sids$BIR74, x=nc.sids$x, y=nc.sids$y)
> 
> dist<-(sids$x-sids$x[1])^2+(sids$y-sids$y[1])^2
> index<-order(dist)
> #Compute the statistic around the first county
> kullnagar.stat(sids[index,], fractpop=.5)
value  size 
    1     1 
> 
> 
> 
> cleanEx(); ..nameEx <- "moranI.boot"
> 
> ### * moranI.boot
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: moranI.boot
> ### Title: Generate bootstrap replicates of Moran's I autocorrelation
> ###   statistic
> ### Aliases: moranI.boot moranI.pboot
> ### Keywords: spatial
> 
> ### ** Examples
> 
> library(spdep)
> data(nc.sids)
> col.W <- nb2listw(ncCR85.nb, zero.policy=TRUE)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74)) 
> 
> niter<-100
> 
> #Permutation model
> moran.boot<-boot(sids, statistic=moranI.boot, R=niter, listw=col.W, 
+         n=length(ncCR85.nb), S0=Szero(col.W) )
> plot(moran.boot)#Display results
> 
> #Multinomial model
> moran.mboot<-boot(sids, statistic=moranI.pboot, sim="parametric", 
+         ran.gen=multinom.sim,  R=niter, listw=col.W,n=length(ncCR85.nb), 
+         S0=Szero(col.W) )
> plot(moran.mboot)#Display results
> 
> #Poisson model
> moran.pboot<-boot(sids, statistic=moranI.pboot, sim="parametric", 
+         ran.gen=poisson.sim,  R=niter, listw=col.W,n=length(ncCR85.nb),
+         S0=Szero(col.W) )
>                 
> plot(moran.pboot)#Display results
> 
> #Poisson-Gamma model
> moran.pgboot<-boot(sids, statistic=moranI.pboot, sim="parametric", 
+         ran.gen=negbin.sim, R=niter,  listw=col.W,n=length(ncCR85.nb),
+         S0=Szero(col.W) )
>                 
> plot(moran.pgboot)#Display results
> 
> 
> 
> cleanEx(); ..nameEx <- "moranI.stat"
> 
> ### * moranI.stat
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: moranI
> ### Title: Compute Moran's I autocorrelation statistic
> ### Aliases: moranI.stat
> ### Keywords: spatial
> 
> ### ** Examples
> 
> library(spdep)
> data(nc.sids)
> col.W <- nb2listw(ncCR85.nb, zero.policy=TRUE)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74) )
> 
> moranI.stat(data=sids, listw=col.W, n=length(ncCR85.nb), S0=Szero(col.W) )
[1] 0.198295
> 
> moranI.stat(data=sids, applyto="SMR", listw=col.W, n=length(ncCR85.nb), 
+         S0=Szero(col.W) )
[1] 0.2385172
> 
> 
> 
> cleanEx(); ..nameEx <- "observed.sim"
> 
> ### * observed.sim
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: observed.sim
> ### Title: Randomly generate observed cases from different statistical
> ###   distributions
> ### Aliases: multinom.sim poisson.sim negbin.sim
> ### Keywords: distribution
> 
> ### ** Examples
> 
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> sids<-cbind(sids, x=nc.sids$x, y=nc.sids$y)
> 
> #Carry out simulations
> datasim<-multinom.sim(sids, mle=calculate.mle(sids, model="multinomal") )
> 
> #Estimators for Poisson distribution
> datasim<-poisson.sim(sids, mle=calculate.mle(sids, model="poisson") )
> 
> #Estimators for Negative Binomial distribution
> datasim<-negbin.sim(sids, mle=calculate.mle(sids, model="negbin") )
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "opgam"
> 
> ### * opgam
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: opgam
> ### Title: Openshaw's GAM
> ### Aliases: opgam opgam.intern
> ### Keywords: spatial
> 
> ### ** Examples
> 
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> sids<-cbind(sids, x=nc.sids$x, y=nc.sids$y)
> 
> #GAM using the centroids of the areas in data
> sidsgam<-opgam(data=sids,  radius=30, step=10, alpha=.002)
> 
> #Plot centroids
> plot(sids$x, sids$y, xlab="Easting", ylab="Northing")
> #Plot points marked as clusters
> points(sidsgam$x, sidsgam$y, col="red", pch="*")
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "pottwhitt.boot"
> 
> ### * pottwhitt.boot
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: pottwhitt.boot
> ### Title: Bootstrap replicates of Potthoff-Whittinghill's statistic
> ### Aliases: pottwhitt.boot pottwhitt.pboot
> ### Keywords: htest
> 
> ### ** Examples
> 
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> sids<-cbind(sids, x=nc.sids$x, y=nc.sids$y)
> 
> niter<-100
> 
> #Permutation model
> pw.boot<-boot(sids, statistic=pottwhitt.boot, R=niter)
> plot(pw.boot)#Plot results
> 
> #Multinomial model
> pw.mboot<-boot(sids, statistic=pottwhitt.pboot, sim="parametric", ran.gen=multinom.sim,  R=niter)
> plot(pw.mboot)#Plot results
> 
> #Poisson model
> pw.pboot<-boot(sids, statistic=pottwhitt.pboot, sim="parametric", ran.gen=poisson.sim,  R=niter)
> plot(pw.pboot)#Plot results
> 
> #Poisson-Gamma model
> pw.pgboot<-boot(sids, statistic=pottwhitt.pboot, sim="parametric", ran.gen=negbin.sim, R=niter)
> plot(pw.pgboot)#Plot results
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "pottwhitt.stat"
> 
> ### * pottwhitt.stat
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: pottwhitt.stat
> ### Title: Compute Potthoff-Whittinghill's statistic
> ### Aliases: pottwhitt.stat
> ### Keywords: htest
> 
> ### ** Examples
> 
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> sids<-cbind(sids, x=nc.sids$x, y=nc.sids$y)
> 
> pottwhitt.stat(sids)
$T
[1] 527848.8

$asintmean
[1] 444222

$asintvat
[1] 87955956

$pvalue
[1] 0

> 
> 
> 
> cleanEx(); ..nameEx <- "rmultin"
> 
> ### * rmultin
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: rmultin
> ### Title: Generate random observations from a multinomial distribution
> ### Aliases: rmultin
> ### Keywords: distribution
> 
> ### ** Examples
> 
> for(i in 1:10)
+         print(rmultin(10, c(1/3, 1/3, 1/3) ))
[1] 2 4 4
[1] 4 5 1
[1] 2 6 2
[1] 6 2 2
[1] 4 1 5
[1] 2 3 5
[1] 4 3 3
[1] 4 3 3
[1] 4 6 0
[1] 3 5 2
> 
> 
> 
> cleanEx(); ..nameEx <- "stone.boot"
> 
> ### * stone.boot
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: stone.boot
> ### Title: Generate boostrap replicates of Stone's statistic
> ### Aliases: stone.boot stone.pboot
> ### Keywords: spatial
> 
> ### ** Examples
> 
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74))
> sids<-cbind(sids, x=nc.sids$x, y=nc.sids$y)
> 
> niter<-100
> 
> #All Tests are performed around county 78.
> 
> #Permutation  model
> st.perboot<-boot(sids, statistic=stone.boot, R=niter, region=78)
> plot(st.perboot)#Display results
> 
> #Multinomial model
> st.mboot<-boot(sids, statistic=stone.pboot, sim="parametric", ran.gen=multinom.sim,  R=niter, region=78)
> plot(st.mboot)#Display results
> 
> #Poisson model
> st.pboot<-boot(sids, statistic=stone.pboot, sim="parametric", ran.gen=poisson.sim,  R=niter, region=78)
> plot(st.pboot)#Display results
> 
> #Poisson-Gamma model
> st.pgboot<-boot(sids, statistic=stone.pboot, sim="parametric", ran.gen=negbin.sim, R=niter, region=78)
> plot(st.pgboot)#Display results
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "tango.boot"
> 
> ### * tango.boot
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: tango.boot
> ### Title: Generate bootstrap replicated of Tango's statistic
> ### Aliases: tango.boot tango.pboot
> ### Keywords: spatial
> 
> ### ** Examples
> 
> library(boot)
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74) )
> sids<-cbind(sids, x=nc.sids$x, y=nc.sids$y)
> 
> #Calculate neighbours based on distance
> coords<-as.matrix(sids[,c("x", "y")])
> 
> dlist<-dnearneigh(coords, 0, Inf)
> dlist<-include.self(dlist)
> dlist.d<-nbdists(dlist, coords)
> 
> #Calculate weights. They are globally standardised but it doesn't
> #change significance.
> col.W.tango<-nb2listw(dlist, glist=lapply(dlist.d, function(x) {exp(-x)}),
+         style="C")
>         
> niter<-100
> 
> #Permutation model
> tn.boot<-boot(sids, statistic=tango.boot, R=niter, listw=col.W.tango, 
+         zero.policy=TRUE)
> plot(tn.boot)#Display results
> 
> #Multinomial model
> tn.mboot<-boot(sids, statistic=tango.pboot, sim="parametric", 
+         ran.gen=multinom.sim,  R=niter, listw=col.W.tango, zero.policy=TRUE)
>                 
> plot(tn.mboot)#Display results
> 
> #Poisson model
> tn.pboot<-boot(sids, statistic=tango.pboot, sim="parametric", 
+         ran.gen=poisson.sim,  R=niter, listw=col.W.tango, zero.policy=TRUE)
>                 
> plot(tn.pboot)#Display results
> 
> #Poisson-Gamma model
> tn.pgboot<-boot(sids, statistic=tango.pboot, sim="parametric", 
+         ran.gen=negbin.sim, R=niter, listw=col.W.tango, zero.policy=TRUE)
> plot(tn.pgboot)#Display results
> 
> 
> 
> cleanEx(); ..nameEx <- "tango.stat"
> 
> ### * tango.stat
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: tango.stat
> ### Title: Compute Tango's statistic for general clustering
> ### Aliases: tango.stat
> ### Keywords: spatial
> 
> ### ** Examples
> 
> library(spdep)
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74) )
> sids<-cbind(sids, x=nc.sids$x, y=nc.sids$y)
> 
> #Calculate neighbours based on distance
> coords<-as.matrix(sids[,c("x", "y")])
> 
> dlist<-dnearneigh(coords, 0, Inf)
> dlist<-include.self(dlist)
> dlist.d<-nbdists(dlist, coords)
> 
> #Calculate weights. They are globally standardised but it doesn't
> #change significance.
> col.W.tango<-nb2listw(dlist, glist=lapply(dlist.d, function(x) {exp(-x)}),
+         style="C")
> 
> niter<-100
> 
> #use exp(-D) as closeness matrix
> tango.stat(sids, col.W.tango, zero.policy=TRUE)
            [,1]
[1,] 0.003844675
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "whittermore.boot"
> 
> ### * whittermore.boot
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: whittermore.boot
> ### Title: Generate bootstrap replicates of Whittermore's statistic
> ### Aliases: whittermore.boot whittermore.pboot
> ### Keywords: spatial
> 
> ### ** Examples
> 
> library(boot)
> library(spdep)
> 
> data(nc.sids)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74) )
> sids<-cbind(sids, x=nc.sids$x, y=nc.sids$y)
> 
> #Calculate neighbours based on distance
> coords<-as.matrix(sids[,c("x", "y")])
> 
> dlist<-dnearneigh(coords, 0, Inf)
> dlist<-include.self(dlist)
> dlist.d<-nbdists(dlist, coords)
> 
> #Calculate weights. They are globally standardised but it doesn't
> #change significance.
> col.W.whitt<-nb2listw(dlist, glist=dlist.d, style="C")
> 
> niter<-100
> 
> #Permutation model
> wt.boot<-boot(sids, statistic=whittermore.boot, R=niter, listw=col.W.whitt,
+         zero.policy=TRUE)
> plot(wt.boot)#Display results
> 
> #Multinomial model
> wt.mboot<-boot(sids, statistic=whittermore.pboot, sim="parametric", 
+         ran.gen=multinom.sim,  R=niter,  listw=col.W.whitt, zero.policy=TRUE)
>                 
> plot(wt.mboot)#Display results
> 
> #Poisson model
> wt.pboot<-boot(sids, statistic=whittermore.pboot, sim="parametric", 
+         ran.gen=poisson.sim,  R=niter,  listw=col.W.whitt, zero.policy=TRUE)
>                 
> plot(wt.pboot)#Display results
> 
> #Poisson-Gamma model
> wt.pgboot<-boot(sids, statistic=whittermore.pboot, sim="parametric", 
+         ran.gen=negbin.sim, R=niter, listw=col.W.whitt, zero.policy=TRUE)
> plot(wt.pgboot)#Display results
> 
> 
> 
> cleanEx(); ..nameEx <- "whittermore.stat"
> 
> ### * whittermore.stat
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: whittermore.stat
> ### Title: Compute Whittermore's statistic
> ### Aliases: whittermore.stat
> ### Keywords: spatial
> 
> ### ** Examples
> 
> library(spdep)
> data(nc.sids)
> col.W <- nb2listw(ncCR85.nb, zero.policy=TRUE)
> 
> sids<-data.frame(Observed=nc.sids$SID74)
> sids<-cbind(sids, Expected=nc.sids$BIR74*sum(nc.sids$SID74)/sum(nc.sids$BIR74) )
> sids<-cbind(sids, x=nc.sids$x, y=nc.sids$y)
> 
> #Calculate neighbours based on distance
> coords<-as.matrix(sids[,c("x", "y")])
> 
> dlist<-dnearneigh(coords, 0, Inf)
> dlist<-include.self(dlist)
> dlist.d<-nbdists(dlist, coords)
> 
> #Calculate weights. They are globally standardised but it doesn't
> #change significance.
> col.W.whitt<-nb2listw(dlist, glist=dlist.d, style="C")
> 
> whittermore.stat(sids, col.W.whitt, zero.policy=TRUE)
            [,1]
[1,] 0.007928723
> 
> 
> 
> ### * <FOOTER>
> ###
> cat("Time elapsed: ", proc.time() - get("ptime", env = .CheckExEnv),"\n")
Time elapsed:  18.37 0.26 18.74 0 0 
> grDevices::dev.off()
null device 
          1 
> ###
> ### Local variables: ***
> ### mode: outline-minor ***
> ### outline-regexp: "\\(> \\)?### [*]+" ***
> ### End: ***
> quit('no')