R：将向量划分为区间，并测试哪个整数落入哪个区间答案

【问题标题】：R: Dividing vector into intervals, and test which integer falls into what intervalR：将向量划分为区间，并测试哪个整数落入哪个区间
【发布时间】：2016-05-07 08:54:25
【问题描述】：

我有 23 条染色体及其长度

chromosome    length
1             249250621
2             243199373
3             198022430
4             191154276
5             180915260
6             171115067 
..            .........
Y             59373566

对于每条染色体，我想创建 5000 个大小相等的箱/区间。

Chr1:
bin_number    start        end
1             1            49850
2             49851        99700
....          .....        .....
5000          249200771    249250621

为此，我尝试同时使用“cut”和“cut2”。 “cut2”无法处理染色体的长度和崩溃，而 cut 为每个单独的位置提供了间隔（249250621 个间隔！）。

cut2(1:249250621, g=5000, onlycuts = TRUE)

cut(1:249250621, breaks=5000)

当我有间隔时，我想分配每个 50.000 个变体属于哪个 bin/interval。

我的数据（染色体 1）：

variant            chromosome    position
1:20000_G/A        1             20000
1:30000_C/CCCCT    1             30000
1:60000_G/T        1             60000
..............     ..            .......

我想要什么：

variant            chromosome    position    bin_number
1:20000_G/A        1             20000       1
1:30000_C/CCCCT    1             30000       1
1:60000_G/T        1             60000       2
..............     ..            .......     ...

我将不胜感激有关将我的染色体分成间隔的方法的任何建议。当我有区间时，我需要能够快速测试变体属于哪个区间的方法。

【问题讨论】：

您在寻找round(seq(1, 249250621, length.out = 5000))（第一个染色体）吗？我使用了round，因为 249250621 不是 5000 的整数倍。
用你的方法剩下的会怎样？它似乎创建了 49860 的间隔，而实际大小应该是 49850.12（或 49851 用于 621 个 bin，而 49850 用于其余 4379 个 bin）。
我以为你不需要小数（这就是我使用round() 的原因）。只需使用seq(1, 249250621, length.out = 5000) 即可获得准确的间隔
给定一个“染色体”的“长度”和它的“变体”的“位置”——即像ff = function(len, pos) findInterval(pos, seq(1, len, length.out = 5000))这样的函数——你想要各自的间隔作为输出（@987654331 @)?您在寻找更具体或更不同的东西吗？

标签： r range intervals cut

【解决方案1】：

如果我理解你的算法，你将每个染色体分成 10000 个箱子。因此，您将获得每个范围的不同尺寸。我曾经应用这个算法来创建独立于染色体的相同大小的范围。

chrSizes <- data.frame(chromosome = c("1", "2", "3", "4", "5", "6", "7", "8", "9", "10", "11", "12", "13", "14", "15", "16", "17", "18", "19", "20", "21", "22", "X", "Y"), 
                       length = c("249250621","243199373", "198022430", "191154276", "180915260", "171115067", "159138663", "146364022", "141213431", "135534747", "135006516", "133851895", "115169878", "107349540", "102531392", "90354753", "81195210", "78077248", "59128983", "63025520", "48129895", "51304566", "155270560", "59373566"), 
                       stringsAsFactors = FALSE)

sizerange <- 5000000
lastranges <- NA
h <- 0

for (i in 1:24) 
{
  thelast <- 1
  bynum <- format(sizerange, scientific = FALSE)
  chrlist <- c(paste0(chrSizes$chromosome[i],":1-",bynum))
  biggest <- chrSizes$length[i]
  while(thelast < biggest)
  {
    bynum1 <- format(as.numeric(bynum)+1, scientific = FALSE)
    bynum2 <- format(as.numeric(bynum1)+sizerange-1, scientific = FALSE)
    berria <- paste0(paste0(chrSizes$chromosome[i],":",bynum1,"-",as.character(bynum2)))
    chrlist <- c(chrlist,berria)
    thelast <- as.numeric(bynum2)+sizerange
    bynum <- format(as.numeric(bynum)+sizerange, scientific = FALSE)
  }
  azkenreg <- paste0(paste0(chrSizes$chromosome[i],":",bynum,"-",as.character(biggest)))
  chrlist <- c(chrlist,azkenreg)
  lastranges <- c(lastranges,chrlist)
}

lastranges <- lastranges[-1]

df <- data.frame(lastranges)
write.table(df,file = "fastacontigs_splited_bysize2.txt",quote = FALSE, row.names = FALSE, col.names = FALSE)

在这种情况下，结果是：

1:1-5000000
1:5000001-10000000
1:10000001-15000000
1:15000000-249250621
2:1-5000000
2:5000001-10000000
2:10000001-15000000
2:15000000-243199373
3:1-5000000
3:5000001-10000000
3:10000001-15000000
3:15000000-198022430
4:1-5000000
4:5000001-10000000
4:10000001-15000000
4:15000000-191154276
5:1-5000000
5:5000001-10000000
5:10000001-15000000
5:15000000-180915260
6:1-5000000
6:5000001-10000000
6:10000001-15000000
6:15000000-171115067
7:1-5000000
7:5000001-10000000
7:10000001-15000000
7:15000000-159138663
8:1-5000000
8:5000001-10000000
8:10000001-15000000
8:15000000-146364022
9:1-5000000
9:5000001-10000000
9:10000001-15000000
9:15000000-141213431
10:1-5000000
10:5000001-10000000
10:10000001-15000000
10:15000000-135534747
11:1-5000000
11:5000001-10000000
11:10000001-15000000
11:15000000-135006516
12:1-5000000
12:5000001-10000000
12:10000001-15000000
12:15000000-133851895
13:1-5000000
13:5000001-10000000
13:10000001-15000000
13:15000000-115169878
14:1-5000000
14:5000001-10000000
14:10000001-15000000
14:15000000-107349540
15:1-5000000
15:5000001-10000000
15:10000001-15000000
15:15000000-102531392
16:1-5000000
16:5000001-10000000
16:10000001-15000000
16:15000001-20000000
16:20000001-25000000
16:25000001-30000000
16:30000001-35000000
16:35000001-40000000
16:40000001-45000000
16:45000001-50000000
16:50000001-55000000
16:55000001-60000000
16:60000001-65000000
16:65000001-70000000
16:70000001-75000000
16:75000001-80000000
16:80000001-85000000
16:85000000-90354753
17:1-5000000
17:5000001-10000000
17:10000001-15000000
17:15000001-20000000
17:20000001-25000000
17:25000001-30000000
17:30000001-35000000
17:35000001-40000000
17:40000001-45000000
17:45000001-50000000
17:50000001-55000000
17:55000001-60000000
17:60000001-65000000
17:65000001-70000000
17:70000001-75000000
17:75000000-81195210
18:1-5000000
18:5000001-10000000
18:10000001-15000000
18:15000001-20000000
18:20000001-25000000
18:25000001-30000000
18:30000001-35000000
18:35000001-40000000
18:40000001-45000000
18:45000001-50000000
18:50000001-55000000
18:55000001-60000000
18:60000001-65000000
18:65000000-78077248
19:1-5000000
19:5000001-10000000
19:10000001-15000000
19:15000001-20000000
19:20000001-25000000
19:25000001-30000000
19:30000001-35000000
19:35000001-40000000
19:40000001-45000000
19:45000000-59128983
20:1-5000000
20:5000001-10000000
20:10000001-15000000
20:15000001-20000000
20:20000001-25000000
20:25000001-30000000
20:30000001-35000000
20:35000001-40000000
20:40000001-45000000
20:45000001-50000000
20:50000001-55000000
20:55000000-63025520
21:1-5000000
21:5000001-10000000
21:10000001-15000000
21:15000001-20000000
21:20000001-25000000
21:25000001-30000000
21:30000001-35000000
21:35000000-48129895
22:1-5000000
22:5000001-10000000
22:10000001-15000000
22:15000001-20000000
22:20000001-25000000
22:25000001-30000000
22:30000001-35000000
22:35000001-40000000
22:40000001-45000000
22:45000000-51304566
X:1-5000000
X:5000001-10000000
X:10000001-15000000
X:15000000-155270560
Y:1-5000000
Y:5000001-10000000
Y:10000001-15000000
Y:15000001-20000000
Y:20000001-25000000
Y:25000001-30000000
Y:30000001-35000000
Y:35000001-40000000
Y:40000001-45000000
Y:45000000-59373566

【讨论】：

【解决方案2】：

如果 bin 范围是恒定的，则此方法有效：

mydata <- data.frame(position = c(20000, 30000, 60000, 
                              49850, 49851, 99700, 99701))
mydata$bin <- ceiling(mydata$position / 49850)

更一般地说，如果 bin 范围不是恒定的，但您已经定义了切点，则可以使用 cut，将其指定为 breaks。

cutpoints <- c(0, 49850, 99700, 149550)
mydata$bin2 <- cut(mydata$position, breaks = cutpoints)

您可以稍作调整来编辑标签。

mydata$bin3 <- cut(mydata$position, breaks = cutpoints,
               labels = seq(length(cutpoints)-1))

【讨论】：

【解决方案3】：

感谢您的意见。我选择使用一个简单的循环来创建间隔，以确保间隔具有所需的大小。

我创建了一个包含染色体大小的 data.frame

chrSizes <- data.frame(chromosome = c("1", "2", "3", "4", "5", "6", "7", "8", "9", "10", "11", "12", "13", "14", "15", "16", "17", "18", "19", "20", "21", "22", "X", "Y"), length = c("249250621","243199373", "198022430", "191154276", "180915260", "171115067", "159138663", "146364022", "141213431", "135534747", "135006516", "133851895", "115169878", "107349540", "102531392", "90354753", "81195210", "78077248", "59128983", "63025520", "48129895", "51304566", "155270560", "59373566"), stringsAsFactors = FALSE)

然后我通过找到精确的块大小然后向下舍入来循环创建间隔的每个染色体。然后将余数用于在前多个间隔中添加一个。

numOfBins <- 10000
chrBinList <- list()
for (i in 1:24) {
  chrBins <- c()
  chrLength <- as.numeric(chrSizes[i,2])
  chunkSize <- floor(chrLength/numOfBins)
  remainder <- chrLength %% chunkSize
  counter <- 1

  # Adding remainder to the first intervals
  for (j in 1:(remainder-1)) {
    chrBins <- c(chrBins, counter)
    counter <- counter + chunkSize + 1
    chrBins <- c(chrBins, counter)
  }

  # Adding normal sized chunks to remaining intervals
  for (k in remainder:numOfBins) {
    chrBins <- c(chrBins, counter)
    counter <- counter + chunkSize
    chrBins <- c(chrBins, counter)
  }

  # Creating a data.frame with intervals
  interval.df <- as.data.frame(matrix(chrBins,ncol = 2, byrow = TRUE))
  colnames(interval.df) <- c("start", "end")

  # Adding to list
  chrBinList[[chrSizes[i,1]]] <- interval.df
}

至于测试值是否属于不同的 bin，我使用 apply 提出了一个缓慢的解决方案。但是，我目前正在研究并行应用功能。

【讨论】：