如果大多数基因列表中存在，如何保留基因？

Question

我有9个基因列表，每个基因的长度为2000个基因。我想保留基因在6个或更多列表中。我不确定如何指定这一点，我一直在使用Intersect函数。

任何帮助都将受到赞赏。

Answer 1

首先，使用少量辅助功能重新创建数据：

# gene generating function as a trigram of lower case letters
gene <- function(...) {
  paste(sample(letters, 3), collapse = "")
}

# creating lists of genes
gene_lists <- lapply(seq(9), function(x) sapply(seq(2000), gene))

然后提取唯一元素：

# getting unique genes
unique_genes <- unique(unlist(gene_lists))
length(unique_genes)
[1] 10661

在这里，我们可以检查合成数据是否具有一定的冗余：

# checking if there are enough redundant genes
stopifnot(length(unique_genes) < length(unlist(gene_lists)))

然后在计数发生时迭代在唯一的基因和列表上：

# iterating over unique_genes
gene_occurence <- sapply(unique_genes, function(gene) {
  # iterating over lists
  # sum counts the total number of occurence
    sum(sapply(gene_lists, function(x) { gene %in% x }))
  })
length(gene_occurence)
[1] 10661
table(gene_occurence)
   1    2    3    4    5    6 
6017 3330 1050  231   31    2 

然后获取常见基因：

limit <- 6
common_genes <- unique_genes[which(gene_occurence >= limit)]
common_genes
[1] "ngu" "het"