added error message for nanc=1

jhavsmith · Feb 2, 2018 · 2aded49 · 2aded49
1 parent 72539ca
commit 2aded49
Show file tree

Hide file tree

Showing 2 changed files with 22 additions and 9 deletions.
diff --git a/R/functions.R b/R/functions.R
@@ -238,11 +238,18 @@ prep.func = function(input.list,params) {
     }
     # Now we divide the data up into halves. The model runs from left to right starting
     # from the selected site, so we have to flip the left side around.
-    left.sample  = sample[,sel.site:1]
-    left.cont    = cont.sample[,sel.site:1]
-    left.pos     = abs(positions[sel.site:1]-positions[sel.site])+1
+    if (is.vector(cont.sample)) {
+        print("You need a minimum of 2 haplotypes in the reference panel")
+        stopifnot(!is.vector(cont.sample))
+        break
+    }
+    if (!is.vector(cont.sample)) { 
+        left.cont    = cont.sample[,sel.site:1]
+        right.cont   = cont.sample[,sel.site:ncol(cont.sample)] 
+    }
+    left.sample  = sample[,sel.site:1]  
     right.sample = sample[,sel.site:ncol(sample)]
-    right.cont   = cont.sample[,sel.site:ncol(cont.sample)]
+    left.pos     = abs(positions[sel.site:1]-positions[sel.site])+1  
     right.pos    = (positions[sel.site:ncol(sample)]-positions[sel.site])+1
     # The bedfile can tell us which (if any) positions were not actually sequenced in the vcf.
     # We need to know this because invariant sites may instead just be unobserved (or unsequenced) sites.

diff --git a/R/functions.R~ b/R/functions.R~
@@ -143,7 +143,6 @@ get.vcfdata.func = function(params) {
 
 # This is another data processing function before running the MCMC.
 prep.func = function(input.list,params) {
-    recover()
     # Here we define the variables we want from the lists 'input.list' and 'params'.
     psel          = params$pos
     nanc          = params$nanc
@@ -239,11 +238,18 @@ prep.func = function(input.list,params) {
     }
     # Now we divide the data up into halves. The model runs from left to right starting
     # from the selected site, so we have to flip the left side around.
-    left.sample  = sample[,sel.site:1]
-    left.cont    = cont.sample[,sel.site:1]
-    left.pos     = abs(positions[sel.site:1]-positions[sel.site])+1
+    if (is.vector(cont.sample)) {
+        print("You need a minimum of 2 haplotypes in the reference panel")
+        stopifnot(!is.vector(cont.sample))
+        break
+    }
+    if (!is.vector(cont.sample)) { 
+        left.cont    = cont.sample[,sel.site:1]
+        right.cont   = cont.sample[,sel.site:ncol(cont.sample)] 
+    }
+    left.sample  = sample[,sel.site:1]  
     right.sample = sample[,sel.site:ncol(sample)]
-    right.cont   = cont.sample[,sel.site:ncol(cont.sample)]
+    left.pos     = abs(positions[sel.site:1]-positions[sel.site])+1  
     right.pos    = (positions[sel.site:ncol(sample)]-positions[sel.site])+1
     # The bedfile can tell us which (if any) positions were not actually sequenced in the vcf.
     # We need to know this because invariant sites may instead just be unobserved (or unsequenced) sites.