提问人：SWK 提问时间：9/7/2023 最后编辑：SWK 更新时间：9/8/2023 访问量：79
计算多个范围的唯一整数数，以计算爆炸结果的水平覆盖率 [关闭]
count number of unique integers from multiple ranges to calculate horizontal coverage of blast result [closed]
问：
闭。这个问题正在寻求有关书籍、工具、软件库等的建议。它不符合 Stack Overflow 准则。它目前不接受答案。
我们不允许提出有关书籍、工具、软件库等建议的问题。您可以编辑问题，以便用事实和引文来回答。
2个月前关闭。
当我尝试将 DNA 序列读数映射/对齐到蛋白质数据库时，我找不到许多提供水平参考/主题序列覆盖率的映射工具，但是 blastx/diamond 确实对齐了每个查询读取的参考基因的开始和结束位置。因此，我需要计算来自多个范围的唯一整数的数量 [subjectStart..subjectEnd] 来获取比对序列的长度。
blastx/diamond 输出被修剪为 3 列，用于参考基因 ID、开始和结束，用于修复反向序列。每个行 ID 的范围可能重叠或重复，也可能没有，这些重复的数字只需要计算一次。awk -F"\t" '{if ($2 <= $3) {print $0} else {print $1"\t"$3"\t"$2}}'
输入
ID Start End
A 1 50
A 2 45
A 25 150
A 50 150
A 155 200
A 205 300
B 5 50
B 61 70
B 81 100
C 1 500`
所需的输出。
ID count
A 292
B 76
C 500
我请求帮助以获得本 Ask 第一版中所需的输出，我可以立即从许多专家那里获得完美的解决方案!! 以下对我来说非常有效，我学到了很多东西。我衷心感谢大家!!
Ruby 的 @dawg
ruby -lane 'BEGIN{h=Hash.new { |hash, key| hash[key] = Set.new() }}
h[$F[0]].merge(($F[1].to_i..$F[2].to_i)) if $.>1
END{
puts "ID\tCount"
h.each{|k,v| puts "#{k}\t#{v.length}"}
}
' Input.file > Output.file
Perl 的 @zdim
perl -MData::Dumper -MList::Util=uniq -wnE'
($id, $beg, $end) = split;
next if not $beg or $beg =~ /[^0-9]/ or not $end or $end =~ /[^0-9]/;
push @{$res{$id}}, $beg..$end;
}{
for (keys %res) { $res{$_} = uniq sort { $a <=> $b } @{$res{$_}} };
say Dumper \%res
' Input.file > Output.file
Awk 由 @pmf
awk 'NR>1 {s[$1] += $3 - ($2 <= b[$1] ? ($3 <= b[$1] ? $3 : b[$1]) + 1 : $2) + 1; b[$1] = b[$1] <= $3 ? $3 : b[$1]} END {OFS="\t"; print "ID", "count"; for (i in s) {print i, s[i]}}' Input.file > Output.file
Awk 由 @dawg
awk '
FNR>1{for(i=$2;i<=$3;i++) ss[$1 "|" i]}
END{
print "ID", "Count"
for (e in ss) {
split(e,idx,"\|")
cnt[idx[1]]++
}
for (e in cnt) print e, cnt[e]
}
' OFS="\t" Input.file > Output.file
答：
0赞
zdim
9/7/2023
#1
单程
perl -MData::Dumper -MList::Util=uniq -wnE'
($id, $beg, $end) = split;
next if not $beg or $beg =~ /[^0-9]/ or not $end or $end =~ /[^0-9]/;
push @{$res{$id}}, $beg..$end;
}{
for (keys %res) { @{$res{$_}} = uniq sort { $a <=> $b } @{$res{$_}} };
say Dumper \%res
' data.txt
此命令行程序（“单行”）是一个演示，请在普通程序中重写。
该语法标记块的开始 -- 一旦处理了文件中的所有行，该代码就会运行。可以在这里写。}{ENDEND { ... }
这将为每个 ID 生成唯一且经过排序的整数列表。如果特别需要计数，只需将列表分配给标量即可
for (keys %res) { $res{$_} = uniq sort { $a <=> $b } @{$res{$_}} };
笔记
检查字段是否为数字是在这里使用正则表达式完成的，因为它们是假定的整数。否则，最好从 Scalar：：Util 使用
looks_like_number在较新的 Perls 中，可以使用后缀取消引用语法来完成取消引用
push $res{$_}->@*, $beg..$end;使用 Key：：Sort 进行排序要容易得多（尤其是在更复杂的情况下）
use Key::Sort qw(usort); # use "isort" if negative ints are possible foreach my $id (sort keys %res) { $res{$id}->@* = uniq usort $res{$id}->@*; }或者，获取每个 ID 的计数而不是数字
$res{$id} = uniq usort $res{$id}->@*;还要注意库的方法，如果不是这样，这里会很好
_inplaceuniq我首先假设它本身在内部使用排序（因为“保留唯一元素的顺序”，这会更快）;但是，如果不排序，那么首先过滤掉重复项，然后进行排序会更快。
sortuniquniquniquniq
评论
1赞
zdim
9/7/2023
@SWK 自发布以来更新了很多。我这样做的目的是你使用“演示”命令行程序（在命令行上复制粘贴，并且应该按原样工作）来编写一个更好的、完整的程序。如果您有任何疑问，请告诉我
0赞
SWK
9/7/2023
我们可以添加命令来输出计数吗？非常感谢!!
1赞
zdim
9/7/2023
@SWK我做到了——在最后一次编辑中（在我发表评论之后）。请再看一遍。
1赞
dawg
9/7/2023
#2
鉴于：
cat file
ID Start End
A 1 50
A 2 45
A 25 150
A 50 150
A 155 200
A 205 300
B 5 50
B 61 70
B 81 100
C 1 500
这里有一个 Ruby 来做到这一点：
ruby -lane 'BEGIN{h=Hash.new { |hash, key| hash[key] = Set.new() }}
h[$F[0]].merge(($F[1].to_i..$F[2].to_i)) if $.>1
END{
puts "ID\tCount"
h.each{|k,v| puts "#{k}\t#{v.length}"}
}
' file
或者这个awk：
awk '
FNR>1{for(i=$2;i<=$3;i++) ss[$1 "|" i]}
END{
print "ID", "Count"
for (e in ss) {
split(e,idx,"\|")
cnt[idx[1]]++
}
for (e in cnt) print e, cnt[e]
}
' OFS="\t" file
打印（但 awk 可能是无序的......
ID Count
A 292
B 76
C 500
从评论来看，B 50 5 而不是 B 5 50，我们可以让它工作吗？
ruby -lane 'BEGIN{h=Hash.new { |hash, key| hash[key] = Set.new() }}
v1,v2=[$F[1],$F[2]].map(&:to_i)
h[$F[0]].merge(v1..v2) if $.>1 && v2>v1
END{
puts "ID\tCount"
h.each{|k,v| puts "#{k}\t#{v.length}"}
}
' file
不过，这应该不是必需的。运算符创建一个范围，该范围将添加到唯一值集中。范围为 null 值，因此不会添加任何内容。..50..5
下面是 IRB（Ruby 交互式 shell）中的测试：
irb(main):043:0> s=Set.new()
=> #<Set: {}>
irb(main):044:0> s.merge(50..5)
=> #<Set: {}>
irb(main):045:0> s.merge(1..5)
=> #<Set: {1, 2, 3, 4, 5}>
所以你可以看到没有添加任何内容。
评论
0赞
SWK
9/7/2023
对于 a 或一些顺序相反的行，例如 B 50 5 而不是 B 5 50，我们可以让它工作吗？非常感谢!!
1赞
dawg
9/7/2023
应该如何解释这样的反线？这些行只是跳过了吗？
0赞
SWK
9/7/2023
是的，这些都被跳过了。可以使用awk print修改数据文件，我想知道如何在Ruby中操作。非常感谢!!
1赞
dawg
9/7/2023
它应该是一个 NoOP，因为是一个负范围，因此是一个空列表。无论如何，我都会在编辑中添加代码。50..5
评论
awkawk 'NR>1 {s[$1] += $3 - ($2 <= b[$1] ? ($3 <= b[$1] ? $3 : b[$1]) + 1 : $2) + 1; b[$1] = b[$1] <= $3 ? $3 : b[$1]} END {OFS="\t"; print "ID", "count"; for (i in s) {print i, s[i]}}'