Boulder::Genbank(3pm)

πŸ“¦
version 1.30 Β· 4 versions
, in debian bookworm, debian trixie, ubuntu noble, ubuntu resolute, ubuntu stonking

🌸 Boulder::Genbank 🌸

part of the User Contributed Perl Documentation, bestie πŸ’–


✨ NAME ✨

Boulder::Genbank β€” grab Genbank data records and hand them to you as parsed Boulder Stones. it's giving "biotech but make it cute" πŸ§¬πŸ’…


πŸ¦‹ SYNOPSIS πŸ¦‹

use Boulder::Genbank

# network access via Entrez
 $gb = Boulder::Genbank->newFh( qw(M57939 M28274 L36028) );

 while ($data = <$gb>) {
     print $data->Accession;

     @introns = $data->features->Intron;
     print "There are ",scalar(@introns)," introns.\n";
     $dna = $data->Sequence;
     print "The dna is ",length($dna)," bp long.\n";

     my @features = $data->features(-type=>[ qw(Exon Source Satellite) ],
                                    -pos=>[90,310] );
     foreach (@features) {
        print $_->Type,"\n";
        print $_->Position,"\n";
        print $_->Gene,"\n";
    }
  }

# another syntax
$gb = new Boulder::Genbank(-accessor=>'Entrez',
                           -fetch => [qw/M57939 M28274 L36028/]);

# local access via Yank
$gb = new Boulder::Genbank(-accessor=>'Yank',
                           -fetch=>[qw/M57939 M28274 L36028/]);
while (my $s = $gb->get) {
   # etc.
}

# parse a file of Genbank records
$gb = new Boulder::Genbank(-accessor=>'File',
                           -fetch => '/usr/local/db/gbpri3.seq');
while (my $s = $gb->get) {
   # etc.
}

# parse flatfile records yourself
open (GB,"/usr/local/db/gbpri3.seq");
local $/ = "//\n";
while (<GB>) {
   my $s = Boulder::Genbank->parse($_);
   # etc.
}

πŸ’– DESCRIPTION πŸ’–

Boulder::Genbank is your girl for retrieval and parsing of NCBI Genbank-format records. she hands the entries back to you in Stone format so you can just... reach in and grab whatever field or value you want, no stress ✨ Boulder::Genbank is a descendent of Boulder::Stream, so it gives you a stream-like interface over a whole series of Stone objects. it's a whole vibe.

>> IMPORTANT NOTE << As of January 2002, NCBI has changed their Batch Entrez interface. I have modified Boulder::Genbank so as to use a "demo" interface, which fixes things, but this isn't guaranteed in the long run. I have written to NCBI, and they may fix this -- or they may not. 😭 >> IMPORTANT NOTE <<

Genbank access is served up by three different accessors, and together they cover both remote and local Genbank databases 🫢 When you create a new Boulder::Genbank stream you pick one of the three accessors, plus whatever accessor-specific parameters control which entries get fetched. the three accessors, bestie, are:

  • Entrez β€” this one gives you access to NetEntrez, pulling the most recent Genbank info straight from NCBI's Web site 🌐 The parameters you pass are either a series of Genbank accession numbers, or an Entrez query (peep http://www.ncbi.nlm.nih.gov/Entrez/linking.html). If you hand it a list of accession numbers, the stream returns a series of stones corresponding to those numbers. Otherwise, if you gave it an Entrez query, the entries come back in the order returned by Entez. (that's their spelling, not mine πŸ’€)
  • File β€” this one gets you local Genbank entries by reading from a flat file (usually one of the .seq files you can download from NCBI's Web site πŸ“). The stream returns a Stone for each entry in the file, starting from the top of the file and working downward. The parameter here is the path to the local file.
  • Yank β€” local Genbank entries via Will Fitzhugh's Yank program. Yank gives fast indexed access to a Genbank flat file using the accession number as the key ⚑ The parameter you pass is a list of accession numbers. Stones come back in the order you requested. By default the yank binary lives in /usr/local/bin/yank. To support other locations, you may define the environment variable YANK to contain the full path.

You can also parse a single Genbank entry out of a text string sitting in a scalar variable, and get a Stone object back πŸ’…

πŸŽ€ Boulder::Genbank methods πŸŽ€

This section lists the public methods that the Boulder::Genbank class makes available.

new()

# Network fetch via Entrez, with accession numbers
$gb=new Boulder::Genbank(-accessor  =>  'Entrez',
                         -fetch     =>  [qw/M57939 M28274 L36028/]);

# Same, but shorter and uses -> operator
$gb = Boulder::Genbank->new qw(M57939 M28274 L36028);

# Network fetch via Entrez, with a query

# Network fetch via Entrez, with a query
$query = 'Homo sapiens[Organism] AND EST[Keyword]';
$gb=new Boulder::Genbank(-accessor  =>  'Entrez',
                         -fetch     =>  $query);

# Local fetch via Yank, with accession numbers
$gb=new Boulder::Genbank(-accessor  =>  'Yank',
                         -fetch     =>  [qw/M57939 M28274 L36028/]);

# Local fetch via File
$gb=new Boulder::Genbank(-accessor  =>  'File',
                         -fetch     =>  '/usr/local/genbank/gbpri3.seq');

The new() method creates a new Boulder::Genbank stream on the accessor you provide. The three possible accessors are Entrez, Yank and File. If it works, the method returns the stream object πŸŽ‰ Otherwise it returns undef πŸ˜”

new() takes the following arguments:

        -accessor       Name of the accessor to use
        -fetch          Parameters to pass to the accessor
        -proxy          Path to an HTTP proxy, used when using
                         the Entrez accessor over a firewall.

Specify the accessor with the -accessor argument. If you don't specify it, it defaults to Entrez.

-fetch is an accessor-specific argument. The possibilities are:

  • For Entrez, the -fetch argument may point to a scalar, in which case it's interpreted as an Entrez query string. See http://www.ncbi.nlm.nih.gov/Entrez/linking.html for a description of the query syntax. Alternatively, -fetch may point to an array reference, in which case it's interpreted as a list of accession numbers to retrieve. If -fetch points to a hash, it's interpreted as extended information. See "Extended Entrez Parameters" below.
  • For Yank, the -fetch argument must point to an array reference containing the accession numbers to retrieve.
  • For File, the -fetch argument must point to a string-valued scalar, which will be interpreted as the path to the file to read Genbank entries from.

For Entrez (and Entrez only, iykyk) Boulder::Genbank lets you use a shortcut syntax where you hand new() a list of accession numbers:

  $gb = new Boulder::Genbank('M57939','M28274','L36028');

newFh()

This works like new(), but returns a filehandle 🫢 To recover each GenBank record read from the filehandle with the <> operator:

  $fh = Boulder::GenBank->newFh('M57939','M28274','L36028');
  while ($record = <$fh>) {
     print $record->asString;
  }

get()

The get() method is inherited from Boulder::Stream, and simply returns the next parsed Genbank Stone, or undef if there's nothing more to fetch. It has the same semantics as the parent class, including the ability to restrict access to certain top-level tags πŸ”’

The object returned is a Stone::GB_Sequence object, which is a descendent of Stone.

put()

The put() method is inherited from the parent Boulder::Stream class, and will write the passed Stone to standard output in Boulder format. Which means: it is currently not possible to write a Boulder::Genbank object back into Genbank flatfile form. 😭

🌈 Extended Entrez Parameters 🌈

The Entrez accessor recognizes extended parameters that let you customize the search ✨ Instead of passing a query string scalar or a list of accession numbers as the -fetch argument, pass a hash reference. That hashref should contain one or more of these keys:

  • -query β€” The Entrez query to process.
  • -accession β€” The list of accession numbers to fetch, as an array ref.
  • -db β€” The database to search. This is a single-letter database code selected from the following list:
  m  MEDLINE
  p  Protein
  n  Nucleotide
  s  Popset
  • -proxy β€” An HTTP proxy to use. For example:
   -proxy => http://www.firewall.com:9000

If you think you need this, get the correct URL from your system administrator. seriously, ask them, don't guess πŸ™ˆ

As an example, here's how to search for ESTs from Oryza sativa that have been entered or modified since 1999:

  my $gb = new Boulder::Genbank( -accessor=>Entrez,
                                 -query=>'Oryza sativa[Organism] AND EST[Keyword] AND 1999[MDAT]',
                                 -db   => 'n'
                                });

πŸ‘‘ METHODS DEFINED BY THE GENBANK STONE OBJECT πŸ‘‘

Every record returned from the Boulder::Genbank stream defines a set of methods that correspond to features and other fields in the Genbank flat file record πŸ’… Stone::GB_Sequence gives the full details, but they're listed here for reference:

⭐ $length = $entry->length

Get the length of the sequence.

⭐ $start = $entry->start

Get the start position of the sequence, currently always "1".

⭐ $end = $entry->end

Get the end position of the sequence, currently always the same as the length.

⭐ @feature_list = $entry->features(-pos=>[50,450],-type=>['CDS','Exon'])

features() will search the entry feature list for those features that meet certain criteria πŸ•΅οΈβ€β™€οΈ The criteria are specified using the -pos and/or -type argument names, as shown below.

  • -pos β€” Provide a position or range of positions which the feature must overlap. A single position is specified in this way:
   -pos => 1500;         # feature must overlap postion 1500

or a range of positions in this way:

   -pos => [1000,1500];  # 1000 to 1500 inclusive

If no criteria are provided, then features() returns all the features, and is equivalent to calling the Features() accessor.

  • -type, -types β€” Filter the list of features by type or a set of types. Matches are case-insensitive, so "exon", "Exon" and "EXON" are all equivalent. You may call with a single type as in:
   -type => 'Exon'

or with a list of types, as in

   -types => ['Exon','CDS']

The names "-type" and "-types" can be used interchangeably.

⭐ $seqObj = $entry->bioSeq;

Returns a Bio::Seq object from the Bioperl project 🧬 Dies with an error message unless the Bio::Seq module is installed.


✨ OUTPUT TAGS ✨

okay so bestie, the tags that come back from the parsing operation? they're straight-up lifted from the NCBI ASN.1 schema πŸ’–. and for consistency's sake they get normalized β€” so the first letter is Capitalized and every letter after that is lowercase. no cap. this section is an abbreviated list of the most useful / most common tags, so like, it's giving greatest hits, not the whole discography. for the full listing you gotta go read "The NCBI Data Model" by James Ostell and Jonathan Kans in "Bioinformatics: A Practical Guide to the Analysis of Genes and Proteins" (Eds. A. Baxevanis and F. Ouellette), pp 121-144. πŸ’…


🌸 Top-Level Tags 🌸

these are the tags that show up at the top level of the parsed Genbank entry. the main characters, fr fr.

Accession

the accession number of this entry πŸŽ€. but because of the vagaries of the Genbank data model (their words, not mine), an entry can have multiple accession numbers β€” like, say, after a merging operation. so Accession may therefore be a multi-valued tag. she contains multitudes.

Example:

my $accessionNo = $s->Accession;

Authors

the list of authors, exactly as they appear on the AUTHORS line of the Genbank record. and nope, no attempt is made to parse them into individual authors. they stay a group. periodt.

Basecount

the nucleotide basecount for the entry 🧬. it's presented as a Boulder Stone with keys a, c, t and g. Example:

my $A = $s->Basecount->A;
my $C = $s->Basecount->C;
my $G = $s->Basecount->G;
my $T = $s->Basecount->T;
print "GC content is ",($G+$C)/($A+$C+$G+$T),"\n";

Blob

the entire flatfile record as an unparsed chunk of text β€” a "blob" 😭. this is the handy way of reassembling the record so a human can actually look at it. iconic for inspection.

Comment

the COMMENT line from the Genbank record. simple. clean. we love her.

Definition

the DEFINITION line from the Genbank record, unmodified. no notes, no edits, we don't touch it.

Features

the FEATURES table πŸ”₯. this one is a complex stone object with multiple subtags, so she gets her own section down below. see "The Features Tag" for the details, bestie.

Journal

the JOURNAL line from the Genbank record, unmodified. 🫢

Keywords

the KEYWORDS line from the Genbank record, unmodified. and no, no attempt is made to parse the keywords into separate values β€” they stay bundled.

Example:

my $keywords = $s->Keywords

Locus

the LOCUS line from the Genbank record. it is not further parsed. that's it, that's the whole tag.

Medline, Nid

references to other database accession numbers πŸ’œ.

Organism

the taxonomic name of the organism this entry was derived from πŸ¦‹. this line is taken from the Genbank entry unmodified. if you want the tea on their taxonomic syntax, go read the NCBI data model documentation.

Reference

the REFERENCE line from the Genbank entry. and there are often multiple Reference lines, so heads up. Example:

my @references = $s->Reference;

Sequence

the DNA or RNA sequence of the entry 🌈. it's presented as a single lower-case string, with all base numbers and formatting characters removed. just the goods.

Source

the entry's SOURCE field; often giving clues on how the sequencing was performed πŸ‘€.

Title

the TITLE field from the paper describing this entry, if any. if there's no paper, there's no title, and that's valid.


πŸŽ€ The Features Tag πŸŽ€

so the Features tag points to a Stone record that contains multiple subtags. each subtag is the name of a feature, and that feature points, in turn, to a Stone that describes the feature's location and other attributes. the full list of features is beyond this document (iykyk), but here are the ones you see most often:

  • Cds β€” a CDS
  • Intron β€” an intron
  • Exon β€” an exon
  • Gene β€” a gene
  • Mrna β€” an mRNA
  • Polya_site β€” a putative polyadenylation signal
  • Repeat_unit β€” a repetitive region
  • Source β€” More information about the organism and cell type the sequence was derived from
  • Satellite β€” a microsatellite (dinucleotide repeat)

each feature will contain one or more of the following subtags, bestie:

DB_xref

a cross-reference to another database, in the form DB_NAME:accession_number. see the NCBI Web site for a description of these cross references. πŸ’«

Evidence

the evidence for this feature β€” either "experimental" or "predicted". those are the two options. no third vibe.

Gene

if the feature involves a gene, this will be the gene's name (or one of its names). this subtag is often seen in "Gene" and Cds features. πŸ‘‘

Example:

foreach ($s->Features->Cds) {
   my $gene = $_->Gene;
   my $position = $_->Position;
   Print "Gene $gene ($position)\n";
}

Map

if the feature is mapped, this provides a map position, usually as a cytogenetic band. πŸ—ΊοΈ

Note

a grab-back for various text notes. literally the junk drawer of the feature, and we respect her for it.

Number

when multiple features of this type occur, this field is used to number them. ordinarily this field is not needed, because Boulder::Genbank preserves the order of features. so like, you usually don't need her, but she's there. πŸ’…

Organism

if the feature is Source, this provides the source organism. 🌱

Position

the position of this feature, usually expresed as a range (1970..1975). (yes, "expresed", that's the original's typo, we're leaving it, it's part of the lore πŸŽ€)

Product

the protein product of the feature, if applicable, as a text string. 🧬

Translation

the protein translation of the feature, if applicable. the glow-up from sequence to protein, honestly.


πŸ’– SEE ALSO πŸ’–

Boulder, Boulder::Blast


πŸ‘‘ AUTHOR πŸ‘‘

Lincoln Stein <lstein@cshl.org>.

Copyright (c) 1997-2000 Lincoln D. Stein

This library is free software; you can redistribute it and/or modify it under the same terms as Perl itself. See DISCLAIMER.txt for disclaimers of warranty.


⭐ EXAMPLE GENBANK OBJECT ⭐

okay bestie, here's the moment. the following is an excerpt from a moderately complex Genbank Stone 🫢. heads up: the Sequence line and several other long lines have been truncated for readability. so don't come for me about the ..., they're intentional.

 Authors=Spritz,R.A., Strunk,K., Surowy,C.S.O., Hoch,S., Barton,D.E. and Francke,U.
 Authors=Spritz,R.A., Strunk,K., Surowy,C.S. and Mohrenweiser,H.W.
 Locus=HUMRNP7011   2155 bp    DNA             PRI       03-JUL-1991
 Accession=M57939
 Accession=J04772
 Accession=M57733
 Keywords=ribonucleoprotein antigen.
 Sequence=aagcttttccaggcagtgcgagatagaggagcgcttgagaaggcaggttttgcagcagacggcagtgacagcccag...
 Definition=Human small nuclear ribonucleoprotein (U1-70K) gene, exon 10 and 11.
 Journal=Nucleic Acids Res. 15, 10373-10391 (1987)
 Journal=Genomics 8, 371-379 (1990)
 Nid=g337441
 Medline=88096573
 Medline=91065657
 Features={
   Polya_site={
     Evidence=experimental
     Position=1989
     Gene=U1-70K
   }
   Polya_site={
     Position=1990
     Gene=U1-70K
   }
   Polya_site={
     Evidence=experimental
     Position=1992
     Gene=U1-70K
   }
   Polya_site={
     Evidence=experimental
     Position=1998
     Gene=U1-70K
   }
   Source={
     Organism=Homo sapiens
     Db_xref=taxon:9606
     Position=1..2155
     Map=19q13.3
   }
   Cds={
     Codon_start=1
     Product=ribonucleoprotein antigen
     Db_xref=PID:g337445
     Position=join(M57929:329..475,M57930:183..245,M57930:358..412, ...
     Gene=U1-70K
     Translation=MTQFLPPNLLALFAPRDPIPYLPPLEKLPHEKHHNQPYCGIAPYIREFEDPRDAPPPTR...
   }
   Cds={
     Codon_start=1
     Product=ribonucleoprotein antigen
     Db_xref=PID:g337444
     Evidence=experimental
     Position=join(M57929:329..475,M57930:183..245,M57930:358..412, ...
     Gene=U1-70K
     Translation=MTQFLPPNLLALFAPRDPIPYLPPLEKLPHEKHHNQPYCGIAPYIREFEDPR...
   }
   Polya_signal={
     Position=1970..1975
     Note=putative
     Gene=U1-70K
   }
   Intron={
     Evidence=experimental
     Position=1100..1208
     Gene=U1-70K
   }
   Intron={
     Number=10
     Evidence=experimental
     Position=1100..1181
     Gene=U1-70K
   }
   Intron={
     Number=9
     Evidence=experimental
     Position=order(M57937:702..921,1..1011)
     Note=2.1 kb gap
     Gene=U1-70K
   }
   Intron={
     Position=order(M57935:272..406,M57936:1..284,M57937:1..599, <1..>1208)
     Gene=U1-70K
   }
   Intron={
     Evidence=experimental
     Position=order(M57935:284..406,M57936:1..284,M57937:1..599, <1..>1208)
     Note=first gap-0.14 kb, second gap-0.62 kb
     Gene=U1-70K
   }
   Intron={
     Number=8
     Evidence=experimental
     Position=order(M57935:272..406,M57936:1..284,M57937:1..599, <1..>1181)
     Note=first gap-0.14 kb, second gap-0.62 kb
     Gene=U1-70K
   }
   Exon={
     Number=10
     Evidence=experimental
     Position=1012..1099
     Gene=U1-70K
   }
   Exon={
     Number=11
     Evidence=experimental
     Position=1182..(1989.1998)
     Gene=U1-70K
   }
   Exon={
     Evidence=experimental
     Position=1209..(1989.1998)
     Gene=U1-70K
   }
   Mrna={
     Product=ribonucleoprotein antigen
     Position=join(M57928:358..668,M57929:319..475,M57930:183..245, ...
     Gene=U1-70K
   }
   Mrna={
     Product=ribonucleoprotein antigen
     Citation=[2]
     Evidence=experimental
     Position=join(M57928:358..668,M57929:319..475,M57930:183..245, ...
     Gene=U1-70K
   }
   Gene={
     Position=join(M57928:207..719,M57929:1..562,M57930:1..577, ...
     Gene=U1-70K
   }
 }
 Reference=1  (sites)
 Reference=2  (bases 1 to 2155)
 =

and that's the whole page, bestie, we made it πŸŽ‰πŸ’– this one came to you courtesy of perl v5.34.0, dated 2022-06-08 β€” so if you're reading this from the future, hi, you're valid, the docs are just vintage. go parse something cute. slay πŸ’…βœ¨

πŸ“¦ where this page comes from