allpy: afed1fd8920c allpy/graph.py

allpy

view allpy/graph.py @ 1091:afed1fd8920c

Added backreferences to `Seqeunce`s from `Monomer`s (closes #49) WARNING! Please note that `Sequence` API almost changed entirely! WARNING! This commit immediately obsoletes classmethods `Monomer.from_code*`, `Monomer.from_name` and `Sequence.from_monomers`. Turns out, python can not pickle sets/dicts which have keys, which inderecly reference the set/dict itself: http://bugs.python.org/issue9269 -- which is excatly what we have in abundance after this change. To allow pickling added `__getstate__` to `Monomer` to return all attributes, except `sequence` and `__setstate__` to `Sequence`, which runs through all monomers and returns the `sequence` attribute back to where it belongs. WARNING! This MAY result in unexpected behaviour in some cases. (Which should be rare enough).

author	Daniil Alexeyevsky <dendik@kodomo.fbb.msu.ru>
date	Sat, 02 Jun 2012 19:33:42 +0400
parents	53690e470eff
children

line source

1 from datetime import datetime, timedelta

2 from copy import copy

4 import config

6 class TimeoutError(Exception):

7 pass

9 class Graph(dict):

10 """ Undirected weighted graph

12 graph[vertex1][vertex2] = weight

13 * vertex1 != vertex2 (no such edges)

14 * weight is float in (0, 1] or 1 (if all edges are equal)

15 * symmetrical

17 Features:

18 * this graph cannot contain vertices without any edges

20 >>> config.debug = True

21 >>> g = Graph()

22 >>> g.set_edge(1, 2)

23 >>> g.fast_cliques()

24 Fast algorithm started

25 [frozenset([1, 2])]

26 >>> g.set_edge(2, 3)

27 >>> g.set_edge(1, 3)

28 >>> g.fast_cliques()

29 Fast algorithm started

30 [frozenset([1, 2, 3])]

31 >>> g.bron_kerbosch()

32 Bron and Kerbosch algorithm started

33 [frozenset([1, 2, 3])]

34 >>> g.cliques()

35 Bron and Kerbosch algorithm started

36 [frozenset([1, 2, 3])]

38 """

40 def copy(self):

41 cls = self.__class__

42 graph = cls()

43 for x, neighbours_dict in self.items():

44 for y, weight in neighbours_dict.items():

45 graph.set_edge(x, y, weight)

46 return graph

48 def set_edge(self, v1, v2, weight=1):

49 assert v1 is not v2

50 self.setdefault(v1, {})

51 self[v1][v2] = weight

52 self.setdefault(v2, {})

53 self[v2][v1] = weight

55 def drop_vertex(self, vertex):

56 """ Remove vertex and all involved edges """

57 for v in self[vertex]:

58 del self[v][vertex]

59 if not self[v]:

60 del self[v]

61 del self[vertex]

63 def vertex_weight(self, vertex):

64 """ Returns sum of weights of all connections of this vertex """

65 return sum(self[vertex].values())

67 def vertex_index(self, vertex):

68 return len(self[vertex])

70 def add_vertex(self, vertex, parent_graph):

71 """ Add vertex and corresponding edges from parent_graph

73 Added edges should be contained in self graph

74 (takes care of hanging edges)

75 """

76 for v in parent_graph[vertex]:

77 if v in self:

78 self.set_edge(vertex, v, parent_graph[vertex][v])

80 def bron_kerbosh(self, timeout=-1, minsize=1):

81 """ Deprecated. Use bron_kerbosch() instead """

82 return self.bron_kerbosch(timeout, minsize)

84 def bron_kerbosch(self, timeout=-1, minsize=1):

85 """ Bron and Kerbosch algorithm implementation

87 returns list of cliques

88 clique is frozenset

89 if timeout=-1, it means infinity

90 if timeout has happened, raises TimeoutError

91 """

92 if timeout == 0:

93 raise TimeoutError()

94 if config.debug:

95 print 'Bron and Kerbosch algorithm started'

96 cliques = []

98 class ExtendCall(object):

99 def __init__(call, candidates, used, compsub):

100 call.candidates = candidates

101 call.used = used

102 call.compsub = compsub

103

104 def get_v_from_candidates(call):

105 """ return some v from candidates """

106 return iter(call.candidates).next()

107

108 def used_candidates(call):

109 """ if v from used connected with all from candidates """

110 for u in call.used:

111 for c in call.candidates:

112 if not self.get(u,{}).get(c,None):

113 break

114 else:

115 return True

116 return False

117

118 def __call__(call):

119 """ process call and return list of recursive calls """

120 while call.candidates and not call.used_candidates():

121 v = call.get_v_from_candidates()

122 call.compsub.add(v)

123 new_candidates = set([i for i in call.candidates

124 if i is not v and self.get(i,{}).get(v,None)])

125 new_used = set([i for i in call.used

126 if i is not v and self.get(i,{}).get(v,None)])

127 if (not new_candidates) and (not new_used):

128 if len(call.compsub) >= minsize:

129 cliques.append(copy(frozenset(call.compsub)))

130 else:

131 yield ExtendCall(new_candidates,

132 new_used, copy(call.compsub))

133 call.compsub.remove(v)

134 call.candidates.remove(v)

135 call.used.add(v)

136

137 stop_time = datetime.now() + timedelta(0, timeout)

138 stack = [ExtendCall(copy(set(self)), set(), set())()]

139 while stack:

140 if timeout != -1 and datetime.now() > stop_time:

141 raise TimeoutError()

142 top = stack[-1]

143 try:

144 call = top.next()

145 stack.append(call())

146 except StopIteration:

147 stack.pop()

148 return cliques

149

150 def drop_bad_vertices(self):

151 """ drop vertices untill self becomes a clique """

152 while len(self) >= 1:

153 indexes = dict((v, self.vertex_index(v)) for v in self)

154 min_index = min(indexes.values())

155 bad_vertices = [v for (v, i) in indexes.items()

156 if i == min_index]

157 if len(bad_vertices) == len(self):

158 break # clique

159 weights = dict([(self.vertex_weight(v),v) for v in bad_vertices])

160 min_weight = min(weights.keys())

161 vertex = weights[min_weight]

162 self.drop_vertex(vertex)

163

164 def expand_to_max_clique(self, parent_graph):

165 """ add vertices untill self becomes a nonextendible clique """

166 candidates = set()

167 for vertex in set(parent_graph) - set(self):

168 neighbours = set(parent_graph[vertex])

169 if set(self).issubset(neighbours):

170 candidates.add(vertex)

171 self.weights = {}

172 for v in candidates:

173 weight = parent_graph.vertex_weight(v)

174 self.weights[v] = weight

175 while True:

176 if not candidates:

177 break

178 _rev = dict((v, k) for (k, v) in self.weights.items())

179 vertex = _rev[max(_rev.keys())]

180 self.add_vertex(vertex, parent_graph)

181 for v in copy(candidates):

182 if v not in parent_graph[vertex]:

183 candidates.discard(v)

184 del self.weights[v]

185

186 def fast_cliques(self, minsize=1):

187 """ returns list of cliques

188

189 clique is frozenset

190 """

191 if config.debug:

192 print 'Fast algorithm started'

193 cliques = []

194 while True:

195 graph = self.copy()

196 for clique in cliques:

197 for v in clique:

198 if v in graph:

199 graph.drop_vertex(v)

200 if not graph:

201 break

202 graph.drop_bad_vertices()

203 graph.expand_to_max_clique(self)

204 clique = copy(frozenset(graph))

205 if len(clique) >= minsize:

206 cliques.append(clique)

207 else:

208 break

209 return cliques

210

211 def connected_components(self):

212 """ return list of connected components, remove single vertices """

213 used = set()

214 graphs = []

215 candidates = set(self)

216 while candidates:

217 Graph = self.__class__

218 graph = Graph()

219 queue = set()

220 queue.add(candidates.pop())

221 while queue:

222 v = queue.pop()

223 candidates.discard(v)

224 for v1 in self.get(v,{}):

225 if v1 in candidates:

226 graph.set_edge(v, v1, self[v][v1])

227 queue.add(v1)

228 graphs.append(graph)

229 return graphs

230

231 def cliques(self, timeout=-1, minsize=1):

232 """ returns length-sorted list of cliques

233

234 clique is frozenset

235

236 try to execute bron_kerbosch

237 if it raises TimeoutError, executes fast_cliques

238 """

239 cliques = []

240 for graph in self.connected_components():

241 if len(graph) >= minsize:

242 for v, neighbours_dict in graph.items():

243 if len(neighbours_dict) < minsize - 1 and v in graph:

244 graph.drop_vertex(v)

245 try:

246 cliques += graph.bron_kerbosch(timeout, minsize)

247 except TimeoutError:

248 cliques += graph.fast_cliques(minsize)

249 cliques.sort(key=lambda clique: len(clique), reverse=True)

250 return cliques

251

252 if __name__ == "__main__":

253 import doctest

254 doctest.testmod()

255