faozhi has asked for the wisdom of the Perl Monks concerning the following question:

Hi, I have a file containing names of contigs eg. 
contig_188
contig_222
contig_19

[download]
I also have another genbank file which looks like 
//
LOCUS       A08_contig188            628 bp    DNA     linear   UNK
DEFINITION  Contig A08_contig131 from balbla
ACCESSION   unknown
FEATURES             Location/Qualifiers
     source          1..628
                     /mol_type="genomic DNA"
                     /db_xref="taxon: 6666666"
                     /genome_md5=""
                     /project="luke_6666666"
                     /genome_id="6666666.23054"
                     /organism="blabla"
     CDS             complement(144..614)
                     /db_xref="GO:0003723"
                     /db_xref="GO:0005525"
                     /db_xref="GO:0006605"
                     /db_xref="GO:0006614"
                     /db_xref="GO:0016020"
                     /translation="MGDIVSLVEKAAQDLGEENLKKAEENLKKGQFSME
+DYLSQLRQM
                     KKMGGIEGIMSFMPGVSKIKSQMDSAGIDESIITKNEAIILSMTKKERE
+NPKIIDGSR
                     KKRIANGSGSDPATINKLLKQFKMMSEMMKKMSKGNLKGMSDKGIPPEL
+FNQLK"
                     /product="Signal recognition particle, subunit Ff
+h SRP54
                     (TC 3.A.5.1.1)"
BASE COUNT      168 a    132 c     76 g    252 t
ORIGIN
        1 tagtagaggg ttaaaaaaac caaccttctc tataaatctt ccgtcccttg caaaa
+cgact
       61 gtcggcaaca acaaccttat aaacaggacg tttttttgtt ccacctcttg ataat
+cttaa
      121 ttttaacatt tcttctcctt tttttacttt aattggttaa atagctctgg cggta
+taccc
      181 ttgtcagaca tacctttcag atttccttta gacatctttt tcatcatttc agaca
+tcatt
      241 ttaaattgct taagcaattt attgatagtg gctggatctg agccagagcc attag
+caatt
      301 cttttttttc tagaaccatc tattattttt gggttctctc tttctttttt tgtca
+ttgat
      361 aaaataatag cttcattttt agtaattata ctttcgtcaa ttcccgctga gtcca
+tttga
      421 gatttaattt tagaaactcc tggcataaaa gacataattc cctcaatgcc accca
+ttttt
      481 ttcatctgtc ttaattgaga caaataatct tccatagaga attgaccttt tttta
+aattt
      541 tcctctgctt ttttaagatt ttcttctcct aaatcttgag ccgccttttc tacaa
+gggat
      601 acaatatccc ccatacctag tattctgt
//
LOCUS       A08_contig131            628 bp    DNA     linear   UNK
DEFINITION  Contig A08_contig131 from blabla
ACCESSION   unknown
FEATURES             Location/Qualifiers
     source          1..628
                     /mol_type="genomic DNA"
                     /db_xref="taxon: 6666666"
                     /genome_md5=""
                     /project="luke_6666666"
                     /genome_id="6666666.23054"
                     /organism="blabla"
     CDS             complement(144..614)
                     /db_xref="GO:0003723"
                     /db_xref="GO:0005525"
                     /db_xref="GO:0006605"
                     /db_xref="GO:0006614"
                     /db_xref="GO:0016020"
                     /translation="MGDIVSLVEKAAQDLGEENLKKAEENLKKGQFSME
+DYLSQLRQM
                     KKMGGIEGIMSFMPGVSKIKSQMDSAGIDESIITKNEAIILSMTKKERE
+NPKIIDGSR
                     KKRIANGSGSDPATINKLLKQFKMMSEMMKKMSKGNLKGMSDKGIPPEL
+FNQLK"
                     /product="Signal recognition particle, subunit Ff
+h SRP54
                     (TC 3.A.5.1.1)"
BASE COUNT      168 a    132 c     76 g    252 t
ORIGIN
        1 tagtagaggg ttaaaaaaac caaccttctc tataaatctt ccgtcccttg caaaa
+cgact
       61 gtcggcaaca acaaccttat aaacaggacg tttttttgtt ccacctcttg ataat
+cttaa
      121 ttttaacatt tcttctcctt tttttacttt aattggttaa atagctctgg cggta
+taccc
      181 ttgtcagaca tacctttcag atttccttta gacatctttt tcatcatttc agaca
+tcatt
      241 ttaaattgct taagcaattt attgatagtg gctggatctg agccagagcc attag
+caatt
      301 cttttttttc tagaaccatc tattattttt gggttctctc tttctttttt tgtca
+ttgat
      361 aaaataatag cttcattttt agtaattata ctttcgtcaa ttcccgctga gtcca
+tttga
      421 gatttaattt tagaaactcc tggcataaaa gacataattc cctcaatgcc accca
+ttttt
      481 ttcatctgtc ttaattgaga caaataatct tccatagaga attgaccttt tttta
+aattt
      541 tcctctgctt ttttaagatt ttcttctcct aaatcttgag ccgccttttc tacaa
+gggat
      601 acaatatccc ccatacctag tattctgt
//

[download]
I would like to grep only the genbank contig information which matches the first file (list of contig names) and write it into a new file.. Anyone can help?

Replies are listed 'Best First'.
Re: How do i extract only contigs of interest
by choroba (Cardinal) on Sep 09, 2015 at 11:03 UTC
    Searching for information taken from one file in another file is quite a common task. The usual solution is to hash the identifiers, then iterate over the second file and check the hash for the existence of the identifier.

    In your case, you have to change the identifiers, as the underscore in contig names is not present in the genbank: 

    #!/usr/bin/perl
use warnings;
use strict;

my %contig;
open my $CNTG, '<', 'contig_names' or die $!;
while (<$CNTG>) {
    chomp;
    s/_//;  # Remove the underscore, as it doesn't appear in the genba
+nk.
    $contig{$_} = 1;
}

local $/ = '//';  # Read the file one record a time.
open my $GB, '<', 'genbank' or die $!;
while (<$GB>) {
    if (my ($id) = /LOCUS \s+ \S+ _ (contig [0-9]+)/x) {
        print if $contig{$id};
    }
}

    [download]
    لսႽ† ᥲᥒ⚪⟊Ⴙᘓᖇ Ꮅᘓᖇ⎱ Ⴙᥲ𝇋ƙᘓᖇ
[reply]
[d/l]
      It is working now: this is my final code 
      #!/usr/bin/perl
use warnings;
use strict;

my %contig;
open (my $CNTG, '<', $ARGV[0]) or die $!;
while (<$CNTG>) {
    chomp;
#    s///;  # Remove the underscore, as it doesn't appear in the genba
+nk.
    $contig{$_} = 1;
#    print $_;
}

local $/ = '//';  # Read the file one record a time.
open (my $GB, '<', $ARGV[1]) or die $!;
while (<$GB>) {
    if (my ($id) = /LOCUS \s+ \S+ _ (contig [0-9]+)/x) {
        print if $contig{$id};
    }
}
~

      [download]
      THANKS!
[reply]
[d/l]
Re: How do i extract only contigs of interest
by Corion (Patriarch) on Sep 09, 2015 at 10:52 UTC

    You have shown us your input data, but we haven't seen the code you already have written yet. Please show us your code and where in your code you have problems. Most likely, filling a hash with the names of your contigs and then extracting only those parts for which an entry in the hash exists is a good approach.

[reply]

Back to Seekers of Perl Wisdom