Код Perl для обратной комплементарности последовательности ДНК

Я пытался написать код на perl, чтобы получить обратное дополнение быстрых последовательностей ДНК в формате файла .fna. файл sequence02C.fna содержит сотни последовательностей ДНК:

>adbca3e
TGCTCCCCACGCTTGCCTCTCCAGTACTCAACCAAAGCAGTCTCTAGAAAAACAGTTTCCAACGCAATACGATGGAATTCCACTTCCCAAATATCTC
>4c2a958
TCCCCACGCTTTCGCGCTTCAGCGTCAGTATCTGTCCAGTGAGCTGACTTCTCCATCGGCATTCCTACACAGTACTCTAGAAAAACAGTTTCTGCTC
>0639b5b
TCGCGCCTCAGTGTCCAACGCAATACGAGTTGCAGACCAGGACACATGGAATTCCACTTCCCTCTCCAGTACTCAACCAAAGCAGTCTCTAGAAAAG

Я использовал следующую команду, которая может открыть файл и сделать реверс, но не показывает идентификатор последовательности (например: >adbca3e) в выводе.

Код:

#!/usr/local/perl

open (NS, "sequence02C.fna");
while (<NS>) {
    if ($_ =~ tr/ATGC/TACG/) {print $_;}
}

выходной файл является только дополнением к последовательности, но не обратным. Кроме того, он не содержит идентификаторов последовательностей ›adbca3e.

Может ли кто-нибудь предложить соответствующий код, чтобы сделать эту обратную комплементарность этой последовательности сразу и получить результат в выходной файл?


perl reverse complement dna-sequence
person mashuk    schedule 25.04.2021    source источник
comment
Может быть, вам следует объяснить разницу между дополнительным и обратным.   -  person TLP    schedule 25.04.2021
comment
Чтобы напечатать все строки, просто удалите условие if.   -  person TLP    schedule 25.04.2021

Ответы (2)


arrow_upward
2
arrow_downward

Вы печатаете только те строки, которые содержат A, T, G или C. Вы хотите напечатать каждую строку, поэтому печать не должна быть условной.

#!/usr/local/perl

use strict;               # Always
use warnings;             # Always

while (<>) {
    if (/^>/) {           # Only modify lines starting with ">".
       tr/ATGC/TACG/;
       $_ = reverse($_);  # You didn't reverse.
    }

    print;                # Print undconditionally.
}

(tr/// и print по умолчанию используют $_.)

Обратите внимание, я не открывал файл. Вы можете использовать программу следующим образом:

perl program.pl sequence02C.fna >sequence02C_revcomp.fna

or

perl -i~ program.pl sequence02C.fna

Последний изменяет файл на месте. (Будьте осторожны! Сначала проверьте его. Однако он делает резервную копию.)

person ikegami    schedule 25.04.2021
comment
Спасибо @ikegami - person mashuk; 26.04.2021

arrow_upward
1
arrow_downward

Вы говорите, что у вас есть программа для реверса, но она дает только комплементарность. Может быть, это очень очевидное описание для вас, но оно не очень понятно для меня.

Если под реверсом вы подразумеваете печать строки в обратном направлении, просто используйте функцию reverse. Дополнительным, я предполагаю, является использование соответствующих азотистых оснований, для чего и предназначена ваша транслитерация tr/ATGC/TACG/.

Чтобы исправить непечатаемые идентификаторы, просто удалите условие if в операторе печати.

Что бы я сделал, так это просто использовал оператор алмаза для небольшой программы, подобной этой:

use strict;
use warnings;
use feature 'say';

while (<>) {
    chomp;
    unless (/^>/) {
        tr/ATGC/TACG/;            # transliterate non-ids
        my $reverse = reverse;    # reverse $_
        say $reverse;             # do something with $reverse
    }
    say;          # print current line
}

Затем вы можете использовать эту программу следующим образом:

$ perl program.pl sequence02C.fna > output.txt
person TLP    schedule 25.04.2021
comment
Спасибо @TLP за ответ с кодами. Я использовал следующую программу, только с добавлением команды open к предложенной программе: open (NS, sequence02C.fna); в то время как (‹NS›) { жевать; но обнаружил вывод ошибки: [siddiquee@intron2 realtask]$ синтаксическая ошибка perl rc.pl в строке 13 rc.pl, рядом с моим глобальным символом $reverse требует явного имени пакета в строке 13 rc.pl. Глобальный символ $reverse требует явного имени пакета в строке 14 rc.pl. Выполнение rc.pl прервано из-за ошибок компиляции. - person mashuk; 25.04.2021
comment
Глобальный символ .. это когда вы не объявляете переменную с помощью my. В вашем случае похоже, что вы пропустили окончание строки точкой с запятой, так что my было повреждено. Однако ошибка не видна в коде, который вы показываете в комментарии. - person TLP; 26.04.2021
comment
@mashuk На самом деле в моем коде отсутствовала точка с запятой, добавил ее сейчас. - person TLP; 26.04.2021
comment
Большое спасибо! @TLP - person mashuk; 26.04.2021