Skip to content

Commit

Permalink
clean up tests etc
Browse files Browse the repository at this point in the history
  • Loading branch information
bbolker committed Oct 28, 2024
1 parent 56668d2 commit a355ff3
Show file tree
Hide file tree
Showing 5 changed files with 46 additions and 60 deletions.
5 changes: 5 additions & 0 deletions NEWS.md
Original file line number Diff line number Diff line change
@@ -1,3 +1,8 @@
# CHANGES IN broom.mixed VERSION 0.2.9.6

- CRAN maintenance release
- `stanreg` tidiers should work for models without random effects

# CHANGES IN broom.mixed VERSION 0.2.9.5

## NEW FEATURES
Expand Down
73 changes: 27 additions & 46 deletions R/VarCorr_tidiers.R
Original file line number Diff line number Diff line change
@@ -1,34 +1,35 @@
## https://github.com/bbolker/broom.mixed/issues/152, @phargarten2

#' Tidying VarCorr of Mixed Effect Models
#' @param object a VarCorr object
#' @param scales see \code{\link{tidy.lme4}}
#' @param x a VarCorr object
#' @param scales see \code{\link{tidy.merMod}}
#' @inheritParams tidy.merMod
#' @importFrom dplyr as_tibble
#' @export
tidy.VarCorr.lme <- function(
object,
x,
## effects = c("ran_pars", "fixed"), #effects are always "ran_pars"
scales = c("sdcor", "vcov"),
conf.int = FALSE,
conf.level = 0.95,
...) {

vcov <- var1 <- NULL ## NSE/R CMD check
grp <- vcov <- var1 <- var2 <- sdcor <- estimate <- s <- NULL ## NSE/R CMD check
check_dots(...)
scales <- match.arg(scales)
if(inherits(object, "VarCorr.lme")) {

if(inherits(x, "VarCorr.lme")) {
## Need to convert nlme object to a proper tibble
object <- convert_VarCorr.lme(object)
x <- convert_VarCorr.lme(x)
}

if (conf.int) {
cli::cli_alert_info("Can't compute confidence intervals for random effect parameters in a crossed model.")
message("Can't compute confidence intervals from VarCorr")
}

if (scales == "vcov") {
ests_random <- as_tibble(object) |>
rename(group = grp, term = var1, estimate_var = vcov) |>
ests_random <- as_tibble(x) |>
rename(group = grp, term = var1, estimate = vcov) |>
mutate(
term = case_when(
!is.na(term) & !is.na(var2) ~ paste0("cov_", term, "_", var2),
Expand All @@ -38,25 +39,11 @@ tidy.VarCorr.lme <- function(
) |>
mutate(effect = "ran_pars", .before = "group")

total_var <- (ests_random
|> filter(is.na(var2))
|> summarize(s = sum(estimate_var))
|> pull(s)
)

ests_random <- ests_random |>
mutate(prop_var = case_when(
!stringr::str_detect(term, "cov_") ~ estimate_var/total_var,
TRUE ~ NA_real_
)
) |>
select(-var2, -sdcor)

##Check
## print(as_tibble(object))

} else if(scales == "sdcor") {
ests_random <- as_tibble(object) |>
ests_random <- as_tibble(x) |>
rename(group = grp, term = var1, estimate = sdcor) |>
mutate(
term = case_when(
Expand All @@ -72,47 +59,41 @@ tidy.VarCorr.lme <- function(
total_var <- ests_random |> filter(is.na(var2)) |> summarize(s = sum(estimate^2)) |> pull(s)

ests_random <- ests_random |>
mutate(prop_var = case_when(
!stringr::str_detect(term, "cor_") ~ estimate^2/total_var,
TRUE ~ NA_real_
)
) |>
## mutate(prop_var = case_when(
## !stringr::str_detect(term, "cor_") ~ estimate^2/total_var,
## TRUE ~ NA_real_
## )
## ) |>
select(-var2, -vcov)

}

return(ests_random)
}

# lmm.lme4 <- lme4::lmer(Reaction ~ Days + (Days | Subject), sleepstudy)
#as.tibble(lme4::VarCorr(lmm.lme4))
# grp var1 var2 vcov sdcor
# <chr> <chr> <chr> <dbl> <dbl>
# 1 Subject (Intercept) NA 612. 24.7
# 2 Subject Days NA 35.1 5.92
# 3 Subject (Intercept) Days 9.60 0.0656
# 4 Residual NA NA 655. 25.6


#' This function converts VarCorr on a nlme object to a tibble
#' This function converts VarCorr on a nlme x to a tibble
#' i.e. VarCorr.merMod -> tibble of VarCorr.lme
#' @noRd
#' @param object A variance-correlation component of nlme::lme object.
#' @param x A variance-correlation component of nlme::lme object.
#' @return A useful (?) tibble
convert_VarCorr.lme <- function(object){
A <- as.matrix(object)
convert_VarCorr.lme <- function(x) {

sdcor <- NULL

A <- as.matrix(x)
row.residual <- stringr::str_which("Residual", rownames(A))

##Unsure how to generalize this
corr.A <- A[-row.residual, -(1:2), drop = TRUE]
corr.A <- A[-row.residual, "Corr", drop = TRUE]
t.corr <- tibble(var = names(corr.A), corr = as.vector(corr.A))
corr <- tibble(
var1 = t.corr[1, "var", drop = TRUE],
var2 = t.corr[2, "var", drop = TRUE],
sdcor = t.corr[2, "corr", drop = TRUE]
)

tib <- tibble(grp = NA_character_, var1 = rownames(A), var2 = NA, vcov= A[ ,"Variance"], sdcor = A[ , "StdDev"]) |>
tib <- tibble(grp = NA_character_, var1 = rownames(A), var2 = NA,
vcov= A[ ,"Variance"], sdcor = A[ , "StdDev"]) |>
bind_rows(corr) |>
mutate(
grp = case_when(
Expand Down
7 changes: 4 additions & 3 deletions R/nlme_tidiers.R
Original file line number Diff line number Diff line change
Expand Up @@ -143,10 +143,11 @@ tidy.lme <- function(x, effects = c("var_model", "ran_pars", "fixed"),
stop(sprintf("unrecognized ran_pars scale %s", sQuote(rscale)))
}
nonlin <- inherits(x, "nlme")

grplen <- attr(x$modelStruct$reStruct, "plen")
multilevel <- length(grplen) > 1
## FIXME: work on multilevel models
if (nonlin) {
warning("ran_pars not yet implemented for nonlinear models")
if (nonlin || multilevel) {
warning("ran_pars not yet implemented for nonlinear or multilevel models")
ret <- dplyr::tibble()
} else {
ret <- tidy_varcov(x, rscale = rscale, conf.int = conf.int, conf.level = conf.level)
Expand Down
18 changes: 9 additions & 9 deletions man/tidy.VarCorr.lme.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

3 changes: 1 addition & 2 deletions tests/testthat/test-VarCorr.R
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,6 @@ if (requireNamespace("lme4", quietly = TRUE) &&
requireNamespace("nlme", quietly=TRUE)) {
library(lme4)
library(nlme)
devtools::load_all()
data("sleepstudy", package="lme4")
lmm.nlme <- lme(Reaction ~ Days, random=~ Days|Subject, sleepstudy)
# > VarCorr(lmm.nlme)
Expand All @@ -22,6 +21,6 @@ if (requireNamespace("lme4", quietly = TRUE) &&
random = list(Subject = ~ Days, Group = ~ 1),
data = sleepstudy2
)
tidy(lmm2)
## tidy(lmm2)

}

0 comments on commit a355ff3

Please sign in to comment.