From 128cca11a60961252c7dbc6a70bea3f7ca75585e Mon Sep 17 00:00:00 2001 From: Phillip Alday Date: Tue, 18 Mar 2025 16:03:40 -0500 Subject: [PATCH 1/4] LikelihoodRatioTest extends StatsAPI.HypothesisTest --- src/MixedModels.jl | 3 ++- src/likelihoodratiotest.jl | 21 ++++++++++++++++++++- test/likelihoodratiotest.jl | 10 +++++++--- 3 files changed, 29 insertions(+), 5 deletions(-) diff --git a/src/MixedModels.jl b/src/MixedModels.jl index d47fb92a5..4a49c657f 100644 --- a/src/MixedModels.jl +++ b/src/MixedModels.jl @@ -31,7 +31,7 @@ using StaticArrays: StaticArrays, SVector using Statistics: Statistics, mean, quantile, std using StatsAPI: StatsAPI, aic, aicc, bic, coef, coefnames, coeftable, confint, deviance using StatsAPI: dof, dof_residual, fit, fit!, fitted, isfitted, islinear, leverage -using StatsAPI: loglikelihood, meanresponse, modelmatrix, nobs, predict, r2, residuals +using StatsAPI: loglikelihood, meanresponse, modelmatrix, nobs, pvalue, predict, r2, residuals using StatsAPI: response, responsename, stderror, vcov, weights using StatsBase: StatsBase, CoefTable, model_response, summarystats using StatsFuns: log2π, normccdf @@ -129,6 +129,7 @@ export @formula, profileσ, profilesigma, profilevc, + pvalue, pwrss, ranef, raneftables, diff --git a/src/likelihoodratiotest.jl b/src/likelihoodratiotest.jl index 4e2bf7287..5781e4546 100644 --- a/src/likelihoodratiotest.jl +++ b/src/likelihoodratiotest.jl @@ -15,7 +15,7 @@ Results of MixedModels.likelihoodratiotest * `pvalues` """ -struct LikelihoodRatioTest +struct LikelihoodRatioTest <: StatsAPI.HypothesisTest formulas::AbstractVector{String} models::NamedTuple{(:dof, :deviance)} tests::NamedTuple{(:dofdiff, :deviancediff, :pvalues)} @@ -26,6 +26,25 @@ function Base.propertynames(lrt::LikelihoodRatioTest, private::Bool=false) return (:deviance, :formulas, :models, :pvalues, :tests) end +""" + StatsAPI.pvalue(lrt::LikelihoodRatioTest) + +Extract the p-value associated with the a likelihood ratio test. + +For `LikelihoodRatioTest`s containing more one model comparison, i.e. more than two models, +this throws an error because it is unclear which p-value is desired. + +To get p-values for multiple tests, use `lrt.pvalues`. +""" +function StatsAPI.pvalue(lrt::LikelihoodRatioTest) + pvalues = lrt.pvalues + if length(pvalues) > 1 + throw(ArgumentError("Cannot extract **only one** p-value from a multiple test result.")) + end + + return only(lrt.pvalues) +end + function Base.getproperty(lrt::LikelihoodRatioTest, s::Symbol) if s == :dof lrt.models.dof diff --git a/test/likelihoodratiotest.jl b/test/likelihoodratiotest.jl index 28f066122..131cf5ab4 100644 --- a/test/likelihoodratiotest.jl +++ b/test/likelihoodratiotest.jl @@ -58,12 +58,13 @@ end @test [deviance(fm0), deviance(fm1)] == lrt.deviance @test deviance(fm0) - deviance(fm1) == only(lrt.tests.deviancediff) @test only(lrt.tests.dofdiff) == 1 - @test sum(map(length,lrt.tests)) == 3 - @test sum(map(length,lrt.pvalues)) == 1 - @test sum(map(length,lrt.models)) == 4 + @test sum(length,lrt.tests) == 3 + @test sum(length,lrt.pvalues) == 1 + @test sum(length,lrt.models) == 4 @test length(lrt.formulae) == 2 show(IOBuffer(),lrt); @test :pvalues in propertynames(lrt) + @test only(lrt.pvalues) == pvalue(lrt) lrt = likelihoodratiotest(lm1,fm1) @test lrt.deviance ≈ likelihoodratiotest(lm1.model,fm1).deviance @@ -73,6 +74,9 @@ end @test occursin("-2 logLik", shown) @test !occursin("deviance", shown) + lrt = likelihoodratiotest(lm0, fm0, fm1) + @test_throws ArgumentError pvalue(lrt) + # non nested FE between non-mixed and mixed @test_throws ArgumentError likelihoodratiotest(lm1, fm0) From 4544d6a264c0c12db01da7b9f7150f4b2a88212b Mon Sep 17 00:00:00 2001 From: Phillip Alday Date: Tue, 18 Mar 2025 16:06:45 -0500 Subject: [PATCH 2/4] JuliaFormatter Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com> --- src/MixedModels.jl | 3 ++- src/likelihoodratiotest.jl | 6 +++++- 2 files changed, 7 insertions(+), 2 deletions(-) diff --git a/src/MixedModels.jl b/src/MixedModels.jl index 4a49c657f..d0ec6df98 100644 --- a/src/MixedModels.jl +++ b/src/MixedModels.jl @@ -31,7 +31,8 @@ using StaticArrays: StaticArrays, SVector using Statistics: Statistics, mean, quantile, std using StatsAPI: StatsAPI, aic, aicc, bic, coef, coefnames, coeftable, confint, deviance using StatsAPI: dof, dof_residual, fit, fit!, fitted, isfitted, islinear, leverage -using StatsAPI: loglikelihood, meanresponse, modelmatrix, nobs, pvalue, predict, r2, residuals +using StatsAPI: + loglikelihood, meanresponse, modelmatrix, nobs, pvalue, predict, r2, residuals using StatsAPI: response, responsename, stderror, vcov, weights using StatsBase: StatsBase, CoefTable, model_response, summarystats using StatsFuns: log2π, normccdf diff --git a/src/likelihoodratiotest.jl b/src/likelihoodratiotest.jl index 5781e4546..2809fa500 100644 --- a/src/likelihoodratiotest.jl +++ b/src/likelihoodratiotest.jl @@ -39,7 +39,11 @@ To get p-values for multiple tests, use `lrt.pvalues`. function StatsAPI.pvalue(lrt::LikelihoodRatioTest) pvalues = lrt.pvalues if length(pvalues) > 1 - throw(ArgumentError("Cannot extract **only one** p-value from a multiple test result.")) + throw( + ArgumentError( + "Cannot extract **only one** p-value from a multiple test result." + ), + ) end return only(lrt.pvalues) From 41bb65271152b69c469b00579ca989cd0b0a688f Mon Sep 17 00:00:00 2001 From: Phillip Alday Date: Tue, 18 Mar 2025 16:26:07 -0500 Subject: [PATCH 3/4] version bump + news --- NEWS.md | 7 ++++++- Project.toml | 2 +- 2 files changed, 7 insertions(+), 2 deletions(-) diff --git a/NEWS.md b/NEWS.md index 0582a7cc5..049757e33 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,3 +1,7 @@ +MixedModels v4.33.0 Release Notes +============================== +- `LikelihoodRatioTest` now extends `StatsAPI.HypothesisTest` and provides a method for `StatsAPI.pvalue`. [#814] + MixedModels v4.32.0 Release Notes ============================== - Added `lmm` and `glmm` as convenience wrappers for `fit(LinearMixedModel, ...)` and `fit(GeneralizedLinearMixedModel, ...)` respectively [#810] @@ -19,7 +23,7 @@ MixedModels v4.29.1 Release Notes MixedModels v4.29.0 Release Notes ============================== - Testbed for experimental support for using PRIMA as an optimization backend introduced via the experimental `prfit!` function. [#799] -- Julia compat bound raised to current 1.10, i.e. current LTS. [#799] +- Julia compat bound raised to current 1.10, i.e. current LTS. [#799] MixedModels v4.28.0 Release Notes ============================== @@ -612,3 +616,4 @@ Package dependencies [#801]: https://github.com/JuliaStats/MixedModels.jl/issues/801 [#802]: https://github.com/JuliaStats/MixedModels.jl/issues/802 [#810]: https://github.com/JuliaStats/MixedModels.jl/issues/810 +[#814]: https://github.com/JuliaStats/MixedModels.jl/issues/814 diff --git a/Project.toml b/Project.toml index 407e0b7e6..26dbea47b 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "MixedModels" uuid = "ff71e718-51f3-5ec2-a782-8ffcbfa3c316" author = ["Phillip Alday ", "Douglas Bates "] -version = "4.32.0" +version = "4.33.0" [deps] Arrow = "69666777-d1a9-59fb-9406-91d4454c9d45" From b91b75a3255d7ffc61a7d80b435e34e7bbd9d79e Mon Sep 17 00:00:00 2001 From: Phillip Alday Date: Wed, 19 Mar 2025 16:15:14 -0500 Subject: [PATCH 4/4] Apply suggestions from code review Co-authored-by: Alex Arslan --- src/likelihoodratiotest.jl | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/likelihoodratiotest.jl b/src/likelihoodratiotest.jl index 2809fa500..bbd6246d6 100644 --- a/src/likelihoodratiotest.jl +++ b/src/likelihoodratiotest.jl @@ -27,11 +27,11 @@ function Base.propertynames(lrt::LikelihoodRatioTest, private::Bool=false) end """ - StatsAPI.pvalue(lrt::LikelihoodRatioTest) + pvalue(lrt::LikelihoodRatioTest) -Extract the p-value associated with the a likelihood ratio test. +Extract the p-value associated with a likelihood ratio test. -For `LikelihoodRatioTest`s containing more one model comparison, i.e. more than two models, +For `LikelihoodRatioTest`s containing more than one model comparison, i.e. more than two models, this throws an error because it is unclear which p-value is desired. To get p-values for multiple tests, use `lrt.pvalues`. @@ -46,7 +46,7 @@ function StatsAPI.pvalue(lrt::LikelihoodRatioTest) ) end - return only(lrt.pvalues) + return only(pvalues) end function Base.getproperty(lrt::LikelihoodRatioTest, s::Symbol)