GMOD_B...AD_GFF3(1p)

📦
version 1.31 · 3 versions
, in debian bookworm, debian trixie, ubuntu jammy, ubuntu noble

💖 GMOD_B...AD_GFF3 💖

User Contributed Perl Documentation · section 1p

ok bestie, sit down, hydrate, and let me tell you about the script that bulk-loads GFF3 files straight into a chado database 🧬✨

🌸 NAME 🌸

$0 — Bulk loads gff3 files into a chado database. no cap, that's the whole job 💅

✨ SYNOPSIS ✨

  % $0 [options]
  % cat <gff-file> | $0 [options]

two ways in: hand it a file, or pipe a file into it like the icon you are 🫶

🎀 OPTIONS 🎀

  • --gfffile — the file holding your GFF3 (it's optional, because you can totally just read from stdin instead 💌)
  • --fastafile — the Fasta file to pull sequence from
  • --organism — the organism this data belongs to. (pro tip: pass the literal value 'fromdata' and it'll read the organism out of the GFF's organism=xxx instead 👀)
  • --dbprofile — database config profile name
  • --dbname — database name
  • --dbuser — database user name
  • --dbpass — database password
  • --dbhost — database host
  • --dbport — database port
  • --analysis — signals that this GFF data came from a computational analysis
  • --noload — build the bulk load files, but do not actually load them. look, don't touch 🫣
  • --nosequence — don't load sequence, even if it's sitting right there in the file
  • --notransact — don't wrap the database load in a single transaction
  • --drop_indexes — drop the indexes on the affected tables before the load starts, then recreate them once it's done. (heads up: generally this does not help performance, so don't get your hopes up 💀)
  • --validate — validate SOFA terms before even attempting the insert. this can make script startup slow, and it's off by default
  • --ontology — give directions for handling misc Ontology_terms
  • --skip_vacuum — skip vacuuming the tables after the inserts (this is the default)
  • --no_skip_vaccum — don't skip vacuuming the tables (yes, that spelling is on purpose, leave it alone 🙈)
  • --inserts — print INSERT statements instead of COPY FROM STDIN
  • --noexon — don't convert CDS features into exons (but polypeptide features still get created, don't worry)
  • --recreate_cache — causes the uniquename cache to be recreated
  • --remove_lock — remove the lock so a new process can run
  • --save_tmpfiles — keep the temp files used for loading the database
  • --random_tmp_dir — use a randomly generated tmp dir (the default is to use the current directory)
  • --no_target_syn — by default the loader drops the targetId into the synonyms list of the feature. this flag deactivates that.
  • --unique_target — trust the unicity of the target IDs. IDs are case sensitive. by default, the uniquename of a new target is 'TargetId_PrimaryKey'. with this flag it becomes 'TargetId'. and on top of that, the Name of the created target will be its TargetId instead of the feature's Name.
  • --dbxref — use either the first Dbxref annotation as the primary dbxref (the one that goes into feature.dbxref_id), or — if you hand it an optional argument — use the first dbxref whose database part (ie, everything before the ':') matches the pattern you supplied.
  • --delete — instead of inserting features, use the GFF lines to delete features, exactly as if the CRUD=delete-all option were set on every line (peep 'Deletes and updates via GFF below'). the loader will ask for confirmation before continuing.
  • --delete_i_really_mean_it — same as --delete, except it does not ask for confirmation. you said what you said 😭
  • --fp_cv — name of the feature property controlled vocabulary (defaults to 'feature_property').
  • --noaddfpcv — by default the loader adds GFF attribute types as new feature_property cv terms whenever they're missing. this flag deactivates that.
** dgg note: should rename this flag: --[no]autoupdate for Chado tables cvterm, cv, db, organism, analysis ... **
  • --manual — detailed manual pages
  • --custom_adapter — use a custom subclass adaptor for Bio::GMOD::DB::Adapter. hand over the path to the adapter as an argument.
  • --private_schema — load the data into a non-public schema.
  • --use_public_cv — when loading into a non-public schema, load any cv and cvterm data into the public schema instead
  • --end_sql — SQL code to execute after the data load is all done
  • --allow_external_parent — allow Parent tags to point at IDs that live outside the current GFF file

ok so, important housekeeping 🧹 — every argument that starts with db, plus organism, can be supplied by default from Bio::GMOD::Config, which got installed back when 'make install' was run. and also: dbprofile and all the other db* options are mutually exclusive — if you supply dbprofile, do not supply any other db* options, because they simply won't be used. pick a lane bestie 💅

🦋 DESCRIPTION 🦋

the GFF in your datafile must be version 3, because GFF3 has tighter control over the spec and actually uses controlled vocabulary. which means: the names of feature types have to be exactly the ones in the Sequence Ontology Feature Annotation (SOFA) — not the synonyms, not the accession numbers. (SO accession numbers might get supported in future versions of this script, so like… watch this space ⭐)

also note that the ##sequence-region directive is not supported as a way of declaring a reference sequence for a GFF3 file. why? because ##sequence-region simply isn't expressive enough to say what kind of thing the sequence is — is it a chromosome? a contig? an arm? a vibe? nobody knows 🫠. if your GFF file uses ##sequence-region this way, you gotta convert it into a full GFF3 line. for example, if you have this line:

  ##sequence-region chrI 1 9999999

then it should be converted into a GFF3 line like this:

  chrI  .       chromosome      1       9999999 .       .       .       ID=chrI

🌈 How GFF3 is stored in chado 🌈

here's the summary of how GFF3 data ends up living in chado, column by column:

Column 1 (reference sequence)

the reference sequence for the feature becomes the srcfeature_id of the feature in the featureloc table for that feature. that featureloc generally gets a rank of zero if there are other locations associated with this feature (for instance, for a match feature); the other locations get featureloc.rank values greater than zero.

Column 2 (source)

the source is stored as a dbxref. the chado instance must have an entry in the db table named 'GFF_source'. the script then creates a dbxref entry for the feature's source and associates it to the feature through the feature_dbxref table.

Column 3 (type)

the cvterm.cvterm_id of the SOFA type goes into feature.type_id.

Column 4 (start)

the value of start minus 1 is stored in featureloc.fmin (one gets subtracted because chado uses interbase coordinates, while GFF uses base coordinates — math is math 🤓).

Column 5 (end)

the value of end is stored in featureloc.fmax.

Column 6 (score)

the score goes into one of the score columns in the analysisfeature table. the default is analysisfeature.significance. (there's a whole section below on analysis results if you want the tea ☕)

Column 7 (strand)

the strand is stored in featureloc.strand.

Column 8 (phase)

the phase is stored in featureloc.phase. heads up: there's currently a problem with the chado schema for the case of single exons having different phases in different transcripts. if your data has exactly that case, go complain to gmod-schema@lists.sourceforge.net so people can find ways to address it. sorry bestie 😭

Column 9 (group)

here is where the magic happens ✨

Assigning feature.name, feature.uniquename

the values of feature.name and feature.uniquename are assigned according to these simple rules:

  • if there's an ID tag, that gets used as feature.uniquename
  • otherwise, it's assigned a uniquename equal to 'auto' concatenated with the feature_id.
  • if there's a Name tag, its value is set to feature.name; otherwise it's null.

note that these rules are much more simple than the ones Bio::DB::GFF uses, and may need to be revisited. it's giving minimal effort, honestly 💀

Assigning feature_relationship entries

all Parent-tagged features get feature_relationship entries of 'part_of' pointing at their parent features. Derived_from tags get 'derived_from' relationships. and important: parent features must appear in the file before any features use Parent or Derived_from tags referring to that feature. order matters, no skipping ahead 🚫

Alias tags

Alias values are stored in the synonym table, under both synonym.name and synonym.synonym_sgml, and are linked to the feature through the feature_synonym table.

Dbxref tags

Dbxref values must be of the form 'db_name:accession', where db_name must have an entry in the db table with a value of db.name equal to 'DB:db_name'. several database names were preinstalled with the database back when 'make prepdb' was run. run 'SELECT name FROM db' to find out which databases are already available. new dbxref entries get created in the dbxref table, and dbxrefs are linked to features through the feature_dbxref table.

Gap tags

currently mostly ignored — the value is stored as a featureprop, but otherwise it's not used yet. it's giving placeholder 🫥

Note tags

the values are stored as featureprop entries for the feature.

Any custom (ie, lowercase-first) tags

custom tags are supported. if the tag doesn't already exist in the cvterm table, it gets created. the value gets stored with its associated cvterm in the featureprop table. so go off, make your own tags 💅

Ontology_term

when Ontology_term tags are used, items from the Gene Ontology and Sequence Ontology get processed automatically whenever the standard DB:accession format is used (e.g. GO:0001234). to use other ontology terms, you must specify the mapping of the DB identifiers in the GFF file and the names of the ontologies in the cv table as a comma separated list of tag=value pairs. for example, to use plant and cell ontology terms, you'd supply this on the command line:

  --ontology 'PO=plant ontology,CL=cell ontology'

where 'plant ontology' and 'cell ontology' are the names in the cv table, exactly as they appear. spelling counts, bestie 🔍

Target tags

okay so Target tags are lowkey extra to process properly. For them to work you need two source features already chillin' in the database first:

  • the 'primary' source feature (your chromosome or contig bestie) 🧬
  • the 'subject' from the similarity analysis — think an EST, a cDNA, or a syntenic chromosome

If the subject feature isn't there? The loader will try to conjure up a placeholder feature object in its place. ✨

And if you've got a fasta file that contains the subject, you can use the perl script gmod_fasta2gff3.pl (it ships with this distribution 💅) to build a GFF3 file that's suitable for loading into chado before you load your analysis results. Order matters, bestie.

CDS and UTR features

The way CDS features are represented in Chado is as an intersection of a transcript's exons and the transcript's polypeptide feature. Wild, right? 🦋

So to allow proper translation of a GFF3 file's CDS features, this loader will:

  • convert CDS and UTR feature lines to corresponding exon features (and add a featureprop note saying the exon was inferred from a GFF3 CDS and/or UTR line)
  • create a polypeptide feature that spans the genomic region from the start of translation to the stop

If your GFF3 file contains both exon and CDS/UTR features, then you'll want to suppress the creation of the exon features and instead only want a polypeptide feature to be created. To do this, use the --noexon option. In this case, the CDS and UTR features will still be converted to exon features as described above. No cap. 👑

Note that in the case where your GFF file contains CDS and/or UTR features that do not belong to 'central dogma' genes (that is, that have a gene, transcript and CDS/exon features), none of the above will happen and the features will be stored as is.

✨ NOTES ✨

Loading fasta file

When the --fastafile is provided with an argument that is the path to a file containing fasta sequence, the loader will attempt to update the feature table with the sequence provided. 🎀

Heads up! The ID provided in the fasta description line must exactly match what is in the feature table uniquename field. Be careful if it's possible that the uniquename of the feature was changed to ensure uniqueness when it was loaded from the original GFF. 😭

Also note that when loading sequence from a fasta file, loading GFF from standard in is disabled. Sorry for any inconvenience. 🙈

##sequence-region

This script does not use sequence-region directives for anything. If it represents a feature that needs to be inserted into the database, it should be represented with a full GFF line. This includes the reference sequence for the features if it is not already in the database, like a chromosome. For example, this:

##sequence-region chr1 1      213456789

should change to this:

chr1  UCSC    chromosome      1       213456789       .       .       .       ID=chr1

Transactions

This application will, by default, try to load all of the data at once as a single transcation. This is safer from the database's point of view, since if anything bad happens during the load, the transaction will be rolled back and the database will be untouched. 💖

The problem occurs if there are many (say, greater than a 2-300,000) rows in the GFF file. When that is the case, doing the load as a single transcation can result in the machine running out of memory and killing processes. 💀 If --notranscat is provided on the commandline, each table will be loaded as a separate transaction.

SQL INSERTs versus COPY FROM

This bulk loader was originally designed to use the PostgreSQL COPY FROM syntax for bulk loading of data. However, as mentioned in the 'Transactions' section, memory issues can sometimes interfere with such bulk loads. In another effort to circumvent this issue, the bulk loader has been modified to optionally create INSERT statements instead of the COPY FROM statements. 🔥

INSERT statements will load much more slowly than COPY FROM statements, but as they load and commit individually, they are more more likely to complete successfully. As an indication of the speed differences involved, loading yeast GFF3 annotations (about 16K rows), it takes about 5 times longer using INSERTs versus COPY on my laptop. 😮‍💨

Deletes and updates via GFF

There is rudimentary support for modifying the features in an existing database via GFF. Currently, there is only support for deleting. 🗑️

In order to delete, the GFF line must have a custom tag in the ninth column, CRUD (for Create, Replace, Update and Delete) and have a recognized value. Currently the two recognized values are CRUD=delete and CRUD=delete-all.

IMPORTANT NOTE: Using the delete operations have the potential of creating orphan features (eg, exons whose gene has been deleted). You should be careful to make sure that doesn't happen. 😭 Included in this distribution is a PostgreSQL trigger (written in plpgsql) that will delete all orphan children recursively, so if a gene is deleted, all transcripts, exons and polypeptides that belong to that gene will be deleted too. See the file modules/sequence/functions/delete-trigger.plpgsql for more information.

delete

The delete option will delete one and only one feature for which the name, type and organism match what is in the GFF line with what is in the database. Note that feature.uniquename are not considered, nor are the coordinates presented in the GFF file. This is so that updates via GFF can be done on the coordinants. 👑

If there is more than one feature for which the name, type and organism match, the loader will print an error message and stop. If there are no features that match the name, type and organism, the loader will print a warning message and continue.

delete-all

The delete-all option works similarly to the delete option, except that it will delete all features that match the name, type and organism in the GFF line (as opposed to allowing only one feature to be deleted). If there are no features that match, the loader will print a warning message and continue.

The run lock

The bulk loader is not a multiuser application. If two separate bulk load processes try to load data into the database at the same time, at least one and possibly all loads will fail. 😵 To keep this from happening, the bulk loader places a lock in the database to prevent other gmod_bulk_load_gff3.pl processes from running at the same time. 🔒

When the application exits normally, this lock will be removed, but if it crashes for some reason, the lock will not be removed. To remove the lock from the command line, provide the flag --remove_lock. Note that if the loader crashed necessitating the removal of the lock, you also may need to rebuild the uniquename cache (see the next section). 🫶

The uniquename cache

The loader uses the chado database to create a table that caches feature_ids, uniquenames, type_ids, and organism_ids of the features that exist in the database at the time the load starts and the features that will be added when the load is complete. 🧠

If it is possilbe that new features have been added via some method that is not this loader (eg, Apollo edits or loads with XORT) or if a previous load using this loader was aborted, then you should supply the --recreate_cache option to make sure the cache is fresh. ✨

Sequence

By default, if there is sequence in the GFF file, it will be loaded into the residues column in the feature table row that corresponds to that feature. By supplying the --nosequence option, the sequence will be skipped. You might want to do this if you have very large sequences, which can be difficult to load. In this context, "very large" means more than 200MB. 🐘

Also note that for sequences to load properly, the GFF file must have the ##FASTA directive (it is required for proper parsing by Bio::FeatureIO), and the ID of the feature must be exactly the same as the name of the sequence following the > in the fasta section.

The ORGANISM table

This script assumes that the organism table is populated with information about your organism. If you are unsure if that is the case, you can execute this command from the psql command-line:

select * from organism;

If you do not see your organism listed, execute this command to insert it:

insert into organism (abbreviation, genus, species, common_name)
              values ('H.sapiens', 'Homo','sapiens','Human');

substituting in the appropriate values for your organism. 🧬

Parents/children order

Parents must come before children in the GFF file. That's it, that's the rule. 💅

Analysis

If you are loading analysis results (ie, blat results, gene predictions), you should specify the -a flag. If no arguments are supplied with the -a, then the loader will assume that the results belong to an analysis set with a name that is the concatenation of the source (column 2) and the method (column 3) with an underscore in between. Otherwise, the argument provided with -a will be taken as the name of the analysis set. Either way, the analysis set must already be in the analysis table. 🚨

The easist way to do this is to insert it directly in the psql shell:

INSERT INTO analysis (name, program, programversion)
             VALUES  ('genscan 2005-2-28','genscan','5.4');

There are other columns in the analysis table that are optional; see the schema documentation and \d analysis in psql for more information.

Chado has four possible columns for storing the score in the GFF score column; please use whichever is most appropriate and identifiy it with --score_col flag (significance is the default). Note that the name of the column can be shortened to one letter. ✨ If you have more than one score associated with each feature, you can put the other scores in the ninth column as a tag=value pair, like 'identity=99', and the bulk loader will put it in the featureprop table (provided there is a cvterm for identity; see the section above concerning custom tags). Available options are:

  • significance (default)
  • identity
  • normscore
  • rawscore

okay so here's the tea on something that is not here yet but is very much planned 💭 — a planned addtion to the functionality of handling analysis results is to allow "mixed" GFF files, where some lines are analysis results and some are literally not. and on top of that, one will be able to supply lists of types (optionally with sources) and their associated entry in the analysis table. the format will probably be tag value pairs, like this:

--analysis match:Rice_est=rice_est_blast, \
           match:Maize_cDNA=maize_cdna_blast, \
           mRNA=genscan_prediction,exon=genscan_prediction
bestie this is a planned addition, so don't come for me if it's not in your version yet 🙈 it's giving roadmap energy.

💞 Grouping features by ID 💞

so the GFF3 specification says you can group features like CDSes and match_parts together by letting them share the same ID. cute concept, love that for them! but this loader does not support that method of grouping. nope. not happening. 🙅

instead, the parent feature must be explicitly created before the parts, and then the parts have to refer back to the parent using the Parent tag. that's the only way this loader is gonna know who belongs to who. periodt. 💅

🚨 External Parent IDs 🚨

so the GFF3 specification states that IDs are only valid within a single GFF file, which means you can't have Parent tags that refer to IDs in another file. those are just the rules, babe. ✨

but — by specificifying the allow_external_parent flag, you can relax this restriction. 👀 and now, the part where I hold your hand, because there is always a warning part:

a word of warning however: if the parent feature's uniquename/ID was modified during loading (to make it unique), this functionality won't work, as it won't be able to find the original feature correctly. 😭 and actually? it may be worse than not working — it may attach child features to the wrong parent. that is literally why it is a bad idea to use this functionality! please use with caution. fr fr, I am begging. 🙏

👑 AUTHORS 👑

Allen Day <allenday@ucla.edu>, Scott Cain <scain@cpan.org>

Copyright (c) 2011

This library is free software; you can redistribute it and/or modify it under the same terms as Perl itself.

and that's the whole page, bestie 🫶 this lil' masterpiece came to you courtesy of perl v5.30.0, dated 2019-12-05, and it's the GMOD_B...AD_GFF3(1p) manpage 💖 go forth and load those features, you absolute icon 👑🌈🦋

📦 where this page comes from