jinseob2kim / jstable
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 2
  if (class(lmer.coef)[1] == "numeric"){
16 2
    lmer.coef = t(data.frame(lmer.coef))
17
  }
18 2
  pv = 2*(1-pnorm(abs(lmer.coef[,3])))
19 2
  if (family == "binomial"){
20 2
    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 2
    return(cbind(OR, pv))
22 2
  } else if (family %in% c(NULL, "gaussian")){
23 2
    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 2
    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 2
  sl = summary(lmerMod.obj)
56 2
  fixef = sl$coefficients[-1,]
57
  
58 2
  mdata = lmerMod.obj@frame
59 2
  forms = as.character(sl$call[[2]])
60 2
  y = forms[2]
61 2
  xfr = strsplit(forms[3]," \\+ ")[[1]]
62 2
  xr = xfr[grepl("\\|", xfr)]
63 2
  xf = xfr[!grepl("\\|", xfr)]
64 2
  family.lmer = ifelse(is.null(sl$family), "gaussian", sl$family)
65 2
  uni.res = ""
66 2
  if (family.lmer == "gaussian"){
67 2
    unis = lapply(xf, function(x){summary(lmer(as.formula(paste(y, "~", x," + ", paste(xr, collapse =" + "), sep="")), data = mdata))$coefficients})
68 2
    unis2 = Reduce(rbind, unis)
69 2
    uni.res <- unis2[rownames(unis2) != "(Intercept)",]
70
  } else{
71 2
    unis = lapply(xf, function(x){summary(glmer(as.formula(paste(y, "~", x," + ", paste(xr, collapse =" + "), sep="")), data = mdata, family = family.lmer))$coefficients})
72 2
    unis2 = Reduce(rbind, unis)
73 2
    uni.res <- unis2[rownames(unis2) != "(Intercept)",]
74
  }
75
  
76 2
  if(length(xf) == 1){
77 2
    fix.all = lmerExp(uni.res, family=family.lmer, dec = dec)
78 2
    family.label = colnames(fix.all)[1]
79 2
    colnames(fix.all) = c(paste(family.label, "(95%CI)",sep = ""), "P value")
80 2
    rownames(fix.all) = rownames(sl$coefficients)[-1]
81
  } else{
82 2
    fix.all = cbind(lmerExp(uni.res, family=family.lmer, dec = dec), lmerExp(fixef, family=family.lmer, dec = dec))
83 2
    family.label = colnames(fix.all)[1]
84 2
    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 2
    rownames(fix.all) = rownames(sl$coefficients)[-1]
86
  }
87
  
88 2
  ranef = data.frame(sl$varcor)[,c(1,4)]
89 2
  ranef.out = cbind(as.character(round(ranef[,2], dec)), matrix(NA, nrow(ranef), ncol(fix.all) - 1 ))
90 2
  rownames(ranef.out) = ranef[,1]
91 2
  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 2
  ranef.na = rbind(rep(NA, ncol(fix.all)), ranef.out)
99 2
  rownames(ranef.na)[1] = "Random effects"
100 2
  colnames(ranef.na) = colnames(fix.all)
101
  
102
  
103
  ## rownames
104 2
  fix.all.list = lapply(xf, function(x){fix.all[grepl(x, rownames(fix.all)),]})      
105 2
  varnum.mfac = which(lapply(fix.all.list, length) > ncol(fix.all))
106 2
  lapply(varnum.mfac, function(x){fix.all.list[[x]] <<- rbind(rep(NA, ncol(fix.all)), fix.all.list[[x]])})
107 2
  fix.all.unlist = Reduce(rbind, fix.all.list)
108
  
109 2
  rn.list = lapply(xf, function(x){rownames(fix.all)[grepl(x, rownames(fix.all))]})
110 2
  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 2
  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 2
  if (class(fix.all.unlist)[1] == "character"){
114 2
    fix.all.unlist = t(data.frame(fix.all.unlist))
115
  }
116 2
  rownames(fix.all.unlist) = unlist(rn.list)
117 2
  tb.df = as.data.frame(rbind(fix.all.unlist, ranef.na))
118
  
119
  
120
  ## metric
121 2
  no.grp = sl$ngrps
122 2
  no.obs = nrow(mdata)
123 2
  ll = round(sl$logLik[[1]], dec)
124 2
  aic = round(sl$AICtab, dec)[1]
125 2
  met = c(NA, no.grp, no.obs, ll, aic)
126 2
  met.mat = cbind(met, matrix(NA, length(met), ncol(fix.all) - 1 ))
127 2
  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 2
  met.mat[, 1] = as.character(met.mat[, 1])
129 2
  colnames(met.mat) = colnames(tb.df)
130 2
  tb.lmerMod = rbind(tb.df, met.mat)
131 2
  lapply(seq(2, ncol(fix.all), by =2), function(x){tb.lmerMod[, x] <<- as.numeric(as.vector(tb.lmerMod[, x]))})
132
  
133
  ## caption
134 2
  cap.lmerMod = paste(sl$methTitle, " : ", y," ~ ", forms[3], sep="")
135 2
  return(list(table = tb.lmerMod, caption = cap.lmerMod))
136
}

Read our documentation on viewing source code .

Loading