Easy Rejection Sampling in msBayes

msBayes is a suite of tools that support a Approximate Bayesian Computation approach for phylogeographic analysis, estimation and hypothesis testing. The following Python script facilitates the execution of the rejection sampling step of the analysis, by wrapping the invocation to the "msReject" program, automatically identifying the summary statistic columns in the data files and composing the correct command arguments. The script depends crucially on a header row of column names as the first line of the files, and assumes that all prior/parameter column names begin with the string "PRI.". #! /usr/bin/env python ############################################################################### ## Copyright 2009 Jeet Sukumaran. ## ## This program is free software; you can redistribute it and/or modify ## it under the terms of the GNU General Public License as published by ## the Free Software Foundation; either version 3 of the License, or ## (at your option) any later version. ## ## This program is distributed in the hope that it will be useful, ## but WITHOUT ANY WARRANTY; without even the implied warranty of ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ## GNU General Public License for more details. ## ## You should have received a copy of the GNU General Public License along ## with this program. If not, see . ## ############################################################################### """ Performs accept/reject sampling on a file of simulations from priors to produce set of samples the posterior distribution. (1) Identifies summary statistic columns in a given data file (based on column names that do not begin with "PRI.". (2) Calls msReject with the appriopriate parameters. """ from optparse import OptionGroup from optparse import OptionParser import subprocess import sys import os _prog_usage = '%prog [options] ' _prog_version = 'SAMPLE-POSTERIORS Version 1.0' _prog_description = """Performs accept/reject sampling on a file of simulations from priors to produce set of samples the posterior distribution.""" _prog_author = 'Jeet Sukumaran' _prog_copyright = 'Copyright (C) 2009 Jeet Sukumaran.' def parse_fields(row, separator="\t", include_empty=False): fields = row.strip().replace("\n", "").split(separator) if not include_empty: fields = [f for f in fields if f != ""] return fields def main(): """ Main CLI handler. """ parser = OptionParser(usage=_prog_usage, add_help_option=True, version=_prog_version, description=_prog_description) parser.add_option('-a', '--application', action='store', dest='app_path', type='string', # also 'float', 'string' etc. default="msReject", metavar='MSREJECT-PATH', help='path to msReject (default="%default")') parser.add_option('-t', '--tolerance', action='store', dest='tolerance', type='float', # also 'float', 'string' etc. default=0.02, metavar='TOLERANCE', help='sampling tolerance (default=%default)') parser.add_option('-s', '--sep', action='store', dest='separator', type='string', # also 'float', 'string' etc. default='\t', metavar='COLUMN-SEPARATOR', help='character to use as column separator (default=)') parser.add_option('-n', '--no-execute', action='store_true', dest="dry_run", default=False, help='print command to standard output, but do not actually run') parser.add_option('-q', '--quiet', action='store_true', dest="quiet", default=False, help='run without any supplemental messages') (opts, args) = parser.parse_args() if len(args)