[adegenet-commits] r1137 - pkg/R

[noreply at r-forge.r-project.org]

Author: jombart
Date: 2013-06-04 12:14:25 +0200 (Tue, 04 Jun 2013)
New Revision: 1137

Modified:
   pkg/R/dapc.R
Log:
added warning to xvalDapc

Modified: pkg/R/dapc.R
===================================================================
--- pkg/R/dapc.R	2013-05-28 19:58:19 UTC (rev 1136)
+++ pkg/R/dapc.R	2013-06-04 10:14:25 UTC (rev 1137)
@@ -88,8 +88,8 @@
         n.da <- as.integer(readLines(n = 1))
     }
 
-    n.da <- min(n.da, length(levels(grp))-1, n.pca) # can't be more than K-1 disc. func., or more than n.pca
-    n.da <- round(n.da)
+    ##n.da <- min(n.da, length(levels(grp))-1, n.pca) # can't be more than K-1 disc. func., or more than n.pca
+    n.da <- round(min(n.da, lda.dim)) # can't be more than K-1 disc. func., or more than n.pca
     predX <- predict(ldaX, dimen=n.da)
 
 
@@ -1020,9 +1020,11 @@
 
     ## FUNCTION GETTING THE % OF ACCURATE PREDICTION FOR ONE NUMBER OF PCA PCs ##
     ## n.pca is a number of retained PCA PCs
+    VOID.GRP <- FALSE # will be TRUE if empty group happened
     get.prop.pred <- function(n.pca){
         f1 <- function(){
             toKeep <- sample(1:N, N.training)
+            if(!(all(table(grp[toKeep])>0) & all(table(grp[-toKeep])>0))) VOID.GRP <<- TRUE
             temp.pca <- pcaX
             temp.pca$li <- temp.pca$li[toKeep,,drop=FALSE]
             temp.dapc <- suppressWarnings(dapc(x[toKeep,,drop=FALSE], grp[toKeep], n.pca=n.pca, n.da=n.da, dudi=temp.pca))
@@ -1031,7 +1033,7 @@
                 out <- mean(temp.pred$assign==grp[-toKeep])
             }
             if(result=="groupMean"){
-                out <- mean(tapply(temp.pred$assign==grp[-toKeep], grp[-toKeep], mean))
+                out <- mean(tapply(temp.pred$assign==grp[-toKeep], grp[-toKeep], mean), na.rm=TRUE)
             }
             return(out)
         }
@@ -1041,6 +1043,7 @@
 
     ## GET %SUCCESSFUL OF ACCURATE PREDICTION FOR ALL VALUES ##
     res.all <- unlist(lapply(n.pca, get.prop.pred))
+    if(VOID.GRP) warning("At least one group was absent from the training / validating sets.\nTry using smaller training sets.")
     res <- data.frame(n.pca=rep(n.pca, each=n.rep), success=res.all)
     return(res)
 } # end xvalDapc.data.frame