silico.biotoul.fr
 

M1 BBS Graphes - Projets

From silico.biotoul.fr

(Difference between revisions)
Jump to: navigation, search
m (Gene Ontology)
m (Documentation pour GeneOntology.py)
Line 44: Line 44:
= Documentation pour <tt>GeneOntology.py</tt> =
= Documentation pour <tt>GeneOntology.py</tt> =
-
<big>'''<tt>loadOBO(filename)</tt>'''</big>
+
<big>'''<tt>load_OBO(filename)</tt>'''</big>
<source lang='python'>
<source lang='python'>
-
def loadOBO(filename):
+
def load_OBO(filename):
# parse OBO file to create a DAG
# parse OBO file to create a DAG
# obsolete terms are discarded
# obsolete terms are discarded
# only is_a and part_of relationships are loaded
# only is_a and part_of relationships are loaded
 +
# [Term]
 +
# id: GO:0000028
 +
# name: ribosomal small subunit assembly
 +
# namespace: biological_process
 +
# def: "The aggregation, arrangement and bonding together of constituent RNAs and proteins to form the small ribosomal subunit." [GOC:jl]
 +
# subset: gosubset_prok
 +
# synonym: "30S ribosomal subunit assembly" NARROW [GOC:mah]
 +
# synonym: "40S ribosomal subunit assembly" NARROW [GOC:mah]
 +
# is_a: GO:0022618 ! ribonucleoprotein complex assembly
 +
# relationship: part_of GO:0042255 ! ribosome assembly
 +
# relationship: part_of GO:0042274 ! ribosomal small subunit biogenesis
def parseTerm(lines):
def parseTerm(lines):
# search for obsolete
# search for obsolete
for l in lines:
for l in lines:
-
if l.startswith('is_obsolete: true'):
+
if l.startswith('is_obsolete: true'): return
-
return
+
# otherwise create node
# otherwise create node
id = lines.pop(0)[4:].rstrip()
id = lines.pop(0)[4:].rstrip()
Line 62: Line 72:
for line in lines:
for line in lines:
# attributes (name, namespace, def)
# attributes (name, namespace, def)
-
if line.startswith('name:'): term['name'] = line[6:]
+
if line.startswith('name: '): term['name'] = line[6:]
-
elif line.startswith('namespace:'): term['namespace'] = line[11:]
+
elif line.startswith('namespace: '): term['namespace'] = line[11:]
-
elif line.startswith('def:'): term['def'] = line[5:]
+
elif line.startswith('def: '): term['def'] = line[5:]
 +
elif line.startswith('alt_id: '): g['alt_id'][ line[8:] ] = id # alternate ids
# relationships
# relationships
elif line.startswith('is_a:'): # is_a
elif line.startswith('is_a:'): # is_a
parent = line[6:line.index('!')].rstrip()
parent = line[6:line.index('!')].rstrip()
-
e = gr.add_egde(g,id, parent)
+
e = gr.add_edge(g,id, parent)
e['type'] = 'is_a'
e['type'] = 'is_a'
elif line.startswith('relationship: part_of '): # part_of
elif line.startswith('relationship: part_of '): # part_of
line = line[line.index('GO:'):]
line = line[line.index('GO:'):]
dest = line[:line.index(' ')]
dest = line[:line.index(' ')]
-
e = gr.add_egde(g,id, dest)
+
e = gr.add_edge(g,id, dest)
e['type'] = 'part_of'
e['type'] = 'part_of'
#
#
-
g=gr.createGraph(directed=True, weighted=False)
+
g=gr.create_graph(directed=True, weighted=False)
 +
g['alt_id'] = {} # alternate GO ids
with open(filename) as f:  
with open(filename) as f:  
line = f.readline().rstrip()
line = f.readline().rstrip()
Line 102: Line 114:
line = f.readline()
line = f.readline()
return g
return g
-
 
</source>
</source>
Line 136: Line 147:
</source>
</source>
-
<big>'''<tt>loadGOA</tt>'''</big>
+
<big>'''<tt>load_GOA</tt>'''</big>
<source lang='python'>
<source lang='python'>
-
def loadGOA(go, filename):
+
def load_GOA(go, filename):
# !gaf-version: 2.1
# !gaf-version: 2.1
# !GO-version: http://purl.obolibrary.org/obo/go/releases/2016-10-29/go.owl
# !GO-version: http://purl.obolibrary.org/obo/go/releases/2016-10-29/go.owl
Line 146: Line 157:
# UniProtKB  P00393  ndh  NOT  GO:0005737  PMID:6784762    IDA              C  NADH dehydrogenase              DHNA_ECOLI|ndh|JW1095|b1109  protein taxon:83333  20100621  EcoliWiki
# UniProtKB  P00393  ndh  NOT  GO:0005737  PMID:6784762    IDA              C  NADH dehydrogenase              DHNA_ECOLI|ndh|JW1095|b1109  protein taxon:83333  20100621  EcoliWiki
#    0        1      2  3      4            5          6        7      8            9                              10
#    0        1      2  3      4            5          6        7      8            9                              10
-
#             id    name        go_id              evidence-codes                    desc                          aliases
+
#             id    name        go_id              evidence-codes                    desc                          aliases
names = {}
names = {}
go['names'] = names
go['names'] = names
Line 155: Line 166:
cols = line.rstrip().split('\t')
cols = line.rstrip().split('\t')
id = cols[1]
id = cols[1]
-
if id not in go['nodes']:
 
-
g = gr.add_node(go,id)
 
-
g['id'] = id
 
-
g['type'] = 'GeneProduct'
 
-
names[cols[2]] = id
 
-
gp = go['nodes'][id]
 
-
gp['name'] = cols[2]
 
-
gp['desc'] = cols[9]
 
-
gp['aliases'] = cols[10]
 
go_id = cols[4]
go_id = cols[4]
-
if go_id not in go['nodes']:
+
if go_id not in go['nodes']: # GOTerm not found search alternate ids
-
#~ print('Attaching a gene product to non existing GO Term (%s)' % (go_id))
+
if go_id in go['alt_id']: # success
-
gr.add_node(go,go_id)
+
go_id = go['alt_id'][go_id] # replace term
-
go_term = go['nodes'][go_id]
+
else: # warn user
-
e = gr.add_egde(go, id, go_id)
+
print('Warning: could not attach a gene product (%s) to a non existing GO Term (%s)' % (id, go_id))
-
e['type'] = 'annotation'
+
if go_id in go['nodes']:
-
if 'evidence-codes' not in e: e['evidence-codes'] = []
+
# create node for gene product if not already present
-
e['evidence-codes'].append( cols[6] )
+
if id not in go['nodes']:
 +
g = gr.add_node(go,id)
 +
g['id'] = id
 +
g['type'] = 'GeneProduct'
 +
names[cols[2]] = id
 +
# create or update gene product attributes
 +
gp = go['nodes'][id]
 +
gp['name'] = cols[2]
 +
gp['desc'] = cols[9]
 +
gp['aliases'] = cols[10]
 +
# attach gene product to GOTerm
 +
go_term = go['nodes'][go_id]
 +
e = gr.add_edge(go, id, go_id)
 +
e['type'] = 'annotation'
 +
if 'evidence-codes' not in e: e['evidence-codes'] = []
 +
e['evidence-codes'].append( cols[6] )
 +
else: # go_id or alt_id not found in GOTerms
 +
print('Error: could not attach a gene product (%s) to non existing GO Term (%s)' % (id, go_id))
line = f.readline()
line = f.readline()
</source>
</source>

Revision as of 09:02, 5 December 2017

Bibliothèque Python

Une partie du projet consiste à terminer la bibliothèque python entamée au cours des TP.

La liste des méthodes à implémenter est la suivante :

  • TP1: dfs, is_acyclic, topological_sort
  • TP2: bfs, BellmanFord, FloydWarshall, FloydWarshallPath, diameter

Une attention particulière sera portée à la qualité du code et de ses commentaires. Un script de tests/validations devra être fourni (ou bien intégré directement dans la bibliothèque).

Gene Ontology

La deuxième partie du projet consiste à étendre la bibliothèque python afin de fournir des utilitaires pour la Gene Ontology. Ses principales fonctionnalités seront :

  • le chargement du graphe représentant la Gene Ontology
  • le chargement des associations gene product - GO Term

Une fois ces étapes réalisées, les méthodes à implémenter sont :

  • détermination du plus long chemin afin d'obtenir la profondeur maximale des trois ontologies (biological process, molecular function et cellular component)
  • obtention des gene products directement associés à un GO Term
  • obtention des GO Terms directement associés à un gene product
  • obtention des gene products associés à un GO Term ou à un de ses descendants
  • pour un gene product, l'obtention des GO Terms associés en incluant les termes ancêtres

Pour le chargement du graphe de la Gene Ontology ainsi que le chargement des annotations spécifiques à un organisme, il vous est fourni en exemple un fichier GeneOntology.py proposant les fonctions load_OBO pour le chargement du graphe de la Gene Ontology et load_GOA pour le chargement des annotations et leur ajout au graphe précédent. Vous pouvez vous servir de GeneOntology.py comme base de départ pour cette partie du projet, auquel cas vous ajouterez donc des fonctions.

Pour le chargement du graphe, il s'agit d'implémenter la méthode load_OBO qui charge le format OBO (version 1.2, cf. OBO_1.2). Toutes les spécifications du format ne sont pas à respecter. On ne s'intéressera qu'aux tags : id, name, namespace, def, is_a, relationship. Les termes indiqués obsolete sont à ignorer. Le fichier à charger sera celui nommé go-basic dans la partie téléchargement de http://geneontology.org.

Pour le chargement des annotations, le format à prendre en charge s'appelle gaf (version 2.1). Vous pourrez travailler sur le(s) fichier(s) de votre choix à sélectionner dans la partie Annotations de la section téléchargement de http://geneontology.org. Un ensemble de proteomes est disponible sur le FTP de l'EBI : ftp://ftp.ebi.ac.uk/pub/databases/GO/goa, notamment le répertoire proteomes, avec le README qui vous indique à quoi correspondent les colonnes des fichiers. Celles (au minimum) à conserver sont :

  • 2) DB_Object_ID = identifiant unique (id)
  • 3) DB_Object_Symbol = nom du produit du gène (name)
  • 7) Evidence Code = qualité de l'annotation (evidence-code)
  • 11) DB_Object_Synonym = autres identifiants pour ce produit de gène (aliases).

Dossier à rendre sous la forme d'une archive au format zip ou tar.gz

Le code de la bibliothèque annoté et assorti de jeux de tests (sans le fichier obo qui est le même pour tous).

Un rapport synthétique au format PDF comprenant a minima les sections :

  • analyse : analyse du contexte, des besoins et des fonctionnalités à fournir
  • conception : choix d'une représentation et algorithmes (avec leur complexité) pour réaliser les fonctionnalités
  • réalisation : choix techniques
  • bilan et perspectives : est-ce que les besoins sont satisfaits ? qu'est-ce que l'on pourrait améliorer

Pour ces différentes sections, des schémas peuvent venir appuyer votre discours.

Documentation pour GeneOntology.py

load_OBO(filename)

def load_OBO(filename):
	# parse OBO file to create a DAG
	# obsolete terms are discarded
	# only is_a and part_of relationships are loaded
# [Term]
# id: GO:0000028
# name: ribosomal small subunit assembly
# namespace: biological_process
# def: "The aggregation, arrangement and bonding together of constituent RNAs and proteins to form the small ribosomal subunit." [GOC:jl]
# subset: gosubset_prok
# synonym: "30S ribosomal subunit assembly" NARROW [GOC:mah]
# synonym: "40S ribosomal subunit assembly" NARROW [GOC:mah]
# is_a: GO:0022618 ! ribonucleoprotein complex assembly
# relationship: part_of GO:0042255 ! ribosome assembly
# relationship: part_of GO:0042274 ! ribosomal small subunit biogenesis
	def parseTerm(lines):
		# search for obsolete
		for l in lines:
			if l.startswith('is_obsolete: true'): return
		# otherwise create node
		id = lines.pop(0)[4:].rstrip()
		term = gr.add_node(g,id)
		term['id'] = id
		term['type'] = 'GOTerm'
		for line in lines:
			# attributes (name, namespace, def)
			if line.startswith('name: '): term['name'] = line[6:]
			elif line.startswith('namespace: '): term['namespace'] = line[11:]
			elif line.startswith('def: '): term['def'] = line[5:]
			elif line.startswith('alt_id: '): g['alt_id'][ line[8:] ] = id # alternate ids
			# relationships
			elif line.startswith('is_a:'): # is_a
				parent = line[6:line.index('!')].rstrip()
				e = gr.add_edge(g,id, parent)
				e['type'] = 'is_a'
			elif line.startswith('relationship: part_of '): # part_of
				line = line[line.index('GO:'):]
				dest = line[:line.index(' ')]
				e = gr.add_edge(g,id, dest)
				e['type'] = 'part_of'
	#
	g=gr.create_graph(directed=True, weighted=False)
	g['alt_id'] = {} # alternate GO ids
	with open(filename) as f: 
		line = f.readline().rstrip()
		# skip header to reach 1st Term
		while not line.startswith('[Term]'): 
			line = f.readline().rstrip()
		buff = []
		line = f.readline()
		stop = False
		while line and not stop:
			# buffer lines until the next Term is found
			line = line.rstrip()
			# new Term
			if line.startswith('[Term]'):
				# next Term found: create corresponding node and edges in parseTerm and empty buffer
				parseTerm(buff)
				buff=[]
			# last Term
			elif line.startswith('[Typedef]'):
				parseTerm(buff)
				stop=True
			# or append to buffer
			else:
				buff.append(line)
			line = f.readline()
	return g

Exemple pour le GOTerm suivant du fichier .obo

[Term]
id: GO:0000028
name: ribosomal small subunit assembly
namespace: biological_process
def: "The aggregation, arrangement and bonding together of constituent RNAs and proteins to form the small ribosomal subunit." [GOC:jl]
subset: gosubset_prok
synonym: "30S ribosomal subunit assembly" NARROW [GOC:mah]
synonym: "40S ribosomal subunit assembly" NARROW [GOC:mah]
is_a: GO:0022618 ! ribonucleoprotein complex assembly
relationship: part_of GO:0042255 ! ribosome assembly
relationship: part_of GO:0042274 ! ribosomal small subunit biogenesis

Dans le script python, un sommet correspondant sera créé :

 go['nodes']['GO:0000028'] = {
   'id'        : 'GO:0000028',
   'name'      : 'ribosomal small subunit assembly',
   'def'       : 'The aggregation, arrangement and bonding together of constituent RNAs and proteins to form the small ribosomal subunit." [GOC:jl]',
   'namespace' : 'biological_process',
   'type'      : 'GOTerm'
 }

Et les arcs suivants :

 go['edges']['GO:0000028']['GO:0022618'] = { 'type': 'is_a' }
 go['edges']['GO:0000028']['GO:0042255'] = { 'type': 'part_of' }
 go['edges']['GO:0000028']['GO:0042274'] = { 'type': 'part_of' }

load_GOA

def load_GOA(go, filename):
# !gaf-version: 2.1
# !GO-version: http://purl.obolibrary.org/obo/go/releases/2016-10-29/go.owl
# UniProtKB  A5A605  ykfM      GO:0006974  PMID:20128927   IMP              P  Uncharacterized protein YkfM    YKFM_ECOLI|ykfM|b4586         protein taxon:83333  20100901  EcoCyc
# UniProtKB  A5A605  ykfM      GO:0016020  GO_REF:0000037  IEA              C  Uncharacterized protein YkfM    YKFM_ECOLI|ykfM|b4586         protein taxon:83333  20161029  UniProt
# UniProtKB  P00448  sodA      GO:0004784  GO_REF:0000003  IEA  EC:1.15.1.1 F  Superoxide dismutase [Mn]       SODM_ECOLI|sodA|JW3879|b3908  protein taxon:83333  20161029  UniProt
# UniProtKB  P00393  ndh  NOT  GO:0005737  PMID:6784762    IDA              C  NADH dehydrogenase              DHNA_ECOLI|ndh|JW1095|b1109   protein taxon:83333  20100621  EcoliWiki
#     0        1       2   3       4             5          6        7      8             9                              10
#              id    name        go_id               evidence-codes                     desc                           aliases
	names = {}
	go['names'] = names
	with open(filename) as f: 
		line = f.readline()
		while line:
			if not line.startswith('!'):
				cols = line.rstrip().split('\t')
				id = cols[1]
				go_id = cols[4]
				if go_id not in go['nodes']: # GOTerm not found search alternate ids
					if go_id in go['alt_id']: # success
						go_id = go['alt_id'][go_id] # replace term
					else: # warn user
						print('Warning: could not attach a gene product (%s) to a non existing GO Term (%s)' % (id, go_id))
				if go_id in go['nodes']:
					# create node for gene product if not already present
					if id not in go['nodes']:
						g = gr.add_node(go,id)
						g['id'] = id
						g['type'] = 'GeneProduct'
						names[cols[2]] = id
					# create or update gene product attributes
					gp = go['nodes'][id]
					gp['name'] = cols[2]
					gp['desc'] = cols[9]
					gp['aliases'] = cols[10]
					# attach gene product to GOTerm
					go_term = go['nodes'][go_id]
					e = gr.add_edge(go, id, go_id)
					e['type'] = 'annotation'
					if 'evidence-codes' not in e: e['evidence-codes'] = []
					e['evidence-codes'].append( cols[6] )
				else: # go_id or alt_id not found in GOTerms
					print('Error: could not attach a gene product (%s) to non existing GO Term (%s)' % (id, go_id))
			line = f.readline()

Exemple pour les lignes suivantes du fichier .goa d'Escherichia coli K12

UniProtKB       A5A605  ykfM            GO:0006974      PMID:20128927   IMP                                     P       Uncharacterized protein YkfM    YKFM_ECOLI|ykfM|b4586   protein taxon:83333     20100901        EcoCyc          
UniProtKB       A5A605  ykfM            GO:0016020      GO_REF:0000037  IEA     UniProtKB-KW:KW-0472            C       Uncharacterized protein YkfM    YKFM_ECOLI|ykfM|b4586   protein taxon:83333     20161029        UniProt         
UniProtKB       A5A605  ykfM            GO:0016020      GO_REF:0000039  IEA     UniProtKB-SubCell:SL-0162       C       Uncharacterized protein YkfM    YKFM_ECOLI|ykfM|b4586   protein taxon:83333     20161029        UniProt         
UniProtKB       A5A605  ykfM            GO:0016021      GO_REF:0000037  IEA     UniProtKB-KW:KW-0812            C       Uncharacterized protein YkfM    YKFM_ECOLI|ykfM|b4586   protein taxon:83333     20161029        UniProt         
UniProtKB       A5A605  ykfM            GO:0046677      PMID:20128927   IMP                                     P       Uncharacterized protein YkfM    YKFM_ECOLI|ykfM|b4586   protein taxon:83333     20100901        EcoCyc          

Le sommet suivant sera créé

 go['nodes']['A5A605'] = {
   'id'      : 'A5A605',
   'name'    : 'ykfM',
   'desc'    : 'Uncharacterized protein YkfM',
   'aliases' : 'YKFM_ECOLI|ykfM|b4586',
   'type'    : 'GeneProduct'
 }

et les arcs suivants :

 go['nodes']['A5A605'] = {
   'GO:0006974': { 'evidence-codes': ['IMP'],        'type': 'annotation'},
   'GO:0016020': { 'evidence-codes': ['IEA', 'IEA'], 'type': 'annotation'},
   'GO:0016021': { 'evidence-codes': ['IEA'],        'type': 'annotation'},
   'GO:0046677': { 'evidence-codes': ['IMP'],        'type': 'annotation'}}