Package Bio :: Package PopGen :: Package SimCoal :: Module Template
[hide private]
[frames] | no frames]

Source Code for Module Bio.PopGen.SimCoal.Template

  1  # Copyright 2007 by Tiago Antao.  All rights reserved. 
  2  # This code is part of the Biopython distribution and governed by its 
  3  # license.  Please see the LICENSE file that should have been included 
  4  # as part of this package. 
  5   
  6  from __future__ import print_function 
  7   
  8  from os import sep 
  9  import re 
 10  from functools import reduce 
 11   
 12  from Bio.PopGen.SimCoal import builtin_tpl_dir 
 13   
 14   
15 -def exec_template(template):
16 executed_template = template 17 match = re.search('!!!(.*?)!!!', executed_template, re.MULTILINE) 18 #while len(match.groups())>0: 19 while match: 20 exec_result = str(eval(match.groups()[0])) 21 executed_template = executed_template.replace( 22 '!!!' + match.groups()[0] + '!!!', 23 exec_result, 1) 24 match = re.search('!!!(.*?)!!!', executed_template, re.MULTILINE) 25 #match = patt.matcher(String(executed_template)) 26 return executed_template
27 28
29 -def process_para(in_string, out_file_prefix, para_list, curr_values):
30 if (para_list == []): 31 template = in_string 32 f_name = out_file_prefix 33 #f_name += '_' + str(total_size) 34 for tup in curr_values: 35 name, val = tup 36 f_name += '_' + str(val) 37 #reg = re.compile('\?' + name, re.MULTILINE) 38 #template = re.sub(reg, str(val), template) 39 template = template.replace('?'+name, str(val)) 40 with open(f_name + '.par', 'w') as f: 41 #executed_template = template 42 executed_template = exec_template(template) 43 clean_template = executed_template.replace('\r\n', '\n').replace('\n\n', '\n') 44 f.write(clean_template) 45 return [f_name] 46 else: 47 name, rng = para_list[0] 48 fnames = [] 49 for val in rng: 50 new_values = [(name, val)] 51 new_values.extend(curr_values) 52 more_names = process_para(in_string, out_file_prefix, para_list[1:], new_values) 53 fnames.extend(more_names) 54 return fnames
55 56
57 -def dupe(motif, times):
58 ret_str = '' 59 for i in range(1, times + 1): 60 ret_str += motif + '\r\n' 61 return ret_str
62 63
64 -def get_xy_from_matrix(x_max, y_max, pos):
65 y = (pos-1) / x_max 66 x = (pos-1) % x_max 67 return x, y
68 69
70 -def get_step_2d(x_max, y_max, x, y, mig):
71 my_x, my_y = get_xy_from_matrix(x_max, y_max, y) 72 other_x, other_y = get_xy_from_matrix(x_max, y_max, x) 73 74 if (my_x-other_x)**2 + (my_y-other_y)**2 == 1: 75 return str(mig) + ' ' 76 else: 77 return '0 '
78 79
80 -def generate_ssm2d_mat(x_max, y_max, mig):
81 mig_mat = '' 82 for x in range(1, x_max*y_max + 1): 83 for y in range(1, x_max*y_max + 1): 84 mig_mat += get_step_2d(x_max, y_max, x, y, mig) 85 mig_mat += "\r\n" 86 return mig_mat
87 88
89 -def generate_island_mat(total_size, mig):
90 mig_mat = '' 91 for x in range(1, total_size + 1): 92 for y in range(1, total_size + 1): 93 if (x == y): 94 mig_mat += '0 ' 95 else: 96 mig_mat += '!!!' + str(mig) + '!!! ' 97 mig_mat += "\r\n" 98 return mig_mat
99 100
101 -def generate_null_mat(total_size):
102 null_mat = '' 103 for x in range(1, total_size + 1): 104 for y in range(1, total_size + 1): 105 null_mat += '0 ' 106 null_mat += '\r\n' 107 return null_mat
108 109
110 -def generate_join_events(t, total_size, join_size, orig_size):
111 events = '' 112 for i in range(1, total_size-1): 113 events += str(t) + ' ' + str(i) + ' 0 1 1 0 1\r\n' 114 events += str(t) + ' ' + str(total_size-1) + ' 0 1 ' + str(1.0*total_size*join_size/orig_size) + ' 0 1\r\n' 115 return events
116 117
118 -def no_processor(in_string):
119 return in_string
120 121
122 -def process_text(in_string, out_file_prefix, para_list, curr_values, 123 specific_processor):
124 text = specific_processor(in_string) 125 return process_para(text, out_file_prefix, para_list, [])
126 127 128 #def prepare_dir(): 129 # try: 130 # mkdir(sep.join([Config.dataDir, 'SimCoal'])) #Should exist, but... 131 # except OSError: 132 # pass #Its ok if already exists 133 # try: 134 # mkdir(sep.join([Config.dataDir, 'SimCoal', 'runs'])) 135 # except OSError: 136 # pass #Its ok if already exists 137 138 139 #sep is because of jython
140 -def generate_model(par_stream, out_prefix, params, 141 specific_processor = no_processor, out_dir = '.'):
142 #prepare_dir() 143 text = par_stream.read() 144 out_file_prefix = sep.join([out_dir, out_prefix]) 145 return process_text(text, out_file_prefix, params, [], specific_processor)
146 147
148 -def get_demography_template(stream, model, tp_dir = None):
149 ''' 150 Gets a demograpy template. 151 152 Most probably this model needs to be sent to GenCases. 153 154 stream - Writable stream. 155 param - Template file. 156 tp_dir - Directory where to find the template, if None 157 use an internal template 158 ''' 159 if tp_dir is None: 160 #Internal Template 161 filename = sep.join([builtin_tpl_dir, model + '.par']) 162 else: 163 #External template 164 filename = sep.join([tp_dir, model + '.par']) 165 with open(filename, 'r') as f: 166 l = f.readline() 167 while l!='': 168 stream.write(l) 169 l = f.readline()
170 171
172 -def _gen_loci(stream, loci):
173 stream.write('//Number of contiguous linkage blocks in chromosome\n') 174 stream.write(str(len(loci)) + '\n') 175 stream.write('//Per Block: Data type, No. of loci, Recombination rate to the right-side locus, plus optional parameters\n') 176 for locus in loci: 177 stream.write(' '.join([locus[0]] + 178 [str(x) for x in list(locus[1])]) + '\n')
179 180
181 -def get_chr_template(stream, chrs):
182 ''' 183 Writes a Simcoal2 loci template part. 184 185 stream - Writable stream. 186 chr - Chromosome list. 187 188 Current loci list: 189 [(chr_repeats,[(marker, (params))])] 190 chr_repeats --> Number of chromosome repeats 191 marker --> 'SNP', 'DNA', 'RFLP', 'MICROSAT' 192 params --> Simcoal2 parameters for markers (list of floats 193 or ints - if to be processed by generate_model) 194 ''' 195 num_chrs = reduce(lambda x, y: x + y[0], chrs, 0) 196 stream.write('//Number of independent (unlinked) chromosomes, and "chromosome structure" flag: 0 for identical structure across chromosomes, and 1 for different structures on different chromosomes.\n') 197 if len(chrs) > 1 or num_chrs == 1: 198 stream.write(str(num_chrs) + ' 1\n') 199 else: 200 stream.write(str(num_chrs) + ' 0\n') 201 for chr in chrs: 202 repeats = chr[0] 203 loci = chr[1] 204 if len(chrs) == 1: 205 _gen_loci(stream, loci) 206 else: 207 for i in range(repeats): 208 _gen_loci(stream, loci)
209 210
211 -def generate_simcoal_from_template(model, chrs, params, out_dir = '.', tp_dir=None):
212 ''' 213 Writes a complete SimCoal2 template file. 214 215 This joins together get_demography_template and get_chr_template, 216 which are feed into generate_model 217 Please check the three functions for parameters (model from 218 get_demography_template, chrs from get_chr_template and 219 params from generate_model). 220 ''' 221 with open(out_dir + sep + 'tmp.par', 'w') as stream: 222 get_demography_template(stream, model, tp_dir) 223 get_chr_template(stream, chrs) 224 #with open(out_dir + sep + 'tmp.par', 'r') as par_stream: 225 #print par_stream.read() 226 with open(out_dir + sep + 'tmp.par', 'r') as par_stream: 227 generate_model(par_stream, model, params, out_dir = out_dir)
228