如何使用R的循环来获取聚合矩阵中除许多子集矩阵之外的部分？
How to use R's loop to get the part of the aggregate matrix other than the many subset matrices?
我有两个纯数字矩阵，它们的格式是一样的：第一列是组名，第二列是起始编号，第三列是结束编号。每个组都是一条染色体的名称。所有数据都是从 fasta 数据转换而来的。我正在尝试通过以整个基因组（table1）为基础来生成带有白名单（table2）的黑名单（table3）。
例如（下表中*空格代表列变化）：
table1的一行：
scahrs1_999 1 12
#表示'scahrs1_999'组的总数字组是1,2,3,4,5,6,7,8,9,10,11,12,13,14.
table2 的一些行：
scahrs1_999 2 3
scahrs1_999 6 8
scahrs1_999 11 12
#表示'scahrs1_999'组中数字的子集为'2, 3', '5, 6, 7', '10, 11'。
我要得到的结果是一组数字，在总集合中不包含任何子集，并且表现与table2的连续数字的子集相同。即：
表 3（结果）：
scahrs1_999 1
scahrs1_999 4 5
scahrs1_999 9 10
scahrs1_999 13 14
也排除只有 1 个数字的子集。即：
table3（最终结果）：
scahrs1_999 4 5
scahrs1_999 9 10
scahrs1_999 13 14
如下图，我有几个和上面类似的组'scahrs1_999'。显然，我不能一一计算。我知道 R 可以用 'loop program' 遍历每个组并得到相应的结果。但是我的编程能力不能胜任这个复杂的工作。
> Table1
1 scahrs1_1001 1 81142
2 scahrs1_1002 1 62661
3 scahrs1_1003 1 104891
4 scahrs1_1004 1 99296
5 scahrs1_1005 1 30919
6 scahrs1_1006 1 97599
7 scahrs1_1008 1 97078
8 scahrs1_1009 1 96958
9 scahrs1_1010 1 45677
> Table2
1 scahrs1_1001 1 753
2 scahrs1_1001 14931 15932
3 scahrs1_1001 17007 18008
4 scahrs1_1001 21211 22212
5 scahrs1_1001 40908 41909
6 scahrs1_1001 63233 64234
7 scahrs1_1001 76009 77010
8 scahrs1_1002 1068 2069
9 scahrs1_1002 12992 13993
10 scahrs1_1002 40448 41449
11 scahrs1_1003 2227 3228
12 scahrs1_1003 18453 19454
13 scahrs1_1003 28679 29680
14 scahrs1_1003 41161 42162
15 scahrs1_1003 41735 42736
16 scahrs1_1003 43040 44041
17 scahrs1_1003 64416 65417
18 scahrs1_1003 71219 72220
19 scahrs1_1003 96090 97091
20 scahrs1_1003 97306 98307
21 scahrs1_1004 1554 2555
22 scahrs1_1004 29086 30087
23 scahrs1_1004 44100 45101
24 scahrs1_1004 47799 48800
25 scahrs1_1004 59550 60551
26 scahrs1_1004 69356 70357
27 scahrs1_1004 71809 72810
28 scahrs1_1004 84272 85273
29 scahrs1_1004 89034 90035
30 scahrs1_1004 98627 99628
31 scahrs1_1005 6695 7696
32 scahrs1_1005 30160 31161
33 scahrs1_1006 298 1299
34 scahrs1_1006 70134 71135
35 scahrs1_1006 93750 94751
36 scahrs1_1008 3859 4860
37 scahrs1_1008 5575 6576
38 scahrs1_1008 7072 8073
39 scahrs1_1008 9342 10343
40 scahrs1_1008 11814 12815
41 scahrs1_1008 15290 16291
42 scahrs1_1008 40167 41168
43 scahrs1_1008 42890 43891
44 scahrs1_1008 44806 45807
45 scahrs1_1008 74442 75443
46 scahrs1_1008 82112 83113
47 scahrs1_1008 93766 94767
48 scahrs1_1008 95233 96234
49 scahrs1_1009 8000 9001
50 scahrs1_1009 37369 38370
51 scahrs1_1009 53086 54087
52 scahrs1_1009 83722 84723
53 scahrs1_1009 90044 91045
54 scahrs1_1010 11341 12342
55 scahrs1_1010 33500 34501
56 scahrs1_1010 34931 35932
57 scahrs1_1010 37937 38938
看起来有些数据是混乱的，但我相信这个 data.table 解决方案将 return 所需的 Table3 （稍微修改发布的数据集以更正错误）：
library(data.table)
Table1 <- data.table(group = c("scahrs1_1", "scahrs1_10", "scahrs1_100", "scahrs1_1000", "scahrs1_1001", "scahrs1_1002", "scahrs1_1003", "scahrs1_1004", "scahrs1_1005", "scahrs1_1006", "scahrs1_1008"), idxStart = 1L, idXEnd = c(1870329L, 925472L, 187969291L, 99113L, 81142L, 62661L, 104891L, 99296L, 30919L, 97599L, 97078L))
Table2 <- data.table(group = c(rep("scahrs1_1", 11), rep("scahrs1_10", 2)), idxStart = c(8158L, 17916L, 18644L, 31439L, 37022L, 62954L, 123548L, 129802L, 135683L, 135942L, 172435L, 22999L, 39260L), idxEnd = c(9159L, 18917L, 19645L, 32440L, 38023L, 63955L, 124549L, 130803L, 136684L, 136943L, 173436L, 24000L, 40261L))
fBetween <- function(grp, idxStart, idxEnd) {
idxRange <- unlist(Table1[group == grp, 2:3])
if (idxStart[1] == idxRange[1]) {
if (last(idxEnd) == idxRange[2]) {
list(idxStart = first(idxEnd, -1) + 1, idxEnd = last(idxStart, -1) - 1)
} else {
list(idxStart = idxEnd + 1, idxEnd = c(last(idxStart, -1) - 1, idxRange[2]))
}
} else {
if (last(idxEnd) == idxRange[2]) {
list(idxStart = c(idxRange[1], first(idxEnd, -1) + 1), idxEnd = idxStart - 1)
} else {
list(idxStart = c(idxRange[1], idxEnd + 1), idxEnd = c(idxStart - 1, idxRange[2]))
}
}
}
Table3 <- setorder(Table2, group, idxStart)[, fBetween(first(group), idxStart, idxEnd), group][idxStart != idxEnd]
Table3
#> group idxStart idxEnd
#> 1: scahrs1_1 1 8157
#> 2: scahrs1_1 9160 17915
#> 3: scahrs1_1 18918 18643
#> 4: scahrs1_1 19646 31438
#> 5: scahrs1_1 32441 37021
#> 6: scahrs1_1 38024 62953
#> 7: scahrs1_1 63956 123547
#> 8: scahrs1_1 124550 129801
#> 9: scahrs1_1 130804 135682
#> 10: scahrs1_1 136685 135941
#> 11: scahrs1_1 136944 172434
#> 12: scahrs1_1 173437 1870329
#> 13: scahrs1_10 1 22998
#> 14: scahrs1_10 24001 39259
#> 15: scahrs1_10 40262 925472
将两者分组后即可应用此功能data.tables
f <- function(s,e,is,ie) {
res = t(rbindlist(list(c(
list(c(s[1],is[1]-1)),
lapply(2:length(is), \(i) c(ie[i-1]+1, is[i]-1)),
list(c(ie[length(ie)]+1, e[length(e)]))
))))
list(start = res[,1], end = res[,2])
}
用法（使用 Table1 和 Table2，@jblood94 提供）
Table1[Table2, on="group"][, f(idxStart, idxEnd, i.idxStart,i.idxEnd), by=group][start<end]
输出：
group start end
1: scahrs1_1 1 8157
2: scahrs1_1 9160 17915
3: scahrs1_1 19646 31438
4: scahrs1_1 32441 37021
5: scahrs1_1 38024 62953
6: scahrs1_1 63956 123547
7: scahrs1_1 124550 129801
8: scahrs1_1 130804 135682
9: scahrs1_1 136944 172434
10: scahrs1_1 173437 1870329
11: scahrs1_10 1 22998
12: scahrs1_10 24001 39259
13: scahrs1_10 40262 925472
我有两个纯数字矩阵，它们的格式是一样的：第一列是组名，第二列是起始编号，第三列是结束编号。每个组都是一条染色体的名称。所有数据都是从 fasta 数据转换而来的。我正在尝试通过以整个基因组（table1）为基础来生成带有白名单（table2）的黑名单（table3）。
例如（下表中*空格代表列变化）： table1的一行：
scahrs1_999 1 12
#表示'scahrs1_999'组的总数字组是1,2,3,4,5,6,7,8,9,10,11,12,13,14.
table2 的一些行：
scahrs1_999 2 3
scahrs1_999 6 8
scahrs1_999 11 12
#表示'scahrs1_999'组中数字的子集为'2, 3', '5, 6, 7', '10, 11'。
我要得到的结果是一组数字，在总集合中不包含任何子集，并且表现与table2的连续数字的子集相同。即：
表 3（结果）：
scahrs1_999 1
scahrs1_999 4 5
scahrs1_999 9 10
scahrs1_999 13 14
也排除只有 1 个数字的子集。即：
table3（最终结果）：
scahrs1_999 4 5
scahrs1_999 9 10
scahrs1_999 13 14
如下图，我有几个和上面类似的组'scahrs1_999'。显然，我不能一一计算。我知道 R 可以用 'loop program' 遍历每个组并得到相应的结果。但是我的编程能力不能胜任这个复杂的工作。
> Table1
1 scahrs1_1001 1 81142
2 scahrs1_1002 1 62661
3 scahrs1_1003 1 104891
4 scahrs1_1004 1 99296
5 scahrs1_1005 1 30919
6 scahrs1_1006 1 97599
7 scahrs1_1008 1 97078
8 scahrs1_1009 1 96958
9 scahrs1_1010 1 45677
> Table2
1 scahrs1_1001 1 753
2 scahrs1_1001 14931 15932
3 scahrs1_1001 17007 18008
4 scahrs1_1001 21211 22212
5 scahrs1_1001 40908 41909
6 scahrs1_1001 63233 64234
7 scahrs1_1001 76009 77010
8 scahrs1_1002 1068 2069
9 scahrs1_1002 12992 13993
10 scahrs1_1002 40448 41449
11 scahrs1_1003 2227 3228
12 scahrs1_1003 18453 19454
13 scahrs1_1003 28679 29680
14 scahrs1_1003 41161 42162
15 scahrs1_1003 41735 42736
16 scahrs1_1003 43040 44041
17 scahrs1_1003 64416 65417
18 scahrs1_1003 71219 72220
19 scahrs1_1003 96090 97091
20 scahrs1_1003 97306 98307
21 scahrs1_1004 1554 2555
22 scahrs1_1004 29086 30087
23 scahrs1_1004 44100 45101
24 scahrs1_1004 47799 48800
25 scahrs1_1004 59550 60551
26 scahrs1_1004 69356 70357
27 scahrs1_1004 71809 72810
28 scahrs1_1004 84272 85273
29 scahrs1_1004 89034 90035
30 scahrs1_1004 98627 99628
31 scahrs1_1005 6695 7696
32 scahrs1_1005 30160 31161
33 scahrs1_1006 298 1299
34 scahrs1_1006 70134 71135
35 scahrs1_1006 93750 94751
36 scahrs1_1008 3859 4860
37 scahrs1_1008 5575 6576
38 scahrs1_1008 7072 8073
39 scahrs1_1008 9342 10343
40 scahrs1_1008 11814 12815
41 scahrs1_1008 15290 16291
42 scahrs1_1008 40167 41168
43 scahrs1_1008 42890 43891
44 scahrs1_1008 44806 45807
45 scahrs1_1008 74442 75443
46 scahrs1_1008 82112 83113
47 scahrs1_1008 93766 94767
48 scahrs1_1008 95233 96234
49 scahrs1_1009 8000 9001
50 scahrs1_1009 37369 38370
51 scahrs1_1009 53086 54087
52 scahrs1_1009 83722 84723
53 scahrs1_1009 90044 91045
54 scahrs1_1010 11341 12342
55 scahrs1_1010 33500 34501
56 scahrs1_1010 34931 35932
57 scahrs1_1010 37937 38938
看起来有些数据是混乱的，但我相信这个 data.table 解决方案将 return 所需的 Table3 （稍微修改发布的数据集以更正错误）：
library(data.table)
Table1 <- data.table(group = c("scahrs1_1", "scahrs1_10", "scahrs1_100", "scahrs1_1000", "scahrs1_1001", "scahrs1_1002", "scahrs1_1003", "scahrs1_1004", "scahrs1_1005", "scahrs1_1006", "scahrs1_1008"), idxStart = 1L, idXEnd = c(1870329L, 925472L, 187969291L, 99113L, 81142L, 62661L, 104891L, 99296L, 30919L, 97599L, 97078L))
Table2 <- data.table(group = c(rep("scahrs1_1", 11), rep("scahrs1_10", 2)), idxStart = c(8158L, 17916L, 18644L, 31439L, 37022L, 62954L, 123548L, 129802L, 135683L, 135942L, 172435L, 22999L, 39260L), idxEnd = c(9159L, 18917L, 19645L, 32440L, 38023L, 63955L, 124549L, 130803L, 136684L, 136943L, 173436L, 24000L, 40261L))
fBetween <- function(grp, idxStart, idxEnd) {
idxRange <- unlist(Table1[group == grp, 2:3])
if (idxStart[1] == idxRange[1]) {
if (last(idxEnd) == idxRange[2]) {
list(idxStart = first(idxEnd, -1) + 1, idxEnd = last(idxStart, -1) - 1)
} else {
list(idxStart = idxEnd + 1, idxEnd = c(last(idxStart, -1) - 1, idxRange[2]))
}
} else {
if (last(idxEnd) == idxRange[2]) {
list(idxStart = c(idxRange[1], first(idxEnd, -1) + 1), idxEnd = idxStart - 1)
} else {
list(idxStart = c(idxRange[1], idxEnd + 1), idxEnd = c(idxStart - 1, idxRange[2]))
}
}
}
Table3 <- setorder(Table2, group, idxStart)[, fBetween(first(group), idxStart, idxEnd), group][idxStart != idxEnd]
Table3
#> group idxStart idxEnd
#> 1: scahrs1_1 1 8157
#> 2: scahrs1_1 9160 17915
#> 3: scahrs1_1 18918 18643
#> 4: scahrs1_1 19646 31438
#> 5: scahrs1_1 32441 37021
#> 6: scahrs1_1 38024 62953
#> 7: scahrs1_1 63956 123547
#> 8: scahrs1_1 124550 129801
#> 9: scahrs1_1 130804 135682
#> 10: scahrs1_1 136685 135941
#> 11: scahrs1_1 136944 172434
#> 12: scahrs1_1 173437 1870329
#> 13: scahrs1_10 1 22998
#> 14: scahrs1_10 24001 39259
#> 15: scahrs1_10 40262 925472
将两者分组后即可应用此功能data.tables
f <- function(s,e,is,ie) {
res = t(rbindlist(list(c(
list(c(s[1],is[1]-1)),
lapply(2:length(is), \(i) c(ie[i-1]+1, is[i]-1)),
list(c(ie[length(ie)]+1, e[length(e)]))
))))
list(start = res[,1], end = res[,2])
}
用法（使用 Table1 和 Table2，@jblood94 提供）
Table1[Table2, on="group"][, f(idxStart, idxEnd, i.idxStart,i.idxEnd), by=group][start<end]
输出：
group start end
1: scahrs1_1 1 8157
2: scahrs1_1 9160 17915
3: scahrs1_1 19646 31438
4: scahrs1_1 32441 37021
5: scahrs1_1 38024 62953
6: scahrs1_1 63956 123547
7: scahrs1_1 124550 129801
8: scahrs1_1 130804 135682
9: scahrs1_1 136944 172434
10: scahrs1_1 173437 1870329
11: scahrs1_10 1 22998
12: scahrs1_10 24001 39259
13: scahrs1_10 40262 925472