【问题标题】:loop through a file in R and select values based on a range循环遍历 R 中的文件并根据范围选择值
【发布时间】:2015-01-25 17:06:53
【问题描述】:

我有一个名为“genes”的文件,它有 4300 行,看起来像

Gene_id  chr  start   stop              
GeneA chr1  10  1000                 
GeneB chr1  2300  7000                     
GeneC chr1 10000 13577                

还有另一个名为“bases”的文件(约 100,000 行),看起来像

Chr Bases          
chr1 160           
chr1 157             
chr1 8500           
chr1 2200               

我想生成一个文件来保存每个基因的起始和终止范围内的碱基

所以输出看起来像

Chr Bases             
chr1 160            
chr1 157                

我已经尝试过这个功能,但它只给了我四次第一个条目:

methC <- apply(bases,1,function(a){
my_bases <- bases[bases[1]==genes$chr & bases[2]>=genes$start & bases[2]<=genes$stop,]
result <- my_bases[,2]
return(result)
})

>methC
# 160 160 160

所以我缺少基数 157,而 160 重复了 4 次。

如果我使用

b <- bases[which(bases[1]==genes$chr & bases[2]>=genes$start & bases[2]<=genes$stop),]
> b
 #  Chr Bases
#chr1   160

我仍然缺少 157,但可能是由于订单。

但是,如果我尝试使用具有此“哪个”功能的真实且更大的文件,我会返回一个空的 data.frame

> b
Chr        Base       
<0 rows> (or 0-length row.names)

,这就是为什么我认为函数会更好地处理大型数据集。

【问题讨论】:

标签: r


【解决方案1】:

我会使用数据表库:

library("data.table")

# Read the data from file
genes <- data.table(read.csv("genes.csv"))
bases <- data.table(read.csv("bases.csv"))

# Define the columns that are used to join the two tables
setkey(genes, chr)
setkey(bases, Chr)

# The first line joins the two tables, and the second line
# filters the result to keep only those that match and defines
# the columns that are to be shown in the output. Note the need
# to take care with capitalisation and the reserved word 'stop'
genes[bases, allow.cartesian = TRUE][
  Bases >= start & Bases <= `stop`, list(Chr = chr, Bases)]

产生这个结果:

    Chr Bases
1: chr1   160
2: chr1   157

【讨论】:

  • 您可以使用fread() 读取文件。和foverlaps() 执行重叠范围连接。
猜你喜欢
  • 1970-01-01
  • 2023-03-18
  • 2015-10-17
  • 1970-01-01
  • 2017-03-26
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 2012-12-28
相关资源
最近更新 更多