Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

use insight::get_variances, add overdispersion() and zerocount() funs #24

Merged
merged 17 commits into from
Mar 29, 2019
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
11 changes: 4 additions & 7 deletions DESCRIPTION
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,7 @@ Authors@R: c(
role = c("aut", "cre"),
comment = c(ORCID = "0000-0001-5375-9967"))
)
Maintainer: The package maintainer <yourself@somewhere.net>
Maintainer: Daniel Lüdecke <d.luedecke@uke.de>
Description: More about what it does (maybe more than one line)
Use four spaces when indenting paragraphs within the Description.
Imports:
Expand All @@ -23,18 +23,15 @@ Imports:
Suggests:
knitr,
lme4,
nlme,
Matrix,
rmarkdown,
testthat,
covr,
rstanarm,
rstantools,
loo,
brms,
circus
Remotes:
easystats/insight,
easystats/bayestestR,
easystats/circus
brms
License: GPL-3
Encoding: UTF-8
LazyData: true
Expand Down
27 changes: 23 additions & 4 deletions NAMESPACE
Original file line number Diff line number Diff line change
@@ -1,19 +1,31 @@
# Generated by roxygen2: do not edit by hand

S3method(check_overdispersion,glm)
S3method(check_overdispersion,glmmTMB)
S3method(check_overdispersion,merMod)
S3method(check_overdispersion,negbin)
S3method(check_singularity,MixMod)
S3method(check_singularity,default)
S3method(check_singularity,glmmTMB)
S3method(check_singularity,lme)
S3method(check_singularity,merMod)
S3method(icc,brmsfit)
S3method(icc,default)
S3method(model_performance,brmsfit)
S3method(model_performance,glm)
S3method(model_performance,lm)
S3method(model_performance,merMod)
S3method(model_performance,stanreg)
S3method(r2,MixMod)
S3method(r2,brmsfit)
S3method(r2,clm)
S3method(r2,clm2)
S3method(r2,glm)
S3method(r2,glmmTMB)
S3method(r2,lm)
S3method(r2,lme)
S3method(r2,merMod)
S3method(r2,mixed)
S3method(r2,mlogit)
S3method(r2,multinom)
S3method(r2,plm)
Expand All @@ -39,6 +51,10 @@ S3method(r2_nagelkerke,glm)
S3method(r2_nagelkerke,multinom)
S3method(r2_nagelkerke,polr)
S3method(r2_nagelkerke,vglm)
export(check_convergence)
export(check_overdispersion)
export(check_singularity)
export(check_zeroinflation)
export(cronbachs_alpha)
export(error_rate)
export(icc)
Expand All @@ -49,7 +65,6 @@ export(item_split_half)
export(looic)
export(model_performance)
export(mse)
export(null_model)
export(r2)
export(r2_bayes)
export(r2_coxnell)
Expand All @@ -61,18 +76,22 @@ export(r2_tjur)
export(rmse)
export(rse)
importFrom(bayestestR,hdi)
importFrom(bayestestR,map_estimate)
importFrom(insight,find_response)
importFrom(insight,get_data)
importFrom(insight,get_response)
importFrom(insight,get_variance)
importFrom(insight,model_info)
importFrom(insight,n_obs)
importFrom(stats,AIC)
importFrom(stats,BIC)
importFrom(stats,binomial)
importFrom(stats,coef)
importFrom(stats,cor)
importFrom(stats,df.residual)
importFrom(stats,dpois)
importFrom(stats,family)
importFrom(stats,fitted)
importFrom(stats,formula)
importFrom(stats,lm)
importFrom(stats,logLik)
importFrom(stats,mad)
importFrom(stats,median)
Expand All @@ -81,7 +100,7 @@ importFrom(stats,nobs)
importFrom(stats,pchisq)
importFrom(stats,predict)
importFrom(stats,predict.glm)
importFrom(stats,reformulate)
importFrom(stats,quantile)
importFrom(stats,residuals)
importFrom(stats,sd)
importFrom(stats,setNames)
Expand Down
66 changes: 0 additions & 66 deletions R/WIP_r2_linmix.R

This file was deleted.

57 changes: 57 additions & 0 deletions R/check_convergence.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,57 @@
#' @title Convergence test for mixed effects models
#' @name check_convergence
#'
#' @description \code{check_convergence()} provides an alternative convergence
#' test for \code{\link[lme4]{merMod}}-objects.
#'
#' @param x A \code{merMod}-object.
#' @param tolerance Indicates up to which value the convergence result is
#' accepted. The smaller \code{tolerance} is, the stricter the test
#' will be.
#'
#' @return \code{TRUE} if convergence is fine and \code{FALSE} if convergence
#' is suspicious. Additionally, the convergence value is returned as return
#' value's name.
#'
#' @details \code{check_convergence()} provides an alternative convergence test for
#' \code{\link[lme4]{merMod}}-objects, as discussed
#' \href{https://github.com/lme4/lme4/issues/120}{here}
#' and suggested by Ben Bolker in
#' \href{https://github.com/lme4/lme4/issues/120#issuecomment-39920269}{this comment}.
#'
#' @examples
#' library(lme4)
#' data(cbpp)
#' set.seed(1)
#' cbpp$x <- rnorm(nrow(cbpp))
#' cbpp$x2 <- runif(nrow(cbpp))
#'
#' model <- glmer(
#' cbind(incidence, size - incidence) ~ period + x + x2 + (1 + x | herd),
#' data = cbpp,
#' family = binomial()
#' )
#'
#' check_convergence(model)
#'
#' @export
check_convergence <- function(x, tolerance = 0.001) {
# check for package availability
if (!requireNamespace("Matrix", quietly = TRUE))
stop("Package `Matrix` needed for this function to work. Please install it.", call. = FALSE)

# is 'x' an lmer object?
if (!inherits(x, "merMod"))
warning("`x` must be a `merMod` object.", call. = F)


relgrad <- with(x@optinfo$derivs, Matrix::solve(Hessian, gradient))

# copy logical value, TRUE if convergence is OK
retval <- max(abs(relgrad)) < tolerance
# copy convergence value
names(retval) <- max(abs(relgrad))

# return result
retval
}
97 changes: 97 additions & 0 deletions R/check_overdispersion.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,97 @@
#' @title Check overdispersion of GL(M)M's
#' @name check_overdispersion
#'
#' @description \code{check_overdispersion()} checks generalized linear (mixed) models
#' for overdispersion.
#'
#' @param x Fitted model of class \code{merMod}, \code{glmmTMB}, \code{glm},
#' or \code{glm.nb} (package \pkg{MASS}).
#' @param ... Currently not used.
#'
#' @return Results from the on the overdispersion test, like chi-squared statistics,
#' p-value or dispersion ratio.
#'
#' @details A p-value < .05 indicates overdispersion. For \code{merMod}- and
#' \code{glmmTMB}-objects, \code{check_overdispersion()} is based on the code in
#' the \href{http://bbolker.github.io/mixedmodels-misc/glmmFAQ.html}{GLMM FAQ},
#' section \emph{How can I deal with overdispersion in GLMMs?}. Note that
#' this function only returns an \emph{approximate} estimate of an
#' overdispersion parameter, and is probably inaccurate for zero-inflated
#' mixed models (fitted with \code{glmmTMB}). The same code is also used to
#' check overdispersion for negative binomial models.
#' \cr \cr
#' For Poisson-models, the overdispersion test is based on the code
#' from \cite{Gelman and Hill (2007), page 115}.
#'
#' @references Bolker B et al. (2017): \href{http://bbolker.github.io/mixedmodels-misc/glmmFAQ.html}{GLMM FAQ.}
#' \cr \cr
#' Gelman A, Hill J (2007) Data Analysis Using Regression and Multilevel/Hierarchical Models. Cambridge, New York: Cambridge University Press
#'
#' @examples
#'
#' @export
check_overdispersion <- function(x, ...) {
UseMethod("check_overdispersion")
}


#' @importFrom insight get_response
#' @importFrom stats fitted nobs coef pchisq
#' @export
check_overdispersion.glm <- function(x, ...) {
# check if we have poisson
if (!stats::family(x)$family %in% c("poisson", "quasipoisson"))
stop("Model must be from Poisson-family.", call. = F)

yhat <- stats::fitted(x)

n <- stats::nobs(x)
k <- length(stats::coef(x))

zi <- (insight::get_response(x) - yhat) / sqrt(yhat)
chisq <- sum(zi^2)
ratio <- chisq / (n - k)
p.value <- stats::pchisq(chisq, df = n - k, lower.tail = FALSE)

list(
chisq_statistic = chisq,
dispersion_ratio = ratio,
residual_df = n - k,
p_value = p.value
)
}


#' @export
check_overdispersion.negbin <- function(x, ...) {
check_overdispersion.lme4(x)
}


#' @export
check_overdispersion.merMod <- function(x, ...) {
check_overdispersion.lme4(x)
}


#' @export
check_overdispersion.glmmTMB <- function(x, ...) {
check_overdispersion.lme4(x)
}


#' @importFrom stats df.residual residuals pchisq
check_overdispersion.lme4 <- function(x) {
rdf <- stats::df.residual(x)
rp <- stats::residuals(x, type = "pearson")
Pearson.chisq <- sum(rp^2)
prat <- Pearson.chisq / rdf
pval <- stats::pchisq(Pearson.chisq, df = rdf, lower.tail = FALSE)

list(
chisq = Pearson.chisq,
ratio = prat,
rdf = rdf,
p = pval
)
}