check for NAs in alt, y, xFull and xNull

R/runvb.R

@@ -1,8 +1,5 @@
 packageName="eagle"
 
-xToListList = function(x) {
-    lapply( x, function(g) lapply( as.list(1:nrow(g)), function(h) g[h,] ))
-}
 
 # alt: list (over exonic SNPs) of alternative read counts
 # n: list (over exonics SNPs) of total read counts
@@ -93,6 +90,13 @@ eagle.helper = function(alt,n,xFull,xNull,s){
         s$normalised.depth=scale(log10(unlist(lapply(n,mean))))
 
   stopifnot( !is.nan(s$normalised.depth))
+
+  if (any(sapply(alt,function(g) any(is.na(g))))) stop("NAs not allowed in alt, please remove these elements.")
+  if (any(sapply(n,function(g) any(is.na(g))))) stop("NAs not allowed in n, please remove these elements.")
+  if (any(sapply(xFull,function(g) any(is.na(g))))) stop("NAs not allowed in xFull, please remove these elements.")
+  if (any(sapply(xNull,function(g) any(is.na(g))))) stop("NAs not allowed in xNull, please remove these elements.")
+
+  xToListList = function(x) lapply( x, function(g) lapply( as.list(1:nrow(g)), function(h) g[h,] ))
 
   xFullList=xToListList(xFull)
   xNullList=xToListList(xNull)

man/eagle.helper.Rd

@@ -16,7 +16,7 @@ eagle.helper(alt,n,xFull,xNull,s)
   \item{n}{ list (over exonics SNPs) of total read counts}
   \item{xFull}{list of design matrices for the alternative hypothesis (e.g. including environment)}
   \item{xNull}{list of design matrices for the null hypothesis} 
-  \item{settings}{list of settings. Create using \code{\link{eagle.settings}}() and then customize. }
+  \item{s}{list of settings. Create using \code{\link{eagle.settings}}() and then customize. }
 }
 \details{
 %%  ~~ If necessary, more details than the description above ~~

src/vbglmm.cpp

@@ -61,7 +61,7 @@ double fit_gamma(double s){
     shape=shape- (log(shape)-::Rf_digamma(shape)-s)/(1.0/shape - ::Rf_trigamma(shape)); 
     counter++;
     if (counter > 1000){
-      cerr << "Warning: fit_gamma(" << s << ") not converging, returning " << init_shape << endl; 
+      Rcerr << "Warning: fit_gamma(" << s << ") not converging, returning " << init_shape << endl; 
       return init_shape;
     }
   }
@@ -538,7 +538,7 @@ class Vbglmm {
 	if (debug && ((lb+1.0e-3)<old_lb)) Rcout << "Warning: lb got worse after optimizing reprep, old:" << old_lb << " new:" << lb << endl; 
       }
 
-      if (!isfinite(random_effect_precision)) { cerr << "err " << expected_err << " n " << total_terms << endl; ::Rf_error("Random effect precision is nan"); }
+      if (!isfinite(random_effect_precision)) { Rcerr << "err " << expected_err << " n " << total_terms << endl; ::Rf_error("Random effect precision is nan"); }
 
       if (debug && ((lb+1.0e-3)<old_lb)) Rcout << "Warning: lb got worse after optimizing random effect var, old:" << old_lb << " new:" << lb << endl; 
     } else {
@@ -750,7 +750,7 @@ class Vbglmm {
 	LDLT<MatrixXd> llt(xxX);
 	double det=MatrixXd(llt.matrixL()).diagonal().prod(); 
 	det = det * det; 
-	if (det < 1e-8) cerr << "Warning: X is ill conditioned for locus " << locus_index << endl; 
+	if (det < 1e-8) Rcerr << "Warning: X is ill conditioned for locus " << locus_index << endl; 
 	choleskys.push_back(llt); 
       } else {
 	LDLT<MatrixXd> lltnull;