基于温度阈值向量的温度数据汇总

Question

我有一个包含每日平均气温数据的数据框架，结构如下：

 'data.frame':  4666 obs. of  6 variables:
 $ Site : chr  "EB" "FFCE" "IB" "FFCE" ...
 $ Date : Date, format: "2013-01-01" "2013-01-01" "2013-01-01" "2014-01-01" ... 
 $ Day  : int  1 1 1 1 1 1 1 1 1 1 ...
 $ Year : int  2013 2013 2013 2014 2014 2014 2014 2015 2015 2015 ...
 $ Month: int  1 1 1 1 1 1 1 1 1 1 ...
 $ Temp : num  28.5 28.3 28.3 27 27.8 ...

我试图制作一个简表，它只是汇总了每个站点超过特定温度阈值(例如25摄氏度、26摄氏度)的天数。我可以像这样用dplyr手动实现这一点-

Days_above = Site_Daily_average %>% 
  group_by(Year, Site) %>% 
  summarise("23" = sum(Temp > 23), "24" = sum(Temp > 24),"25"= sum(Temp > 
25), "26"= sum(Temp > 26),  "27"= sum(Temp > 27), "28"= sum(Temp > 28), "29" 
= sum(Temp > 29),"30"= sum(Temp > 30), "31" = sum(Temp > 31), "ABOVE 
THRESHOLD" = sum(Temp > maxthreshold))%>% as.data.frame()  

它产生了这样一个表：

   Year Site  23  24  25  26  27  28  29 30 31 ABOVE THRESHOLD
1  2012   EB 142 142 142  91  64  22   0  0  0               0
2  2012 FFCE 238 238 238 210 119  64   0  0  0               0
3  2012   IB 238 238 238 218 138  87   1  0  0               0
4  2013   EB 115 115 115 115 115 109  44  0  0               0
5  2013 FFCE 223 223 216 197 148 114  94  0  0               0
6  2013   IB 365 365 365 348 299 194 135  3  0               0

..。

但是，正如您所看到的，代码相当冗长。我遇到的问题是为一系列的温度阈值产生同样的输出，即Tempclasses = Seq(16,32,0.25)。

如您所见，手动输入需要很长时间。我觉得这是一个非常简单的计算，应该使用dplyr来识别序列向量中的每个变量，执行这个函数并以完整的表格式生成一个输出。对不起，如果这是不清楚，因为我是相对较新的R，欢迎任何建议，谢谢。

Answer 1

下面是一种tidyverse方法，同样使用mtcars作为说明：

library(tidyverse)

mtcars %>% 
  mutate(threshold = cut(mpg, 
                         breaks=seq(10, max(mtcars$mpg)+10, 5), 
                         labels=seq(10, max(mtcars$mpg)+5, 5))) %>% 
  group_by(cyl, threshold) %>% 
  tally %>% 
  ungroup %>% 
  complete(threshold, nesting(cyl), fill=list(n=0)) %>% 
  arrange(desc(threshold)) %>% 
  group_by(cyl) %>% 
  mutate(N_above = cumsum(n)) %>% 
  select(-n) %>% 
  arrange(cyl, threshold)

threshold cyl N\_above 1 10 4 11 2 15 4 11 3 20 4 11 4 25 4 6 5 30 4 4 6 35 4 0 7 10 6 7 8 15 6 7 9 20 6 3 10 25 6 0 11 30 6 0 12 35 6 0 13 10 8 14 14 15 8 8 15 20 8 0 16 25 8 0 17 30 8 0 18 35 8 0

如果您希望最终数据采用宽格式，请在末尾添加一个spread，然后删除arrange。

... %>%
select(-n) %>% 
spread(threshold, N_above)

亚洲联盟10 15 20 25 30 35 1 11 11 6 4 0 2 6 7 3 0 0 3 8 14 8 0 0 0

Answer 2

正如@dww评论的那样，我们可以使用cut获得所需的格式。我已经在mtcars dataset上尝试过了，我们为mpg列创建了从10到35之间的范围，步骤为5。

df <- mtcars
df$group <- cut(df$mpg, seq(10, 35, 5))

然后我们按照cyl进行分组，并使用table来计算它们中有多少落在各自的桶中。

table(df$cyl, df$group)

#  (10,15] (15,20] (20,25] (25,30] (30,35]
#4       0       0       5       2       4
#6       0       4       3       0       0
#8       6       8       0       0       0

现在，如果某个值大于10，那么它也大于15，因此(15，20)桶中的数字也应该包括(10,15)桶中的数字，而(20，15)桶中的数字应该包括前一个数字。因此，我们需要对这个表进行逐行累加。

t(apply(table(df$cyl, df$group), 1, cumsum))

#   (10,15] (15,20] (20,25] (25,30] (30,35]
# 4       0       0       5       7      11
# 6       0       4       7       7       7
# 8       6      14      14      14      14

对于你的情况，代码会

Site_Daily_average$group <- cut(Site_Daily_average$Temp, seq(16,32,0.25))

#and then do table to get required answer.
t(apply(table(Site_Daily_average$Year,Site_Daily_average$Site, 
              Site_Daily_average$group), 1, cumsum)