source: orange-bioinformatics/obiKEGG.py @ 1348:3642b9d48d66

Revision 1348:3642b9d48d66, 41.3 KB checked in by markotoplak, 3 years ago (diff)

Removed dependency on orngServerFiles.createPath

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