根据截止将区域分成更小的区域
Break region into smaller regions based on cutoff
我假设这是一个有点简单的编程问题,但我一直在努力解决它。主要是因为我不知道该用什么词,也许吧?
给定一组"ranges"(形式为1-一组数字,2-IRanges,或3-GenomicRanges),我想把它分成一组更小的范围。
示例开头:
Chr Start End
1 1 10000
2 1 5000
休息时间示例:2000
新数据集:
Chr Start End
1 1 2000
1 2001 4000
1 4001 6000
1 6001 8000
1 8001 10000
2 1 2000
2 2001 4000
2 4001 5000
我在 R 中做这个。我知道我可以简单地用 seq
生成这些,但我希望能够基于 list/df 个区域而不是每次我有新的区域列表时都必须手动执行此操作。
这是我使用 seq 制作的示例:
给定 22 条染色体,遍历它们并将每条染色体分解成碎片
# initialize df
Regions <- data.frame(Chromosome = c(), Start = c(), End = c())
# for each row, do the following
for(i in 1:nrow(Chromosomes)){
# create a sequence from the minimum start to the max end by some value
breks <- seq(min(Chromosomes$Start[Chromosomes$Chromosome == i]), max(Chromosomes$End[Chromosomes$Chromosome == i]), by=2000000)
# put this into a dataframe
database <- data.frame(Chromosome = i, Start = breks, End = c(breks[2:length(breks)]-1, max(Chromosomes$End[Chromosomes$Chromosome == i])))
# bind with what we already have
Regions <- rbind(Regions, database)
rm(database)
}
这很好用,我想知道包中是否已经内置了一些东西可以作为单行代码执行此操作,或者更灵活,因为它有其局限性。
使用 R / Bioconductor package GenomicRanges,这是您的初始范围
library(GenomicRanges)
rngs = GRanges(1:2, IRanges(1, c(10000, 5000)))
然后创建一个跨基因组的滑动 window,首先作为列表生成(每个染色体一组图块),然后根据您在问题中的格式取消列出
> windows = slidingWindows(rngs, width=2000, step=2000)
> unlist(windows)
GRanges object with 8 ranges and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] 1 [ 1, 2000] *
[2] 1 [2001, 4000] *
[3] 1 [4001, 6000] *
[4] 1 [6001, 8000] *
[5] 1 [8001, 10000] *
[6] 2 [ 1, 2000] *
[7] 2 [2001, 4000] *
[8] 2 [4001, 5000] *
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
用as(df, "GRanges")
或as(unlist(tiles), "data.frame")
.
从/强制到data.frame
在 ?"slidingWindows,GenomicRanges-method"
寻求帮助(制表符完成是你的朋友,?"slidingW<tab>
)。
尴尬的是,这似乎只在 GenomicRanges(v. 1.25.93?)的 'devel' version 中实现; tile
做了类似的事情,但在跨越 GRanges 的宽度时将范围的宽度四舍五入为大致相等。这是一个穷人的版本
windows <- function(gr, width, withMcols=FALSE) {
starts <- Map(seq, start(rngs), end(rngs), by=width)
ends <- Map(function(starts, len) c(tail(starts, -1) - 1L, len),
starts, end(gr))
seq <- rep(seqnames(gr), lengths(starts))
strand <- rep(strand(gr), lengths(starts))
result <- GRanges(seq, IRanges(unlist(starts), unlist(ends)), strand)
seqinfo(result) <- seqinfo(gr)
if (withMcols) {
idx <- rep(seq_len(nrow(gr)), lengths(starts))
mcols(result) = mcols(gr)[idx,,drop=FALSE]
}
result
}
调用为
> windows(rngs, 2000)
如果该方法有用,请考虑在 Bioconductor 上提出后续问题 support site。
我假设这是一个有点简单的编程问题,但我一直在努力解决它。主要是因为我不知道该用什么词,也许吧?
给定一组"ranges"(形式为1-一组数字,2-IRanges,或3-GenomicRanges),我想把它分成一组更小的范围。
示例开头:
Chr Start End
1 1 10000
2 1 5000
休息时间示例:2000
新数据集:
Chr Start End
1 1 2000
1 2001 4000
1 4001 6000
1 6001 8000
1 8001 10000
2 1 2000
2 2001 4000
2 4001 5000
我在 R 中做这个。我知道我可以简单地用 seq
生成这些,但我希望能够基于 list/df 个区域而不是每次我有新的区域列表时都必须手动执行此操作。
这是我使用 seq 制作的示例:
给定 22 条染色体,遍历它们并将每条染色体分解成碎片
# initialize df
Regions <- data.frame(Chromosome = c(), Start = c(), End = c())
# for each row, do the following
for(i in 1:nrow(Chromosomes)){
# create a sequence from the minimum start to the max end by some value
breks <- seq(min(Chromosomes$Start[Chromosomes$Chromosome == i]), max(Chromosomes$End[Chromosomes$Chromosome == i]), by=2000000)
# put this into a dataframe
database <- data.frame(Chromosome = i, Start = breks, End = c(breks[2:length(breks)]-1, max(Chromosomes$End[Chromosomes$Chromosome == i])))
# bind with what we already have
Regions <- rbind(Regions, database)
rm(database)
}
这很好用,我想知道包中是否已经内置了一些东西可以作为单行代码执行此操作,或者更灵活,因为它有其局限性。
使用 R / Bioconductor package GenomicRanges,这是您的初始范围
library(GenomicRanges)
rngs = GRanges(1:2, IRanges(1, c(10000, 5000)))
然后创建一个跨基因组的滑动 window,首先作为列表生成(每个染色体一组图块),然后根据您在问题中的格式取消列出
> windows = slidingWindows(rngs, width=2000, step=2000)
> unlist(windows)
GRanges object with 8 ranges and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] 1 [ 1, 2000] *
[2] 1 [2001, 4000] *
[3] 1 [4001, 6000] *
[4] 1 [6001, 8000] *
[5] 1 [8001, 10000] *
[6] 2 [ 1, 2000] *
[7] 2 [2001, 4000] *
[8] 2 [4001, 5000] *
-------
seqinfo: 2 sequences from an unspecified genome; no seqlengths
用as(df, "GRanges")
或as(unlist(tiles), "data.frame")
.
在 ?"slidingWindows,GenomicRanges-method"
寻求帮助(制表符完成是你的朋友,?"slidingW<tab>
)。
尴尬的是,这似乎只在 GenomicRanges(v. 1.25.93?)的 'devel' version 中实现; tile
做了类似的事情,但在跨越 GRanges 的宽度时将范围的宽度四舍五入为大致相等。这是一个穷人的版本
windows <- function(gr, width, withMcols=FALSE) {
starts <- Map(seq, start(rngs), end(rngs), by=width)
ends <- Map(function(starts, len) c(tail(starts, -1) - 1L, len),
starts, end(gr))
seq <- rep(seqnames(gr), lengths(starts))
strand <- rep(strand(gr), lengths(starts))
result <- GRanges(seq, IRanges(unlist(starts), unlist(ends)), strand)
seqinfo(result) <- seqinfo(gr)
if (withMcols) {
idx <- rep(seq_len(nrow(gr)), lengths(starts))
mcols(result) = mcols(gr)[idx,,drop=FALSE]
}
result
}
调用为
> windows(rngs, 2000)
如果该方法有用,请考虑在 Bioconductor 上提出后续问题 support site。