-
Notifications
You must be signed in to change notification settings - Fork 0
/
Admix.py
62 lines (44 loc) · 1.97 KB
/
Admix.py
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
57
58
59
60
61
# Tests for admixture times and proportions.
from __future__ import division, print_function
import numpy as np
import sys, getopt
import Admix_estimate, snp_data
from parse import parse_pops
###########################################################################
def parse_options():
"""
data: Root of genotype data in eigenstrat format, i.e. root{.geno .snp .ind}
1,2, and 3: Corresponding to populations A,B and C
inbred: Comma separated list of pops that might be inbred
"""
options ={ "data":"", "1":"", "2":"", "3":"", "inbred":[],
"out":"out"}
try:
opts, args = getopt.getopt(sys.argv[1:], "d:1:2:3:i:o:", # @UnusedVariable
["data", "1", "2", "3", "inbred", "out"])
except Exception as err:
print(str(err))
sys.exit()
for o, a in opts:
if o in ["-d","--data"]: options["data"] = a
elif o in ["-o","--out"]: options["out"] = a
elif o in ["-i","--inbred"]: options["inbred"] = parse_pops(a)
elif o in ["-1","--1"]: options["1"] = a
elif o in ["-2","--2"]: options["2"] = a
elif o in ["-3","--3"]: options["3"] = a
options["pops"]=[options[x] for x in ["1", "2", "3"]]
print("found options:", file=sys.stderr)
print(options, file=sys.stderr)
return options
###########################################################################
def main(options):
# Load population data - TODO: Move all this to an eigenstrat class
data=snp_data.eigenstrat_data(options["data"], options["pops"],
True, options["inbred"], 0)
test=Admix_estimate.patterson_estimator(data, np.array([options[x] for x in ["1", "2", "3"]]))
test.estimate_admix()
return
###########################################################################
if __name__=="__main__":
options=parse_options()
main(options)