How set up spatstat::densityHeat()?

classic Classic list List threaded Threaded
4 messages Options
Reply | Threaded
Open this post in threaded view
|

How set up spatstat::densityHeat()?

Michal Kvasnička
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.*



*Masaryk University | School of Economics and Administration*
Department of Economics
A: Lipová 41a | 602 00 Brno
E: [hidden email] | W: www.econ.muni.cz/~qasar

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: 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

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: How set up spatstat::densityHeat()?

Michal Kvasnička
Dear Rolf.

Thanks a lot for the suggestions! I'll check it out.

Best wishes to you and all spatstat team,
*Ing. Michal Kvasnička, Ph.D.*



*Masaryk University | School of Economics and Administration*
Department of Economics
A: Lipová 41a | 602 00 Brno
E: [hidden email] | W: www.econ.muni.cz/~qasar


pá 17. 7. 2020 v 7:34 odesílatel Rolf Turner <[hidden email]>
napsal:

>
> 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
>

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: How set up spatstat::densityHeat()?

Rolf Turner

On 17/07/20 9:04 pm, Michal Kvasnička wrote:

> Dear Rolf.
>
> Thanks a lot for the suggestions! I'll check it out.
>
> Best wishes to you and all spatstat team,
> *Ing. Michal Kvasnička, Ph.D.*


Good luck with your endeavours.  I'd like to point out a slight glitch
in the advice that I proffered, explicitly my suggestion that you use

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

to install the latest version of spatstat.  The code above was copied
and pasted from my .Rhistory and I neglected to  change the "lib=.Rlib",
which is peculiar to my personal R set-up.

You should replace ".Rlib" with a text string specifying the library
where you keep the packages that you install.

Apologies for any confusion that may have resulted.

cheers,

Rolf

<SNIP>

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

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo