Old setup
Jump to navigation
Jump to search
vue d'ensemble
Marche a suivre pour la mise en place d'une base de donnée pour l'anaylse comparative d'un groupe de génomes. Les données sont stockées dans une base de donnée MysSQL. Les génomes annotés sont importé à l'aide du shéma biosql (http://biopython.org/wiki/BioSQL). Ce shéma permet de stocker facilement des données biologiques mais n'est pas adapté à des analyses comparatives. J'ai donc ajouté un ensemble de bases de données et tables à la base de donnée par moi-même.
- l'ensemble des tables crées par mes scripts s'appellent _<nom de la base de donnée>.
- il est possible d'explorer ces données via une interface web développée avec django.
bases de données
biosqldb (cf. http://biopython.org/wiki/BioSQL )
- ensemble de tables qui permet de stocker des données biologiques
key tables pour l'utilisation des données comparatives
- biodatabase: contient la liste des bases de données importée. Une base de donnée est un ensemble de bioentries (et dans notre cas précis, un ensemble de genbank record contenant des genomes/plasmides). Chaque database à un identifiant propre (biodatabase_id).
- bioentry: contient la liste des bioentries de chaque biodatabase. Les chromosomes et plasmides d'un même génome ont un même taxon_id reporté dans cette table. Si des plasmides et le chromosome d'un même organisme n'ont pas le même taxon_id, c'est qu'il y a eu un problème lors de l'import des genbank, et donc un problème lié au formatage des genbank (cf script gbk_check.py utilisé pour vérifier le format des genbank).
tables ajoutée au shéma originel
- genomes_info_<biodb name>: contient quelques statistiques sur les génomes importés (longueur, GC,…)
- orthology_detail\_<biodb name>: contient l'ensemble des "protein encoding CDS" de la base de donnée <biodb>, avec leur groupe d'orthologue correspondant et quelques infos générales (TM domains, SP, start, stop, gene name, aa sequence,…).
- interpro\_<biodb>: annotation interproscan de l'ensemble des protéines de la base de donnée
COG
- contient des tables générales pour la description des COGs + une table par <biodb> avec le COG du best hit de chaque protéine.
enzyme
- contient les données PRIAM et GHOSTKOALA (KEGG) avec des tables générales avec l'annotation de tous les EC, KO, pathways, modules et les différentes tables permettant de savoir quels KO/EC font partie de quels modules/pathways kegg.
- contient une table/biodb avec les EC priam et les KO
comparative\_tables
nterpro\_pMK190\_1:- ensemble de tables permettant de faire les analyses comparatives.
- orthology\_<biodb>: matrice avec le count dans chaque <taxon\_id> de tous les groupes d'orthologues.
- Pfam\_<biodb>: matrice avec le count dans chaque <taxon\_id> de tous les domaines Pfam.
- EC\_<biodb>: matrice avec le count dans chaque <taxon\_id> de tous le EC numbers.
- COG\_<biodb>: matrice avec le count dans chaque <taxon\_id> de tous les COGs.
- interpro\_<biodb>: matrice avec le count dans chaque <taxon\_id> de toutes les entrées interpro.
- cogs\_profiles\_euclidian\_distance\_<biodb>: distance euclidienne entre pairs de profile présence/absence de cogs.
- identity\_closest\_homolog\_<biodb>: table avec pour chaque locus\_tag, l'identité de l'homologue le plus proche. Redondant avec l'information contenue dans la base de donnéée orth\_<biodb>, mais beaucoup plus rapide d'accès. Utilisé pour générer les circos.
- shared\_orthogroups\_average\_identity\_<biodb>: identity moyenne et médianne entre toutes les protéines d'un groupe d'orthologue.
- shared\_orthogroups\_<biodb>: nombre d'orthologues partagés entre pairs de taxons.
orth\_<biodb name>
- contient une table/groupe d'orthologue qui est une matrice d'itentité entre toutes les sequences aa. Identité calculée sur la base des alignements.
biosqldb\_phylogenies
- contient une table par base de donnée avec l'ensemble des phylogénies des groupes d'orthologues au format newick
blastnr
plast\_nr2sqltable\_v2.py
- contient les résultats des blast vs nr de chaque accession (chromosomes et plasmides sont séparés, une table/accession)
- structure un peu complexe car un hit peut avoir plusieurs hsp,…
- table blastnr\_taxonomy contient la taxonomie détaillée des hits. Si on importe un nouveau set de résultats blastnr, les taxon pas encore présent dans cette table seront ajoutés (récupérés sur le ncbi avec eutils).
installation
installation dépendances
sudo apt-get install mysql-server
sudo apt-get install memcached
sudo apt-get install python-django
sudo pip install simplejson
sudo apt-get install python-pip python-dev build-essential
sudo pip install dill
sudo apt-get install zlib1g-dev
sudo apt-get install libmemcached-dev
sudo pip install pylibmc
sudo pip install django-pylibmc
sudo pip install django-sorting-bootstrap
sudo pip install gunicorn
sudo pip install 'django-autocomplete-light>=2.0.0a1'
sudo apt-get install python-biopython
sudo apt-get install python-biopython-doc
sudo apt-get install python-biopython-sql
sudo pip install pandas
sudo pip install django-bootstrap3 sudo pip install blah
sudo pip install BeautifulSoup
sudo pip install simplejson
sudo apt-get install python-qt4
sudo pip install -U pip setuptools
sudo pip install matplotlib
sudo apt-get install python-tk
sudo apt-get install libxml2-dev libxslt1-dev
sudo pip install lxml
sudo apt-get install python-tk
sudo pip install matplotlib
sudo pip install ete2
installation base de donnée mysql
- installer mysql
- Créer un utilisateur. Les identifiants de l'utilisateur root de Mysql doivent etre fournis à l'application Django. J'utilise estrella3 comme password root pour accéder à cette base de donnée, si tu met le meme ça nous évitera d'avoir des problèmes de synchronisation de nos fichiers, qu'en dis-tu?
loader la base de donnée minimale
- contient le shema biosql
- contient des tables pour les COG, KO,…
TODO faire les liste de ces tables et comment les loader si necessaire
# sauvegarde # mysqldump -u root -p --all-databases > biodb_sheme_db.sql # relaod mysql -u root -pestrella3 < biodb_sheme_db.sql
installation django
répertoires bitbucket
les cloner et les ajouter au $PATH et au $PYTHONPATH
- utilitaires plein de petits utilitaires dont certains sont nécessaires à la base de donnée
- biosqldb nombreux scripts pour interagir avec la base de donnée mysql
- plots scripts pour déssiner les régions génomiques, circos,…
- chlamdb contient les différents scripts pythons pour l'application django, les templates pour les pages html, les fichiers statiques tels que scripts js et les pages de style css
- orthomcl execution d'orthomcl sur le cluster et parsing du ficher avec les groupes
- annotation_pipeline pour le clustering et l'annotation des protéines (remplace orthomcl depuis début 2016)
installation memcached (mise en memoire de données pour gagner en réactivité)
- Installer memcached
augmenter la taille minimale des objets stockés
- éditer /etc/memcached.conf
# Increase limit -I 128M
- redémarrer: /etc/init.d/memcached restart
configuration
configuration de l'application django
- créer le fichier custom\_tags.py dans templatetags
# adapter le path sudo emacs /usr/lib/python2.7/dist-packages/django/templatetags/custom_tags.py
# le fichier doit contenir: from django import template register = template.Library() @register.filter def keyvalue(dict, key): return dict[key]Configurer l'emplacement de l'aplication et du dossier template
Copier le fichier settings.py de chlamdb\_instructions et le mettre dans chlamdb/settings/settings.py
Editer ensuite celui-ci
TEMPLATE_DIRS = ( "/home/trestan/work/dev/django/chlamydia/templates/", ) PROJECT_DIR = "/home/trestan/work/dev/django/chlamydia/chlamdb/"configuration de l'accès mysql
Ne changer que si tu as un autre mot de passe qu'estrella3
biosqldb/mysqldb_plot_genomic_feature.py ligne 208
conn = MySQLdb.connect(host="localhost", # your host, usually localhost user="root", # your username passwd="estrella3", # your password db="orth_%s" % biodb_name) # name of the data base
biosqldb/manipulate_mysqldb.py ligne 51
t server = BioSeqDatabase.open_database(driver="MySQLdb", user="root", passwd = "estrella3", host = "localhost", db="biosqldb")création des utilisateurs
- les informations concernant les utilisateur sont dans une base de donnée sqlite. Cette base de donnée ne peut pas se synchroniser, il faut donc:
- Initialiser
python manage.py syncdb
- Créer un super utilisateur:
python manage.py createsuperuser
Créer de nouveaux utilisateurs
- lancer le shell django (dans le dossier de l'application django)
python manage.py shell
- loader le module nécessaire et créer un utilisateur
from django.contrib.auth.models import User user = User.objects.create_user('john', 'lennon@thebeatles.com', 'johnpassword') user.save()
lancement d'une aplication django
Se fait à l'aide du script manage.py situé dans le répertoire chlamdb.
python manage.py runserver --nothreading 0.0.0.0:8000
Va démarrer un serveur. Pour accèder à la page d'accueil, aller à l'adresse: http://localhost:8000/chlamdb/login
utilisateur: crib
pw: estrella3
mise en place d'une base de donnée
téléchargement des génomes
# c complete genomes only taxid2genomes.py -t <taxon id> -c
- si bcp de génomes d'une meme espèce, il vaut mieux mettre le taxid d'une souche particulière
formatage gbff files
mv G*/*gbff.gz . gunzip *gbff # à executer avec l'ensemble des fichiers gbff à include dans la base de donnée gbk_check.py -i *gbff
gbk_check.py génomes complets
- écrit un nouveau fichier par record (pour les génomes complets, écrits des fichiers séparés pour le chromosomes et le(s) plasmides) + un ficher avec tous les records associés à un taxid.
- modifie la description pour ne garder que genre-espèce-souche dans la déscription du génome
- supprime le taxid des records: un taxid local sera attribué à chaque souche lors de l'import des gbk dans mysql (chromosome(s) et plasmide(s) auront le même taxid)
gbk_check.py draft genomes
- filtration contigs < 1000
- concaténastion des séquences avec 200\*N entre chaque contig
- si pas annoté: réannotation avec prokka
- ATTENTION: il est nécessaire d'exécuter gbk_check.py -i *gbff avec l'ensemble des fichiers en une seule fois si certains génomes doivent etre réannotés, car en cas de réannotation, il faut attribuer un locus tag à chaque génome (basé sur genre-espèce-souche). Le script inclu un controle qui évite que 2 génomes aient des locus\_tag identiques.
- Il vaut la peine d'ouvrir les gbk est voir si la description nécessite des changement. gbk_check.py enlève quelques info inutile (Chlamydia trachomatis D, nil), mais peut passer à coté de certaines choses. Noms parfois incomplets. Description utilisée par la suite comme label pour chaque génome.
- IL EST IMPORTANT DE CHECKER SUR MYSQL SI CHAQUE GENOME A UN TAXON\_ID UNIQUE
select taxon_id from bioentry where biodatabase_id="140" group by taxon_id having count(*)>1;
extraction des séquences
ATTENTION: bien vérifier que chaque CDS a bien un locus tag, par ex avec le script:
for i in `ls *gbk`;do echo $i; add_locus_tag.py -i $i -c;done
- si problèmes de ce coté là, y remédier. Le script add_locus_tag.py peut remédier à 2 situations. Crée un nouveau fichier avec extension _locus.gbk
- ajout de locus tag pour tous les CDS (définir un prefix qui sera incrémenté)
add_locus_tag.py -i <file>.gbk -l <prefix>
- utilise gene comme locus tag uniquement si aucun locus n'est identifié (cas observé avec le génome de C. felis)
add_locus_tag.py -i <file>.gbk -g
création d'un fichier faa par gbk (inclu seq des chromosome + plasmide(s))
- ignore les loci labellés comme pseudogènes
- crée des fasta avec comme entete ">locus\_tag description": il vaut mieux travailler avec les locus\_tag, car les protein\_id ne sont pas uniques (2 proteins identiques entre 2 génomes peuvent avoir le meme protein\_id)
TODO chhanger pour les protein id, car plusieurs protéines peuvent avoir un même locus tag. Va nécessiter pas mal d'adaptations.
for i in `ls *gbk`; do gbk2faa.py -i -f $i;done
import des génomes dans la base de donnée (schéma biosql)
# va créer une base de donnée si elle n'existe pas load_gbk_into_mysql_db.py -d <database name> -g *merged.gbk
- NOTE: possibilité d'effacer une base de donnée avec ce meme script: load_gbk_into_mysql_db.py -d <database name> -r
annotation des génomes
il faut d'abord faire des cluster de sequences homologues
# sur vital-it (script local en dévelopment) bsub -o out_main.txt -e err_main.txt ~/annotation_pipeline/orthofinder_vital.py -i faa
- crée un ficher avec les clusters
annotation des séquences (aa)
- input: 1 fichier faa par génome (avec seq chr + plasmides), ficher de sortie d'orthofinder (OrthologousGroups.txt)
# -i INPUT_FASTA [INPUT_FASTA ...], --input_fasta INPUT_FASTA [INPUT_FASTA ...] # input fasta file # -m MCL_FILE, --mcl_file MCL_FILE # mcl_file # -s, --interpro execute interproscan # -p, --priam execute PRIAM # -t, --phylo execute RAXML for phylogenetic tree reconstruction of # all orthogroups with >2 sequences # -r, --rerun try to rerun jobs killed by LFS system, increasing # missing ressources # -l, --run_local Run locally # -c, --cogs Blast COGs database (using BlastP) # -n, --blast_refseq Blast ncbi refseq database (using PLASTp) # -w, --blast_swissprot # Blast ncbi nr database (using PLASTp) # -a, --get_align get multi fasta & alignments from orthofinder results # (need fasta nucl or aa as input, with locus tag as # header) # -gc, --get_core get single copy core groups # -gct CORE_TOLERANCE_MISSING, --core_tolerance_missing CORE_TOLERANCE_MISSING # number of missing data accepted to include in core # (default=0, only single copy conserved in all genomes) # -nu, --nucl use nucleotide sequences (alignment only) # soumettre à la queue longue car peut nécessiter plus de 24 heures si certains sous-jobs se font killer pour des raisons de durée ou de ressource (automatiquement resoumis) bsub -q long -o annot_out.txt -e annot_err.txt ~/annotation_pipeline/annotation_pipeline.py -i *faa -m Ortogroups.txt -s -p -c -n -w -r -a -t -gc -gct 0
loading des résultats dans la base de donnée
- interproscan
- cave at créer la database interpro et ses tables correspondantes
cd interpro_results interpro2sqltable.py -i interpro_*/*tsv -d <database_name>
or
cd interpro_results for i in `ls ./interpro_*/*tsv`; do interpro2sqltable.py -i $i -d <database_name>; done
- groupes d'orthologues
# utiliser le meme nom de base de donnée que précédement mysqldb_load_mcl_output.py -d <mysqldb name> -o -m Orthogroups.txt # CAVE THE -O IS IMPORTANT !!!
- COGs
# import des résultats blast dans la base de donnée, génération des tables correspondates # coverage_cutoff=50, identity_cutoff=25 blastCOG2COG_classification.py -i blast_*txt -d <database_name>
- KO (groupes d'orthologues de la base de donnée KEGG)
- on peut les obtenir en soumettant les séquences à leur server (http://www.kegg.jp/ghostkoala/). On peut soumettre toutes les séquences d'un seul coup, en utilisant le fichier all.faa.
priam_genomeEnzyme2sqltable.py -k user_ko_all.txt -d <database_name>
or
priam_genomeEnzyme2sqltable.py -k user_*.txt -d <database_name>
- EC priam (même script)
priam_genomeEnzyme2sqltable.py -i PRIAM_*/PRIAM_*/ANNOTATION/sequenceECs.txt -d <database_name>
- ajout table d'identité de tous les groupes d'orthologues
- calcul identité parwise ajout à la base de donnée orth\_<base de donnee>
cd aa_alignments orthogroup_identity_db.py -a *faa -d <database name>
- ajout des phylogénies de tous les groupes d'orthologues
load_phylogenies_into_mysql_db.py -d <database_name> -t RAxML_result.*
reference phylogeny
- utilisée pour visualiser la distribution des homologues dans le groupe d'espèces étudié
get single copy core groups fasta alignment
- les entetes des séquences vont contenir les taxons ids
# sauve les alignement dans le dossier d'execution mysqldb_load_mcl_output.py -d chlamydia_11_15 -c .
# alignements avec MAFFT for i in `ls *txt`; do mafft --maxiterate 1000 --quiet $i > ${i/.txt/_mafft.faa}; done# conctenation des alignement sur la base des entetes (taxon id du génome) concat_align.py -f group*mafft.faa
phylogenetic reconstruction
# remove spaces in header (TODO: modifier script qui génère les fasta pour enlever cet espace) sed -i 's/ //' msa.fa # phylogenetic reconstruction FastTree -gamma -spr 4 -mlacc 2 -slownni msa.fa > chlamydiales_fasttree.tree
rename nodes (optional, only if you want to replace taxon ids b genome labels for figtree display)
phylo_convert_leaf_labels.py -i core_chlamydia.tree -a -d chlamydia_04_16
ajout phylogénie de référence
- insertion d'une phylégénie de référence (par ex basée dur core genes concaténés, 16S,…)
login mysql
mysql -uroot -pestrell3 biosqldb
commade sql ajout phylogénie
- nécessite de savoir l'identifiant de la base de donnée. Dans la base biosqldb Mysql:
select * from biodatabase;
- relever le database id et l'utiliser dans la commande suivant, avec l'arbre en format newick
insert into reference_phylogeny (biodatabase_id, tree) values(<database id>, "(1035343:0.21057,(301:0.20823,(1172027:0.00064,(1069693:0.00064,1172028:0.00101)0.493:0.00089)1.000:0.17718)0.920:0.03258,(((62:0.00127,(1069694:0.00144,313:0.00159)0.739:0.00109)1.000:0.04219,(1143376:0.00064,330758:0.00061)1.000:0.03251)1.000:0.16344,((308:0.00064,(59:0.00064,(60:0.00064,307:0.00064)0.783:0.00064)0.904:0.00064)1.000:0.23067,(66:0.25811,((414781:0.18184,1137444:0.14988)1.000:0.11694,((1279496:0.48960,(64:0.27784,1279497:0.36339)0.925:0.05753)1.000:0.13235,((1021545:0.72605,((134:0.70304,((133:0.13884,132:0.21859)1.000:0.23610,(291:0.41268,(129:0.26749,130:0.27988)1.000:0.13232)0.999:0.08764)1.000:0.16114)1.000:0.24198,(262606:0.61870,(1265182:0.51047,((128:0.37760,(126:0.12204,(1276848:0.16315,616422:0.11342)0.975:0.05046)1.000:0.22251)1.000:0.26239,((125:0.09276,896740:0.08499)1.000:0.43896,(1279493:0.41975,(127:0.33659,131:0.27939)1.000:0.15519)0.978:0.09479)0.999:0.08432)0.999:0.09551)0.699:0.06018)0.999:0.07929)1.000:0.14003)1.000:0.39573,((48:0.00994,(283:0.02641,55:0.01940)0.444:0.00519)1.000:0.08775,((804807:0.05632,((314:0.03612,886707:0.03610)1.000:0.05981,((45:0.01097,52:0.00608)1.000:0.02279,(46:0.03300,47:0.03460)0.428:0.00903)1.000:0.03179)1.000:0.02932)0.944:0.01454,(49:0.10299,87925:0.08127)1.000:0.03259)0.866:0.03656)1.000:0.38582)0.734:0.07865)1.000:0.12121)0.299:0.02953)1.000:0.04305)0.549:0.02263)0.548:0.01948);")
operon Door2
table obtenue sur la base de donnée door2
>select * from DOOR2_operon_accessions limit 10; +---------+-----------+ | door_id | accession | +---------+-----------+ | 13 | NC_009930 | | 111 | NC_009928 | | 1101 | NC_009931 | | 1190 | NC_009929 | | 1740 | NC_009932 | | 2209 | NC_009926 | | 2477 | NC_009933 | | 2636 | NC_009925 | | 3287 | NC_009927 | | 3396 | NC_009934 |
récuperation de tous les operons disponibles pour les accessions d'une base de donnée
- utilise les anciens locus tag, donc pas toujours complet
- necessite table DOOR2_operon_accessions
- extraite de la page qui liste tous les génomes sur door
- permet de savoir quels sont des génomes disponibles sur door et des télécharger les data disponibles
- souvent basé sur les vieux locus tag et non les locus tags refseq
+---------+-----------+ | door_id | accession | +---------+-----------+ | 13 | NC_009930 | | 111 | NC_009928 | | 1101 | NC_009931 | | 1190 | NC_009929 | | 1740 | NC_009932 | | 2209 | NC_009926 | | 2477 | NC_009933 | | 2636 | NC_009925 | | 3287 | NC_009927 |
operon_door2.py -d chlamydia_04_16
create seqfeature\_id2old\_locus\_tag
- old locus tags from genbank
./create_locus2old_locus_table.py -d chlamydia_04_16
create seqfeature\_id2locus\_tag
- le but est de remplacer tous les identifiants par des identifiants numériques pour gagner en vitesse
- à terme: changer orthology\_detail, pathways,COG,KO,…
get_locus2seqfeature_table.py -d chlamydia_04_16
get uniprot data
- search if protein is on uniprot
- get annotation status, annotation score and annotation
- get go annotation
- get all cross references to other databases
- some refseq sequences are not on uniprot (the ones differing from initial annotation)
accession2db_links.py -d chlamydia_04_16
get go annotation
- from interpro annotation using:
- and from uniprot, which includes annotation derived from unirules, hamap,…
- basic idea: show go from uniprot when available, or interpro/COG when not available
classification: available from go sql shema
- molecular function
- cellular component
- biological process
update table identity\_closest\_homolog2\_chlamydia\_04\_16
- on remplace les locus tags par des seqfeature\_id
- script comparative_tables.py est déjà modifié
add table with locus of interest
columns of the table:
- category
- gene
- locus\_tag
- description
- reference
+-----------------------------------------------------+------+-------------+--------------------------------------+---------------------------------------------------------------------------------+ | category | gene | locus_tag | description | reference | +-----------------------------------------------------+------+-------------+--------------------------------------+---------------------------------------------------------------------------------+ | Components of the Parachlamydia chemosensory system | Tar | PUV_RS13085 | Methyl-accepting chemotaxis protein | Unity in variety--the pan-genome of the Chlamydiae. Collingro et al. MBE (2011) | | Components of the Parachlamydia chemosensory system | Tar | PUV_RS05455 | Methyl-accepting chemotaxis protein | Unity in variety--the pan-genome of the Chlamydiae. Collingro et al. MBE (2011) | +-----------------------------------------------------+------+-------------+--------------------------------------+---------------------------------------------------------------------------------|
- on peut laisser des champs libres, il faut juste le bon nombre de colonnes
load tab delimited file
load_annotation_table.py -t table -d chlamydia_04_16
Générer les tables "comparative\_tables"
comparative_tables.py -d <database_name> #and options
-Charger les fichier interpro \*.html dans assets/
TODO générer la table "genomes"
-Table "identity best hit"
orthofinder_blast2pairwise_BBH.py -s SequenceIDs.txt -b Blast*.txt -d <biodb_name>