如何从fasta文件中删除重复项，但根据标头保留每个组至少一个
如何从fasta文件中删除重复项，但根据标头保留每个组至少一个
提问于 2020-07-25 19:52:49
我有一个多快件文件，如下所示：
(所有序列均大于100 and，不止一条直线，长度相同)
>Lineage1_samplenameA
CGCTTCAACGGAATGGATCTACGTTACAGCCTGCATAAAGAAAACGGAGTTGCCGAGGACGAAAGCGACTTTAGGTTCTGTCCGTTGTCTTTGGCGGAAA
>Lineage2_samplenameB
AAATTCAACGGAATGGATCTACGTTACAGCCTGCATAAAGAAAACGGAGTTGCCGAGGACGAAAGCGACTTTAGGTTCTGTCCGTTGTCTTTGGCGGAAG
>Lineage3_samplenameC
CGCTTCAACGGAATGGATCTACGTTACAGCCTGCATAAAGAAAACGGAGTTGCCGAGGACGAAAGCGACTTTAGGTTCTGTCCGTTGTCTTTGGCGGAAA
>Lineage3_samplenameD
CGCTTCAACGGAATGGATCTACGTTACAGCCTGCATAAAGAAAACGGAGTTGCCGAGGACGAAAGCGACTTTAGGTTCTGTCCGTTGTCTTTGGCGGAAA我需要移除副本，但至少要保持每个谱系的顺序。因此，在这个简单的示例(注意，samplenameA、C和D是相同的)中，我只想删除samplenameD或samplenameC，而不是两者都删除。最后，我希望获得与原始文件相同的头信息。
示例输出：
>Lineage1_samplenameA
CGCTTCAACGGAATGGATCTACGTTACAGCCTGCATAAAGAAAACGGAGTTGCCGAGGACGAAAGCGACTTTAGGTTCTGTCCGTTGTCTTTGGCGGAAA
>Lineage2_samplenameB
AAATTCAACGGAATGGATCTACGTTACAGCCTGCATAAAGAAAACGGAGTTGCCGAGGACGAAAGCGACTTTAGGTTCTGTCCGTTGTCTTTGGCGGAAG
>Lineage3_samplenameC
CGCTTCAACGGAATGGATCTACGTTACAGCCTGCATAAAGAAAACGGAGTTGCCGAGGACGAAAGCGACTTTAGGTTCTGTCCGTTGTCTTTGGCGGAAA我找到了一种移除重复物品的方法。多亏了皮埃尔·林登鲍姆。
sed -e '/^>/s/$/@/' -e 's/^>/#/'
file.fasta |\
tr -d '\n' | tr "#" "\n" | tr "@"
"\t" |\
sort -u -t ' ' -f -k 2,2 |\
sed -e 's/^/>/' -e 's/\t/\n/'在上面的示例中运行此操作将得到以下结果：
>Lineage1_samplenameA
CGCTTCAACGGAATGGATCTACGTTACAGCCTGCATAAAGAAAACGGAGTTGCCGAGGACGAAAGCGACTTTAGGTTCTGTCCGTTGTCTTTGGCGGAAA
>Lineage2_samplenameB
AAATTCAACGGAATGGATCTACGTTACAGCCTGCATAAAGAAAACGGAGTTGCCGAGGACGAAAGCGACTTTAGGTTCTGTCCGTTGTCTTTGGCGGAAG因此，->丢失了3序列
现在，我只是在寻找一个快速的解决方案，以删除重复，但保持至少一个序列的每个谱系基础上的fasta头。
我对脚本很陌生..。任何在bash/python/R中的想法都是受欢迎的。
谢谢！
回答 1
Stack Overflow用户
回答已采纳
发布于 2020-07-25 21:31:45
在这种情况下，我可以看到两个相对好的替代方案。( A)查看现有工具(如Biopython库或FASTX工具包)。我认为他们都有很好的命令来完成这里的大部分工作，所以学习它们也许是值得的。或者，B写你自己的。在这种情况下，您可能需要尝试(我将继续使用python)：
循环遍历文件，逐行，并将沿袭/序列数据添加到字典中。我建议用这个序列作为钥匙。这样，您就可以很容易地知道您是否已经遇到了这个密钥。
myfasta = {}
if myfasta[sequence]:
myfasta[sequence].append(lineage_id)
else:
myfasta[sequence] = [lineage_id]这样，您的键(序列)将保存具有相同序列的lineage_ids列表。请注意，此解决方案的烦人之处在于遍历文件，将沿袭id与序列分离，说明可能扩展到多行的序列，等等。
在此之后，您可以遍历字典，并只使用字典中列表中的第一个lineage_id将序列写入文件。
页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接：
https://stackoverflow.com/questions/63092913
复制相关文章
相似问题
腾讯云开发者
