0

I have raster bricks of temperature and rainfall consisting of one layer for each month of the year. The raster brick looks like this:

> tmx
class      : RasterBrick
dimensions : 3875, 8017, 31065875, 240  (nrow, ncol, ncell, nlayers)
resolution : 0.04490319, 0.04490319  (x, y)
extent     : -180.0393, 179.9495, -90.04088, 83.95898  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
source     : /scratch/bntjoa002/chelsa_cmip5_ts/tasmax/CHELSAcmip5ts_tasmax_ACCESS1-3_rcp8                                      5_2070-2089.nc
names      : X2070.01.15, X2070.02.15, X2070.03.15, X2070.04.15, X2070.05.15, X2070.06.15,                                       X2070.07.15, X2070.08.15, X2070.09.15, X2070.10.15, X2070.11.15, X2070.12.15, X2071.01.15                                      , X2071.02.15, X2071.03.15, ...
Date       : 2070-01-15, 2089-12-15 (min, max)
varname    : air_temperature

I have been using the following code to get the average monthly values for a specified period (e.g. resulting in 12 averaged layers, one for each month, over a X year period):

new_tmx <- subset(tmx, which(getZ(tmx) >= as.Date("2070-01-15") & getZ(tmx) <= as.Date("2089-12-15"))) 
indices <- format(as.Date(names(new_tmx), format = "X%Y.%m.%d"), format = "%m")
indices <- as.numeric(indices)
Monthlytmx<- stackApply(new_tmx, indices, fun = mean)

I do this separately for the maximum temperature, minimum temperature, and precipitation bricks.

What I need now is to do the same process, but instead to output one raster layer for each year consisting of the averaged values for one-year periods (i.e. for 2070, 2071, 2072...). I've only done very simple loops in R, and I'm not sure how to deal with the time/date function in this example?

1
  • Take a look at the rts package that has functions that do exactly what you are asking (ie., apply.monthly). The package is specifically designed for raster time series analysis. Oct 22 '20 at 17:54
2

I think what you need to do is create a new index from the year in the date - you can use the "format" function on a date to get the year out:

> indices = format(getZ(tmx),"%Y")
> indices
  [1] "2070" "2070" "2070" "2070" "2070" "2070" "2070" "2070" "2070" "2070"
 [11] "2070" "2070" "2071" "2071" "2071" "2071" "2071" "2071" "2071" "2071"
 [21] "2071" "2071" "2071" "2071" "2072" "2072" "2072" "2072" "2072" "2072"
...

then stackApply with that:

> Yearlytmx<- stackApply(new_tmx, indices, fun = mean)
> Yearlytmx
class      : RasterBrick 
dimensions : 3, 4, 12, 20  (nrow, ncol, ncell, nlayers)
resolution : 0.25, 0.3333333  (x, y)
extent     : 0, 1, 0, 1  (xmin, xmax, ymin, ymax)
crs        : NA 
source     : memory
names      : index_2070, index_2071, index_2072, index_2073, index_2074, index_2075, index_2076, index_2077, index_2078, index_2079, index_2080, index_2081, index_2082, index_2083, index_2084, ... 
min values :  0.2871428,  0.3901444,  0.3483226,  0.2906245,  0.2947116,  0.3599598,  0.3503530,  0.3295298,  0.4077774,  0.3867792,  0.3916475,  0.3714683,  0.3746881,  0.3817867,  0.3694289, ... 
max values :  0.6998091,  0.6197922,  0.6100803,  0.5943091,  0.6323317,  0.6812757,  0.7028588,  0.6113572,  0.5682294,  0.6368843,  0.6010818,  0.6163263,  0.6534441,  0.5922601,  0.5980651, ... 

That brick now has 20 layers, which should be one for each year.

You can then add the year as a Z to the data:

> Yearlytmx = setZ(Yearlytmx,unique(indices),"Year")

Or use the "names", which have an "index_" attached to them:

> names(Yearlytmx)
 [1] "index_2070" "index_2071" "index_2072" "index_2073" "index_2074"
 [6] "index_2075" "index_2076" "index_2077" "index_2078" "index_2079"
[11] "index_2080" "index_2081" "index_2082" "index_2083" "index_2084"
[16] "index_2085" "index_2086" "index_2087" "index_2088" "index_2089"
3
  • This is very helpful, but I've realised I've gone about this in the wrong way. For the next step, I need to input the values into the biovars() function, which requires that each of tmx, tmn, and prec have 12 layers (for each month). So actually I need to index month+year (e.g. 012070, 022070, etc.) for each of tmn, tmx, and pre and then I need to loop the biovars function to create a biovars raster for each year using the monthly values. Would a for loop or lapply be better for this? Oct 23 '20 at 5:53
  • Actually maybe I should simply create a stack of 12 layers for each year and then loop the biovars function through the stacks. That seems the simplest. Oct 23 '20 at 6:20
  • 1
    So have I answered the question you asked, but that's not the question you meant to ask? If so, suggest you mark this as an answer and start a new question - you could edit your question but then I'd have to delete my answer or edit it if I get round to re-answering the new question...
    – Spacedman
    Oct 23 '20 at 8:31

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.