+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+Iberia01 - A new gridded dataset of daily precipitation and temperatures over Iberia.
+S. Herrera, R.M. Cardoso, P.M.M. Soares, F. Espírito-Santo, J.M. Gutiérrez
+27/05/2019
+
+
+
+Abstract
+Iberia01 is a gridded dataset of daily precipitation and mean, maximum and minimum daily temperatures over the Iberian Peninsula covering the period 1971-2015 and reaching a spatial resolution of 0.1º. This dataset has been developed based on previous works done for Spain and Portugal. The present annex includes the code, in R, to reproduce some of the figures shown in the paper Herrera et al. 2019 submitted to Earth System Science Data (ESSD). +
+
+R PACKAGE REQUIREMENTS
+The following R packages (R Core Team 2017) are required to follow the examples. Note the last Session Information Section of this manual, where the specific package versions used to build these examples are indicated.
+
+
+Installing the adequate version
+All the climate4R packages are currently under active development. Thus, to ensure the reproducibility of these examples, it is recommend to install the same versions of the packages that produced this document, even though in most cases more recent versions exist. The package versions used are indicated in the Session Information Section at the end of this document. The function install_github from package devtools (Wickham and Chang 2016) is recommended to install a particular package version from the GitHub repositories. For, instance, to install transformeR v0.0.14, it suffices with pointing to the specific version tag after the @ symbol.
Thus, in order to install the proper versions of the packages of the climate4R bundle needed to run these examples:
+devtools::install_github(c("SantanderMetGroup/transformeR@v0.0.14",
+ "SantanderMetGroup/loadeR@v1.0.9")rm(list = ls())
+options(java.parameters = "-Xmx8000m")
+library(loadeR)
+
+R package for data loading
+The R package loadeR will perform the data loading task, including authentication against the Santander User Data Gateway (UDG) server (see the installation instructions) and enabling the access to any remote dataset via OPeNDAP (more details in the loadeR’s wiki page), as well as creating and accesing datasets locally.
Prior to data access through the Santander User Data Gateway, authentication is required to access the UDG. The authentication is performed in one step:
+loginUDG(username = "jDoe", password = "*****")
+
+R package for data post-processing and plotting
+In addition, the R package transformeR (Bedia and Iturbide 2017) enables climate data transformation (plotting, aggregation, subsetting, PCA/EOF analysis …). Further details on its capabilities and installation instructions are available in the wiki page. It is seamlessly integrated with the data structures provided by loadeR.ECOMS.
library(transformeR)
+library(convertR)
+library(visualizeR)
+
+Additional packages
+
+
+Functions for Extreme Value Distributions (evd):
+In order to illustrate the analysis included in Herrera et al. 2019 the evd package should be installed to obtain the Generalized Extreme Value (GEV) distribution.
+install.packages("evd")library(evd)In order to extend the analysis shown in this document to other indices, the reader can consider the following R packages to easily estimate other climate and drought indices.
+
+
+R package for Climate Indices calculation:
+climate4R.indices is the package to compute several indices within the climate4R framework, therefore is seamlessly integrated with the climate4R data structures, and provides support for parallel computing.
devtools::install_github(c("SantanderMetGroup/transformeR", "SantanderMetGroup/climate4R.indices"))library(climate4R.indices)
+indexShow()
+
+
+?indexGrid # see the examples
+
+R package for Drought Assessment:
+drought4R is a drought assessment package for the climate4R Framework for climate data access and analysis. This package is a wrapper to the R package SPEI (Beguería and Serrano 2017), tailored to its usage with large climate model datasets within climate4R.
devtools::install_github("SantanderMetGroup/drought4R")library(drought4R)
+
+Colour Palettes:
+The package RColorBrewer is used to replicate the spectral color palette used in the paper in the correlation maps. Next, the palette veri.colors is defined, that will be used in the verification maps:
library(RColorBrewer)
+cols <- brewer.pal(n = 11, name = "Spectral")
+veri.colors <- colorRampPalette(rev(cols))
+
+E-OBS v17.0: Climate Indices
+First, we obtain the climate indices defined in the publication for the Iberia01 gridded dataset. To this aim, the dataset and geographical domain should be defined:
+eobsNcML <- "http://opendap.knmi.nl/knmi/thredds/dodsC/e-obs_0.25regular/rr_0.25deg_reg_v17.0.nc"
+di <- dataInventory(eobsNcML)
+## Variables:
+names(di)## [1] "rr"## Defining the temporal and geographical domain:
+latLim <- c(34, 44)
+lonLim <- c(-10, 6)
+years <- c(1971:2010)
+
+Daily Precipitation Climatology:
+
pr <- loadGridData(eobsNcML, "rr", dictionary = FALSE, lonLim = lonLim, latLim = latLim,
+ season = c(1:12), years = years, time = "DD", aggr.d = "sum")
+setGridUnits(pr, "mm", var = "rr")
+
+spatialPlot(climatology(pr, clim.fun = list(FUN = "mean", na.rm = TRUE)), backdrop.theme = "countries", at = seq(-0.5,6,0.5),
+ col.regions = veri.colors(51), colorkey = list(space = "bottom"), main = "Precipitation (1971-2010)")
+
+50-Year Return Value of the Annual Maximum Daily Precipitation:
+
pr.maxYear <- aggregateGrid(grid = pr, aggr.m = list(FUN = "max", na.rm = TRUE), aggr.y = list(FUN = "max", na.rm = TRUE))
+
+spatialPlot(climatology(pr.maxYear, clim.fun = list(FUN = "mean", na.rm = TRUE)), backdrop.theme = "countries", at = seq(-5,75,5),
+ col.regions = veri.colors(51), colorkey = list(space = "bottom"), main = "Max. Daily Precipitation (1971-2010)")Once the time series of annual maximum daily precipitation have been obtained, the GEV-analysis for each gridbox can be estimated:
+
pr.rv50y <- climatology(pr, clim.fun = list(FUN = "mean", na.rm = TRUE))
+
+for (i in c(1:dim(pr.maxYear$Data)[2])){
+ for (j in c(1:dim(pr.maxYear$Data)[3])){
+ auxData <- pr.maxYear$Data[,i,j]
+ auxData[which(is.infinite(auxData))] <- NA
+ if (any(!is.na(auxData))){
+ auxGEV <- fgev(auxData)
+ if ((auxGEV$estimate[3] - auxGEV$std.err[3] < 0) & (0 < auxGEV$estimate[3] + auxGEV$std.err[3])){
+ auxGEV <- fgev(auxData, shape = 0)
+ auxRV <- qgev(1-1/50, loc = auxGEV$estimate[1], scale = auxGEV$estimate[2], shape = 0)
+ }else{
+ auxRV <- qgev(1-1/50, loc = auxGEV$estimate[1], scale = auxGEV$estimate[2], shape = auxGEV$estimate[3])
+ }
+ pr.rv50y$Data[1,i,j] <- as.numeric(auxRV)
+ }
+ }
+}
+spatialPlot(pr.rv50y, backdrop.theme = "countries", at = seq(-10,200,10),
+ col.regions = veri.colors(51), colorkey = list(space = "bottom"), main = "50-Year Return Value Precipitation (1971-2010)")
+
+Wet-Day Frequency:
+
indWet <- which(pr$Data > 1)
+indDry <- which(pr$Data <= 1)
+pr$Data[indWet] <- 1
+pr$Data[indDry] <- 0
+spatialPlot(climatology(pr, clim.fun = list(FUN = "mean", na.rm = TRUE)), backdrop.theme = "countries", at = seq(-0.05,0.5,0.05),
+ col.regions = veri.colors(51), colorkey = list(space = "bottom"), main = "Wet-Day Frequency (1971-2010)")
+
+Daily Mean Temperature Climatology:
+
eobsNcML <- "http://opendap.knmi.nl/knmi/thredds/dodsC/e-obs_0.25regular/tg_0.25deg_reg_v17.0.nc"
+tg <- loadGridData(eobsNcML, "tg", dictionary = FALSE, lonLim = lonLim, latLim = latLim,
+ season = c(1:12), years = years, time = "DD", aggr.d = "sum")
+
+spatialPlot(climatology(tg, clim.fun = list(FUN = "mean", na.rm = TRUE)), backdrop.theme = "countries", at = seq(5,20,1),
+ col.regions = veri.colors(51), colorkey = list(space = "bottom"), main = "Mean Temperature (1971-2010)")
+
+50-Year Return Value of the Annual Maximum Mean Temperature:
+
tg.maxYear <- aggregateGrid(grid = tg, aggr.m = list(FUN = "max", na.rm = TRUE), aggr.y = list(FUN = "max", na.rm = TRUE))
+spatialPlot(climatology(tg.maxYear, clim.fun = list(FUN = "mean", na.rm = TRUE)), backdrop.theme = "countries", at = seq(15,40,2.5),
+ col.regions = veri.colors(51), colorkey = list(space = "bottom"), main = "Max. Daily Mean Temperature (1971-2010)")Once the time series of annual maximum daily precipitation have been obtained, the GEV-analysis for each gridbox can be estimated:
+
tg.rv50y <- climatology(tg, clim.fun = list(FUN = "mean", na.rm = TRUE))
+for (i in c(1:dim(tg.maxYear$Data)[2])){
+ for (j in c(1:dim(tg.maxYear$Data)[3])){
+ auxData <- tg.maxYear$Data[,i,j]
+ auxData[which(is.infinite(auxData))] <- NA
+ if (any(!is.na(auxData))){
+ auxGEV <- fgev(auxData[which(!is.na(auxData))], std.err = FALSE)
+ auxRV <- qgev(1-1/50, loc = auxGEV$estimate[1], scale = auxGEV$estimate[2], shape = auxGEV$estimate[3])
+ tg.rv50y$Data[1,i,j] <- as.numeric(auxRV)
+ }
+ }
+}
+spatialPlot(tg.rv50y, backdrop.theme = "countries", at = seq(25,50,2.5),
+ col.regions = veri.colors(51), colorkey = list(space = "bottom"), main = "50-Year Return Value Temperature (1971-2010)")
+
+Iberia01: Climate Indices
+In this section we show the equivalent code to make the same analysis considering the dataset Iberia01. First, we obtain the climate indices defined in the publication for the Iberia01 gridded dataset. To this aim, the dataset and geographical domain should be defined:
IberiaNcML <- "http://meteo.unican.es/tds5/dodsC/Iberia01/Iberia01_v1.0_010reg_aa_3d.ncml"
+di <- dataInventory(IberiaNcML)
+## Variables:
+names(di)
+## Defining the temporal and geographical domain:
+latLim <- c(34, 44)
+lonLim <- c(-10, 6)
+years <- c(1971:2010)
+
+Daily Precipitation Climatology:
+
pr <- loadGridData(IberiaNcML, "pr", dictionary = FALSE, lonLim = lonLim, latLim = latLim,
+ season = c(1:12), years = years, time = "DD", aggr.d = "sum")
+setGridUnits(pr, "mm", var = "pr")
+spatialPlot(climatology(pr, clim.fun = list(FUN = "mean", na.rm = TRUE)), backdrop.theme = "countries", at = seq(-0.5,6,0.5),
+ col.regions = veri.colors(51), colorkey = list(space = "bottom"), main = "Precipitation (1971-2010)")
+
+50-Year Return Value of the Annual Maximum Daily Precipitation:
+
pr.maxYear <- aggregateGrid(grid = pr, aggr.m = list(FUN = "max", na.rm = TRUE), aggr.y = list(FUN = "max", na.rm = TRUE))
+spatialPlot(climatology(pr.maxYear, clim.fun = list(FUN = "mean", na.rm = TRUE)), backdrop.theme = "countries", at = seq(-5,75,5),
+ col.regions = veri.colors(51), colorkey = list(space = "bottom"), main = "Max. Daily Precipitation (1971-2010)")Once the time series of annual maximum daily precipitation have been obtained, the GEV-analysis for each gridbox can be estimated:
+
pr.rv50y <- climatology(pr, clim.fun = list(FUN = "mean", na.rm = TRUE))
+
+for (i in c(1:dim(pr.maxYear$Data)[2])){
+ for (j in c(1:dim(pr.maxYear$Data)[3])){
+ auxData <- pr.maxYear$Data[,i,j]
+ auxData[which(is.infinite(auxData))] <- NA
+ if (any(!is.na(auxData))){
+ auxGEV <- fgev(auxData)
+ if ((auxGEV$estimate[3] - auxGEV$std.err[3] < 0) & (0 < auxGEV$estimate[3] + auxGEV$std.err[3])){
+ auxGEV <- fgev(auxData, shape = 0)
+ auxRV <- qgev(1-1/50, loc = auxGEV$estimate[1], scale = auxGEV$estimate[2], shape = 0)
+ }else{
+ auxRV <- qgev(1-1/50, loc = auxGEV$estimate[1], scale = auxGEV$estimate[2], shape = auxGEV$estimate[3])
+ }
+ pr.rv50y$Data[1,i,j] <- as.numeric(auxRV)
+ }
+ }
+}
+spatialPlot(pr.rv50y, backdrop.theme = "countries", at = seq(-10,200,10),
+ col.regions = veri.colors(51), colorkey = list(space = "bottom"), main = "50-Year Return Value Precipitation (1971-2010)")
+
+Wet-Day Frequency:
+
indWet <- which(pr$Data > 1)
+indDry <- which(pr$Data <= 1)
+pr$Data[indWet] <- 1
+pr$Data[indDry] <- 0
+spatialPlot(climatology(pr, clim.fun = list(FUN = "mean", na.rm = TRUE)), backdrop.theme = "countries", at = seq(-0.05,0.5,0.05),
+ col.regions = veri.colors(51), colorkey = list(space = "bottom"), main = "Wet-Day Frequency (1971-2010)")
+
+Daily Mean Temperature Climatology:
+
tg <- loadGridData(eobsNcML, "tas", dictionary = FALSE, lonLim = lonLim, latLim = latLim,
+ season = c(1:12), years = years, time = "DD", aggr.d = "sum")
+
+spatialPlot(climatology(tg, clim.fun = list(FUN = "mean", na.rm = TRUE)), backdrop.theme = "countries", at = seq(5,20,1),
+ col.regions = veri.colors(51), colorkey = list(space = "bottom"), main = "Mean Temperature (1971-2010)")
+
+50-Year Return Value of the Annual Maximum Mean Temperature:
+
tg.maxYear <- aggregateGrid(grid = tg, aggr.m = list(FUN = "max", na.rm = TRUE), aggr.y = list(FUN = "max", na.rm = TRUE))
+spatialPlot(climatology(tg.maxYear, clim.fun = list(FUN = "mean", na.rm = TRUE)), backdrop.theme = "countries", at = seq(15,40,2.5),
+ col.regions = veri.colors(51), colorkey = list(space = "bottom"), main = "Max. Daily Mean Temperature (1971-2010)")Once the time series of annual maximum daily precipitation have been obtained, the GEV-analysis for each gridbox can be estimated:
+
tg.rv50y <- climatology(tg, clim.fun = list(FUN = "mean", na.rm = TRUE))
+for (i in c(1:dim(tg.maxYear$Data)[2])){
+ for (j in c(1:dim(tg.maxYear$Data)[3])){
+ auxData <- tg.maxYear$Data[,i,j]
+ auxData[which(is.infinite(auxData))] <- NA
+ if (any(!is.na(auxData))){
+ auxGEV <- fgev(auxData[which(!is.na(auxData))], std.err = FALSE)
+ auxRV <- qgev(1-1/50, loc = auxGEV$estimate[1], scale = auxGEV$estimate[2], shape = auxGEV$estimate[3])
+ tg.rv50y$Data[1,i,j] <- as.numeric(auxRV)
+ }
+ }
+}
+spatialPlot(tg.rv50y, backdrop.theme = "countries", at = seq(25,50,2.5),
+ col.regions = veri.colors(51), colorkey = list(space = "bottom"), main = "50-Year Return Value Temperature (1971-2010)")
+
+REFERENCES
+-
+
- Bedia, J. and Iturbide, M. (2017). transformeR: Climate data post-processing. R package version 0.0.9. https://github.com/SantanderMetGroup/transformeR/wiki. +
- Beguería, S., Vicente‐Serrano, S.M., Reig, F., Latorre, B., (2014): Standardized precipitation evapotranspiration index (SPEI) revisited: parameter fitting, evapotranspiration models, tools, datasets and drought monitoring. Int. J. Climatol. 34, 3001–3023. https://doi.org/10.1002/joc.3887 +
- Frías M.D., Iturbide M., Manzanas R., Bedia J., Fernández J., Herrera S., Cofiño A.S., Gutiérrez J.M., (2018): An R package to visualize and communicate uncertainty in seasonal climate prediction, Environmental Modelling & Software, 99, 101-110, https://doi.org/10.1016/j.envsoft.2017.09.008. +
- Gutiérrez J. M., Bedia J., Iturbide M., Herrera S., Manzanas R., Baño-Medina J., Frías M. D., San-Martín D., Fernández J., Cofiño A. S., (2018): Climate Research Reproducibility with the Climate4R R-based framework; 8 th International Workshop on Climate Information. Available on-line here: http://www.meteo.unican.es/work/climate4r/climate4r_CI2018.pdf +
- Masih, I., Maskey, S., Mussá, F.E.F., Trambauer, P., (2014): A review of droughts on the African continent: a geospatial and long-term perspective. Hydrol. Earth Syst. Sci. 18, 3635–3649. https://doi.org/10.5194/hess-18-3635-2014 +
- MeteoSwiss (2017). easyVerification: Ensemble Forecast Verification for Large Data Sets. R package version 0.4.1. https://CRAN.R-project.org/package=easyVerification +
- R Core Team (2017). R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. URL https://www.R-project.org/. +
- Schwalm, C.R., Anderegg, W.R.L., Michalak, A.M., Fisher, J.B., Biondi, F., Koch, G., Litvak, M., Ogle, K., Shaw, J.D., Wolf, A., Huntzinger, D.N., Schaefer, K., Cook, R., Wei, Y., Fang, Y., Hayes, D., Huang, M., Jain, A., Tian, H., (2017): Global patterns of drought recovery. Nature 548, 202–205. https://doi.org/10.1038/nature23021 +
- Siegert, S. (2017). SpecsVerification: Forecast Verification Routines for Ensemble Forecasts of Weather and Climate. R package version 0.5-2. https://CRAN.R-project.org/package=SpecsVerification +
- Vicente-Serrano S.M., Beguería S., López-Moreno J.I., (2010): A Multi-scalar drought index sensitive to global warming: The Standardized Precipitation Evapotranspiration Index - SPEI. Journal of Climate 23(7), 1696-1718, DOI: 10.1175/2009JCLI2909.1.– http://digital.csic.es/handle/10261/22405. +
- Vicente-Serrano S.M., Beguería S., López-Moreno J.I., Angulo M., El Kenawy A. (2010) A global 0.5° gridded dataset (1901-2006) of a multiscalar drought index considering the joint effects of precipitation and temperature. Journal of Hydrometeorology 11(4), 1033-1043, DOI: 10.1175/2010JHM1224.1. http://digital.csic.es/handle/10261/23906. Dataset SPEIbase v.2.5 available at: http://hdl.handle.net/10261/153475 (DOI: http://dx.doi.org/10.20350/digitalCSIC/8508) +
- Wickham, H. and Chang, W. (2016). devtools: Tools to Make Developing R Packages Easier. R package version 1.12.0. https://CRAN.R-project.org/package=devtools +
–
+
+
+
+
+
+SESSION INFORMATION AND PACKAGE VERSIONS
+print(sessionInfo(package = c("loadeR", "convertR",
+ "drought4R", "RColorBrewer",
+ "transformeR", "visualizeR")))## R version 3.5.2 (2018-12-20)
+## Platform: x86_64-pc-linux-gnu (64-bit)
+## Running under: Ubuntu 14.04.3 LTS
+##
+## Matrix products: default
+## BLAS: /usr/lib/atlas-base/atlas/libblas.so.3.0
+## LAPACK: /usr/lib/lapack/liblapack.so.3.0
+##
+## locale:
+## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
+## [3] LC_TIME=es_ES.UTF-8 LC_COLLATE=en_US.UTF-8
+## [5] LC_MONETARY=es_ES.UTF-8 LC_MESSAGES=en_US.UTF-8
+## [7] LC_PAPER=es_ES.UTF-8 LC_NAME=es_ES.UTF-8
+## [9] LC_ADDRESS=es_ES.UTF-8 LC_TELEPHONE=es_ES.UTF-8
+## [11] LC_MEASUREMENT=es_ES.UTF-8 LC_IDENTIFICATION=es_ES.UTF-8
+##
+## attached base packages:
+## character(0)
+##
+## other attached packages:
+## [1] loadeR_1.4.10 convertR_0.1.2 drought4R_0.2.0
+## [4] RColorBrewer_1.1-2 transformeR_1.4.8 visualizeR_1.3.1
+##
+## loaded via a namespace (and not attached):
+## [1] maps_3.3.0 SpecsVerification_0.5-2
+## [3] dotCall64_1.0-0 sm_2.2-5.6
+## [5] assertthat_0.2.0 sp_1.3-1
+## [7] latticeExtra_0.6-28 grDevices_3.5.2
+## [9] yaml_2.2.0 pillar_1.3.1
+## [11] lattice_0.20-38 glue_1.3.1
+## [13] base_3.5.2 RcppEigen_0.3.3.5.0
+## [15] digest_0.6.18 lmomco_2.3.2
+## [17] colorspace_1.4-0 htmltools_0.3.6
+## [19] Matrix_1.2-15 plyr_1.8.4
+## [21] pkgconfig_2.0.2 raster_2.8-19
+## [23] padr_0.4.1 purrr_0.3.1
+## [25] scales_1.0.0 CircStats_0.2-6
+## [27] dtw_1.20-1 tibble_2.0.1
+## [29] proxy_0.4-23 datasets_3.5.2
+## [31] ggplot2_3.1.0 Lmoments_1.2-3
+## [33] verification_1.42 pbapply_1.4-0
+## [35] lazyeval_0.2.1 SPEI_1.7
+## [37] magrittr_1.5 crayon_1.3.4
+## [39] easyVerification_0.4.4 evaluate_0.13
+## [41] methods_3.5.2 MASS_7.3-51.1
+## [43] utils_3.5.2 tools_3.5.2
+## [45] data.table_1.12.0 stringr_1.4.0
+## [47] munsell_0.5.0 stats_3.5.2
+## [49] climate4R.indices_0.0.0 akima_0.6-2
+## [51] compiler_3.5.2 mapplots_1.5.1
+## [53] evd_2.3-3 rlang_0.3.1
+## [55] grid_3.5.2 RCurl_1.95-4.12
+## [57] rstudioapi_0.9.0 graphics_3.5.2
+## [59] goftest_1.1-1 spam_2.2-2
+## [61] bitops_1.0-6 tcltk_3.5.2
+## [63] rmarkdown_1.11 boot_1.3-20
+## [65] vioplot_0.3.0 gtable_0.2.0
+## [67] codetools_0.2-16 abind_1.4-5
+## [69] R6_2.4.0 zoo_1.8-4
+## [71] knitr_1.22 dplyr_0.8.0.1
+## [73] loadeR.java_1.1.1 udunits2_0.13
+## [75] rJava_0.9-10 stringi_1.4.3
+## [77] parallel_3.5.2 Rcpp_1.0.1
+## [79] fields_9.7 tidyselect_0.2.5
+## [81] xfun_0.5