-
Notifications
You must be signed in to change notification settings - Fork 0
/
helpers.py
439 lines (405 loc) · 16.9 KB
/
helpers.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
import random
import os
import pickle
import itertools
import pymnet as pn
import networkx
import matplotlib.pyplot as plt
import matplotlib
import numpy as np
import math
import collections
class persistent(object):
def __init__(self,function):
self.function=function
def __call__(self,*args,**kwargs):
if 'persistent_file' in kwargs:
persistent_file=kwargs.pop('persistent_file')
if os.path.exists(persistent_file): #read from the file
with open(persistent_file,'rb') as f:
result = pickle.load(f)
else: #create it, write it
result = self.function(*args,**kwargs)
with open(persistent_file,'wb') as f:
pickle.dump(result,f)
return result
else:
return self.function(*args,**kwargs)
def er_multilayer_any_aspects(l=[10,4,4],p=0.05):
# generate all nodelayers and then add random edges between them
M = pn.MultilayerNetwork(aspects=len(l)-1,directed=False,fullyInterconnected=True)
for ii,layers_in_aspect in enumerate(l):
for jj in range(layers_in_aspect):
M.add_layer(layer=jj,aspect=ii)
nls = list(M.iter_node_layers())
for ii in range(len(nls)):
for jj in range(ii+1,len(nls)):
if random.random() < p:
M[nls[ii]][nls[jj]] = 1
return M
def er_multilayer_any_aspects_deg_1_or_greater(l=[5,5,5],p=0.05):
# only add nodelayers which have at least one edge, net not fully interconnected
M = pn.MultilayerNetwork(aspects=len(l)-1,directed=False,fullyInterconnected=False)
possible_nls = list(itertools.product(*(list(range(elem_number)) for elem_number in l)))
for ii in range(len(possible_nls)):
for jj in range(ii+1,len(possible_nls)):
if random.random() < p:
M[possible_nls[ii]][possible_nls[jj]] = 1
return M
def geo_multilayer_any_aspects(l,edges_in_layers,edges_between_layers):
# TODO: currently does not work with 1 aspect, since layers are ints instead of tuples
# tempfix: use 1 as number of elementary layers in second aspect
#
# l : number of elementary layers in each aspect
# edges_in_layers : int with approximate number of edges within layers
# edges_between_layers : int with approximate number of edges between each pair of layers
M = pn.MultilayerNetwork(aspects=len(l)-1,directed=False,fullyInterconnected=True)
n = l[0]
pos = None
for ii,layers_in_aspect in enumerate(l):
for jj in range(layers_in_aspect):
M.add_layer(layer=jj,aspect=ii)
layers = list(M.iter_layers()) # for getting all pairs later
for layer in layers:
netX,pos = get_single_geo_instance(n,edges_in_layers,pos)
for node in netX.nodes:
M.add_node(node,layer=layer)
for e in netX.edges:
# add edges between nodelayers inside layer
node1 = e[0]
node2 = e[1]
medge = [node1, node2, *layer]
M[tuple(medge)] = 1
for layer_pair in itertools.combinations(layers, 2):
netX,pos = get_single_geo_instance(n,edges_between_layers,pos)
for e in netX.edges:
# add edges between nodelayers between layers
iledge = []
node1 = e[0]
node2 = e[1]
iledge.append(node1)
iledge.append(node2)
for ii in range(len(layer_pair[0])):
iledge.append(layer_pair[0][ii])
iledge.append(layer_pair[1][ii])
M[tuple(iledge)] = 1
return M
def geo_er_mixed_multilayer_any_aspects(l,edges_in_layers,edges_between_layers,p_er_intra,p_er_inter):
# geo model with er model on top (each layer and each pair of layers)
#
# TODO: currently does not work with 1 aspect, since layers are ints instead of tuples
# tempfix: use 1 as number of elementary layers in second aspect
#
# l : number of elementary layers in each aspect
# edges_in_layers : int with approximate number of edges within layers
# edges_between_layers : int with approximate number of edges between each pair of layers
M = pn.MultilayerNetwork(aspects=len(l)-1,directed=False,fullyInterconnected=True)
n = l[0]
pos = None
for ii,layers_in_aspect in enumerate(l):
for jj in range(layers_in_aspect):
M.add_layer(layer=jj,aspect=ii)
layers = list(M.iter_layers()) # for getting all pairs later
for layer in layers:
netX,pos = get_single_geo_instance(n,edges_in_layers,pos)
for node in netX.nodes:
M.add_node(node,layer=layer)
for e in netX.edges:
# add edges between nodelayers inside layer
node1 = e[0]
node2 = e[1]
medge = [node1, node2, *layer]
M[tuple(medge)] = 1
# add additional er edges
all_nodes = list(netX.nodes)
for node_pair in itertools.combinations(all_nodes,2):
if random.random() < p_er_intra:
meredge = [node_pair[0], node_pair[1], *layer]
M[tuple(meredge)] = 1
for layer_pair in itertools.combinations(layers, 2):
netX,pos = get_single_geo_instance(n,edges_between_layers,pos)
for e in netX.edges:
# add edges between nodelayers between layers
iledge = []
node1 = e[0]
node2 = e[1]
iledge.append(node1)
iledge.append(node2)
for ii in range(len(layer_pair[0])):
iledge.append(layer_pair[0][ii])
iledge.append(layer_pair[1][ii])
M[tuple(iledge)] = 1
# add additional er edges
# consider all possible edge pairs, i.e.
# bipartite er where number of edges is n^2
# also including self-edges
all_nodes = list(netX.nodes)
for node1 in all_nodes:
for node2 in all_nodes:
if random.random() < p_er_inter:
ileredge = []
ileredge.append(node1)
ileredge.append(node2)
for ii in range(len(layer_pair[0])):
ileredge.append(layer_pair[0][ii])
ileredge.append(layer_pair[1][ii])
M[tuple(ileredge)] = 1
return M
def geo_er_mixed_multilayer_any_aspects_mean_deg_parametrization(l,d_geo_intra,d_geo_inter,d_er_intra,d_er_inter):
# geo model with er model on top (each layer and each pair of layers)
#
# TODO: currently does not work with 1 aspect, since layers are ints instead of tuples
# tempfix: use 1 as number of elementary layers in second aspect
#
# l : number of elementary layers in each aspect
# edges_in_layers : int with approximate number of edges within layers
# edges_between_layers : int with approximate number of edges between each pair of layers
#
# transform to other parametrization:
# geo
edges_in_layers = int(np.floor((d_geo_intra*l[0])/2.0))
edges_between_layers = int(np.floor((l[0]*d_geo_inter)/(np.product(l[1:])-1)))
# er
# testing with 0
p_er_intra = d_er_intra/(float(l[0]-1))
p_er_inter = d_er_inter/(float(np.product(l[1:])-1)*l[0])
#
M = pn.MultilayerNetwork(aspects=len(l)-1,directed=False,fullyInterconnected=True)
n = l[0]
pos = None
for ii,layers_in_aspect in enumerate(l):
for jj in range(layers_in_aspect):
M.add_layer(layer=jj,aspect=ii)
layers = list(M.iter_layers()) # for getting all pairs later
for layer in layers:
netX,pos = get_single_geo_instance(n,edges_in_layers,pos)
for node in netX.nodes:
M.add_node(node,layer=layer)
for e in netX.edges:
# add edges between nodelayers inside layer
node1 = e[0]
node2 = e[1]
medge = [node1, node2, *layer]
M[tuple(medge)] = 1
# add additional er edges
all_nodes = list(netX.nodes)
for node_pair in itertools.combinations(all_nodes,2):
if random.random() < p_er_intra:
meredge = [node_pair[0], node_pair[1], *layer]
M[tuple(meredge)] = 1
for layer_pair in itertools.combinations(layers, 2):
netX,pos = get_single_geo_instance(n,edges_between_layers,pos)
for e in netX.edges:
# add edges between nodelayers between layers
iledge = []
node1 = e[0]
node2 = e[1]
iledge.append(node1)
iledge.append(node2)
for ii in range(len(layer_pair[0])):
iledge.append(layer_pair[0][ii])
iledge.append(layer_pair[1][ii])
M[tuple(iledge)] = 1
# add additional er edges
# consider all possible edge pairs, i.e.
# bipartite er where number of edges is n^2
# also including self-edges
all_nodes = list(netX.nodes)
for node1 in all_nodes:
for node2 in all_nodes:
if random.random() < p_er_inter:
ileredge = []
ileredge.append(node1)
ileredge.append(node2)
for ii in range(len(layer_pair[0])):
ileredge.append(layer_pair[0][ii])
ileredge.append(layer_pair[1][ii])
M[tuple(ileredge)] = 1
return M
def get_single_geo_instance(n,n_edges,pos):
# slightly overgenerates the number of edges
r = math.sqrt(2.2 * float(n_edges) / ((n - 1.0) * n) / math.pi)
netX = networkx.soft_random_geometric_graph(n, r, pos=pos)
pos = networkx.get_node_attributes(netX, 'pos')
return netX, pos
def get_edge_information(M):
edge_info = dict()
for layer in M.iter_layers():
if isinstance(layer,int):
edge_info[(layer,)] = 0
else:
edge_info[layer] = 0
for layer_pair in itertools.combinations(M.iter_layers(),2):
unordered_pair = frozenset(layer_pair)
edge_info[unordered_pair] = 0
for e in M.edges:
layer_part = e[2:-1]
layer1 = tuple(layer_part[0::2])
layer2 = tuple(layer_part[1::2])
if layer1 == layer2:
edge_info[layer1] = edge_info[layer1] + 1
else:
edge_info[frozenset((layer1,layer2))] = edge_info[frozenset((layer1,layer2))] + 1
return edge_info
def load_edgelist(fname):
M = pn.MultiplexNetwork(couplings='categorical',directed=False,fullyInterconnected=False)
with open(fname,'r') as f:
for line in f:
layer,node1,node2,weight = line.split()
if node1 != node2:
M[node1,layer][node2,layer] = 1
return M
def save_edgelist_cpp_format(M,savename,include_isolated_nls=False):
# pymnet edge notation: a0, b0, a1, b1, ..., ak, bk
# where aN is elem layer of node-layer a for aspect N
# saving notation: a0 a1 a2 ... b0 b1 b2 ...
# if include_isolated_nls:
# all edges are written first and lone nodelayers after them (one per line)
all_nodelayers_with_edges = set()
with open(savename,'w') as f:
for e in M.edges:
nl1 = e[0:(M.aspects+1)*2:2]
nl2 = e[1:(M.aspects+1)*2:2]
all_nodelayers_with_edges.add(nl1)
all_nodelayers_with_edges.add(nl2)
f.write(' '.join([str(x) for x in nl1+nl2])+'\n')
if include_isolated_nls:
for nl in M.iter_node_layers():
if nl not in all_nodelayers_with_edges:
f.write(' '.join([str(y) for y in nl])+'\n')
def load_edgelist_in_cpp_format(savename):
# only works when first line is an edge (has two nodelayers)
aspect_flag = False
with open(savename,'r') as f:
for line in f:
line_list = [int(x) for x in line.strip().split()]
if not aspect_flag:
n_aspects = int(len(line_list)/2)-1
M = pn.MultilayerNetwork(aspects=n_aspects,fullyInterconnected=False)
aspect_flag = True
nl1 = tuple(line_list[:n_aspects+1])
nl2 = tuple(line_list[n_aspects+1:])
if nl2:
M[nl1][nl2] = 1
else:
M.add_node(nl1[0],layer=nl1[1:])
return M
def heatmap(data, row_labels, col_labels, xlabel, ylabel, title='', ax=None,
cbar_kw={}, cbarlabel="", **kwargs):
"""
Create a heatmap from a numpy array and two lists of labels.
Parameters
----------
data
A 2D numpy array of shape (M, N).
row_labels
A list or array of length M with the labels for the rows.
col_labels
A list or array of length N with the labels for the columns.
ax
A `matplotlib.axes.Axes` instance to which the heatmap is plotted. If
not provided, use current axes or create a new one. Optional.
cbar_kw
A dictionary with arguments to `matplotlib.Figure.colorbar`. Optional.
cbarlabel
The label for the colorbar. Optional.
**kwargs
All other arguments are forwarded to `imshow`.
"""
if not ax:
ax = plt.gca()
# Plot the heatmap
im = ax.imshow(data, **kwargs)
# Create colorbar
cbar = ax.figure.colorbar(im, ax=ax, **cbar_kw)
cbar.ax.set_ylabel(cbarlabel, rotation=-90, va="bottom")
cbar.outline.set_visible(False)
# Show all ticks and label them with the respective list entries.
ax.set_xticks(np.arange(data.shape[1]))
ax.set_xticklabels(col_labels)
ax.set_yticks(np.arange(data.shape[0]))
ax.set_yticklabels(row_labels)
# Let the horizontal axes labeling appear on top.
#ax.tick_params(top=True, bottom=False,
# labeltop=True, labelbottom=False)
# Rotate the tick labels and set their alignment.
plt.setp(ax.get_xticklabels(), rotation=30, ha="right",
rotation_mode="anchor")
# Turn spines off and create white grid.
#ax.spines[:].set_visible(False)
ax.set_frame_on(False)
ax.set_xlabel(xlabel)
ax.set_ylabel(ylabel)
if title:
ax.set_title(title)
ax.set_xticks(np.arange(data.shape[1]+1)-.5, minor=True)
ax.set_yticks(np.arange(data.shape[0]+1)-.5, minor=True)
ax.grid(which="minor", color="w", linestyle='-', linewidth=3)
ax.tick_params(which="minor", bottom=False, left=False)
return im, cbar
def annotate_heatmap(im, data=None, valfmt="{x:.2f}",
textcolors=("black", "white"),
threshold=None, **textkw):
"""
A function to annotate a heatmap.
Parameters
----------
im
The AxesImage to be labeled.
data
Data used to annotate. If None, the image's data is used. Optional.
valfmt
The format of the annotations inside the heatmap. This should either
use the string format method, e.g. "$ {x:.2f}", or be a
`matplotlib.ticker.Formatter`. Optional.
textcolors
A pair of colors. The first is used for values below a threshold,
the second for those above. Optional.
threshold
Value in data units according to which the colors from textcolors are
applied. If None (the default) uses the middle of the colormap as
separation. Optional.
**kwargs
All other arguments are forwarded to each call to `text` used to create
the text labels.
"""
if not isinstance(data, (list, np.ndarray)):
data = im.get_array()
# Normalize the threshold to the images color range.
if threshold is not None:
threshold = im.norm(threshold)
else:
threshold = im.norm(data.max())/2.
# Set default alignment to center, but allow it to be
# overwritten by textkw.
kw = dict(horizontalalignment="center",
verticalalignment="center")
kw.update(textkw)
# Get the formatter in case a string is supplied
if isinstance(valfmt, str):
valfmt = matplotlib.ticker.StrMethodFormatter(valfmt)
# Loop over the data and create a `Text` for each "pixel".
# Change the text's color depending on the data.
texts = []
for i in range(data.shape[0]):
for j in range(data.shape[1]):
kw.update(color=textcolors[int(im.norm(data[i, j]) > threshold)])
text = im.axes.text(j, i, valfmt(data[i, j], None), **kw)
texts.append(text)
return texts
def reformat_ppi_file(inputfile,outputfile):
with open(inputfile,'r') as f:
with open(outputfile,'w') as g:
nodes_by_layer = collections.defaultdict(set)
for line in f:
layer, node1, node2, edgeweight = line.split()
nodes_by_layer[layer] = nodes_by_layer[layer].union((node1,node2))
reformatted_line = ' '.join((node1, layer, node2, layer))
g.write(reformatted_line+'\n')
# interlayer edges
for layer1,layer2 in itertools.combinations(nodes_by_layer.keys(),2):
for node in nodes_by_layer[layer1].intersection(nodes_by_layer[layer2]):
reformatted_line = ' '.join((node,layer1,node,layer2))
g.write(reformatted_line+'\n')