Чтение последовательностей .fasta для извлечения нуклеотидных данных, а затем запись в файл с ограничением табуляции

Прежде чем продолжить, я подумал, что порекомендую читателям ознакомиться с моими предыдущими проблемами с Perl, поскольку я был новичком во всем этом.

Это были мои сообщения за последние несколько дней в хронологическом порядке:

  1. Как сделать Я усредняю ​​значения столбцов из данных, разделенных табуляцией ... (решено)
  2. Почему я не вижу вычисленных результатов в моем выводе файл? (решено)
  3. Использование файла .fasta для вычисления относительного содержания последовательностей < / а>

Как я уже говорил выше, благодаря помощи некоторых из вас мне удалось выяснить первые два запроса, и я действительно извлек из этого урок. Я искренне благодарен. Для человека, который ничего об этом не знает и все еще чувствует, что не знает, помощь была практически находкой.

Последний вопрос остается нерешенным, и это продолжение. Я просмотрел некоторые из рекомендуемых текстов, но поскольку я пытаюсь закончить это до понедельника, я не уверен, что я что-то полностью упустил. В любом случае, я попытался выполнить эту задачу.

Просто чтобы вы знали, задача состоит в том, чтобы открыть и прочитать файл .fasta (я думаю, что наконец-то кое-что неплохо придумал, аллилуйя!), прочитать каждую последовательность, вычислить относительное содержание нуклеотидов G + C, а затем записать в файл TABDelimited имена генов и их соответствующее содержание G + C.

Несмотря на то, что я попытался это сделать, я знаю, что я далеко не готов выполнить программу, чтобы обеспечить результаты, которые мне нужны, поэтому я снова обращаюсь к вам, ребята, за некоторыми советами. , или примеры того, как это сделать. Как и в случае с моими предыдущими решенными запросами, я хотел бы, чтобы они были в стиле, аналогичном тому, в котором я их уже делал, - даже если это может быть не самый удобный / эффективный способ. Это просто позволяет мне знать, что я делаю на каждом этапе пути, даже если мне кажется, что я спамлю это!

В любом случае, файл .fasta читается примерно так:

>label
sequence
>label
sequence
>label
sequence

Я не уверен, как открыть файл .fasta, поэтому я не уверен, какие метки к каким относятся, но я знаю, что гены должны быть помечены как gag, pol или env. Нужно ли мне открывать файл .fasta, чтобы знать, что я делаю, или я могу сделать это «вслепую», выбрав указанный выше формат?

Это может быть совершенно очевидно, но я все еще борюсь со всем этим. Я чувствую, что уже должен был это понять!

Во всяком случае, текущий код у меня следующий:

#!/usr/bin/perl -w
# This script reads several sequences and computes the relative content of G+C of each sequence.
use strict; 

my $infile = "Lab1_seq.fasta";                               # This is the file path
open INFILE, $infile or die "Can't open $infile: $!";        # This opens file, but if file isn't there it mentions this will not open
my $outfile = "Lab1_SeqOutput.txt";             # This is the file's output
open OUTFILE, ">$outfile" or die "Cannot open $outfile: $!"; # This opens the output file, otherwise it mentions this will not open

my $sequence = ();  # This sequence variable stores the sequences from the .fasta file
my $GC = 0;         # This variable checks for G + C content

my $line;                             # This reads the input file one-line-at-a-time
while ($line = <INFILE>) {
    chomp $line;                      # This removes "\n" at the end of each line (this is invisible)

    foreach my $line ($infile) {
        if($line = ~/^\s*$/) {         # This finds lines with whitespaces from the beginning to the ending of the sequence. Removes blank line.
            next;
        } elsif($line = ~/^\s*#/) {        # This finds lines with spaces before the hash character. Removes .fasta comment
            next; 
        } elsif($line = ~/^>/) {           # This finds lines with the '>' symbol at beginning of label. Removes .fasta label
            next;
        } else {
            $sequence = $line;
        }
    }
    {
        $sequence =~ s/\s//g;               # Whitespace characters are removed
        return $sequence;
    }

Я не уверен, что здесь что-то правильно, но его выполнение оставило у меня синтаксическую ошибку ar line 35 (за последней строкой, и, следовательно, там ничего нет!). Об этом говорится на "EOF". Это все, что я могу отметить. В противном случае я пытаюсь выяснить, как вычислить количество нуклеотидов G + C в каждой из последовательностей, а затем правильно свести это в таблицу в выходном файле .txt. Я считаю, что именно это подразумевается под файлом TABDelimited?

В любом случае, прошу прощения, если этот запрос кажется слишком длинным, «тупым» или повторяющимся, но, сказав это, я не смог найти никакой информации, непосредственно относящейся к этому, поэтому ваша помощь будет очень признательна, а объяснения для каждого шага, если возможно !!

Добрый.


person PkC    schedule 17.03.2012    source источник


Ответы (1)


Прямо в конце у вас есть дополнительная скоба. Это должно работать:

#!/usr/bin/perl -w
# This script reads several sequences and computes the relative content of G+C of each sequence.

use strict; 

my $infile = "Lab1_seq.fasta";                               # This is the file path
open INFILE, $infile or die "Can't open $infile: $!";        # This opens file, but if file isn't there it mentions this will not open
my $outfile = "Lab1_SeqOutput.txt";             # This is the file's output
open OUTFILE, ">$outfile" or die "Cannot open $outfile: $!"; # This opens the output file, otherwise it mentions this will not open

my $sequence = ();  # This sequence variable stores the sequences from the .fasta file
my $GC = 0;         # This variable checks for G + C content

my $line;                             # This reads the input file one-line-at-a-time

while ($line = <INFILE>) {
    chomp $line;                      # This removes "\n" at the end of each line (this is invisible)

    if($line =~ /^\s*$/) {         # This finds lines with whitespaces from the beginning to the ending of the sequence. Removes blank line.
        next;

    } elsif($line =~ /^\s*#/) {        # This finds lines with spaces before the hash character. Removes .fasta comment
        next; 
    } elsif($line =~ /^>/) {           # This finds lines with the '>' symbol at beginning of label. Removes .fasta label
        next;
    } else {
        $sequence = $line;
    }

    $sequence =~ s/\s//g;               # Whitespace characters are removed
    print OUTFILE $sequence;
}

Также я редактировал вашу обратную строку. Возврат выйдет из вашего цикла. Я подозреваю, что вы хотите распечатать это в файл, поэтому я это сделал. Возможно, вам потребуется сначала выполнить некоторые дополнительные преобразования, чтобы преобразовать его в формат, разделенный табуляцией.

person Ilion    schedule 17.03.2012
comment
Вроде было выдано много ошибок. Интересно, не потому ли, что он неполный. В настоящее время я пытаюсь прочитать, как вычислить количество нуклеотидов в последовательностях. Я предполагаю, что для этого потребуется еще один цикл, но терминология все еще вне меня. - person PkC; 17.03.2012
comment
Извините, я пропустил проблему с строкой foreach. Цикл WHILE уже перебирает каждую строку из $ infile, поэтому цикл foreach для этого не нужен. Я скорректировал приведенный выше код, и это должно решить проблему. - person Ilion; 17.03.2012
comment
Hrm ваше регулярное выражение, ищущее комментарий, вероятно, тоже вызывает проблемы. Возможно, вам понадобится шанс } elsif($line =~ /^\s*#/) { на } elsif($line =~ qr(^\s*#/)) {. qr создает регулярное выражение в кавычках. - person Ilion; 17.03.2012
comment
Woohoo !! Теперь код печатает в текстовый файл всю последовательность без пробелов. Единственная проблема в том, что я не знаю, где последовательности начинались или заканчивались, поэтому я не уверен, какие последовательности применимы к каждому гену. Хотя кодон остановки / запуска должен дать мне указание. Принимая это во внимание, как мне изменить / добавить код, чтобы вычислить количество G + C в последовательностях, а затем распечатать их в файле с разделителями табуляции? - person PkC; 17.03.2012
comment
Что ж, мои познания в молекулярной биологии невелики, поэтому я не знаю, как проводить расчеты. В предоставленном вами коде есть переменная $ GC, но она нигде не используется. Рискну предположить, что после строки $sequence =~ s/\s//g; вы хотите как-то вычислить G + C последовательности $. Получив это, вы можете распечатать строку с помощью print OUTFILE "geneName\t$GC\n". \ T представляет собой TAB, а \ n - это NEWLINE. Если вам нужно просмотреть файлы в Windows, вы можете вместо этого использовать \ r \ n. - person Ilion; 18.03.2012