【问题标题】:How to subset consecutive rows if they meet a condition如果满足条件,如何对连续行进行子集化
【发布时间】:2015-12-09 15:08:02
【问题描述】:

我正在使用 R 来分析一些包含每日最高和最低温度值的时间序列(1951-2013)。数据结构如下:

YEAR MONTH  DAY     MAX    MIN
1985     1    1    22.8    9.4
1985     1    2    28.6   11.7
1985     1    3    24.7   12.2
1985     1    4    17.2    8.0
1985     1    5    17.9    7.6
1985     1    6    17.7    8.1

我需要根据以下定义找到热浪的频率:连续三天或更长时间,每日最高和最低温度超过研究中所有天的最高和最低温度的 90%期间。

基本上,当最高和最低温度超过阈值时,我想对那些连续天(三天或更长时间)进行子集化。输出将是这样的:

YEAR MONTH   DAY     MAX     MIN
1989     7    18    45.0    23.5
1989     7    19    44.2    26.1
1989     7    20    44.7    24.4
1989     7    21    44.6    29.5
1989     7    24    44.4    31.6
1989     7    25    44.2    26.7
1989     7    26    44.5    25.0
1989     7    28    44.8    26.0
1989     7    29    44.8    24.6
1989     8    19    45.0    24.3
1989     8    20    44.8    26.0
1989     8    21    44.4    24.0
1989     8    22    45.2    25.0

我已尝试以下方法将我的完整数据集子集为仅超过 90% 温度的天数:

HW<- subset(Mydata, Mydata$MAX >= (quantile(Mydata$MAX,.9)) &
                    Mydata$MIN >= (quantile(Mydata$MIN,.9)))

但是,我陷入了如何仅对满足条件的连续天进行子集化的方法。

【问题讨论】:

    标签: r subset


    【解决方案1】:

    data.table 的方法与@jlhoward 的方法略有不同(使用相同的数据):

    library(data.table)
    
    setDT(df)
    df[, hotday := +(MAX>=44.5 & MIN>=24.5)
       ][, hw.length := with(rle(hotday), rep(lengths,lengths))
         ][hotday == 0, hw.length := 0]
    

    这会生成一个带有热波长变量 (hw.length) 的数据表,而不是特定热波长的 TRUE/FALSE 变量:

    > df
        YEAR MONTH DAY  MAX  MIN hotday hw.length
     1: 1989     7  18 45.0 23.5      0         0
     2: 1989     7  19 44.2 26.1      0         0
     3: 1989     7  20 44.7 24.4      0         0
     4: 1989     7  21 44.6 29.5      1         1
     5: 1989     7  22 44.4 31.6      0         0
     6: 1989     7  23 44.2 26.7      0         0
     7: 1989     7  24 44.5 25.0      1         3
     8: 1989     7  25 44.8 26.0      1         3
     9: 1989     7  26 44.8 24.6      1         3
    10: 1989     7  27 45.0 24.3      0         0
    11: 1989     7  28 44.8 26.0      1         1
    12: 1989     7  29 44.4 24.0      0         0
    13: 1989     7  30 45.2 25.0      1         1
    

    【讨论】:

      【解决方案2】:

      我可能在这里遗漏了一些东西,但我没有看到预先设置子集的意义。如果您有每天的数据,按时间顺序排列,您可以使用运行长度编码(请参阅rle(...) 函数的文档)。

      在此示例中,我们创建了一个人工数据集并将“热浪”定义为 MAX >= 44.5 和 MIN >= 24.5。那么:

      # example data set
      df <- data.frame(YEAR=1989, MONTH=7, DAY=18:30, 
                       MAX=c(45, 44.2, 44.7, 44.6, 44.4, 44.2, 44.5, 44.8, 44.8, 45, 44.8, 44.4, 45.2),
                       MIN=c(23.5, 26.1, 24.4, 29.5, 31.6, 26.7, 25, 26, 24.6, 24.3, 26, 24, 25))
      
      r <- with(with(df, rle(MAX>=44.5 & MIN>=24.5)),rep(lengths,lengths))
      df$heat.wave <- with(df,MAX>=44.5&MIN>=24.5) & (r>2)
      df
      #    YEAR MONTH DAY  MAX  MIN heat.wave
      # 1  1989     7  18 45.0 23.5     FALSE
      # 2  1989     7  19 44.2 26.1     FALSE
      # 3  1989     7  20 44.7 24.4     FALSE
      # 4  1989     7  21 44.6 29.5     FALSE
      # 5  1989     7  22 44.4 31.6     FALSE
      # 6  1989     7  23 44.2 26.7     FALSE
      # 7  1989     7  24 44.5 25.0      TRUE
      # 8  1989     7  25 44.8 26.0      TRUE
      # 9  1989     7  26 44.8 24.6      TRUE
      # 10 1989     7  27 45.0 24.3     FALSE
      # 11 1989     7  28 44.8 26.0     FALSE
      # 12 1989     7  29 44.4 24.0     FALSE
      # 13 1989     7  30 45.2 25.0     FALSE
      

      这将创建一个列 heat.wave,如果当天有热浪,则为 TRUE。如果您只需要提取 hw 天数,请使用

      df[df$heat.wave,]
      #   YEAR MONTH DAY  MAX  MIN heat.wave
      # 7 1989     7  24 44.5 25.0      TRUE
      # 8 1989     7  25 44.8 26.0      TRUE
      # 9 1989     7  26 44.8 24.6      TRUE
      

      【讨论】:

        【解决方案3】:

        您的问题实际上归结为在您的子集数据集中查找连续 3 天以上的分组,删除所有剩余数据。

        让我们考虑一个示例,我们希望保留一些行并删除其他行:

        dat <- data.frame(year = 1989, month=c(6, 7, 7, 7, 7, 7, 8, 8, 8, 10, 10), day=c(12, 11, 12, 13, 14, 21, 5, 6, 7, 12, 13))
        dat
        #    year month day
        # 1  1989     6  12
        # 2  1989     7  11
        # 3  1989     7  12
        # 4  1989     7  13
        # 5  1989     7  14
        # 6  1989     7  21
        # 7  1989     8   5
        # 8  1989     8   6
        # 9  1989     8   7
        # 10 1989    10  12
        # 11 1989    10  13
        

        我已排除温度数据,因为我假设我们已经使用您问题中的代码将超过 90% 的天数进行了子集化。

        在此数据集中,7 月有 4 天的热浪,8 月有 3 天的热浪。第一步是将数据转换为日期对象并计算连续观察之间的天数(我假设数据已经按天排序):

        dates <- as.Date(paste(dat$year, dat$month, dat$day, sep="-"))
        (dd <- as.numeric(difftime(tail(dates, -1), head(dates, -1), units="days")))
        # [1] 29  1  1  1  7 15  1  1 66  1
        

        我们已经接近了,因为现在我们可以看到有 1 天的多个日期间隔的时间段 - 这些是我们想要抓取的时间段。我们可以使用rle 函数来分析数字 1 的游程,只保留长度为 2 或以上的游程:

        (valid.gap <- with(rle(dd == 1), rep(values & lengths >= 2, lengths)))
        # [1] FALSE  TRUE  TRUE  TRUE FALSE FALSE  TRUE  TRUE FALSE FALSE
        

        最后,我们可以将数据集子集为热浪中 1 天日期间隔两侧的日期:

        dat[c(FALSE, valid.gap) | c(valid.gap, FALSE),]
        #   year month day
        # 2 1989     7  11
        # 3 1989     7  12
        # 4 1989     7  13
        # 5 1989     7  14
        # 7 1989     8   5
        # 8 1989     8   6
        # 9 1989     8   7
        

        【讨论】:

          【解决方案4】:

          一种简单的方法,不是完全矢量化的..

          # play data
          year <- c("1960")
          month <- c(rep(1,30), rep(2,30), rep(3,30))
          day <- rep(1:30,3)
          maxT <- round(runif(90, 20, 22),1)
          minT <- round(runif(90, 10, 12),1)
          
          df <- data.frame(year, month, day, maxT, minT)
          
          # target and tricky data...
          df[1:3, 4] <- 30
          df[1:4, 5] <- 14
          df[10:13, 4] <- 30
          df[10:11, 5] <- 14
          
          # limits
          df$maxTope <- df$maxT - quantile(df$maxT,0.9)
          df$minTope <- df$minT - quantile(df$minT,0.9)
          
          # define heat day
          df$heat <- ifelse(df$maxTope > 0 & df$minTope >0, 1, 0)
          
          # count heat day2
          for(i in 2:dim(df)[1]){ 
              df$count[1] <- ifelse(df$heat[1] == 1, 1, 0)
              df$count[i] <- ifelse(df$heat[i] == 1, df$count[i-1]+1, 0)
          }
          
          # select last day of heat wave (and show the number of days in $count)
          df[which(df$count >= 3),]
          

          【讨论】:

            【解决方案5】:

            这是一个快速的小解决方案:

            is_High_Temp <- ((quantile(Mydata$MAX,.9)) &
                                Mydata$MIN >= (quantile(Mydata$MIN,.9)))
            start_of_a_series <- c(T,is_High_Temp[-1] != is_High_Temp[-length(x)]) # this is the tricky part
            series_number <- cumsum(start_of_a_series) 
            series_length <- ave(series_number,series_number,FUN=length())
            is_heat_wave  <-  series_length >= 3 & is_High_Temp 
            

            【讨论】:

              【解决方案6】:

              使用 dplyr 的解决方案,同样使用 rle()

              library(dplyr)
              
              cond <- expr(MAX >= 44.5 & MIN >= 24.5)
              
              df %>% 
                mutate(heatwave = 
                         rep(rle(!!cond)$values & rle(!!cond)$lengths >= 3, 
                             rle(!!cond)$lengths)) %>%
                filter(heatwave)
              
              #>   YEAR MONTH DAY  MAX  MIN heatwave
              #> 1 1989     7  24 44.5 25.0     TRUE
              #> 2 1989     7  25 44.8 26.0     TRUE
              #> 3 1989     7  26 44.8 24.6     TRUE
              

              reprex package (v0.3.0) 于 2020-05-16 创建

              数据

              #devtools::install_github("alistaire47/read.so")
              df <- read.so::read.so("YEAR MONTH   DAY     MAX     MIN
              1989     7    18    45.0    23.5
              1989     7    19    44.2    26.1
              1989     7    20    44.7    24.4
              1989     7    21    44.6    29.5
              1989     7    24    44.4    31.6
              1989     7    25    44.2    26.7
              1989     7    26    44.5    25.0
              1989     7    28    44.8    26.0
              1989     7    29    44.8    24.6
              1989     8    19    45.0    24.3
              1989     8    20    44.8    26.0
              1989     8    21    44.4    24.0
              1989     8    22    45.2    25.0")
              

              【讨论】:

                猜你喜欢
                • 1970-01-01
                • 1970-01-01
                • 1970-01-01
                • 1970-01-01
                • 1970-01-01
                • 2020-07-26
                • 2022-07-22
                • 1970-01-01
                • 2017-06-23
                相关资源
                最近更新 更多