From 33da62c9fe9aebe1f28d89bcb58ea2c0390db00d Mon Sep 17 00:00:00 2001 From: Hui Lan Date: Fri, 6 Dec 2019 15:16:11 +0800 Subject: draw_subnetwork.py: add a known network (G0) from the thermomorphogenesis paper. Fixed a bug. Now I close figure (plt.close()) before creating a new one, to avoid that the current figure is drawn on top of the old one. -Hui --- Code/draw_subnetwork.py | 161 ++++++++++++++++++++++++++++++++++++------------ 1 file changed, 120 insertions(+), 41 deletions(-) diff --git a/Code/draw_subnetwork.py b/Code/draw_subnetwork.py index 79ab5cd..d05f2d4 100644 --- a/Code/draw_subnetwork.py +++ b/Code/draw_subnetwork.py @@ -1,7 +1,10 @@ # Usage: python draw_subnetwork.py edges.txt -# Purpose: draw a sub-network given a list of genes. +# Purpose: draw a sub-network given a list of genes from thermomorphogenesis paper. +# The paper Molecular and genetic control of plant thermomorphogenesis. https://doi.org/10.1038/nplants.2015.190 +# # Created on 5 December 2019 by Hui Lan (lanhui@zjnu.edu.cn) + import os, sys import networkx as nx import pylab as plt @@ -9,6 +12,7 @@ import glob import math from networkx.algorithms.distance_measures import diameter, eccentricity + def build_network_from_file(edge_fname, gene_lst, gene_dict): G = nx.DiGraph() @@ -46,6 +50,33 @@ def build_network_from_file(edge_fname, gene_lst, gene_dict): return G + +def build_thermomorphogenesis_network(gene_lst, gene_dict): + ''' Edges from thermo.png in my asus laptop. ''' + G = nx.DiGraph() + + for g in gene_lst: + G.add_node(gene_dict[g]) + + pairs = [('PIF4', 'PAR1'), ('PIF4', 'PRE1'), ('PIF4', 'YUC8'), ('PIF4', 'TAA1'), ('PIF4', 'IAA4'), ('PIF4', 'SAUR21'), + ('HY5', 'PAR1'), ('HY5', 'PRE1'), ('HY5', 'YUC8'), ('HY5', 'TAA1'), ('HY5', 'IAA4'), ('HY5', 'SAUR21'), + ('FCA', 'PAR1'), ('FCA', 'PRE1'), ('FCA', 'YUC8'), ('FCA', 'TAA1'), ('FCA', 'IAA4'), ('FCA', 'SAUR21'), + ('BZR1', 'PAR1'), ('BZR1', 'PRE1'), ('BZR1', 'YUC8'), ('BZR1', 'TAA1'), ('BZR1', 'IAA4'), ('BZR1', 'SAUR21'), + ('ARF6', 'PAR1'), ('ARF6', 'PRE1'), ('ARF6', 'YUC8'), ('ARF6', 'TAA1'), ('ARF6', 'IAA4'), ('ARF6', 'SAUR21'), + ('PAR1', 'IBH1'), ('PRE1', 'IBH1'), + ('PAR1', 'PIF4'), ('PRE1', 'PIF4'), + ('IBH1', 'HBI1'), + ('IAA4', 'ARF6'), + ('ARF6', 'SAUR21') + ] # see paper Molecular and genetic control of plant thermomorphogenesis. https://doi.org/10.1038/nplants.2015.190 + + for (g2, g1) in pairs: # g2 is source, g1 is target + G.add_edge(g2, g1, weight=2) # g2 is source, and g1 is target + + + return G + + def compute_total_edge_weight(edges, G): total = 0 for e in edges: @@ -71,18 +102,22 @@ def draw_graph(G, fname): nx.draw_networkx_edges(G,pos,edgelist=esmall,width=1,alpha=0.1,edge_color='k',style='dashed') nx.draw_networkx_labels(G,pos,font_size=8,font_color='k',font_family='sans-serif') plt.axis('off') - plt.savefig(fname) # save as png + plt.savefig(fname) + plt.close() #plt.show() # display def better_date(s): + ''' Add a dash between year and month, and a dash between month and day.''' if len(s) == 8: return '-'.join([s[:4], s[4:6], s[6:]]) else: return s + def draw_graph2(G, fname, date): - pos=nx.circular_layout(G) + + pos = nx.circular_layout(G) all_edges = [] all_widths = [] @@ -91,11 +126,29 @@ def draw_graph2(G, fname, date): all_widths.append(math.sqrt(d['weight'])) nx.draw_networkx_nodes(G,pos,alpha=0.05) - nx.draw_networkx_edges(G,pos,edgelist=all_edges,width=all_widths,alpha=0.05,edge_color='k',style='dashed') + nx.draw_networkx_edges(G,pos,edgelist=all_edges,width=all_widths,alpha=0.2,edge_color='k',style='dashed') nx.draw_networkx_labels(G,pos,font_size=11,font_color='b',font_family='sans-serif') plt.axis('off') plt.title(better_date(date)) - plt.savefig(fname) # save as png + plt.savefig(fname) + plt.close() + #plt.show() # display + + +def draw_graph3(G, fname): + + pos = nx.circular_layout(G) + all_edges = [] + + for (u, v, d) in G.edges(data=True): + all_edges.append((u, v)) + + nx.draw_networkx_nodes(G,pos,alpha=0.05) + nx.draw_networkx_edges(G,pos,edgelist=all_edges,alpha=0.2,edge_color='k',style='dashed') + nx.draw_networkx_labels(G,pos,font_size=11,font_color='b',font_family='sans-serif') + plt.axis('off') + plt.savefig(fname) + plt.close() # it is important to close the plot before creating another one. Otherwise, plots will overlap. #plt.show() # display @@ -165,6 +218,7 @@ thermomorphogenesis_genes_small = [ 'AT5G11260', #HY5 'AT2G42870', #PAR1 'AT5G39860', #PRE1 + 'AT5G43700', #IAA4 'AT4G16280', #FCA 'AT2G43060', #IBH1 'AT2G18300', #HBI1 @@ -172,18 +226,19 @@ thermomorphogenesis_genes_small = [ 'AT1G70560', #TAA1 'AT1G30330', #ARF6 'AT1G19850', #ARF5 - 'AT1G75080', #BZR1 - 'AT2G25930', #ELF3 - 'AT2G40080', #ELF4 - 'AT3G46640', #LUX + 'AT5G01830', #SAUR21 + 'AT1G75080' #BZR1 +# 'AT2G25930', #ELF3 +# 'AT2G40080', #ELF4 +# 'AT3G46640', #LUX ] - gene_dict = { 'AT2G43010':'PIF4', 'AT5G11260':'HY5', 'AT2G42870':'PAR1', 'AT5G39860':'PRE1', + 'AT5G43700':'IAA4', 'AT4G16280':'FCA', 'AT2G43060':'IBH1', 'AT2G18300':'HBI1', @@ -191,13 +246,18 @@ gene_dict = { 'AT1G70560':'TAA1', 'AT1G30330':'ARF6', 'AT1G19850':'ARF5', + 'AT5G01830':'SAUR21', 'AT1G75080':'BZR1', 'AT2G25930':'ELF3', 'AT2G40080':'ELF4', 'AT3G46640':'LUX' } -print('Make sub graph...') +print('Make sub graphs ...') + +G0 = build_thermomorphogenesis_network(thermomorphogenesis_genes_small, gene_dict) +print('Number of edges in the paper thermomorphogenesis is %d' % (len(G0.edges()))) +draw_graph2(G0, '../Data/temp/graph-%s.pdf' % ('20160101'), '') graph_lst = [] graph_names = [] @@ -206,43 +266,62 @@ for fname in sorted(glob.glob('../Analysis/edges.txt.2019*')): continue print(fname) graph_names.append(fname) - graph_lst.append(build_network_from_file(fname, thermomorphogenesis_genes_small, gene_dict)) + G = build_network_from_file(fname, thermomorphogenesis_genes_small, gene_dict) + graph_lst.append(G) -#G1205 = build_network_from_file('../Analysis/edges.txt.20191205', thermomorphogenesis_genes) -#G1203 = build_network_from_file('../Analysis/edges.txt.20191203', thermomorphogenesis_genes) -#G1126 = build_network_from_file('../Analysis/edges.txt.20191126', thermomorphogenesis_genes) -#G1108 = build_network_from_file('../Analysis/edges.txt.20191108', thermomorphogenesis_genes) +for i in range(len(graph_lst)): + G = graph_lst[i] + print('Graph from %s' % (graph_names[i])) + e = G.edges(data=True) + n = len(e) + print('Number of edges is %d' % (n)) + print('------------------------------------------------------------------------') + draw_graph2(G, '../Data/temp/graph-%s.pdf' % (graph_names[i].split('.')[-1]), graph_names[i].split('.')[-1]) -for i in range(len(graph_lst)-1): - G1 = graph_lst[i] - G2 = graph_lst[i+1] +print('Compute network differences ...') +Gdiff1 = nx.difference(G0, G) +draw_graph3(Gdiff1, '../Data/temp/graph-in-G0-not-in-G.pdf') - print(nx.is_isomorphic(G1, G2)) - print('Graph 1 from %s' % (graph_names[i])) - e1 = G1.edges(data=True) - n1 = len(e1) - tw1 = compute_total_edge_weight(e1, G1) - print('Number of edges is %d' % (n1)) - if n1 > 0: - print('Total edge association strength is %4.2f (avg=%4.2f).' % (tw1, tw1/n1)) - #print(e1) +Gdiff2 = nx.difference(G, G0) +draw_graph3(Gdiff2, '../Data/temp/graph-in-G-not-in-G0.pdf') - print('Graph 2 from %s' % (graph_names[i+1])) - e2 = G2.edges(data=True) - n2 = len(e2) - tw2 = compute_total_edge_weight(e2, G2) - print('Number of edges is %d' % (n2)) - if n2 > 0: - print('Total edge association strength is %4.2f (avg=%4.2f).' % (tw2, tw2/n2)) - #print(e2) +print('Compute network intersection ...') +Gdiff1 = nx.intersection(G0, G) +draw_graph3(Gdiff1, '../Data/temp/graph-in-G0-and-in-G.pdf') - ged = nx.algorithms.similarity.graph_edit_distance(G1, G2) - print(ged) +print('Compute edit distance ...') +ged = nx.algorithms.similarity.graph_edit_distance(G0, G) # sometimes take very long time to finish +print('Edit distance is %4.0f' % (ged)) - print('------------------------------------------------------------------------') +# for i in range(len(graph_lst)-1): +# G1 = graph_lst[i] +# G2 = graph_lst[i+1] - draw_graph2(G1, '../Data/temp/graph-%s.png' % (graph_names[i].split('.')[-1]), graph_names[i].split('.')[-1]) - +# print(nx.is_isomorphic(G1, G2)) +# print('Graph 1 from %s' % (graph_names[i])) +# e1 = G1.edges(data=True) +# n1 = len(e1) +# tw1 = compute_total_edge_weight(e1, G1) +# print('Number of edges is %d' % (n1)) +# if n1 > 0: +# print('Total edge association strength is %4.2f (avg=%4.2f).' % (tw1, tw1/n1)) +# #print(e1) + +# print('Graph 2 from %s' % (graph_names[i+1])) +# e2 = G2.edges(data=True) +# n2 = len(e2) +# tw2 = compute_total_edge_weight(e2, G2) +# print('Number of edges is %d' % (n2)) +# if n2 > 0: +# print('Total edge association strength is %4.2f (avg=%4.2f).' % (tw2, tw2/n2)) +# #print(e2) + +# ged = nx.algorithms.similarity.graph_edit_distance(G1, G2) +# print(ged) + +# print('------------------------------------------------------------------------') + +# draw_graph2(G1, '../Data/temp/graph-%s.png' % (graph_names[i].split('.')[-1]), graph_names[i].split('.')[-1]) -- cgit v1.2.1