fix display and confint of bootstrap for models without dispersion parameter (#861)

* fix display and confint of bootstrap for models without dispersion parameter

* NEWS

* patch not minor

(cherry picked from commit ec7b986)

Author: palday

NEWS.md

@@ -1,3 +1,34 @@
+MixedModels v5.0.2 Release Notes
+==============================
+- The default display and `confint` methods for bootstrap results from models without dispersion parameters has been fixed. [#861]
+
+MixedModels v5.0.1 Release Notes
+==============================
+- Fixes a method error with `Grouping()` contrasts when used with recent CategoricalArray releases. [#860]
+
+MixedModels v5.0.0 Release Notes
+==============================
+- Optimization is now performed _without constraints_. In a post-fitting step, the Cholesky factor is canonicalized to have non-negative diagonal elements. [#840]
+- The default optimizer has changed to NLopt's implementation of NEWUOA where possible. NLopt's implementation fails on 1-dimensional problems, so in the case of a single, scalar random effect, BOBYQA is used instead. In the future, the default optimizer backend will likely change to PRIMA and NLopt support will be moved to an extension. Blocking this change in backend is an issue with PRIMA.jl when running in VSCode's built-in REPL on Linux. [#840]
+- [BREAKING] Support for constrained optimization has been completely removed, i.e. the field `lowerbd` has been removed from `OptSummary`. [#849]
+- [BREAKING] The deprecated `use_threads` kwarg has been dropped from `parametricbootstrap`. It had been a no-op since v4.10.0. [#841]
+- [BREAKING] The deprecated `hide_progress` kwarg has been dropped from `parametricbootstrap`. It had been replaced by `progress` since v4.22.0. [#841]
+- [BREAKING] A fitlog is always kept -- the deprecated keyword argument `thin` has been removed as has the `fitlog` keyword argument. [#850]
+- The fitlog is now stored as Tables.jl-compatible column table. [#850]
+- Internal code around the default optimizer has been restructured. In particular, the NLopt backend has been moved to a submodule, which will make it easier to move it to an extension if we promote another backend to the default. [#853]
+- Internal code around optimization in profiling has been restructuring so that fitting done during calls to `profile` respect the `backend` and `optimizer` settings. [#853]
+- The `prfit!` convenience function has been removed. [#853]
+- The `dataset` and `datasets` functions have been removed. They are now housed in `MixedModelsDatasets`.[#854]
+- The local implementation of `fulldummy` and the nesting syntax has been removed and a dependency on RegressionFormulae.jl for their implementation has been added. [#855]
+- One argument `predict(::GeneralizedLinearMixedModel)`, i.e. prediction on the original data, now supports the `type` keyword argument. [#856]
+- `isnested(A::ReMat, B::ReMat)` is now a method of `StatsModels.isnested`.[#858]
+- [BREAKING ]`likelihoodratiotest` has been reworked to be a thin wrapper around `StatsModels.lrtest`. The historical difference in behavior in terms of nesting checks created some confusion. Users advanced enough to create models with non-obvious nesting are assumed to be advanced enough to manually compute the likelihood ratio test. The function `likelihoodratiotest` and associated `LikelihoodRatioTest` type (now with a type parameter for number of models) has been kept to enable printing of test results with model formulae. Most users should not notice a difference in behavior, but the display has been slightly changed and the internal field structure has changed.[#858]
+- Failures to fit a spline in profiling now generates a more helpful error message. [#857]
+
+MixedModels v4.48.1 Release Notes
+==============================
+- [BACKPORT] The default display and `confint` methods for bootstrap results from models without dispersion parameters has been fixed. [#861]
+
 MixedModels v4.38.0 Release Notes
 ==============================
 - Experimental support for evaluating `FiniteDiff.finite_difference_gradient` and `FiniteDiff.finite_difference_hessian of the objective of a fitted `LinearMixedModel`. [#842]
@@ -655,5 +686,16 @@ Package dependencies
 [#828]: https://github.com/JuliaStats/MixedModels.jl/issues/828
 [#829]: https://github.com/JuliaStats/MixedModels.jl/issues/829
 [#836]: https://github.com/JuliaStats/MixedModels.jl/issues/836
+[#840]: https://github.com/JuliaStats/MixedModels.jl/issues/840
 [#841]: https://github.com/JuliaStats/MixedModels.jl/issues/841
 [#842]: https://github.com/JuliaStats/MixedModels.jl/issues/842
+[#849]: https://github.com/JuliaStats/MixedModels.jl/issues/849
+[#850]: https://github.com/JuliaStats/MixedModels.jl/issues/850
+[#853]: https://github.com/JuliaStats/MixedModels.jl/issues/853
+[#854]: https://github.com/JuliaStats/MixedModels.jl/issues/854
+[#855]: https://github.com/JuliaStats/MixedModels.jl/issues/855
+[#856]: https://github.com/JuliaStats/MixedModels.jl/issues/856
+[#857]: https://github.com/JuliaStats/MixedModels.jl/issues/857
+[#858]: https://github.com/JuliaStats/MixedModels.jl/issues/858
+[#860]: https://github.com/JuliaStats/MixedModels.jl/issues/860
+[#861]: https://github.com/JuliaStats/MixedModels.jl/issues/861

Project.toml

@@ -1,7 +1,7 @@
 name = "MixedModels"
 uuid = "ff71e718-51f3-5ec2-a782-8ffcbfa3c316"
 author = ["Phillip Alday <me@phillipalday.com>", "Douglas Bates <dmbates@gmail.com>"]
-version = "4.38.0"
+version = "4.38.1"
 
 [deps]
 Arrow = "69666777-d1a9-59fb-9406-91d4454c9d45"

src/bootstrap.jl

@@ -368,14 +368,16 @@ function StatsBase.confint(
     tbl = Table(bsamp.tbl)
     lower = T[]
     upper = T[]
-    v = similar(tbl.σ)
-    par = sort!(
-        collect(
-            filter(
-                k -> !(startswith(string(k), 'θ') || string(k) == "obj"), propertynames(tbl)
-            ),
-        ),
-    )
+    v = similar(tbl.σ, T)
+    par = filter(collect(propertynames(tbl))) do k
+        k = string(k)
+        # σ is missing in models without a dispersion parameter
+        if k == "σ" && Missing <: eltype(tbl.σ)
+            return false
+        end
+        return !startswith(k, 'θ') && k != "obj"
+    end
+    sort!(par)
     tails = [(1 - level) / 2, (1 + level) / 2]
     for p in par
         if method === :shortest
@@ -643,6 +645,8 @@ push! `σ` times the row lengths (σs) and the inner products of normalized rows
 """
 function σρ!(v::AbstractVector{<:Union{T,Missing}}, t::LowerTriangular, σ) where {T}
     dat = t.data
+    # for models without a dispersion parameter, σ is missing, but for the math below we can treat it as 1
+    σ = coalesce(σ, one(T))
     for i in axes(dat, 1)
         ssqr = zero(T)
         for j in 1:i

test/bootstrap.jl

@@ -200,15 +200,13 @@ end
         bs = parametricbootstrap(StableRNG(42), 100, gm0; progress=false)
         # make sure we're not copying
         @test length(bs.lowerbd) == length(gm0.θ)
-        bsci = filter!(:type => ==("β"), DataFrame(shortestcovint(bs)))
-        ciwidth = 2 .* stderror(gm0)
-        waldci = DataFrame(coef=fixefnames(gm0),
-                           lower=fixef(gm0) .- ciwidth,
-                           upper=fixef(gm0) .+ ciwidth)
+        bsci = confint(bs)
+        waldci = confint(gm0)
 
         # coarse tolerances because we're not doing many bootstrap samples
-        @test all(isapprox.(bsci.lower, waldci.lower; atol=0.5))
-        @test all(isapprox.(bsci.upper, waldci.upper; atol=0.5))
+        # end-1 because the bootstrap CIs include the variance component
+        @test all(isapprox.(collect(bsci.lower)[1:end-1], collect(waldci.lower); atol=0.5))
+        @test all(isapprox.(collect(bsci.upper)[1:end-1], collect(waldci.upper); atol=0.5))
 
         σbar = mean(MixedModels.tidyσs(bs)) do x; x.σ end
         @test σbar ≈ 0.56 atol=0.1