source: orange-bioinformatics/_bioinformatics/obiKEGG.py @ 1682:2e2c671e4300

Revision 1682:2e2c671e4300, 41.8 KB checked in by Lan Zagar <lan.zagar@…>, 23 months ago (diff)

Replaced split(\n) with splitlines (closes #1160).

Line 
1"""
2=======
3obiKEGG
4=======
5obiKEGG module is an interface to `Kyoto Encyclopedia of Genes and Genomes
6<http://www.genome.jp/kegg/>`_ that allows easy access to KEGG pathway
7and genes data.
8
9Example ::
10   
11    >>> genes = KEGGOrganism("human")
12    >>> genes
13"""
14
15from __future__ import absolute_import
16
17try:
18    import Image, ImageDraw, ImageMath
19except:
20    pass
21
22import cStringIO
23import math
24import time
25import os, sys, tarfile
26import re
27
28from cPickle import load, loads, dump
29from collections import defaultdict
30
31import xml.dom.minidom as minidom
32
33from functools import partial, wraps
34import urllib2
35import cPickle
36
37from Orange.orng import orngEnviron, orngServerFiles
38
39from . import obiData, obiProb, obiTaxonomy
40
41DEFAULT_DATABASE_PATH = orngServerFiles.localpath("KEGG")
42KEGG_FTP_PATH = "ftp://ftp.genome.jp/pub/kegg/"
43
44_join = lambda *args: os.path.join(DEFAULT_DATABASE_PATH, *args)
45
46def special_tags(domain, filename):
47    return dict([tuple(tag.split(":", 1)) for tag in orngServerFiles.info(domain, filename)["tags"] \
48                 if tag.startswith("#") and ":" in tag])
49def loads(domain, filename, version=None):
50    """ A function decorator for a function that loads data from server files.
51    Makes sure the file is downloaded and of the right version
52    """
53    def decorator(func):
54        @wraps(func)
55        def wrapper(*args, **kwargs):
56            try:
57                ## First try downloading the required server file 
58                file = filename if isinstance(filename, basestring) else filename(*args, **kwargs)
59                orngServerFiles.localpath_download(domain, file)
60                if version is not None and special_tags(domain, file).get("#version", str(version)) < str(version):
61                    orngServerFiles.update(domain, file) 
62            except Exception, ex:
63                print ex, args, kwargs
64            return func(*args, **kwargs)
65        return wrapper
66    return decorator
67           
68def deprecated(first_arg, text=None):
69    import warnings
70    def wrapper(*args, **kwargs):
71        warnings.warn_explicit(("Call to deprecated function %s. " % first_arg.__name__) + (text or ""), 
72                               category=DeprecationWarning,
73                               filename=first_arg.func_code.co_filename,
74                               lineno=first_arg.func_code.co_firstlineno + 1
75                      )
76#        warnings.warn(("Call to deprecated function %s." % first_arg.__name__) + (text or ""),
77#                               category=DeprecationWarning,;
78#                               stacklevel =  -1 #-2 if text is None else -3
79#                        )
80        return first_arg(*args, **kwargs)
81    if isinstance(first_arg, basestring): #deprecated called with a string as first argument
82        return partial(deprecated, text=first_arg)
83    else:
84        return wraps(first_arg)(wrapper)
85       
86def cached(func):
87    @wraps(func)
88    def wrapper(self, *args, **kwargs):
89        sig = args + tuple(sorted(kwargs.items()))
90        if sig in wrapper:
91            pass
92    return wrapper
93           
94def cached_property(func):
95    @property
96    @wraps(func)
97    def wrapper_func(self):
98        cache_name = "_" + func.__name__
99        if not hasattr(self, cache_name):
100            setattr(self, cache_name, func(self))
101        return getattr(self, cache_name)
102    return wrapper_func
103
104def cached_method(func, cache_name="_cached_method_cache", store=None):
105    def wrapper(self, *args, **kwargs):
106        sig = (func.__name__,) + args + tuple(sorted(kwargs.items()))
107        if not hasattr(self, cache_name):
108            setattr(self, cache_name, store() if store is not None else {})
109        if sig not in getattr(self, cache_name):
110            getattr(self, cache_name)[sig] = func(self, *args, **kwargs)
111        return getattr(self, cache_name)[sig]
112    return wrapper
113
114class pickle_cache(dict):
115    def __init__(self, filename, sync=True):
116        dict.__init__(self)
117        self.filename = filename
118        if os.path.exists(self.filename):
119            try:
120                self.update(cPickle.load(open(self.filename, "rb")))
121            except Exception, ex:
122                print self.filename, "loading failed", ex
123               
124    def __setitem__(self, key, value):
125        super(pickle_cache, self).__setitem__(key, value)
126        self.sync()
127       
128    def __delitem__(self, key):
129        super(pickle_cache, self).__delitem__(key)
130        self.sync()
131       
132    def sync(self):
133#        print "sync", self, self.filename
134        cPickle.dump(dict(self.items()), open(self.filename, "wb"), 2)
135       
136def persistent_cached_method(func):
137    return cached_method(func, cache_name="_persistent_cached_method_cache")
138
139def persistent_cached_class(cls, store=pickle_cache):
140    def cache_loader(self):
141        if not hasattr(self, "_cache_filename"):
142            setattr(self, "_cache_filename", getattr(self, "filename") + ".cache")
143#        print "creating cache", self._cache_filename
144        try:
145#            print os.stat(self.filename).st_mtime, os.stat(self._cache_filename).st_mtime
146            if os.stat(self.filename).st_mtime > os.stat(self._cache_filename).st_mtime:
147#                print "cache obsolete"
148                os.remove(self._cache_filename)
149        except Exception, ex:
150            pass
151        return store(self._cache_filename)
152    cls._persistent_cached_method_cache = cached_property(cache_loader)
153    return cls
154
155class cache_manager(object):
156    """ takes an fuction that returns a dict like object (cache)
157    """
158    def __init__(self, opener):
159        self.opener = opener
160        self.name = opener.__name__
161       
162    def __call__(self, *args, **kwargs):
163        return self.opener(self, *args, **kwargs)
164           
165    def __get__(self, instance, type_=None):
166        """ Return an cache store from instance """
167        ## TODO: check if instance is a type
168        if getattr(instance, self.name):
169            return getattr(instance, self.name)
170        return self.opener.__get__(instance, type)
171   
172    def cached_method(manager, wrapped):
173        def wrapper(self, *args, **kwargs):
174            sig = (wrapped.__name__,) + args + tuple(sorted(kwargs.items()))
175            cache = getattr(self, manager.__name__)
176            if sig not in cache:
177                cache[sig] = func(self, *args, **kwargs)
178            return cache[sig]
179
180def proxy_dict_decorator(cls, proxy_name):
181    for name in ["__len__", "__contains__", "__getitem__", "__iter__", "get", "has_key", "items", 
182                 "iteritems", "iterkeys", "itervalues", "keys", "values", "update"]:
183        def proxy(func):
184            def f(self, *args, **kwargs):
185                return func(getattr(self, proxy_name), *args, **kwargs)
186            return f
187        setattr(cls, name, proxy(getattr(dict, name)))
188    return cls
189
190def _create_path_for_file(target): 
191    try:
192        os.makedirs(os.path.dirname(target))
193    except OSError:
194        pass
195
196def downloader(func):
197    def wrapper(*args, **kwargs):
198        import shutil
199        def download(src, dst):
200            if isinstance(dst, basestring):
201                _create_path_for_file(dst)
202                dst = open(dst, "wb")
203            shutil.copyfileobj(src, dst)
204        jobs = func(*args, **kwargs)
205        for src, dst in [jobs] if type(jobs) == tuple else jobs:
206            download(src, dst)
207    return wrapper
208
209from Orange.orng.orngMisc import with_gc_disabled
210
211SLINE, MLINE, LINELIST = range(1, 4)
212
213class _field_mutators(object):
214    @classmethod
215    def deindent(cls, text, field):
216        pad = len(text) - len(text[len(field):].lstrip())
217        return "\n".join([line[pad:] for line in text.splitlines()])
218   
219    @classmethod
220    def to_list(cls, text, field):
221        return [line for line in text.splitlines() if line.strip()]
222   
223    @classmethod
224    def to_id_list(cls, text, field):
225        return [tuple(line.split(": ", 1)) for line in cls.to_list(cls.deindent(text, field), field)]
226   
227    @classmethod
228    def to_dict(cls, text, field):
229        text = cls.deindent(text, field).replace("\n ", " ")
230        lines = [t.split(": ",1) for t in text.splitlines()]
231        return dict([(key, [t.strip() for t in value.split(" ") if t]) for key, value in lines])
232   
233    @classmethod
234    def to_ids(cls, text, field):
235        lines = [line.split(" ", 2)[:2] for line in cls.deindent(text, field).splitlines() if ":" in line.split(" ", 1)[0]]
236        return [db.lower() + id for db, id in lines] 
237#        ids = reduce(lambda col, line:1, text.splitlines(), [])
238        return ids
239    default = deindent
240#    entry = classmethod(lambda cls, text, field: cls.deindent(text, field).split()[0])
241    orthology = to_id_list
242    pathway = classmethod(lambda cls, text, field: [t[1].split()[0] for t in cls.to_id_list(text, field)])
243    codon_usage = classmethod(lambda cls, text, field: text.replace(field, " " * len(field)))
244    structure = to_id_list
245    motif = to_dict
246    pathway = to_ids
247    orthology = to_ids
248   
249    @classmethod
250    def genes(cls, text, field):
251        genes = [line.split(": ", 1) for line in cls.deindent(text, field).replace("\n    ", "").splitlines()]
252        return dict([(org.lower(), [gene.split("(")[0] for gene in genes.split(" ")]) for org, genes in genes])
253   
254    @classmethod
255    def dblinks(cls, text, field):
256        links = [line.split(": ", 1) for line in cls.deindent(text, field).replace("\n   ", "").splitlines()]
257        return dict([(db.lower(), [name for name in links.split(" ")]) for db, links in links])
258   
259    dblinks = to_dict
260
261def entry_decorate(cls):
262    reserved_map = {"class": "class_"}
263    for field in cls.FIELDS:
264        mutator = getattr(_field_mutators, field.lower(), _field_mutators.default)
265        def construct(field, mutator):
266            def get(self):
267                text = self.get(field)
268                return mutator(text, field) if text else None
269            return get
270        if not hasattr(cls, field.lower()):
271            setattr(cls, reserved_map.get(field.lower(), field.lower()), property(construct(field, mutator)))
272        else:
273            import warnings
274            warnings.warn(str(cls)+ "already contains " + field.lower())
275    return cls
276
277def borg_class(cls):
278    borg__init__ = cls.__init__
279    cls._borg_instance = None
280    def __init__(self, *args, **kwargs):
281        if not cls._borg_instance:
282            borg__init__(self, *args, **kwargs)
283            cls._borg_instance = args + tuple(sorted(kwargs.items())), self
284        elif cls._borg_instance[0] == args + tuple(sorted(kwargs.items())):
285            self.__dict__ = cls._borg_instance[1].__dict__
286        else:
287            print "Warning. Atempting to create an individual instance of '%s'.\nWe are borg. You will be assimilated. Resistance is futile." % cls.__name__
288            self.__dict__ = cls._borg_instance[1].__dict__
289    cls.__init__ = __init__
290    return cls
291
292def defaultpath(func):
293#    @wraps
294    def wrapper(*args, **kwargs):
295        if "path" not in kwargs:
296            kwargs["path"] = DEFAULT_DATABASE_PATH
297        return func(*args, **kwargs)
298    return wrapper
299             
300       
301class KEGGDBEntry(object):
302    FIELDS = []
303    MULTIPLE_FIELDS = []
304    def __init__(self, entrytext, index=None):
305        self.entrytext = entrytext
306        if not self.FIELDS:
307            self.FIELDS = set([line.split()[0] for line in entrytext.splitlines() if line.strip() and not line.startswith(" ")])
308        self._index(index)
309       
310    def _indices(self, field):
311        text = self.entrytext
312        fieldlen = len(field)
313        return reduce(lambda indices,s: indices.append(text.find(field, indices[-1] + fieldlen)) or indices, range(text.count(field)), [-fieldlen])[1:]
314         
315    def _index(self, index=None):
316        if index is not None:
317            index, self.partitions = index
318            self.index = dict(zip(self.FIELDS, index))
319            return 
320        index = dict([(field, self._indices(field)) for field in self.FIELDS])
321        sorted_indices = sorted(reduce(set.union, index.values(), set())) + [-1]
322        self.partitions = sorted_indices #[(sorted_indices[i], sorted_indices[i+1]) for i in range(len(sorted_indices) - 1)]
323        self.index = dict([(key, [sorted_indices.index(ind) for ind in indices]) for key, indices in index.items()])
324       
325    def section_iter(self, field):
326        for ind in self.index[field]:
327            yield self.entrytext[self.partitions[ind]: self.partitions[ind+1]]
328       
329    def get(self, key, default=None):
330        if key not in self.FIELDS:
331            return default
332        sections = [section for section in self.section_iter(key)]
333        if key in self.MULTIPLE_FIELDS:
334            return sections
335        elif sections:
336            return sections[0]
337        else: 
338            return default
339       
340    def entry_key(self):
341        return self.entry.split()[0]
342       
343class KEGGDataBase(object):
344    Entry = KEGGDBEntry
345    VERSION = "v1.0"
346    def __init__(self, file=None):
347        self.load(file)
348       
349    @with_gc_disabled
350    def load(self, file=None):
351        if file is None:
352            file = self.FILENAME % {"path":DEFAULT_DATABASE_PATH}
353        self.filename = file
354        data = open(file, "rb").read().split("///\n")
355        file_timestamp = os.stat(self.filename).st_mtime
356#        print file_timestamp
357        build_index = False
358#        if os.path.exists(self.filename + ".index"):
359        try:
360            version, timestamp, index = cPickle.load(open(self.filename + ".index", "rb"))
361            assert(version == self.VERSION and timestamp == file_timestamp)
362        except Exception, ex:
363#            print "error", ex
364            build_index = True
365            index = [None] * len(data)
366        self.entrys = [self.Entry(entry, ind) for entry, ind in zip(data, index) if entry.strip()]
367        if build_index:
368#            self.entrys = [self.Entry(entry) for entry in data if entry.strip()]
369            index = [([e.index[field] for field in e.FIELDS], e.partitions) for e in self.entrys]
370            cPickle.dump((self.VERSION, file_timestamp, index), open(self.filename + ".index", "wb"), cPickle.HIGHEST_PROTOCOL)
371       
372        self.entry_dict = dict([(entry.entry_key(), entry) for entry in self.entrys])
373       
374proxy_dict_decorator(KEGGDataBase, "entry_dict")
375       
376class KEGGGeneEntry(KEGGDBEntry):
377    FIELDS = ["ENTRY", "NAME", "DEFINITION","ORTHOLOGY", "PATHWAY", "CLASS", "POSITION",
378              "MOTIF", "DBLINKS", "STRUCTURE", "CODON_USAGE", "AASEQ", "NTSEQ"]
379   
380    def aliases(self):
381        return [self.entry_key()] + (self.name.split(",") if self.name else []) + [link[1][0] for link in self.dblinks.items() if self.dblinks]
382   
383    @property
384    def alt_names(self):
385        return self.aliases()
386       
387entry_decorate(KEGGGeneEntry)
388   
389class KEGGGenes(KEGGDataBase):
390    Entry = KEGGGeneEntry
391    FILENAME = "genes/organisms/%(org_code)s/%(org_name)s.ent"
392   
393    @loads("KEGG", lambda self, org: "kegg_genes_%s.tar.gz" % org)
394    def load(self, org):
395        super(KEGGGenes, self).load(_join(self.filename(org)))
396        self.entry_dict = dict([(org + ":" + key, value) for key, value in self.entry_dict.items()])
397        self.org_code = org
398       
399    def gene_aliases(self):
400        aliases = {}
401        for entry in self.entrys:
402            aliases.update(dict.fromkeys(entry.aliases(), self.org_code + ":" + entry.entry_key()))
403        return aliases
404           
405    @classmethod
406    def filename(cls, org):
407        return cls.FILENAME % {"org_code":org, "org_name":KEGGGenome()[org].name.split(",")[1].strip()}
408     
409    @classmethod
410    @downloader
411    def download(cls, org):
412#        local_dir = orngServerFiles.localpath("KEGG") if local_dir is None else local_dir
413        filename = cls.filename(org)
414        return (urllib2.urlopen(KEGG_FTP_PATH + filename), 
415                _join(filename))
416   
417class KEGGGenomeEntry(KEGGDBEntry):
418    FIELDS = ["ENTRY", "NAME", "DEFINITION", "ANNOTATION", "TAXONOMY", "DATA_SOURCE",
419              "ORIGINAL_DB", "CHROMOSOME", "STATISTICS", "REFERENCE"]
420    MULTIPLE_FIELDS = ["REFERENCE"]
421   
422    def org_code(self):
423        if self.name is not None:
424            return self.name.split(",")[0]
425        else:
426            return self.entry.split()[0]
427       
428    def entry_key(self):
429        return self.org_code()
430   
431entry_decorate(KEGGGenomeEntry)
432       
433class KEGGGenome(KEGGDataBase):
434    Entry = KEGGGenomeEntry
435    VERSION = "v1.1"
436    TAXID_MAP = {"4932" : "559292"}
437    def name(self):
438        return super(KEGGGeneome, self).name()
439   
440    @loads("KEGG", "kegg_genome.tar.gz", version=VERSION)
441    def load(self, file=None):
442        filename = _join("genes", "genome")
443        super(KEGGGenome, self).load(filename)
444   
445    @classmethod
446    def common_organisms(cls):
447#        genome = KEGGGenome()
448#        id_map = {"562":"511145", "2104":"272634", "5833":"36329", "4896":"284812", "11103":None, "4754":None, "4577":None}
449#        return [genome.get(id_map.get(taxid, taxid)).entry_key() for taxid in obiTaxonomy.common_taxids() if id_map.get(taxid, taxid) is not None]
450        return ['ath', 'bta', 'cel', 'cre', 'dre', 'ddi', 'dme', 'eco', 'hsa', 'mmu', 'mpn', 'osa',
451                'pfa', 'rno', 'sce', 'spo', 'zma', 'xla']
452       
453    @classmethod
454    def essential_organisms(cls):
455        genome = KEGGGenome()
456       
457        return [genome.get(cls.TAXID_MAP.get(taxid, taxid)).entry_key() for taxid in obiTaxonomy.essential_taxids()]
458           
459    def get(self, key, *args):
460        if key in self.TAXID_MAP:
461            key = self.TAXID_MAP[key]
462        if key not in self:
463            keys = self.search(key)
464            keys = [k for k in keys if key in self[k].name or key in self[k].taxonomy]
465            key = keys.pop(0) if keys else key
466        return super(KEGGGenome, self).get(key, *args)
467   
468    def search(self, string, relevance=False):
469        if string in self.TAXID_MAP:
470            string = self.TAXID_MAP[string]
471           
472        def match(entry, string):
473            rel = 0
474            if string in entry.entrytext:
475                weight_f = [(entry.definition or "", 2), (entry.name or "", 4), (entry.taxonomy or "", 8), (entry.entry_key(), 16)]
476                rel += 1 + sum([w for text, w in weight_f if string in text]) + (16 if entry.entry_key() in self.common_organisms() else 0)
477            return rel
478        matched = sorted(zip([match(entry, string) for entry in self.entrys], self.entrys), reverse=True)
479        return [(entry.entry_key(), rel) if relevance else entry.entry_key() for rel, entry in matched if rel]
480           
481    @classmethod
482    @downloader
483    def download(cls, file=None):
484        return (urllib2.urlopen(KEGG_FTP_PATH + "/genes/genome"), 
485                _join("genes", "genome") if file is None else file)
486       
487borg_class(KEGGGenome)
488       
489class KEGGCompoundEntry(KEGGDBEntry):
490    FIELDS = ["ENTRY", "NAME", "FORMULA", "MASS", "REMARK", "REACTION", "PATHWAY",
491              "ENZYME", "DBLINKS", "ATOM", "BOND"]
492    def entry_key(self):
493        return "cpd:" + super(KEGGCompoundEntry, self).entry_key()
494   
495entry_decorate(KEGGCompoundEntry)
496   
497class KEGGCompounds(KEGGDataBase):
498    Entry = KEGGCompoundEntry
499    FILENAME = "%(path)s/ligand/compound/compound"
500   
501    @loads("KEGG", "kegg_ligand.tar.gz")
502    def load(self, file=None):
503        super(KEGGCompounds, self).load(file)
504       
505    @classmethod
506    @downloader
507    def download(cls, file=None):
508        return (urllib2.urlopen(KEGG_FTP_PATH + "/ligand/compound/compound"),
509                _join("ligand", "compound", "compound") if file is None else file)
510   
511borg_class(KEGGCompounds)
512   
513class KEGGEnzymeEntry(KEGGDBEntry):
514    FIELDS = ["ENTRY", "NAME", "CLASS", "SYSNAME", "REACTION", "ALL_REAC", "SUBSTRATE",
515              "PRODUCT", "COFACTOR", "COMMENT", "REFERENCE", "PATHWAY", "ORTHOLOGY", 
516              "GENES", "STRUCTURES", "DBLINKS"]
517    MULTIPLE_FIELDS = ["REFERENCE"]
518    def entry_key(self):
519        return "EC:" + self.entry.split()[1]
520   
521entry_decorate(KEGGEnzymeEntry)
522
523class KEGGEnzymes(KEGGDataBase):
524    Entry = KEGGEnzymeEntry
525    FILENAME = "%(path)s/ligand/enzyme/enzyme"
526
527    @loads("KEGG", "kegg_ligand.tar.gz")
528    def load(self, file=None):
529        super(KEGGEnzymes, self).load(file)
530   
531    @classmethod
532    @downloader
533    def download(cls, file=None):
534        return (urllib2.urlopen(KEGG_FTP_PATH + "ligand/enzyme/enzyme"),
535                _join("ligand", "enzyme", "enzyme"))
536   
537borg_class(KEGGEnzymes)
538
539class KEGGReactionEntry(KEGGDBEntry):
540    FIELDS = ["ENTRY", "NAME", "DEFINITION", "EQUATION", "COMMENT", "RPAIR", "PATHWAY", 
541              "ENZYME", "ORTHOLOGY"]
542    def entry_key(self):
543        return "rn:" + super(KEGGReactionEntry, self).entry_key()
544   
545entry_decorate(KEGGReactionEntry)
546
547class KEGGReactions(KEGGDataBase):
548    Entry = KEGGReactionEntry
549    FILENAME = "%(path)s/ligand/reaction/reaction"
550   
551    @loads("KEGG", "kegg_ligand.tar.gz")
552    def load(self, file=None):
553        super(KEGGReactions, self).load(file)
554       
555    @classmethod
556    @downloader
557    def download(cls, file=None):
558        return (urllib2.urlopen(KEGG_FTP_PATH + "ligand/reaction/reaction"),
559                _join("ligand", "reaction", "reaction"))
560       
561borg_class(KEGGReactions)
562
563class KEGGKOEntry(KEGGDBEntry):
564    FIELDS = ["ENTRY", "NAME", "DEFINITION", "CLASS", "DBLINKS", "GENES"]
565   
566    def entry_key(self):
567        return "ko:" + super(KEGGKOEntry, self).entry_key()
568   
569entry_decorate(KEGGKOEntry)
570   
571class KEGGKO(KEGGDataBase):
572    Entry = KEGGKOEntry
573    FILENAME = "%(path)s/genes/ko"
574   
575    @loads("KEGG", "kegg_orthology.tar.gz")
576    def load(self, file=None):
577        super(KEGGKO, self).load(file)
578       
579    @classmethod
580    @downloader
581    def download(cls, file=None):
582        return (urllib2.urlopen(KEGG_FTP_PATH + "genes/ko"),
583                _join("genes", "ko"))
584       
585borg_class(KEGGKO)
586
587class KEGGBriteEntry(object):
588    _search_re = {"ids": re.compile('(?P<ids>\[.*:.*\])'),
589                  "title": re.compile(r'(<[Bb]>)?(?P<title>\b[a-zA-Z0-9_/\s,;:.+=\-\[\]{}\(\)]+?)(?(1)</[Bb]>)$'),
590                  "links": re.compile('(?P<links><a href=".+?">.*?</a>)')}
591    def __init__(self, line, entrys=None):
592        self.entrys = [] #entrys if entrys is not None else []
593        self.line = line[1:].strip()
594        for name, re in self._search_re.items():
595            search = re.search(self.line)
596            setattr(self, name, search.group(name) if search else None)
597#        self.identifiers = self.groups.get("ids")
598#        self.title = self.groups. # TODO: parse line (html, kegg identifiers ...)
599
600    def __iter__(self):
601        return iter(self.entrys)
602
603class KEGGBrite(KEGGBriteEntry):
604    VERSION = "v1.0"
605    def __init__(self, brite_id):
606        super(KEGGBrite, self).__init__("")
607        self.load(brite_id)
608       
609    @loads("KEGG", "kegg_brite.tar.gz")
610    def load(self, brite_id):
611        file = self.filename(brite_id)
612        lines = open(file, "rb").read().split("\n!\n")[1].splitlines()
613        def collect(lines, depth, collection):
614            while lines:
615                line = lines[0]
616                if line.startswith("#"):
617                    lines.pop(0)
618                elif line.startswith(depth) and len(line.strip()) > 1:
619                    collection.append(KEGGBriteEntry(lines.pop(0))) 
620                elif line[0] > depth:
621                    collect(lines, line[0], collection[-1].entrys)
622                elif line[0] < depth:
623                    return
624                else:
625                    lines.pop(0)
626                       
627        collect([line for line in lines if not line.startswith("#") and len(line) > 1], "A", self.entrys)
628   
629    @classmethod
630    @defaultpath
631    def filename(cls, brite_id, path=DEFAULT_DATABASE_PATH):
632        if path is None:
633            path = "%(path)s"
634        if brite_id.startswith("br"):
635            return ("%(path)s/brite/br/" + brite_id + ".keg") % dict(path=path)
636        elif brite_id.startswith("ko"):
637            return ("%(path)s/brite/ko/" + brite_id + ".keg") % dict(path=path)
638        else:
639            org = brite[:-5]
640            return ("%(path)s/brite/organisms/" + org + "/" + brite_id + ".keg") % dict(path=path)
641       
642    @classmethod
643    @downloader
644    def download(cls, brite_id):
645        filename = cls.filename(brite_id, path="")
646        return (urllib2.urlopen(cls.filename(brite_id, path=KEGG_FTP_PATH.rstrip("/"))),
647                cls.filename(brite_id))
648
649class KEGGOrganism(object):
650    VERSION = "v2.0"
651    DOMAIN = "KEGG"
652    def __init__(self, org, genematcher=None):
653        self.load(org)
654        self.genematcher = genematcher
655       
656#    @loads("KEGG", lambda self, org: "kegg_genes_%s.tar.gz" % self.organism_name_search(org))
657    def load(self, org):
658        org = self.organism_name_search(org)
659        self.org_code = org
660   
661    @property
662#    @deprecated("Use org_code instead")
663    def org(self):
664        return self.org_code
665   
666    @cached_property
667    def genes(self):
668        return KEGGGenes(self.org_code)
669       
670    @cached_property
671    def gene_aliases(self):
672        return self.genes.gene_aliases()
673   
674    def pathways(self, with_ids=None):
675        return [pathway for pathway, values in KEGGPathway.list(self.org_code).items() if all([id in values for id in (with_ids or [])])]
676   
677#    @deprecated("Use KEGGOrganism.pathways instead")
678    def list_pathways(self):
679        return self.pathways()
680   
681#    @deprecated
682    def get_linked_pathways(self, pathway_id):
683        return KEGGPathway(pathway_id).genes()
684   
685    def enzymes(self, genes=None):
686        enzymes = KEGGEnzymes()
687        return [enzyme.entry_key() for enzyme in enzymes.itervalues() if enzyme.genes and self.org_code in enzyme.genes]
688   
689    def get_enriched_pathways(self, genes, reference=None, prob=obiProb.Binomial(), callback=None):
690        """Return a dictionary with enriched pathways ids as keys and (list_of_genes, p_value, num_of_reference_genes) tuples as items."""
691        allPathways = defaultdict(lambda :[[], 1.0, []])
692        from Orange.orng import orngMisc
693        milestones = orngMisc.progressBarMilestones(len(genes), 100)
694        for i, gene in enumerate(genes):
695            pathways = self.pathways([gene])
696            for pathway in pathways:
697                allPathways[pathway][0].append(gene)
698            if callback and i in milestones:
699                callback(i*100.0/len(genes))
700        reference = set(reference if reference is not None else self.genes.keys())
701        for p_id, entry in allPathways.items():
702            entry[2].extend(reference.intersection(KEGGPathway(p_id).genes()))
703            entry[1] = prob.p_value(len(entry[0]), len(reference), len(entry[2]), len(genes))
704        return dict([(pid, (genes, p, len(ref))) for pid, (genes, p, ref) in allPathways.items()])
705   
706    def get_genes_by_enzyme(self, enzyme):
707        enzyme = KEGGEnzymes()[enzyme]
708        return enzyme.genes.get(self.org_code, []) if enzyme.genes else []
709   
710    def get_genes_by_pathway(self, pathway_id):
711        return KEGGPathway(pathway_id).genes()
712   
713    def get_enzymes_by_pathway(self, pathway_id):
714        return KEGGPathway(pathway_id).enzymes()
715   
716    def get_compounds_by_pathway(self, pathway_id):
717        return KEGGPathway(pathway_id).compounds()
718   
719    def get_pathways_by_genes(self, gene_ids):
720        gene_ids = set(gene_ids)
721        pathways = [self.genes[id].pathway for id in gene_ids if self.genes[id].pathway]
722        pathways = reduce(set.union, pathways, set())
723        return [id for id in pathways if gene_ids.issubset(KEGGPathway(id).genes())] 
724   
725    def get_pathways_by_enzymes(self, enzyme_ids):
726        enzyme_ids = set(enzyme_ids)
727        pathways = [KEGGEnzymes()[id].pathway for id in enzyme_ids]
728        pathwats = reduce(set.union, pathways, set())
729        return [id for id in pathways if enzyme_ids.issubset(KEGGPathway(id).enzymes())]
730   
731    def get_pathways_by_compounds(self, compound_ids):
732        compound_ids = set(compound_ids)
733        pathways = [KEGGCompounds()[id].pathway for id in compound_ids]
734        pathwats = reduce(set.union, pathways, set())
735        return [id for id in pathways if compound_ids.issubset(KEGGPathway(id).compounds())]
736   
737    def get_enzymes_by_compound(self, compound_id):
738        return KEGGCompound()[compound_id].enzyme
739   
740    def get_enzymes_by_gene(self, gene_id):
741        return self.genes[gene_id].enzymes
742   
743    def get_compounds_by_enzyme(self, enzyme_id):
744        return self._enzymes_to_compounds.get(enzyme_id)
745   
746    def get_unique_gene_ids(self, genes, caseSensitive=True):
747        """Return a tuple with three elements. The first is a dictionary mapping from unique gene
748        ids to gene names in genes, the second is a list of conflicting gene names and the third is a list
749        of unknown genes.
750        """
751        unique, conflicting, unknown = {}, [], []
752        for gene in genes:
753            names = self.genematcher.match(gene)
754            if len(names) == 1:
755                unique[names[0]] = gene
756            elif len(names) == 0:
757                unknown.append(gene)
758            else:
759                conflicting.append(gene)
760        return unique, conflicting, unknown
761   
762    @property
763    @cached_method
764    def _enzymes_to_compounds(self):
765        dd = {}
766        for val in KEGGCompounds().values():
767            dd.update(dict.fromkeys(val.enzymes, val.entry_key()))
768        return dd
769   
770    def _set_genematcher(self, genematcher):
771        setattr(self, "_genematcher", genematcher)
772       
773    def _get_genematcher(self):
774        if getattr(self, "_genematcher", None) == None:
775            from . import obiGene
776            if self.org_code == "ddi":
777                self._genematcher = obiGene.matcher([obiGene.GMKEGG(self.org_code), obiGene.GMDicty(),
778                                                     [obiGene.GMKEGG(self.org_code), obiGene.GMDicty()]])
779            else:
780                self._genematcher = obiGene.matcher([obiGene.GMKEGG(self.org_code)])
781            self._genematcher.set_targets(self.genes.keys())
782        return self._genematcher
783    genematcher = property(_get_genematcher, _set_genematcher)
784   
785    def get_genes(self):
786        return self.genes
787   
788    @classmethod
789    def organism_name_search(cls, name):
790        genome = KEGGGenome()
791        if name not in genome:
792            ids = genome.search(name)
793            if not ids:
794                from . import obiTaxonomy
795                ids = obiTaxonomy.search(name)
796                ids = [id for id in ids if genome.search(id)]
797            name = ids.pop(0) if ids else name
798           
799        try:
800            return KEGGGenome()[name].entry_key()
801        except KeyError:
802            raise ValueError("Organims with name='%s' not found in KEGG." % name)
803       
804    @classmethod
805    def organism_version(cls, name):
806        name = cls.organism_name_search(name)
807        orngServerFiles.localpath_download("KEGG", "kegg_genes_%s.tar.gz" % name)
808        return orngServerFiles.info("KEGG", "kegg_genes_%s.tar.gz" % name)["datetime"]
809             
810class KEGGPathway(object):
811    PNG_FILENAME = "%(path)s/pathway/organisms/%(org)s/%(org)s%(map_id)s.png"
812    KGML_FILENAME = "%(path)s/xml/kgml/metabolic/organisms/%(org)s/%(org)s%(map_id)s.xml"
813    VERSION = "v2.0"
814   
815    class entry(object):
816        def __init__(self, dom_element):
817            self.__dict__.update(dom_element.attributes.items())
818            self.graphics = ()
819            self.components = []
820            self.graphics = dict(dom_element.getElementsByTagName("graphics")[0].attributes.items())
821            self.components = [node.getAttribute("id") for node in dom_element.getElementsByTagName("component")]
822    class reaction(object):
823        def __init__(self, dom_element):
824            self.__dict__.update(dom_element.attributes.items())
825            self.substrates = [node.getAttribute("name") for node in dom_element.getElementsByTagName("substrate")]
826            self.products = [node.getAttribute("name") for node in dom_element.getElementsByTagName("product")]
827           
828    class relation(object):
829        def __init__(self, dom_element):
830            self.__dict__.update(dom_element.attributes.items())
831            self.subtypes = [node.attributes.items() for node in dom_element.getElementsByTagName("subtype")]
832       
833    def __init__(self, file):
834        if not os.path.exists(file):
835            path, org, map_id = self.split_pathway_id(file)
836#            file = self.KGML_FILENAME % dict(path=DEFAULT_DATABASE_PATH, org=org, map_id=id)
837            file = self.filename_kgml(org, map_id)
838        if not os.path.exists(file) and os.path.exists(file.replace("metabolic", "non-metabolic")):
839            file = file.replace("metabolic", "non-metabolic")
840        self.filename = file
841#        self.load(file)
842       
843    @persistent_cached_method
844    def pathway_attributes(self):
845        return dict(self.pathway_dom().attributes.items())
846   
847    @property
848    def name(self):
849        return self.pathway_attributes().get("name")
850   
851    @property
852    def org(self):
853        return self.pathway_attributes().get("org")
854   
855    @property
856    def number(self):
857        return self.pathway_attributes().get("number")
858   
859    @property   
860    def title(self):
861        return self.pathway_attributes().get("title")
862   
863    @property
864    def image(self):
865        return self.pathway_attributes().get("image")
866   
867    @property
868    def link(self):
869        return self.pathway_attributes().get("link")
870   
871    @cached_method
872    @loads("KEGG", lambda self: "kegg_pathways_%s.tar.gz" % self.filename.split("/")[-1][:-9])
873    def pathway_dom(self):
874        return minidom.parse(self.filename).getElementsByTagName("pathway")[0]
875   
876    @cached_method
877    def entrys(self):
878        return [self.entry(e) for e in self.pathway_dom().getElementsByTagName("entry")]
879   
880    @cached_method
881    def reactions(self):
882        return [self.reaction(e) for e in self.pathway_dom().getElementsByTagName("reaction")]
883   
884    @cached_method
885    def relations(self):
886        return [self.relation(e) for e in self.pathway_dom().getElementsByTagName("relation")]
887   
888    @classmethod
889    def split_pathway_id(cls, id):
890        path, id = id.split(":") if ":" in id else ("path", id)
891        org, id = id[:-5], id[-5:]
892        return path, org, id 
893       
894    def __iter__(self):
895        """ Iterate over all elements in the pathway
896        """
897        return iter(self.all_elements())
898   
899    def __contains__(self, element):
900        """ Retrurn true if element in the pathway
901        """
902        return element in self.all_elements()
903   
904    def __getitem__(self, key):
905        return
906   
907    @cached_method
908    def all_elements(self):
909        """ Return all elements
910        """
911        return reduce(list.__add__, [self.genes(), self.compounds(), self.enzmes(), self.reactions()], [])
912   
913    def _get_entrys_by_type(self, type):
914        return reduce(set.union, [entry.name.split() for entry in self.entrys() if entry.type == type], set())
915   
916    @persistent_cached_method
917    def genes(self):
918        """ Return all genes on the pathway
919        """
920        return self._get_entrys_by_type("gene")
921   
922    @persistent_cached_method
923    def compounds(self):
924        """ Return all compounds on the pathway
925        """
926        return self._get_entrys_by_type("compound")
927   
928    @persistent_cached_method
929    def enzymes(self):
930        """ Return all enzymes on the pathway
931        """
932        return self._get_entrys_by_type("enzyme")
933   
934    @persistent_cached_method
935    def orthologs(self):
936        """ Return all orthologs on the pathway
937        """
938        return self._get_entrys_by_type("ortholog")
939   
940    @persistent_cached_method
941    def maps(self):
942        """ Return all linked maps on the pathway
943        """
944        return self._get_entrys_by_type("map")
945   
946    @persistent_cached_method
947    def groups(self):
948        """ Return all groups on the pathway
949        """
950        return self._get_entrys_by_type("ortholog")
951
952    def get_image(self):
953        """ Return an image of the pathway
954        """
955        return self.filename_png(self.org, self.number) % dict(path=DEFAULT_DATABASE_PATH)
956   
957    @persistent_cached_method
958    def get_bounding_box_dict(self):
959        return dict([(element.id, element.graphics) for element in self.entrys() if element.graphics])
960   
961    @persistent_cached_method
962    def graphics(self, item):
963        return [entry.graphics for entry in self.entrys() if item in entry.name and entry.graphics]
964       
965    @classmethod
966    @partial(cached_method, cache_name="_cls_cached_method_cache_")
967    @loads("KEGG", lambda cls, org: "kegg_pathways_%s.tar.gz" % org)
968    def list(cls, org):
969        file = cls.directory_png(org) + org + ".list"
970        data = [line.split() for line in open(file, "rb").read().splitlines()]
971        return reduce(lambda dict, line: dict[line[0]].update(line[1:]) or dict, data, defaultdict(set))
972   
973    @classmethod
974    @defaultpath
975    def filename_kgml(cls, org, map_id, path=DEFAULT_DATABASE_PATH):
976        path = "%(path)s" if path is None else path
977        if org in ["map", "ec"]:
978            return "%(path)s/xml/kgml/metabolic/ec/ec%(map_id)s.xml" % dict(org=org, map_id=map_id, path=path)
979        elif org == "ko":
980            return "%(path)s/xml/kgml/metabolic/ko/ko%(map_id)s.xml" % dict(org=org, map_id=map_id, path=path)
981        else:
982            return "%(path)s/xml/kgml/metabolic/organisms/%(org)s/%(org)s%(map_id)s.xml" % dict(org=org, map_id=map_id, path=path)
983       
984    @classmethod
985    @defaultpath
986    def filename_png(cls, org, map_id, path=DEFAULT_DATABASE_PATH):
987        path = "%(path)s" if path is None else path
988        if org in ["map", "ec", "ko", "rn"]:
989            return "%(path)s/pathway/%(org)s/%(org)s%(map_id)s.png" % dict(org=org, map_id=map_id, path=path)
990        else:
991            return "%(path)s/pathway/organisms/%(org)s/%(org)s%(map_id)s.png" % dict(org=org, map_id=map_id, path=path)
992       
993    @classmethod
994    @defaultpath
995    def directory_png(cls, org, path=DEFAULT_DATABASE_PATH):
996        return cls.filename_png(org, "", path=path).rsplit("/", 1)[0] + "/" 
997   
998    @classmethod
999    @defaultpath
1000    def directory_kgml(cls, org, path=DEFAULT_DATABASE_PATH):
1001        return cls.filename_kgml(org, "", path=path).rsplit("/", 1)[0] + "/"
1002       
1003    @classmethod
1004    def pathways(cls, org):
1005        file = cls.directory_png(org) + org + ".list"
1006        pathways = [line.split()[0] for line in open(file, "rb").read().splitlines() if line.strip()]
1007        return sorted(set(pathways))
1008       
1009    @classmethod
1010    def download(cls, pathway, target=None):
1011        """ Download the pathway (xml and png files) to target directory
1012        """
1013        org = pathway[:3]
1014        import urllib2
1015        xml_data = urllib2.urlopen("ftp://ftp.genome.jp/pub/kegg/xml/kgml/metabolic/organisms/%s/%s.xml" % (org, pathway)).read()
1016        open(os.path.join(target, pathway + ".xml"), "wb").write(xml_data)
1017        png_data = urllib2.urlopen("ftp://ftp.genome.jp/pub/kegg/pathway/organisms/%s/%s.png" % (org, pathway)).read()
1018        open(os.path.join(target, pathway + ".png"), "wb").write(png_data)
1019       
1020    @classmethod
1021    @downloader
1022    def download_pathways(cls, org):
1023        png_path = cls.directory_png(org, path=None).replace("%(path)s/", "")
1024        xml_path = cls.directory_kgml(org, path=None).replace("%(path)s/", "")
1025        data = urllib2.urlopen("ftp://ftp.genome.jp/pub/kegg/" + png_path + org + ".list").read()
1026        pathways = sorted(set([line.split()[0][5:] for line in data.splitlines() if line]))
1027        from .obiData import FtpDownloader
1028        ftp = FtpDownloader("ftp.genome.jp", _join(), "pub/kegg/", numOfThreads=15)
1029        ftp.massRetrieve([png_path + pathway + ".png" for pathway in pathways])
1030           
1031        if org != "map" :
1032            pathways = ftp.listdir(xml_path)
1033            pathways = [p for p in pathways if p.endswith(".xml")]
1034            ftp.massRetrieve([xml_path + pathway for pathway in pathways])
1035            if org not in ["ec", "rn"]: # Add non-metabolic pathways to the download list
1036                xml_path = xml_path.replace("metabolic", "non-metabolic")
1037                pathways = ftp.listdir(xml_path)
1038                pathways = [p for p in pathways if p.endswith(".xml")] 
1039                ftp.massRetrieve([xml_path + pathway for pathway in pathways])
1040        ftp.massRetrieve([png_path + org + ".list"])
1041        return []
1042   
1043persistent_cached_class(KEGGPathway)
1044       
1045def organism_name_search(name):
1046    return KEGGOrganism.organism_name_search(name)
1047
1048def pathways(org):
1049    return KEGGPathway.list(org)
1050
1051def organisms():
1052    return KEGGOrganism.organisms()
1053
1054def to_taxid(name):
1055    genome = KEGGGenome()
1056    names = genome.search(name)
1057    if genome[names[0]].taxonomy:
1058        return re.findall(r"TAX:\s*(\d+)", genome[names[0]].taxonomy)[0]
1059    else:
1060        return None #Should this raise an error?
1061
1062def from_taxid(taxid):
1063    return KEGGGenome().search(taxid)[0]
1064   
1065def test():
1066    p = KEGGPathway("sce00010.xml")
1067    print p.genes
1068    print p.reactions
1069    print p.compounds
1070    print p.image
1071    g = KEGGGenome()
1072    org = KEGGOrganism("Homo sapiens")
1073    print list(org.genes)[:10]
1074    org.gene_aliases
1075    print org.pathways(with_ids=org.genes.keys()[:5])
1076    print org.enzymes()
1077    print org.enriched_pathways(org.genes.keys()[:10])
1078    print org.genematcher
1079
1080if __name__=="__main__":
1081    test()
1082#    org1 = KEGGOrganism("ddi")
1083#    org2 = KEGGOrganism("ddi")
1084#    org2.api = KEGGInterface()
1085#    tests = [("get_genes", ()),
1086#             ("get_genes_by_enzyme", ("ec:1.1.1.1",)),
1087#             ("get_genes_by_pathway", ("path:ddi00010",)),
1088#             ("get_pathways_by_genes", (["ddi:DDB_0191256"],)),
1089#             ("get_pathways_by_enzymes", (["ec:1.1.1.1"],)),
1090#             ("get_pathways_by_compounds", (["cpd:C00001"],)),
1091#             ("get_linked_pathways", ("path:ddi00010",)),
1092#             ("list_pathways", ()),
1093#             ("get_compounds_by_enzyme", ("ec:1.1.1.1",)),
1094#             ("get_compounds_by_pathway", ("path:ddi00010",)),
1095#             ("get_enzymes_by_compound", ("cpd:C00001",)),
1096#             ("get_enzymes_by_pathway", ("path:ddi00010",)),
1097#             ("get_enzymes_by_gene", ("ddi:DDB_0191256",))]
1098#    for name, args in tests:
1099#        s1 = set(getattr(org1, name)(*args))
1100#        s2 = set(getattr(org2, name)(*args))
1101#        if s1 and s2:
1102#            print name
1103#            print s1-s2
1104#            print s2-s1
1105#        else:
1106#            print name
1107#            print "both empty"
Note: See TracBrowser for help on using the repository browser.