Calling on Spatstat

library(spatstat)
## Loading required package: spatstat.data
## Loading required package: spatstat.univar
## spatstat.univar 3.1-4
## Loading required package: spatstat.geom
## spatstat.geom 3.5-0
## Loading required package: spatstat.random
## spatstat.random 3.4-1
## Loading required package: spatstat.explore
## Loading required package: nlme
## spatstat.explore 3.5-2
## Loading required package: spatstat.model
## Loading required package: rpart
## spatstat.model 3.4-0
## Loading required package: spatstat.linnet
## spatstat.linnet 3.3-1
## 
## spatstat 3.4-0 
## For an introduction to spatstat, type 'beginner'

Plots

First random point plot:

Other ways to write this code would have been:

random1 = rpoispp(500) plot (random1)

or

rpoispp(500) |> plot()

Point Density and Denisty Contours

pp <- rpoispp(500)
plot(density(pp))
contour(density(pp), add=T)

plot(density(pp))
plot(pp,add=T)

Point Density Plots

add=T adds this new data to the same image, rather than creating a new plot.

Quadrat Count

plot(quadratcount(pp,nx=10,ny=10))

The numbers in the squares represent the number of times an “event” occurs in them (the number of times a random point will be in that quadrat).

Quadrat Tests

quadrat.test(rpoispp(500),nx=20,ny=20)
## Warning: Some expected counts are small; chi^2 approximation may be inaccurate
## 
##  Chi-squared test of CSR using quadrat counts
## 
## data:  rpoispp(500)
## X2 = 401.18, df = 399, p-value = 0.92
## alternative hypothesis: two.sided
## 
## Quadrats: 20 by 20 grid of tiles
quadrat.test(rpoispp(50),nx=10,ny=10)
## Warning: Some expected counts are small; chi^2 approximation may be inaccurate
## 
##  Chi-squared test of CSR using quadrat counts
## 
## data:  rpoispp(50)
## X2 = 86.111, df = 99, p-value = 0.3622
## alternative hypothesis: two.sided
## 
## Quadrats: 10 by 10 grid of tiles
quadrat.test(rpoispp(10000),nx=20,ny=40)
## 
##  Chi-squared test of CSR using quadrat counts
## 
## data:  rpoispp(10000)
## X2 = 780.08, df = 799, p-value = 0.6452
## alternative hypothesis: two.sided
## 
## Quadrats: 20 by 40 grid of tiles
quadrat.test(rpoispp(500),nx=4,ny=4)
## 
##  Chi-squared test of CSR using quadrat counts
## 
## data:  rpoispp(500)
## X2 = 11.919, df = 15, p-value = 0.6297
## alternative hypothesis: two.sided
## 
## Quadrats: 4 by 4 grid of tiles
Table Comparing P-values
Number of Points Number of Cells Col3 Col4
500 20x20 0.362 It is likely that this pattern could be observed under conditions of CSR.
50 10x10 0.7622 It is even more likely that this pattern could be observed under conditions of CSR.
10000 20x40 0.4829 It is also likely this pattern could be observed under conditions of CSR.
500 4x4 0.4114 Similarly, it is likely this pattern could be observed under CSR conditions.

1st Order Trend

pp <- rpoispp(function(x,y) {200*x + 200*y})
quadrat.test(pp, nx=8, ny=8)
## Warning: Some expected counts are small; chi^2 approximation may be inaccurate
## 
##  Chi-squared test of CSR using quadrat counts
## 
## data:  pp
## X2 = 79.41, df = 63, p-value = 0.1587
## alternative hypothesis: two.sided
## 
## Quadrats: 8 by 8 grid of tiles
plot(density(pp))
plot(pp, pch=1, add=TRUE)

First p-value = 0.1601 Second p-value = 0.03271 Third p-value - 0.2255 All three of these are larger than 0.05 and therefore likely they could be observed under CSR conditions. However, the second p-value is lower, making it slightly less likely that this result would be observed under CSR conditions. None of these p-values provide enough evidence to reject the null hypothesis.

pp <- rpoispp(function(x,y) {200*x + 200*y})
quadrat.test(pp, nx=16, ny=16)
## Warning: Some expected counts are small; chi^2 approximation may be inaccurate
## 
##  Chi-squared test of CSR using quadrat counts
## 
## data:  pp
## X2 = 329.77, df = 255, p-value = 0.002196
## alternative hypothesis: two.sided
## 
## Quadrats: 16 by 16 grid of tiles
plot(density(pp))
plot(pp, pch=1, add=TRUE)

p-value = 0.7271

pp <- rpoispp(function(x,y) {200*x + 200*y})
quadrat.test(pp, nx=32, ny=32)
## Warning: Some expected counts are small; chi^2 approximation may be inaccurate
## 
##  Chi-squared test of CSR using quadrat counts
## 
## data:  pp
## X2 = 991.83, df = 1023, p-value = 0.4956
## alternative hypothesis: two.sided
## 
## Quadrats: 32 by 32 grid of tiles
plot(density(pp))
plot(pp, pch=1, add=TRUE)

p-value = 0.2561. I expected the larger the area, the great the p-value to be, though this is not exactly what happened. This fluctuation seemingly reflects the domed shape of significance, showing extremes on either end.

Lansing Tree Data

plot(split(lansing))

plot(split(lansing)$hickory)

summary(lansing)
## Marked planar point pattern:  2251 points
## Average intensity 2251 points per square unit (one unit = 924 feet)
## 
## *Pattern contains duplicated points*
## 
## Coordinates are given to 3 decimal places
## i.e. rounded to the nearest multiple of 0.001 units (one unit = 924 feet)
## 
## Multitype:
##          frequency proportion intensity
## blackoak       135 0.05997335       135
## hickory        703 0.31230560       703
## maple          514 0.22834300       514
## misc           105 0.04664594       105
## redoak         346 0.15370950       346
## whiteoak       448 0.19902270       448
## 
## Window: rectangle = [0, 1] x [0, 1] units
## Window area = 1 square unit
## Unit of length: 924 feet

The area in square feet of the study region is 924 square feet. The total number of trees is 2,251 trees. The most abundant tree is the hickory, and the least abundant tree is the black oak.

Different Species

plot(split(lansing)$hickory)

plot(split(lansing)$maple)

plot(split(lansing)$blackoak)

The distribution of hickory trees looks aggregated. The distribution of maple trees looks even more aggregated. There seems to be a higher concentration on the bottom half and going into the middle. I would maybe still expect the p-value to be less than 0.05, but perhaps a lower value. The distribution of black oak trees looks fairly random, perhaps a bit of clustering around the type.

Quadrat Tests for Trees

quadrat.test(split(lansing)$hickory, nx=20, ny=20)
## Warning: Some expected counts are small; chi^2 approximation may be inaccurate
## 
##  Chi-squared test of CSR using quadrat counts
## 
## data:  split(lansing)$hickory
## X2 = 766.7, df = 399, p-value < 2.2e-16
## alternative hypothesis: two.sided
## 
## Quadrats: 20 by 20 grid of tiles

p-value is 2.2e-16 because this is a small p-value, it provides greater evidence that we can reject a null hypothesis hickory trees are likely not random

quadrat.test(split(lansing)$maple, nx=20, ny=20)
## Warning: Some expected counts are small; chi^2 approximation may be inaccurate
## 
##  Chi-squared test of CSR using quadrat counts
## 
## data:  split(lansing)$maple
## X2 = 818.3, df = 399, p-value < 2.2e-16
## alternative hypothesis: two.sided
## 
## Quadrats: 20 by 20 grid of tiles

p-value is 2.2e-16, also extremely small, likely not random

quadrat.test(split(lansing)$blackoak, x=20, ny=20)
## Warning: : unrecognised argument 'x' was ignored
## Warning: Some expected counts are small; chi^2 approximation may be inaccurate
## 
##  Chi-squared test of CSR using quadrat counts
## 
## data:  split(lansing)$blackoak
## X2 = 258.33, df = 99, p-value = 7.226e-16
## alternative hypothesis: two.sided
## 
## Quadrats: 5 by 20 grid of tiles

p-value is 2.2e-16, all equally small For all three of the above commands, I got a message saying “chi^2 approximation may be inaccurate)

##Quadrat Variations for Hickory

quadrat.test(split(lansing)$hickory, nx=40, ny=40)
## Warning: Some expected counts are small; chi^2 approximation may be inaccurate
## 
##  Chi-squared test of CSR using quadrat counts
## 
## data:  split(lansing)$hickory
## X2 = 2103.3, df = 1599, p-value = 4.19e-16
## alternative hypothesis: two.sided
## 
## Quadrats: 40 by 40 grid of tiles

p-value = 4.19e-16

quadrat.test(split(lansing)$hickory, nx=100, ny=100)
## Warning: Some expected counts are small; chi^2 approximation may be inaccurate
## 
##  Chi-squared test of CSR using quadrat counts
## 
## data:  split(lansing)$hickory
## X2 = 10492, df = 9999, p-value = 0.0005939
## alternative hypothesis: two.sided
## 
## Quadrats: 100 by 100 grid of tiles

p-value = 0.0005939

quadrat.test(split(lansing)$hickory, nx=1000, ny=1000)
## Warning: Some expected counts are small; chi^2 approximation may be inaccurate
## 
##  Chi-squared test of CSR using quadrat counts
## 
## data:  split(lansing)$hickory
## X2 = 1002142, df = 1e+06, p-value = 0.1299
## alternative hypothesis: two.sided
## 
## Quadrats: 1000 by 1000 grid of tiles

p-value = 0.1299 As the quadrat numbers get bigger and bigger, the p-value gets smaller and smaller. This reveals that there may not be sufficient evidence to reject the null hypothesis (hickory trees may, in fact, be random)

Nearest Neighbor difference

plot(apply(nndist(rpoispp(500), k=1:500), 2, FUN = mean))

Alternative ways of achieving this: neighbor1 = rpoispp(500) plot(apply(nndist(neighbor1, k=1:500),2, FUN = mean))

or

rpoispp(500) |> nndist(k = 1:500) |> apply(2, mean) |> plot()

plot(apply(nndist(rpoispp(500), k=1:100), 2, FUN = mean),
     xlab = "Neighbor Order (k)",
     ylab = "Average Nearest Neighbor Distance",
     type = "b",
     pch = 19)
nn_distance <- nndist(split(lansing)$maple, k = 1:100)
ann_values <- colMeans(nn_distance)
points(1:length(ann_values), ann_values,
       pch = 5)

This plot is showing both pure random points (CSR) and the maple tree distribution. We can see the maple tree’s distance from it’s neighbors increases slower than pure random points. Because of this, we can say the maple trees are likely not totally random, and may be aggregated or concentrated.

plot(apply(nndist(rpoispp(42), k=1:42), 2, FUN = mean),
     type = "b",
     pch = 19)
nn_distance <- nndist(cells, k=1:42)
ann_values <- colMeans(nn_distance)
points(1:length(ann_values), ann_values,
       pch = 5)

I removed the x and y labels because the plot would not fit otherwise. This plot suggests that the cell data is even more aggregated, because the neighbors are all very close to each other, much closer than the random points predict.

G env

plot(Gest(rpoispp(500)))

G_env <- envelope(rpoispp(500), Gest, nsim = 95, alpha  = 0.05)
## Generating 95 simulations of CSR  ...
## 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20,
## 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40,
## 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60,
## 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80,
## 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 
## 95.
## 
## Done.
plot(G_env)

##Hickory G en

Hickory_G_env <- envelope(split(lansing)$hickory, Gest, nsim = 95, alpha  = 0.05)
## Generating 95 simulations of CSR  ...
## 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20,
## 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40,
## 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60,
## 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80,
## 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 
## 95.
## 
## Done.
plot(Hickory_G_env)

The data wanders above the null expectation (though not by much), suggesting that the hickory trees are clumped.

##Black oak G env

Blackoak_G_env <- envelope(split(lansing)$blackoak, Gest, nsim = 95, alpha  = 0.05)
## Generating 95 simulations of CSR  ...
## 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20,
## 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40,
## 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60,
## 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80,
## 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 
## 95.
## 
## Done.
plot(Blackoak_G_env)

Similarly, the blackoak data wanders above the random points, also indicating clumping.

##Cells G env

Cells_G_env <- envelope(cells, Gest, nsim = 95, alpha  = 0.05)
## Generating 95 simulations of CSR  ...
## 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20,
## 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40,
## 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60,
## 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80,
## 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 
## 95.
## 
## Done.
plot(Cells_G_env)

This shows the cell data veering below the random line, which would indicate more even spacing because the distances are decreasing slower.

##Longleaf G env

data(longleaf)
Longleaf_G_env <- envelope(longleaf, Gest, nsim = 95, alpha  = 0.05)
## Generating 95 simulations of CSR  ...
## 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20,
## 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40,
## 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60,
## 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80,
## 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 
## 95.
## 
## Done.
plot(Longleaf_G_env)

The data wanders above the confidence band, suggesting clumping.