R/2kyears_temp.R

In cortinah/hockeystick: Download and Visualize Essential Climate Change Data

#' Download and plot essential climate data
#'
#' Retrieves the Common Era Global Surface Temperature Reconstructions.
#' Source is PAGES2k Consortium and NOAA National Centers for Environmental Information.
#' \url{https://www.ncei.noaa.gov/access/paleo-search/study/26872}
#'
#' @name get_temp2k
#' @param use_cache (boolean) Return cached data if available, defaults to TRUE. Use FALSE to fetch updated data.
#' @param write_cache (boolean) Write data to cache, defaults to FALSE. Use TRUE to write data to cache for later use. Can also be set using options(hs_write_cache=TRUE)
#'
#' @return Invisibly returns a tibble with filtered and unfiltered temperature reconstructions and Cowtan & Way instrumental temperatures.
#'
#' `get_temp2k` invisibly returns a tibble with the PAGES2k Consortium temperature reconstruction (years 1-2000 CE) and instrumental record (years 1850-2017 CE).  Temperatures represent deviations from the 1961-1990 mean.
#'
#' The returned object includes annual average temperature anomalies as well as filtered anomalies using a 31-year Butterworth filter.
#' Reconstructions use seven different statistical methods that draw from a global collection of temperature-sensitive palaeoclimate records.
#'
#' Methodology described in PAGES2k (2019) \url{https://www.nature.com/articles/s41561-019-0400-0}
#'
#' @importFrom readr read_csv
#' @importFrom lubridate ymd
#' @importFrom utils download.file
#'
#' @examples
#' \donttest{
#' # Fetch temp anomaly from cache if available:
#' anomaly <- get_temp2k()
#' #
#' # Force cache refresh:
#' anomaly <- get_temp2k(use_cache=FALSE)
#' #
#' # Review cache contents and last update dates:
#' hockeystick_cache_details()
#' #
#' # Plot output using package's built-in ggplot2 settings
#' plot_temp2k(anomaly) }
#'
#' @author Hernando Cortina, \email{hch@@alum.mit.edu}
#' @references
#' \itemize{
#' \item PAGES2k Common Era Surface Temperature Reconstructions. \url{https://www.ncei.noaa.gov/access/paleo-search/study/26872}
#' \item PAGES 2k Consortium., Neukom, R., Barboza, L.A. et al. Consistent multidecadal variability in global temperature reconstructions and simulations over the Common Era. \emph{Nat. Geosci.} 12, 643–649 (2019).  \doi{10.1038/s41561-019-0400-0}
#'  }
#'
#' @export
get_temp2k <- function(use_cache = TRUE, write_cache = getOption("hs_write_cache")) {

hs_path <- tools::R_user_dir("hockeystick","cache")

if (use_cache) {
  if (file.exists(file.path(hs_path,'temp2k.rds'))) return(invisible(readRDS((file.path(hs_path,'temp2k.rds')))))   }

file_url <- 'https://www.ncei.noaa.gov/pub/data/paleo/pages2k/neukom2019temp/recons/Full_ensemble_median_and_95pct_range.txt'
connected <- .isConnected(file_url)
if (!connected) {message("Retrieving remote data requires internet connectivity."); return(invisible(NULL))}

dl <- tempfile()
download.file(file_url, dl)
temp2k <- suppressMessages( read_table(dl, col_names = FALSE, skip = 5,
                                         col_types = 'innnnnnnn') )

colnames(temp2k) <- c('year', 'instrumental', 'ensemble_median', 'ensemble_2.5', 'ensemble_97.5', 'filtered_instrumental',
                       'filtered_ensemble_median', 'filtered_ensemble_2.5', 'filtered_ensemble_97.5')

temp2k_l <- pivot_longer(temp2k, -year, names_to='measure')

temp2k_l <- filter(temp2k_l, measure %in% c('instrumental', 'ensemble_median', 'filtered_instrumental', 'filtered_ensemble_median'))

dir.create(hs_path, showWarnings = FALSE, recursive = TRUE)
if (write_cache) saveRDS(temp2k_l, file.path(hs_path, 'temp2k.rds'))

invisible(temp2k_l)
}


#' Download and plot essential climate data
#'
#' Plots the Common Era 2000-year global temperature anomaly retrieved using `get_temp2k()` with ggplot2. The output ggplot2 object may be further modified.
#'
#'
#' @name plot_temp2k
#' @param dataset Name of the tibble generated by `get_temp2k`
#' @param instrumental (boolean) Include the Cowtan & Way instrumental temperatures through 2017. Defaults to TRUE.
#' @param filtered (boolean) Use the filtered temperatures provided by PAGES2k Consortium. Temperatures filtered using a 31-year Butterworth filter. Defaults to TRUE.
#' @param print (boolean) Display temperature anomaly ggplot2 chart, defaults to TRUE. Use FALSE to not display chart.
#'
#' @return Invisibly returns a ggplot2 object with temperature anomaly chart
#'
#' @details `plot_temp2k` invisibly returns a ggplot2 object with a pre-defined temperature anomaly chart using data from `get_temp2k`.
#' By default the chart is also displayed. Users may further modify the output ggplot2 chart.
#'
#' @import ggplot2
#'
#' @examples
#' \donttest{
#' # Fetch temperature anomaly:
#' anomaly <- get_temp2k()
#' #
#' # Plot output using package's built-in ggplot2 defaults
#' plot_temp2k(anomaly)
#'
#' # Or just call plot_temp2k(), which defaults to get_temp2k() dataset
#' plot_temp2k()
#'
#' p <- plot_temp2k(anomaly, print = FALSE)
#' # Modify plot such as: p + ggplot2::labs(title='Temperature Increase in the Common Era') }
#'
#' @author Hernando Cortina, \email{hch@@alum.mit.edu}
#'
#' @export

plot_temp2k <- function(dataset = get_temp2k(), instrumental = TRUE,
                        filtered = TRUE, print = TRUE) {

  if (is.null(dataset)) return(invisible(NULL))

  if (filtered) dataset <- filter(dataset, substr(measure,1, 8) == "filtered") else
                dataset <- filter(dataset, substr(measure,1, 8) != "filtered")

  if (!instrumental) dataset <- filter(dataset, measure != "filtered_instrumental" & measure != "instrumental")

plot <- ggplot(dataset, aes(x=year, y=value, color=measure)) +geom_line(aes(color=measure), linewidth=1.05) + theme_bw(base_size=12) +  scale_y_continuous(n.breaks = 10, minor_breaks = NULL) +
  scale_color_manual(name=NULL, values=c('dodgerblue2','firebrick1')) +theme(legend.position = c(0.3, 0.85),legend.background=element_blank()) +
  labs(title='Global Common Era Temperature Reconstruction', subtitle='Global surface temperature relative to 1961-1990 mean', x='Year C.E.',
       y='Temperature Anomaly (C\U00B0)', caption='Source: PAGES2k Consortium, NOAA\nhttps://www.ncei.noaa.gov/access/paleo-search/study/26872')

if (print) suppressWarnings( print(plot) )
invisible(plot)
}