[Genabel-commits] r1220 - in pkg/GenABEL: . R

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

Author: yurii
Date: 2013-05-16 19:19:12 +0200 (Thu, 16 May 2013)
New Revision: 1220

Modified:
   pkg/GenABEL/CHANGES.LOG
   pkg/GenABEL/R/polygenic.R
Log:
negative eigenvalues set to 'very small' (1e-16) instead of using abs

Modified: pkg/GenABEL/CHANGES.LOG
===================================================================
--- pkg/GenABEL/CHANGES.LOG	2013-05-16 12:49:21 UTC (rev 1219)
+++ pkg/GenABEL/CHANGES.LOG	2013-05-16 17:19:12 UTC (rev 1220)
@@ -1,6 +1,8 @@
 ***  v. 1.7-6
 
 (2013.05.16)
+Negative eigenvalues set to 'very small' (1e-16) instead using abs 
+in 'polygenic': this is more justified 
 Fixed too long lines in Rd files.
 
 (2013.05.15)

Modified: pkg/GenABEL/R/polygenic.R
===================================================================
--- pkg/GenABEL/R/polygenic.R	2013-05-16 12:49:21 UTC (rev 1219)
+++ pkg/GenABEL/R/polygenic.R	2013-05-16 17:19:12 UTC (rev 1220)
@@ -295,9 +295,10 @@
 			eigres <- eigenOfRel
 		else
 			eigres <- eigen(relmat,symmetric=TRUE)
-		if (any(eigres$values<0)) {
-			#eigres$values <- abs(eigres$values)
-			warning("some eigenvalues <=0, taking ABS for det; try option llfun='polylik'")
+		if (any(eigres$values<1e-16)) {
+			eigres$values[eigres$values<1e-16] <- 1e-16
+            msg <- paste("some eigenvalues close/less than 1e-16, setting them to 1e-16\nyou can also try option llfun='polylik' instead")
+			warning(msg)
 		}
 	} else stop("cannot be here...")
 	if (!quiet) {
@@ -568,7 +569,7 @@
 #	print(proc.time()-time0)
 # end old variant
 # new variant
-	a <- sum(log(abs(es))) # logarithm of determinant of sigma as sum of logarithm of eigenvalues
+	a <- sum(log(es)) # logarithm of determinant of sigma as sum of logarithm of eigenvalues
 # end new variant
 	b <- (t(resY) %*% InvSigma_x_residualY)
 # this is -2*lnLikelihood