Summarize daily data which lacks explicit grouping variable (month)

Question

I have dataframe that has 6000 locations. For each location, I have 36 years daily data of rainfall in wide format.

A sample data:

      set.seed(123)

      mat <- matrix(round(rnorm(6000*36*365), digits = 2),nrow = 6000*36, ncol = 365)
      dat <- data.table(mat)
      names(dat) <- rep(paste0("d_",1:365))

      dat$loc.id <- rep(1:6000, each = 36)
      dat$year <- rep(1980:2015, times = 6000)                     

What I want to do is for each location, generate the long term average rainfall for each month. For e.g. for loc.id = 1, mean rainfall in Jan, Feb, March... Dec.

Let' say this data is called df which is a data table

    library(dplyr)

Here's what I did:

    loc.list <- unique(dat$loc.id)
      my.list <- list() # a list to store results 

      ptm <- proc.time()

      for(i in seq_along(loc.list)){

          n <- loc.list[i]
          df1 <- dat[dat$loc.id == n,]
          df2 <- gather(df1, day, rain, -year)   # this melts the data in long format

          df3 <- df2 %>% mutate(day = gsub("d_","", day)) %>% # since the day column was in "d_1" format, I converted into integer (1,2,3..365)
                         mutate(day = as.numeric(as.character(day))) %>%  # ensure that day column is numeric. For some reasonson, some NA.s appear.
                         arrange(year,day) %>% # ensure that they are arranged in order
                         mutate(month = strptime(paste(year, day), format = "%Y %j")$mon + 1) %>% # assing each day to a month
                         group_by(year,month) %>%  # group by year and month
                         summarise(month.rain = sum(rain)) %>% # calculate for each location, year and month, total rainfall
                         group_by(month) %>% # group by month
                         summarise(month.mean = round(mean(month.rain), digits = 2)) #  calculate for each month, the long term mean

          my.list[[i]] <- df3
          }
      proc.time() - ptm

      user  system elapsed 
      1036.17    0.20 1040.68

I wanted to ask if there are more efficient and faster way to achieve this task

Answer 1

Another data.table alternative:

# change column names to month, grabbed from 365 dates of a non-leap year
setnames(dat, c(format(as.Date("2017-01-01") + 0:364, "%b"),
                "loc.id", "year"))

# melt to long format
d <- melt(dat, id.vars = c("loc.id", "year"),
          variable.name = "month", value.name = "rain")

# calculate mean rain by location and month
d2 <- d[ , .(mean_rain = mean(rain)), by = .(loc, month)]

---

This seems ~7 times faster than the answer by caw5cs. The result by Martin Morgan is in a different format though, which prevents a direct comparison of timings.

---

If you rather have unique column names in 'dat', you may use %b_%d (month-day) instead of %b only. Then use substr in by to grab the month part:

# change column names to month_day, using 365 dates of a non-leap year
setnames(dat, c(format(as.Date("2017-01-01") + 0:364, "%b_%d"),
                "loc.id", "year"))

# melt to long format
d <- melt(dat, id.vars = c("loc.id", "year"),
          variable.name = "month_day", value.name = "rain")

# calculate mean rain by location and month
d2 <- d[ , .(mean_rain = mean(rain)), by = .(loc.id, month = substr(month_day, 1, 3))]

Answer 2

Use the cryptically named rowsum() to sum daily rainfall at each site, over all years

loc.id = rep(1:6000, each = 36)
daily.by.loc = rowsum(mat, loc.id)

and use the same trick on the transposed matrix to sum by month (since there are 365 columns leap years must be ignored).

month = factor(
    months(as.Date(0:364, origin="1970-01-01")),
    levels = month.name
)
loc.by.month = rowsum(t(daily.by.loc), month)

Calculate the average by dividing by number of observations; R's column-major matrix representation and recycling rules apply. Transpose so the orientation is the same as the data.

days.per.month = tabulate(month)
ans = t(loc.by.month / (36 * days.per.month))

The result is a 6000 x 12 matrix

> dim(ans)
[1] 6000   12
> head(ans, 3)
      January     February       March       April         May         June
1  0.01554659  0.002043651 -0.02950717 -0.02700926 0.003521505 -0.011268519
2  0.04953405  0.032926587 -0.04959677  0.02808333 0.022051971  0.009768519
3 -0.01125448 -0.023343254 -0.02672939  0.04012963 0.018530466  0.035583333
          July       August   September     October    November    December
1  0.009874552 -0.030824373 -0.04958333 -0.03366487 -0.07390741 -0.07899642
2 -0.011630824 -0.003369176 -0.00100000 -0.00594086 -0.02817593 -0.01161290
3  0.031810036  0.059641577 -0.01109259  0.04646953 -0.01601852  0.03103943

in less than a second.

Answer 3

Grossly misread the question the first time, oops! Seems to be working as intended this time.

library(data.table)
set.seed(123)

mat <- matrix(round(rnorm(6000*36*365), digits = 2),nrow = 6000*36, ncol = 365)
dat <- data.table(mat)
names(dat) <- rep(paste0("d_",1:365))

dat$loc.id <- rep(1:6000, each = 36)
dat$year <- rep(1980:2015, times = 6000)


system.time({

# convert to long format with month # as column name
date_cols <- colnames(dat)[1:365]
setnames(dat, date_cols, as.character(1:365))
dat.long <- melt(dat, measure.vars=as.character(1:365), variable="day", value="rainfall")

# R date starts at 0 for Jan 1, so we offset the day by 1
dat.long[, day := as.numeric(day) - 1]
setkey(dat.long, year, day)

# Make table for merging year/day/month
months <- CJ(year=1980:2015, day=0:365)
months[, date := as.Date(day, origin=paste(year, "-01-01", sep=""))]
months[, month := tstrsplit(date, "-")[2]]
setkey(months, year, day)

# Merge tables to get month column
dat.merge <- merge(dat.long, months)



# aggregate by location an dmonth
dat.ag <- dat.merge[, list(mean_rainfall = mean(rainfall)), by=list(loc.id, month)]
})

Yielding

  user  system elapsed
14.420   4.205  18.626

> dat.ag
       loc.id month mean_rainfall
    1:      1    01   0.015546595
    2:      2    01   0.049534050
    3:      3    01  -0.011254480
    4:      4    01  -0.019453405
    5:      5    01   0.005860215
   ---
71996:   5996    12   0.027407407
71997:   5997    12   0.020334237
71998:   5998    12   0.043360434
71999:   5999    12  -0.006856369
72000:   6000    12   0.040542005