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