Search code examples
rclassificationgiskernel-densityspatstat

Point pattern classification with spatstat: how to choose the right bandwidth?


I'm still trying to find the best way to classify bivariate point patterns:

Point pattern classification with spatstat: what am I doing wrong?

I now analysed 110 samples of my dataset using @Adrian's suggestion with sigma=bw.diggle (as I wanted an automatic bandwidth selection). f is a "resource selection function" (RSF) which describes the relationship between the intensity of the Cancer point process and the covariate (here kernel density of Immune):

Cancer <- split(cells)[["tumor"]]
Immune <- split(cells)[["bcell"]]
Dimmune <- density(Immune,sigma=bw.diggle)
f <- rhohat(Cancer, Dimmune)

I am in doubt about some results I've got. A dozen of rho-functions looked weird (disrupted, single peak). After changing to default sigma=NULL or sigma=bw.scott (which are smoother) the functions became "better" (see examples below). I also experimented with the following manipulations:

  cells # bivariate point pattern with marks "tumor" and "bcell"
  o.marks<-cells$marks  # original marks
    #A) randomly re-assign original marks
  a.marks <- sample(cells$marks)
    #B) replace marks randomly with a 50/50 proportion
  b.marks<-as.factor(sample(c("tumor","bcell"), replace=TRUE, size=length(o.marks))) 
    #C) random (homogenious?) pattern with the original number of points 
  randt<-runifpoint(npoints(subset(cells,marks=="tumor")),win=cells$window) 
  randb<-runifpoint(npoints(subset(cells,marks=="bcell")),win=cells$window)
  cells<-superimpose(tumor=randt,bcell=randb)
    #D) tumor points are associated with bcell points (is "clustered" a right term?)
  Cancer<-rpoint(npoints(subset(cells,marks=="tumor")),Dimmune,win=cells$window)
    #E) tumor points are segregated from bcell points
  reversedD<-Dimmune
  density.scale.v<-sort(unique((as.vector(Dimmune$v)[!is.na(as.vector(Dimmune$v))]))) # density scale
  density.scale.v.rev<-rev(density.scale.v)# reversed density scale
  new.image.v<-Dimmune$v
  # Loop over matrix
  for(row in 1:nrow(Dimmune$v)) {
    for(col in 1:ncol(Dimmune$v)) {
      if (is.na(Dimmune$v[row, col])==TRUE){next}
      number<-which(density.scale.v==Dimmune$v[row, col])
      new.image.v[row, col]<-density.scale.v.rev[number]}
  }
  reversedD$v<-new.image.v # reversed density
  Cancer<-rpoint(npoints(subset(cells,marks=="tumor")),reversedD,win=cells$window)

A better way to generate inverse density heatmaps is given by @Adrian in his post below.

I could not generate rpoint patterns for the bw.diggle density as it produced negative numbers.Thus I replaced the negatives Dimmune$v[which(Dimmune$v<0)]<-0 and could run rpoint then. As @Adrian explained in the post below, this is normal and can be solved easier by using a density.ppp option positive=TRUE.

I first used bw.diggle, because hopskel.test indicarted "clustering" for all my patterns. Now I'm going to use bw.scott for my analysis but can this decision be somehow justified? Is there a better method besides "RSF-function is looking weird"?

some examples:

sample10: sample10 bw.diggle vs bw.scott sample20: sample20 bw.diggle vs bw.scott sample110: sample110 bw.diggle vs bw.scott


Solution

  • That is a lot of questions!

    Please try to ask only one question per post.

    But here are some answers to your technical questions about spatstat.

    Negative values: The help for density.ppp explains that small negative values can occur because of numerical effects. To force the density values to be non-negative, use the argument positive=TRUE in the call to density.ppp. For example density(Immune, bw.diggle, positive=TRUE).

    Reversed image: to reverse the ordering of values in an image Z you can use the following code:

    V <- Z
    A <- order(Z[])
    V[][A] <- Z[][rev(A)]
    

    Then V is the order-reversed image.