Извлечь заданную пользователем последовательность из обратной цепи файла FASTA с помощью samtools
У меня есть список регионов с начальной и конечной точками.
Я использовал команду samtools faidx ref.fa <region>. Эта команда дала мне последовательность прямой цепи для этой области.
В руководстве по Samtools есть опция для извлечения обратной нити, но я не мог понять, как ее использовать.
Кто-нибудь знает, как запустить эту команду для обратной цепи в samtools?
Мои регионы похожи на:
LG2:124522-124572 (Forward)
LG3:250022-250072 (Reverse)
LG29:4822278-4822318 (Reverse)
LG12:2,595,915-2,596,240 (Forward)
LG16:5,405,500-5,405,828 (Reverse)
Ответы 2
Как вы заметили, samtools имеет опцию --reverse-complement (или -i) для вывода последовательности из обратной цепи.
Насколько мне известно, samtools не поддерживает обозначение региона, которое позволяет указывать нить.
Быстрое решение - разделить файл региона на прямое и обратное и дважды запустить samtools.
Приведенные ниже шаги довольно подробны, поэтому шаги понятны. Это довольно просто очистить, например, с помощью замены процесса в bash.
# Separate the strand regions.
# Use grep and sed twice, or awk (below).
grep -F '(Forward)' regions.txt | sed 's/ (Forward)//' > forward-regions.txt
grep -F '(Reverse)' regions.txt | sed 's/ (Reverse)//' > reverse-regions.txt
# Above as an awk one-liner.
awk '{ strand=($2 == "(Forward)") ? "forward" : "reverse"; print $1 > strand"-regions.txt" }' regions.txt
# Run samtools, marking the strand as +/- in the FASTA output.
samtools faidx ref.fa -r forward-regions.txt --mark-strand sign -o forward-sequences.fa
samtools faidx ref.fa -r reverse-regions.txt --mark-strand sign -o reverse-sequences.fa --reverse-complement
# Combine the FASTA output to a single file.
cat forward-sequences.fa reverse-sequences.fa > sequences.fa
rm forward-sequences.fa reverse-sequences.fa
просто хочу упомянуть, что вам, вероятно, потребуется обновить samtools до последней версии, если вы столкнулись с проблемой. В моем случае samtools V1.2 не работал, а V1.10 работал.