Fix bug in method.CC_nloglik when computing coefficients

ecpolley · Aug 6, 2014 · 76bfdc6 · 76bfdc6
1 parent fc143ef
commit 76bfdc6
Showing 1 changed file with 9 additions and 7 deletions.
diff --git a/R/method.R b/R/method.R
@@ -199,33 +199,35 @@ method.CC_nloglik <- function() {
     plogis(trimLogit(predY, trim = control$trimLogit) %*% matrix(coef))
   }
   computeCoef = function(Z, Y, libraryNames, obsWeights, control, verbose, ...) {
-    cvRisk <- apply(Z, 2, function(x) { -sum(2 * obsWeights * ifelse(Y, log(x), log(1-x))) } )
+    logitZ = trimLogit(Z, control$trimLogit)
+    cvRisk <- apply(logitZ, 2, function(x) -sum(2 * obsWeights * 
+                                       ifelse(Y, plogis(x, log.p=TRUE),
+                                                 plogis(x, log.p=TRUE, lower.tail=FALSE))))
     names(cvRisk) <- libraryNames
     obj_and_grad <- function(y,x, w=NULL) {
         y <- y
         x <- x
       function(beta) {
-        xB <- x %*% cbind(beta) 
+        xB <- x %*% cbind(beta)
         loglik <- y * plogis(xB, log.p=TRUE) + (1-y) * plogis(xB, log.p=TRUE, lower.tail=FALSE)
         if (!is.null(w)) loglik <- loglik * w
         obj <- -2 * sum(loglik)
-
         p <- plogis(xB)
         grad <- if (is.null(w)) 2 * crossprod(x, cbind(p - y))
         else 2 * crossprod(x, w*cbind(p - y))
         list(objective=obj, gradient=grad)
       }
-    }    
+    }
 
     r <- nloptr(x0=rep(1/ncol(Z), ncol(Z)),
-            eval_f=obj_and_grad(Y, Z),
+            eval_f=obj_and_grad(Y, logitZ),
             lb=rep(0, ncol(Z)),
             ub=rep(1, ncol(Z)),
             eval_g_eq = function(beta) (sum(beta)-1),
             eval_jac_g_eq = function(beta) rep(1, length(beta)),
             opts=list("algorithm"="NLOPT_LD_SLSQP","ftol_rel"=1.0e-8))
-    if (r$status < 1 || r$status > 4) { 
-      warning(r$messate)
+    if (r$status < 1 || r$status > 4) {
+      warning(r$message)
     }
     coef <- r$solution
     if (any(is.na(coef))) {