source: orange-bioinformatics/obiGeneMania.py @ 1527:66a207e0663b

Revision 1527:66a207e0663b, 32.7 KB checked in by ales_erjavec <ales.erjavec@…>, 2 years ago (diff)

Now using Orange.network.Graph.

Line 
1""" Interface to retrieve gene networks from GeneMANIA server
2
3Example::
4   
5    >>> conn = Connection("http://localhost:8080/genemania")
6    >>> net = conn.retrieve(org="3702", ["PHYB", "ELF3", 'COP1", "SPA1", "FUS9"])
7    >>> net.save("graph.net")
8    >>> net.retrieve(org="3702", genes=["CIP1"], m="bp", r=100).save("CIP1.net")
9
10   
11
12"""
13
14import urllib2
15import urllib
16import re
17import posixpath
18from xml.dom import minidom
19
20import orange
21import orngNetwork
22
23DEFAULT_SERVER = "http://193.2.72.57:8080/genemania"
24
25_TAX_ID_2_INDEX = {"3702": 1,
26                   "6239": 2,
27                   "7227": 3,
28                   "9606": 4,
29                   "10090": 5,
30                   "4932": 6
31                   }
32class Connection(object):
33    _RE_TOKEN = re.compile(r'<li\s+id\s*=\s*"menu_save"\s*token\s*=\s*"([0-9]+)"><label>Save</label>')
34    _RE_NETWORK_TAB = re.compile(r'^<div\s*id\s*=\s*"networks_tab"\s*class\s*=\s*"tab">*?^</div>', re.MULTILINE)
35    _RE_NETWORK_GROUP_NAMES = re.compile(r'<div\s*class\s*=\s*"network_name"\s*id\s*=\s*"networkGroupName([0-9]+)"\s*>\s*([a-zA-Z0-9_\- ]+)\s*</div>')
36    _RE_NETWORK_NAMES = re.compile(r'<div\s*class\s*=\s*"network_name"\s*id\s*=\s*"networkName([0-9]+)"\s*>\s*([a-zA-Z0-9_\- ]+)\s*</div>')
37
38    def __init__(self, address=DEFAULT_SERVER):
39        """ Construct a Connection instance for GeneMANIA server at `address`
40       
41        :param address: URL address of GeneMANIA server
42        :type address: str
43        """
44        self.address = address
45                 
46       
47    def retrieveXML(self, org="9606", genes=[], m="automatic", r=10, token=None):
48        """ Same as `retrieve` but return the network as an xml string
49        """
50        if token is None:
51            page = self.retrieveHtmlPage(org, genes, m, r)
52#            query = self._queryPage(org, genes, m, r)
53#            stream = urllib2.urlopen(query)
54#            page = stream.read()
55            match = self._RE_TOKEN.findall(page)
56       
57            if match:
58                token = match[0]
59            else:
60                raise ValueError("Invalid query. %s" % self._queryPage(org, genes, m, r))
61       
62        query = self._queryGraph(token)
63        stream = urllib2.urlopen(query)
64        graph = stream.read()
65        self._graph = graph
66        return graph
67   
68    def retrieveHtmlPage(self, org="9606", genes=[], m="automatic", r=10):
69        """ Retrieve the HTML page (contains token to retrieve the graph, network descriptions ...)"
70        """
71        query = self._queryPage(org, genes, m, r)
72        stream = urllib2.urlopen(query)
73        page = stream.read()
74        self._page = page
75        return page
76   
77    def validate(self, org, genes):
78        """ Validate gene names for organism. Return a two
79        tuple, one with known and one with unknown genes
80        """
81       
82        organism = _TAX_ID_2_INDEX.get(org, 1)
83        genes = "; ".join(genes)
84        data = urllib.urlencode([("organism", str(organism)), ("genes", genes)])
85        validatorUrl = posixpath.join(self.address, "validator")
86        stream = urllib2.urlopen(validatorUrl, data)
87        response = stream.read()
88        dom = minidom.parseString(response)
89        return parseValidationResponse(dom)
90       
91       
92       
93    def _queryPage(self, org, genes, m, r):
94        return posixpath.join(self.address, "link?o=%s&g=%s&m=%s&r=%i" % (org, "|".join(genes), m, r)) 
95   
96   
97    def _queryGraph(self, token):
98        return posixpath.join(self.address, "pages/graph.xhtml?token=%s" % token)
99   
100   
101    def retrieve(self, org, genes, m="automatic", r=10):
102        """ Retrieve orngNetwork.Network instance representing the network for
103        the query, (See
104        `http://193.2.72.57:8080/genemania/pages/help.jsf#section/link`_ for
105        more details)
106       
107        :param org: NCBI taxonomy identifier (A. thaliana=3702, C. elegans=6239,
108                    D. melanogaster=7227, H. sapiens=9606, M. musculus=10090
109                    S. cerevisiae=4932)
110        :type org: str
111       
112        :param genes: query genes
113        :type genes: list
114       
115        :param m: network combining method; must be one of the following:
116                    * "automatic_relevance": Assigned based on query genes
117                    * "automatic": Automatically selected weighting method
118                       (Default)
119                    * "bp": biological process based
120                    * "mf": molecular function based
121                    * "cc": cellular component based
122                    * "average": equal by data type
123                    * "average_category: equal by network
124        :type m: str
125       
126        :param r: the number of results generated by GeneMANIA (must be in
127                  range 1..100
128        :type r: int
129        """
130        xml = self.retrieveXML(org, genes, m, r)
131        dom = minidom.parseString(xml)
132        graph = parse(dom)
133        return graph
134   
135       
136   
137def parse(DOM):
138    """ Parse the graph DOM as returned from geneMANIA server and return
139    an orngNetwork.Network instance
140    """
141    nodes = DOM.getElementsByTagName("node")
142    edges = DOM.getElementsByTagName("edge")
143    from collections import defaultdict
144    graphNodes = {}
145    graphEdges = defaultdict(list)
146   
147    def parseAttributes(element):
148        return dict([(key, value) for key, value in element.attributes.items()])
149   
150    def parseText(element):
151        text = u""
152        for el in element.childNodes:
153            if isinstance(el, minidom.Text):
154                text += el.wholeText
155        return text
156               
157    def parseData(node):
158        data = node.getElementsByTagName("data")
159        parsed = {}
160        for el in data:
161            attrs = parseAttributes(el)
162            key = attrs["key"]
163            parsed[key] = parseText(el)
164        return parsed
165   
166    for node in nodes:
167        attrs = parseAttributes(node)
168        id = attrs["id"]
169        data = parseData(node)
170        graphNodes[id] = data
171   
172    for edge in edges:
173        attrs = parseAttributes(edge)
174        source, target = attrs["source"], attrs["target"]
175        data = parseData(edge)
176        graphEdges[source, target].append(data)
177       
178    allData = reduce(list.__add__, graphEdges.values(), [])
179    edgeTypes = set([int(data["networkGroupId"]) for data in allData])
180    groupId2int = dict(zip(edgeTypes, range(len(edgeTypes))))
181    groupId2groupCode = dict([(int(data["networkGroupId"]), str(data["networkGroupCode"])) for data in allData])
182    graphNode2nodeNumber = dict(zip(graphNodes, range(len(graphNodes))))
183   
184    import Orange
185    graph = Orange.network.Graph()
186    for id, data in graphNodes.items():
187        graph.add_node(graphNode2nodeNumber[id],
188                       original_id=str(id),
189                       symbol=data["symbol"],
190                       score=float(data["score"]))
191         
192    graph.add_nodes_from(sorted(graphNode2nodeNumber.values()))
193   
194    edgeWeights = []
195    for (source, target), edge_data in graphEdges.items():
196        edgesDefined = [None] * len(edgeTypes)
197        for data in edge_data:
198            networkGroupId = int(data["networkGroupId"])
199            edgeInd = groupId2int[networkGroupId]
200            edgesDefined[edgeInd] = float(data["weight"])
201            graph.add_edge(graphNode2nodeNumber[source], 
202                           graphNode2nodeNumber[target],
203                           weight=float(data["weight"]),
204                           networkGroupId=networkGroupId)
205           
206        edgesDefined = [0 if w is None else w for w in edgesDefined]
207        edgeWeights.append(edgesDefined)
208       
209       
210    nodedomain = orange.Domain([orange.StringVariable("label"),
211                                orange.StringVariable("id"),
212                                orange.FloatVariable("score"),
213                                orange.StringVariable("symbol"),
214                                orange.StringVariable("go"),
215                                orange.EnumVariable("source", values=["true", "false"])], None)
216   
217    edgedomain = orange.Domain([orange.FloatVariable("u"),
218                                orange.FloatVariable("v")] +\
219                               [orange.FloatVariable("weight_%s" % groupId2groupCode[id]) for id in edgeTypes],
220                               None)
221   
222    node_items = graphNodes.items()
223    node_items = sorted(node_items, key=lambda t: graphNode2nodeNumber[t[0]])
224   
225    nodeitems = orange.ExampleTable(nodedomain,
226                  [[str(node["symbol"]), str(id), float(node["score"]),
227                    str(node["symbol"]), str(node["go"]), str(node["source"])]\
228                     for id, node in node_items])
229   
230    edgeitems = orange.ExampleTable(edgedomain,
231                  [[str(graphNode2nodeNumber[source] + 1), 
232                    str(graphNode2nodeNumber[target] + 1)] + weights \
233                   for ((source, target), _), weights in zip(graphEdges.items(), edgeWeights)])
234       
235    graph.set_items(nodeitems)
236    graph.set_links(edgeitems)
237   
238    return graph
239
240def parseValidationResponse(dom):
241    def getData(node):
242        data = []
243        for c in node.childNodes:
244            if c.nodeType == node.TEXT_NODE:
245                data.append(c.data)
246               
247        return " ".join([d.strip() for d in data])
248       
249    def getStrings(node):
250        strings = []
251        for string in node.getElementsByTagName("string"):
252            strings.append(getData(string))
253        return strings
254    errorCode = dom.getElementsByTagName("errorCode")[0]
255    errorCode = getData(errorCode)
256    invalidSymbols = getStrings(dom.getElementsByTagName("invalidSymbols")[0])
257    geneIds = getStrings(dom.getElementsByTagName("geneIds")[0])
258   
259    return errorCode, invalidSymbols, geneIds
260   
261
262from HTMLParser import HTMLParser
263
264class NetworkGroup(object):
265    """ Network group descriptor
266    """
267    def __init__(self):
268        self.weight = ""
269        self.networks = []
270        self.name = ""
271        self.id = ""
272
273
274class Network(object):
275    """ Source network descriptor
276    """
277   
278    def __init__(self):
279        self.weight = ""
280        self.name = ""
281        self.id = ""
282        self.description = ""
283       
284       
285class _NetworkTabParser(HTMLParser):
286    """ Parses the "Network" tab from the GeneMANIA HTML pages
287    """
288    _RE_GROUP_ID = re.compile(r"networkGroup(\d+)")
289    _RE_GROUP_WEIGHT_ID = re.compile(r"networkGroupWeight(\d+)")
290    _RE_GROUP_NAME_ID = re.compile(r"networkGroupName(\d+)")
291   
292    _RE_NETWORK_ID = re.compile(r"network(\d+)")
293    _RE_NETWORK_WEIGHT_ID = re.compile(r"networkWeight(\d+)")
294    _RE_NETWORK_NAME_ID = re.compile(r"networkName(\d+)")
295    _RE_NETWORK_DESCRIPTION_ID = re.compile("networkDescription(\d+)")
296   
297   
298    def __init__(self, *args, **kwargs):
299        HTMLParser.__init__(self)
300        self.networkGroups = []
301        self.networks = {}
302       
303        self.currentGroup = None
304        self.currentNetwork = None
305       
306        self.data_handler = None
307       
308    def handle_start_group(self, tag, attrs):
309        """ Handle '<li class=... id="networkGroup%i">'
310        """
311        self.currentGroup = NetworkGroup()
312        self.currentGroup.id = attrs.get("id")
313       
314        self.networkGroups.append(self.currentGroup)
315       
316       
317    def handle_start_group_weight(self, tag, attrs):
318        """ Handle '<span tooltip="..." id="networkGroupWeight%i">'
319        """
320        self.data_handler = self.handle_group_weight_data
321       
322    def handle_group_weight_data(self, data):
323        self.currentGroup.weight += data
324       
325    def handle_end_group_weight(self, tag):
326        self.data_handler = None
327       
328    def handle_start_group_name(self, tag, attrs):
329        """ Handle '<div class="network_name" id="networkGroupName%i">'
330        """
331        self.data_handler = self.handle_group_name_data
332       
333    def handle_group_name_data(self, data):
334        self.currentGroup.name += data
335       
336    def handle_start_network(self, tag, attrs):
337        """ Handle '<li class="checktree_network" id="network%i">'
338        """
339        self.currentNetwork = Network()
340        self.currentNetwork.id = attrs.get("id")
341       
342        self.currentGroup.networks.append(self.currentNetwork)
343       
344    def handle_start_network_weight(self, tag, attrs):
345        """ Handle '<span tooltip="..." id="networkWeight%i">'
346        """
347        self.data_handler = self.handle_network_weight_data
348       
349    def handle_network_weight_data(self, data):
350        self.currentNetwork.weight += data
351       
352    def handle_start_network_name(self, tag, attrs):
353        """ Handle '<div class="network_name" id="networkName%i">'
354        """
355        self.data_handler = self.handle_network_name_data
356       
357    def handle_network_name_data(self, data):
358        self.currentNetwork.name += data
359       
360    def handle_start_network_description(self, tag, attrs):
361        """ Handle '<div class="text" id="networkDescription%i">'
362        """
363        self.data_handler = self.handle_network_description_data
364       
365    def handle_network_description_data(self, data):
366        self.currentNetwork.description += data
367       
368    def handle_data(self, data):
369        if self.data_handler:
370            self.data_handler(data)
371   
372    def handle_starttag(self, tag, attrs):
373        attrs = dict(attrs)
374        if tag == "li" and self._RE_GROUP_ID.search(attrs.get("id", "")):
375            self.handle_start_group(tag, attrs)
376        elif tag == "span" and self._RE_GROUP_WEIGHT_ID.search(attrs.get("id", "")):
377            self.handle_start_group_weight(tag, attrs)
378        elif tag == "div" and self._RE_GROUP_NAME_ID.search(attrs.get("id", "")):
379            self.handle_start_group_name(tag, attrs)
380        elif tag == "li" and self._RE_NETWORK_ID.search(attrs.get("id", "")):
381            self.handle_start_network(tag, attrs)
382        elif tag == "span" and self._RE_NETWORK_WEIGHT_ID.search(attrs.get("id", "")):
383            self.handle_start_network_weight(tag, attrs)
384        elif tag == "div" and self._RE_NETWORK_NAME_ID.search(attrs.get("id", "")):
385            self.handle_start_network_name(tag, attrs)
386        elif tag == "div" and self._RE_NETWORK_DESCRIPTION_ID.search(attrs.get("id", "")):
387            self.handle_start_network_description(tag, attrs)
388        else:
389            HTMLParser.handle_starttag(self, tag, attrs)
390           
391    def handle_endtag(self, tag):
392        self.data_handler = None
393           
394
395def parsePage(html):
396    parser = _NetworkTabParser()
397    parser.feed(html)
398    return parser.networkGroups
399   
400
401def retrieve(org=None, genes=[], m="automatic", r=10):
402    """ A helper function, same as Connection().retrive(*args, **kwargs)
403    """
404    return Connection().retrieve(org, genes, m, r)
405
406
407"""
408======================
409PPI Database interface
410======================
411
412"""
413
414
415import sqlite3
416import csv
417import os
418import posixpath
419
420from contextlib import contextmanager
421import StringIO
422
423@contextmanager
424def finishing(obj):
425    """ Calls obj.finish() on context exit.
426    """
427    yield obj
428    obj.finish()
429
430def guess_size(fileobj):
431    try:
432        if isinstance(fileobj, file):
433            return os.fstat(fileobj.fileno()).st_size
434        elif isinstance(fileobj, StringIO.StringIO):
435            pos = fileobj.tell()
436            fileobj.seek(0, 2)
437            length = fileobj.tell() - pos
438            fileobj.seek(pos, 0)
439            return length
440        elif isinstance(fileobj, urllib.addinfourl):
441            length = fileobj.headers.get("content-length", None)
442            return length
443    except Exception, ex:
444        pass
445
446
447def copyfileobj(src, dst, buffer=2**10, content_len=None, progress=None):
448    count = 0
449    if content_len is None:
450        content_len = guess_size(src) or sys.maxint
451    while True:
452        data = src.read(buffer)
453        dst.write(data)
454        count += len(data)
455        if progress:
456            progress(100.0 * count / content_len)
457        if not data:
458            break
459           
460           
461def wget(url, directory=".", dst_obj=None, progress=None):
462    """
463    .. todo:: Move to Orange.misc
464   
465    """
466    stream = urllib2.urlopen(url)
467    length = stream.headers.get("content-length", None)
468    if length is None:
469        length = sys.maxint
470    else:
471        length = int(length)
472   
473    basename = posixpath.basename(url)
474       
475    if dst_obj is None:
476        dst_obj = open(os.path.join(directory, basename), "wb")
477   
478    if progress == True:
479        from Orange.misc import ConsoleProgressBar
480        progress = ConsoleProgressBar("Downloading %r." % basename)
481        with finishing(progress):
482            copyfileobj(stream, dst_obj, buffer=2**10, content_len=length,
483                        progress=progress)
484    else:
485        copyfileobj(stream, dst_obj, buffer=2**10, content_len=length,
486                    progress=progress)
487   
488import obiPPI
489import orngServerFiles
490
491import obiTaxonomy
492from collections import namedtuple
493from operator import itemgetter
494from Orange.misc import lru_cache
495
496GENE_MANIA_INTERACTION_FIELDS = \
497    ["gene_a", "gene_b", "weight", "network_name",
498     "network_group", "source", "pubmed_id"]
499     
500GeneManiaInteraction = namedtuple("GeneManiaInteraction",
501                                  field_names=GENE_MANIA_INTERACTION_FIELDS
502                                 )
503
504import weakref
505class Internalizer(object):
506    """ A class that acts as the python builtin function ``intern``,
507    for as long as it is alive.
508   
509    .. note:: This is for memory optimization only, it does not affect
510        dict lookup speed.
511   
512    """
513    def __init__(self):
514        self._intern_dict = {}
515       
516    def __call__(self, obj):
517        return self._intern_dict.setdefault(obj, obj)
518   
519class GeneManiaDatabase(obiPPI.PPIDatabase):
520    DOMAIN = "PPI"
521    SERVER_FILE = "gene-mania-{taxid}.sqlite"
522   
523    TAXID2NAME = ""
524   
525    # DB schema
526    SCHEMA = """
527    table: `genes`
528        - `internal_id`: int (pk)
529        - `gene_name`: text (preferred name)
530       
531    table: `synonyms`:
532        - `internal_id: int (foreign key `genes.internal_id`)
533        - `synonym`: text
534        - `source_id`: int
535       
536    table: `source`:
537        - `source_id`: int
538        - `source_name`: text
539       
540    table: `links`:
541        - `gene_a`: int (fk `genes.internal_key`)
542        - `gene_b`: int (fk `genes.internal_key`)
543        - `network_id`: (fk `networks.network_id`)
544        - `weight`: real
545       
546    table: `networks`:
547        - `network_id`: int
548        - `network_name`: text
549        - `network_group`: text
550        - `source`: text
551        - `pubmed_id`: text
552       
553    view: `links_annotated`:
554        - `gene_name_a`
555        - `gene_name_b`
556        - `network_name`
557        - `network_group`
558        - `weight`
559       
560    """
561   
562    def __init__(self, taxid):
563        self.taxid = taxid
564       
565    @classmethod
566    def common_taxids(self):
567        return ["3702", "6239", "7227", "9606", "10090", "10116", "4932"]
568   
569    def organisms(self):
570        """ Return all organism taxids contained in this database.
571       
572        .. note:: a single taxid is returned (the one used at
573            instance initialization)   
574       
575        """
576        return [self.taxid]
577   
578    def ids(self, taxid=None):
579        """ Return all primary ids for `taxid`.
580        """
581        if taxid is None:
582            taxids = self.organisms()
583            return reduce(list.__add__, map(self.ids, taxids), [])
584       
585        con = self._db(taxid)
586        cur = con.execute("""\
587            SELECT gene_name FROM genes
588            """)
589        return map(itemgetter(0), cur)
590       
591    def synonyms(self, id):
592        """ Return a list of synonyms for primary `id`.
593        """
594        con = self._db(self.taxid)
595        cur = con.execute("""\
596            SELECT synonyms.synonym
597            FROM synonyms NATURAL LEFT JOIN genes
598            WHERE genes.gene_name=?
599            """, (id,))
600        return map(itemgetter(0), cur)
601       
602    def all_edges(self, taxid=None):
603        """ Return a list of all edges.
604        """
605        con = self._db(self.taxid)
606        cur = con.execute("""
607            SELECT links.gene_a, links.gene_b, links.weight
608            FROM links""")
609        id_to_name = self._gene_id_to_name()
610        return [(id_to_name[r[0]], id_to_name[r[1]], r[2]) \
611                for r in cur]
612       
613    def all_edges_annotated(self, taxid=None):
614        """ Return a list of all edges with all available annotations
615        """
616        con = self._db(self.taxid)
617        cur = con.execute("""\
618            SELECT links.gene_a, links.gene_b, links.weight, links.network_id
619            FROM links""")
620        gene_to_name = self._gene_id_to_name()
621        network_to_description = self._network_id_to_description()
622        res = []
623        for gene_a, gene_b, w, n_id in cur:
624            n_desc = network_to_description[n_id]
625           
626            res.append(GeneManiaInteraction(gene_to_name[gene_a],
627                            gene_to_name[gene_b], w, *n_desc))
628        return res
629       
630    def edges(self, id1):
631        """ Return all edges for primary id `id1`.
632        """       
633        con = self._db(self.taxid)
634        cur = con.execute("""\
635            SELECT genes1.gene_name, genes2.gene_name, links.weight
636            FROM genes AS genes1 
637                JOIN links
638                    ON genes1.internal_id=links.gene_a
639                JOIN genes AS genes2
640                    ON genes2.internal_id=links.gene_b
641            WHERE genes1.gene_name=?
642            """, (id1,))
643        res = cur.fetchall()
644        cur = con.execute("""\
645            SELECT genes1.gene_name, genes2.gene_name, links.weight
646            FROM genes AS genes1 
647                JOIN  links
648                    ON genes1.internal_id=links.gene_a
649                JOIN genes AS genes2
650                    ON genes2.internal_id=links.gene_b
651            WHERE genes2.gene_name=?
652            """, (id1,))
653        res += cur.fetchall()
654       
655        return res
656   
657    def edges_annotated(self, id=None):
658        """ Return a list of annotated edges for primary `id`
659        """
660        con = self._db(self.taxid)
661        cur = con.execute("""\
662            SELECT genes1.gene_name, genes2.gene_name, links.weight,
663                   networks.network_name, networks.network_group,
664                   networks.source, networks.pubmed_id
665            FROM genes AS genes1
666                JOIN  links
667                    ON genes1.internal_id=links.gene_a
668                JOIN genes AS genes2
669                    ON genes2.internal_id=links.gene_b
670                NATURAL JOIN networks
671            WHERE genes1.gene_name=?
672            """, (id,))
673        res = cur.fetchall()
674        cur = con.execute("""\
675            SELECT genes1.gene_name, genes2.gene_name, links.weight,
676                   networks.network_name, networks.network_group,
677                   networks.source, networks.pubmed_id
678            FROM genes AS genes1
679                JOIN links
680                    ON genes1.internal_id=links.gene_a
681                JOIN genes AS genes2
682                    ON genes2.internal_id=links.gene_b
683                NATURAL JOIN networks
684            WHERE genes2.gene_name=?
685            """, (id,))
686        res += cur.fetchall()
687        return [GeneManiaInteraction(*r) for r in res]
688   
689    def search_id(self, name, taxid=None):
690        """ Search the database for gene name. Return a list of matching
691        primary ids. Use `taxid` to limit the results to a single organism.
692       
693        """
694        con = self._db(self.taxid)
695        cur = con.execute("""\
696            SELECT genes.gene_name
697            FROM genes NATURAL JOIN synonyms
698            WHERE synonyms.synonym=?
699            """, (name,))
700        return map(itemgetter(0), cur)
701       
702    def _db(self, taxid=None):
703        """ Return an open sqlite3.Connection object. 
704        """
705        taxid = taxid or self.taxid
706        filename = orngServerFiles.localpath_download("PPI",
707                            self.SERVER_FILE.format(taxid=taxid))
708        if not os.path.exists(filename):
709            raise ValueError("Database is missing.")
710       
711        return sqlite3.connect(filename)
712   
713    @lru_cache(maxsize=1)
714    def _gene_id_to_name(self):
715        """ Return a dictionary mapping internal gene ids to
716        primary gene identifiers.
717       
718        """
719        con = self._db(self.taxid)
720        cur = con.execute("SELECT * FROM genes")
721        return dict(cur)
722   
723    @lru_cache(maxsize=1)
724    def _network_id_to_description(self):
725        """ Return a dictionary mapping internal network ids
726        to (name, group, source, pubmed id).
727         
728        """
729        con = self._db(self.taxid)
730        cur = con.execute("SELECT * FROM networks")
731        return dict((t[0], t[1:]) for t in cur)
732   
733    #####################################
734    # Data download and DB initialization
735    #####################################
736     
737    @classmethod
738    def download_data(cls, taxid=None, progress_callback=None):
739        """ Download the data for ``taxid`` from the GeneMANIA
740        website and initialize the local database.
741       
742        """
743        import tarfile
744       
745        baseurl = "http://genemania.org/data/current/"
746        directory = orngServerFiles.localpath("PPI")
747        if taxid is None:
748            taxid = cls.common_taxids()
749       
750        if isinstance(taxid, (list, tuple)):
751            taxids = taxid
752        else:
753            taxids = [taxid]
754        for taxid in taxids:
755            name = obiTaxonomy.name(taxid)
756            name = name.replace(" ", "_")
757           
758            if progress_callback is None:
759                progress = True #orngServerFiles.ConsoleProgressBar("Downloading %r." % filename)
760            else:
761                progress = progress_callback
762           
763            filename = name + ".tgz"
764            url = baseurl + "networks/" + filename   
765            wget(url, directory=directory, progress=progress)
766           
767            tgz_filename = os.path.join(directory, filename)   
768            tgz = tarfile.open(tgz_filename)
769            tgz.extractall(directory)
770           
771            filename = name + ".COMBINED.tgz"
772            url = baseurl + "precombined/" + filename
773            wget(url, directory=directory, progress=progress)
774           
775            tgz_filename = os.path.join(directory, filename)
776            tgz = tarfile.open(tgz_filename)
777            tgz.extractall(directory)
778       
779            cls.init_db([taxid])
780       
781    @classmethod
782    def init_db(cls, taxid=None):
783        """ Init the local data base.
784        """
785        from functools import partial
786        directory = orngServerFiles.localpath("PPI")
787        pjoin = partial(os.path.join, directory)
788        if taxid is None:
789            taxid = cls.common_taxids()
790           
791        if isinstance(taxid, (list, tuple)):
792            for tid in taxid:
793                cls.init_db(tid)
794            return
795               
796        if not isinstance(taxid, basestring):
797            raise ValueError("wrong taxid")
798           
799#        taxid = taxids
800        name = obiTaxonomy.name(taxid).replace(" ", "_")
801        networks = csv.reader(open(pjoin(name, "networks.txt")), delimiter="\t")
802        networks.next() # Header
803        networks = list(networks)
804       
805        database = pjoin(cls.SERVER_FILE.format(taxid=taxid))
806        with sqlite3.connect(database) as con:
807            con.execute("""DROP TABLE IF EXISTS genes""")
808            con.execute("""DROP TABLE IF EXISTS synonyms""")
809            con.execute("""DROP TABLE IF EXISTS source""")
810            con.execute("""DROP TABLE IF EXISTS links""")
811            con.execute("""DROP TABLE IF EXISTS networks""")
812           
813            con.execute("""DROP INDEX IF EXISTS genes_index""")
814            con.execute("""DROP INDEX IF EXISTS links_index_a""")
815            con.execute("""DROP INDEX IF EXISTS links_index_b""")
816           
817            con.execute("""\
818                CREATE TABLE networks
819                    (network_id INTEGER PRIMARY KEY ASC AUTOINCREMENT,
820                     network_name TEXT,
821                     network_group TEXT,
822                     source TEXT,
823                     pubmed_id TEXT
824                    )""")
825           
826            con.executemany("""\
827                INSERT INTO networks
828                VALUES (?, ?, ?, ?, ?)""", [(i, r[2], r[1], r[3], r[4]) \
829                                        for i, r in enumerate(networks)])
830           
831            con.execute("""\
832                CREATE TABLE genes
833                    (internal_id INTEGER PRIMARY KEY ASC AUTOINCREMENT,
834                     gene_name TEXT
835                    )""")
836           
837            identifiers = csv.reader(open(pjoin(name, "identifier_mappings.txt"), "rb"),
838                                    delimiter="\t")
839            identifiers.next() # skip header
840            identifiers = list(identifiers)
841            genes = sorted(set(r[0] for r in identifiers))
842            sources = sorted(set(r[2] for r in identifiers))
843           
844            con.executemany("""\
845                INSERT INTO genes
846                VALUES (?, ?)""", enumerate(genes))
847           
848            con.execute("""\
849            CREATE TABLE source
850                (source_id INTEGER PRIMARY KEY ASC AUTOINCREMENT,
851                 source_name TEXT
852                )""")
853           
854            con.executemany("""\
855                INSERT INTO source
856                VALUES (?, ?)""", enumerate(sources))
857           
858            con.execute("""\
859                CREATE TABLE synonyms
860                    (internal_id INTEGER REFERENCES genes (internal_id),
861                     synonym TEXT,
862                     source_id INT REFERENCES source (source_id)
863                    )""")
864           
865            gene_to_id = dict((g, i) for i, g in enumerate(genes))
866            source_to_id = dict((s, i) for i, s in enumerate(sources))
867            con.executemany("""\
868                INSERT INTO synonyms
869                VALUES (?, ?, ?)""", [(gene_to_id[r[0]], r[1], source_to_id[r[2]])\
870                                       for r in identifiers])
871           
872            con.execute("""\
873                CREATE TABLE links
874                    (gene_a INTEGER REFERENCES genes (internal_id),
875                     gene_b INTEGER REFERENCES genes (internal_id),
876                     network_id INTEGER REFERENCES networks (network_id),
877                     weight REAL
878                     -- , PRIMARY KEY (gene_a, gene_b, network_id)
879                    )""")
880           
881            for i, (filename, group, _, _, _) in enumerate(networks):
882                nf  = open(pjoin(name, filename), "rb")
883                interactions = csv.reader(nf, delimiter="\t")
884                interactions.next() # skip header
885                con.executemany("""\
886                    INSERT INTO links
887                    VALUES (?, ?, ?, ?)""",
888                    [(gene_to_id[r[0]], gene_to_id[r[1]], i, float(r[2])) \
889                     for r in interactions]
890                )
891               
892            # Add special combined network entry
893            combined_id = len(networks)
894            con.execute("""\
895                INSERT INTO networks
896                VALUES (?, ?, ?, ?, ?)""", 
897                (combined_id, "BP_COMBINING", "COMBINED", "GeneMANIA", ""))
898           
899            # Add the combined network links.
900            combined = open(pjoin(name + ".COMBINED", "COMBINED.DEFAULT_NETWORKS.BP_COMBINING.txt"), "rb")
901            combined = csv.reader(combined, delimiter="\t")
902            combined.next()
903            con.executemany("""\
904                INSERT INTO links
905                VALUES (?, ?, ?, ?)""",
906                    ((gene_to_id[r[0]], gene_to_id[r[1]], combined_id, float(r[2])) \
907                     for r in combined))
908           
909           
910            con.execute("""
911                CREATE VIEW IF NOT EXISTS links_annotated
912                AS SELECT genes1.gene_name AS gene_name_a,
913                          genes2.gene_name AS gene_name_b,
914                          links.weight,
915                          networks.network_name,
916                          networks.network_group,
917                          networks.source,
918                          networks.pubmed_id
919                   FROM  genes AS genes1
920                        JOIN links
921                              ON genes1.internal_id=links.gene_a
922                        JOIN genes AS genes2
923                              ON links.gene_b=genes2.internal_id
924                        JOIN networks
925                              ON links.network_id=networks.network_id
926                    """)
927           
928           
929            con.execute("""\
930                CREATE INDEX IF NOT EXISTS genes_index ON genes (gene_name)
931                """)
932            con.execute("""\
933                CREATE INDEX IF NOT EXISTS links_index_a ON links (gene_a)
934                """)
935            con.execute("""\
936                 CREATE INDEX IF NOT EXISTS links_index_b ON links (gene_b)
937                """)
938       
939           
940if __name__ == "__main__":
941    retrieve("9606", [ 'MRE11A', 'RAD51', 'MLH1', 'MSH2', 'DMC1', 'RAD51AP1', 'RAD50', 'MSH6', 'XRCC3', 'PCNA', 'XRCC2' ])
Note: See TracBrowser for help on using the repository browser.