Add HAB encoding

2 files changed, 34 insertions, 10 deletions

diff --git a/wqflask/wqflask/my_pylmm/pyLMM/convertlmm.py b/wqflask/wqflask/my_pylmm/pyLMM/convertlmm.py
index f35cc47c..8a1f03ad 100644
--- a/wqflask/wqflask/my_pylmm/pyLMM/convertlmm.py
+++ b/wqflask/wqflask/my_pylmm/pyLMM/convertlmm.py
@@ -148,20 +148,26 @@ if options.geno:
     wrln("# Genotype format version 1.0")
     wrln("# Individuals = "+str(num_inds))
     wrln("# Phenotypes = "+str(len(num_snps)))
-    for i in range(len(num_snps)):
+    for i in range(num_inds):
         wr("\t"+str(i+1))
-        wr("\n")
+    wr("\n")
 
     m = []
     def out(i,x):
-        wr(str(i)+"\t")
-        wr("\t".join(x))
-        wr("\n")
+        # wr(str(i)+"\t")
+        # wr("\t".join(x))
+        # wr("\n")
+        m.append(x)
         
-    snps = plink.readbed(options.geno+'.bed',num_inds, ('A','H','B','NA'), out)
-    # for i in range(count):
-    #     wr("\t".join(phenos[i]))
-    #     wr("\n")
+    snps = plink.readbed(options.geno+'.bed',num_inds, ('A','H','B','-'), out)
+
+    msg("Write transposed genotype matrix")
+    for g in range(len(num_snps)):
+        wr(str(g+1)+"\t")
+        for i in range(num_inds):
+            wr(m[g][i])
+        wr("\n")
+            
     msg(str(count)+" geno lines written (with "+str(snps)+" snps)")
    
 msg("Converting done")
diff --git a/wqflask/wqflask/my_pylmm/pyLMM/runlmm.py b/wqflask/wqflask/my_pylmm/pyLMM/runlmm.py
index 27afd6f4..4398926f 100644
--- a/wqflask/wqflask/my_pylmm/pyLMM/runlmm.py
+++ b/wqflask/wqflask/my_pylmm/pyLMM/runlmm.py
@@ -44,6 +44,8 @@ parser.add_option("--kinship",dest="kinship",
                   help="Kinship file format 1.0")
 parser.add_option("--pheno",dest="pheno",
                   help="Phenotype file format 1.0")
+parser.add_option("--geno",dest="geno",
+                  help="Genotype file format 1.0")
 parser.add_option("-q", "--quiet",
                   action="store_false", dest="verbose", default=True,
                   help="don't print status messages to stdout")
@@ -84,4 +86,20 @@ if options.pheno:
             Y1.append(ns) # <--- slow
     Y = np.array(Y1)
 
-print Y
+if options.geno:
+    G1 = []
+    print options.geno
+    with open(options.geno,'r') as tsvin:
+        assert(tsvin.readline().strip() == "# Genotype format version 1.0")
+        tsvin.readline()
+        tsvin.readline()
+        tsvin.readline()
+        tsvin.readline()
+        tsv = csv.reader(tsvin, delimiter='\t')
+        for row in tsv:
+            print(row)
+            # ns = np.genfromtxt(row[1:])
+            G1.append(ns) # <--- slow
+    G = np.array(G1)
+
+print G