从有条件的组中选择一行？

Question

我有一个基因行数据集，其中基因也分组。我希望根据几个条件为每组选择 1 个基因到一个新的数据框中：

1. 如果组内其他人的得分差异>0.02，则选择得分最高的基因
2. 如果组中基因之间的得分差异direct_count的基因
3. 如果direct_count相同，选择secondary_count最高的基因
4. 如果一切都相同，请选择两个基因。

我一直在尝试在这里使用类似的问题，但在设置这么多条件时，我无法让其他示例适用于我的代码。

我的数据如下：

  Group Gene      Score     direct_count   secondary_count 
    1   AQP11    0.5566507       4               5
    1   CLNS1A   0.2811747       0               2
    1   RSF1     0.5469924       3               6
    2   CFDP1    0.4186066       1               2
    2   CHST6    0.4295135       1               3
    3   ACE      0.634           1               1
    3   NOS2     0.6345          1               1
    4   Gene1    0.1            10               20
    4   Gene2    0.68            3                1
    4   Gene3    0.7             0                1

每组基因的输出选择：

 Group Gene      Score     direct_count   secondary_count 
    1   AQP11    0.5566507       4               5       #highest direct_count
    2   CHST6    0.4295135       1               3       #highest secondary_count after matching direct_count
    3   ACE      0.634           1               1       #ACE and NOS2 have matching counts
    3   NOS2     0.6345          1               1

我目前正在尝试将dplyr::group_by() 与 if 语句一起使用。

输入数据：

structure(list(Group = c(1L, 1L, 1L, 2L, 2L, 3L, 3L, 4L, 4L, 4L), Gene = c("AQP11", 
"CLNS1A", "RSF1", "CFDP1", "CHST6", "ACE", "NOS2", "Gene1", "Gene2", "Gene3"), Score = c(0.5566507, 
0.2811747, 0.5269924, 0.4186066, 0.4295135, 0.634, 0.6345, 0.1, 0.68, 0.7), direct_count = c(4L, 
0L, 3L, 1L, 1L, 1L, 1L, 10L, 3L, 0L ), secondary_count = c(5L, 2L, 6L, 2L, 
3L, 1L, 1L, 20L, 1L, 1L)), row.names = c(NA, -10L), class = c("data.table", 
"data.frame"))

编辑： 包括 sessioninfo 并且还想在我的真实数据中注意一些行的 direct_count 和 secondary_count 不适用。

> sessionInfo()
R version 4.0.2 (2020-06-22)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 18362)

Matrix products: default

locale:
[1] LC_COLLATE=English_United Kingdom.1252 
[2] LC_CTYPE=English_United Kingdom.1252   
[3] LC_MONETARY=English_United Kingdom.1252
[4] LC_NUMERIC=C                           
[5] LC_TIME=English_United Kingdom.1252    

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] forcats_0.5.0     stringr_1.4.0     purrr_0.3.4       readr_1.4.0      
 [5] tibble_3.0.4      ggplot2_3.3.2     tidyverse_1.3.0   tidyr_1.1.2      
 [9] dplyr_1.0.2       data.table_1.13.2

loaded via a namespace (and not attached):
 [1] Rcpp_1.0.5       cellranger_1.1.0 pillar_1.4.6     compiler_4.0.2  
 [5] dbplyr_1.4.4     tools_4.0.2      jsonlite_1.7.1   lubridate_1.7.9 
 [9] lifecycle_0.2.0  gtable_0.3.0     pkgconfig_2.0.3  rlang_0.4.8     
[13] reprex_0.3.0     cli_2.1.0        DBI_1.1.0        rstudioapi_0.11 
[17] haven_2.3.1      withr_2.3.0      xml2_1.3.2       httr_1.4.2      
[21] fs_1.5.0         generics_0.0.2   vctrs_0.3.4      gtools_3.8.2    
[25] hms_0.5.3        grid_4.0.2       tidyselect_1.1.0 glue_1.4.1      
[29] R6_2.4.1         fansi_0.4.1      readxl_1.3.1     modelr_0.1.8    
[33] blob_1.2.1       magrittr_1.5     backports_1.1.10 scales_1.1.1    
[37] ellipsis_0.3.1   rvest_0.3.6      assertthat_0.2.1 colorspace_1.4-1
[41] stringi_1.5.3    munsell_0.5.0    broom_0.7.2      crayon_1.3.4  

真实数据选择的编辑问题：

structure(list(Group = c(2L, 2L, 2L, 2L, 2L), Gene = c("CFDP1", 
"CHST6", "RNU6-758P", "Gene1", "TMEM170A"), Score = c(0.551740109920502, 
0.598918557167053, 0.564491391181946, 0.567291617393494, 0.616708278656006
), direct_count = c(1, 1, 0, 0, 0), secondary_count = c(62, 
6, 1, 1, 2)), row.names = c(NA, -5L), class = c("data.table", 
"data.frame"), .internal.selfref = <pointer: 0x00000183dc6b1ef0>)

从这个组中选择Gene1，而实际上应该是CHST6，我找不到原因。

数据如下：

    
   Group Gene         Score      direct_count      secondary_count
1   2    CFDP1        0.5517401        1                  62
2   2    CHST6        0.5989186        1                   6
3   2    RNU6-758P    0.5644914        0                   1
4   2    Gene1        0.5672916        0                   1
5   2    TMEM170A     0.6167083        0                   2

CHST6 在所有基因中的direct_count 最高，

Answer 1

来自 tidyverse 的group_by 和 filter 是你的朋友。

library(dplyr)
library(tidyr)


df <- structure(list(Group = c(1L, 1L, 1L, 2L, 2L, 3L, 3L), 
Gene = c("AQP11", "CLNS1A", "RSF1", "CFDP1", "CHST6", "ACE", "NOS2"), Score = c(0.5566507,0.2811747, 0.5269924, 0.4186066, 0.4295135, 0.634, 0.6345), direct_count = c(4L, 0L, 3L, 1L, 1L, 1L, 1L), secondary_count = c(5L, 2L, 6L, 2L, 
3L, 1L, 1L)), row.names = c(NA, -7L), class = c("data.table", "data.frame"))
                 
                 
new_df <- df %>%
  #first condition
  group_by(Group) %>%
  mutate(max_score_difference = abs(max(Score)-min(Score))) %>%
  filter((max_score_difference > 0.02 & Score == max(Score)) | max_score_difference < 0.02) %>%
  # second condition
  filter(max_score_difference > 0.02 | (max_score_difference < 0.02 & direct_count == max(direct_count))) %>%
  # third condition
  filter(max_score_difference > 0.02 | (max_score_difference < 0.02 & secondary_count == max(secondary_count))) %>%
  ungroup() %>%
  #fourth condition met by max statements in filters above
  select(-max_score_difference) %>%
  data.frame()

print(new_df)

Answer 2

您的数据：

> db
# A tibble: 7 x 5
  Group `p Gene` Score direct_count secondary_count
  <dbl> <chr>    <dbl>        <dbl>           <dbl>
1     1 AQP11    0.557            4               5
2     1 CLNS1A   0.281            0               2
3     1 RSF1     0.547            3               6
4     2 CFDP1    0.419            1               2
5     2 CHST6    0.430            1               3
6     3 ACE      0.634            1               1
7     3 NOS2     0.634            1               1

现在我们首先编写一个函数，这样对于只有一个组的 df 就可以做任何事情：

gene_selection = function (df) {
  if (dim(df)[1]==1) {
    return(df)
  }
  else {
    df=arrange(df,-Score)
    if((df$Score[1]-df$Score[2])>0.2) {
      return(df[1,])
    }
    else{
      if (df$direct_count[1]!=df$direct_count[2]) {
        return(df[which.max(df$direct_count[1:2]),])
      }
      else {
        if (df$secondary_count[1]!=df$secondary_count[2]) {
          return(df[which.max(df$secondary_count[1:2]),])
        }
        else {
          return(df[1:2,])
        }
      }
    }
  }  
}

现在使用group_map在所有组上实现这个功能：

> db%>%
+     mutate(Group2=Group) %>%
+     group_by(Group2) %>%
+     group_map(~gene_selection(.)) %>%
+     bind_rows()
# A tibble: 4 x 5
  Group `p Gene` Score direct_count secondary_count
  <dbl> <chr>    <dbl>        <dbl>           <dbl>
1     1 AQP11    0.557            4               5
2     2 CHST6    0.430            1               3
3     3 NOS2     0.634            1               1
4     3 ACE      0.634            1               1

Answer 3

这是一个data.table 方法，但如果有一个更大的数据集来验证这个方法会更好。

library(data.table)

dt[,
   {
     if (.N == 1L) 
       ind = 1L
     else {
      o = sort(Score, decreasing = TRUE, index.return = TRUE)
      x = o$x
      ix = o$ix
      if (x[1L] - x[2L] > 0.05) 
        ind = ix[1L]
      else {
        search_inds = ix[which(x[1L] - x <= 0.05)]
        direct_count_sub = direct_count[search_inds]
        wm_dc = which.max(direct_count_sub)
        
        if (sum((dc_max <- direct_count_sub[wm_dc] == direct_count_sub)) == 1L)
          ind = search_inds[wm_dc]
        else {
          secondary_count_sub = secondary_count[search_inds]
          wm_sc = which.max(secondary_count_sub)
          if (sum((sc_max <- secondary_count_sub[wm_sc] == secondary_count_sub)) == 1L)
            ind = search_inds[wm_sc]
          else
            ind = search_inds[dc_max & sc_max]
        }
      }
     }
     .SD[ind]
   },
   by = Group]

##   Group   Gene     Score direct_count secondary_count
##   <int> <char>     <num>        <int>           <int>
##1:     1  AQP11 0.5566507            4               5
##2:     2  CHST6 0.4295135            1               3
##3:     3    ACE 0.6340000            1               1
##4:     3   NOS2 0.6345000            1               1

注意，data.table 受基地影响很大。我们主要按组检查您的 4 个标准，然后根据哪种情况为真对组进行细分。