[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