1 #' @title lmerExp: transform the unit of coefficients (internal function) 2 #' @description Transform the unit of coefficients to "Coeff", "OR" or "RR" 3 #' @param lmer.coef lmer coefficients. 4 #' @param family Family: "gaussian", "binomial", "poisson", "quasipoisson", etc..., Default: 'binomial' 5 #' @param dec Decimal point 6 #' @return The transforemed coefficients(95% CI), p-value 7 #' @details DETAILS 8 #' @examples 9 #' #EXAMPLE1 10 #' @rdname lmerExp 11 #' @importFrom stats pnorm 12 13 14 lmerExp = function(lmer.coef, family ="binomial", dec){ 15 1 if (class(lmer.coef)[1] == "numeric"){ 16 1 lmer.coef = t(data.frame(lmer.coef)) 17 } 18 1 pv = 2*(1-pnorm(abs(lmer.coef[,3]))) 19 1 if (family == "binomial"){ 20 1 OR = paste(round(exp(lmer.coef[,1]), dec), " (", round(exp(lmer.coef[,1] - 1.96*lmer.coef[,2]), dec), ",", round(exp(lmer.coef[,1] + 1.96*lmer.coef[,2]), dec),")", sep="") 21 1 return(cbind(OR, pv)) 22 1 } else if (family %in% c(NULL, "gaussian")){ 23 1 coeff = paste(round(lmer.coef[,1], dec), " (", round(lmer.coef[,1] - 1.96*lmer.coef[,2], dec), ",", round(lmer.coef[,1] + 1.96*lmer.coef[,2], dec), ")", sep="") 24 1 return(cbind(coeff, pv)) 25 0 } else if (family %in% c("poisson", "quasipoisson")){ 26 0 RR = paste(round(exp(lmer.coef[,1]), dec), " (", round(exp(lmer.coef[,1] - 1.96*lmer.coef[,2]), dec), ",", round(exp(lmer.coef[,1] + 1.96*lmer.coef[,2]), dec),")", sep="") 27 0 return(cbind(RR, pv)) 28 } 29 } 30 31 32 33 34 35 36 37 #' @title lmer.display: table for "lmerMod" or "glmerMod" object (lme4 package) 38 #' @description Make mixed effect model results from "lmerMod" or "glmerMod" object (lme4 package) 39 #' @param lmerMod.obj "lmerMod" or "glmerMod" object 40 #' @param dec Decimal, Default: 2 41 #' @param ci.ranef Show confidence interval of random effects?, Default: F 42 #' @return Table: fixed & random effect 43 #' @details DETAILS 44 #' @examples 45 #' library(geepack) 46 #' data(dietox) 47 #' dietox\$Cu <- as.factor(dietox\$Cu) 48 #' l1 <- lme4::lmer(Weight ~ Time + Cu + (1|Pig) + (1|Evit), data = dietox) 49 #' lmer.display(l1) 50 #' @rdname lmer.display 51 #' @export 52 #' @importFrom lme4 lmer glmer confint.merMod 53 54 lmer.display = function(lmerMod.obj, dec = 2, ci.ranef = F){ 55 1 sl = summary(lmerMod.obj) 56 1 fixef = sl\$coefficients[-1,] 57 58 1 mdata = lmerMod.obj@frame 59 1 forms = as.character(sl\$call[[2]]) 60 1 y = forms[2] 61 1 xfr = strsplit(forms[3]," \\+ ")[[1]] 62 1 xr = xfr[grepl("\\|", xfr)] 63 1 xf = xfr[!grepl("\\|", xfr)] 64 1 family.lmer = ifelse(is.null(sl\$family), "gaussian", sl\$family) 65 1 uni.res = "" 66 1 if (family.lmer == "gaussian"){ 67 1 unis = lapply(xf, function(x){summary(lmer(as.formula(paste(y, "~", x," + ", paste(xr, collapse =" + "), sep="")), data = mdata))\$coefficients}) 68 1 unis2 = Reduce(rbind, unis) 69 1 uni.res <- unis2[rownames(unis2) != "(Intercept)",] 70 } else{ 71 1 unis = lapply(xf, function(x){summary(glmer(as.formula(paste(y, "~", x," + ", paste(xr, collapse =" + "), sep="")), data = mdata, family = family.lmer))\$coefficients}) 72 1 unis2 = Reduce(rbind, unis) 73 1 uni.res <- unis2[rownames(unis2) != "(Intercept)",] 74 } 75 76 1 if(length(xf) == 1){ 77 1 fix.all = lmerExp(uni.res, family=family.lmer, dec = dec) 78 1 family.label = colnames(fix.all)[1] 79 1 colnames(fix.all) = c(paste(family.label, "(95%CI)",sep = ""), "P value") 80 1 rownames(fix.all) = rownames(sl\$coefficients)[-1] 81 } else{ 82 1 fix.all = cbind(lmerExp(uni.res, family=family.lmer, dec = dec), lmerExp(fixef, family=family.lmer, dec = dec)) 83 1 family.label = colnames(fix.all)[1] 84 1 colnames(fix.all) = c(paste("crude ", family.label, "(95%CI)",sep = ""), "crude P value", paste("adj. ", family.label, "(95%CI)",sep = ""), "adj. P value") 85 1 rownames(fix.all) = rownames(sl\$coefficients)[-1] 86 } 87 88 1 ranef = data.frame(sl\$varcor)[,c(1,4)] 89 1 ranef.out = cbind(as.character(round(ranef[,2], dec)), matrix(NA, nrow(ranef), ncol(fix.all) - 1 )) 90 1 rownames(ranef.out) = ranef[,1] 91 1 if (ci.ranef){ 92 0 ranef.ci = confint.merMod(lmerMod.obj, parm = 1:nrow(ranef), oldNames = F)^2 93 0 ranef.paste = paste(round(ranef\$vcov, dec)," (", round(ranef.ci[,1], dec), ",", round(ranef.ci[,2], dec),")", sep="") 94 0 ranef.out = cbind(ranef.paste, matrix(NA, nrow(ranef), 3)) 95 0 rownames(ranef.out) = ranef[,1] 96 } 97 98 1 ranef.na = rbind(rep(NA, ncol(fix.all)), ranef.out) 99 1 rownames(ranef.na)[1] = "Random effects" 100 1 colnames(ranef.na) = colnames(fix.all) 101 102 103 ## rownames 104 1 fix.all.list = lapply(xf, function(x){fix.all[grepl(x, rownames(fix.all)),]}) 105 1 varnum.mfac = which(lapply(fix.all.list, length) > ncol(fix.all)) 106 1 lapply(varnum.mfac, function(x){fix.all.list[[x]] <<- rbind(rep(NA, ncol(fix.all)), fix.all.list[[x]])}) 107 1 fix.all.unlist = Reduce(rbind, fix.all.list) 108 109 1 rn.list = lapply(xf, function(x){rownames(fix.all)[grepl(x, rownames(fix.all))]}) 110 1 varnum.2fac = which(lapply(xf, function(x){length(sapply(mdata, levels)[[x]])}) == 2) 111 0 lapply(varnum.2fac, function(x){rn.list[[x]] <<- paste(xf[x], ": ", levels(mdata[, xf[x]])[2], " vs ", levels(mdata[, xf[x]])[1], sep="")}) 112 1 lapply(varnum.mfac, function(x){rn.list[[x]] <<- c(paste(xf[x],": ref.=", levels(mdata[, xf[x]])[1], sep=""), gsub(xf[x]," ", rn.list[[x]]))}) 113 1 if (class(fix.all.unlist)[1] == "character"){ 114 1 fix.all.unlist = t(data.frame(fix.all.unlist)) 115 } 116 1 rownames(fix.all.unlist) = unlist(rn.list) 117 1 tb.df = as.data.frame(rbind(fix.all.unlist, ranef.na)) 118 119 120 ## metric 121 1 no.grp = sl\$ngrps 122 1 no.obs = nrow(mdata) 123 1 ll = round(sl\$logLik[[1]], dec) 124 1 aic = round(sl\$AICtab, dec)[1] 125 1 met = c(NA, no.grp, no.obs, ll, aic) 126 1 met.mat = cbind(met, matrix(NA, length(met), ncol(fix.all) - 1 )) 127 1 rownames(met.mat) = c("Metrics", paste("No. of groups (", rownames(met.mat)[2:(length(no.grp) + 1)], ")", sep=""), "No. of observations", "Log-likelihood", "AIC value") 128 1 met.mat[, 1] = as.character(met.mat[, 1]) 129 1 colnames(met.mat) = colnames(tb.df) 130 1 tb.lmerMod = rbind(tb.df, met.mat) 131 1 lapply(seq(2, ncol(fix.all), by =2), function(x){tb.lmerMod[, x] <<- as.numeric(as.vector(tb.lmerMod[, x]))}) 132 133 ## caption 134 1 cap.lmerMod = paste(sl\$methTitle, " : ", y," ~ ", forms[3], sep="") 135 1 return(list(table = tb.lmerMod, caption = cap.lmerMod)) 136 }

Read our documentation on viewing source code .