[R-sig-Geo] Using pixel image as spatstat's intensity surface

Allar Haav allar.haav at gmail.com
Tue Jul 23 03:23:02 CEST 2013


Hi all,

I am currently trying to use spatstat's inhomogeneous cluster analysis 
methods (mostly Linhom), but ran into a problem. Namely, I have a pixel 
image (type "im") with values ranging from 0 to 1 indicating point 
probabilities. When creating and plotting random points with it, this 
intensity surface works as it should. I can also use it in regular 
envelope functions without issues. However, any attempt to use the 
intensity surface with Linhom function results in empty plot. I 
suspected that this might have something to do with the integral of the 
surface, therefore made a crude transformation so that it would be 
approximately the number of points under study (see below), but of 
course it didn't change anything. Generally I still feel that the crux 
lies in the intensity surface's values, but I can't quite put my finger 
on it.

A simplified code I used is as follows:

area <- as.owin(area_polygon)
sites <- as.ppp(as(readShapeSpatial('sites.shp')))
sites$window <- area
model <- as.im(read.asciigrid('model.asc'))
sites.n <- sites$n
integral <- summary(model)$integral

# This creates a nice plot with simulation envelopes
plot(envelope(sites, Lest, nsim=19, 
simulate=expression(rpoispp(eval.im(model*(sites.n/integral ))))))

# These produce only empty plots
plot(envelope(sites, Linhom, nsim=19, lambda=model ))
plot(envelope(sites, Linhom, nsim=19, lambda=eval.im(model 
*(sites.n/integral ))))

# A further check without envelope
Linhom(sites, lambda=model)
# This resulted in:
"Error in Kwtsum(dIJ, bI, wIJ, b, w = reciplambda, breaks) : Weights in 
K-function were infinite or NA"

That last error message left me a bit perplexed. Indeed, I have NA 
values in the dataset as the study area is not perfectly rectangular. 
But this should be a non-issue anyway.
Frankly I'm running out of ideas. I've been thinking of alternative 
approaches, e.g. creating random points from the intensity surface and 
fitting a model based on it using "ppm", perhaps using some covariates. 
Yet this seems actually an awful idea considering the probable loss of 
accuracy and the fact that basically the model already exists.

As a side note, a bit of googling on the issue revealed an older post 
about inhomogeneous process analysis 
(http://r-sig-geo.2731867.n2.nabble.com/memory-allocation-problem-with-envelope-function-in-spatstat-td5029218.html). 
There a command equivalent to my first envelope plot was used: 
Lin=envelope(SPECIES.ppp,Linhom,nsim=5,simulate=expression(rpoispp(lambda)), 
! + correction="border")
That made me think that I might be overthinking the problem and I've 
already found the solution with first try. However I am still skeptical 
as the difference between Lest and Linhom does not lie merely in 
simulated patterns - weights are also to be considered. If so, why was 
this solution used?

Any suggestions or ideas? All sort of feedback would be appreciated.

Regards,
Allar Haav



More information about the R-sig-Geo mailing list