[R-sig-Geo] Loosmore and Ford (2006) Goodness-of-fit (GoF) test for CSR
marcelino.delacruz at upm.es
marcelino.delacruz at upm.es
Fri Mar 2 08:18:11 CET 2012
Hi Aurelie,
Please find below my own version of the Loosmore and Ford GOF test.
You only need to compute an envelope object with the argument
savefuns=TRUE, and then run Lof.gof on the result.
Example:
>data(nztrees)
>nztrees.env <- envelope(nztrees, Kest, r=seq(0, 25, by=0.5),savefuns=TRUE)
> LF.gof(nztrees.env)
$u
[1] 126421.9
$p
[1] 0.29
LF.gof gives you the result of the summary statistic (u) and its p-value
(p)
The arguments of LF.gof are "X" = the envelope object and "rmin" and
"rmax" that define the limits of integration to compute the summary
statistics.
By the way, maybe I'm wrong but it seems to me that you are trying to
test the random labelling hypothesis instead of CSR.
If this is true, you should go this way:
# first create a list of simulated patterns
nsim <- 99
GROT <- superimpose(GR, OT)
sim <-list()
for (i in 1:nsim) {sim[[i]] = rlabel(GROT); sim[[i]] <-
sim[[i]][sim[[i]]$marks=="GR"]}
# then compute the envelopes
GR.env <- envelope (GR, Kest, simulate = sim)
# then compute the GOF test
LF.gof(GR.env)
HTH,
Marcelino
#############################################################################
#
# Loosmore and Ford (2006) GOF test
#
#############################################################################
LF.gof<- function(X, rmin=NULL, rmax=NULL){
require(spatstat)
gmin <- NULL
gmax <- NULL
r <- as.data.frame(attributes(X)$simfuns)[,1]
H1 <- X$obs
Hi <- as.data.frame(attributes(X)$simfuns)[,-1]
#Join observed and simulated functions
H <- cbind(H1,Hi)
# select the range of integrable results
if (!is.null(rmin)) gmin <- r >= rmin
if (!is.null(rmax)) gmax <- r <= rmax
if (!is.null(gmin) & !is.null(gmax)){ simfuns <- simfuns[gmin &
gmax,]; Hobs <- Hobs[gmin & gmax]}
if (!is.null(gmin) & is.null(gmax)) {simfuns <- simfuns[gmin,]; Hobs
<- Hobs[gmin]}
if (is.null(gmin) & !is.null(gmax)) {simfuns <- simfuns[gmax,]; Hobs
<- Hobs[gmax]}
s <- dim(H)[2]
# compute summary statistics for the observed (i.e., first) pattern
Hmeans <- rowSums (H[,-1], na.rm = FALSE)/(s-1)
u <- sum((H[,1]-Hmeans)^2)
# compute summary statistics for the simulated patterns
for ( i in 2:s) {
Himeans <- rowSums (H[,-i])/(s-1)
u <- c(u, sum((H[,i]-Himeans)^2))
}
p <- 1-((rank(u)[1]-1)/(s))
return(list(u=u[1], p=p))
}
###############################################################################
Con fecha 1/3/2012, "Aurelie C.Godin" <godina at dal.ca> escribió:
>Dear all,
>
>I am trying to use Loosmore and Ford (2006) GoF test for CSR. I was
>wondering if anyone here is familiar with this test and could help me coding
>it with my data. Or perhaps suggest something else.
>
>The authors provided an R source code and user's guide for the
>implementation of their statistical test described in their paper. Available
>at: http://www.esapubs.org/archive/ecol/E087/120/suppl-1.htm#anchorAuthors
>http://www.esapubs.org/archive/ecol/E087/120/suppl-1.htm#anchorAuthors
>
>There are multiple steps and I am not quiet sure if I need to run them all
>or not... but regardless, I get an error message at step 1 (see below)!
>
> "GR" is my point pattern (with my marks of specific interest) and it is a
>subset of "OT". Basically, I want to test the Ho : that my marks (GR) are
>distributed randomly within the non-random spatial distribution of my event
>(OT).
>
>http://r-sig-geo.2731867.n2.nabble.com/file/n7332928/OTGR.RData OTGR.RData
>
>> OT
> planar point pattern: 2691 points
>window: polygonal boundary
>enclosing rectangle: [-80, -57.6343] x [66.25, 78.16666] units
>
>> GR
> planar point pattern: 563 points
>window: polygonal boundary
>enclosing rectangle: [-80, -57.6343] x [66.25, 78.16666] units
>
>Loosmore & Ford GoF
>
>source("http://www.esapubs.org/archive/ecol/E087/120/cedl_scr.r")
>library(spatstat)
># STEP 1 in Loosmore & Ford R code
># gen.ext.pp.list() generates an external CSR point pattern list.
>#Note, this function overwrites any externally named variables `my.pp.list'
>and `pp.indx' if they already exist.
># input variables:
># - nsim.pp is the number of simulated point patterns in this list
># - npts is the number of points in each pattern
># - type is the pattern type to be generated
>
>gen.ext.pp.list <- function(99, GR$n, type="CSR")
>{
> # create a blank list at the top level environment
> my.pp.list <<- c() # should be equivalent to <<- NULL
> # now popoulate it with the appropriate patterns.
> for (i in 1:99)
> if (type=="CSR")
> my.pp.list[[i]] <<- runifpoint(GR$n)
> else # for now, just have an inhibition option
> my.pp.list[[i]] <<- rSSI(r=0.0045, n=GR$n,
> giveup=10000)
>
> # lastly, create the top level index to the pattern list
> pp.indx <<- 1
> # nothing to return because the list is created externally.
> return(0)
>}
>
>>Error: object 'type' not found
>
>Thank you very much in advance,
>Best,
>
>
>
>
>
>-----
>Aurelie Cosandey-Godin
>Ph.D. Student, Department of Biology, Dalhousie University
>Industrial Graduate Fellow, WWF-Canada
>Email: godina at dal.ca | Web: wormlab.biology.dal.ca
>
>--
>View this message in context: http://r-sig-geo.2731867.n2.nabble.com/Loosmore-and-Ford-2006-Goodness-of-fit-GoF-test-for-CSR-tp7332928p7332928.html
>Sent from the R-sig-geo mailing list archive at Nabble.com.
>
>_______________________________________________
>R-sig-Geo mailing list
>R-sig-Geo at r-project.org
>https://stat.ethz.ch/mailman/listinfo/r-sig-geo
More information about the R-sig-Geo
mailing list