时间戳

首页 > 解决方案 > 将基因组区域转换为 R 数据框或 GenomicRanges 对象中的基因组位置

r - 将基因组区域转换为 R 数据框或 GenomicRanges 对象中的基因组位置

问题描述

我有一个数据框，其中包含一些基因组区间及其在几个样本中的相应覆盖率：

       sample1   sample2   sample3
1:1-3    30        NA         NA
1:1-4    NA        40         35
1:4-5    35        NA         NA
1:5-7    NA        50         50
1:6-7    60        NA         NA

我想获得相同的数据框，但对于基因组位置：

       sample1    sample2     sample3
1:1      30         40          35
1:2      30         40          35
1:3      30         40          35
1:4      35         40          35
1:5      35         50          50
1:6      60         50          50
1:7      60         50          50

你知道我怎么能得到这个吗？（我也尝试过在 GenomicRanges 对象中转换数据框，但我仍然不知道该怎么做）

标签： rgenomicranges

解决方案

我读到的你的数据：

tab = structure(list(sample1 = c(30L, NA, 35L, NA, 60L), sample2 = c(NA, 
40L, NA, 50L, NA), sample3 = c(NA, 35L, NA, 50L, NA)), class = "data.frame", row.names = c("1:1-3", 
"1:1-4", "1:4-5", "1:5-7", "1:6-7"))

这取决于您的数据集有多大，因此这是基于 GenomicRange 的解决方案：

library(GenomicRanges)
gr = GRanges(rownames(tab))
seq_range = range(gr)
W = width(seq_range)

COV = lapply(tab,function(i){
  i[is.na(i)] = 0
  coverage(gr,weight=i,width=W)
})

cov_samples = sapply(COV,function(i)as.matrix(i[seq_range]))
cov_samples
         sample1 sample2 sample3
[1,]      30      40      35
[2,]      30      40      35
[3,]      30      40      35
[4,]      35      40      35
[5,]      35      50      50
[6,]      60      50      50
[7,]      60      50      50

现在我们将它与坐标结合起来：

final = data.frame(
seqnames=rep(as.character(seqnames(seq_range)),W),
pos =  unlist(lapply(W,seq,from=1)),
cov_samples)

  seqnames pos sample1 sample2 sample3
1        1   1      30      40      35
2        1   2      30      40      35
3        1   3      30      40      35
4        1   4      35      40      35
5        1   5      35      50      50
6        1   6      60      50      50
7        1   7      60      50      50

推荐阅读