source: orange-bioinformatics/_bioinformatics/obiTaxonomy.py @ 1828:857938f13079

Revision 1828:857938f13079, 15.1 KB checked in by Flashpoint <vid.flashpoint@…>, 10 months ago (diff)

Added an organisms shortname() function to the obiTaxonomy.py to add short names to geneset tags

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