[R-sig-Geo] [FORGED] Simple Ripley's CRS test for market point patters

[Rolf Turner]

Thanks to Marcelino for chiming in on this.  I should have responded 
earlier, but was "busy on a personal matter".

To add to what Marcelino said:

(1) Your post in HTML was horribly garbled and a struggle to read. 
*PLEASE* set your email client so that it does *NOT* post in HTML when 
posting to this list.

(2) A very minor point:  Your construction of "syn.ppp" was 
unnecessarily complicated and convoluted.  You could simply do:

     syn.ppp  <- ppp(x=xp,y=yp,window=W,marks=area)

(3) You appear to be confused as to the distinction between "marks" and 
"covariates".   These are superficially similar but conceptually 
completely different.  What you are dealing with is *marks*.  There are
no covariates involved.  See [1], p. 147.

(4) Your calls to envelope() are mal-formed; the expression "area >= 0" 
and "area < 25" will be taken as the values of "nrank" and ... who knows 
what?  Did you not notice the warning messages you got about there being 
something wrong with "nrank"?

You are being hopelessly naïve in expecting envelope() to interpret 
"area >= 0" and "area < 25" in the way you want it to interpret them.
The spatstat package does not read minds.

Marcelino has shown you how to make a proper call.

(5) Since you are interested in categorising "area" into groups, rather 
than being interested in the *numerical* value of area, you should do 
the categorisation explicitly:

acat <- cut(area,breaks=c(-Inf,25,55,Inf),right=FALSE,
             labels=c("s","m","l")
# right=FALSE since you want "area < 25" rather than
# "area <= 25" --- although this makes no difference for the
# area values actually used.

syn.ppp <- ppp(x=xp,y=yp,marks=acat)

(6) It is not clear to me what your "research question" is.  Do you want 
to test whether the intensities differ between the area categories? 
Unless my thinking is hopelessly confused, this has nothing to do with 
(or at least does not require the use of) the envelope() function:

f1 <- ppm(syn.ppp) # Same (constant) intensity for each area category.
f2 <- ppm(syn.ppp ~ marks) # Allows different (constant) intensity
                            # for each area category.
anova(f1,f2,test="Chi")

This test (not surprisingly) rejects the hypothesis that the intensities 
are the same; p-value = 0.0020.

If this is not the hypothesis that you are interested in, please clarify 
your thinking and your question, and get back to us.

cheers,

Rolf Turner

-- 
Honorary Research Fellow
Department of Statistics
University of Auckland
Phone: +64-9-373-7599 ext. 88276

[1] Spatial Point Patterns: Methodology and Applications with R
1st Edition, Adrian Baddeley, Ege Rubak, Rolf Turner
Chapman and Hall/CRC, 2015

P. S.  It is of course (!!!) highly recommended that you spend some time 
reading the book cited above, if you are going to work in this area.

R. T.

On 23/07/19 9:36 AM, Alexandre Santos via R-sig-Geo wrote:
> Dear R-Sig-Geo Members,     I"ve like to find any simple way for apply CRS test for market point patters, for this I try to create a script below:
> #Packages require(spatstat)require(sp)
> 
> # Create some points that represents ant nests in UTMxp<-c(371278.588,371250.722,371272.618,371328.421,371349.974,371311.95,371296.265,371406.46,371411.551,371329.041,371338.081,371334.182,371333.756,371299.818,371254.374,371193.673,371172.836,371173.803,371153.73,371165.051,371140.417,371168.279,371166.367,371180.575,371132.664,371129.791,371232.919,371208.502,371289.462,371207.595,371219.008,371139.921,371133.215,371061.467,371053.69,371099.897,371108.782,371112.52,371114.241,371176.236,371159.185,371159.291,371158.552,370978.252,371120.03,371116.993)
> yp<-c(8246507.94,8246493.176,8246465.974,8246464.413,8246403.465,8246386.098,8246432.144,8246394.827,8246366.201,8246337.626,8246311.125,8246300.039,8246299.594,8246298.072,8246379.351,8246431.998,8246423.913,8246423.476,8246431.658,8246418.226,8246400.161,8246396.891,8246394.225,8246400.391,8246370.244,8246367.019,8246311.075,8246255.174,8246255.085,8246226.514,8246215.847,8246337.316,8246330.197,8246311.197,8246304.183,8246239.282,8246239.887,8246241.678,8246240.361,8246167.364,8246171.581,8246171.803,8246169.807,8246293.57,8246183.194,8246189.926)
> # Then I create the size of each nest - my covariate used as marked processarea<-c(117,30,4,341,15,160,35,280,108,168,63,143,2,48,182,42,88,56,27,156,288,45,49,234,72,270,91,40,304,56,35,4,56.7,9,4.6,105,133,135,23.92,190,12.9,15.2,192.78,104,255,24)
> 
> # Make a countour - only as exerciseW <- convexhull.xy(xp,yp)
> #Create a ppp objectp_xy<-cbind(xp,yp)syn.ppp<-ppp(x=coordinates(p_xy)[,1],y=coordinates(p_xy)[,2],window=W, marks=area)syn.ppp <- as.ppp(syn.ppp)plot(syn.ppp, main=" ")
> First, I've like to study CSR of market point process (my hypothesis is that different size have a same spatial distribution) when area >= 0, area < 25 and area >=25, area < 55, for this I make:
> # Area 0-25env.formi1<-envelope(syn.ppp,nsim=99,fun=Kest, area >= 0, area < 25)plot(env.formi1,lwd=list(3,1,1,1), main="")
> # Area 25-55env.formi2<-envelope(syn.ppp,nsim=99,fun=Kest, area >=25, area < 55)plot(env.formi2,lwd=list(3,1,1,1), main="")
> My first problem is that I have the same plot in both conditions and I don't know why.
> Second, if I try to estimate the market intensity observed pattern
> est.int <- density(syn.ppp)est_xy <-  rmpoispp(est.int)plot(est_xy, main=" ")
> My output is only my points position without marked area in my ppp object created.
> My question is what is the problem with this Ripley's reduced second moment function approach? There are any way for study my point process when the area is a covariate of my point process?
> Thanks in advanced
> Alexandre