Pergunta

Basicamente, um GenBank arquivo consiste no gene entradas (anunciado pela 'gene', seguido por seu correspondente 'CDS' entrada (somente uma por gene), como os dois que eu mostrar aqui abaixo.Eu gostaria de obter locus_tag vs produto em uma guia, delimitado por duas colunas do arquivo.'gene' e 'CDS' são sempre precedido e seguido por espaços.

Uma pergunta anterior sugeriu um script.

O problema é que parece que porque o "produto" tem, por vezes, o carácter '/' dentro de seu nome, seus conflitos com este script, que, até onde posso entender, é usando '/' como separador de campo para armazenar informações em uma matriz?

Eu gostaria de resolver este problema, modificar esse script ou a construção de outro.

perl -nE'
  BEGIN{ ($/, $") = ("CDS", "\t") }
  say "@r[0,1]" if @r= m!/(?:locus_tag|product)="(.+?)"!g and @r>1
' file


 gene            complement(8972..9094)
                 /locus_tag="HAPS_0004"
                 /db_xref="GeneID:7278619"
 CDS             complement(8972..9094)
                 /locus_tag="HAPS_0004"
                 /codon_start=1
                 /transl_table=11
                 /product="hypothetical protein"
                 /protein_id="YP_002474657.1"
                 /db_xref="GI:219870282"
                 /db_xref="GeneID:7278619"
                 /translation="MYYKALAHFLPTLSTMQNILSKSPLSLDFRLLFLAFIDKR"
 gene            68..637
                 /locus_tag="HPNK_00040"
 CDS             68..637
                 /locus_tag="HPNK_00040"
                 /codon_start=1
                 /transl_table=11
                 /product="NinG recombination protein/bacteriophage lambda
                 NinG family protein"
                 /protein_id="CRESA:HPNK_00040"
                 /translation="MIKPKVKKRKCKCCGGEFKSADSFRKWCSAECGVKLAKIAQEKA
                 RQKAIEKRNREERAKIKATRERLKSRSEWLKDAQAIFNEYIRLRDKDEPCISCRRFHQ
                 GQYHAGHYRTVKAMPELRFNEDNVHKQCSACNNHLSGNITEYRINLVRKIGAERVEAL
                 ESYHPPVKWSVEDCKEIIKTYRAKIKELK"
Foi útil?

Solução

Como o exemplo do GenBank arquivo incompleto, fui on-line para encontrar um arquivo de exemplo que pode ser usado em um exemplo, e eu encontrei este arquivo.

Usando este código e o Bio::GenBankParser módulo, foi analisada a adivinhar o que as partes da estrutura que foram depois.Neste caso, "recursos" que continha uma locus_tag campo e um product de campo.

use strict;
use warnings;
use feature 'say';
use Bio::GenBankParser;

my $file = shift;
my $parser = Bio::GenBankParser->new( file => $file );
while ( my $seq = $parser->next_seq ) {
    my $feat = $seq->{'FEATURES'};
    for my $f (@$feat) {
        my $tag = $f->{'feature'}{'locus_tag'};
        my $prod = $f->{'feature'}{'product'};
        if (defined $tag and defined $prod) {
            say join "\t", $tag, $prod;
        }
    }
}

Uso:

perl script.pl input.txt > output.txt

Saída:

MG_001  DNA polymerase III, beta subunit
MG_470  CobQ/CobB/MinD/ParA nucleotide binding domain-containing protein

A saída do seu one-liner para a mesma entrada seria:

MG_001  DNA polymerase III, beta subunit
MG_470  CobQ/CobB/MinD/ParA nucleotide binding
                     domain-containing protein

Supondo, claro, que você adicionar o /s modificador para a regex para dar conta de várias linhas de entradas (o que leeduhem apontado nos comentários):

m!/(?:locus_tag|product)="(.+?)"!sg
#                                ^---- this

Outras dicas

Depois de ter lido o seu duplicado pergunta http://www.biostars.org/p/94164/ (por favor, não fazer double post como este), aqui está um mínimo de Biopython resposta:

import sys
from Bio import SeqIO
filename = sys.argv[1] # Takes first command line argument input filename
for record in SeqIO.parse(filename, "genbank"):
    for feature in record.features:
        if feature.type == "CDS":
            locus_tag = feature.qualifiers.get("locus_tag", ["???"])[0]
            product = feature.qualifiers.get("product", ["???"])[0]
            print("%s\t%s" % (locus_tag, product))

Com pequenas mudanças você pode escrever isso para um ficheiro em vez disso.

Licenciado em: CC-BY-SA com atribuição
Não afiliado a StackOverflow
scroll top