R/globaldiag.R

In stopp: Spatio-Temporal Point Pattern Methods, Model Fitting, Diagnostics, Simulation, Local Tests

#' Global diagnostics of a spatio-temporal point process first-order intensity
#'
#' @description
#' This function performs global diagnostics of a model fitted for the
#' first-order intensity of a  spatio-temporal point pattern, by returning
#' the inhomogeneous K-function weighted by the
#' provided intensity to diagnose, its theoretical value,
#' and their difference.
#'
#' @details
#' If applied to a \code{stp} object, it resorts to  the
#' spatio-temporal inhomogeneous K-function (Gabriel and Diggle, 2009)
#'  documented by the function
#'  \link{STIKhat}  of the \code{stpp} package (Gabriel et al, 2013).
#'
#'  If applied to a \code{stlp} object, it uses the
#' spatio-temporal inhomogeneous K-function on a linear network (Moradi and Mateu, 2020)
#'  documented by the function
#'  \link{STLKinhom}  of the \code{stlnpp} package (Moradi et al., 2020).
#'
#'
#' @param x A \code{stp} object
#' @param intensity A vector of intensity values, of the same length as the number
#' of point in \code{x}
#'
#' @return A list of class \code{globaldiag}, containing
#' \describe{
#' \item{\code{x}}{The observed point pattern}
#' \item{\code{dist}}{The spatial ranges of the K-function}
#' \item{\code{times}}{The temporal ranges of the K-function}
#' \item{\code{est}}{The estimated K-function weighted by the intensity function in input}
#' \item{\code{theo}}{The theoretical K-function}
#' \item{\code{diffK}}{The difference between the estimated and the theoretical K-functions}
#' \item{\code{squared.diff}}{The sum of the squared differences between the estimated and the theoretical K-functions}
#' }
#'
#'
#'
#' @export
#'
#' @author Nicoletta D'Angelo and Giada Adelfio
#'
#' @seealso
#' \link{plot.globaldiag}, \link{print.globaldiag},
#' \link{summary.globaldiag}
#'
#'
#'
#'
#' @examples
#'
#' set.seed(2)
#' inh <- rstpp(lambda = function(x, y, t, a) {exp(a[1] + a[2]*x)}, par = c(.3, 6))
#' 
#' mod1 <- stppm(inh, formula = ~ 1)
#' mod2 <- stppm(inh, formula = ~ x)
#' 
#' g1 <- globaldiag(inh, mod1$l)
#' g2 <- globaldiag(inh, mod2$l)
#' 
#'
#'
#'
#' @references
#' Adelfio, G., Siino, M., Mateu, J., and Rodríguez-Cortés, F. J. (2020). Some properties of local weighted second-order statistics for spatio-temporal point processes. Stochastic Environmental Research and Risk Assessment, 34(1), 149-168.
#'
#' D’Angelo, N., Adelfio, G.  and Mateu, J. (2022) Local inhomogeneous second-order characteristics for spatio-temporal point processes on linear networks. Stat Papers. https://doi.org/10.1007/s00362-022-01338-4
#'
#' Gabriel, E., and Diggle, P. J. (2009). Second‐order analysis of inhomogeneous spatio‐temporal point process data. Statistica Neerlandica, 63(1), 43-51.
#'
#' Gabriel, E., Rowlingson, B. S., & Diggle, P. J. (2013). stpp: An R Package for Plotting, Simulating and Analyzing Spatio-Temporal Point Patterns. Journal of Statistical Software, 53(2), 1–29. https://doi.org/10.18637/jss.v053.i02
#'
#' Moradi M, Cronie O, and Mateu J (2020). stlnpp: Spatio-temporal analysis of point patterns on linear networks.
#'
#' Moradi, M. M., and Mateu, J. (2020). First-and second-order characteristics of spatio-temporal point processes on linear networks. Journal of Computational and Graphical Statistics, 29(3), 432-443.
#'
#'
globaldiag <- function(x, intensity){
  
  if (!inherits(x, c("stp", "stlp"))) stop("x should be either from class stp or stlp")
  
  nx <- nrow(x$df)
  
  if (!is.numeric(intensity)) {
    stop("intensity should be a numeric vector")
  } else {
    if(length(intensity) != nx) {
      stop("intensity should be the same length as the number of points in x")
    }
  } 

  if(any(class(x) == "stp")) {
    Kfunct <- stpp::STIKhat(xyt = as.stpp(x), lambda = intensity)

    est <- Kfunct$Khat
    theo <- Kfunct$Ktheo

    dist <- Kfunct$dist
    times <- Kfunct$times
  } else if(any(class(x) == "stlp")){
    Kfunct <- stlnpp::STLKinhom(x = as.stlpp(x), lambda = intensity)

    est <- Kfunct$Kinhom
    theo <- Kfunct$Ktheo

    dist <- Kfunct$r
    times <- Kfunct$t
  }

  diffK <- est - theo
  squared.diff <- sum(diffK ^ 2 / (est))
  
  out <- list(x = x, dist = dist, times = times, est = est, theo = theo,
              diffK = diffK, squared.diff = squared.diff)
  
  class(out) <- "globaldiag"
  
  return(out)

}