scores.grid.time: Scores for gridded reference data with a varying time...
In amber: Automated Model Benchmarking R Package

Description Usage Arguments Value Examples

This function compares model output against remote-sensing based reference data that vary in time. The performance of a model is expressed through scores that range from zero to one, where increasing values imply better performance. These scores are computed in five steps: (i) computation of a statistical metric, (ii) nondimensionalization, (iii) conversion to unit interval, (iv) spatial integration, and (v) averaging scores computed from different statistical metrics. The latter includes the bias, root-mean-square error, phase shift, inter-annual variability, and spatial distribution. The corresponding equations are documented in amber-package.

scores.grid.time(long.name, nc.mod, nc.ref, mod.id, ref.id, unit.conv.mod,
  unit.conv.ref, variable.unit, score.weights = c(1, 2, 1, 1, 1),
  outlier.factor = 1000, irregular = FALSE,
  my.projection = "+proj=ob_tran +o_proj=longlat +o_lon_p=83. +o_lat_p=42.5 +lon_0=263.",
  numCores = 2, timeInt = "month", outputDir = FALSE, myLevel = 1,
  variable.name = FALSE, phaseMinMax = "phaseMax")

`long.name`	A string that gives the full name of the variable, e.g. 'Gross primary productivity'
`nc.mod`	A string that gives the path and name of the netcdf file that contains the model output, e.g. '/home/model_gpp.nc'
`nc.ref`	A string that gives the path and name of the netcdf file that contains the reference data output, e.g. '/home/reference_gpp.nc'
`mod.id`	A string that identifies the source of the reference data set, e.g. 'CLASSIC'
`ref.id`	A string that identifies the source of the reference data set, e.g. 'MODIS'
`unit.conv.mod`	A number that is used as a factor to convert the unit of the model data, e.g. 86400
`unit.conv.ref`	A number that is used as a factor to convert the unit of the reference data, e.g. 86400
`variable.unit`	A string that gives the final units using LaTeX notation, e.g. 'gC m$^-2$ day$^-1$'
`score.weights`	R object that gives the weights of each score (S_{bias}, S_{rmse}, S_{phase}, S_{iav}, S_{dist}) that are used for computing the overall score, e.g. c(1,2,1,1,1)
`outlier.factor`	A number that is used to define outliers, e.g. 10. Plotting raster objects that contain extreme outliers lead to figures where most grid cells are presented by a single color since the color legend covers the entire range of values. To avoid this, the user may define outliers that will be masked out and marked with a red dot. Outliers are all values that exceed the interquartile range multiplied by the outlier factor defined here.
`irregular`	Logical. If TRUE the data is converted from an irregular to a regular grid. Default is FALSE.
`my.projection`	A string that defines the projection of the irregular grid
`numCores`	An integer that defines the number of cores, e.g. 2
`timeInt`	A string that gives the time interval of the model data, e.g. 'month' or 'year'
`outputDir`	A string that gives the output directory, e.g. '/home/project/study'. The output will only be written if the user specifies an output directory.
`myLevel`	A number that determines what level of the output netCDF file to use. This is relevant for files with multiple levels, which applies to soil data. By default, myLevel is set to 1.
`variable.name`	A string with the variable name, e.g. 'GPP'. If FALSE, the variable name stored in the NetCDF file will be used instead. Default is FALSE.
`phaseMinMax`	A string (either 'phaseMax' or 'phaseMin') that determines whether to assess the seasonal peak as a maximum or a minimum. The latter may be appropriate for variables that tend to be negative, such as net longwave radiation or net ecosystem exchange.

(1) A list that contains three elements. The first element is a raster stack with model data (mean, mod.mean; standard deviation; interannual-variability, mod.iav; monthly mean climatology; month of annual cycle maximum, mod.max.month), the reference data (mean, ref.mean; standard deviation; interannual-variability, ref.iav; monthly mean climatology; month of annual cycle maximum, ref.max.month), statistical metrics (bias, bias; centralized root mean square error, crmse; time difference of the annual cycle maximum, phase), and scores (bias score, bias.score; root mean square error score, rmse.score; inter-annual variability score iav.score; annual cycle score (phase.score). The second and third element of the list are spatial point data frames that give the model and reference outliers, respectively. Most of the content of the list can be plotted using plotGrid. The only exception is monthly mean climatology.

(2) NetCDF files for each of the statistical variables listed above.

(3) Three text files: (i) score values and (ii) score inputs averaged across the entire grid, and (iii) score values for each individual grid cell.

library(amber)
library(classInt)
library(doParallel)
library(foreach)
library(Hmisc)
library(latex2exp)
library(ncdf4)
library(parallel)
library(raster)
library(rgdal)
library(rgeos)
library(scico)
library(sp)
library(stats)
library(utils)
library(viridis)
library(xtable)

# (1) Global plots on a regular grid
long.name <- 'Gross primary productivity'
nc.mod <- system.file('extdata/modelRegular', 'gpp_monthly.nc', package = 'amber')
nc.ref <- system.file('extdata/referenceRegular', 'gpp_GBAF_128x64.nc', package = 'amber')
mod.id <- 'CLASSIC' # define a model experiment ID
ref.id <- 'GBAF' # give reference dataset a name
unit.conv.mod <- 86400*1000 # optional unit conversion for model data
unit.conv.ref <- 86400*1000 # optional unit conversion for reference data
variable.unit <- 'gC m$^{-2}$ day$^{-1}$' # unit after conversion (LaTeX notation)

plot.me <- scores.grid.time(long.name, nc.mod, nc.ref, mod.id, ref.id, unit.conv.mod,
unit.conv.ref, variable.unit)
plotGrid(long.name, plot.me)

# (2) Regional plots on a rotated grid
long.name <- 'Gross primary productivity'
nc.mod <- system.file('extdata/modelRotated', 'gpp_monthly.nc', package = 'amber')
nc.ref <- system.file('extdata/referenceRotated', 'gpp_GBAF_rotated.nc', package = 'amber')
mod.id <- 'CLASSIC' # define a model experiment ID
ref.id <- 'GBAF' # give reference dataset a name
unit.conv.mod <- 86400*1000 # optional unit conversion for model data
unit.conv.ref <- 86400*1000 # optional unit conversion for reference data
variable.unit <- 'gC m$^{-2}$ day$^{-1}$' # unit after conversion (LaTeX notation)
my.projection <- '+proj=ob_tran +o_proj=longlat +o_lon_p=83. +o_lat_p=42.5 +lon_0=263.'

plot.me <- scores.grid.time(long.name, nc.mod, nc.ref, mod.id, ref.id, unit.conv.mod,
unit.conv.ref, variable.unit, score.weights = c(1, 2, 1, 1, 1), outlier.factor = 1000,
irregular = TRUE, my.projection = my.projection, numCores = 2, timeInt = 'month',
outputDir = FALSE, myLevel = 1, variable.name = FALSE)

# Plot results:
irregular <- TRUE # data is on an irregular grid
my.projection <- '+proj=ob_tran +o_proj=longlat +o_lon_p=83. +o_lat_p=42.5 +lon_0=263.'
# shp.filename <- system.file('extdata/ne_50m_admin_0_countries/ne_50m_admin_0_countries.shp',
#  package = 'amber')
shp.filename <- system.file("extdata/ne_110m_land/ne_110m_land.shp", package = "amber")
my.xlim <- c(-171, 0) # longitude range that you wish to plot
my.ylim <- c(32, 78) # latitude range that you wish to plot
plot.width <- 7 # plot width
plot.height <- 3.8 # plot height

plotGrid(long.name, plot.me, irregular, my.projection,
shp.filename, my.xlim, my.ylim, plot.width, plot.height)
 #donttest