warning for GLMM with dispersion parameter (#373)* warning for GLMM with dispersion parameter(cherry picked from commit 66b8863)

Author: palday

src/generalizedlinearmixedmodel.jl

@@ -233,9 +233,9 @@ GeneralizedLinearMixedModel(f::FormulaTerm, tbl,
 GeneralizedLinearMixedModel(f::FormulaTerm, tbl::Tables.ColumnTable,
                             d::Normal,
                             l::IdentityLink;
-                            wts = [],
-                            offset = [],
-                            contrasts = Dict{Symbol,Any}()) =
+                            wts = wts,
+                            offset = offset,
+                            contrasts = contrasts) =
     throw(ArgumentError("use LinearMixedModel for Normal distribution with IdentityLink"))
 function GeneralizedLinearMixedModel(f::FormulaTerm, tbl::Tables.ColumnTable,
                                      d::Distribution,
@@ -247,7 +247,14 @@ function GeneralizedLinearMixedModel(f::FormulaTerm, tbl::Tables.ColumnTable,
         d = Bernoulli()
     end
     (isa(d, Normal) && isa(l, IdentityLink)) &&
-        throw(ArgumentError("use LinearMixedModel for Normal distribution with IdentityLink"))
+    throw(ArgumentError("use LinearMixedModel for Normal distribution with IdentityLink"))
+
+    if !any(isa(d, dist) for dist in (Bernoulli, Binomial, Poisson))
+        @warn """Results for families with a dispersion parameter are not reliable. 
+                 It is best to avoid trying to fit such models in MixedModels until 
+                 the authors get a better understanding of those cases."""
+    end
+
     LMM = LinearMixedModel(f, tbl, contrasts = contrasts; wts = wts)
     y = copy(LMM.y)
         # the sqrtwts field must be the correct length and type but we don't know those

test/Project.toml

@@ -1,6 +1,7 @@
 [deps]
 BlockArrays = "8e7c35d0-a365-5155-bbbb-fb81a777f24e"
 DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
+GLM = "38e38edf-8417-5370-95a0-9cbb8c7f171a"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 NamedArrays = "86f7a689-2022-50b4-a561-43c23ac3c673"
 RData = "df47a6cb-8c03-5eed-afd8-b6050d6c41da"

test/pirls.jl

@@ -1,5 +1,7 @@
 using DataFrames, LinearAlgebra, MixedModels, RData, Test
 
+using GLM: SqrtLink
+
 if !@isdefined(dat) || !isa(dat, Dict{Symbol, DataFrame})
     const dat = Dict(Symbol(k) => v for (k, v) in
         load(joinpath(dirname(pathof(MixedModels)), "..", "test", "dat.rda")))
@@ -76,3 +78,13 @@ end
     #@test isapprox(sum(x -> sum(abs2, x), gm4.u), 196.8695297987013, atol=0.1)
     #@test isapprox(sum(gm4.resp.devresid), 220.92685781326136, atol=0.1)
 end
+
+@testset "dispersion" begin
+
+    form = @formula(Y ~ 1 + U + (1 + U | G))
+
+    @test_logs (:warn, r"dispersion parameter") GeneralizedLinearMixedModel(form, dat[:sleepstudy], Gamma())
+    @test_logs (:warn, r"dispersion parameter") GeneralizedLinearMixedModel(form, dat[:sleepstudy], InverseGaussian())
+    @test_logs (:warn, r"dispersion parameter") GeneralizedLinearMixedModel(form, dat[:sleepstudy], Normal(), SqrtLink())
+
+end