substrRight <- function(x, n=5){
substr(x, nchar(x)-n+1, nchar(x)) ### Thank you Andrie (http://stackoverflow.com/questions/7963898/extracting-the-last-n-characters-from-a-string-in-r)
}
substrRightStar <- function(x, n=6){
substr(x, nchar(x)-n+1, nchar(x)-1) ### Thank you Andrie (http://stackoverflow.com/questions/7963898/extracting-the-last-n-characters-from-a-string-in-r)
}
#' @export
print.GEEmediate <- function(x, digits = max(options()$digits - 4, 3),...)
{
cat("Call:\n")
print(x$call)
cat("\n")
if(x$pres=="tog")
{
cat("Coefficients:\n")
coeffs <- x$GEEfit$coefficients
cat("\n")
coeffs.names <- names(coeffs)
covmat <- x$GEEfit$robust.variance
sd.err <- sqrt(diag(covmat))
zvalue <- coeffs/sd.err
pvalue <- 2*pnorm(abs(zvalue),lower.tail = F)
coef.table <- cbind(coeffs, sd.err, zvalue, pvalue)
dimnames(coef.table) <- list(coeffs.names,
c("Estimate", "Std. Error", "z value", "Pr(>|z|)"))
printCoefmat(coef.table, digits = digits)
} else if (x$pres=="sep") {
cat("------------------------------------------------------------------------------")
cat("\n------------------------------------------------------------------------------")
cat("\nMarginal Model (Model without the Mediator):\n")
coeffs <- x$GEEfit$coefficients
coeffs.names <- names(coeffs)
#stars <- which(sapply(coeffs.names, substrRight)==".star" | sapply(coeffs.names, substrRightStar)==".star")
#no.stars <- which(sapply(coeffs.names, substrRight)!=".star" & sapply(coeffs.names, substrRightStar)!=".star")
stars <- which(sapply(coeffs.names, function(x) grepl(".star", x)))
no.stars <- which(!sapply(coeffs.names, function(x) grepl(".star", x)))
covmat <- x$GEEfit$robust.variance
sd.err <- sqrt(diag(covmat))
zvalue <- coeffs/sd.err
pvalue <- 2*pnorm(abs(zvalue),lower.tail = F)
coef.table <- cbind(coeffs, sd.err, zvalue, pvalue)
dimnames(coef.table) <- list(names(coeffs),
c("Estimate", "Std. Error", "z value", "Pr(>|z|)"))
cond.table <- coef.table[no.stars,]
rownames(cond.table)[1] <- "(Intercept)"
marg.table <- coef.table[stars,]
rownames(marg.table)[1] <- "(Intercept)"
rownames(marg.table)[-1] <- sapply(rownames(marg.table)[-1], function(x) gsub(pattern = ".star", "", x))
#cat(row.names(marg.table), "\n")
stats::printCoefmat(marg.table, digits = digits)
cat("\n------------------------------------------------------------------------------")
cat("\n------------------------------------------------------------------------------")
cat("\nConditional Model (Model with the Mediator):\n \n")
stats::printCoefmat(cond.table, digits = digits)
}
cat("\n---------------------------")
cat("\nNatural Indirect Effect: ", format(x$nie, digits = digits),
"\np=", format(x$nie.pval, digits = 2),
" for ", paste0(x$alter), " test for mediation \n", sep = "")
cat("Confidence Interval = [", format(x$nie.ci[1],digits = digits),",",format(x$nie.ci[2],digits = digits),"]", sep = "")
cat("\n---------------------------")
cat("\nNatural Direct Effect:", format(x$nde, digits = digits))
cat("\n---------------------------")
if(x$pm >=0 & x$pm < 1)
{
cat("\nMediation Proportion:", format(100*x$pm,digits = 3),"%",
"\np=", format.pval(x$pm.pval, digits = 2),
" for one-sided test for mediation \n", sep = "")
cat("Confidence Interval = [", format(100*x$pm.ci[1],digits = 3),"%",",",format(100*x$pm.ci[2],digits = 3),"%","]", sep = "")
} else
{
cat("\nMediation Proportion:", format(100*x$pm,digits = 3),"%")
}
cat("\n---------------------------")
}
#' @export
predict.GEEmediate <- function (object, newdata = NULL, model.pred = c("cond", "marg"),type = c("link", "response", "terms"),
se.fit = FALSE, dispersion = NULL, terms = NULL, na.action = na.pass, ...)
{
model.pred <- match.arg(model.pred)
type <- match.arg(type)
if (se.fit==T) {warning("se.fit=T is currently not supported for GEEmediate")}
if (type=="terms") {
warning("type='terms' is currently not supported for GEEmediate, using type='response' intead")
type <- "response"
}
gee.object <- object$GEEfit
if(!missing(newdata))
{
newdf <- newdata[gee.object$xnames[substrRight(gee.object$xnames)!=".star" & substrRight(gee.object$xnames)!="INT"]]
dupl.df.new <- DupliData(df = newdf, mediator = object$call[[4]], surv = F)
dupl.df.new <- dupl.df.new[gee.object$xnames]
eta <- gee.object$linear.predictors
if (type=="response") {
pred <- make.link(tolower(gee.object$model$link))$linkinv(eta)
} else {
pred <- eta
}
if (model.pred=="cond") {out <- pred[1:(gee.object$nobs/2)]}
if (model.pred=="marg") {out <- pred[(gee.object$nobs/2+1):gee.object$nobs]}
} else {
pred <- predict.glm(object = gee.object, type = type, se.fit = F, na.action = na.action, ...)
if (model.pred=="cond") {out <- pred[1:(gee.object$nobs/2)]}
if (model.pred=="marg") {out <- pred[(gee.object$nobs/2+1):gee.object$nobs]}}
out
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.