[GIS] How to connect the point pattern shapefile with the polygon window in R

rshapefilespatstat

I would like to do some spatial statistics in R. I have a point layer called '1983' that represents the distribution of one species, and a polygon layer called 'SA-R' that represents the study area.

I followed the steps from the pdf "Handling shapefiles in the spatstat package":

library(spatstat)

spatstat 1.36-0 (nickname: ‘Intense Scrutiny’)
For an introduction to spatstat, type ‘beginner’

library(sp)

library(maptools)

Checking rgeos availability: FALSE
Note: when rgeos is not available, polygon geometry computations in maptools depend on gpclib,
which has a restricted licence. It is disabled by default;
to enable gpclib, type gpclibPermit()

Year1 <- readShapeSpatial("1983.shp")

Studyarea <- readShapeSpatial("SA-R.shp")

Field name: ‘opp m²’ changed to: ‘opp.m²’
Field name: ‘opp ha’ changed to: ‘opp.ha’

as.owin(Studyarea)

window: polygonal boundary
enclosing rectangle: [226675.38, 234961.33] x [21233.66, 32693.92] units

as.ppp(Year1)

marked planar point pattern: 7 points
Mark variables: ID, Year
window: rectangle = [227669.54, 229542.88] x [22341.009, 30989.048] units

Data1 <- Year1[Studyarea]

Error in x[j] : invalid subscript type 'S4'

I can plot each of my shapefiles in R, but do not manage to plot my point layer within my study area polygon. Does someone know what's wrong with my last command?
I have already checked, and both of my shapefile layers are in the same projection system.

Your syntax works until this step:

w -> as.owin(Studyarea)
Error in as.owin.default(Studyarea) : Can't interpret W as a window

For some reason, it did work like this:

library(spatstat)
library(maptools)
library(sp)
Year1 <- readShapeSpatial("1983.shp")
as.ppp(Year1)

marked planar point pattern: 7 points
Mark variables: ID, Year
window: rectangle = [227669.54, 229542.88] x [22341.009, 30989.048] units

Studyarea <- readShapeSpatial("SA-R.shp")

Field name: ‘opp m²’ changed to: ‘opp.m²’
Field name: ‘opp ha’ changed to: ‘opp.ha’

w <- as.owin(Studyarea)

Year1 <- ppp(coordinates(Year1)[,1],coordinates(Year1)[,2], window=w, marks=Year1@data$MyVar)

plot(Year1)

Best Answer

First, I would recommend using rgdal for reading/writing shapefiles. It looks like you are trying to directly assign your owin object to your ppp object. This will not work in the way you are thinking. You can, in fact, specify the owin window object when creating the ppp object. This is detailed in the as.ppp help which can be accessed using: ?as.ppp

Is this a marked process problem (values associated with the points)? If so then you need to change the way you are coercing the point data, using ppp directly, so you can define the marks. 

# Add required packages 
require(sp)
require(rgdal)
require(spatstat)  

# Set working directory and read shapefiles
setwd("C:/MyData")
  Year1 <- readOGR(getwd(), "1983")
  Studyarea <- readOGR(getwd(), "SA-R")

# Coerce study area to owin object
w <- as.owin(Studyarea)

# Coerce points to ppp object and pass the object the study area window 
Year1 <- as.ppp(Year1, W=w)

# Alternate method with marks
Year1 <- ppp(coordinates(Year1)[,1],coordinates(Year1)[,2], window=w, marks=Year1@data$MyVar) 

# Plot points with study area window
plot(Year1)
Related Question