[R-sig-Geo] How set up spatstat::densityHeat()?

[Rolf Turner]

On 14/07/20 1:39 am, Michal Kvasnička wrote:

> Hello.
> 
> I have to find hotspots of road accidents. It seems that
> spatstat::densityHeat() implements McSwiggan et al. (2017) algorithm I want
> to use. I plan to follow path of DrHotNet package:
> 
> 1. estimate the density function
> 2. lixelize the network
> 3. calculate the density in the middle point of each lixel from (2) using
> function estimated in (1)
> 4. find clusters of lixels from (3) with sufficiently high density
> 
> However, I don't know how to set the function's parameters.
> 
> I cannot use the implicit setting (finespacing=TRUE) because it suggests
> tens of billions iterations. (I use OSM road maps of Czech districts
> converted to Krovak projection. My testing district has about 2,700 km of
> roads with an average length of road segments of about 30 meters due to
> road curvatures. This is not the largest district I have to handle---Prague
> is much larger.)
> 
> I can set finespacing = FALSE and then finetune the parameters myself.
> However, I can't find any documentation. Having skimmed the source code,
> I'm trying something like this:
> 
> density_function <- densityHeat(accident_lpp, sigma, finespacing = FALSE,
> eps = 10, fun = TRUE)
> 
> This way I get the function I need. However, I'm not certain what I'm
> doing. It seems that the densityHeat() uses some grid, based on lixelized
> network. My questions are:
> 
> * Is it sufficient to set the eps parameter or do I need to set other
> parameters too? Which ones?
> * Is the parameter in meters (my map is projected with meters as its units)?
> * If eps = 10 means grid of 10x10 meters, what precisely does it mean? How
> is the network lixelized (and why)? How does it affect the estimated
> function (not the linim object)?
> * What is the (economically/substantially) reasonable value of eps (and
> other necessary parameters)?
> * Is there any documentation to the function? (I've read chapter 17 of the
> spatstat book but found nothing there. The function help page doesn't cover
> this either.)
> * Does it any difference if I lixelize the network myself before running
> densityHeat()?
> 
> Many thanks for your advice. (And an apology for such a long question.)
> 
> Best wishes,
> *Ing. Michal Kvasnička, Ph.D.*

Dear Michal,

I am replying on behalf of the spatstat authors (and not claiming any 
personal competence in the area of point processes on networks).

My co-authors suggest that you:

    * upgrade to the latest development version of spatstat.  You
      could use

      > remotes::install_github("spatstat/spatstat",lib=.Rlib)

      to do this,

    * look at help(densityHeat), in particular the
      section "Discretisation and Error Messages" which explains
      how the discretisation is performed,

    * use densityfun.lpp instead of density.lpp if you want a
      function rather than a pixel image.

Note that if finespacing=FALSE is selected, the discretisation is 
determined by the default pixel resolution of the window containing the 
network. The arguments eps and dimyx would be passed to as.mask to 
determine pixel resolution in the usual way (their interpretation is 
explained in the help for as.mask()).

Note also that density.lpp() calls densityHeat() which calls the 
internal algorithm. Similarly densityfun.lpp() calls the same internal 
algorithm. Consequently they all handle discretisation arguments the 
same way.

The spatstat authors know nothing about the package DrHotNet and do not 
necessarily endorse whatever they are doing. Definitely *do not* 
'lixellate' (discretise) the network before applying density(), because 
this will screw up the discretisation procedure.

We are very sceptical about the validity of applying clustering to the 
elements of the discretisation.

HTH

cheers,

Rolf Turner

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