`variogramST()` error when using tibble

[vincenzocoia]

Summary of the problem

The gamma (sample variogram) estimate output by variogramST() is NA when the data input is an STIDF() object that itself has a data object entered as a tibble -- a special type of data frame that is commonly used these days with the popularization of the tidyverse. The documentation of STIDF()'s data argument says that it should be a data frame -- which a tibble is.

Summary of the culprit

The reason behind the NAs is because of the line

gamma[j,i] <- 0.5*mean((data[,,colnames(m)[1]]@data[indSp[,1],1] - data[,,colnames(m)[1]]@data[indSp[,1]+indSp[,2],1])^2,
                             na.rm=TRUE)

in the variogramST() function. It's subsetting the data@data object using [<vector>, 1], which returns a vector when data@data is a plain data frame, and remains a tibble when data@data is a tibble.

Summary of suggested fixes

• The simplest fix is to update the documentation of STIDF()'s data argument, specifying that it has to be a regular data frame, not a tibble data frame.
• Better would be to add a line to STIDF(): data <- as.data.frame(data) (this strips away the special tibble class, if present).
• Best would be to add , drop = TRUE to the [] indexing.

Reproducible Example of the Computation Error

Load dependencies:

library(tidyverse)
library(sp)
#> The legacy packages maptools, rgdal, and rgeos, underpinning the sp package,
#> which was just loaded, will retire in October 2023.
#> Please refer to R-spatial evolution reports for details, especially
#> https://r-spatial.org/r/2023/05/15/evolution4.html.
#> It may be desirable to make the sf package available;
#> package maintainers should consider adding sf to Suggests:.
#> The sp package is now running under evolution status 2
#>      (status 2 uses the sf package in place of rgdal)
library(gstat)
library(spacetime)
library(lubridate)
#> 
#> Attaching package: 'lubridate'
#> The following objects are masked from 'package:base':
#> 
#>     date, intersect, setdiff, union

Generate temperature data at different (x, y) locations across 4 days

n <- 20
set.seed(2)
full <- tibble(x = runif(n), 
               y = runif(n), 
               time = ymd("2023-07-07") + 
                 ddays(sample(1:4, size = n, replace = TRUE)),
               temperature = rnorm(n))

ggplot(full, aes(x, y, colour = temperature)) +
  facet_wrap(~ time) +
  geom_point() +
  theme_bw()

Turn that data into (1) a spatial dataset, and (2) a dataset of the dependent variable (temperature)

## Make spatial data
xy <- select(full, x, y)
coordinates(xy) <- ~ x + y

## Make dataset of dependent variable (temperature).
## NOTE that class is a tibble -- a special type of data frame
data <- select(full, temperature)
class(data)
#> [1] "tbl_df"     "tbl"        "data.frame"

Make spatiotemporal (STIDF()) data and corresponding sample variograms:

1. Using data as a tibble.
2. Using data as a data frame only.

stdat1 <- STIDF(xy, full$time, data = data)
stdat2 <- STIDF(xy, full$time, data = as.data.frame(data))

vv1 <- variogramST(temperature ~ 1, data = stdat1, tunit = "days", tlags = 0:4, 
                   progress = FALSE) 
#> Warning in mean.default((data[, , colnames(m)[1]]@data[indSp[, 1], 1] - :
#> argument is not numeric or logical: returning NA

(snip 40 warnings like this)

vv2 <- variogramST(temperature ~ 1, data = stdat2, tunit = "days", tlags = 0:4,
                   progress = FALSE) 

Sample variogram is NA when data is a tibble:

vv1$gamma
#>  [1] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
#> [26] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
#> [51] NA NA NA NA NA NA NA NA NA NA
vv2$gamma
#>  [1]         NA         NA         NA         NA 0.01638317 0.30720071
#>  [7]         NA 0.41096178 0.50923723 0.25442567 1.46648330 0.04606630
#> [13] 0.41838934 0.23442360 0.05686060         NA 0.13897625 1.40272114
#> [19] 0.15111853 0.37774027 0.90697259 1.83911154 4.01188360 2.63673603
#> [25] 4.18618535 1.55124883 1.24112682 0.21521518 0.40065672 0.65988774
#> [31]         NA 1.03942542         NA         NA 1.74519573 3.34520046
#> [37] 0.61374727 0.79421057 0.42605686 0.70225243 1.99483642 1.10960392
#> [43] 0.05751237 1.83829579 0.89214221         NA         NA         NA
#> [49]         NA         NA 0.01868505         NA         NA         NA
#> [55]         NA         NA 0.21618180 0.03618512         NA 0.54229038

^{Created on 2023-07-07 with reprex v2.0.2}