通过合并范围中的R - 应用循环

我在这里发表问题：Matched Range Merge in R有关合并在一个文件中落入一个范围在第二个文件基于多项两个文件。到目前为止，我还没有成功将代码拼凑在一起来实现这一点。我遇到的问题是我正在使用的代码逐行比较文件。这是一个问题，因为1）一个文件比另一个文件长得多，并且2）我需要较短文件中的行扫描长文件中的每个范围对 - 不仅仅是同一行中的范围。通过合并范围中的R - 应用循环

我一直在使用原始问题中发布的函数，我觉得应该有一种方法可以将它应用于比较第一个文件中的每一行与第二个文件中的每一行的更一般的循环，但我还没有弄明白。如果有人有任何建议，我将不胜感激。

****编辑。

数据的性质是这样的：每个范围内不一定是唯一的，但大多数都是。它们的尺寸也不相同，有些完全落入其他尺寸。因此，findInterval会产生错误，因为无法对范围进行排序以便处于“非降序”顺序。

以下是前6行，每行数据帧：

file1test <- data.frame(SNP=c("rs2343", "rs211", "rs754", "rs854", "rs343", "rs626"), BP=c(860269, 369640, 861822, 367934, 706940, 717244)) 


file2 <- data.frame(Gene=c("E613", "E92", "E49", "E3543", "E11", "E233"), BP_start=c(367640, 621059, 721320, 860260, 861322, 879584), BP_end = c(368634, 622053, 722513, 879955, 879533, 894689))

所以，你可以看到，5日线的范围位于4号线的范围内，并且从第一文件中的两个单核苷酸多态性在第四行的范围内，但只有一行落在第二行的范围内。

第一个文件，其中包含的SNP，只有〜400行。但是包含范围的第二个文件大约有20K。我想作为输出产生的是一个数据框，其中包含来自第一个文件（SNP）的行，并且BP位于第二个文件中的BP范围内。如果SNP分为两个范围，那么它会出现两次等

来源

2012-08-09 mfk534

Was [this]（http://stackoverflow.com/q/10474590/324364）问题根本没有帮助？当你找到解决方案时，你确实看到了，不是吗？ – joran 2012-08-09 21:42:34

谢谢，我确实在我最初的搜索中查看了这个问题，但是我的范围并不一定是唯一的，这带来了额外的问题。 – mfk534 2012-08-09 22:57:35

好的。然后，您必须提供一个可重复的示例，以清楚表明'findInterval'不起作用的方式来说明您的数据特征。（当然这也解释了你想如何合并才能工作，因为如果范围不是唯一的，合并的定义并不明确。 – joran 2012-08-09 22:59:39

的GenomicRanges包在Bioconductor的是专为这种类型的操作。阅读与你的数据，例如，read.delim使

con <- textConnection("SNP  BP 
rs064 12292 
rs319 345367 
rs285 700042") 
snps <- read.delim(con, head=TRUE, sep="") 

con <- textConnection("Gene BP_start BP_end 
E613 345344  363401 
E92  694501  705408 
E49  362370  368340") ## missing trailing digit on BP_end?? 
genes <- read.delim(con, head=TRUE, sep="")

然后创建“IRanges”出来的每个

library(IRanges) 
isnps <- with(snps, IRanges(BP, width=1, names=SNP)) 
igenes <- with(genes, IRanges(BP_start, BP_end, names=Gene)

（注意坐标系，IRanges预计开始和结束被列入在范围内;同样，当end = start - 1时end> =开始期望0宽度范围）。然后找到的SNP重叠的基因（“主体”）

olaps <- findOverlaps(isnps, igenes)

两个SNP的重叠

> queryHits(olaps) 
[1] 2 3

（在IRanges术语“查询”）和它们重叠在第一和第二基因

> subjectHits(olaps) 
[1] 1 2

如果查询重叠的多个基因，它会一直重复queryHits（反之亦然）。然后，您可以合并您的数据帧作为

> cbind(snps[queryHits(olaps),], genes[subjectHits(olaps),]) 
    SNP  BP Gene BP_start BP_end 
2 rs319 345367 E613 345344 363401 
3 rs285 700042 E92 694501 705408

通常基因和SNP有染色体链（“+”，“ - ”，或“*”，以表明链并不重要）的信息，和你” d想要在这些背景下做重叠;而不是创建“IRanges”的情况下，你会创建“农庄”（基因组范围）以及随后的簿记将采取对您

library(GenomicRanges) 
isnps <- 
    with(snps, GRanges("chrA", IRanges(BP, width=1, names=SNP), "*") 
igenes <- 
    with(genes, GRanges("chrA", IRanges(BP_start, BP_end, names=Gene), "+"))

来源

2012-08-09 23:27:47

Dr. Morgan - 感谢您向我介绍'GenomicRanges'和Bioconductor。非常酷这些工具在我的工作中将会非常有用，我很感激！ – mfk534 2012-08-10 00:29:03

我相信你要求的是conditional join护理。它们在SQL中很容易，并且sqldf包使用SQL查询R中的数据帧变得很容易。

只需选择一个版本，具体取决于您希望如何处理不匹配的SNP。

内加入版本：

> sqldf("select * from file1test f1 inner join file2 f2 
+  on (f1.BP > f2.BP_start and f1.BP<= f2.BP_end) ")

输出：

 SNP  BP Gene BP_start BP_end 
1 rs2343 860269 E3543 860260 879955 
2 rs754 861822 E3543 860260 879955 
3 rs754 861822 E11 861322 879533 
4 rs854 367934 E613 367640 368634 
>

左加入版本：

> sqldf("select * from file1test f1 left join file2 f2 
+  on (f1.BP > f2.BP_start and f1.BP<= f2.BP_end) ")

输出：

 SNP  BP Gene BP_start BP_end 
1 rs2343 860269 E3543 860260 879955 
2 rs211 369640 <NA>  NA  NA 
3 rs754 861822 E3543 860260 879955 
4 rs754 861822 E11 861322 879533 
5 rs854 367934 E613 367640 368634 
6 rs343 706940 <NA>  NA  NA 
7 rs626 717244 <NA>  NA  NA 
>

注意，你可能需要您放置=，如果它是相当重要的BP将下降对于其中英国石油公司BP_start或BP_end完全匹配的情况下哪个组要小心。

来源

2012-08-10 04:09:21

感谢您介绍sqldf，Tommy！ – mfk534 2012-08-10 16:36:36

通过合并范围中的R - 应用循环

回答

相关问题