diff --git a/NEWS.md b/NEWS.md index 3df72616b..5b90612b5 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,8 +1,6 @@ # Items for next release -* New vignette on AB testing. (#409) - -* `stan_jm()` gains an offset term for the longitudinal submodel. (#415, @pamelanluna) +### Bug fixes * Fix bug where `loo()` with `k_threshold` argument specified would error if the model formula was a string instead of a formula object. (#454) @@ -16,6 +14,13 @@ models fit with `stan_polr()`. (#450) * Fix bug where `stan_glmer()` would error if `prior_aux=NULL`. (#482) +### New features + +* New vignette on AB testing. (#409) + +* `stan_jm()` gains an offset term for the longitudinal submodel. (#415, @pamelanluna) + +* Effective number of parameters are computed for K-fold CV not just LOO CV. (#462) # rstanarm 2.21.1 diff --git a/R/loo-kfold.R b/R/loo-kfold.R index e1d7a1ca8..019df7646 100644 --- a/R/loo-kfold.R +++ b/R/loo-kfold.R @@ -263,15 +263,22 @@ kfold.stanreg <- elpds <- rep(NA, length(elpds_unord)) elpds[obs_order] <- elpds_unord - pointwise <- cbind(elpd_kfold = elpds, p_kfold = NA, kfoldic = -2 * elpds) + # for computing effective number of parameters + ll_full <- log_lik(x) + lpds <- apply(ll_full, 2, log_mean_exp) + ps <- lpds - elpds + + pointwise <- cbind(elpd_kfold = elpds, p_kfold = ps, kfoldic = -2 * elpds) est <- colSums(pointwise) se_est <- sqrt(N * apply(pointwise, 2, var)) out <- list( estimates = cbind(Estimate = est, SE = se_est), pointwise = pointwise, - elpd_kfold = est[1], - se_elpd_kfold = se_est[1] + elpd_kfold = est[["elpd_kfold"]], + se_elpd_kfold = se_est[["elpd_kfold"]], + p_kfold = est[["p_kfold"]], + se_p_kfold = se_est[["p_kfold"]] ) rownames(out$estimates) <- colnames(pointwise) diff --git a/tests/testthat/test_loo.R b/tests/testthat/test_loo.R index a080f3f85..9044cfc28 100644 --- a/tests/testthat/test_loo.R +++ b/tests/testthat/test_loo.R @@ -166,12 +166,14 @@ test_that("kfold works on some examples", { SW(kf <- kfold(fit_gaus, 4)) SW(kf2 <- kfold(example_model, 2)) - expect_named(kf, c("estimates", "pointwise", "elpd_kfold", "se_elpd_kfold")) - expect_named(kf2, c("estimates", "pointwise", "elpd_kfold", "se_elpd_kfold")) + expect_named(kf, c("estimates", "pointwise", "elpd_kfold", "se_elpd_kfold", "p_kfold", "se_p_kfold")) + expect_named(kf2, c("estimates", "pointwise", "elpd_kfold", "se_elpd_kfold", "p_kfold", "se_p_kfold")) expect_named(attributes(kf), c("names", "class", "K", "dims", "model_name", "discrete", "yhash", "formula")) expect_named(attributes(kf2), c("names", "class", "K", "dims", "model_name", "discrete", "yhash", "formula")) expect_s3_class(kf, c("kfold", "loo")) expect_s3_class(kf2, c("kfold", "loo")) + expect_false(is.na(kf$p_kfold)) + expect_false(is.na(kf2$p_kfold)) SW(kf <- kfold(fit_gaus, K = 2, save_fits = TRUE))