Working with CFtime

Climate change models and calendars

Around the world, many climate change models are being developed (100+) under the umbrella of the World Climate Research Programme to assess the rate of climate change. Published data is generally publicly available to download for research and other (non-commercial) purposes through partner organizations in the Earth Systems Grid Federation.

The data are all formatted to comply with the CF Metadata Conventions, a set of standards to support standardization among research groups and published data sets. These conventions greatly facilitate use and analysis of the climate projections because standard processing work flows (should) work across the various data sets.

On the flip side, the CF Metadata Conventions needs to cater to a wide range of modeling requirements and that means that some of the areas covered by the standards are more complex than might be assumed. One of those areas is the temporal dimension of the data sets. The CF Metadata Conventions supports no less than nine different calendar definitions, that, upon analysis, fall into five distinct calendars (from the perspective of computation of climate projections):

• standard or gregorian: The international civil calendar that is in common use in many countries around the world, adopted by edict of Pope Gregory XIII in 1582 and in effect from 15 October of that year. The proleptic_gregorian calendar is the same as the gregorian calendar, but with validity extended to periods prior to 1582-10-15.
• julian: Adopted in the year 45 BCE, every fourth year is a leap year. Originally, the julian calendar did not have a monotonically increasing year assigned to it and there are indeed several julian calendars in use around the world today with different years assigned to them. Common interpretation is currently that the year is the same as that of the standard calendar. The julian calendar is currently 13 days behind the gregorian calendar.
• 365_day or noleap: No years have a leap day.
• 366_day or all_leap: All years have a leap day.
• 360_day: Every year has 12 months of 30 days each.

The three latter calendars are specific to the CF Metadata Conventions to reduce computational complexities of working with dates. These three, and the julian calendar, are not compliant with the standard POSIXt date/time facilities in R and using standard date/time procedures would quickly lead to problems. In the below code snippet, the date of 1949-12-01 is the datum from which other dates are calculated. When adding 43,289 days to this datum for a data set that uses the 360_day calendar, that should yield a date some 120 years after the datum:

# POSIXt calculations on a standard calendar - INCORRECT
as.Date("1949-12-01") + 43289
#> [1] "2068-06-08"

# CFtime calculation on a "360_day" calendar - CORRECT
# See below examples for details on the two functions
as_timestamp(CFtime("days since 1949-12-01", "360_day", 43289))
#> [1] "2070-02-30"

Using standard POSIXt calculations gives a result that is about 21 months off from the correct date - obviously an undesirable situation. This example is far from artificial: 1949-12-01 is the datum for all CORDEX data, covering the period 1951 - 2005 for historical experiments and the period 2006 - 2100 for RCP experiments (with some deviation between data sets), and several models used in the CORDEX set use the 360_day calendar. The 365_day or noleap calendar deviates by about 1 day every 4 years (disregarding centurial years), or about 24 days in a century. The 366_day or all_leap calendar deviates by about 3 days every 4 years, or about 76 days in a century.

The CFtime package deals with the complexity of the different calendars allowed by the CF Metadata Conventions. It properly formats dates and times (even oddball dates like 2070-02-30) and it can generate calendar-aware factors for further processing of the data.

Using CFtime to deal with calendars

Data sets that are compliant with the CF Metadata Conventions always include a datum, a specific point in time in reference to a specified calendar, from which other points in time are calculated by adding a specified offset of a certain unit. This approach is encapsulated in the CFtime package by the S4 class CFtime.

# Create a CF time object from a definition string, a calendar and some offsets
cf <- CFtime("days since 1949-12-01", "360_day", 19830:90029)
cf
#> CF datum of origin:
#>   Origin  : 1949-12-01 00:00:00
#>   Units   : days
#>   Calendar: 360_day
#> CF time series:
#>   Elements: [2005-01-01 .. 2199-12-30] (average of 1.000000 days between 70200 elements)
#>   Bounds  : not set

The CFtime() function takes a datum description (which is actually a unit - “days” - in reference to a datum - “1949-12-01”), a calendar description, and a vector of offsets from that datum. Once a CFtime instance is created its datum and calendar cannot be changed anymore. Offsets may be added.

In practice, these parameters will be taken from the data set of interest. CF Metadata Conventions require data sets to be in the NetCDF format, with all metadata describing the data set included in a single file, including the mandatory “Conventions” global attribute which should have a string identifying the version of the CF Metadata Conventions that this file adheres to (among possible others). Not surprisingly, all the pieces of interest are contained in the mandatory time dimension of the file. The process then becomes as follows, for a CMIP6 file of daily precipitation:

# Opening a data file that is included with the package and showing some attributes.
# Usually you would `list.files()` on a directory of your choice.
fn <- list.files(path = system.file("extdata", package = "CFtime"), full.names = TRUE)[1]
nc <- nc_open(fn)
attrs <- ncatt_get(nc, "")
attrs$title
#> [1] "NOAA GFDL GFDL-ESM4 model output prepared for CMIP6 update of RCP4.5 based on SSP2"

# "Conventions" global attribute must have a string like "CF-1.*" for this package to work reliably
attrs$Conventions
#> [1] "CF-1.7 CMIP-6.0 UGRID-1.0"

# Create the CFtime instance from the metadata in the file.
cf <- CFtime(nc$dim$time$units, 
             nc$dim$time$calendar, 
             nc$dim$time$vals)
cf
#> CF datum of origin:
#>   Origin  : 1850-01-01 00:00:00
#>   Units   : days
#>   Calendar: noleap
#> CF time series:
#>   Elements: [2015-01-01 12:00:00 .. 2099-12-31 12:00:00] (average of 1.000000 days between 31025 elements)
#>   Bounds  : not set

You can see from the global attribute “Conventions” that the file adheres to the CF Metadata Conventions, among others. According to the CF conventions, units and calendar are required attributes of the time dimension in the NetCDF file, and nc$dim$time$vals are the offset values, or dimnames() in R terms, for the time dimension of the data.

The above example (and others in this vignette) use the ncdf4 package. If you are using the RNetCDF package, checking for CF conventions and then creating a CFtime instance goes like this:

library(RNetCDF)
nc <- open.nc(fn)
att.get.nc(nc, -1, "Conventions")
#> [1] "CF-1.7 CMIP-6.0 UGRID-1.0"
cf <- CFtime(att.get.nc(nc, "time", "units"), 
             att.get.nc(nc, "time", "calendar"), 
             var.get.nc(nc, "time"))
cf
#> CF datum of origin:
#>   Origin  : 1850-01-01 00:00:00
#>   Units   : days
#>   Calendar: noleap
#> CF time series:
#>   Elements: [2015-01-01 12:00:00 .. 2099-12-31 12:00:00] (average of 1.000000 days between 31025 elements)
#>   Bounds  : not set

The corresponding character representations of the time series can be easily generated:

dates <- as_timestamp(cf, format = "date")
dates[1:10]
#>  [1] "2015-01-01" "2015-01-02" "2015-01-03" "2015-01-04" "2015-01-05"
#>  [6] "2015-01-06" "2015-01-07" "2015-01-08" "2015-01-09" "2015-01-10"

…as well as the full range of the time series:

range(cf)
#> [1] "2015-01-01 12:00:00" "2099-12-31 12:00:00"

Note that in this latter case, if any of the timestamps in the time series have a time that is other than 00:00:00 then the time of the extremes of the time series is also displayed. This is a common occurrence because the CF Metadata Conventions prescribe that the middle of the time period (month, day, etc) is recorded, which for months with 31 days would be something like 2005-01-15T12:00:00.

Supporting processing of climate projection data

When working with high resolution climate projection data, typically at a “day” resolution, one of the processing steps would be to aggregate the data to some lower resolution such as a dekad (10-day period), a month or a meteorological season, and then compute a derivative value such as the dekadal sum of precipitation, monthly minimum/maximum daily temperature, or seasonal average daily short-wave irradiance.

It is also possible to create factors for multiple “epochs” in one go. This greatly reduces programming effort if you want to calculate anomalies over multiple future periods. A complete example is provided in the vignette “Processing climate projection data”.

It is easy to generate the factors that you need once you have a CFtime instance prepared:

# Create a dekad factor for the whole `cf` time series that was created above
f_k <- CFfactor(cf, "dekad")
str(f_k)
#>  Factor w/ 3060 levels "2015D01","2015D02",..: 1 1 1 1 1 1 1 1 1 1 ...
#>  - attr(*, "epoch")= int -1
#>  - attr(*, "period")= chr "dekad"
#>  - attr(*, "CFtime")=Formal class 'CFtime' [package "CFtime"] with 4 slots
#>   .. ..@ datum     :Formal class 'CFdatum' [package "CFtime"] with 5 slots
#>   .. .. .. ..@ definition: chr "days since 1850-01-01"
#>   .. .. .. ..@ unit      : int 4
#>   .. .. .. ..@ origin    :'data.frame':  1 obs. of  8 variables:
#>   .. .. .. .. ..$ year  : int 1850
#>   .. .. .. .. ..$ month : num 1
#>   .. .. .. .. ..$ day   : num 1
#>   .. .. .. .. ..$ hour  : num 0
#>   .. .. .. .. ..$ minute: num 0
#>   .. .. .. .. ..$ second: num 0
#>   .. .. .. .. ..$ tz    : chr "+0000"
#>   .. .. .. .. ..$ offset: num 0
#>   .. .. .. ..@ calendar  : chr "noleap"
#>   .. .. .. ..@ cal_id    : int 4
#>   .. ..@ resolution: num 14.9
#>   .. ..@ offsets   : num [1:3060] 60230 60240 60250 60261 60271 ...
#>   .. ..@ bounds    : logi TRUE

# Create monthly factors for a baseline epoch and early, mid and late 21st century epochs
baseline <- CFfactor(cf, epoch = 1991:2020)
future <- CFfactor(cf, epoch = list(early = 2021:2040, mid = 2041:2060, late = 2061:2080))
str(future)
#> List of 3
#>  $ early: Factor w/ 12 levels "01","02","03",..: NA NA NA NA NA NA NA NA NA NA ...
#>   ..- attr(*, "epoch")= int 20
#>   ..- attr(*, "period")= chr "month"
#>  $ mid  : Factor w/ 12 levels "01","02","03",..: NA NA NA NA NA NA NA NA NA NA ...
#>   ..- attr(*, "epoch")= int 20
#>   ..- attr(*, "period")= chr "month"
#>  $ late : Factor w/ 12 levels "01","02","03",..: NA NA NA NA NA NA NA NA NA NA ...
#>   ..- attr(*, "epoch")= int 20
#>   ..- attr(*, "period")= chr "month"

For the “epoch” version, there are two interesting things to note here:

• The epochs do not have to coincide with the boundaries of the time series. In the example above, the time series starts in 2015, while the baseline epoch is from 1991. Obviously, the number of time steps from the time series that then fall within this epoch will then be reduced.
• The factor is always of the same length as the time series, with NA values where the time series values are not falling in the epoch. This ensures that the factor is compatible with the data set which the time series describes, such that functions like tapply() will not throw an error.

There are six periods defined for CFfactor():

• year, to summarize data to yearly timescales
• season, the meteorological seasons. Note that the month of December will be added to the months of January and February of the following year, so the date “2020-12-01” yields the factor value “2021S1”.
• quarter, the standard quarters of the year.
• month, monthly summaries, the default period.
• dekad, 10-day period. Each month is subdivided in dekads as follows: (1) days 01 - 10; (2) days 11 - 20; (3) remainder of the month.
• day, to summarize sub-daily data.

  new_time <- attr(f_k, "CFtime")
  new_time
#> CF datum of origin:
#>   Origin  : 1850-01-01 00:00:00
#>   Units   : days
#>   Calendar: noleap
#> CF time series:
#>   Elements: [1974-12-26 12:00:00 .. 2099-12-16 00:00:00] (average of 14.911572 days between 3060 elements)
#>   Bounds  : regular and consecutive

# Is the time series complete?
is_complete(cf)
#> [1] TRUE

# How many time units fit in a factor level?
CFfactor_units(cf, baseline)
#>  [1] 31 28 31 30 31 30 31 31 30 31 30 31

# What's the absolute and relative coverage of our time series
CFfactor_coverage(cf, baseline, "absolute")
#>  [1] 186 168 186 180 186 180 186 186 180 186 180 186
CFfactor_coverage(cf, baseline, "relative")
#>  [1] 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2

# 4 years of data on a `365_day` calendar, keep 80% of values
n <- 365 * 4
cov <- 0.8
offsets <- sample(0:(n-1), n * cov)

cf <- CFtime("days since 2020-01-01", "365_day", offsets)
cf
#> CF datum of origin:
#>   Origin  : 2020-01-01 00:00:00
#>   Units   : days
#>   Calendar: 365_day
#> CF time series:
#>   Elements: [2020-01-01 .. 2023-12-31] (average of 1.250214 days between 1168 elements)
#>   Bounds  : not set
# Note that there are about 1.25 days between observations

mon <- CFfactor(cf, "month")
CFfactor_coverage(cf, mon, "absolute")
#>  [1] 26 25 28 27 26 25 26 21 24 23 24 25 26 26 24 26 25 23 26 24 22 21 20 19 26
#> [26] 18 21 25 27 24 24 26 22 27 22 25 26 22 27 27 27 20 27 28 23 23 27 22
CFfactor_coverage(cf, mon, "relative")
#>  [1] 0.8387097 0.8928571 0.9032258 0.9000000 0.8387097 0.8333333 0.8387097
#>  [8] 0.6774194 0.8000000 0.7419355 0.8000000 0.8064516 0.8387097 0.9285714
#> [15] 0.7741935 0.8666667 0.8064516 0.7666667 0.8387097 0.7741935 0.7333333
#> [22] 0.6774194 0.6666667 0.6129032 0.8387097 0.6428571 0.6774194 0.8333333
#> [29] 0.8709677 0.8000000 0.7741935 0.8387097 0.7333333 0.8709677 0.7333333
#> [36] 0.8064516 0.8387097 0.7857143 0.8709677 0.9000000 0.8709677 0.6666667
#> [43] 0.8709677 0.9032258 0.7666667 0.7419355 0.9000000 0.7096774

CFtime and POSIXt

The CF Metadata Conventions support nine different calendars. Of these, only the standard, gregorian and proleptic_gregorian calendars are fully compatible with POSIXt. The other calendars have varying degrees of discrepancies:

• julian: Every fouth year is a leap year. Dates like 2100-02-29 and 2200-02-29 are valid.
• 365_day or noleap: No leap year exists. 2020-02-29 does not occur.
• 366_day or all_leap: All years are leap years.
• 360_day: All months have 30 days in every year. This means that 31 January, March, May, July, August, October and December never occur, while 29 and 30 February occur in every year.

Converting time series using these incompatible calendars to Dates is likely to produce problems. This is most pronounced for the 360_day calendar:

# Days in January and February
cf <- CFtime("days since 2023-01-01", "360_day", 0:59)
cf_days <- as_timestamp(cf, "date")
as.Date(cf_days)
#>  [1] "2023-01-01" "2023-01-02" "2023-01-03" "2023-01-04" "2023-01-05"
#>  [6] "2023-01-06" "2023-01-07" "2023-01-08" "2023-01-09" "2023-01-10"
#> [11] "2023-01-11" "2023-01-12" "2023-01-13" "2023-01-14" "2023-01-15"
#> [16] "2023-01-16" "2023-01-17" "2023-01-18" "2023-01-19" "2023-01-20"
#> [21] "2023-01-21" "2023-01-22" "2023-01-23" "2023-01-24" "2023-01-25"
#> [26] "2023-01-26" "2023-01-27" "2023-01-28" "2023-01-29" "2023-01-30"
#> [31] "2023-02-01" "2023-02-02" "2023-02-03" "2023-02-04" "2023-02-05"
#> [36] "2023-02-06" "2023-02-07" "2023-02-08" "2023-02-09" "2023-02-10"
#> [41] "2023-02-11" "2023-02-12" "2023-02-13" "2023-02-14" "2023-02-15"
#> [46] "2023-02-16" "2023-02-17" "2023-02-18" "2023-02-19" "2023-02-20"
#> [51] "2023-02-21" "2023-02-22" "2023-02-23" "2023-02-24" "2023-02-25"
#> [56] "2023-02-26" "2023-02-27" "2023-02-28" NA           NA

31 January is missing from the vector of Dates because the 360_day calendar does not include it and 29 and 30 February are NAs because POSIXt rejects them. This will produce problems later on when processing your data.

The general advice is therefore: do not convert CFtime objects to Date objects unless you are sure that the CFtime object uses a POSIXt-compatible calendar.

Final observations

• This package is intended to facilitate processing of climate projection data. It is not a full implementation of the CF Metadata Conventions “time” component.
• In parsing and deparsing of offsets and timestamps, data is rounded to 3 digits of precision of the unit of the datum. When using a description of time that is very different from the datum unit, this may lead to some loss of precision due to rounding errors. For instance, if milli-second precision is required, use a unit of “seconds”. The authors have no knowledge of published climate projection data that requires milli-second precision so for the intended use of the package this issue is marginal.