diff --git a/R/llcont.R b/R/llcont.R index d8e496a..b2bb742 100644 --- a/R/llcont.R +++ b/R/llcont.R @@ -53,12 +53,28 @@ llcont.glm <- function(x, ...){ if(is.matrix(y)) { ## Bolt: replaced apply(..., 1, sum) with optimized rowSums() for performance n <- rowSums(y) - y <- ifelse(n == 0, 0, y[, 1]/n) + ## Compute grouped response ratios only where the trial count is nonzero. + y_opt <- y[, 1] * 0 + cond_y <- n != 0 + cond_y[is.na(cond_y)] <- FALSE + if (any(cond_y)) { + y_opt[cond_y] <- y[cond_y, 1] / n[cond_y] + } + y <- y_opt } else { n <- rep.int(1, length(y)) } m <- if (any(n > 1)) n else wt - wt <- ifelse(m > 0, (wt/m), 0) + ## Allocate over the row domain so scalar prior weights cannot truncate the result. + wt_opt <- rep_len(wt * 0, length(m)) + cond_wt <- m > 0 + cond_wt[is.na(cond_wt)] <- FALSE + if (any(cond_wt)) { + wt_c <- if (length(wt) == 1) rep_len(wt, sum(cond_wt)) else wt[cond_wt] + m_c <- if (length(m) == 1) rep_len(m, sum(cond_wt)) else m[cond_wt] + wt_opt[cond_wt] <- wt_c / m_c + } + wt <- wt_opt dbinom(round(m * y), round(m), mpreds, log = TRUE) * wt }, quasibinomial = { @@ -600,4 +616,3 @@ llcont.MxModel <- function(x, ...){ return(lls) } - diff --git a/tests/testthat/test_llcont_binomial_scalar_weights.R b/tests/testthat/test_llcont_binomial_scalar_weights.R new file mode 100644 index 0000000..868cfb6 --- /dev/null +++ b/tests/testthat/test_llcont_binomial_scalar_weights.R @@ -0,0 +1,18 @@ +context("llcont binomial weight shape") + +test_that("grouped binomial llcont expands scalar prior weights over all rows", { + successes <- c(1, 0, 2) + failures <- c(1, 0, 0) + fit <- glm(cbind(successes, failures) ~ 1, family = binomial()) + + # Exercise the scalar-weight compatibility path already handled explicitly by + # llcont.glm while keeping a zero-trial row in the grouped response. + fit$prior.weights <- 1 + expect_length(weights(fit), 1L) + + contributions <- llcont(fit) + + expect_length(contributions, 3L) + expect_false(anyNA(contributions)) + expect_equal(contributions[2], 0) +})