-
Notifications
You must be signed in to change notification settings - Fork 0
Expand file tree
/
Copy pathmain_config_snp_parallel_automated_mut_rate.py
More file actions
56 lines (45 loc) · 1.83 KB
/
Copy pathmain_config_snp_parallel_automated_mut_rate.py
File metadata and controls
56 lines (45 loc) · 1.83 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
import subprocess
import argparse
import yaml
def main():
parser = argparse.ArgumentParser(description="Run commands in series with automatic mutation rate configuration.")
parser.add_argument("-t", required=True, help="Type parameter for the second command.")
parser.add_argument("-f", required=True, help="File path for the input probability matrix.")
parser.add_argument("-p", type=float, required=True, help="Probability parameter for the first command.")
parser.add_argument("-m", type=int, default=1000, help="Number of histories used to estimate mut rate.")
parser.add_argument("-l", type=int, required=True, help="Length parameter for the second command.")
args = parser.parse_args()
# Run the first command
print("\n")
print("AUTOMATED ESTIMATION OF MUTATION RATE\n")
print("Using gaussian approximation of tree length distribution to obtain mutation rate..\n")
cmd1 = [
"python", "scr/main_est_mutation_rate_oneshot.py",
"-a", "gauss",
"-p", str(args.p),
"-m", str(args.m),
"-t", args.t,
"-f", args.f,
"-rs", "1"
]
result1 = subprocess.run(cmd1, capture_output=True, text=True, check=True)
# Parse the output of the first command
output_lines = result1.stdout.strip().split("\n")
mut_rate, error_bar = map(float, output_lines[-1].split())
print(f"Estimated mutation rate: {mut_rate}\n")
print(f"Estimate error in the fraction of invariant sites: {error_bar}\n")
# Run the second command
print("GENERATING THE SNP CONFIGURATION")
cmd2 = [
"caffeinate", "recopops", "snps",
"-t", args.t,
"-f", args.f,
"-rs", "1",
"-u", str(mut_rate),
"-v", "0",
"-l", str(args.l)
]
subprocess.run(cmd2, check=True)
print("\n")
if __name__ == "__main__":
main()