我有以下尝试计数的简单脚本
SAM/BAM file 中用 "CB:Z"编码的标签:
samtools view -h small.bam | grep "CB:Z:" |
sed 's/.*CB:Z:\([ACGT]*\).*/\1/' |
sort |
uniq -c |
awk '{print $2 " " $1}'
通常它需要处理 4000 万行。该代码需要大约 1 小时才能完成。这一行
sed 's/.*CB:Z:\([ACGT]*\).*/\1/'
非常耗时。我怎样才能加快速度?
我使用正则表达式的原因是“CB”标签列位置
不固定。有时在第 20 列,有时在第 21 列。
示例 BAM 文件可以在 HERE 中找到。
更新
完整 4000 万行文件的速度比较:
我的初始代码:
real 21m47.088s
user 26m51.148s
sys 1m27.912s
James Brown 与 AWK:real 1m28.898s
user 2m41.336s
sys 0m6.864s
詹姆斯布朗与 MAWK:real 1m10.642s
user 1m41.196s
sys 0m6.484s
最佳答案
另一个 awk,很像 @tripleee 的,我假设:
$ samtools view -h small.bam | awk '
match($0,/CB:Z:[ACGT]*/) { # use match for the regex match
a[substr($0,RSTART+5,RLENGTH-5)]++ # len(CB:z:)==5, hence +-5
}
END {
for(i in a)
print i,a[i] # sample output,tweak it to your liking
}'
示例输出:...
TCTTAATCGTCC 175
GGGAAGGCCTAA 190
TCGGCCGATCGG 32
GACTTCCAAGCC 76
CCGCGGCATCGG 36
TAGCGATCGTGG 125
...
注意 :您的 sed 's/.*CB:Z:...
匹配最后一个实例,而我的 awk 'match($0,/CB:Z:[ACGT]*/)...
匹配第一个实例。注意 2 :在评论中引用 @Sundeep:- - 使用
LC_ALL=C mawk '..'
将提供更好的速度。
关于linux - 如何加速在非常大的单单元 BAM 文件上使用 Regex 的 sed,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/65000162/