Получить кодирующую аминокислоту, когда в последовательности ДНК есть определенный образец

Вот некоторый код, который должен хотя бы помочь вам начать. Например, вы можете запустить так:

./retrieve_coding_aa.pl file.fa ATAGTA

Содержание retrieve_coding_aa.pl:

#!/usr/bin/perl 

use strict;
use warnings;

use File::Basename;
use Bio::SeqIO;
use Bio::Tools::CodonTable;
use Data::Dumper;

my $pattern = $ARGV[1];

my $fasta = Bio::SeqIO->new ( -file => $ARGV[0], -format => 'fasta');

while (my $seq = $fasta->next_seq ) {

    my $pos = 0;

    my %counts;

    for (split /($pattern)/ => $seq->seq) {

        if ($_ eq $pattern) {

            my $dist = $pos % 3;

            unless ($dist == 0) {

                my $num = 3 - $dist;

                s/.{$num}//;

                chop until length () % 3 == 0;
            }

            my $table = Bio::Tools::CodonTable->new();

            $counts{$_}++ for split (//, $table->translate($_));
        }

        $pos += length;
    }

    print $seq->display_id() . ":\n";

    map {

        print "$_ => $counts{$_}\n"
    }
    sort {

        $counts{$a} <=> $counts{$b}
    }
    keys %counts;

    print "\n";
}

Вот результаты с использованием примера ввода:

sequence1:
S => 1

sequence2:
V => 1
I => 1

Класс Bio::Tools::CodonTable также поддерживает нестандартные таблицы использования кодонов. Вы можете изменить таблицу, используя указатель id. Например:

$table = Bio::Tools::CodonTable->new( -id => 5 );

or:

$table->id(5);

Для получения дополнительной информации, в том числе о том, как исследовать эти таблицы, см. документацию здесь: http://metacpan.org/pod/Bio%3a%3aTools%3a%3aCodonTable

Я буду придерживаться первой версии того, что вы хотели, потому что дополнение только еще больше запутало меня. (кадр?) Я нашел ATAGTA только один раз в sequence2, но я предполагаю, что вам также нужны зеркальные изображения / обратная последовательность, которая в этом случае будет ATGATA. Ну, мой скрипт этого не делает, так что вам придется прописать это дважды в файле input_sequences, но я думаю, это не должно быть проблемой.

Я работаю с файлом, подобным вашему, который я называю «dna.txt», и с файлом входных последовательностей, называемым «input_seq.txt». Файл результатов представляет собой список шаблонов и их вхождений в файле dna.txt (включая результаты перекрытия, но его можно настроить на отсутствие перекрытия, как описано в awk).

input_seq.txt:

GC
ATA
ATAGTA
ATGATA

ДНК.txt:

>sequence1
ATGGCGCATAGTAATGC
>sequence2
ATGATAGTAATGCGCGC

результаты.txt:

GC,6
ATA,2
ATAGTA,2
ATGATA,1

Код awk вызывает другой awk (но один из них простой). Вы должны запустить "./match_patterns.awk input_seq.txt", чтобы получить сгенерированный файл результатов.:

*match_patterns.awk:*

#! /bin/awk -f
{return_value= system("awk -vsubval="$1" -f test.awk dna.txt")}

test.awk:

#! /bin/awk -f
{string=$0
do
{
where = match(string, subval)
# code is for overlapping matches (i.e ATA matches twice in ATATAC)
# for non-overlapping replace +1 by +RLENGTH in following line
if (RSTART!=0){count++; string=substr(string,RSTART+1)}
}
while (RSTART != 0)
}
END{print subval","count >> "results.txt"}

Файлы должны быть все в одном каталоге.

Удачи!