1

I am currently exploring three shapefiles, each with point data, and all confined to the same window. I also have them in ppp format which I've used to create kernel density maps.

plot(density.ppp(smktppp, 0.5, edge=T), main="Supermarket Density")
plot(density.ppp(tptppp, 0.5, edge=T), main="Transport Density")
plot(density.ppp(farmppp, 0.5, edge=T), main="Urban Farm Density")

I would like to overlay these plots, using map algebra, or fuzzy logic, etc, to create one output map showing the density of the three combined. How would I go about doing this in R?

Ege Rubak
  • 4,347
  • 1
  • 10
  • 18
skyandbuildings
  • 9
  • 1
  • Welcome to StackOverflow! Please read the info about [how to ask a good question](http://stackoverflow.com/help/how-to-ask) and how to give a [reproducible example](http://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example/5963610). This will make it much easier for others to help you. – Jaap Feb 25 '17 at 17:57
  • Part of asking a good question is to make life easier for the answerers. That includes putting in the needed library calls to load the needed functions and putting in some data to test and demonstrate. – IRTFM Feb 25 '17 at 19:52

2 Answers2

0

If you simply want to estimate the overall density (usually called intensity since it doesn't integrate to one) of points disregarding whether it is "Supermarket", "Transport" or "Urban Farm" you just combine all the points and do as before:

library(spatstat)
combined <- superimpose(smktppp, tptppp, farmppp)
plot(density(combined), main="Density of all points.")

Of course you can choose the smoothing bandwidth to be 0.5 as before or any other value you like.

You can also do normal algebra with the raster images produced by density.ppp (object of class im). If you saved these as smktim, tptim and farmim you can do something like rslt <- smktim + tptim + farmim to get the sum of the three estimates.

Ege Rubak
  • 4,347
  • 1
  • 10
  • 18
0

Alternatively, if you combine the three point patterns into a single object

library(spatstat)
X <- superimpose(Supermarket=smktppp, Transport=tptppp, Farm=farmppp)

then you can display both the original data and the intensity estimates, with or without identifying the kind of point:

# original data:
plot(X)  # single plot with 3 different plot characters for 3 types 
plot(split(X)) # three plot panels, one for each type of point
plot(unmark(X)) # single plot without distinguishing types of points

# intensity images:
plot(density(X, 0.5)) # single plot: intensity regardless of type
plot(density(split(X), 0.5)) # three panels: intensity for each type
plot(relrisk(X), 0.5) # three panels: relative probabilities of each type

See the spatstat book for details.

Adrian Baddeley
  • 1,956
  • 1
  • 8
  • 7