تحليل ملف GenBank:الحصول على علامة الموقع مقابل المنتج

Question

في الأساس، يتكون ملف GenBank من إدخالات جينية (يتم الإعلان عنها بواسطة "gene" متبوعة بإدخال "CDS" المقابل لها (واحد فقط لكل جين) مثل المدخلين اللذين أعرضهما هنا أدناه.أرغب في الحصول على locus_tag vs Product في ملف مكون من عمودين مفصول بعلامات جدولة.دائمًا ما يسبق "gene" و"CDS" وتليهما مسافات.

اقترح سؤال سابق نصًا.

تكمن المشكلة في أنه نظرًا لأن "المنتج" يحتوي أحيانًا على حرف "/" داخل اسمه، فإنه يتعارض مع هذا البرنامج النصي، وبقدر ما أستطيع أن أفهم، يستخدم "/" كفاصل حقل لتخزين المعلومات في ملف مجموعة مصفوفة؟

أرغب في حل هذه المشكلة، إما بتعديل هذا البرنامج النصي أو إنشاء برنامج آخر.

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"

Answer 1

نظرًا لأن نموذج ملف GenBank الخاص بك لم يكن مكتملًا، فقد ذهبت عبر الإنترنت للعثور على نموذج ملف يمكن استخدامه في مثال، ووجدت هذا الملف.

باستخدام هذا الرمز و Bio::GenBankParser الوحدة النمطية، تم تحليلها لتخمين أجزاء الهيكل التي كنت تبحث عنها.في هذه الحالة، "الميزات" التي تحتوي على كل من أ locus_tag الميدان و أ product مجال.

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;
        }
    }
}

الاستخدام:

perl script.pl input.txt > output.txt

انتاج:

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

سيكون الإخراج من سطر واحد لنفس الإدخال:

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

على افتراض بالطبع أنك قمت بإضافة /s معدّل إلى regex لحساب الإدخالات متعددة الأسطر (والتي leeduhem تمت الإشارة إليه في التعليقات):

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

Answer 2

Having read your duplicated question http://www.biostars.org/p/94164/ (please don't double post like this), here's a minimal Biopython answer:

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))

With minor changes you can write this out to a file instead.