2017-04-07 6 views
2

私は以下のようなデータベースを持っています。行の間の区間の交差を比較して取得する

pos1<-c(5,15,25,40,80,5,18,22,38,84,5,16,50,92,31,50,20,30,50,70,27,50,60,50,90,20,40) 
pos2<-c(10,17,30,42,90,10,20,24,42,87,10,19,52,100,40,70,25,32,60,90,30,60,71,60,100,25,50) 
chr<-c(1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,2,2,2,2) 
n<-c(25,65,78,56,35,78,58,98,14,25,65,85,98,74,20,36,48,98,52,69,21,47,53,10,12,37,82) 
pop<-c("A","A","A","A","A","B","B","B","B","C","C","C","C","C","D","D","A","A","A","A","B","B","B","C","C","D","D") 
data<-data.frame(pos1,pos2,chr,pop,n) 

位置1および位置2は、各chrおよび母集団の間隔の開始点および終了点を設計しました。私の意図は、ポップA、B、C(Dではない)と各母集団ごとに一意の間隔との交点を求めることです。

ので、独特の間隔で私は次のような成果data.frameを持っているでしょう:

pos1.u<-c(25,50,92,20,30,27,90) 
pos2.u<-c(30,52,100,25,32,30,100) 
chr.u<-c(1,1,1,2,2,2,2) 
pop.u<-c("A","B","C","A","A","B","C") 
n.u<-c(78,98,74,48,98,21,12) 
data.u<-data.frame(pos1.u,pos2.u,chr.u,pop.u,n.u) 

そして、それらの3つの集団の間で、次のようなdata.frameを交差区間ため

pos1.c<-c(5,15,40,80,5,38,85,5,16,50,70,50,60,50) 
pos2.c<-c(10,17,42,90,10,42,87,10,19,60,90,60,71,60) 
chr.c<-c(1,1,1,1,1,1,1,1,1,2,2,2,2,2) 
pop.c<-c("A","A","A","A","B","B","B","C","C","A","A","B","B","C") 
n.c<-c(25,65,56,35,78,14,25,65,85,52,69,47,53,10) 
data.c<-data.frame(pos1.c,pos2.c,chr.c,pop.c,n.c) 

正確にこれを行うスクリプトを書く方法がわかりません。お手伝いできますか?

+0

「これらの3つの集団の交差」とはどういう意味ですか? A、B、Cで発生するpos1、pos2、およびchrの組み合わせは、5,10、および1と50,60,2という2つしかありません。 – ulfelder

+0

完全な交差を有するセグメントである。しかし、私は重複するすべてのセグメントに興味があります。おそらく、私は交差よりも重複して使うべきです...ごめんなさい。重複しているすべてのセグメントと、重複していないすべてのセグメントを見つけたいと思っています。あなたの質問をありがとう!あなたが私をさらに助けてくれることを願っています... – Cisco

+0

"重複"とは、特定の組み合わせ 'chr'と' pop'に 'pos1'で始まり' pos2'で終わるシーケンスの一部が'chr'の値は同じですが、' pop'の値は違っていますか? – ulfelder

答えて

1

私は次のコードはあなたが求めていることをしますが、それはあなたとは異なる結果を生むと思います。私が考えている不一致は、開かれた間隔と閉じた間隔の定義にあると思います。以下はエンドポイントが含まれていないと仮定しているが、これはあなたが意味するものではないと思われる(そうでなければ、(15,18)と(17,19)はオーバーラップするとはみなされない。 。したがって、以下のオープン/クローズド定義を調整する必要があるかもしれません。

pos1<-c(5,15,25,40,80,5,18,22,38,84,5,16,50,92,31,50,20,30,50,70,27,50,60,50,90,20,40) 
pos2<-c(10,17,30,42,90,10,20,24,42,87,10,19,52,100,40,70,25,32,60,90,30,60,71,60,100,25,50) 
chr<-c(1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,2,2,2,2) 
n<-c(25,65,78,56,35,78,58,98,14,25,65,85,98,74,20,36,48,98,52,69,21,47,53,10,12,37,82) 
pop<-c("A","A","A","A","A","B","B","B","B","C","C","C","C","C","D","D","A","A","A","A","B","B","B","C","C","D","D") 
data<-data.frame(pos1,pos2,chr,pop,n,stringsAsFactors = FALSE) 

library(intervals) 
data<-data[data$pop!="D",] #remove irrelevant D entries 
rownames(data) <- seq_len(nrow(data)) #reset rownames to allow for removed Ds 

#set ints as a list of intervals (as required by intervals package) 
ints <- tapply(1:nrow(data),data$pop,function(v) 
     Intervals(as.matrix(data[v,c("pos1","pos2")]), 
     closed=c(FALSE,FALSE), #this is where you adjust open/closed lower and upper ends of the intervals - TRUE means end value included 
     type="Z")) #Z is integers 
pops <- unique(data$pop) #unique values of pop 
popidx <- lapply(pops,function(x) which(data$pop==x)) #list of indices of these values in data 
names(popidx) <- pops 

#sets is a df of all pairwise combinations to check 
sets <- expand.grid(pops,pops,stringsAsFactors = FALSE) 
sets <- sets[sets$Var1!=sets$Var2,] 

olap <- lapply(1:nrow(sets),function(i) 
     interval_overlap(ints[[sets$Var1[i]]],ints[[sets$Var2[i]]])) #list of overlaps 
olap <- lapply(1:nrow(sets),function(i) { 
    df<-as.data.frame(olap[[i]],stringsAsFactors=FALSE) 
    df$pos1 <- as.numeric(rownames(df)) 
    df$pos2 <- sapply(1:nrow(df),function(j) popidx[[sets$Var2[i]]][df[j,1][[1]][1]]) 
    return(df)}) #tidy up as dfs, with correct indices in data (rather than in ints) 
olap <- do.call(rbind,olap)[,-1] #join dataframes 
olap$olaps <- !is.na(olap$pos2) #identify those with overlaps 

#group by unique pos1 and identify max and min no of overlaps with other groups 
olap <- data.frame(minoverlap=tapply(olap$olaps,olap$pos1,min),maxoverlap=tapply(olap$olaps,olap$pos1,max)) 
olap$rowno <- as.numeric(rownames(olap)) 

uniques <- data[olap$rowno[olap$maxoverlap==0],] #intervals appearing in just one pop 
commons <- data[olap$rowno[olap$minoverlap>0],] #intervals with an overlap in all other pops 
+0

よろしくお願いします!それは完全に動作します! – Cisco

+0

それは安堵です!しかし、よりエレガントな方法が必要であると確信しています。 –

+0

Andrew。質問があります。私はこのスクリプトをもう少し修正して使用しています。そして、変数 "Chr"が考慮されていない可能性があることがわかりました。同じ「Chr」(この例では1と2)内の間隔を比較するのに役立つことを願っています。私は自分自身を説明していますか? – Cisco