bash 脚本 - 使用 sed 中的模式列表删除子字符串
bash script - Use patterns list in sed to remove substrings
我有这个文件 (adapters.txt),其中包含模式列表:
cactctttccctacacgacgctcttccg
cactctttccctacacgacgctcttccgaatcta
cactctttccctacacgacgctcttccgaatctaatt
cactctttccctacacgacgctcttccgaatctaatta
cactctttccctacacgacgctcttccgaatctag
cactctttccctacacgacgctcttccgaatctagc
cactctttccctacacgacgctcttccgacctcattcc
cactctttccctacacgacgctcttccgacctcattcccaccctcttccg
cactctttccctacacgacgctcttccgatc
cactctttccctacacgacgctcttccgatccaatt
cactctttccctacacgacgctcttccgatttagc
cactctttccctacacgacgctcttccgatttagct
cactctttccctacacgacgctcttccgatttcattc
cactctttccctacacgacgctcttccgatttcattcttcccc
cactctttccctacacgacgctcttccgattttatttc
cactctttccctacacgacgctcttccggatcta
cactctttccctacacgacgctcttccggatctaatt
cactctttccctacacgacgctcttccggatctaattc
cactctttccctacacgacgctcttccggatctaattca
cactctttccctacacgacgctcttccggatctagctt
cactctttccctacacgacgctcttccggttcta
cactctttccctacacgacgctttccgatcta
cactctttccctacacgacgctttccgatctaattc
cactctttccctacacgacgtcttccgatctaattctggaccatagtgcaatgt
cactctttccctacacgcgctcttccgatcta
cactctttccctacacgcgctcttccgatctaattcg
cactctttccctacacgcgctcttccgatctaattcgg
cactctttccctacacgcgctcttccgatctaattcggcgg
cactctttccctacacgcgctcttccgatctagct
cactctttccctaccgacgctcttccgatcta
cactctttccctacacgacg
我需要从“sequences.fasta”文件中找到并删除这些模式:
>seq01
cactctttccctacacgacgctcttccgWANTEDSEQUENCE
>seq01
cactctttccctacacgacgctcttccgaatctaWANTEDSEQUENCE
>seq03
cactctttccctacacgacgctcttccgaatctaattWANTEDSEQUENCE
>seq04
cactctttccctacacgacgctcttccgaatctaattaWANTEDSEQUENCE
>seq05
cactctttccctacacgcgctcttccgatctaattcggWANTEDSEQUENCE
>seq06
cactctttccctacacgcgctcttccgatctaattcggcggWANTEDSEQUENCE
>seq07
cactctttccctacacgcgctcttccgatctagctWANTEDSEQUENCE
>seq08
cactctttccctaccgacgctcttccgatctaWANTEDSEQUENCE
所以想要的输出应该是:
>seq01
WANTEDSEQUENCE
>seq02
WANTEDSEQUENCE
>seq03
WANTEDSEQUENCE
>seq04
WANTEDSEQUENCE
>seq05
WANTEDSEQUENCE
>seq06
WANTEDSEQUENCE
>seq07
WANTEDSEQUENCE
>seq08
WANTEDSEQUENCE
(只是为了示例,我使用了“WANTEDSEQUENCE”而不是真正的序列)
我尝试了以下方法(以及一些变体。我也尝试了 while read
):
ADAPS=($(cat adapters.txt))
FASTA="sequences.fasta"
for ADAP in "${ADAPS[@]}";
do
sed "s/${ADAP}//g" "${FASTA}" > output.fasta
done
但是我得到了这个:
>seq01
ctcttccgWANTEDSEQUENCE
>seq01
ctcttccgaatctaWANTEDSEQUENCE
>seq03
ctcttccgaatctaattWANTEDSEQUENCE
>seq04
ctcttccgaatctaattaWANTEDSEQUENCE
>seq05
cactctttccctacacgcgctcttccgatctaattcggWANTEDSEQUENCE
>seq06
cactctttccctacacgcgctcttccgatctaattcggcggWANTEDSEQUENCE
>seq07
cactctttccctacacgcgctcttccgatctagctWANTEDSEQUENCE
>seq08
cactctttccctaccgacgctcttccgatctaWANTEDSEQUENCE
我该如何解决这个问题?
Sort adapters.txt
in reverse order by its line length,从其输出创建一个 sed
脚本,并将其与 bash 的命令替换 <(...)
一起使用,并应用第二个 sed
它到 sequences.fasta
:
sed -f <(awk '{ print length, [=10=] }' adapters.txt | sort -rn | cut -d" " -f2- | sed -E 's/(.*)/s|&||/') sequences.fasta
输出:
>seq01
WANTEDSEQUENCE
>seq01
WANTEDSEQUENCE
>seq03
WANTEDSEQUENCE
>seq04
WANTEDSEQUENCE
>seq05
WANTEDSEQUENCE
>seq06
WANTEDSEQUENCE
>seq07
WANTEDSEQUENCE
>seq08
WANTEDSEQUENCE
adapters.txt
的排序是必要的,因为它包含同一文件中其他字符串的子字符串。
多行和文件中的相同代码:
awk '{ print length, [=11=] }' adapters.txt | sort -rn | cut -d" " -f2- > adapters_sorted.txt
sed -E 's/(.*)/s|&||/' adapters_sorted.txt > sed.script
sed -f sed.script sequences.fasta
这是一个 POSIX awk 解决方案:
$ awk 'NR==FNR{seq[FNR]=[=10=]; x=FNR; next}
{for(i=1; i<=x; i++) if ([=10=] ~ "^" seq[i]) {sub(seq[i],""); print [=10=]; next}
print}
' <(awk '{ print length()"\t"[=10=]}' adapters.txt | sort -nr | cut -f2) sequences.fasta
>seq01
WANTEDSEQUENCE
>seq01
WANTEDSEQUENCE
>seq03
WANTEDSEQUENCE
>seq04
WANTEDSEQUENCE
>seq05
WANTEDSEQUENCE
>seq06
WANTEDSEQUENCE
>seq07
WANTEDSEQUENCE
>seq08
WANTEDSEQUENCE
或gawk
您可以在其中对序列从长到短进行内部排序:
$ gawk 'BEGIN{PROCINFO["sorted_in"] = "@val_num_desc"}
NR==FNR { seq[[=11=]] = length([=11=]); next }
{for (e in seq) if([=11=]~"^" e) {sub(e,""); print [=11=]; next}
print}
' adapters.txt sequences.fasta
# same output
使用 GNU awk sorted_in:
$ cat tst.awk
NR==FNR {
adapters2lengths[] = length()
next
}
!/^>/ {
PROCINFO["sorted_in"] = "@val_num_desc"
for (adapter in adapters2lengths) {
if ( index([=10=],adapter) == 1 ) {
[=10=] = substr([=10=],adapters2lengths[adapter]+1)
break
}
}
}
{ print }
.
$ awk -f tst.awk adapters.txt sequences.fasta
>seq01
WANTEDSEQUENCE
>seq01
WANTEDSEQUENCE
>seq03
WANTEDSEQUENCE
>seq04
WANTEDSEQUENCE
>seq05
WANTEDSEQUENCE
>seq06
WANTEDSEQUENCE
>seq07
WANTEDSEQUENCE
>seq08
WANTEDSEQUENCE
这个和@dawg 的 gawk 解决方案在功能上的区别在于,这个解决方案只进行字符串比较,而他们的解决方案进行正则表达式比较——只有当你的“adapters.txt”文件包含正则表达式元字符时才重要,所有否则等于我只是更喜欢使用字符串,除非我需要正则表达式。
使用 GNU ed
和 bash
。
#!/usr/bin/env bash
ed -s sequences.fasta < <(
printf '%s\n' '1,$-1s/$/\|/' '1,$j' 's/^/,s\//' 's/$/\/\//' '$a' ,p Q . ,p Q |
ed -s adapters.txt
)
我有这个文件 (adapters.txt),其中包含模式列表:
cactctttccctacacgacgctcttccg
cactctttccctacacgacgctcttccgaatcta
cactctttccctacacgacgctcttccgaatctaatt
cactctttccctacacgacgctcttccgaatctaatta
cactctttccctacacgacgctcttccgaatctag
cactctttccctacacgacgctcttccgaatctagc
cactctttccctacacgacgctcttccgacctcattcc
cactctttccctacacgacgctcttccgacctcattcccaccctcttccg
cactctttccctacacgacgctcttccgatc
cactctttccctacacgacgctcttccgatccaatt
cactctttccctacacgacgctcttccgatttagc
cactctttccctacacgacgctcttccgatttagct
cactctttccctacacgacgctcttccgatttcattc
cactctttccctacacgacgctcttccgatttcattcttcccc
cactctttccctacacgacgctcttccgattttatttc
cactctttccctacacgacgctcttccggatcta
cactctttccctacacgacgctcttccggatctaatt
cactctttccctacacgacgctcttccggatctaattc
cactctttccctacacgacgctcttccggatctaattca
cactctttccctacacgacgctcttccggatctagctt
cactctttccctacacgacgctcttccggttcta
cactctttccctacacgacgctttccgatcta
cactctttccctacacgacgctttccgatctaattc
cactctttccctacacgacgtcttccgatctaattctggaccatagtgcaatgt
cactctttccctacacgcgctcttccgatcta
cactctttccctacacgcgctcttccgatctaattcg
cactctttccctacacgcgctcttccgatctaattcgg
cactctttccctacacgcgctcttccgatctaattcggcgg
cactctttccctacacgcgctcttccgatctagct
cactctttccctaccgacgctcttccgatcta
cactctttccctacacgacg
我需要从“sequences.fasta”文件中找到并删除这些模式:
>seq01
cactctttccctacacgacgctcttccgWANTEDSEQUENCE
>seq01
cactctttccctacacgacgctcttccgaatctaWANTEDSEQUENCE
>seq03
cactctttccctacacgacgctcttccgaatctaattWANTEDSEQUENCE
>seq04
cactctttccctacacgacgctcttccgaatctaattaWANTEDSEQUENCE
>seq05
cactctttccctacacgcgctcttccgatctaattcggWANTEDSEQUENCE
>seq06
cactctttccctacacgcgctcttccgatctaattcggcggWANTEDSEQUENCE
>seq07
cactctttccctacacgcgctcttccgatctagctWANTEDSEQUENCE
>seq08
cactctttccctaccgacgctcttccgatctaWANTEDSEQUENCE
所以想要的输出应该是:
>seq01
WANTEDSEQUENCE
>seq02
WANTEDSEQUENCE
>seq03
WANTEDSEQUENCE
>seq04
WANTEDSEQUENCE
>seq05
WANTEDSEQUENCE
>seq06
WANTEDSEQUENCE
>seq07
WANTEDSEQUENCE
>seq08
WANTEDSEQUENCE
(只是为了示例,我使用了“WANTEDSEQUENCE”而不是真正的序列)
我尝试了以下方法(以及一些变体。我也尝试了 while read
):
ADAPS=($(cat adapters.txt))
FASTA="sequences.fasta"
for ADAP in "${ADAPS[@]}";
do
sed "s/${ADAP}//g" "${FASTA}" > output.fasta
done
但是我得到了这个:
>seq01
ctcttccgWANTEDSEQUENCE
>seq01
ctcttccgaatctaWANTEDSEQUENCE
>seq03
ctcttccgaatctaattWANTEDSEQUENCE
>seq04
ctcttccgaatctaattaWANTEDSEQUENCE
>seq05
cactctttccctacacgcgctcttccgatctaattcggWANTEDSEQUENCE
>seq06
cactctttccctacacgcgctcttccgatctaattcggcggWANTEDSEQUENCE
>seq07
cactctttccctacacgcgctcttccgatctagctWANTEDSEQUENCE
>seq08
cactctttccctaccgacgctcttccgatctaWANTEDSEQUENCE
我该如何解决这个问题?
Sort adapters.txt
in reverse order by its line length,从其输出创建一个 sed
脚本,并将其与 bash 的命令替换 <(...)
一起使用,并应用第二个 sed
它到 sequences.fasta
:
sed -f <(awk '{ print length, [=10=] }' adapters.txt | sort -rn | cut -d" " -f2- | sed -E 's/(.*)/s|&||/') sequences.fasta
输出:
>seq01 WANTEDSEQUENCE >seq01 WANTEDSEQUENCE >seq03 WANTEDSEQUENCE >seq04 WANTEDSEQUENCE >seq05 WANTEDSEQUENCE >seq06 WANTEDSEQUENCE >seq07 WANTEDSEQUENCE >seq08 WANTEDSEQUENCE
adapters.txt
的排序是必要的,因为它包含同一文件中其他字符串的子字符串。
多行和文件中的相同代码:
awk '{ print length, [=11=] }' adapters.txt | sort -rn | cut -d" " -f2- > adapters_sorted.txt
sed -E 's/(.*)/s|&||/' adapters_sorted.txt > sed.script
sed -f sed.script sequences.fasta
这是一个 POSIX awk 解决方案:
$ awk 'NR==FNR{seq[FNR]=[=10=]; x=FNR; next}
{for(i=1; i<=x; i++) if ([=10=] ~ "^" seq[i]) {sub(seq[i],""); print [=10=]; next}
print}
' <(awk '{ print length()"\t"[=10=]}' adapters.txt | sort -nr | cut -f2) sequences.fasta
>seq01
WANTEDSEQUENCE
>seq01
WANTEDSEQUENCE
>seq03
WANTEDSEQUENCE
>seq04
WANTEDSEQUENCE
>seq05
WANTEDSEQUENCE
>seq06
WANTEDSEQUENCE
>seq07
WANTEDSEQUENCE
>seq08
WANTEDSEQUENCE
或gawk
您可以在其中对序列从长到短进行内部排序:
$ gawk 'BEGIN{PROCINFO["sorted_in"] = "@val_num_desc"}
NR==FNR { seq[[=11=]] = length([=11=]); next }
{for (e in seq) if([=11=]~"^" e) {sub(e,""); print [=11=]; next}
print}
' adapters.txt sequences.fasta
# same output
使用 GNU awk sorted_in:
$ cat tst.awk
NR==FNR {
adapters2lengths[] = length()
next
}
!/^>/ {
PROCINFO["sorted_in"] = "@val_num_desc"
for (adapter in adapters2lengths) {
if ( index([=10=],adapter) == 1 ) {
[=10=] = substr([=10=],adapters2lengths[adapter]+1)
break
}
}
}
{ print }
.
$ awk -f tst.awk adapters.txt sequences.fasta
>seq01
WANTEDSEQUENCE
>seq01
WANTEDSEQUENCE
>seq03
WANTEDSEQUENCE
>seq04
WANTEDSEQUENCE
>seq05
WANTEDSEQUENCE
>seq06
WANTEDSEQUENCE
>seq07
WANTEDSEQUENCE
>seq08
WANTEDSEQUENCE
这个和@dawg 的 gawk 解决方案在功能上的区别在于,这个解决方案只进行字符串比较,而他们的解决方案进行正则表达式比较——只有当你的“adapters.txt”文件包含正则表达式元字符时才重要,所有否则等于我只是更喜欢使用字符串,除非我需要正则表达式。
使用 GNU ed
和 bash
。
#!/usr/bin/env bash
ed -s sequences.fasta < <(
printf '%s\n' '1,$-1s/$/\|/' '1,$j' 's/^/,s\//' 's/$/\/\//' '$a' ,p Q . ,p Q |
ed -s adapters.txt
)