R/ggplot.R

In fbachl/inlabru: Bayesian Latent Gaussian Modelling using INLA and Extensions

#' @describeIn inlabru-deprecated
#' Plot a map using extent of a spatial object
#'
#' `r lifecycle::badge("deprecated")` This function is deprecated as `ggmap`
#' isn't supported.
#'
#' Used `ggmap::get_map()` to query map services like Google Maps for a region
#' centered around the spatial object provided. Then calls `ggmap()` to plot the
#' map.
#'
#' This function required the `ggmap` package.
#' @export
gmap <- function(...) {
  lifecycle::deprecate_stop(
    when = "2.11.1.9002",
    what = "gmap()",
    details = "ggmap isn't supported anymore."
  )
  # data <- fm_transform(data, crs = fm_crs("longlat_globe"))
  # df <- cbind(sp::coordinates(data), data@data)
  #
  # # Figure out a sensible bounding box (range of data plus 30%)
  # extend <- 2.5
  # crange <- apply(sp::coordinates(data), MARGIN = 2, range)
  # lonlim <- (extend * (crange[, 1] - mean(crange[, 1]))) + mean(crange[, 1])
  # latlim <- (extend * (crange[, 2] - mean(crange[, 2]))) + mean(crange[, 2])
  #
  # # Create map
  # requireNamespace("ggmap")
  # myMap <-
  #   ggmap::get_map(
  #     c(lonlim[1], latlim[1], lonlim[2], latlim[2]),
  #     ...
  #   )
  #
  # # Return map
  # ggmap::ggmap(myMap)
}

#' ggplot2 geomes for inlabru related objects
#'
#' @description
#' gg is a generic function for generating geomes from various kinds of spatial
#' objects, e.g. Spatial* data, meshes, Raster objects and inla/inlabru
#' predictions. The function invokes particular methods which depend on the
#' [class] of the first argument.
#'
#' @name gg
#' @export
#' @param data an object for which to generate a geom.
#' @param ... Arguments passed on to the geom method.
#' @return The form of the value returned by gg depends on the class of its
#'   argument. See the documentation of the particular methods for details of
#'   what is produced by that method.
#'
#' @family geomes for inla and inlabru predictions
#' @family geomes for spatial data
#' @family geomes for meshes
#' @family geomes for Raster data
#'
#' @examples
#' if (require("ggplot2", quietly = TRUE)) {
#'   # Load Gorilla data
#'
#'   gorillas <- inlabru::gorillas_sf
#'
#'   # Invoke ggplot and add geomes for the Gorilla nests and the survey
#'   # boundary
#'
#'   ggplot() +
#'     gg(gorillas$boundary) +
#'     gg(gorillas$nests)
#' }
gg <- function(data, ...) {
  UseMethod("gg")
}


#' @describeIn inlabru-deprecated
#' `r lifecycle::badge("deprecated")` This function is deprecated as `ggmap`
#' isn't supported.
#'
#' ggplot geom for spatial data
#'
#' gm is a wrapper for the [gg] method. It will take the first argument and
#' transform its coordinate system to latitude and longitude. Thereafter, [gg]
#' is called using the transformed data and the arguments provided via `...`. gm
#' is intended to replace gg whenever the data is supposed to be plotted over a
#' spatial map generated by [gmap], which only works if the coordinate system is
#' latitude/longitude.
gm <- function(...) {
  lifecycle::deprecate_warn(
    when = "2.11.1.9002",
    what = "gm()",
    with = "gg()",
    details = "ggmap isn't supported anymore."
  )
  gg(..., crs = fm_CRS("+proj=longlat"))
}


#' Geom for matrix
#'
#' Creates a tile geom for plotting a matrix
#'
#' Requires the `ggplot2` package.
#'
#' @name gg.matrix
#' @export
#' @param data A `matrix` object.
#' @param mapping a set of aesthetic mappings created by `aes`. These are passed
#'   on to `geom_tile`.
#' @param ... Arguments passed on to `geom_tile`.
#' @return A `geom_tile` with reversed y scale.
#' @family geomes for inla and inlabru predictions
#'
#' @importFrom rlang .data
#' @examples
#' if (require("ggplot2", quietly = TRUE)) {
#'   A <- matrix(runif(100), nrow = 10)
#'   ggplot() +
#'     gg(A)
#' }
gg.matrix <- function(data, mapping = NULL, ...) {
  requireNamespace("ggplot2")
  A <- as.matrix(data)
  grd <- expand.grid(row = seq_len(nrow(A)), column = seq_len(ncol(A)))
  df <- data.frame(value = as.vector(A), grd)

  dmap <- ggplot2::aes(x = .data$column, y = .data$row, fill = .data$value)

  if (!is.null(mapping)) {
    dmap <- modifyList(dmap, mapping)
  }

  ggp <- c(
    ggplot2::geom_tile(dmap, data = df, ...),
    ggplot2::scale_y_reverse()
  )
}


#' Geom for data.frame
#'
#' This geom constructor will simply call [gg.bru_prediction()] for the data
#' provided.
#'
#' Requires the `ggplot2` package.
#'
#' @name gg.data.frame
#' @export
#' @param ... Arguments passed on to [gg.bru_prediction()].
#' @return Concatenation of a `geom_line` value and optionally a `geom_ribbon`
#'   value.
#' @family geomes for inla and inlabru predictions
#' @example inst/examples/gg.prediction.R

gg.data.frame <- function(...) {
  gg.bru_prediction(...)
}


#' Geom for predictions
#'
#' @description
#'
#' This geom serves to visualize `prediction` objects which usually results from
#' a call to [predict.bru()]. Predictions objects provide summary statistics
#' (mean, median, sd, ...) for one or more random variables. For single
#' variables (or if requested so by setting `bar = TRUE`), a boxplot-style geom
#' is constructed to show the statistics. For multivariate predictions the mean
#' of each variable (y-axis) is plotted against the row number of the variable
#' in the prediction data frame (x-axis) using `geom_line`. In addition, a
#' `geom_ribbon` is used to show the confidence interval.
#'
#' Note: `gg.bru_prediction` also understands the format of INLA-style posterior
#' summaries, e.g. `fit$summary.fixed` for an inla object `fit`
#'
#' Requires the `ggplot2` package.
#'
#' @name gg.bru_prediction
#' @export
#' @importFrom utils modifyList
#' @param data A prediction object, usually the result of a [predict.bru()]
#'   call.
#' @param mapping a set of aesthetic mappings created by `aes`. These are passed
#'   on to `geom_line`.
#' @param ribbon If TRUE, plot a ribbon around the line based on the smallest
#'   and largest quantiles present in the data, found by matching names starting
#'   with `q` and followed by a numerical value.  `inla()`-style
#'   `numeric+"quant"` names are converted to inlabru style before matching.
#' @param alpha The ribbons numeric alpha (transparency) level in `[0,1]`.
#' @param bar If TRUE plot boxplot-style summary for each variable.
#' @param \dots Arguments passed on to `geom_line`.
#' @return Concatenation of a `geom_line` value and optionally a `geom_ribbon`
#'   value.
#' @family geomes for inla and inlabru predictions
#' @example inst/examples/gg.prediction.R

gg.bru_prediction <- function(data,
                              mapping = NULL,
                              ribbon = TRUE,
                              alpha = NULL,
                              bar = FALSE,
                              ...) {
  if (inherits(data, c("Spatial", "SpatRaster", "RasterLayer", "sf"))) {
    return(NextMethod())
  }

  if (is.null(alpha)) {
    alpha <- 0.3
  }

  if (is.list(data) && !is.data.frame(data)) {
    data <- lapply(names(data), function(x) {
      if (NROW(data[[x]]) == 1) {
        rownames(data[[x]]) <- x
      } else {
        rownames(data[[x]]) <- paste0(x, ".", seq_len(NROW(data[[x]])))
      }
      data[[x]]
    })
    data <- dplyr::bind_rows(data)
  }

  requireNamespace("ggplot2")
  # Convert from old and inla style names
  if (("median" %in% names(data)) &&
    !("0.5quant" %in% names(data)) &&
    !("q0.5" %in% names(data))) {
    names(data)[names(data) == "median"] <- "q0.5"
  }
  new_quant_names <- list(
    q0.025 = "0.025quant",
    q0.5 = "0.5quant",
    q0.975 = "0.975quant"
  )
  for (quant_name in names(new_quant_names[new_quant_names %in% names(data)])) {
    names(data)[names(data) == new_quant_names[[quant_name]]] <- quant_name
  }
  # Find quantile levels
  quant_names <- names(data)[grepl("^q[01]\\.?[0-9]*$", names(data))]
  if (length(quant_names) > 0) {
    quant_probs <- as.numeric(sub("^q", "", quant_names))
    quant_names <- quant_names[order(quant_probs)]
    quant_probs <- quant_probs[order(quant_probs)]
    lqname <- quant_names[1]
    uqname <- quant_names[length(quant_names)]
  }

  if (bar || (nrow(data) == 1)) {
    sz <- 10 # bar width
    med_sz <- 4 # median marker size

    data <- cbind(
      data.frame(data),
      data.frame(
        variable = rownames(data),
        summary = data$mean,
        sdmax = data$mean + data$sd,
        sdmin = data$mean - data$sd
      )
    )

    geom <- c(
      ggplot2::geom_point(
        data = data,
        mapping = ggplot2::aes(
          x = .data$variable,
          y = .data$summary,
          color = .data$variable
        ),
        shape = 95, size = 0
      )
    )
    if (length(quant_names) > 0) {
      geom <- c(
        geom,
        ggplot2::geom_segment(
          data = data,
          mapping = ggplot2::aes(
            y = .data[[lqname]],
            yend = .data[[uqname]],
            x = .data$variable,
            xend = .data$variable,
            color = .data$variable
          ),
          linewidth = sz
        )
      )
    }

    # Mean and median
    geom <- c(
      geom,
      ggplot2::geom_point(
        data = data,
        mapping = ggplot2::aes(
          x = .data$variable,
          y = .data$mean
        ),
        color = "black", shape = 20, size = med_sz
      )
    )
    if ("q0.5" %in% quant_names) {
      geom <- c(
        geom,
        ggplot2::geom_point(
          data = data,
          mapping = ggplot2::aes(
            x = .data$variable,
            y = .data$q0.5
          ),
          color = "black", shape = 3, size = sz * 2 / 3
        )
      )
    }
    geom <- c(
      geom,
      ggplot2::coord_flip()
    )
  } else {
    if ("pdf" %in% names(data)) {
      y.str <- "pdf"
      ribbon <- FALSE
    } else if ("mean" %in% names(data)) {
      y.str <- "mean"
    } else if ("q0.5" %in% names(data)) {
      y.str <- "q0.5"
    } else if ("median" %in% names(data)) {
      y.str <- "q0.5"
    } else {
      stop(paste0(
        "Prediction has neither mean nor median/q0.5 or pdf as column. ",
        "Don't know what to plot."
      ))
    }

    line.map <- ggplot2::aes(
      x = .data[[names(data)[1]]],
      y = .data[[y.str]]
    )

    if (!is.null(mapping)) {
      line.map <- utils::modifyList(line.map, mapping)
    }

    geom <- ggplot2::geom_line(data = data, line.map, ...)

    if (ribbon && length(quant_names) > 0) {
      ribbon.map <- ggplot2::aes(
        x = .data[[names(data)[1]]],
        ymin = .data[[lqname]],
        ymax = .data[[uqname]]
      )
      # Use line color for ribbon filling
      if ("colour" %in% names(line.map)) {
        ribbon.map <- modifyList(
          ribbon.map,
          ggplot2::aes(fill = .data[[line.map[["colour"]]]])
        )
      }
      geom <- c(
        geom,
        ggplot2::geom_ribbon(data = data, ribbon.map, alpha = alpha)
      )
    }
  }
  geom
}

#' @rdname gg.bru_prediction
#' @export
gg.prediction <- function(data, ...) {
  gg.bru_prediction(data = data, ...)
}

#' Geom for SpatialPoints objects
#'
#' This function coerces the `SpatialPoints` into a `data.frame` and uses
#' `geom_point` to plot the points. Requires the `ggplot2` package.
#'
#' @export
#' @param data A SpatialPoints object.
#' @param mapping Aesthetic mappings created by `aes` used to update the default
#'   mapping. The default mapping is
#'   ```{r eval = FALSE}
#'   ggplot2::aes(
#'     x = .data[[sp::coordnames(data)[1]]],
#'     y = .data[[sp::coordnames(data)[2]]]
#'   )
#'   ```
#' @param crs A `sp::CRS` object defining the coordinate system to project the
#'   data to before plotting.
#' @param ... Arguments passed on to `geom_point`.
#' @return A `geom_point` return value
#' @family geomes for spatial data
#'
#' @example inst/examples/gg.sp.R

gg.SpatialPoints <- function(data, mapping = NULL, crs = NULL, ...) {
  bru_safe_sp(force = TRUE)
  requireNamespace("ggplot2")
  if (!is.null(crs)) {
    data <- fm_transform(data, crs)
  }

  df <- as.data.frame(data)
  cnames <- sp::coordnames(data)
  if (is.null(cnames)) {
    cnames <- c("x", "y")
  }

  dmap <- ggplot2::aes(
    x = .data[[cnames[1]]],
    y = .data[[cnames[2]]]
  )

  if (!is.null(mapping)) {
    dmap <- modifyList(dmap, mapping)
  }

  ggplot2::geom_point(data = df, mapping = dmap, ...)
}

#' Geom helper for sf objects
#'
#' This function uses `geom_sf()`, unless overridden by the geom argument.
#' Requires the `ggplot2` package.
#'
#' @name gg.sf
#' @export
#' @param data An `sf` object.
#' @param mapping Default mapping is `ggplot2::aes(geometry = ...)`,
#' where the geometry name is obtained from `attr(data, "sf_column")`.
#' This is merged with the user supplied mapping.
#' @param geom Either "sf" (default) or "tile". For "tile", uses
#' `geom_tile(..., stat = "sf_coordinates")`, intended for converting point data
#' to grid tiles with the `fill` aesthetic, which is by default set to the first
#' data column.
#' @param ... Arguments passed on to `geom_sf` or `geom_tile`.
#' @return A ggplot return value
#' @family geomes for spatial data
#' @example inst/examples/gg.sf.R

gg.sf <- function(data, mapping = NULL, ..., geom = "sf") {
  requireNamespace("ggplot2")

  sf_column <- attr(data, "sf_column")
  geom <- match.arg(geom, c("sf", "tile"))

  if (identical(geom, "sf")) {
    dmap <- ggplot2::aes(
      geometry = .data[[sf_column]]
    )
    if (!is.null(mapping)) {
      dmap <- modifyList(dmap, mapping)
    }
    result <- ggplot2::geom_sf(data = data, mapping = dmap, ...)
  } else {
    data_column <- setdiff(names(data), sf_column)
    if (length(data_column) > 0) {
      data_column <- data_column[1]
      dmap <- ggplot2::aes(
        geometry = .data[[sf_column]],
        fill = .data[[data_column]]
      )
    } else {
      dmap <- ggplot2::aes(
        geometry = .data[[sf_column]]
      )
    }
    if (!is.null(mapping)) {
      dmap <- modifyList(dmap, mapping)
    }
    result <-
      c(
        ggplot2::geom_tile(
          data = data,
          mapping = dmap,
          stat = "sf_coordinates",
          ...
        ),
        if (fmesher::fm_crs_is_null(fm_crs(data))) {
          ggplot2::coord_sf(default = TRUE)
        } else {
          ggplot2::coord_sf(default = TRUE, crs = fm_crs(data))
        }
      )
  }

  result
}

#' Geom for SpatialLines objects
#'
#' @description
#'
#' Extracts start and end points of the lines and calls `geom_segment` to plot
#' lines between them. Requires the `ggplot2` package.
#'
#' @export
#' @param data A `SpatialLines` or `SpatialLinesDataFrame` object.
#' @param mapping Aesthetic mappings created by `ggplot2::aes` or
#'   `ggplot2::aes_` used to update the default mapping. The default mapping is
#'   ```{r eval = FALSE}
#'   ggplot2::aes(
#'     x = .data[[sp::coordnames(data)[1]]],
#'     y = .data[[sp::coordnames(data)[2]]],
#'     xend = .data[[paste0("end.", sp::coordnames(data)[1])]],
#'     yend = .data[[paste0("end.", sp::coordnames(data)[2])]])
#'   ```
#' @param crs A `CRS` object defining the coordinate system to project the data
#'   to before plotting.
#' @param ... Arguments passed on to `ggplot2::geom_segment`.
#' @return A `geom_segment`` return value.
#' @family geomes for spatial data
#' @example inst/examples/gg.sp.R

gg.SpatialLines <- function(data, mapping = NULL, crs = NULL, ...) {
  bru_safe_sp(force = TRUE)
  requireNamespace("ggplot2")
  if (!is.null(crs)) {
    data <- fm_transform(data, crs)
  }

  qq <- sp::coordinates(data)
  cnames <- sp::coordnames(data)
  if (is.null(cnames)) {
    cnames <- c("x", "y")
  }
  sp <- do.call(rbind, lapply(
    qq,
    function(k) do.call(rbind, lapply(k, function(x) x[1:(nrow(x) - 1), ]))
  ))
  ep <- do.call(rbind, lapply(
    qq,
    function(k) do.call(rbind, lapply(k, function(x) x[2:(nrow(x)), ]))
  ))
  colnames(sp) <- cnames
  colnames(ep) <- paste0("end.", cnames)
  if (!inherits(data, "SpatialLinesDataFrame")) {
    df <- data.frame(cbind(sp, ep))
  } else {
    df <- data.frame(
      cbind(sp, ep),
      data@data[
        rep(
          seq_len(NROW(data)),
          times = vapply(
            qq,
            FUN = function(x) {
              NROW(x[[1]]) - 1
            },
            0
          )
        ), ,
        drop = FALSE
      ]
    )
  }
  dmap <- ggplot2::aes(
    x = .data[[cnames[1]]],
    y = .data[[cnames[2]]],
    xend = .data[[paste0("end.", cnames[1])]],
    yend = .data[[paste0("end.", cnames[2])]]
  )

  if (!is.null(mapping)) {
    dmap <- modifyList(dmap, mapping)
  }

  ggplot2::geom_segment(data = df, mapping = dmap, ...)
}

#' Geom for SpatialPolygons objects
#'
#' Uses the `ggplot2::fortify()` function to turn the `SpatialPolygons` objects
#' into a `data.frame`. Then calls `geom_polygon` to plot the polygons. Requires
#' the `ggplot2` package.
#'
#' @export
#' @param data A `SpatialPolygons` or `SpatialPolygonsDataFrame` object.
#' @param mapping Aesthetic mappings created by `aes` used to update the default
#'   mapping.
#' @param crs A `CRS` object defining the coordinate system to project the data
#'   to before plotting.
#' @param ... Arguments passed on to `geom_sf`.
#' Unless specified by the user,
#' the argument `alpha = 0.2` (alpha level for polygon filling) is added.
#' @return A `geom_sf` object.
#' @details Up to version `2.10.0`, the `ggpolypath` package was used to ensure
#'   proper plotting, since the `ggplot2::geom_polygon` function doesn't always
#'   handle geometries with holes properly. After `2.10.0`, the object is
#'   converted to `sf` format and passed on to [gg.sf()] instead, as `ggplot2`
#'   version `3.4.4` deprecated the intenrally used `ggplot2::fortify()` method
#'   for `SpatialPolygons/DataFrame` objects.
#' @family geomes for spatial data
#' @example inst/examples/gg.sp.R

gg.SpatialPolygons <- function(data, mapping = NULL, crs = NULL, ...) {
  bru_safe_sp(force = TRUE)
  data <- sf::st_as_sf(data)
  if (!is.null(crs)) {
    data <- fm_transform(data, crs, passthrough = TRUE)
  }
  arg <- list(...)
  if ("alpha" %in% names(arg)) {
    gg(data, mapping = mapping, ...)
  } else {
    gg(data, mapping = mapping, ..., alpha = 0.2)
  }
}

#' Geom for SpatialGridDataFrame objects
#'
#' Coerces input `SpatialGridDataFrame` to `SpatialPixelsDataFrame` and calls
#' [gg.SpatialPixelsDataFrame()] to plot it.
#' Requires the `ggplot2` package.
#'
#' @export
#' @param data A SpatialGridDataFrame object.
#' @param ... Arguments passed on to [gg.SpatialPixelsDataFrame()].
#' @return A `geom_tile` value.
#' @family geomes for spatial data
#' @example inst/examples/gg.sp.R

gg.SpatialGridDataFrame <- function(data, ...) {
  bru_safe_sp(force = TRUE)
  data <- as(data, "SpatialPixelsDataFrame")
  gg(data, ...)
}

#' Geom for SpatialPixelsDataFrame objects
#'
#' Coerces input SpatialPixelsDataFrame to data.frame and uses `geom_tile` to
#' plot it.
#' Requires the `ggplot2` package.
#'
#' @export
#' @param data A SpatialPixelsDataFrame object.
#' @param mapping Aesthetic mappings created by `aes` used to update the default
#'   mapping. The default mapping is
#'   ```{r eval=FALSE}
#'   ggplot2::aes(
#'     x = .data[[sp::coordnames(data)[1]]],
#'     y = .data[[sp::coordnames(data)[2]]],
#'     fill = .data[[names(data)[[1]]]]
#'   )
#'   ```
#' @param crs A `sp::CRS` object defining the coordinate system to project the
#'   data to before plotting.
#' @param mask A `sp::SpatialPolygons` object defining the region that is
#'   plotted.
#' @param ... Arguments passed on to `geom_tile`.
#' @return A `geom_tile` return value.
#' @family geomes for spatial data
#' @example inst/examples/gg.sp.R
gg.SpatialPixelsDataFrame <- function(data,
                                      mapping = NULL,
                                      crs = NULL,
                                      mask = NULL, ...) {
  bru_safe_sp(force = TRUE)
  if (!is.null(crs)) {
    data <- fm_transform(data, crs)
  }
  if (!is.null(mask)) {
    data <- data[as.vector(!is.na(sp::over(data, mask))), ]
  }

  df <- as.data.frame(data)

  if (length(names(data)) == 0) {
    dmap <- ggplot2::aes(
      x = .data[[sp::coordnames(data)[1]]],
      y = .data[[sp::coordnames(data)[2]]]
    )
  } else {
    dmap <- ggplot2::aes(
      x = .data[[sp::coordnames(data)[1]]],
      y = .data[[sp::coordnames(data)[2]]],
      fill = .data[[names(data)[[1]]]]
    )
  }

  if (!is.null(mapping)) {
    dmap <- modifyList(dmap, mapping)
  }
  gm <- ggplot2::geom_tile(data = df, mapping = dmap, ...)

  gm
}


#' Geom for SpatialPixels objects
#'
#' Uses `geom_point` to plot the pixel centers.
#' Requires the `ggplot2` package.
#'
#' @export
#' @param data A `sp::SpatialPixels` object.
#' @param ... Arguments passed on to `geom_tile`.
#' @return A `geom_tile` return value.
#' @family geomes for spatial data
#' @examples
#' if (require("ggplot2", quietly = TRUE) &&
#'   requireNamespace("terra", quietly = TRUE) &&
#'   bru_safe_sp()) {
#'   # Load Gorilla data
#'
#'   gcov <- gorillas_sf_gcov()
#'   elev <- terra::as.data.frame(gcov$elevation, xy = TRUE)
#'   pxl <- sf::as_Spatial(sf::st_as_sf(elev, coords = c("x", "y")))
#'
#'   # Turn elevation covariate into SpatialPixels
#'   pxl <- sp::SpatialPixels(pxl)
#'
#'   # Plot the pixel centers
#'   ggplot() +
#'     gg(pxl, size = 0.1)
#' }
gg.SpatialPixels <- function(data, ...) {
  bru_safe_sp(force = TRUE)
  gg(sp::SpatialPoints(data), ...)
}


#' Geom wrapper for SpatRaster objects
#'
#' Convenience wrapper function for `tidyterra::geom_spatraster()`.
#' Requires the `ggplot2` and `tidyterra` packages.
#'
#' @export
#' @param data A SpatRaster object.
#' @param ... Arguments passed on to `geom_spatraster`.
#' @return The output from `geom_spatraster.
#' @family geomes for spatial data
#' @examples
#' if (require("ggplot2", quietly = TRUE) &&
#'   requireNamespace("terra", quietly = TRUE) &&
#'   require("tidyterra", quietly = TRUE)) {
#'   # Load Gorilla covariates
#'
#'   gcov <- gorillas_sf_gcov()
#'
#'   # Plot the pixel centers
#'   ggplot() +
#'     gg(gcov$elevation)
#' }
gg.SpatRaster <- function(data, ...) {
  if (!requireNamespace("tidyterra", quietly = TRUE)) {
    stop(paste0(
      "gg.SpatRaster requires the 'tidyterra' package, ",
      "which is not installed.\n",
      "Please install it and try again!"
    ))
  }
  if (".mask" %in% names(data)) {
    requireNamespace("terra")
    data <-
      terra::mask(
        data,
        mask = data$.mask,
        maskvalues = FALSE,
        updatevalue = NA
      )
  }
  tidyterra::geom_spatraster(data = data, ...)
}





#' Geom for fm_mesh_2d objects
#'
#' @description
#' This function extracts the graph of an [fmesher::fm_mesh_2d] object and uses
#' `geom_line` to visualize the graph's edges. Alternatively, if the `color`
#' argument is provided, interpolates the colors across for a set of
#' SpatialPixels covering the mesh area and calls [gg.SpatialPixelsDataFrame()]
#' to plot the interpolation.
#' Requires the `ggplot2` package.
#'
#' Also see the [fmesher::geom_fm()] method.
#'
#' @export
#'
#' @param data An `fm_mesh_2d` object.
#' @param color A vector of scalar values to fill the mesh with colors.
#' The length of the vector mus correspond to the number of mesh vertices.
#' The alternative name `colour` is also recognised.
#' @param alpha A vector of scalar values setting the alpha value of the colors
#'   provided.
#' @param edge.color Color of the regular mesh edges.
#' @param edge.linewidth Line width for the regular mesh edges. Default 0.25
#' @param interior If TRUE, plot the interior boundaries of the mesh.
#' @param int.color Color used to plot the interior constraint edges.
#' @param int.linewidth Line width for the interior constraint edges. Default
#'   0.5
#' @param exterior If TRUE, plot the exterior boundaries of the mesh.
#' @param ext.color Color used to plot the exterior boundary edges.
#' @param ext.linewidth Line width for the exterior boundary edges. Default 1
#' @param crs A CRS object supported by [fm_transform()] defining the coordinate
#'   system to project the mesh to before plotting.
#' @param nx Number of pixels in x direction (when plotting using the color
#'   parameter).
#' @param ny Number of pixels in y direction (when plotting using the color
#'   parameter).
#' @param mask A `SpatialPolygon` or `sf` polygon defining the region that is
#'   plotted.
#' @param ... ignored arguments (S3 generic compatibility).
#' @return `geom_line` return values or, if the color argument is used, the
#'   values of [gg.SpatialPixelsDataFrame()].
#' @family geomes for meshes
#' @example inst/examples/gg.fm_mesh_2d.R

gg.fm_mesh_2d <- function(data,
                          color = NULL,
                          alpha = NULL,
                          edge.color = "grey",
                          edge.linewidth = 0.25,
                          interior = TRUE,
                          int.color = "blue",
                          int.linewidth = 0.5,
                          exterior = TRUE,
                          ext.color = "black",
                          ext.linewidth = 1,
                          crs = NULL,
                          mask = NULL,
                          nx = 500, ny = 500,
                          ...) {
  requireNamespace("ggplot2")
  if (is.null(color) && ("colour" %in% names(list(...)))) {
    color <- list(...)[["colour"]]
  }
  if (!is.null(color)) {
    if (inherits(mask, "Spatial")) {
      # For backwards compatibility with old plotting code, use Spatial
      # format, and mask only in the gg.Spatial* method.
      px <- fm_pixels(data, dims = c(nx, ny), format = "sp")
      A <- fm_basis(data, px)
      px$color <- as.vector(A %*% color)
      if (!is.null(alpha)) {
        px$alpha <- as.vector(A %*% alpha)
        gg <- gg(px, mask = mask, alpha = "alpha")
      } else {
        gg <- gg(px, mask = mask)
      }
    } else {
      px <- fm_pixels(data, dims = c(nx, ny), mask = mask, format = "sf")
      proj <- fm_evaluator(data, px)
      px$color <- fm_evaluate(proj, field = color)
      if (!is.null(alpha)) {
        if (length(alpha) == 1) {
          px$alpha <- alpha
        } else {
          px$alpha <- fm_evaluate(proj, field = alpha)
        }
        gg <- gg(px,
          ggplot2::aes(fill = .data[["color"]]),
          alpha = px[["alpha"]],
          geom = "tile"
        )
      } else {
        gg <- gg(px,
          ggplot2::aes(fill = .data[["color"]]),
          geom = "tile"
        )
      }
    }

    return(gg)
  }

  if (data$manifold == "S2") {
    stop("Geom not implemented for spherical meshes (manifold = S2)")
  }
  if (!is.null(crs)) {
    data <- fm_transform(data, crs = crs)
  }

  df <- rbind(
    data.frame(
      a = data$loc[data$graph$tv[, 1], c(1, 2)],
      b = data$loc[data$graph$tv[, 2], c(1, 2)]
    ),
    data.frame(
      a = data$loc[data$graph$tv[, 2], c(1, 2)],
      b = data$loc[data$graph$tv[, 3], c(1, 2)]
    ),
    data.frame(
      a = data$loc[data$graph$tv[, 1], c(1, 2)],
      b = data$loc[data$graph$tv[, 3], c(1, 2)]
    )
  )

  colnames(df) <- c("x", "y", "xend", "yend")
  mp <- ggplot2::aes(
    x = .data$x,
    y = .data$y,
    xend = .data$xend,
    yend = .data$yend
  )
  msh <- ggplot2::geom_segment(
    data = df,
    mapping = mp,
    color = edge.color,
    linewidth = edge.linewidth
  )

  # Outer boundary
  if (exterior) {
    df <- data.frame(
      data$loc[data$segm$bnd$idx[, 1], 1:2],
      data$loc[data$segm$bnd$idx[, 2], 1:2]
    )
    colnames(df) <- c("x", "y", "xend", "yend")
    bnd <- ggplot2::geom_segment(
      data = df,
      mapping = mp,
      color = ext.color,
      linewidth = ext.linewidth
    )
  } else {
    bnd <- NULL
  }

  if (interior) {
    # Interior boundary
    df <- data.frame(
      data$loc[data$segm$int$idx[, 1], 1:2],
      data$loc[data$segm$int$idx[, 2], 1:2]
    )
    colnames(df) <- c("x", "y", "xend", "yend")
    if (nrow(df) == 0) {
      int <- NULL
    } else {
      int <- ggplot2::geom_segment(
        data = df,
        mapping = mp,
        color = int.color,
        linewidth = int.linewidth
      )
    }
  } else {
    int <- NULL
  }

  # Return combined geomes
  c(msh, bnd, int)
}


#' Geom for fm_mesh_1d objects
#'
#' This function generates a `geom_point` object showing the knots (vertices)
#' of a 1D mesh.
#' Requires the `ggplot2` package.
#'
#' @export
#' @param data An [fmesher::fm_mesh_1d] object.
#' @param mapping aesthetic mappings created by `aes`. These are passed on to
#'   `geom_point`.
#' @param y Single or vector numeric defining the y-coordinates of the mesh
#'   knots to plot.
#' @param shape Shape of the knot markers.
#' @param ... parameters passed on to `geom_point`.
#' @return An object generated by `geom_point`.
#' @family geomes for meshes
#'
#' @examples
#' \donttest{
#' if (require("fmesher", quietly = TRUE) &&
#'   require("ggplot2", quietly = TRUE)) {
#'   # Create a 1D mesh
#'
#'   mesh <- fm_mesh_1d(seq(0, 10, by = 0.5))
#'
#'   # Plot it
#'
#'   ggplot() +
#'     gg(mesh)
#'
#'   # Plot it using a different shape and size for the mesh nodes
#'
#'   ggplot() +
#'     gg(mesh, shape = "|", size = 5)
#' }
#' }
#'
gg.fm_mesh_1d <- function(data,
                          mapping = ggplot2::aes(.data[["x"]], .data[["y"]]),
                          y = 0, shape = 4, ...) {
  df <- data.frame(x = data$loc, y = y)
  ggplot2::geom_point(data = df, mapping = mapping, shape = shape, ...)
}

#' @describeIn gg.fm_mesh_2d Alias for `gg.fm_mesh_2d`, supporting `inla.mesh`
#' objects.
#' @export
gg.inla.mesh <- gg.fm_mesh_2d

#' @describeIn gg.fm_mesh_1d Alias for `gg.fm_mesh_1d`, supporting
#' `inla.mesh.1d` objects.
#' @export
gg.inla.mesh.1d <- gg.fm_mesh_1d



#' Geom for RasterLayer objects
#'
#' This function takes a RasterLayer object, converts it into a
#' `SpatialPixelsDataFrame` and uses `geom_tile` to plot the data.
#'
#' This function requires the `raster` and `ggplot2` packages.
#'
#' @export
#' @param data A RasterLayer object.
#' @param mapping aesthetic mappings created by `aes`. These are passed on to
#'   `geom_tile`.
#' @param ... Arguments passed on to `geom_tile`.
#' @return An object returned by `geom_tile`
#' @family geomes for Raster data
#'
#' @examples
#' \dontrun{
#' # Some features require the raster and spatstat.data packages.
#' if (require("spatstat.data", quietly = TRUE) &&
#'   require("raster", quietly = TRUE) &&
#'   require("ggplot2", quietly = TRUE)) {
#'   # Load Gorilla data
#'   data("gorillas", package = "spatstat.data", envir = environment())
#'
#'   # Convert elevation covariate to RasterLayer
#'
#'   elev <- as(gorillas.extra$elevation, "RasterLayer")
#'
#'   # Plot the elevation
#'
#'   ggplot() +
#'     gg(elev)
#' }
#' }
gg.RasterLayer <- function(data,
                           mapping = ggplot2::aes(
                             x = .data[["x"]],
                             y = .data[["y"]],
                             fill = .data[["layer"]]
                           ),
                           ...) {
  requireNamespace("ggplot2")
  requireNamespace("raster")
  spdf <- as(data, "SpatialPixelsDataFrame")
  df <- as.data.frame(spdf)
  # head(r.df)
  # g <- ggplot(r.df, ggplot2::aes(x=x, y=y)) +
  #   geom_tile(ggplot2::aes(fill = layer)) +
  #   coord_equal()
  ggplot2::geom_tile(data = df, mapping = mapping, ...)
}

#' Plot method for posterior marginals estimated by bru
#'
#' @description
#'
#' From version `2.11.0`, `plot.bru(x, ...)` calls `plot.inla(x, ...)`
#' from the `INLA` package, unless the first argument after `x` is a
#' `character`, in which case the pre-`2.11.0` behaviour is used, calling
#' `plotmarginal.inla(x, ...)` instead.
#'
#' Requires the `ggplot2` package.
#'
#' @method plot bru
#' @export
#' @param x a fitted [bru()] model.
#' @param ... Options passed on to other methods.
#'
#' @examples
#' \dontrun{
#' if (require("ggplot2", quietly = TRUE)) {
#'   # Generate some data and fit a simple model
#'   input.df <- data.frame(x = cos(1:10))
#'   input.df <- within(
#'     input.df,
#'     y <- 5 + 2 * x + rnorm(length(x), mean = 0, sd = 0.1)
#'   )
#'   fit <- bru(y ~ x, family = "gaussian", data = input.df)
#'   summary(fit)
#'
#'   # Plot the posterior density of the model's x-effect
#'   plot(fit, "x")
#' }
#' }
#'
plot.bru <- function(x, ...) {
  ll <- list(...)
  if (length(ll) == 0) {
    NextMethod("plot")
  } else if (is.character(ll[[1]])) {
    plotmarginal.inla(x, ...)
  } else {
    NextMethod("plot")
  }
}

#' @title Plot prediction using ggplot2
#'
#' @description
#' Generates a base ggplot2 using `ggplot()` and adds a geom for input `x` using
#' [gg].
#'
#' Requires the `ggplot2` package.
#'
#' @name plot.bru_prediction
#' @param x a prediction object.
#' @param y Ignored argument but required for S3 compatibility.
#' @param ... Arguments passed on to [gg.prediction()].
#' @return an object of class `gg`
#' @example inst/examples/gg.prediction.R
#' @export
#' @method plot bru_prediction

plot.bru_prediction <- function(x, y = NULL, ...) {
  requireNamespace("ggplot2")
  ggplot2::ggplot() +
    gg(x, ...)
}

#' @rdname plot.bru_prediction
#' @export
#' @method plot prediction
plot.prediction <- function(x, y = NULL, ...) {
  requireNamespace("ggplot2")
  ggplot2::ggplot() +
    gg(x, ...)
}



#' @title Multiple ggplots on a page.
#'
#' @description
#' Renders multiple ggplots on a single page.
#'
#' @param ... Comma-separated `ggplot` objects.
#' @param plotlist A list of `ggplot` objects - an alternative to the
#'   comma-separated argument above.
#' @param cols Number of columns of plots on the page.
#' @param layout A matrix specifying the layout. If present, `cols` is ignored.
#'   If the layout is something like `matrix(c(1,2,3,3), nrow=2, byrow=TRUE)`,
#'   then plot 1 will go in the upper left, 2 will go in the upper right, and 3
#'   will go all the way across the bottom.
#'
#' @author David L. Borchers \email{dlb@@st-andrews.ac.uk}
#'
#' @source
#' <http://www.cookbook-r.com/Graphs/Multiple_graphs_on_one_page_(ggplot2)/>
#'
#' @examples
#' if (require("ggplot2", quietly = TRUE)) {
#'   df <- data.frame(x = 1:10, y = 1:10, z = 11:20)
#'   pl1 <- ggplot(data = df) +
#'     geom_line(mapping = aes(x, y), color = "red")
#'   pl2 <- ggplot(data = df) +
#'     geom_line(mapping = aes(x, z), color = "blue")
#'   multiplot(pl1, pl2, cols = 2)
#' }
#' @export
#
multiplot <- function(..., plotlist = NULL, cols = 1, layout = NULL) {
  # Make a list from the ... arguments and plotlist
  plots <- c(list(...), plotlist)

  numPlots <- length(plots)

  # If layout is NULL, then use 'cols' to determine layout
  if (is.null(layout)) {
    # Make the panel
    # ncol: Number of columns of plots
    # nrow: Number of rows needed, calculated from # of cols
    layout <- matrix(seq(1, cols * ceiling(numPlots / cols)),
      ncol = cols, nrow = ceiling(numPlots / cols)
    )
  }

  if (numPlots == 1) {
    print(plots[[1]])
  } else {
    # Set up the page
    grid::grid.newpage()
    grid::pushViewport(grid::viewport(
      layout = grid::grid.layout(nrow(layout), ncol(layout))
    ))

    # Make each plot, in the correct location
    for (i in 1:numPlots) {
      # Get the i,j matrix positions of the regions that contain this subplot
      matchidx <- as.data.frame(which(layout == i, arr.ind = TRUE))

      print(plots[[i]], vp = grid::viewport(
        layout.pos.row = matchidx$row,
        layout.pos.col = matchidx$col
      ))
    }
  }
}



# Default color palette for plots using \link{gg}
#
# @export
# @param
# @return Color values as returned by
#   `RColorBrewer::brewer.pal(9, "YlOrRd"`

bru.pal <- function() {
  # library(RColorBrewer)
  # brewer.pal(9, "YlOrRd")
  pal <- c(
    "#FFFFCC",
    "#FFEDA0",
    "#FED976",
    "#FEB24C",
    "#FD8D3C",
    "#FC4E2A",
    "#E31A1C",
    "#BD0026",
    "#800026"
  )
}