Here we will use the vista package to plot temporal residuals from spatiotemporal models
library("knitr") opts_chunk$set(message = FALSE, fig.width = 5.5)
library(vista) library(mgcv) data(pcod)
We'll use a GAM for demonstration purposes, but any other package randomForest
, sdmTMB
, etc. could be used instead.
m <- gam(density ~ 0 + as.factor(year) + s(X,Y,by=as.factor(year)), data=pcod, family=tw()) pcod$pred = predict(m) pcod$resid = residuals(m)
The pred_time
function displays the range of predictions across time slices.
pred_time(df = pcod, time = "year")
Similarly, the resid_time
function is designed to plot residuals across time slices
resid_time(df = pcod, time = "year")
In some cases, the variability of a spatial process may be non-stationary and increasing or decreasing over time. We can visualize this using the
sd_resid_time
function, where ideally there is no trend or outlier.
sd_resid_time(df = pcod, time = "year")
The qq
function is also useful in displaying the quantile residuals (via a call to stats::qqnorm
). By default this will be faceted by time, but faceting can be turned off with by_time = FALSE
.
qq(pcod, time = "year")
The scale_loc
function displays scale - location plots, which may be useful at diagnosing trends in residuals vs predicted values,
scale_loc(df = pcod, time="year")
We can calculate lots of variants of Moran's I. The moran_ts
plot is useful for displaying time series of the Moran-I statistic, calculated separately for each time slice. This can be done on predicted values,
moran_ts(df = pcod, time = "year", response = "pred")
But if a goal is evaluating whether a particular model removes spatial autocorrelation, it may be more useful to calculate this on the residuals,
moran_ts(df = pcod, time = "year", response = "resid")
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.