diff --git a/.github/workflows/build.yml b/.github/workflows/build.yml index 7cc595b1..91c0f4c4 100644 --- a/.github/workflows/build.yml +++ b/.github/workflows/build.yml @@ -17,7 +17,7 @@ jobs: with: python-version: ${{ matrix.python-version }} - name: Add micromamba to system path - uses: mamba-org/provision-with-micromamba@main + uses: mamba-org/setup-micromamba@main with: environment-name: dadi-cli activate-environment: dadi-cli diff --git a/dadi_cli/GenerateFs.py b/dadi_cli/GenerateFs.py index 17693777..36a5e4ae 100644 --- a/dadi_cli/GenerateFs.py +++ b/dadi_cli/GenerateFs.py @@ -94,9 +94,10 @@ def generate_fs( if calc_coverage: import pickle - coverage_dd = {chrom_pos:{'coverage':dd[chrom_pos]['coverage']} for chrom_pos in dd} - print(f"\nSaving coverage dictionary in pickle named:\n{output}.coverage.pickle\n") - pickle.dump(coverage_dd, open(f"{output}.coverage.pickle","wb")) + import dadi.LowPass.LowPass as lp + cov_dist = lp.compute_cov_dist(dd, pop_ids) + print(f"\nSaving coverage distribution data in pickle named:\n{output}.coverage.pickle\n") + pickle.dump(cov_dist, open(f"{output}.coverage.pickle","wb")) if bootstrap is None: fs = dadi.Spectrum.from_data_dict( diff --git a/dadi_cli/InferDFE.py b/dadi_cli/InferDFE.py index f2985265..4f512b3e 100644 --- a/dadi_cli/InferDFE.py +++ b/dadi_cli/InferDFE.py @@ -13,6 +13,7 @@ def infer_dfe( cache2d: dadi.DFE.Cache2D_mod, sele_dist: str, sele_dist2: str, + pdf_file: str, theta: float, p0: list[float], upper_bounds: list[float], @@ -96,9 +97,9 @@ def infer_dfe( func = cache2d.integrate if sele_dist is not None: - sele_dist = get_dadi_pdf(sele_dist) + sele_dist, _ = get_dadi_pdf(sele_dist, pdf_file) if sele_dist2 is not None: - sele_dist2 = get_dadi_pdf(sele_dist2) + sele_dist2, _ = get_dadi_pdf(sele_dist2, pdf_file) if (sele_dist is None) and (sele_dist2 is not None): sele_dist = sele_dist2 @@ -125,8 +126,8 @@ def infer_dfe( func = func_cov(func, cov_args[0], fs.pop_ids, nseq, fs.sample_sizes, Fx=Fx) p0_len = len(p0) - lower_bounds = convert_to_None(lower_bounds, p0_len) - upper_bounds = convert_to_None(upper_bounds, p0_len) + # lower_bounds = convert_to_None(lower_bounds, p0_len) + # upper_bounds = convert_to_None(upper_bounds, p0_len) fixed_params = convert_to_None(fixed_params, p0_len) # Fit a DFE to the data diff --git a/dadi_cli/InferDM.py b/dadi_cli/InferDM.py index 0a3662bb..00104efb 100644 --- a/dadi_cli/InferDM.py +++ b/dadi_cli/InferDM.py @@ -99,9 +99,8 @@ def infer_demography( func_ex = func_cov(func_ex, cov_args[0], fs.pop_ids, nseq, fs.sample_sizes, Fx=Fx) p0_len = len(p0) - lower_bounds = convert_to_None(lower_bounds, p0_len) - upper_bounds = convert_to_None(upper_bounds, p0_len) - fixed_params = convert_to_None(fixed_params, p0_len) + if fixed_params == -1: + fixed_params = convert_to_None(fixed_params, p0_len) p0 = dadi.Misc.perturb_params( p0, fold=1, upper_bound=upper_bounds, lower_bound=lower_bounds @@ -224,9 +223,8 @@ def infer_global_opt( func_ex = func_cov(func_ex, cov_args[0], fs.pop_ids, nseq, fs.sample_sizes, Fx=Fx) p0_len = len(p0) - lower_bounds = convert_to_None(lower_bounds, p0_len) - upper_bounds = convert_to_None(upper_bounds, p0_len) - fixed_params = convert_to_None(fixed_params, p0_len) + if fixed_params == -1: + fixed_params = convert_to_None(fixed_params, p0_len) p0 = dadi.Misc.perturb_params( p0, fold=1, upper_bound=upper_bounds, lower_bound=lower_bounds diff --git a/dadi_cli/Pdfs.py b/dadi_cli/Pdfs.py index a3722236..014fa5de 100644 --- a/dadi_cli/Pdfs.py +++ b/dadi_cli/Pdfs.py @@ -1,7 +1,12 @@ +import dadi, importlib, os, sys import dadi.DFE as DFE +from inspect import getmembers, isfunction -def get_dadi_pdf(pdf: str) -> callable: +def get_dadi_pdf( + pdf: str, + pdf_file: str = None + ) -> tuple[callable, list[str]]: """ Obtains a built-in probability density function for modeling the distribution of fitness effects in dadi. @@ -22,7 +27,27 @@ def get_dadi_pdf(pdf: str) -> callable: If the specified probability density function name is not recognized. """ - if pdf == "beta": + if pdf_file is not None: + # If the user has the model folder in their PATH + #try: + # func = getattr(importlib.import_module(model_file), model_name) + # If the user does not have the model folder in their PATH we add it + # This currently can mess with the User's PATH while running dadi-cli + #except: + # model_file = os.path.abspath(model_file) + # model_path = os.path.dirname(model_file) + # model_file = os.path.basename(model_file) + # model_file = os.path.splitext(model_file)[0] + # sys.path.append(model_path) + # func = getattr(importlib.import_module(model_file), model_name) + try: + module_name = os.path.splitext(os.path.basename(pdf_file))[0] + model_path = os.path.dirname(os.path.abspath(pdf_file)) + sys.path.append(model_path) + func = getattr(importlib.import_module(module_name), pdf) + except ImportError as e: + raise ImportError(f"Failed to import module: {pdf_file}") from e + elif pdf == "beta": func = DFE.PDFs.beta elif pdf == "biv_ind_gamma": func = DFE.PDFs.biv_ind_gamma @@ -41,7 +66,10 @@ def get_dadi_pdf(pdf: str) -> callable: else: raise ValueError("Probability density function " + pdf + " is not available!") - return func + if pdf_file is not None: + return func, func.__param_names__ + else: + return func, None def get_dadi_pdf_params(pdf: str) -> list[str]: diff --git a/dadi_cli/Plot.py b/dadi_cli/Plot.py index 75fee6e7..a4f215e8 100644 --- a/dadi_cli/Plot.py +++ b/dadi_cli/Plot.py @@ -229,6 +229,7 @@ def plot_fitted_dfe( projections: list[int], pdf: str, pdf2: str, + pdf_file: str, cov_args: list, cov_inbreeding: list, output: str, @@ -255,6 +256,8 @@ def plot_fitted_dfe( Name of the 1D probability density function file for modeling the DFE. pdf2 : str Name of the 2D probability density function file for modeling the DFE. + pdf_file : str + Name of file with custom probability density function model(s) in it. output : str Path where the comparison plot will be saved. The file format is inferred from the file extension. vmin : float @@ -275,9 +278,9 @@ def plot_fitted_dfe( fs = dadi.Spectrum.from_file(fs) if pdf != None: - pdf = get_dadi_pdf(pdf) + pdf, _ = get_dadi_pdf(pdf, pdf_file) if pdf2 != None: - pdf2 = get_dadi_pdf(pdf2) + pdf2, _ = get_dadi_pdf(pdf2, pdf_file) if pdf == None: pdf = pdf2 diff --git a/dadi_cli/SimulateFs.py b/dadi_cli/SimulateFs.py index 72e0684e..906b051c 100644 --- a/dadi_cli/SimulateFs.py +++ b/dadi_cli/SimulateFs.py @@ -102,6 +102,7 @@ def simulate_dfe( cache2d: dadi.DFE.Cache2D_mod, sele_dist: str, sele_dist2: str, + pdf_file: str, ratio: float, misid: bool, output: str @@ -143,9 +144,9 @@ def simulate_dfe( func = cache2d.integrate if sele_dist is not None: - sele_dist = get_dadi_pdf(sele_dist) + sele_dist, _ = get_dadi_pdf(sele_dist, pdf_file) if sele_dist2 is not None: - sele_dist2 = get_dadi_pdf(sele_dist2) + sele_dist2, _ = get_dadi_pdf(sele_dist2, pdf_file) if (sele_dist is None) and (sele_dist2 is not None): sele_dist = sele_dist2 diff --git a/dadi_cli/Stat.py b/dadi_cli/Stat.py index 25a06db4..4999e0b4 100644 --- a/dadi_cli/Stat.py +++ b/dadi_cli/Stat.py @@ -6,7 +6,7 @@ from typing import Optional from dadi_cli.Models import get_model from dadi_cli.Pdfs import get_dadi_pdf -from dadi_cli.utilities import get_opts_and_theta, convert_to_None +from dadi_cli.utilities import get_opts_and_theta, convert_to_None, pts_l_func def godambe_stat_demograpy( @@ -49,11 +49,15 @@ def godambe_stat_demograpy( # We want the best fits from the demograpic fit. # We set the second argument to true, since we want # all the parameters from the file. - demo_popt, _, param_names = get_opts_and_theta(demo_popt, post_infer=True) + demo_popt, theta, param_names = get_opts_and_theta(demo_popt, post_infer=True) fixed_params = convert_to_None(fixed_params, len(demo_popt)) free_params = _free_params(demo_popt, fixed_params) fs = dadi.Spectrum.from_file(fs) + if grids is None: + grids = pts_l_func(fs.sample_sizes) fs_boot_files = glob.glob(bootstrap_dir + "/*.fs") + if fs_boot_files == []: + raise ValueError(f"ERROR:\n{bootstrap_dir} is empty\n") all_boot = [] for f in fs_boot_files: boot_fs = dadi.Spectrum.from_file(f) @@ -74,9 +78,10 @@ def gfunc(free_params, ns, pts): uncerts_adj = dadi.Godambe.GIM_uncert( gfunc, grids, all_boot, popt, fs, multinom=True, log=logscale, eps=eps ) - # The uncertainty for theta is predicted, so we slice the - # the uncertainties for just the parameters. - uncerts_adj = uncerts_adj[:-1] + # Add theta into popt + popt = np.append(popt, np.array([theta])) + # param_string = '\t'.join(param_names[1:]) + # f.write(f"Description\t{param_string}\n") f.write( "Estimated 95% uncerts (theta adj), with step size " @@ -118,6 +123,7 @@ def godambe_stat_dfe( cache2d: str, sele_dist: str, sele_dist2: str, + pdf_file: str, output: str, bootstrap_syn_dir: str, bootstrap_non_dir: str, @@ -141,6 +147,8 @@ def godambe_stat_dfe( Name of the 1D PDF for modeling DFE. sele_dist2 : str Name of the 2D PDF for modeling DFE. + pdf_file : str + Name of file with custom probability density function model(s) in it. output : str File path for the output results. bootstrap_syn_dir : str @@ -171,11 +179,17 @@ def godambe_stat_dfe( fs = dadi.Spectrum.from_file(fs) non_fs_files = glob.glob(bootstrap_non_dir + "/*.fs") + # Raise error if directory is empty + if non_fs_files == []: + raise ValueError(f"ERROR:\n{bootstrap_dir} is empty\n") all_non_boot = [] for f in non_fs_files: boot_fs = dadi.Spectrum.from_file(f) all_non_boot.append(boot_fs) syn_fs_files = glob.glob(bootstrap_syn_dir + "/*.fs") + # Raise error if directory is empty + if syn_fs_files == []: + raise ValueError(f"ERROR:\n{bootstrap_dir} is empty\n") all_syn_boot = [] for f in syn_fs_files: boot_fs = dadi.Spectrum.from_file(f) @@ -189,9 +203,9 @@ def godambe_stat_dfe( sfunc = s2.integrate if sele_dist2 is not None: - sele_dist2 = get_dadi_pdf(sele_dist2) + sele_dist2, _ = get_dadi_pdf(sele_dist2, pdf_file) if sele_dist is not None: - sele_dist = get_dadi_pdf(sele_dist) + sele_dist, _ = get_dadi_pdf(sele_dist, pdf_file) else: sele_dist = sele_dist2 diff --git a/dadi_cli/parsers/common_arguments.py b/dadi_cli/parsers/common_arguments.py index 244a70c7..1a1b8573 100644 --- a/dadi_cli/parsers/common_arguments.py +++ b/dadi_cli/parsers/common_arguments.py @@ -39,7 +39,7 @@ def add_bounds_argument(parser: argparse.ArgumentParser) -> None: parser.add_argument( "--lbounds", - type=float, + type=str, nargs="+", required=boundary_req, help="Lower bounds of the optimized parameters.", @@ -47,7 +47,7 @@ def add_bounds_argument(parser: argparse.ArgumentParser) -> None: parser.add_argument( "--ubounds", - type=float, + type=str, nargs="+", required=boundary_req, help="Upper bounds of the optimized parameters.", @@ -191,7 +191,7 @@ def add_constant_argument(parser: argparse.ArgumentParser) -> None: parser.add_argument( "--constants", default=-1, - type=float, + type=str, nargs="+", help="Fixed parameters during the inference or using Godambe analysis. Use -1 to indicate a parameter is NOT fixed. Default: None.", ) @@ -315,6 +315,14 @@ def add_dfe_argument(parser: argparse.ArgumentParser) -> None: type=str, help="2D probability density function for the joint DFE inference. To check available probability density functions, please use `dadi-cli Pdf`.", ) + parser.add_argument( + "--pdf-file", + type=existed_file, + required=False, + dest="pdf_file", + help="Name of python probability density function module file (not including .py) that contains custom probability density functions to use. Default: None.", + ) + # Currently this command is only needed for param_names. diff --git a/dadi_cli/parsers/infer_dfe_parsers.py b/dadi_cli/parsers/infer_dfe_parsers.py index cb6749cc..5fc5322d 100644 --- a/dadi_cli/parsers/infer_dfe_parsers.py +++ b/dadi_cli/parsers/infer_dfe_parsers.py @@ -27,8 +27,8 @@ def _run_infer_dfe(args: argparse.Namespace) -> None: The 1D probability distribution function name. - pdf2d : str or None The 2D probability distribution function name. - - pdf_file : bool - Flag indicating whether a custom PDF file is used. + - pdf_file : str or None + Name of file with custom probability density function model(s) in it. - constants : list or int List of constants for the PDFs or -1 if not using constants. - lbounds : list or None @@ -80,6 +80,11 @@ def _run_infer_dfe(args: argparse.Namespace) -> None: The mixed PDF model if applicable. """ + # Make sure flags are used: + if args.pdf1d == None and args.pdf2d == None: + raise ValueError("Require --pdf1d and/or --pdf2d depending on DFE model") + if args.cache1d == None and args.cache2d == None: + raise ValueError("cache1d --pdf1d and/or --cache2d depending on DFE model") if "://" in args.fs: import urllib.request sfs_fi = open("sfs.fs","w") @@ -87,6 +92,13 @@ def _run_infer_dfe(args: argparse.Namespace) -> None: sfs_fi.write(f.read().decode('utf-8')) sfs_fi.close() args.fs ="sfs.fs" + if args.pdf_file is not None: + if "://" in args.pdf_file: + model_fi = open("dadi_pdfs.py","w") + with urllib.request.urlopen(args.pdf_file) as f: + model_fi.write(f.read().decode('utf-8')) + model_fi.close() + args.pdf_file = "dadi_pdfs" fs = dadi.Spectrum.from_file(args.fs) # Due to development history, much of the code expects a args.misid variable, so create it. @@ -97,6 +109,9 @@ def _run_infer_dfe(args: argparse.Namespace) -> None: make_dir(args.output_prefix) + # Converts str to float and None string to None value + args.lbounds, args.ubounds, args.constants = convert_bounds_and_constants(args.lbounds, args.ubounds, args.constants) + # # Things need to be updated for these to work if None not in [args.pdf1d, args.pdf2d]: pass @@ -107,11 +122,11 @@ def _run_infer_dfe(args: argparse.Namespace) -> None: args.constants, _ = check_pdf_params( args.constants, pdf, "--constant", args.misid ) - if not args.pdf_file and args.lbounds != -1: + if not args.pdf_file and args.lbounds != None: args.lbounds, _ = check_pdf_params( args.lbounds, pdf, "--lbounds", args.misid ) - if not args.pdf_file and args.ubounds != -1: + if not args.pdf_file and args.ubounds != None: args.ubounds, _ = check_pdf_params( args.ubounds, pdf, "--ubounds", args.misid ) @@ -151,7 +166,7 @@ def _run_infer_dfe(args: argparse.Namespace) -> None: else: cache2d = args.cache2d - if not np.all(cache_ns == fs.sample_sizes): + if not np.all(cache_ns == fs.sample_sizes) and args.cov_args == []: raise ValueError('Cache and frequencey spectrum do not have the same sample sizes') if args.maxeval == False: @@ -175,7 +190,12 @@ def _run_infer_dfe(args: argparse.Namespace) -> None: else: for pdf in [args.pdf1d, args.pdf2d]: if pdf != None: - _, param_names = check_pdf_params(args.p0, pdf, "", args.misid) + if args.pdf_file != None: + # try: + _, param_names = get_dadi_pdf(pdf, args.pdf_file) + # except: + else: + _, param_names = check_pdf_params(args.p0, pdf, "", args.misid) param_names = "# Log(likelihood)\t" + "\t".join(param_names) if args.misid: @@ -245,6 +265,7 @@ def _run_infer_dfe(args: argparse.Namespace) -> None: cache2d, args.pdf1d, args.pdf2d, + args.pdf_file, theta, args.p0, args.ubounds, @@ -270,6 +291,7 @@ def _run_infer_dfe(args: argparse.Namespace) -> None: cache2d, args.pdf1d, args.pdf2d, + args.pdf_file, theta, args.p0, args.ubounds, @@ -390,13 +412,7 @@ def add_infer_dfe_parsers(subparsers: argparse.ArgumentParser) -> None: required=True, help="Ratio for the nonsynonymous mutations to the synonymous mutations.", ) - parser.add_argument( - "--pdf-file", - type=str, - required=False, - dest="pdf_file", - help="Name of python probability density function module file (not including .py) that contains custom probability density functions to use. Default: None.", - ) + add_inference_argument(parser) add_delta_ll_argument(parser) add_misid_argument(parser) diff --git a/dadi_cli/parsers/infer_dm_parsers.py b/dadi_cli/parsers/infer_dm_parsers.py index 4fb77ae5..4982e4b8 100644 --- a/dadi_cli/parsers/infer_dm_parsers.py +++ b/dadi_cli/parsers/infer_dm_parsers.py @@ -91,6 +91,9 @@ def _run_infer_dm(args: argparse.Namespace) -> None: make_dir(args.output_prefix) + # Converts str to float and None string to None value + args.lbounds, args.ubounds, args.constants = convert_bounds_and_constants(args.lbounds, args.ubounds, args.constants) + if args.cov_args != []: import pickle args.cov_args[0] = pickle.load(open(args.cov_args[0], 'rb')) @@ -124,9 +127,9 @@ def _run_infer_dm(args: argparse.Namespace) -> None: args.constants = check_params( args.constants, args.model, "--constant", args.misid ) - if not args.model_file and args.lbounds != -1: + if not args.model_file and args.lbounds != None: args.lbounds = check_params(args.lbounds, args.model, "--lbounds", args.misid) - if not args.model_file and args.ubounds != -1: + if not args.model_file and args.ubounds != None: args.ubounds = check_params(args.ubounds, args.model, "--ubounds", args.misid) if args.p0 == -1: diff --git a/dadi_cli/parsers/plot_parsers.py b/dadi_cli/parsers/plot_parsers.py index 27bfed15..cfbdd162 100644 --- a/dadi_cli/parsers/plot_parsers.py +++ b/dadi_cli/parsers/plot_parsers.py @@ -33,6 +33,8 @@ def _run_plot(args: argparse.Namespace) -> None: The 1D probability density function name for DFE plotting. - pdf2 : str, optional The 2D probability density function name for DFE plotting. + - pdf_file : str or None + Name of file with custom probability density function model(s) in it. - cache1d : str, optional Path to a cache file for 1D computation. - cache2d : str, optional @@ -103,6 +105,7 @@ def _run_plot(args: argparse.Namespace) -> None: cache2d=args.cache2d, pdf=args.pdf1d, pdf2=args.pdf2d, + pdf_file=args.pdf_file, cov_args=args.cov_args, cov_inbreeding=args.cov_inbreeding, sele_popt=args.dfe_popt, @@ -192,13 +195,6 @@ def add_plot_parsers(subparsers: argparse.ArgumentParser) -> None: default=0.1, help="Minimum value to be plotted in the frequency spectrum, default: 0.1.", ) - parser.add_argument( - "--ratio", - type=positive_num, - default=False, - required=False, - help="Ratio for the nonsynonymous mutations to the synonymous mutations.", - ) parser.add_argument( "--interactive", default=False, diff --git a/dadi_cli/parsers/simulate_dfe_parsers.py b/dadi_cli/parsers/simulate_dfe_parsers.py index 0f52de31..c5d30ed8 100644 --- a/dadi_cli/parsers/simulate_dfe_parsers.py +++ b/dadi_cli/parsers/simulate_dfe_parsers.py @@ -29,6 +29,8 @@ def _run_simulate_dfe(args: argparse.Namespace) -> None: The 1D probability distribution function name for the DFE. - pdf2d : str, optional The 2D probability distribution function name for the DFE. + - pdf_file : str, optional + Name of file with custom probability density function model(s) in it. - ratio : float Ratio for adjusting the selection parameters, typically used to scale between different types of mutations or fitness effects. @@ -63,6 +65,7 @@ def _run_simulate_dfe(args: argparse.Namespace) -> None: cache2d=cache2d, sele_dist=args.pdf1d, sele_dist2=args.pdf2d, + pdf_file=args.pdf_file, ratio=args.ratio, misid=args.misid, output=args.output, diff --git a/dadi_cli/parsers/stat_dfe_parsers.py b/dadi_cli/parsers/stat_dfe_parsers.py index 339cf303..be9f4ebd 100644 --- a/dadi_cli/parsers/stat_dfe_parsers.py +++ b/dadi_cli/parsers/stat_dfe_parsers.py @@ -28,6 +28,8 @@ def _run_stat_dfe(args: argparse.Namespace) -> None: The 1D probability distribution function used for the DFE. - pdf2d : str The 2D probability distribution function used for the DFE. + - pdf_file : str + Name of file with custom probability density function model(s) in it. - dfe_popt : str Path to the file containing optimized parameter values for the DFE model. - constants : list @@ -67,6 +69,7 @@ def _run_stat_dfe(args: argparse.Namespace) -> None: cache2d=args.cache2d, sele_dist=args.pdf1d, sele_dist2=args.pdf2d, + pdf_file=args.pdf_file, dfe_popt=args.dfe_popt, fixed_params=args.constants, logscale=args.logscale, diff --git a/dadi_cli/parsers/stat_dm_parsers.py b/dadi_cli/parsers/stat_dm_parsers.py index 24b4051a..8e635bf9 100644 --- a/dadi_cli/parsers/stat_dm_parsers.py +++ b/dadi_cli/parsers/stat_dm_parsers.py @@ -106,6 +106,7 @@ def add_stat_dm_parsers(subparsers: argparse.ArgumentParser) -> None: parser.add_argument( "--demo-popt", type=existed_file, + required=True, dest="demo_popt", help="File contains the bestfit demographic parameters, generated by `dadi-cli BestFit`.", ) diff --git a/dadi_cli/utilities.py b/dadi_cli/utilities.py index 01a3887c..d3dc9425 100644 --- a/dadi_cli/utilities.py +++ b/dadi_cli/utilities.py @@ -1,6 +1,7 @@ import dadi, multiprocessing import numpy as np from typing import Optional +import ast def pts_l_func( @@ -78,6 +79,27 @@ def convert_to_None( return inference_input +def convert_bounds_and_constants(lbounds, ubounds, constants): + if lbounds != None: + try: + lbounds = [ast.literal_eval(ele) for ele in lbounds] + except ValueError: + raise ValueError("--lbounds must be int, float, or None") + + if ubounds != None: + try: + ubounds = [ast.literal_eval(ele) for ele in ubounds] + except ValueError: + raise ValueError("--ubounds must be int, float, or None") + + if constants != -1: + try: + constants = [ast.literal_eval(ele) for ele in constants] + except ValueError: + raise ValueError("--constants must be int, float, or None") + + return lbounds, ubounds, constants + def get_opts_and_theta( filename: str, gen_cache: bool = False, diff --git a/docs/index.md b/docs/index.md index 104702d0..98619e54 100644 --- a/docs/index.md +++ b/docs/index.md @@ -46,12 +46,13 @@ To get help information, users can use dadi-cli -h ``` -There are thirteen subcommands in `dadi-cli`: +There are fourteen subcommands in `dadi-cli`: - GenerateFs - GenerateCache - InferDM - InferDFE +- LowPass - BestFit - StatDM - StatDFE diff --git a/docs/userguide/lowpass.md b/docs/userguide/lowpass.md new file mode 100644 index 00000000..33d10beb --- /dev/null +++ b/docs/userguide/lowpass.md @@ -0,0 +1,41 @@ +# dadi.LowPass Integration + +dadi-cli can utilize the [Low-Pass module of dadi](https://dadi.readthedocs.io/en/latest/user-guide/low-pass/). This enables the ability to correct models based on the coverage each sample has. + +## GenerateFs Input + +Before you use Low-Pass you'll want to make sure your VCF has the allele depth (AD) entry. +
+FORMAT
+GT:AD:DP:GQ:PL
+
+The information from the AD entry is used when generating a demographic model. Users can save the information in a pickled dictionary when running `GenerateFs` with the `--calc-coverage`: +```python +# Generate synonymous SFS +dadi-cli GenerateFs --vcf examples/data/mus.syn.subset.vcf.gz --pop-info examples/data/mouse.popfile.txt --pop-ids Mmd_IRA Mmd_FRA --projections 4 8 --polarized --output examples/results/lowpass/mus.syn.fs --calc-coverage +# Generate nonsynonymous SFS +dadi-cli GenerateFs --vcf examples/data/mus.nonsyn.subset.vcf.gz --pop-info examples/data/mouse.popfile.txt --pop-ids Mmd_IRA Mmd_FRA --projections 4 8 --polarized --output examples/results/lowpass/mus.nonsyn.fs --calc-coverage +``` +In addition to the `mus.syn.fs` and `mus.nonsyn.fs` SFS files, `mus.syn.fs.coverage.pickle` and `mus.nonsyn.fs.coverage.pickle` are generated. These files will be used for `InferDM` and `InferDFE` to preform the Low-Pass model correction. + +## InferDM with Low-Pass + +Low-Pass will correct the demographic and DFE models based on the coverage of samples, to use it users will need to use the `--coverage-model` flag, which takes the `.coverage.pickle` file and total number of haplotypes in the data for each population. There were 10 Mmd_IRA samples and 16 Mmd_FRA samples being requested from the `--pop-info` file from `GenerateFs`. So, to add the correction to the demographic mode, users would include `--coverage-model examples/results/lowpass/mus.syn.fs.coverage.pickle 10 16` in their `InferDM` command. +```python +dadi-cli InferDM --fs examples/results/lowpass/mus.syn.fs --model split_mig --p0 4.5 0.8 0.8 0.36 0.01 --lbounds 1e-5 1e-5 1e-5 1e-5 1e-5 --ubounds 10 10 1 10 1 --output-prefix examples/results/lowpass/mus.split_mig.lowpass --optimizations 20 --check-convergence 5 --grids 50 60 70 --coverage-model examples/results/lowpass/mus.syn.fs.coverage.pickle 10 16 +``` + +There are no extra parameters for the Low-Pass model correction, so results can look similar to non-Low-Pass corrected results. + +## InferDFE with Low-Pass + +InferDFE is largely the same as InferDM. When users run `GenerateCache`, user's will want to use the number of haplotypes in the data for `--sample-sizes` instead of the sample sizes in the SFS file. Given that there were 10 Mmd_IRA samples and 16 Mmd_FRA haplotypes users would use `--sample-sizes 10 16` instead of `--sample-sizes 4 8`: + +```python +dadi-cli GenerateCache --model split_mig_sel_single_gamma --sample-sizes 10 16 --grids 50 60 70 --gamma-pts 20 --gamma-bounds 1e-4 100 --demo-popt examples/results/lowpass/mus.split_mig.lowpass.InferDM.bestfits --output examples/results/lowpass/mus.1d.cache +``` + +When users run InferDFE they should use the `--coverage-model` flag with the similar inputs as `InferDM`, except using the `.coverage.pickle` file for the nonsynonymous data. +```python +dadi-cli InferDFE --fs examples/results/lowpass/mus.nonsyn.fs --demo-popt examples/results/lowpass/mus.split_mig.lowpass.InferDM.bestfits --cache1d examples/results/lowpass/mus.1d.cache --pdf1d lognormal --ratio 2.31 --p0 1 1 0.01 --lbounds 1e-5 1e-5 1e-5 --ubounds 30 10 1 --optimizations 20 --check-convergence 5 --output-prefix examples/results/lowpass/mus.lognormal.lowpass --coverage-model examples/results/lowpass/mus.nonsyn.fs.coverage.pickle 10 16 +``` \ No newline at end of file diff --git a/examples/data/mouse.popfile.txt b/examples/data/mouse.popfile.txt new file mode 100644 index 00000000..b4ab2cc7 --- /dev/null +++ b/examples/data/mouse.popfile.txt @@ -0,0 +1,13 @@ +Mmd_IRA1_AH15.variant3 Mmd_IRA +Mmd_IRA2_AH23.variant3 Mmd_IRA +Mmd_IRA5_JR2-F1C.variant3 Mmd_IRA +Mmd_IRA7_JR7-F1C.variant3 Mmd_IRA +Mmd_IRA8_JR8-F1A.variant3 Mmd_IRA +Mmd_FRA1_14.variant Mmd_FRA +Mmd_FRA2_15B.variant Mmd_FRA +Mmd_FRA3_16B.variant Mmd_FRA +Mmd_FRA4_18B.variant Mmd_FRA +Mmd_FRA5_B2C.variant Mmd_FRA +Mmd_FRA6_C1.variant Mmd_FRA +Mmd_FRA7_E1.variant Mmd_FRA +Mmd_FRA8_F1B.variant Mmd_FRA \ No newline at end of file diff --git a/examples/data/mus.nonsyn.subset.vcf.gz b/examples/data/mus.nonsyn.subset.vcf.gz new file mode 100644 index 00000000..c124130c Binary files /dev/null and b/examples/data/mus.nonsyn.subset.vcf.gz differ diff --git a/examples/data/mus.syn.subset.vcf.gz b/examples/data/mus.syn.subset.vcf.gz new file mode 100644 index 00000000..d4be2289 Binary files /dev/null and b/examples/data/mus.syn.subset.vcf.gz differ diff --git a/mkdocs.yml b/mkdocs.yml index 34fdba36..6a4c4fe3 100644 --- a/mkdocs.yml +++ b/mkdocs.yml @@ -33,6 +33,7 @@ nav: - Demographic inference: 'userguide/demog.md' - DFE inference: 'userguide/dfe.md' - Joint DFE inference: 'userguide/jdfe.md' + - LowPass sequence coverage correction: 'userguide/lowpass.md' - Statistical testing: 'userguide/stat.md' - Cloud computing: 'userguide/cloud.md' - Simulation: 'userguide/simulation.md' diff --git a/pyproject.toml b/pyproject.toml index d7fdbfc2..0ed0f1ca 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -4,7 +4,7 @@ build-backend = "hatchling.build" [project] name = "dadi-cli" -version = "0.9.7" +version = "0.9.9" requires-python = ">=3.9" description = "A command line interface for dadi" readme = "README.md" diff --git a/tests/example_data/LowPass-files/cov.fs.coverage.pickle b/tests/example_data/LowPass-files/cov.fs.coverage.pickle new file mode 100644 index 00000000..19561ec2 Binary files /dev/null and b/tests/example_data/LowPass-files/cov.fs.coverage.pickle differ diff --git a/tests/test_InferDFE.py b/tests/test_InferDFE.py index 7451749b..6869f26d 100644 --- a/tests/test_InferDFE.py +++ b/tests/test_InferDFE.py @@ -22,6 +22,7 @@ def test_infer_dfe_code(): cache2d = None sele_dist = "lognormal" sele_dist2 = None + pdf_file = None theta = 1000 * 2.31 p0 = [1, 1] upper_bounds = [10, 10] @@ -39,6 +40,7 @@ def test_infer_dfe_code(): cache2d, sele_dist, sele_dist2, + pdf_file, theta, p0, upper_bounds, @@ -59,6 +61,7 @@ def test_infer_dfe_code(): cache2d=None, sele_dist=sele_dist, sele_dist2=sele_dist2, + pdf_file=pdf_file, theta=theta, p0=p0, upper_bounds=upper_bounds, @@ -80,6 +83,7 @@ def test_infer_coverage_model_dfe_code(): cache2d = None sele_dist = "lognormal" sele_dist2 = None + pdf_file = None theta = 1000 * 2.31 p0 = [1, 1] upper_bounds = [10, 10] @@ -97,6 +101,7 @@ def test_infer_coverage_model_dfe_code(): cache2d, sele_dist, sele_dist2, + pdf_file, theta, p0, upper_bounds, diff --git a/tests/test_Pdfs.py b/tests/test_Pdfs.py index ed620643..c30c2bf4 100644 --- a/tests/test_Pdfs.py +++ b/tests/test_Pdfs.py @@ -11,14 +11,14 @@ def test_get_dadi_pdf(): - assert get_dadi_pdf("beta") == DFE.PDFs.beta - assert get_dadi_pdf("biv_ind_gamma") == DFE.PDFs.biv_ind_gamma - assert get_dadi_pdf("biv_lognormal") == DFE.PDFs.biv_lognormal - assert get_dadi_pdf("exponential") == DFE.PDFs.exponential - assert get_dadi_pdf("gamma") == DFE.PDFs.gamma - assert get_dadi_pdf("lognormal") == DFE.PDFs.lognormal - assert get_dadi_pdf("normal") == DFE.PDFs.normal - assert get_dadi_pdf("mixture") == DFE.mixture + assert get_dadi_pdf("beta")[0] == DFE.PDFs.beta + assert get_dadi_pdf("biv_ind_gamma")[0] == DFE.PDFs.biv_ind_gamma + assert get_dadi_pdf("biv_lognormal")[0] == DFE.PDFs.biv_lognormal + assert get_dadi_pdf("exponential")[0] == DFE.PDFs.exponential + assert get_dadi_pdf("gamma")[0] == DFE.PDFs.gamma + assert get_dadi_pdf("lognormal")[0] == DFE.PDFs.lognormal + assert get_dadi_pdf("normal")[0] == DFE.PDFs.normal + assert get_dadi_pdf("mixture")[0] == DFE.mixture with pytest.raises(Exception) as e_info: get_dadi_pdf("haha") diff --git a/tests/test_Plot.py b/tests/test_Plot.py index 0d6fb729..1adeaf85 100644 --- a/tests/test_Plot.py +++ b/tests/test_Plot.py @@ -206,6 +206,7 @@ def test_plot_fitted_dfe(files): projections=[8, 8], pdf="lognormal", pdf2="biv_lognormal", + pdf_file=None, cov_args=[], cov_inbreeding=[], output="tests/test_results/plot_fitted_dfe_mixture_w_proj.png", @@ -222,6 +223,7 @@ def test_plot_fitted_dfe(files): projections=[8, 8], pdf="lognormal", pdf2="biv_lognormal", + pdf_file=None, cov_args=[pickle.load(open("./tests/example_data/LowPass-files/cov.fs.coverage.pickle", 'rb')), 20, 20], cov_inbreeding=[], output="tests/test_results/plot_fitted_dfe_mixture_w_proj.png", diff --git a/tests/test_SimulateFs.py b/tests/test_SimulateFs.py index 810d3a94..3b1ce5c6 100644 --- a/tests/test_SimulateFs.py +++ b/tests/test_SimulateFs.py @@ -104,6 +104,7 @@ def test_simulate_dfe_code(): cache2d = pickle.load(open("tests/example_data/cache_split_mig_2d.bpkl", "rb")) sele_dist = "lognormal" sele_dist2 = "biv_lognormal" + pdf_file = None theta_ns = 2.31 misid = False output_1d = "tests/test_results/simulate_dfe_split_mig_one_s_lognormal.fs" @@ -114,20 +115,22 @@ def test_simulate_dfe_code(): ) # Test 1D - simulate_dfe(p0_1d, cache1d, None, sele_dist, None, theta_ns, misid, output_1d) + simulate_dfe(p0_1d, cache1d, None, sele_dist, None, pdf_file, theta_ns, misid, output_1d) dadi_cli_fs = dadi.Spectrum.from_file(output_1d) - dadi_fs = cache1d.integrate(p0_1d, None, get_dadi_pdf(sele_dist), theta_ns) + func, _ = get_dadi_pdf(sele_dist, pdf_file) + dadi_fs = cache1d.integrate(p0_1d, None, get_dadi_pdf(sele_dist, pdf_file)[0], theta_ns) assert np.allclose(dadi_cli_fs, dadi_fs) # Test 2D - simulate_dfe(p0_2d, None, cache2d, None, sele_dist2, theta_ns, misid, output_2d) + simulate_dfe(p0_2d, None, cache2d, None, sele_dist2, pdf_file, theta_ns, misid, output_2d) dadi_cli_fs = dadi.Spectrum.from_file(output_2d) - dadi_fs = cache2d.integrate(p0_2d, None, get_dadi_pdf(sele_dist2), theta_ns, None) + func, _ = get_dadi_pdf(sele_dist2, pdf_file) + dadi_fs = cache2d.integrate(p0_2d, None, get_dadi_pdf(sele_dist2, pdf_file)[0], theta_ns, None) assert np.allclose(dadi_cli_fs, dadi_fs) # Test mix simulate_dfe( - p0_mix, cache1d, cache2d, sele_dist, sele_dist2, theta_ns, misid, output_mix + p0_mix, cache1d, cache2d, sele_dist, sele_dist2, pdf_file, theta_ns, misid, output_mix ) dadi_cli_fs = dadi.Spectrum.from_file(output_mix) dadi_fs = DFE.mixture( @@ -135,8 +138,8 @@ def test_simulate_dfe_code(): None, cache1d, cache2d, - get_dadi_pdf(sele_dist), - get_dadi_pdf(sele_dist2), + get_dadi_pdf(sele_dist)[0], + get_dadi_pdf(sele_dist2)[0], theta_ns, None, ) @@ -144,11 +147,11 @@ def test_simulate_dfe_code(): # Test 1D misid simulate_dfe( - p0_1d_misid, cache1d, None, sele_dist, None, theta_ns, True, output_1d_misid + p0_1d_misid, cache1d, None, sele_dist, None, pdf_file, theta_ns, True, output_1d_misid ) dadi_cli_fs = dadi.Spectrum.from_file(output_1d_misid) dadi_fs = dadi.Numerics.make_anc_state_misid_func(cache1d.integrate)( - p0_1d_misid, None, get_dadi_pdf(sele_dist), theta_ns + p0_1d_misid, None, get_dadi_pdf(sele_dist)[0], theta_ns ) assert np.allclose(dadi_cli_fs, dadi_fs) diff --git a/tests/test_Stat.py b/tests/test_Stat.py index f3e0b748..4606409f 100644 --- a/tests/test_Stat.py +++ b/tests/test_Stat.py @@ -70,6 +70,7 @@ def test_StatDFE_code_mix_lognormal(): cache2d = "tests/example_data/cache_split_mig_2d.bpkl" sele_dist1 = "lognormal" sele_dist2 = "biv_lognormal" + pdf_file = None output = "./tests/test_results/example.split_mig.dfe.mixture_lognormal.params.ci" bootstrap_syn_dir = "tests/example_data/split_mig_bootstrap_syn/" bootstrap_non_dir = "tests/example_data/split_mig_bootstrap_non/bootstrap_non_mix/" @@ -83,6 +84,7 @@ def test_StatDFE_code_mix_lognormal(): cache2d, sele_dist1, sele_dist2, + pdf_file, output, bootstrap_syn_dir, bootstrap_non_dir, @@ -99,6 +101,7 @@ def test_StatDFE_code_1d_lognormal(): cache2d = None sele_dist1 = "lognormal" sele_dist2 = None + pdf_file = None output = "./tests/test_results/example.split_mig.dfe.1D_lognormal.params.ci" bootstrap_syn_dir = "tests/example_data/split_mig_bootstrap_syn/" bootstrap_non_dir = "tests/example_data/split_mig_bootstrap_non/bootstrap_non_1d/" @@ -114,6 +117,7 @@ def test_StatDFE_code_1d_lognormal(): cache2d, sele_dist1, sele_dist2, + pdf_file, output, bootstrap_syn_dir, bootstrap_non_dir, @@ -130,6 +134,7 @@ def test_StatDFE_code_2d_lognormal(): cache2d = "tests/example_data/cache_split_mig_2d.bpkl" sele_dist1 = None sele_dist2 = "biv_lognormal" + pdf_file = None output = "./tests/test_results/example.split_mig.dfe.2D_lognormal.params.ci" bootstrap_syn_dir = "tests/example_data/split_mig_bootstrap_syn/" bootstrap_non_dir = "tests/example_data/split_mig_bootstrap_non/bootstrap_non_2d/" @@ -145,6 +150,7 @@ def test_StatDFE_code_2d_lognormal(): cache2d, sele_dist1, sele_dist2, + pdf_file, output, bootstrap_syn_dir, bootstrap_non_dir, diff --git a/tests/test_cli.py b/tests/test_cli.py index e3ef5f19..4ff7f1c7 100644 --- a/tests/test_cli.py +++ b/tests/test_cli.py @@ -1,6 +1,6 @@ import pytest import dadi_cli.__main__ as cli -import os +import os, ast try: @@ -310,8 +310,8 @@ def test_infer_dm_args(model, model_file, nomisid): assert args.grids == grids assert args.nomisid == False assert args.constants == -1 - assert args.lbounds == lbounds - assert args.ubounds == ubounds + assert [ast.literal_eval(ele) for ele in args.lbounds] == lbounds + assert [ast.literal_eval(ele) for ele in args.ubounds] == ubounds assert args.global_optimization == True assert args.seed == seed @@ -334,9 +334,9 @@ def test_infer_dfe_args(): pdf_file = False p0 = [1, 1, 0, 0.001] grids = [120, 140, 160] - ubounds = [5, 5, -1, 0.999] - lbounds = [1e-4, 1e-4, -1, 1e-4] - constants = [-1, -1, 0, -1] + ubounds = [5, 5, None, 0.999] + lbounds = [1e-4, 1e-4, None, 1e-4] + constants = [None, None, 0, None] nomisid = True cuda = False maxeval = 100 @@ -399,11 +399,11 @@ def test_infer_dfe_args(): str(delta_ll), "--nomisid", "--constants", - "-1", "-1", "0", "-1", + "None", "None", "0", "None", "--lbounds", - "1e-4", "1e-4", "-1", "1e-4", + "1e-4", "1e-4", "None", "1e-4", "--ubounds", - "5", "5", "-1", "0.999", + "5", "5", "None", "0.999", "--seed", str(seed), ] @@ -423,9 +423,9 @@ def test_infer_dfe_args(): assert args.delta_ll == delta_ll assert args.gpus == gpus assert args.nomisid == True - assert args.constants == constants - assert args.lbounds == lbounds - assert args.ubounds == ubounds + assert [ast.literal_eval(ele) for ele in args.constants] == constants + assert [ast.literal_eval(ele) for ele in args.lbounds] == lbounds + assert [ast.literal_eval(ele) for ele in args.ubounds] == ubounds assert args.seed == seed assert args.port == 9123 assert args.pdf_file == None @@ -471,7 +471,6 @@ def test_plot_args(): projections = [10, 10] resid_range = 3 vmin = 1e-5 - ratio = 2.31 args = parser.parse_args( [ @@ -502,8 +501,6 @@ def test_plot_args(): str(resid_range), "--vmin", str(vmin), - "--ratio", - str(ratio) ] ) @@ -520,7 +517,6 @@ def test_plot_args(): assert args.projections == projections assert args.resid_range == resid_range assert args.vmin == vmin - assert args.ratio == ratio def test_stat_dm_args(): @@ -578,7 +574,7 @@ def test_stat_dfe_args(): bootstrapping_synonymous_dir = "tests/example_data/split_mig_bootstrap_syn/" bootstrapping_nonsynonymous_dir = "tests/example_data/split_mig_bootstrap_non/bootstrap_non_mix/" dfe_popt = "tests/example_data/example.split_mig.dfe.lognormal_mixture.params.with.misid.InferDFE.bestfits" - constants = [-1, -1, 0, -1, -1] + constants = [None, None, 0, None, None] logscale = False args = parser.parse_args( @@ -597,7 +593,7 @@ def test_stat_dfe_args(): "--output", output, "--constants", - "-1", "-1", "0", "-1", "-1", + "None", "None", "0", "None", "None", "--dfe-popt", dfe_popt, "--bootstrapping-synonymous-dir", @@ -613,7 +609,7 @@ def test_stat_dfe_args(): assert args.pdf1d == pdf1d assert args.pdf2d == pdf2d assert args.output == output - assert args.constants == constants + assert [ast.literal_eval(ele) for ele in args.constants] == constants assert args.dfe_popt == dfe_popt assert args.bootstrapping_syn_dir == bootstrapping_synonymous_dir assert args.bootstrapping_non_dir == bootstrapping_nonsynonymous_dir diff --git a/tests/test_main.py b/tests/test_main.py index e019db35..372d2001 100644 --- a/tests/test_main.py +++ b/tests/test_main.py @@ -141,6 +141,7 @@ def simulate_args(): simulate_args.cache2d = "tests/example_data/cache_split_mig_2d.bpkl" simulate_args.pdf1d = "lognormal" simulate_args.pdf2d = "biv_lognormal" + simulate_args.pdf_file = None simulate_args.ratio = 2.31 simulate_args.nomisid = False simulate_args.output = "tests/test_results/main_simulate_mix_dfe.fs" @@ -157,6 +158,7 @@ def simulate_args(): simulate_args.cache2d = "https://github.com/xin-huang/dadi-cli/blob/master/tests/example_data/cache_split_mig_2d.bpkl?raw=true" simulate_args.pdf1d = "lognormal" simulate_args.pdf2d = "biv_lognormal" + simulate_args.pdf_file = None simulate_args.ratio = 2.31 simulate_args.nomisid = False simulate_args.output = "tests/test_results/main_simulate_mix_dfe_html.fs" @@ -212,8 +214,8 @@ def infer_dm_args(): pytest.model_file = "tests/example_data/example_models" pytest.p0 = [1, 0.5] pytest.grids = [30, 40, 50] - pytest.ubounds = [10, 10] - pytest.lbounds = [1e-3, 1e-3] + pytest.ubounds = ["10", "10"] + pytest.lbounds = ["1e-3", "1e-3"] pytest.constants = -1 pytest.nomisid = True pytest.cov_args = [] @@ -273,9 +275,9 @@ def test_run_infer_dm_convergence(check_convergence, force_convergence, optimiza def test_run_infer_no_custom_model_dm(infer_dm_args): pytest.model = "three_epoch_inbreeding" pytest.model_file = None - pytest.constants = [None, 1, None, None, None] - pytest.ubounds = [10, 10, 10, 10, 1] - pytest.lbounds = [1e-3, 1e-3, 1e-3, 1e-3, 1e-5] + pytest.constants = ["None", "1", "None", "None", "None"] + pytest.ubounds = ["10", "10", "10", "10", "1"] + pytest.lbounds = ["1e-3", "1e-3", "1e-3", "1e-3", "1e-5"] pytest.p0 = -1 pytest.maxeval = False pytest.inbreeding = True @@ -290,8 +292,8 @@ def test_run_infer_dm_misid(infer_dm_args): pytest.nomisid = False pytest.output_prefix = "tests/test_results/main.test.two_epoch.demo_misid.params" pytest.p0 = [1, 0.5, 1e-2] - pytest.ubounds = [10, 10, 0.999] - pytest.lbounds = [1e-3, 1e-3, 1e-4] + pytest.ubounds = ["10", "10", "0.999"] + pytest.lbounds = ["1e-3", "1e-3", "1e-4"] _run_infer_dm(pytest) for ele in glob.glob(pytest.output_prefix+"*"): os.remove(ele) @@ -304,8 +306,8 @@ def test_run_infer_dm_lowpass(infer_dm_args): pytest.model_file = None pytest.output_prefix = "tests/test_results/main.test.split_mig.demo_lowpass.params" pytest.p0 = [1, 1, 0.1, 1] - pytest.ubounds = [10, 10, 1, 10] - pytest.lbounds = [1e-3, 1e-3, 1e-3, 1e-3] + pytest.ubounds = ["10", "10", "1", "10"] + pytest.lbounds = ["1e-3", "1e-3", "1e-3", "1e-3"] _run_infer_dm(pytest) for ele in glob.glob(pytest.output_prefix+"*"): os.remove(ele) @@ -317,8 +319,8 @@ def test_run_infer_dm_global_bestfit(infer_dm_args): pytest.output_prefix = "tests/test_results/main.test.two_epoch.demo_bestfit_p0.params" pytest.nomisid = False pytest.p0 = [1, 0.5, 1e-2] - pytest.ubounds = [10, 10, 0.999] - pytest.lbounds = [1e-3, 1e-3, 1e-4] + pytest.ubounds = ["10", "10", "0.999"] + pytest.lbounds = ["1e-3", "1e-3", "1e-4"] pytest.maxeval = 10 _run_infer_dm(pytest) for ele in glob.glob(pytest.output_prefix+"*"): @@ -343,8 +345,8 @@ def test_run_infer_dm_html(infer_dm_args): pytest.model_file = "https://raw.githubusercontent.com/xin-huang/dadi-cli/master/tests/example_data/example_models.py" pytest.bestfit_p0 = "https://raw.githubusercontent.com/xin-huang/dadi-cli/master/tests/example_data/example.bestfit.two_epoch.demo.params.InferDM.opts.0" pytest.nomisid = False - pytest.ubounds = [10, 10, 0.999] - pytest.lbounds = [1e-3, 1e-3, 1e-4] + pytest.ubounds = ["10", "10", "0.999"] + pytest.lbounds = ["1e-3", "1e-3", "1e-4"] pytest.output_prefix += ".html" _run_infer_dm(pytest) for ele in glob.glob(pytest.output_prefix+"*"): @@ -363,13 +365,13 @@ def infer_dfe_args(): pytest.cache2d = "tests/example_data/cache_split_mig_2d.bpkl" pytest.pdf1d = "lognormal" pytest.pdf2d = "biv_lognormal" + pytest.pdf_file = None pytest.mix_pdf = False pytest.ratio = 2.31 - pytest.pdf_file = False pytest.p0 = [1, 1] pytest.grids = [120, 140, 160] - pytest.ubounds = [10, 10] - pytest.lbounds = [1e-3, 1e-3] + pytest.ubounds = ["10", "10"] + pytest.lbounds = ["1e-3", "1e-3"] pytest.constants = -1 pytest.nomisid = True pytest.cov_args = [] @@ -398,6 +400,7 @@ def test_run_infer_dfe_1d(infer_dfe_args): pytest.cache2d = None pytest.mix_pdf = None pytest.output_prefix += "1d_lognormal_dfe" + print('pdf:',pytest.pdf_file) _run_infer_dfe(pytest) for ele in glob.glob(pytest.output_prefix+"*"): os.remove(ele) @@ -455,8 +458,8 @@ def test_run_infer_dfe_2d_lognormal(infer_dfe_args): pytest.mix_pdf = None pytest.output_prefix += "2d_lognormal_dfe" pytest.p0 = [1, 1, 0.5] - pytest.ubounds = [10, 10, 0.999] - pytest.lbounds = [1e-3, 1e-3, 1e-3] + pytest.ubounds = ["10", "10", "0.999"] + pytest.lbounds = ["1e-3", "1e-3", "1e-3"] _run_infer_dfe(pytest) for ele in glob.glob(pytest.output_prefix+"*"): @@ -468,9 +471,9 @@ def test_run_infer_dfe_mix(infer_dfe_args): pytest.mix_pdf = 'mixture_lognormal' pytest.output_prefix += "mix_lognormal_dfe" pytest.p0 = [1, 1, 0, 0.5] - pytest.ubounds = [10, 10, None, 0.999] - pytest.lbounds = [1e-3, 1e-3, None, 1e-3] - pytest.constants = [None, None, 0, None] + pytest.ubounds = ["10", "10", "None", "0.999"] + pytest.lbounds = ["1e-3", "1e-3", "None", "1e-3"] + pytest.constants = ["None", "None", "0", "None"] _run_infer_dfe(pytest) fids = glob.glob(pytest.output_prefix+"*") @@ -493,9 +496,9 @@ def test_run_infer_dfe_mix_html(infer_dfe_args): pytest.output_prefix += "mix_lognormal_dfe_html" pytest.nomisid = False # pytest.p0 = None - pytest.ubounds = [10, 10, None, 0.999, 0.4] - pytest.lbounds = [1e-3, 1e-3, None, 1e-3, 1e-3] - pytest.constants = [None, None, 0, None, None] + pytest.ubounds = ["10", "10", "None", "0.999", "0.4"] + pytest.lbounds = ["1e-3", "1e-3", "None", "1e-3", "1e-3"] + pytest.constants = ["None", "None", "0", "None", "None"] _run_infer_dfe(pytest) fids = glob.glob(pytest.output_prefix+"*") @@ -561,6 +564,7 @@ def test_run_infer_dm_workqueue(infer_dm_args): ) pytest.output_prefix = "tests/test_results/main.test.two_epoch.demo_wq.params" pytest.work_queue = ['pytest-dadi-cli', 'tests/mypwfile'] + pytest.port = 9123 _run_infer_dm(pytest) factory.kill() @@ -597,6 +601,7 @@ def stat_args(): stat_args.cache2d=None stat_args.pdf1d="lognormal" stat_args.pdf2d=None + stat_args.pdf_file = None stat_args.bootstrapping_syn_dir="tests/example_data/split_mig_bootstrap_syn/" stat_args.bootstrapping_non_dir="tests/example_data/split_mig_bootstrap_non/bootstrap_non_1d/" stat_args.grids=[30, 40, 50] @@ -618,6 +623,7 @@ def plot_args(): # pytest.fs2_1d = "tests/example_data/two_epoch_non.fs" pytest.pdf1d = "lognormal" pytest.pdf2d = "biv_lognormal" + pytest.pdf_file = None pytest.output = "tests/test_results/main_test_plot.pdf" # pytest.demo_popt = "tests/example_data/example.two_epoch.demo.params.InferDM.bestfits" # pytest.fs1d_cache1d = "tests/example_data/cache_two_epoch_1d.bpkl"