Lab 5 RMD
Lindsey Meier
2025-10-11
Question 1
Row 2 columns 2 and 3 have the lowest density of trees, with counts of 49 and 17 respectively.
Q <- quadratcount(bei, nx=6, ny=3)
plot(bei, cex =0.5, pch="+")
plot(Q, add=TRUE, cex=2)
Question 2
The highest point density contour is 0.02.
den <- density(bei)
plot(den)
plot(bei, add=T, cex=0.5)
persp(den)
contour(den)
##Question 3 The lowest density value in the resulting tessellation is 0. The is evidence to reject the null hypothesis of complete spatial randomness, because the p-value is 2.2e-16, very low and indicative of some pattern (clustering) rather than CSR.
mytess <- hextess(bei, 50)
Q <- quadratcount(bei, tess=mytess)
plot(bei, cex=0.5, pch="+")
plot(Q, add=TRUE, cex=1)
mytess <- hextess(bei, 10)
Q <- quadratcount(bei, tess=mytess)
plot(bei, cex=0.5, pch="+")
plot(Q, add=TRUE, cex=1)
mytest <- quadrat.test(bei, tess=mytess)## Warning: Some expected counts are small; chi^2 approximation may be inaccurateQuestion 4 & 5
The breaks for the slope are 1-5, instead of 1-4, since this time we did quintiles.
slope <- bei.extra$grad
elev <- bei.extra$elev
b <- quantile(slope, probs =
(0:4)/4)
slopecut <- cut(slope,
breaks=b, labels=1:4)
V <- tess(image = slopecut)
plot(V)
plot(bei, add=T, pch="+")
quadcount <-
quadratcount(bei, tess=V)
plot(quadcount)
quadrat.test(bei, tess=V)##
## Chi-squared test of CSR using quadrat counts
##
## data: bei
## X2 = 661.84, df = 3, p-value < 2.2e-16
## alternative hypothesis: two.sided
##
## Quadrats: 4 tiles (levels of a pixel image)b <- quantile(slope, probs =
(0:5)/5)
slopecut <- cut(slope,
breaks=b, labels=1:5)
V <- tess(image = slopecut)
plot(V)
plot(bei, add=T, pch="+")
quadcount <-
quadratcount(bei, tess=V)
plot(quadcount)
quadrat.test(bei, tess=V)##
## Chi-squared test of CSR using quadrat counts
##
## data: bei
## X2 = 587, df = 4, p-value < 2.2e-16
## alternative hypothesis: two.sided
##
## Quadrats: 5 tiles (levels of a pixel image)Question 6
1:490, 2:719, 3:861, 4:1174, 5:360
e <- quantile(elev, probs = (0:5)/5)
elevcut <- cut(elev, breaks=e, labels=1:5)
x <- tess(image = elevcut)
plot(x)
plot(bei, add=T, pch="+")
quadcount1 <-
quadratcount(bei, tess=x)
plot(quadcount1)
quadrat.test(bei, tess=x)##
## Chi-squared test of CSR using quadrat counts
##
## data: bei
## X2 = 567.09, df = 4, p-value < 2.2e-16
## alternative hypothesis: two.sided
##
## Quadrats: 5 tiles (levels of a pixel image)quadcount1## tile
## 1 2 3 4 5
## 490 719 861 1174 360Question 7
Both p-values are 2.2e-16, so both slope and elevation could be a contributing cause of the non-random spatial distirbution of the points in the bei dataset.
quadrat.test(bei, tess=x)##
## Chi-squared test of CSR using quadrat counts
##
## data: bei
## X2 = 567.09, df = 4, p-value < 2.2e-16
## alternative hypothesis: two.sided
##
## Quadrats: 5 tiles (levels of a pixel image)quadrat.test(bei, tess=V)##
## Chi-squared test of CSR using quadrat counts
##
## data: bei
## X2 = 587, df = 4, p-value < 2.2e-16
## alternative hypothesis: two.sided
##
## Quadrats: 5 tiles (levels of a pixel image)Question 8 & 9
The dashed line in the plot is used as a reference line, and represents the observed intensity of something occurring under the null hypothesis, and having a line above or below the line will mean the observed intensity is higher or lower. The plot suggests that these trees prefer slope, and their presence along steep and steeper slopes grows, until at a point the graph sharply drops off.
plot(rhohat(bei, slope))
Question 10
Slope and slope:elevation appear to be significant in this model; slope on its own and slope with elevation.
model <- ppm(bei ~ slope * elev)
summary(model)## Point process model
## Fitted to data: bei
## Fitting method: maximum likelihood (Berman-Turner approximation)
## Model was fitted using glm()
## Algorithm converged
## Call:
## ppm.formula(Q = bei ~ slope * elev)
## Edge correction: "border"
## [border correction distance r = 0 ]
## --------------------------------------------------------------------------------
## Quadrature scheme (Berman-Turner) = data + dummy + weights
##
## Data pattern:
## Planar point pattern: 3604 points
## Average intensity 0.00721 points per square metre
## Window: rectangle = [0, 1000] x [0, 500] metres
## Window area = 5e+05 square metres
## Unit of length: 1 metre
##
## Dummy quadrature points:
## 130 x 130 grid of dummy points, plus 4 corner points
## dummy spacing: 7.692308 x 3.846154 metres
##
## Original dummy parameters: =
## Planar point pattern: 16904 points
## Average intensity 0.0338 points per square metre
## Window: rectangle = [0, 1000] x [0, 500] metres
## Window area = 5e+05 square metres
## Unit of length: 1 metre
## Quadrature weights:
## (counting weights based on 130 x 130 array of rectangular tiles)
## All weights:
## range: [1.64, 29.6] total: 5e+05
## Weights on data points:
## range: [1.64, 14.8] total: 41000
## Weights on dummy points:
## range: [1.64, 29.6] total: 459000
## --------------------------------------------------------------------------------
## FITTED :
##
## Nonstationary Poisson process
##
## ---- Intensity: ----
##
## Log intensity: ~slope * elev
## Model depends on external covariates 'slope' and 'elev'
## Covariates provided:
## slope: im
## elev: im
##
## Fitted trend coefficients:
## (Intercept) slope elev slope:elev
## -4.403426761 -36.500458701 -0.007024363 0.292813497
##
## Estimate S.E. CI95.lo CI95.hi Ztest
## (Intercept) -4.403426761 0.614763210 -5.60834051 -3.198513010 ***
## slope -36.500458701 5.261456400 -46.81272375 -26.188193651 ***
## elev -0.007024363 0.004192219 -0.01524096 0.001192235
## slope:elev 0.292813497 0.036257125 0.22175084 0.363876155 ***
## Zval
## (Intercept) -7.162801
## slope -6.937330
## elev -1.675572
## slope:elev 8.076026
##
## ----------- gory details -----
##
## Fitted regular parameters (theta):
## (Intercept) slope elev slope:elev
## -4.403426761 -36.500458701 -0.007024363 0.292813497
##
## Fitted exp(theta):
## (Intercept) slope elev slope:elev
## 1.223534e-02 1.406217e-16 9.930002e-01 1.340193e+00plot(model)
Question 11
The results, the coefficients, are not the same, but very similar. The significance tests are still telling us the similar story as above. Slope is a good predictor or tree density in this dataset and model, slope with elevation is as well.
model <- slrm(bei ~ slope)
ppmmodel <- ppm(bei ~ slope * elev)
slrmmodel <- slrm(bei ~ slope * elev)
ppmmodel## Nonstationary Poisson process
## Fitted to point pattern dataset 'bei'
##
## Log intensity: ~slope * elev
##
## Fitted trend coefficients:
## (Intercept) slope elev slope:elev
## -4.403426761 -36.500458701 -0.007024363 0.292813497
##
## Estimate S.E. CI95.lo CI95.hi Ztest
## (Intercept) -4.403426761 0.614763210 -5.60834051 -3.198513010 ***
## slope -36.500458701 5.261456400 -46.81272375 -26.188193651 ***
## elev -0.007024363 0.004192219 -0.01524096 0.001192235
## slope:elev 0.292813497 0.036257125 0.22175084 0.363876155 ***
## Zval
## (Intercept) -7.162801
## slope -6.937330
## elev -1.675572
## slope:elev 8.076026slrmmodel## Fitted spatial logistic regression model
## Formula: bei ~ slope * elev
## Fitted coefficients:
## (Intercept) slope elev slope:elev
## -4.330416600 -36.301093115 -0.009931764 0.303478908plot(ppmmodel)
plot(slrmmodel)
summary(slrmmodel)## Fitted spatial logistic regression model
## Call: [1] "slrm(bei ~ slope * elev)"
## Formula: bei ~ slope * elev
## Fitted coefficients:
## Estimate S.E. CI95.lo CI95.hi Ztest
## (Intercept) -4.330416600 0.766832930 -5.83338152 -2.827452e+00 ***
## slope -36.301093115 6.687911473 -49.40915873 -2.319303e+01 ***
## elev -0.009931764 0.005245837 -0.02021342 3.498863e-04
## slope:elev 0.303478908 0.046448056 0.21244239 3.945154e-01 ***
## Zval
## (Intercept) -5.647145
## slope -5.427867
## elev -1.893266
## slope:elev 6.533727summary(ppmmodel)## Point process model
## Fitted to data: bei
## Fitting method: maximum likelihood (Berman-Turner approximation)
## Model was fitted using glm()
## Algorithm converged
## Call:
## ppm.formula(Q = bei ~ slope * elev)
## Edge correction: "border"
## [border correction distance r = 0 ]
## --------------------------------------------------------------------------------
## Quadrature scheme (Berman-Turner) = data + dummy + weights
##
## Data pattern:
## Planar point pattern: 3604 points
## Average intensity 0.00721 points per square metre
## Window: rectangle = [0, 1000] x [0, 500] metres
## Window area = 5e+05 square metres
## Unit of length: 1 metre
##
## Dummy quadrature points:
## 130 x 130 grid of dummy points, plus 4 corner points
## dummy spacing: 7.692308 x 3.846154 metres
##
## Original dummy parameters: =
## Planar point pattern: 16904 points
## Average intensity 0.0338 points per square metre
## Window: rectangle = [0, 1000] x [0, 500] metres
## Window area = 5e+05 square metres
## Unit of length: 1 metre
## Quadrature weights:
## (counting weights based on 130 x 130 array of rectangular tiles)
## All weights:
## range: [1.64, 29.6] total: 5e+05
## Weights on data points:
## range: [1.64, 14.8] total: 41000
## Weights on dummy points:
## range: [1.64, 29.6] total: 459000
## --------------------------------------------------------------------------------
## FITTED :
##
## Nonstationary Poisson process
##
## ---- Intensity: ----
##
## Log intensity: ~slope * elev
## Model depends on external covariates 'slope' and 'elev'
## Covariates provided:
## slope: im
## elev: im
##
## Fitted trend coefficients:
## (Intercept) slope elev slope:elev
## -4.403426761 -36.500458701 -0.007024363 0.292813497
##
## Estimate S.E. CI95.lo CI95.hi Ztest
## (Intercept) -4.403426761 0.614763210 -5.60834051 -3.198513010 ***
## slope -36.500458701 5.261456400 -46.81272375 -26.188193651 ***
## elev -0.007024363 0.004192219 -0.01524096 0.001192235
## slope:elev 0.292813497 0.036257125 0.22175084 0.363876155 ***
## Zval
## (Intercept) -7.162801
## slope -6.937330
## elev -1.675572
## slope:elev 8.076026
##
## ----------- gory details -----
##
## Fitted regular parameters (theta):
## (Intercept) slope elev slope:elev
## -4.403426761 -36.500458701 -0.007024363 0.292813497
##
## Fitted exp(theta):
## (Intercept) slope elev slope:elev
## 1.223534e-02 1.406217e-16 9.930002e-01 1.340193e+00Question 12
There are 150 sites in this reduced dataset.
gor_subset <-
gorillas[marks(gorillas)$group
=="minor"&
marks(gorillas)
$season == "rainy",]
gor_subsetdry <-
gorillas[marks(gorillas)$group
=="major"&
marks(gorillas)
$season=="dry",]
gor_subsetdry## Marked planar point pattern: 150 points
## Mark variables: group, season, date
## window: polygonal boundary
## enclosing rectangle: [580457.9, 585934] x [674172.8, 678739.2] metres##Question 13 There is evidence for a non-random distribution across vegetation types for the dry season, primary vegetation type has the highest number of nests.
veg_tess <- tess(image =
gorillas.extra$vegetation)
quadratcount(gor_subset,
tess=veg_tess)## tile
## Disturbed Colonising Grassland Primary Secondary Transition
## 18 1 5 140 5 3quadrat.test(gor_subset,
tess=veg_tess)## Warning: Some expected counts are small; chi^2 approximation may be inaccurate##
## Chi-squared test of CSR using quadrat counts
##
## data: gor_subset
## X2 = 225.47, df = 5, p-value < 2.2e-16
## alternative hypothesis: two.sided
##
## Quadrats: 6 tiles (levels of a pixel image)veg_tess <-(image=
gorillas.extra$vegetation)
quadratcount(gor_subsetdry,
tess=veg_tess)## tile
## Disturbed Colonising Grassland Primary Secondary Transition
## 24 0 2 120 2 2quadrat.test(gor_subsetdry,
tess=veg_tess)## Warning: Some expected counts are small; chi^2 approximation may be inaccurate##
## Chi-squared test of CSR using quadrat counts
##
## data: gor_subsetdry
## X2 = 183.29, df = 5, p-value < 2.2e-16
## alternative hypothesis: two.sided
##
## Quadrats: 6 tiles (levels of a pixel image)Question 14
There is evidence for a non-random distribution across slope type for the dry season, there are 58 nests on ridges, although, valley is not far behind with 55 nests.
veg_tess <- tess(image=
gorillas.extra$slopetype)
quadratcount(gor_subsetdry,
tess=veg_tess)## tile
## Valley Toe Flat Midslope Upper Ridge
## 55 12 0 21 4 58quadrat.test(gor_subsetdry,
tess=veg_tess)## Warning: Some expected counts are small; chi^2 approximation may be inaccurate##
## Chi-squared test of CSR using quadrat counts
##
## data: gor_subsetdry
## X2 = 5.0863, df = 5, p-value = 0.8109
## alternative hypothesis: two.sided
##
## Quadrats: 6 tiles (levels of a pixel image)Question 15
The results strongly reject the null hypothesis of CSR. This means the gorilla nest sites are not distributed randomly or uniformly with respect to the categories of slope and elevation; they show a significant preference or avoidance for certain ranges of these two environmental variables.
#slope
quadrat.test(gor_subsetdry,
tess=V)##
## Chi-squared test of CSR using quadrat counts
##
## data: gor_subsetdry
## X2 = 600, df = 4, p-value < 2.2e-16
## alternative hypothesis: two.sided
##
## Quadrats: 5 tiles (levels of a pixel image)#Elev
quadrat.test(gor_subsetdry,
tess=x)##
## Chi-squared test of CSR using quadrat counts
##
## data: gor_subsetdry
## X2 = 599.63, df = 4, p-value < 2.2e-16
## alternative hypothesis: two.sided
##
## Quadrats: 5 tiles (levels of a pixel image)Question 16
The model predicts that the likelihood of a dry-season nest site being present increases with steeper slopes and higher elevations. The resulting prediction map would show the highest probability of nest sites in areas that combine these two characteristics.
slopegor<-
gorillas.extra$slopeangle
elevgor <-
gorillas.extra$elevation
slrmm <- slrm(gor_subsetdry
~ slopegor*elevgor)
slrmm## Fitted spatial logistic regression model
## Formula: gor_subsetdry ~ slopegor * elevgor
## Fitted coefficients:
## (Intercept) slopegor elevgor slopegor:elevgor
## -2.880176e+01 2.372978e-01 9.322878e-03 -1.236946e-04plot(slrmm)