diff --git a/R/pcm.R b/R/pcm.R
index c65e62b..97b6aa8 100644
--- a/R/pcm.R
+++ b/R/pcm.R
@@ -43,11 +43,13 @@ pcm <- function(df=NULL,items=NULL,grp=NULL,dif.items=NULL,type.dif=NULL,weights
   if (!("id"%in%colnames(df))) {
     stop('ERROR: no column named id provided')
   }
-  if (!is.null(grp) & any(apply(df[df[,grp]==0,items],2,max)<max(df[,items])) | any(apply(df[df[,grp]==1,items],2,max)<max(df[,items])) ) {
+  suppressWarnings(
+  if (!is.null(grp) & any(apply(df[df[,grp]==0,items],2,function(k) max(k,na.rm = T))<max(df[,items],na.rm = T)) | any(apply(df[df[,grp]==1,items],2,function(k) max(k,na.rm = T))<max(df[,items],na.rm=T)) ) {
     if (fit=="ucminf") {
       fit <- "optim"
     }
   }
+  )
   ##### Analysis
   restab.diftype <- NULL
   se.beta <- NULL
@@ -79,12 +81,12 @@ pcm <- function(df=NULL,items=NULL,grp=NULL,dif.items=NULL,type.dif=NULL,weights
     if (is.null(weights)) {
       colnames(df.long) <- c("id","item","resp")
       nbitems <- length(2:(length(colnames(df))))
-      maxmod <- max(df[,2:(length(colnames(df)))])
+      maxmod <- max(df[,2:(length(colnames(df)))],na.rm = T)
       df.long$item <- factor(df.long$item,levels=seq(1,length(colnames(df))-1),ordered = F)
     } else {
       colnames(df.long) <- c("id","weights","item","resp")
       nbitems <- length(2:(length(colnames(df))-1))
-      maxmod <- max(df[,2:(length(colnames(df))-1)])
+      maxmod <- max(df[,2:(length(colnames(df))-1)],na.rm = T)
       df.long$item <- factor(df.long$item,levels=seq(1,length(colnames(df))-2),ordered = F)
     }
     df.long$resp <- factor(df.long$resp,0:maxmod,ordered=T)
@@ -128,12 +130,12 @@ pcm <- function(df=NULL,items=NULL,grp=NULL,dif.items=NULL,type.dif=NULL,weights
       if (is.null(weights)) {
         colnames(df.long) <- c("id","grp","item","resp")
         nbitems <- length(2:(length(colnames(df))-1))
-        maxmod <- max(df[,2:(length(colnames(df))-1)])
+        maxmod <- max(df[,2:(length(colnames(df))-1)],na.rm = T)
         df.long$item <- factor(df.long$item,levels=seq(1,length(colnames(df))-2),ordered = F)
       } else {
         colnames(df.long) <- c("id","grp","weights","item","resp")
         nbitems <- length(2:(length(colnames(df))-2))
-        maxmod <- max(df[,2:(length(colnames(df))-2)])
+        maxmod <- max(df[,2:(length(colnames(df))-2)],na.rm = T)
         df.long$item <- factor(df.long$item,levels=seq(1,length(colnames(df))-3),ordered = F)
       }
       df.long$resp <- factor(df.long$resp,0:maxmod,ordered=T)
@@ -230,12 +232,12 @@ pcm <- function(df=NULL,items=NULL,grp=NULL,dif.items=NULL,type.dif=NULL,weights
       if (is.null(weights)) {
         colnames(df.long) <- c("id","grp","item","resp")
         nbitems <- length(2:(length(colnames(df))-1))
-        maxmod <- max(df[,2:(length(colnames(df))-1)])
+        maxmod <- max(df[,2:(length(colnames(df))-1)],na.rm = T)
         df.long$item <- factor(df.long$item,levels=seq(1,length(colnames(df))-2),ordered = F)
       } else {
         colnames(df.long) <- c("id","grp","weights","item","resp")
         nbitems <- length(2:(length(colnames(df))-2))
-        maxmod <- max(df[,2:(length(colnames(df))-2)])
+        maxmod <- max(df[,2:(length(colnames(df))-2)],na.rm = T)
         df.long$item <- factor(df.long$item,levels=seq(1,length(colnames(df))-3),ordered = F)
       }
       df.long$resp <- factor(df.long$resp,0:maxmod,ordered=T)
@@ -268,13 +270,16 @@ pcm <- function(df=NULL,items=NULL,grp=NULL,dif.items=NULL,type.dif=NULL,weights
     }
 
   }
-  if (method.theta=="eap") {
-    theta <- c(-1*ranef(mod,norm=F)+ ifelse(is.null(grpo),0, ifelse(grp==1,beta,0) ) )
-  } else if (method.theta=="wle") {
-    theta <- PP::PP_gpcm(as.matrix(df[,items]),t(restab),rep(1,length(items)))$resPP$resPP[,1]
-  } else if (method.theta=="mle") {
-    theta <- PP::PP_gpcm(as.matrix(df[,items]),t(restab),rep(1,length(items)),type="mle")$resPP$resPP[,1]
+  if (is.null(grpo)) {
+    if (method.theta=="eap") {
+      theta <- c(-1*ranef(mod,norm=F)+ ifelse(is.null(grpo),0, ifelse(grp==1,beta,0) ) )
+    } else if (method.theta=="wle") {
+      theta <- PP::PP_gpcm(as.matrix(df[,items]),t(restab),rep(1,length(items)))$resPP$resPP[,1]
+    } else if (method.theta=="mle") {
+      theta <- PP::PP_gpcm(as.matrix(df[,items]),t(restab),rep(1,length(items)),type="mle")$resPP$resPP[,1]
+    }
   }
+
   resid <- apply(matrix(1:nbitems,ncol=length(nbitems)),1, function(k) sapply(1:nrow(df), function(j) resi(theta[j],restab[k,],df[j,items[k]],beta=0)))
   colnames(resid) <- items_o