library(sf)
library(spatstat)

The zip file ‘urkiola.zip’ contains two shapefiles. The first urkiola.shp contains point locations of two species of tree (oak and birch) in a Spanish National Park. The second, urkiolaWindow.shp, contains the park boundary as a polygon. Code for reading in these files, and converting them into a point process object (ppp) is given below. Once you have obtained the ppp object:

# Tree locations
urkiola.sf <- st_read("../datafiles/urkiola/urkiola/urkiola.shp")
## Reading layer `urkiola' from data source 
##   `C:\Users\erics\OneDrive\Documents\geog6000\datafiles\urkiola\urkiola\urkiola.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 1245 features and 5 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: 2.1 ymin: 0.1 xmax: 219.1 ymax: 149.9
## CRS:           NA
# Park boundary
urkiola.win <- st_read("../datafiles/urkiola/urkiola/urkiolaWindow.shp")
## Reading layer `urkiolaWindow' from data source 
##   `C:\Users\erics\OneDrive\Documents\geog6000\datafiles\urkiola\urkiola\urkiolaWindow.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 1 feature and 3 fields
## Geometry type: POLYGON
## Dimension:     XY
## Bounding box:  xmin: 0.05 ymin: 0.05 xmax: 219.95 ymax: 149.95
## CRS:           NA
urkiola.win <- as.owin(urkiola.win)
# First extract coordinates
urkiola.crds <- st_coordinates(urkiola.sf)
## Convert to ppp
urkiola.ppp <- ppp(x = urkiola.crds[, 1], y = urkiola.crds[,2], 
                   marks = as.factor(urkiola.sf$tree), window = urkiola.win)
plot(urkiola.ppp)

Give the intensity (spatial density) of the combined set of two tree species by using the summary() function

summary(urkiola.ppp)
## Marked planar point pattern:  1245 points
## Average intensity 0.06564029 points per square unit
## 
## Coordinates are given to 1 decimal place
## i.e. rounded to the nearest multiple of 0.1 units
## 
## Multitype:
##       frequency proportion  intensity
## birch       886  0.7116466 0.04671269
## oak         359  0.2883534 0.01892760
## 
## Window: polygonal boundary
## single connected closed polygon with 44 vertices
## enclosing rectangle: [0.05, 219.95] x [0.05, 149.95] units
##                      (219.9 x 149.9 units)
## Window area = 18967 square units
## Fraction of frame area: 0.575

the intensity is about 0.066 points per square unit

Plot the distribution of the two tree species as separate figures. You will need to use the split() function to access the data by individual species (e.g. split(urkiola.ppp)$oak will return only the location of the oak trees)

urkiola.oak = split(urkiola.ppp)$oak

urkiola.birch = split(urkiola.ppp)$birch

plot(urkiola.oak)

plot(urkiola.birch)

Now make plots of the kernel density of each species, and compare these plots. Do they suggest that the two species co-occur within the park?

plot(density(urkiola.oak, sigma = 10))

plot(density(urkiola.birch, sigma = 10))

It seems like the two species may co-occur, but not completely overlapping. Especially in the bottom-most portion of the window where both species are dense, there looks to be some overlap. However, I am pretty inconclusive with a visual check.

Using Ripley’s K function, examine the individual distributions of the two species for non-random spatial organization. Use the envelope() function to test the obtained Ripley’s K against a hypothesis of Complete Spatial Randomness. Plot the results and state whether or not each species has a random distribution or not, and if not, whether the trees are clustered or have a regular distribution. Note that using summary() on the output of the envelope() function with give the significance level of the Monte Carlo test

oak.kest = Kest(urkiola.oak)
birch.kest = Kest(urkiola.birch)

plot(birch.kest)

plot(oak.kest)

oak.kest.mc = envelope(urkiola.oak, 
                            fun = 'Kest', 
                            nsim = 99, 
                            verbose = FALSE, global = TRUE)

plot(oak.kest.mc, shade = c("hi", "lo"))

summary(oak.kest.mc)
## Simultaneous critical envelopes for K(r)
## and observed value for 'urkiola.oak'
## Obtained from 99 simulations of CSR
## Alternative: two.sided
## Envelopes computed as theoretical curve plus/minus maximum simulated value of 
## maximum absolute deviation
## Significance level of Monte Carlo test: 1/100 = 0.01
## Data: urkiola.oak

The oak trees seem to be spatially clustered; significance is indicated at .01

birch.kest.mc = envelope(urkiola.birch,
                         fun = 'Kest',
                         nsim = 99,
                         verbose = FALSE, global = TRUE)
plot(birch.kest.mc, shade = c("hi", "lo"))

summary(birch.kest.mc)
## Simultaneous critical envelopes for K(r)
## and observed value for 'urkiola.birch'
## Obtained from 99 simulations of CSR
## Alternative: two.sided
## Envelopes computed as theoretical curve plus/minus maximum simulated value of 
## maximum absolute deviation
## Significance level of Monte Carlo test: 1/100 = 0.01
## Data: urkiola.birch

The birch trees also indicate spatial clustering.

Finally, test the spatial dependence between the two distributions to test if the two species truly co-occur or not. Use the envelope() function again, but with the Kcross() function to estimate the cross-dependence. Plot the results and state whether the distributions of the two species show positive, negative or no correlation with each other

birchoak.kc <- envelope(urkiola.ppp, 
                       Kcross, 
                       i = "oak", 
                       j = "birch",
                       nsim = 99, 
                       verbose = FALSE,
                       global = TRUE)

plot(birchoak.kc)

The birch and the oak, by comparing the observed line with the theoretical line, look to be randomly distributed (the combined distribution is random)