source: orange-bioinformatics/_bioinformatics/obiTaxonomy.py @ 1720:354d91b1af9f

Revision 1720:354d91b1af9f, 14.1 KB checked in by markotoplak, 20 months ago (diff)

Fixed update scripts for STRING and Taxonomy.

RevLine 
[1632]1from __future__ import absolute_import, division
[579]2
3from collections import defaultdict
[1632]4import cPickle, os, shutil, sys, StringIO, tarfile, urllib2
5
6from Orange.orng import orngEnviron, orngServerFiles
7
8from . import obiData, obiGenomicsUpdate
[579]9
[707]10# list of common organisms from http://www.ncbi.nlm.nih.gov/Taxonomy
11def common_taxids():
12    """Return taxonomy IDs for common organisms."""
13    return ["3702",  # Arabidopsis thaliana
14            "9913",  # Bos taurus
15            "6239",  # Caenorhabditis elegans
16            "3055",  # Chlamydomonas reinhardtii
17            "7955",  # Danio rerio (zebrafish)
[784]18            "352472", # Dictyostelium discoideum
[707]19            "7227",  # Drosophila melanogaster
20            "562",   # Escherichia coli
21            "11103", # Hepatitis C virus
22            "9606",  # Homo sapiens
23            "10090", # Mus musculus
24            "2104",  # Mycoplasma pneumoniae
25            "4530",  # Oryza sativa
26            "5833",  # Plasmodium falciparum
27            "4754",  # Pneumocystis carinii
28            "10116", # Rattus norvegicus
29            "4932",  # Saccharomyces cerevisiae
30            "4896",  # Schizosaccharomyces pombe
31            "31033", # Takifugu rubripes
32            "8355",  # Xenopus laevis
33            "4577",  # Zea mays
34             ] 
35
36def essential_taxids():
37    """Return taxonomy IDs for organisms that are included in (default) Orange Bioinformatics installation."""
[784]38    return ["352472", # Dictyostelium discoideum
[707]39            "7227",  # Drosophila melanogaster
40            "9606",  # Homo sapiens
41            "10090", # Mus musculus
42            "4932",  # Saccharomyces cerevisiae
43            ] 
44
[579]45default_database_path = os.path.join(orngEnviron.bufferDir, "bigfiles", "Taxonomy")
46
[611]47class MultipleSpeciesException(Exception):
[581]48    pass
49
[611]50class UnknownSpeciesIdentifier(Exception):
51    pass
[579]52
[798]53def pickled_cache(filename=None, dependencies=[], version=1, maxSize=30):
54    """ Return a cache function decorator.
55    """
[804]56    def datetime_info(domain, filename):
57        try:
58            return orngServerFiles.info(domain, filename)["datetime"]
59        except IOError:
60            return orngServerFiles.ServerFiles().info(domain, filename)["datetime"]
61           
[753]62    def cached(func):
63        default_filename = os.path.join(orngEnviron.bufferDir, func.__module__ + "_" + func.__name__ + "_cache.pickle")
64        def f(*args, **kwargs):
[804]65            currentVersion = tuple([datetime_info(domain, file) for domain, file in dependencies]) + (version,)
[753]66            try:
67                cachedVersion, cache = cPickle.load(open(filename or default_filename, "rb"))
[798]68                if cachedVersion != currentVersion:
[753]69                    cache = {}
70            except IOError, er:
71                cacheVersion, cache = "no version", {}
72            allArgs = args + tuple([(key, tuple(value) if type(value) in [set, list] else value)\
73                                     for key, value in kwargs.items()])
74            if allArgs in cache:
75                return cache[allArgs]
76            else:
77                res = func(*args, **kwargs)
78                if len(cache) > maxSize:
79                    del cache[iter(cache).next()]
80                cache[allArgs] = res
[798]81                cPickle.dump((currentVersion, cache), open(filename or default_filename, "wb"), protocol=cPickle.HIGHEST_PROTOCOL)
[753]82                return res
83        return f
84
85    return cached
86
[620]87def cached(func):
88    """Cached one arg method
89    """
90    def f(self, arg):
91        if arg not in self._cache:
92            self._cache[arg] = func(self, arg)
93        return self._cache[arg]
94    f._cache = {}
95    f.__name__ = "Cached " + func.__name__
96    return f
97   
98class TextDB(object):
99    entry_start_string = chr(255)
100    entry_end_string = chr(254)+"\n"
101    entry_separator_string = chr(253)
102   
103    @property
104    def _text_lower(self):
105        if id(self._text) == self._lower_text_id:
106            return self._lower_text
107        else:
108            self._lower_text_id = id(self._text)
109            self._lower_text = self._text.lower()
110            return self._lower_text
111       
112    def __init__(self, file=None, **kwargs):
113        self._text = ""
114        self._lower_text_id = id(self._text) - 1
115        self._cache = {}
116        self.__dict__.update(kwargs)
117       
118        if file != None:
[632]119            self._text = open(file, "rb").read()
[620]120
121    def _find_all(self, string, start=0, text=None, unique=True):
122        text = text if text != None else self._text_lower
123        while True:
124            index = text.find(string, start)
125            if index != -1:
126                yield index
127                if unique:
128                    start = text.find(self.entry_start_string, index + 1)
129                else:
130                    start = index + 1
131            else:
132                raise StopIteration
133
134    def _get_entry_at(self, index, text=None):
135        text = text if text != None else self._text
136        start = text.rfind(self.entry_start_string, 0, index + 1)
137        end = text.find(self.entry_end_string, index)
138        return self._text[start+1:end]
139
140    @cached
141    def get_entry(self, id):
142        try:
143            index = self._find_all(self.entry_start_string + id + self.entry_separator_string).next()
144        except StopIteration:
145            raise KeyError, id
146        return self._get_entry_at(index)
147               
148    def search(self, string):
149        string = string.lower()
150        res = []
151        for idx in self._find_all(string):
152            entry = self._get_entry_at(idx)
153            id , rest = entry.split(self.entry_separator_string, 1)
154            self._cache[id] = entry
155            res.append(id)
156        return res
157
158    def insert(self, entry):
159        self._text += self.entry_start_string + self.entry_separator_string.join(entry) + self.entry_end_string
160
161    def __iter__(self):
162        for idx in self._find_all(self.entry_start_string):
163            entry = self._get_entry_at(idx)
164            if entry:
165                yield entry.split(self.entry_separator_string ,1)[0]
166
167    def __getitem__(self, id):
168        entry = self.get_entry(id)
169        return entry.split(self.entry_separator_string)[1:]
170
171    def __setitem__(self, id, entry):
172        self.insert([id] + list(entry))
173
174    def create(self, filename):
[632]175        f = open(filename, "wb")
[620]176        f.write(self._text)
177        def write(entry):
178            f.write(self.entry_start_string + self.entry_separator_string.join(entry) + self.entry_end_string)
179        return write
180   
[611]181class Taxonomy(object):
[620]182    __shared_state = {"_text":None, "_info":None}
183    def __init__(self):
184        self.__dict__ = self.__shared_state
185        if not self._text:
186            self.Load()
187           
188    def Load(self):
189        try:
190            self._text = TextDB(os.path.join(default_database_path, "ncbi_taxonomy.tar.gz", "ncbi_taxonomy.db"))
191            self._info = TextDB(os.path.join(default_database_path, "ncbi_taxonomy.tar.gz", "ncbi_taxonomy_inf.db"))
192            return
193        except Exception, ex:
[880]194            pass
[620]195        try:
[1693]196            orngServerFiles.download("Taxonomy", "ncbi_taxonomy.tar.gz")
[620]197            self._text = TextDB(os.path.join(default_database_path, "ncbi_taxonomy.tar.gz", "ncbi_taxonomy.db"))
198            self._info = TextDB(os.path.join(default_database_path, "ncbi_taxonomy.tar.gz", "ncbi_taxonomy_inf.db"))
199            return
200        except Exception, ex:
201            raise
[656]202
203    def get_entry(self, id):
204        try:
205            entry = self._text[id]
206        except KeyError:
[753]207            raise UnknownSpeciesIdentifier
[656]208        return entry
[620]209               
210    def search(self, string, onlySpecies=True):
211        res = self._text.search(string)
212        if onlySpecies:
[632]213            res = [r for r in res if "species" in self._text[r][1]]
[620]214        return res
[579]215
[620]216    def __iter__(self):
217        return iter(self._text)
218
219    def __getitem__(self, id):
[656]220        entry = self.get_entry(id)
[620]221        return entry[2] ## item with index 2 is allways scientific name
222
223    def other_names(self, id):
[656]224        entry = self.get_entry(id)
[620]225        info = self._info[id]
[661]226        names = entry[2:] ## index 2 and larger are names
[620]227        return list(zip(names, info))[1:] ## exclude scientific name
[656]228
229    def rank(self, id):
230        entry = self.get_entry(id)
231        return entry[1]
232
233    def parent(self, id):
234        entry = self.get_entry(id)
235        return entry[0]
[661]236
237    def subnodes(self, id, levels=1):
238        res = self._text.search(self._text.entry_separator_string + id + self._text.entry_separator_string)
239        res = [r for r in res if self.get_entry(r)[0] == id]
240        if levels > 1:
241            for r in list(res):
242                res.extend(self.subnodes(r, levels-1))
243        return res
[708]244
245    def taxids(self):
246        return list(self)
[661]247   
[620]248    @staticmethod
249    def ParseTaxdumpFile(file=None, outputdir=None, callback=None):
250        from cStringIO import StringIO
[1720]251        import Orange.utils
[620]252        if file == None:
[1720]253            so = StringIO()
254            Orange.utils.wget("ftp://ftp.ncbi.nih.gov/pub/taxonomy/taxdump.tar.gz", dst_obj=so)
255            file = tarfile.open(None, "r:gz", StringIO(so.getvalue()))
256            so.close()
257        elif type(file) == str:
[620]258            file = tarfile.open(file)
259        names = file.extractfile("names.dmp").readlines()
260        nodes = file.extractfile("nodes.dmp").readlines()
261        namesDict = defaultdict(list)
[611]262        for line in names:
263            if not line.strip():
264                continue
265            line = line.rstrip("\t\n|").split("\t|\t")
266            id, name, unique_name, name_class = line
267            if unique_name:
268                namesDict[id].append((unique_name , name_class))
269            else:
270                namesDict[id].append((name , name_class))
271
[620]272        nodesDict = {}
[611]273        for line in nodes:
274            if not line.strip():
275                continue
276            line = line.split("\t|\t")[:3]
277            id, parent, rank = line
278            nodesDict[id] = (parent, rank)
[620]279       
280        name_class_codes = defaultdict(iter(range(255)).next)
281        name_class_codes["scientific name"] ## Force scientific name to be first
282        if outputdir == None:
283            outputdir = default_database_path
284        text = TextDB().create(os.path.join(outputdir, "ncbi_taxonomy.db"))
285        info = TextDB().create(os.path.join(outputdir, "ncbi_taxonomy_inf.db"))
[1720]286        milestones = set(range(0, len(namesDict), max(int(len(namesDict)/100), 1)))
[620]287        for i, (id, names) in enumerate(namesDict.items()):
288            parent, rank = nodesDict[id]
289            ## id, parent and rank go first
290            entry = [id, parent, rank]
291            ## all names and name class codes pairs follow ordered so scientific name is first
292            names = sorted(names, key=lambda (name, class_): name_class_codes[class_])
293            entry.extend([name for name ,class_ in names])
294            info_entry = [id] + [class_ for name, class_ in names]
295            text(entry)
296            info(info_entry)
297            if callback and i in milestones:
298                callback(i)
[579]299   
300def name(taxid):
[649]301    """ Return the scientific name for organism with taxid.
302    """
[620]303    return Taxonomy()[taxid]
[579]304
305def other_names(taxid):
[649]306    """ Return a list of (name, name_type) tuples but exclude the scientific name.
307    """
[622]308    return  Taxonomy().other_names(taxid)
[579]309
[798]310@pickled_cache(None, [("Taxonomy", "ncbi_taxonomy.tar.gz")], version=1)
[713]311def search(string, onlySpecies=True, exact=False):
[649]312    """ Search the NCBI taxonomy database for an organism
313    Arguments::
314            - *string*      Search string
315            - *onlySpecies* Return only taxids of species (and subspecies)
[713]316            - *exact*       Return only taxids of organism that exactly match the string
[649]317    """
[713]318    ids = Taxonomy().search(string, onlySpecies)
319    if exact:
[753]320        ids = [id for id in ids if string in [name(id)] + [t[0] for t in other_names(id)]]
[713]321    return ids
[620]322
[656]323def lineage(taxid):
324    """ Return a list of taxids ordered from the topmost node (root) to taxid.
325    """
326    tax = Taxonomy()
327    result = [taxid]
328    while True:
329        parent = tax.parent(result[-1])
330        result.append(parent)
331        if tax[parent] == "root" or parent=="1":
332            break
333    result.reverse()
334    return result
335   
[658]336def to_taxid(code, mapTo=None):
[753]337    """ See if the code is a valid code in any database and return a set of its taxids.
[649]338    """
[1716]339    from . import obiKEGG, obiGO
[649]340    results = set()
[1716]341    for test in [obiKEGG.to_taxid, obiGO.to_taxid]:
[649]342        try:
343            r = test(code)
344            if type(r) == set:
345                results.update(r)
346            else:
347                results.add(r)
348        except Exception, ex:
349            pass
350
[658]351    if mapTo:
[661]352        mapped = [[parent_id for parent_id in mapTo if parent_id in lineage(r)].pop() for r in results if any(parent_id in lineage(r) for parent_id in mapTo)]
353        if not mapped:
354            t = Taxonomy()
355            subnodes = dict([(r, t.subnodes(r, 4)) for r in results])
356            mapped = [id for id in mapTo if any(id in subnodes[r] for r in results)]
357        results = mapped
[658]358
[649]359    return results
360
[708]361def taxids():
362    return Taxonomy().taxids()
363
[1632]364from . import obiGenomicsUpdate
[620]365
366class Update(obiGenomicsUpdate.Update):
367    def GetDownloadable(self):
368        return [Update.UpdateTaxonomy]
369   
370    def IsUpdatable(self, func, args):
371        from datetime import datetime
[1632]372        from . import obiData
[620]373        if func == Update.UpdateTaxonomy:
374##            stream = urllib2.urlopen("ftp://ftp.ncbi.nih.gov/pub/taxonomy/taxdump.tar.gz")
375##            date = datetime.strptime(stream.headers.get("Last-Modified"), "%a, %d %b %Y %H:%M:%S %Z")
376            ftp = obiData.FtpWorker("ftp.ncbi.nih.gov")
377            size, date = ftp.statFtp("pub/taxonomy/taxdump.tar.gz")
378            return date > self.GetLastUpdateTime(func, args)
379
380    def UpdateTaxonomy(self):
381        Taxonomy.ParseTaxdumpFile(outputdir=self.local_database_path)
382        import tarfile
383        tFile = tarfile.open(os.path.join(self.local_database_path, "ncbi_taxonomy.tar.gz"), "w:gz")
384        tFile.add(os.path.join(self.local_database_path, "ncbi_taxonomy.db"), "ncbi_taxonomy.db")
385        tFile.add(os.path.join(self.local_database_path, "ncbi_taxonomy_inf.db"), "ncbi_taxonomy_inf.db")
386        tFile.close()
[579]387
388if __name__ == "__main__":
[581]389    ids = search("Homo sapiens")
390    print ids
391    print other_names(ids[0])
[632]392   
Note: See TracBrowser for help on using the repository browser.