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 .

Loading