From 028551bdfdcdefa2c23f48c16373ad66203dc0e1 Mon Sep 17 00:00:00 2001 From: Phillip Alday Date: Wed, 29 May 2024 12:51:33 -0500 Subject: [PATCH 1/3] pkg extension for Vcov --- .github/workflows/documenter.yml | 2 +- .github/workflows/minimum.yml | 2 +- Project.toml | 15 ++++++++++++--- ext/MixedModelsVcovExt.jl | 28 ++++++++++++++++++++++++++++ src/MixedModels.jl | 3 ++- 5 files changed, 44 insertions(+), 6 deletions(-) create mode 100644 ext/MixedModelsVcovExt.jl diff --git a/.github/workflows/documenter.yml b/.github/workflows/documenter.yml index 483fc4b8f..972c81db3 100644 --- a/.github/workflows/documenter.yml +++ b/.github/workflows/documenter.yml @@ -24,7 +24,7 @@ jobs: - uses: actions/checkout@v4 - uses: julia-actions/setup-julia@v2 with: - version: 1.8 + version: 1.9 - uses: julia-actions/cache@v2 with: cache-compiled: "true" diff --git a/.github/workflows/minimum.yml b/.github/workflows/minimum.yml index f7e37a65f..4731ce9db 100644 --- a/.github/workflows/minimum.yml +++ b/.github/workflows/minimum.yml @@ -20,7 +20,7 @@ jobs: runs-on: ${{ matrix.os }} strategy: matrix: - julia-version: [1.8] + julia-version: [1.9] julia-arch: [x64] os: [ubuntu-22.04, macos-11, windows-2019] steps: diff --git a/Project.toml b/Project.toml index 225f9aa54..e20c854b1 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "MixedModels" uuid = "ff71e718-51f3-5ec2-a782-8ffcbfa3c316" author = ["Phillip Alday ", "Douglas Bates ", "Jose Bayoan Santiago Calderon "] -version = "4.24.0" +version = "4.25.0" [deps] Arrow = "69666777-d1a9-59fb-9406-91d4454c9d45" @@ -28,6 +28,13 @@ StatsModels = "3eaba693-59b7-5ba5-a881-562e759f1c8d" StructTypes = "856f2bd8-1eba-4b0a-8007-ebc267875bd4" Tables = "bd369af6-aec1-5ad0-b16a-f7cc5008161c" TypedTables = "9d95f2ec-7b3d-5a63-8d20-e2491e220bb9" +Vcov = "ec2bfdc2-55df-4fc9-b9ae-4958c2cf2486" + +[weakdeps] +Vcov = "ec2bfdc2-55df-4fc9-b9ae-4958c2cf2486" + +[extensions] +MixedModelsVcovExt = "Vcov" [compat] Aqua = "0.8" @@ -61,7 +68,8 @@ Suppressor = "0.2" Tables = "1" Test = "1" TypedTables = "1" -julia = "1.8" +Vcov = "0.7, 0.8" +julia = "1.9" [extras] Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595" @@ -71,6 +79,7 @@ InteractiveUtils = "b77e0a4c-d291-57a0-90e8-8db25a27a240" StableRNGs = "860ef19b-820b-49d6-a774-d7a799459cd3" Suppressor = "fd094767-a336-5f1f-9728-57cf17d0bbfb" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" +Vcov = "ec2bfdc2-55df-4fc9-b9ae-4958c2cf2486" [targets] -test = ["Aqua", "DataFrames", "ExplicitImports", "InteractiveUtils", "StableRNGs", "Suppressor", "Test"] +test = ["Aqua", "DataFrames", "ExplicitImports", "InteractiveUtils", "StableRNGs", "Suppressor", "Test", "Vcov"] diff --git a/ext/MixedModelsVcovExt.jl b/ext/MixedModelsVcovExt.jl new file mode 100644 index 000000000..7c3446d5b --- /dev/null +++ b/ext/MixedModelsVcovExt.jl @@ -0,0 +1,28 @@ +module MixedModelsVcovExt + +using MixedModels +using MixedModels: StatsAPI +using Vcov: Vcov, CovarianceEstimator, VcovData + + +# in 0.8, Vcov uses invcrossmodelmatrix as an accessor function +# for a precomputed inv(crossmodelmatrix(m)) +if pkgversion(Vcov) >= v"0.8" + Vcov.invcrossmodelmatrix(m::MixedModel) = inv(crossmodelmatrix(m)) + + function Vcov.VcovData(m::MixedModel) + cm = crossmodelmatrix(m) + ddf = nobs(m) - round(Int, sum(leverage(m))) + return VcovData(modelmatrix(m), cm, inv(cm), residuals(m), ddf) + end + +end + +# to avoid method ambiguities, we define this separately for each +# relevant CovarianceEstimator +StatsAPI.vcov(m::MixedModel, c::Vcov.SimpleCovariance) = vcov(VcovData(m), c) +StatsAPI.vcov(m::MixedModel, c::Vcov.RobustCovariance) = vcov(VcovData(m), c) + +# Vcov.vcov(m::MixedModel, ::Vcov.SimpleCovariance) = vcov(m) + +end # module diff --git a/src/MixedModels.jl b/src/MixedModels.jl index d03de81cd..a97a962bc 100644 --- a/src/MixedModels.jl +++ b/src/MixedModels.jl @@ -28,7 +28,7 @@ using SparseArrays: SparseArrays, SparseMatrixCSC, SparseVector, dropzeros!, nnz using SparseArrays: nonzeros, nzrange, rowvals, sparse using StaticArrays: StaticArrays, SVector using Statistics: Statistics, mean, quantile, std -using StatsAPI: StatsAPI, aic, aicc, bic, coef, coefnames, coeftable, confint, deviance +using StatsAPI: StatsAPI, aic, aicc, bic, coef, coefnames, coeftable, confint, crossmodelmatrix, deviance using StatsAPI: dof, dof_residual, fit, fit!, fitted, isfitted, islinear, leverage using StatsAPI: loglikelihood, meanresponse, modelmatrix, nobs, predict, r2, residuals using StatsAPI: response, responsename, stderror, vcov, weights @@ -89,6 +89,7 @@ export @formula, condVar, condVartables, confint, + crossmodelmatrix, dataset, datasets, deviance, From bf0a6f1c7089394206b2cd0a637ed089368619c5 Mon Sep 17 00:00:00 2001 From: Phillip Alday Date: Wed, 29 May 2024 13:27:22 -0500 Subject: [PATCH 2/3] just the fixed effects --- ext/MixedModelsVcovExt.jl | 26 ++++++++++++++------------ 1 file changed, 14 insertions(+), 12 deletions(-) diff --git a/ext/MixedModelsVcovExt.jl b/ext/MixedModelsVcovExt.jl index 7c3446d5b..a9a99abda 100644 --- a/ext/MixedModelsVcovExt.jl +++ b/ext/MixedModelsVcovExt.jl @@ -2,18 +2,21 @@ module MixedModelsVcovExt using MixedModels using MixedModels: StatsAPI -using Vcov: Vcov, CovarianceEstimator, VcovData - +using MixedModels.StatsModels: hasintercept -# in 0.8, Vcov uses invcrossmodelmatrix as an accessor function -# for a precomputed inv(crossmodelmatrix(m)) -if pkgversion(Vcov) >= v"0.8" - Vcov.invcrossmodelmatrix(m::MixedModel) = inv(crossmodelmatrix(m)) +using Vcov: Vcov, CovarianceEstimator, VcovData - function Vcov.VcovData(m::MixedModel) - cm = crossmodelmatrix(m) - ddf = nobs(m) - round(Int, sum(leverage(m))) - return VcovData(modelmatrix(m), cm, inv(cm), residuals(m), ddf) +function Vcov.VcovData(m::MixedModel) + cm = crossmodelmatrix(m) + mm = modelmatrix(m) + # from just the fixed effects + ddf = nobs(m) - rank(m) - hasintercept(formula(m)) # dof_residual(m) # nobs(m) - round(Int, sum(leverage(m))) + resid = response(m) - mm * coef(m) # residuals(m) + @static if pkgversion(Vcov) >= v"0.8" + # Vcov 0.8 stores the inverse of crossmodelmatrix + return VcovData(mm, cm, inv(cm), resid, ddf) + else + return VcovData(mm, cm, resid, ddf) end end @@ -22,7 +25,6 @@ end # relevant CovarianceEstimator StatsAPI.vcov(m::MixedModel, c::Vcov.SimpleCovariance) = vcov(VcovData(m), c) StatsAPI.vcov(m::MixedModel, c::Vcov.RobustCovariance) = vcov(VcovData(m), c) - -# Vcov.vcov(m::MixedModel, ::Vcov.SimpleCovariance) = vcov(m) +# StatsAPI.vcov(m::MixedModel, c::Vcov.ClusterCovariance) = vcov(VcovData(m), c) end # module From d7d1a05be851bf2c17061cc6fa0e17fb487ec16c Mon Sep 17 00:00:00 2001 From: Phillip Alday Date: Wed, 29 May 2024 21:07:16 -0500 Subject: [PATCH 3/3] wip --- Project.toml | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index e20c854b1..0f8f670a9 100644 --- a/Project.toml +++ b/Project.toml @@ -31,10 +31,11 @@ TypedTables = "9d95f2ec-7b3d-5a63-8d20-e2491e220bb9" Vcov = "ec2bfdc2-55df-4fc9-b9ae-4958c2cf2486" [weakdeps] +GroupedArrays = "6407cd72-fade-4a84-8a1e-56e431fc1533" Vcov = "ec2bfdc2-55df-4fc9-b9ae-4958c2cf2486" [extensions] -MixedModelsVcovExt = "Vcov" +MixedModelsVcovExt = ["Vcov", "GroupedArrays"] [compat] Aqua = "0.8" @@ -45,6 +46,7 @@ DataFrames = "1" Distributions = "0.21, 0.22, 0.23, 0.24, 0.25" ExplicitImports = "1.3" GLM = "1.8.2" +GroupedArrays = "0.3" InteractiveUtils = "1" JSON3 = "1" LinearAlgebra = "1"