Mercurial > repos > miller-lab > genome_diversity
diff dpmix_plot.py @ 31:a631c2f6d913
Update to Miller Lab devshed revision 3c4110ffacc3
author | Richard Burhans <burhans@bx.psu.edu> |
---|---|
date | Fri, 20 Sep 2013 13:25:27 -0400 |
parents | 8997f2ca8c7a |
children |
line wrap: on
line diff
--- a/dpmix_plot.py Fri Jul 26 12:51:13 2013 -0400 +++ b/dpmix_plot.py Fri Sep 20 13:25:27 2013 -0400 @@ -3,8 +3,10 @@ import os import sys import math + import matplotlib as mpl mpl.use('PDF') +from matplotlib.backends.backend_pdf import PdfPages import matplotlib.pyplot as plt from matplotlib.path import Path import matplotlib.patches as patches @@ -226,18 +228,49 @@ return vals, labels ################################################################################ +################################################################################ +################################################################################ +################################################################################ + +def space_for_legend(plot_params): + space = 0.0 + + legend_states = plot_params['legend_states'] + if legend_states: + ind_space = plot_params['ind_space'] + ind_height = plot_params['ind_height'] + space += len(legend_states) * (ind_space + ind_height) - ind_space + + return space + +################################################################################ + +def space_for_chroms(plot_params, chroms, individuals, data): + space_dict = {} + + chrom_height = plot_params['chrom_height'] + ind_space = plot_params['ind_space'] + ind_height = plot_params['ind_height'] + + for chrom in chroms: + space_dict[chrom] = chrom_height + + individual_count = 0 + for individual in individuals: + if individual in data[chrom]: + individual_count += 1 + + space_dict[chrom] += individual_count * (ind_space + ind_height) + + return space_dict + +################################################################################ def make_dpmix_plot(input_dbkey, input_file, output_file, galaxy_data_index_dir, state2name=None, populations=3): fs_chrom_len = build_chrom_len_dict(input_dbkey, galaxy_data_index_dir) chroms, individuals, data, chrom_len, used_states = parse_input_file(input_file) - if populations == 3: - make_state_rectangle = make_state_rectangle_3pop - elif populations == 2: - make_state_rectangle = make_state_rectangle_2pop - else: - pass - + ## populate chrom_len for chrom in chrom_len.keys(): if chrom in fs_chrom_len: chrom_len[chrom] = fs_chrom_len[chrom] @@ -245,135 +278,177 @@ #check_chroms(chroms, chrom_len, input_dbkey) check_data(data, chrom_len, input_dbkey) - ## units below are inches - top_space = 0.10 - chrom_space = 0.25 - chrom_height = 0.25 - ind_space = 0.10 - ind_height = 0.25 + ## plot parameters + plot_params = { + 'plot_dpi': 300, + 'page_width': 8.50, + 'page_height': 11.00, + 'top_margin': 0.10, + 'bottom_margin': 0.10, + 'chrom_space': 0.25, + 'chrom_height': 0.25, + 'ind_space': 0.10, + 'ind_height': 0.25, + 'legend_space': 0.10 + } - total_height = 0.0 - - ## make a legend - ## only print out states that are + ## in the legend, only print out states that are ## 1) in the data ## - AND - ## 2) in the state2name map - ## here, we only calculate the space needed legend_states = [] if state2name is not None: for state in used_states: if state in state2name: legend_states.append(state) - if legend_states: - total_height += len(legend_states) * (ind_space + ind_height) - total_height += (top_space - ind_space) - at_top = False - else: - at_top = True + plot_params['legend_states'] = legend_states + + ## choose the correct make_state_rectangle method + if populations == 3: + plot_params['rectangle_method'] = make_state_rectangle_3pop + elif populations == 2: + plot_params['rectangle_method'] = make_state_rectangle_2pop + + pdf_pages = PdfPages(output_file) + + ## generate a list of chroms for each page + + needed_for_legend = space_for_legend(plot_params) + needed_for_chroms = space_for_chroms(plot_params, chroms, individuals, data) - for chrom in chroms: - if at_top: - total_height += (top_space + chrom_height) - at_top = False + chrom_space_per_page = plot_params['page_height'] + chrom_space_per_page -= plot_params['top_margin'] + plot_params['bottom_margin'] + chrom_space_per_page -= needed_for_legend + plot_params['legend_space'] + chrom_space_per_page -= plot_params['chrom_space'] + + chroms_left = chroms[:] + pages = [] + + space_left = chrom_space_per_page + chrom_list = [] + + while chroms_left: + chrom = chroms_left.pop(0) + space_needed = needed_for_chroms[chrom] + plot_params['chrom_space'] + if (space_needed > chrom_space_per_page): + print >> sys.stderr, 'Multipage chroms not yet supported' + sys.exit(1) + + ## sometimes 1.9 - 1.9 < 0 (-4.4408920985e-16) + ## so, we make sure it's not more than a millimeter over + if space_left - space_needed > -0.04: + chrom_list.append(chrom) + space_left -= space_needed else: - total_height += (top_space + chrom_space + chrom_height) - - individual_count = 0 - for individual in individuals: - if individual in data[chrom]: - individual_count += 1 - total_height += individual_count * (ind_space + ind_height) + pages.append(chrom_list[:]) + chrom_list = [] + chroms_left.insert(0, chrom) + space_left = chrom_space_per_page - width = 7.5 - height = math.ceil(total_height) - - bottom = 1.0 - - fig = plt.figure(figsize=(width, height)) + ############################################################################ - if legend_states: - at_top = True - for state in sorted(legend_states): - if at_top: - bottom -= (top_space + ind_height)/height - at_top = False - else: - bottom -= (ind_space + ind_height)/height - ## add code here to draw legend - # [left, bottom, width, height] - ax1 = fig.add_axes([0.0, bottom, 0.09, ind_height/height]) - plt.axis('off') - ax1.set_xlim(0, 1) - ax1.set_ylim(0, 1) - for patch in make_state_rectangle(0, 1, state, 'legend', state2name[state]): - ax1.add_patch(patch) + plot_dpi = plot_params['plot_dpi'] + page_width = plot_params['page_width'] + page_height = plot_params['page_height'] + top_margin = plot_params['top_margin'] + ind_space = plot_params['ind_space'] + ind_height = plot_params['ind_height'] + make_state_rectangle = plot_params['rectangle_method'] + legend_space = plot_params['legend_space'] + chrom_space = plot_params['chrom_space'] + chrom_height = plot_params['chrom_height'] + + for page in pages: + fig = plt.figure(figsize=(page_width, page_height), dpi=plot_dpi) + bottom = 1.0 - (top_margin/page_height) + + # print legend + if legend_states: + top = True + for state in sorted(legend_states): + if top: + bottom -= ind_height/page_height + top = False + else: + bottom -= (ind_space + ind_height)/page_height - ax2 = fig.add_axes([0.10, bottom, 0.88, ind_height/height], frame_on=False) - plt.axis('off') - plt.text(0.0, 0.5, state2name[state], fontsize=10, ha='left', va='center') - else: - at_top = True + ax1 = fig.add_axes([0.0, bottom, 0.09, ind_height/page_height]) + plt.axis('off') + ax1.set_xlim(0, 1) + ax1.set_ylim(0, 1) + for patch in make_state_rectangle(0, 1, state, 'legend', state2name[state]): + ax1.add_patch(patch) - for_webb = False + ax2 = fig.add_axes([0.10, bottom, 0.88, ind_height/page_height], frame_on=False) + plt.axis('off') + plt.text(0.0, 0.5, state2name[state], fontsize=10, ha='left', va='center') + + bottom -= legend_space/page_height - for chrom in chroms: - length = chrom_len[chrom] - vals, labels = tick_foo(0, length) + # print chroms + top = True + for chrom in page: + length = chrom_len[chrom] + vals, labels = tick_foo(0, length) - if at_top: - bottom -= (top_space + chrom_height)/height - at_top = False - else: - bottom -= (top_space + chrom_space + chrom_height)/height + if top: + bottom -= chrom_height/page_height + top = False + else: + bottom -= (chrom_space + chrom_height)/page_height - if not for_webb: - ax = fig.add_axes([0.0, bottom, 1.0, chrom_height/height]) + ax = fig.add_axes([0.0, bottom, 1.0, chrom_height/page_height]) plt.axis('off') plt.text(0.5, 0.5, chrom, fontsize=14, ha='center') - individual_count = 0 - for individual in individuals: - if individual in data[chrom]: - individual_count += 1 + individual_count = 0 + for individual in individuals: + if individual in data[chrom]: + individual_count += 1 - i = 0 - for individual in individuals: - if individual in data[chrom]: - i += 1 + i = 0 + for individual in individuals: + if individual in data[chrom]: + i += 1 + bottom -= (ind_space + ind_height)/page_height - bottom -= (ind_space + ind_height)/height - if not for_webb: - # [left, bottom, width, height] - ax1 = fig.add_axes([0.0, bottom, 0.09, ind_height/height]) + ax1 = fig.add_axes([0.0, bottom, 0.09, ind_height/page_height]) plt.axis('off') plt.text(1.0, 0.5, individual, fontsize=10, ha='right', va='center') - # [left, bottom, width, height] - ax2 = fig.add_axes([0.10, bottom, 0.88, ind_height/height], frame_on=False) - ax2.set_xlim(0, length) - ax2.set_ylim(0, 1) - if i != individual_count: - plt.axis('off') - else: - if not for_webb: + + ax2 = fig.add_axes([0.10, bottom, 0.88, ind_height/page_height], frame_on=False) + ax2.set_xlim(0, length) + ax2.set_ylim(0, 1) + + if i != individual_count: + plt.axis('off') + else: ax2.tick_params(top=False, left=False, right=False, labelleft=False) ax2.set_xticks(vals) ax2.set_xticklabels(labels) - else: - plt.axis('off') - for p1, p2, state in sorted(data[chrom][individual]): - #for patch in make_state_rectangle(p1, p2, state, chrom, individual): - for patch in make_state_rectangle(p1, p2, state, chrom, individual): - ax2.add_patch(patch) + + for p1, p2, state in sorted(data[chrom][individual]): + for patch in make_state_rectangle(p1, p2, state, chrom, individual): + ax2.add_patch(patch) - plt.savefig(output_file) + # extend last state to end of chrom + if p2 < length: + for patch in make_state_rectangle(p2, length, state, chrom, individual): + ax2.add_patch(patch) + + + pdf_pages.savefig(fig) + plt.close(fig) + + pdf_pages.close() ################################################################################ if __name__ == '__main__': - input_dbkey, input_file, output_file, galaxy_data_index_dir = sys.argv[1:5] - make_dpmix_plot(input_dbkey, input_file, output_file, galaxy_data_index_dir) + make_dpmix_plot('loxAfr3', 'output.dat', 'output2_files/picture.pdf', '/scratch/galaxy/home/oocyte/galaxy_oocyte/tool-data', state2name={0: 'heterochromatin', 1: 'reference', 2: 'asian'}, populations=2) +# input_dbkey, input_file, output_file, galaxy_data_index_dir = sys.argv[1:5] +# make_dpmix_plot(input_dbkey, input_file, output_file, galaxy_data_index_dir) sys.exit(0) ## notes