BioSQL is a joint effort between the OBF projects (BioPerl, BioJava etc) to support a shared database schema for storing sequence data. In theory, you could load a GenBank file into the database with BioPerl, then using Biopython extract this from the database as a record object with featues - and get more or less the same thing as if you had loaded the GenBank file directly as a SeqRecord using SeqIO.
We have some existing documentation (HTML, PDF) for the Biopython interfaces to BioSQL, covering installing Python database adaptors and basic usage of BioSQL. This is a little old, and I am hoping to use this wiki page to update the above documentation in future.
NOTE - At the time of writing, there are a few problems with BioSQL and Biopython 1.44 which are being fixed in CVS.
This is fairly complicated - partly because there are so many options. For example, you can use a range of different SQL database packages (we'll focus on MySQL), you can have the database on your own computer (the assumption here) or on a separate server, and of course there are usernames and passwords associated with database. And finally the details will also vary depending on your operating system.
This text is based in part on the BioSQL scheme INSTALL instructions, which also covers alternatives to MySQL.
Installing Required Software
You will need to install some database software plus the associated python library so that Biopython can "talk" to the database. In this example we'll talk about the most common choice, MySQL. How you do this will also depend on your operating system, for example on a Debian or Ubuntu Linux machine try this:
sudo apt-get install mysql-common mysql-server python-mysqldb
It will also be important to have perl (to run some of the setup scripts) and cvs (to get some BioSQL files). Again, on a Debian or Ubuntu Linux machine try this:
sudo apt-get install perl cvs
You may find perl is already installed.
Downloading the BioSQL Schema & Scripts
One the software is installed, your next task is to setup a database and import the BioSQL scheme (i.e. setup the relevant tables within the database). If you have CVS installed, then on Linux you can download the latest schema like this (password is 'cvs').
cd ~ mkdir repository cd repository cvs -d :pserver:firstname.lastname@example.org:/home/repository/biopython checkout biosql cvs -d :pserver:email@example.com:/home/repository/biosql checkout biosql-schema cd biosql-schema/sql
If you don't want to use CVS, then download the files via the View CVS web interface. Click the Download tarball link to get a tar.gz file containing all the current CVS file, and then unzip that.
Or, navigate to the relevant schema file for your database and download just that, e.g. biosql-schema/sql/biosqldb-mysql.sql for MySQL. You will also want the NCBI Taxonomy loading perl script, scripts/load_ncbi_taxonomy.pl.
Creating the empty database
Assuming you are using MySQL, the following command line should create a new database on your own computer called bioseqdb, belonging to the root user account:
mysqladmin -u root create bioseqdb
We can then tell MySQL to load the BioSQL scheme we downloaded above (file ~/repository/biosql-schema/sql/biosqldb-mysql.sql if you downloaded the files where we suggested).
cd ~/repository/biosql-schema/sql mysql -u root bioseqdb < biosqldb-mysql.sql
You can have a quick play using the mysql command line tool, for example:
mysql --user=root bioseqdb -e "show tables"
+----------------------------+ | Tables_in_bioseqdb | +----------------------------+ | biodatabase | | bioentry | | bioentry_dbxref | | bioentry_path | | bioentry_qualifier_value | | bioentry_reference | | bioentry_relationship | | biosequence | | comment | | dbxref | | dbxref_qualifier_value | | location | | location_qualifier_value | | ontology | | reference | | seqfeature | | seqfeature_dbxref | | seqfeature_path | | seqfeature_qualifier_value | | seqfeature_relationship | | taxon | | taxon_name | | term | | term_dbxref | | term_path | | term_relationship | | term_relationship_term | | term_synonym | +----------------------------+
Or, to look inside a table:
mysql --user=root bioseqdb -e "select * from bioentry;"
This should return no rows as the table is empty.
Before you start trying to load sequences into the database, it is a good idea to load the NCBI taxonomy database using the scripts/load_ncbi_taxonomy.pl script in the BioSQL package.
The script should be able to download the files it needs from the NCBI taxonomy FTP site automatically:
cd ~/repository/biosql-schema/scripts ./load_ncbi_taxonomy.pl --dbname bioseqdb --driver mysql --dbuser root --download true
There is about 10MB to fetch, so it can take a little while (and doesn't give any feedback while this happens). If you are worried, open a file browser window and check to see it is downloading a file called taxdump.tar.gz to the taxdata subdirectory.
You should see this output at the command prompt - be warned that some of these steps do take a while (especially rebuilding nested set left/right values):
Loading NCBI taxon database in taxdata: ... retrieving all taxon nodes in the database ... reading in taxon nodes from nodes.dmp ... insert / update / delete taxon nodes ... (committing nodes) ... rebuilding nested set left/right values ... reading in taxon names from names.dmp ... deleting old taxon names ... inserting new taxon names ... cleaning up Done.
This might be a good point for a tea break - I didn't time this but it was over ten minutes.
Running the unit test
Because there are so many ways you could have setup your BioSQL database, you my have to tell the unit test a few bits of information. If you have followed these instructions, then the unit test should just work (using CVS, what will be the next release after Biopython 1.44).
If you have done things differently (e.g. PostgreSQL instead of MySQL, or using a different database username and password), then you need to edit the script Tests/setup_BioSQL.py and fill in the following fields:
DBDRIVER = 'MySQLdb' DBTYPE = 'mysql'
and a little lower down,
DBHOST = 'localhost' DBUSER = 'root' DBPASSWD = '' TESTDB = 'biosql_test'
Change these to match your setup. You can then run the unit test as normal, e.g.
python run_tests.py test_BioSQL
Creating a (sub) database
BioSQL lets us define named "sub" databases within the single SQL database (which we called bioseqdb earlier). For this example, lets create a one for some orchid sequences:
from BioSQL import BioSeqDatabase server = BioSeqDatabase.open_database(driver="MySQLdb", user="root", passwd = "", host = "localhost", db="bioseqdb") db = server.new_database("orchids", description="Just for testing") server.adaptor.commit()
The call server.adaptor.commit() tells the database to save the changes so far (commit the SQL transaction). It is up to you to decide when to commit the SQL transaction(s), and/or rollback changes, rather than having Biopython try and decide for you and risk getting it wrong. See Explicit is better than implicit (The Zen of Python) and bug 2395.
There should now be a single row in the biodatabase table for our new orchid database. You can check this at the command line:
mysql --user=root bioseqdb -e "select * from biodatabase;"
Which should give something like this (assuming you haven't done any other testing yet):
+----------------+---------+-----------+------------------+ | biodatabase_id | name | authority | description | +----------------+---------+-----------+------------------+ | 1 | orchids | NULL | Just for testing | +----------------+---------+-----------+------------------+
Now that we have setup an orchids database within our biosqldb MySQL database, lets add some sequences to it.
Loading Sequences into a database
When loading sequences into a BioSQL database with Biopython we have to provide annotated SeqRecord objects. This gives us another excuse to use the SeqIO module! A quick recap on reading in sequences as SeqRecords, based on one of the orchid examples in the Biopython Tutorial:
from Bio import GenBank from Bio import SeqIO handle = GenBank.download_many(['6273291', '6273290', '6273289']) for seq_record in SeqIO.parse(handle, "genbank") : print seq_record.id, seq_record.description[:50] + "..." print "Sequence length %i," % len(seq_record.seq), print "from: %s" % seq_record.annotations['source'] handle.close()
The expected output is below, note we have three records with a total of nine features:
AF191665.1 Opuntia marenae rpl16 gene; chloroplast gene for c... Sequence length 902, 3 features, from: chloroplast Opuntia marenae AF191664.1 Opuntia clavata rpl16 gene; chloroplast gene for c... Sequence length 899, 3 features, from: chloroplast Grusonia clavata AF191663.1 Opuntia bradtiana rpl16 gene; chloroplast gene for... Sequence length 899, 3 features, from: chloroplast Opuntia bradtianaa
Now, instead of printing things on screen, let's add these three records to a new (empty) orchid database:
from Bio import GenBank from Bio import SeqIO from BioSQL import BioSeqDatabase server = BioSeqDatabase.open_database(driver="MySQLdb", user="root", passwd = "", host = "localhost", db="bioseqdb") db = server["orchids"] handle = GenBank.download_many(['6273291', '6273290', '6273289']) db.load(SeqIO.parse(handle, "genbank")) server.adaptor.commit()
Again, you must explicitly call server.adaptor.commit() to record the SQL transaction which is otherwise left pending.
The db.load() function should have returned the number of records loaded (three in this example), and again have a look in the database and you should see new rows in several tables.
The bioentry and biosequence tables should have three new rows:
mysql --user=root bioseqdb -e "select * from bioentry;" mysql --user=root bioseqdb -e "select * from biosequence;"
The should also be nine new features:
mysql --user=root bioseqdb -e "select * from seqfeature;"
Next, we'll try and load these three records back from the database.
Extracting Sequences from the database
This continues from the previous example, where we loaded three records into an orchids database:
from BioSQL import BioSeqDatabase server = BioSeqDatabase.open_database(driver="MySQLdb", user="root", passwd = "", host = "localhost", db="bioseqdb") db = server["orchids"] for identifiers in ['6273291', '6273290', '6273289'] : seq_record = db.lookup(gi=identifiers) print seq_record.id, seq_record.description[:50] + "..." print "Sequence length %i," % len(seq_record.seq)
AF191665.1 Opuntia marenae rpl16 gene; chloroplast gene for c... Sequence length 902 AF191664.1 Opuntia clavata rpl16 gene; chloroplast gene for c... Sequence length 899 AF191663.1 Opuntia bradtiana rpl16 gene; chloroplast gene for... Sequence length 899
Todo - sort out the annotation, e.g. bug 2396.
The objects you get back from BioSQL act like a SeqRecord object with a Seq object as the sequence, but they are not exactly the same. You actually get their BioSQL database equivalent, a DBSeqRecord object with a DBSeq object for the sequence. These will only load the sequence and annotation from the database on demand.
Deleting a (sub) database
As mentioned above, BioSQL lets us define named "sub" databases within the single SQL database (which we called bioseqdb). In the previous example, we created a sub-database for some orchid sequences. The following code will delete the orchid database (and all the records in it):
from BioSQL import BioSeqDatabase server = BioSeqDatabase.open_database(driver="MySQLdb", user="root", passwd = "", host = "localhost", db="bioseqdb") server.remove_database("orchids") server.adaptor.commit()
Again, the call server.adaptor.commit() is a work around for bug 2395.
There should now be one less row in the biodatabase table, check this at the command line:
mysql --user=root bioseqdb -e "select * from biodatabase;"
You can also check that the three orchid sequences have gone from the other tables.
MySQL Tips and Tricks
If you are getting timeout errors, check to see if your SQL server has any orphaned threads.
mysql --user=root bioseqdb -e "SHOW INNODB STATUS\G" | grep "thread id"
And if there are, assuming you are the only person using this database, you might try killing them off using the thread id given by the above command:
mysql --user=root bioseqdb -e "KILL 123;"
Use at your own risk!