|
1
|
|
|
#!/usr/bin/env python |
|
2
|
|
|
"""Ancestors/Descendants.""" |
|
3
|
|
|
|
|
4
|
|
|
import os |
|
5
|
|
|
import sys |
|
6
|
|
|
import timeit |
|
7
|
|
|
import numpy as np |
|
8
|
|
|
from numpy.random import shuffle |
|
9
|
|
|
from scipy import stats |
|
10
|
|
|
|
|
11
|
|
|
from goatools.base import download_go_basic_obo |
|
12
|
|
|
from goatools.obo_parser import GODag |
|
13
|
|
|
from goatools.test_data.godag_timed import prt_hms |
|
14
|
|
|
from goatools.gosubdag.gosubdag import GoSubDag |
|
15
|
|
|
|
|
16
|
|
|
|
|
17
|
|
|
def test_go_pools(): |
|
18
|
|
|
"""Print a comparison of GO terms from different species in two different comparisons.""" |
|
19
|
|
|
objr = _Run() |
|
20
|
|
|
# Check that subset GoSubDags have the same ancestors/descendants as Full GoSubDag |
|
21
|
|
|
for qty in np.random.randint(10, 100, size=10): |
|
22
|
|
|
print("") |
|
23
|
|
|
goids = objr.get_goids_rand(qty) |
|
24
|
|
|
# No relationships loaded; GoSubDag subset equivalent to Full subset? |
|
25
|
|
|
gosubdag_r0 = objr.get_gosubdag_r0(goids) |
|
26
|
|
|
for goid in gosubdag_r0.go2obj: |
|
27
|
|
|
r0_u = gosubdag_r0.rcntobj.go2parents[goid] |
|
28
|
|
|
r0_d = gosubdag_r0.rcntobj.go2descendants[goid] |
|
29
|
|
|
assert r0_u == objr.gosubdag_r0.rcntobj.go2parents[goid] |
|
30
|
|
|
assert r0_d == objr.gosubdag_r0.rcntobj.go2descendants[goid] |
|
31
|
|
|
# All relationships loaded; GoSubDag(r0) vs. GoSubDag(r1) |
|
32
|
|
|
gosubdag_r1 = objr.get_gosubdag_r1(goids) |
|
33
|
|
|
assert gosubdag_r0.go_sources == gosubdag_r1.go_sources |
|
34
|
|
|
assert set(gosubdag_r0.go2obj).issubset(gosubdag_r1.go2obj) |
|
35
|
|
|
cnts = {'r0_u':[], 'r1_u':[], 'r0_d':[], 'r1_d':[]} |
|
36
|
|
|
for goid in gosubdag_r0.go2obj: |
|
37
|
|
|
r0_u = gosubdag_r0.rcntobj.go2parents[goid] |
|
38
|
|
|
r0_d = gosubdag_r0.rcntobj.go2descendants[goid] |
|
39
|
|
|
r1_u = gosubdag_r1.rcntobj.go2parents[goid] |
|
40
|
|
|
r1_d = gosubdag_r1.rcntobj.go2descendants[goid] |
|
41
|
|
|
assert r0_d.issubset(r1_d), "R1({}) R0({})".format(len(r1_d), len(r0_d)) |
|
42
|
|
|
assert r0_u.issubset(r1_u), "R1({}) R0({})".format(len(r1_u), len(r0_u)) |
|
43
|
|
|
cnts['r0_u'].append(len(r0_u)) |
|
44
|
|
|
cnts['r1_u'].append(len(r1_u)) |
|
45
|
|
|
cnts['r0_d'].append(len(r0_d)) |
|
46
|
|
|
cnts['r1_d'].append(len(r1_d)) |
|
47
|
|
|
objr.prt_cnts(cnts) |
|
48
|
|
|
|
|
49
|
|
|
|
|
50
|
|
|
class _Run(object): |
|
51
|
|
|
"""Group entire go-basic.obo""" |
|
52
|
|
|
|
|
53
|
|
|
obo = os.path.join(os.path.dirname(os.path.abspath(__file__)), "../../go-basic.obo") |
|
54
|
|
|
|
|
55
|
|
|
def __init__(self): |
|
56
|
|
|
download_go_basic_obo(self.obo, sys.stdout, loading_bar=None) |
|
57
|
|
|
self.godag_r0 = GODag("go-basic.obo") |
|
58
|
|
|
self.godag_r1 = GODag("go-basic.obo", optional_attrs=set(['relationship'])) |
|
59
|
|
|
self.goids = list(set(o.id for o in self.godag_r0.values())) |
|
60
|
|
|
# GoSubDag (plain) |
|
61
|
|
|
tic = timeit.default_timer() |
|
62
|
|
|
self.gosubdag_r0 = GoSubDag(self.goids, self.godag_r0, prt=None) |
|
63
|
|
|
prt_hms(tic, "GoSubDag r0 {N:4} GOs {S:3} srcs".format( |
|
64
|
|
|
N=len(self.gosubdag_r0.go2obj), S=len(self.gosubdag_r0.go_sources))) |
|
65
|
|
|
# GoSubDag with relationships |
|
66
|
|
|
self.gosubdag_r1 = GoSubDag(self.goids, self.godag_r1, prt=None, relationships=True) |
|
67
|
|
|
prt_hms(tic, "GoSubDag r1 {N:4} GOs {S:3} srcs".format( |
|
68
|
|
|
N=len(self.gosubdag_r1.go2obj), S=len(self.gosubdag_r1.go_sources))) |
|
69
|
|
|
|
|
70
|
|
|
def prt_cnts(self, cnts): |
|
71
|
|
|
"""Compare ancestor/descendant counts with relatives=False/True.""" |
|
72
|
|
|
k2v = {k:self.str_stats(v) for k, v in cnts.items()} |
|
73
|
|
|
print(k2v) |
|
74
|
|
|
|
|
75
|
|
|
@staticmethod |
|
76
|
|
|
def str_stats(vals): |
|
77
|
|
|
"""Print statistics on values.""" |
|
78
|
|
|
ntd = stats.describe(vals) |
|
79
|
|
|
std = int(round(np.sqrt(ntd.variance))) |
|
80
|
|
|
return "({m} {M}) STD={STD:,}".format(m=ntd.minmax[0], M=ntd.minmax[1], STD=std) |
|
81
|
|
|
|
|
82
|
|
|
def get_gosubdag_r0(self, goids): |
|
83
|
|
|
"""Return a GoSubDag with N randomly chosen GO sources.""" |
|
84
|
|
|
tic = timeit.default_timer() |
|
85
|
|
|
gosubdag = GoSubDag(goids, self.godag_r0, relationships=None, |
|
86
|
|
|
#rcntobj=self.gosubdag_r0.rcntobj, |
|
87
|
|
|
prt=None) |
|
88
|
|
|
prt_hms(tic, "GoSubDag r0 {N:4} GOs {S:3} srcs".format( |
|
89
|
|
|
N=len(gosubdag.go2obj), S=len(gosubdag.go_sources))) |
|
90
|
|
|
return gosubdag |
|
91
|
|
|
|
|
92
|
|
|
def get_gosubdag_r1(self, goids): |
|
93
|
|
|
"""Return a GoSubDag with N randomly chosen GO sources.""" |
|
94
|
|
|
tic = timeit.default_timer() |
|
95
|
|
|
gosubdag = GoSubDag(goids, self.godag_r1, relationships=True, |
|
96
|
|
|
#rcntobj=self.gosubdag_r1.rcntobj, |
|
97
|
|
|
prt=None) |
|
98
|
|
|
prt_hms(tic, "GoSubDag r1 {N:4} GOs {S:3} srcs".format( |
|
99
|
|
|
N=len(gosubdag.go2obj), S=len(gosubdag.go_sources))) |
|
100
|
|
|
return gosubdag |
|
101
|
|
|
|
|
102
|
|
|
def get_goids_rand(self, qty): |
|
103
|
|
|
"""Return N randomly chosen GO IDs.""" |
|
104
|
|
|
shuffle(self.goids) |
|
105
|
|
|
return self.goids[:qty] |
|
106
|
|
|
|
|
107
|
|
|
|
|
108
|
|
|
if __name__ == '__main__': |
|
109
|
|
|
test_go_pools() |
|
110
|
|
|
|
tanghaibao /
goatools
