This session covers tools for investigating the dependence of the point process intensity on a covariate. They include nonparametric methods and parametric models. The lecturer’s R script is available here (right click and save).
The bei
dataset gives the locations of trees in a survey area with additional covariate information in a list bei.extra
.
Assign the elevation covariate to a variable elev
by typing
elev <- bei.extra$elev
Assume that the intensity of trees is a function \(\lambda(u) = \rho(e(u))\) where \(e(u)\) is the terrain elevation at location u. Compute a nonparametric estimate of the function \(\rho\) and plot it by
rh <- rhohat(bei, elev)
plot(rh)
Compute the predicted intensity based on this estimate of \(\rho\).
Compute a non-parametric estimate of intensity by kernel smoothing and compare with the predicted intensity above.
Continuing with the dataset bei
conduct both Berman’s Z1 and Z2 tests for dependence on elev
, and plot the results.
Returning to the Japanese Pines data,
Fit the uniform Poisson point process model to the Japanese Pines data
ppm(japanesepines~1)
Read off the fitted intensity. Check that this is the correct value of the maximum likelihood estimate of the intensity.
The japanesepines
dataset is believed to exhibit spatial inhomogeneity.
Plot a kernel smoothed intensity estimate.
Fit the Poisson point process models with loglinear intensity (trend formula ~x+y
) and log-quadratic intensity (trend formula ~polynom(x,y,2)
) to the Japanese Pines data.
extract the fitted coefficients for these models using coef
.
Plot the fitted model intensity (using plot(fit)
)
perform the Likelihood Ratio Test for the null hypothesis of a loglinear intensity against the alternative of a log-quadratic intensity, using anova
.
Generate 10 simulated realisations of the fitted log-quadratic model, and plot them, using plot(simulate(fit, nsim=10))
where fit
is the fitted model.
The update
command can be used to re-fit a point process model using a different model formula.
Type the following commands and interpret the results:
fit0 <- ppm(japanesepines ~ 1)
fit1 <- update(fit0, . ~ x)
fit1
fit2 <- update(fit1, . ~ . + y)
fit2
Now type step(fit2)
and interpret the results.
The bei
dataset gives the locations of trees in a survey area with additional covariate information in a list bei.extra
.
Assign both terrain elevation and slope (gradient) sensible names
elev <- bei.extra$elev
grad <- bei.extra$grad
Fit a Poisson point process model to the data which assumes that the intensity is a loglinear function of terrain slope and elevation
Read off the fitted coefficients and write down the fitted intensity function.
Plot the fitted intensity as a colour image.
extract the estimated variance-covariance matrix of the coefficient estimates, using vcov
.
Compute and plot the standard error of the intensity estimate (see help(predict.ppm)
).
Fit Poisson point process models to the Japanese Pines data, with the following trend formulas. Read off an expression for the fitted intensity function in each case.
Trend formula | Fitted intensity function |
---|---|
~1 |
|
~x |
|
~sin(x) |
|
~x+y |
|
~polynom(x,y,2) |
|
~factor(x < 0.4) |
Make image plots of the fitted intensities for the inhomogeneous models above.