From 692051ddd2bc6a74b97e6d11fea1f9a8e3266be9 Mon Sep 17 00:00:00 2001 From: yufongpeng <54415349+yufongpeng@users.noreply.github.com> Date: Thu, 29 Feb 2024 02:55:14 +0800 Subject: [PATCH 1/3] Update AnovaBase to 0.8.0 --- Project.toml | 4 +-- src/AnovaFixedEffectModels.jl | 4 +-- src/anova.jl | 12 ++++---- src/fit.jl | 5 ++- src/io.jl | 57 ++++++++++++++++++++++++++++++++++- test/runtests.jl | 3 +- 6 files changed, 69 insertions(+), 16 deletions(-) diff --git a/Project.toml b/Project.toml index 713ff0f..086fe1f 100644 --- a/Project.toml +++ b/Project.toml @@ -16,11 +16,11 @@ StatsModels = "3eaba693-59b7-5ba5-a881-562e759f1c8d" Tables = "bd369af6-aec1-5ad0-b16a-f7cc5008161c" [compat] -AnovaBase = "0.7" +AnovaBase = "0.8" Distributions = "0.23, 0.24, 0.25" FixedEffectModels = "1.3, 1.4, 1.5, 1.6, 1.7, 1.8" Reexport = "0.2, 1" StatsBase = "0.33, 0.34" StatsModels = "0.7" Tables = "1.7" -julia = "1.6, 1.7, 1.8, 1.9" +julia = "1.6, 1.7, 1.8, 1.9, 1.10" diff --git a/src/AnovaFixedEffectModels.jl b/src/AnovaFixedEffectModels.jl index 1e2a00c..ace2c68 100644 --- a/src/AnovaFixedEffectModels.jl +++ b/src/AnovaFixedEffectModels.jl @@ -7,10 +7,10 @@ import StatsBase: fit!, fit import StatsModels: TableRegressionModel, vectorize, width, apply_schema, ModelFrame, ModelMatrix, columntable, asgn -using AnovaBase: select_super_interaction, extract_contrasts, canonicalgoodnessoffit, subformula, dof_asgn, lrt_nested, ftest_nested, _diff, _diffn, has_intercept +using AnovaBase: select_super_interaction, extract_contrasts, canonicalgoodnessoffit, subformula, dof_asgn, lrt_nested, ftest_nested, _diff, _diffn, testname import AnovaBase: anova, nestedmodels, anovatable, prednames, predictors, formula using Tables: columntable - +import Base: show export anova_lfe, lfe include("anova.jl") diff --git a/src/anova.jl b/src/anova.jl index e004c36..d0d42a5 100644 --- a/src/anova.jl +++ b/src/anova.jl @@ -66,10 +66,10 @@ function anova(::Type{FTest}, end elseif aovm.type == 2 fstat = ntuple(last(fullasgn) - offset) do fix - select1 = sort!(collect(select_super_interaction(fullpred, fix + offset))) - select2 = setdiff(select1, fix + offset) - select1 = findall(in(select1), fullasgn) - select2 = findall(in(select2), fullasgn) + s1 = sort!(collect(select_super_interaction(fullpred, fix + offset))) + s2 = setdiff(s1, fix + offset) + select1 = findall(in(s1), fullasgn) + select2 = findall(in(s2), fullasgn) (β[select1]' * (varβ[select1, select1] \ β[select1]) - β[select2]' * (varβ[select2, select2] \ β[select2])) / df[fix] end else @@ -83,7 +83,7 @@ function anova(::Type{FTest}, σ² = rss(aovm.model) / dfr devs = @. fstat * σ² * df pvalue = @. ccdf(FDist(df, dfr), abs(fstat)) - AnovaResult{FTest}(aovm, df, devs, fstat, pvalue, NamedTuple()) + AnovaResult(aovm, FTest, df, devs, fstat, pvalue, NamedTuple()) end # ================================================================================================================= @@ -102,7 +102,7 @@ function anova(::Type{FTest}, dev = rss.(models) # check comparable and nested check && @warn "Could not check whether models are nested: results may not be meaningful" - ftest_nested(NestedModels{M}(models), df, dfr, dev, last(dev) / last(dfr)) + ftest_nested(NestedModels(models), df, dfr, dev, last(dev) / last(dfr)) end function anova(::Type{FTest}, aovm::NestedModels{M}) where {M <: FixedEffectModel} diff --git a/src/fit.jl b/src/fit.jl index c2c83ef..15ea9fd 100644 --- a/src/fit.jl +++ b/src/fit.jl @@ -1,6 +1,5 @@ # ========================================================================================================== # Backend funcion -formula(model::T) where {T <: FixedEffectModel} = model.formula predictors(model::T) where {T <: FixedEffectModel} = model.formula_schema.rhs.terms # Variable dispersion @@ -19,10 +18,10 @@ Generate nested models from modeltype, formula and data. The null model will be function nestedmodels(modeltype::Type{FixedEffectModel}, f::FormulaTerm, data; null = true, kwargs...) fullm = lfe(f, data; kwargs...) predterms = predictors(fullm) - has_intercept(predterms) || (length(predterms) > 1 ? (predterms = predterms[2:end]) : throw(ArgumentError("Empty model is given!"))) + hasintercept(predterms) || (length(predterms) > 1 ? (predterms = predterms[2:end]) : throw(ArgumentError("Empty model is given!"))) feterms = filter(isfe, f.rhs) subm = ntuple(length(predterms) - 1) do i lfe(FormulaTerm(f.lhs, (predterms[1:i]..., feterms...)), data; kwargs...) end - NestedModels{FixedEffectModel}(null ? (lfe(FormulaTerm(f.lhs, (ConstantTerm(0), feterms...)), data; kwargs...), subm..., fullm) : (subm..., fullm)) + NestedModels(null ? (lfe(FormulaTerm(f.lhs, (ConstantTerm(0), feterms...)), data; kwargs...), subm..., fullm) : (subm..., fullm)) end diff --git a/src/io.jl b/src/io.jl index d7e17fd..5adfd97 100644 --- a/src/io.jl +++ b/src/io.jl @@ -38,4 +38,59 @@ function anovatable(aov::AnovaResult{<: NestedModels{<: FixedEffectModel, N}, FT ], ["DOF", "ΔDOF", "Res.DOF", "R²", "ΔR²", "R²_within", "ΔR²_within", "Res.SS", "Exp.SS", "F value", "Pr(>|F|)"], rownames, 11, 10) -end \ No newline at end of file +end + +function show(io::IO, anovamodel::FullModel{<: T}) where {T <: FixedEffectModel} + println(io, "FullModel for type $(anovamodel.type) test") + println(io) + println(io, "Predictors:") + println(io, join(prednames(anovamodel), ", ")) + println(io) + println(io, "Formula:") + println(io, anovamodel.model.formula) + println(io) + println(io, "Coefficients:") + show(io, coeftable(anovamodel.model)) +end + +function show(io::IO, anovamodel::NestedModels{M, N}) where {M <: FixedEffectModel, N} + println(io, "NestedModels with $N models") + println(io) + println(io, "Formulas:") + for(id, m) in enumerate(anovamodel.model) + println(io, "Model $id: ", m.formula) + end + println(io) + println(io, "Coefficients:") + show(io, coeftable(first(anovamodel.model))) + println(io) + N > 2 && print(io, " .\n" ^ 3) + show(io, coeftable(last(anovamodel.model))) +end + +# Show function that delegates to anovatable +function show(io::IO, aov::AnovaResult{<: FullModel{<: FixedEffectModel}, T}) where {T <: GoodnessOfFit} + at = anovatable(aov) + println(io, "Analysis of Variance") + println(io) + println(io, "Type $(anova_type(aov)) test / $(testname(T))") + println(io) + println(io, aov.anovamodel.model.formula) + println(io) + println(io, "Table:") + show(io, at) +end + +function show(io::IO, aov::AnovaResult{<: MultiAovModels{<: FixedEffectModel}, T}) where {T <: GoodnessOfFit} + at = anovatable(aov) + println(io,"Analysis of Variance") + println(io) + println(io, "Type $(anova_type(aov)) test / $(testname(T))") + println(io) + for(id, m) in enumerate(aov.anovamodel.model) + println(io, "Model $id: ", m.formula) + end + println(io) + println(io, "Table:") + show(io, at) +end \ No newline at end of file diff --git a/test/runtests.jl b/test/runtests.jl index 0fdff8b..c18a38a 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -52,7 +52,7 @@ isapprox(x::NTuple{N, Float64}, y::NTuple{N, Float64}, atol::NTuple{N, Float64} global aovl1 = AnovaGLM.anova(lm1) global aovf2 = AFE.anova(fem2, type = 3) global aovl2 = AnovaGLM.anova(lm2, type = 3) - global aovfs = AFE.anova(NestedModels{FixedEffectModel}(fem0, fem1)) + global aovfs = AFE.anova(NestedModels(fem0, fem1)) global aovfs2 = AFE.anova(fem0, fem1) global aovls = AnovaGLM.anova(lm0, lm1) @test !(@test_error test_show(aovf1)) @@ -68,7 +68,6 @@ isapprox(x::NTuple{N, Float64}, y::NTuple{N, Float64}, atol::NTuple{N, Float64} df = DataFrame(y = randn(1000), x = rand(1:5, 1000), z = rand(["1", "2"], 1000), t = 1:1000) fems1 = nestedmodels(FixedEffectModel, @formula(y ~ t + fe(z) + fe(x)), df) fems2 = nestedmodels(FixedEffectModel, @formula(y ~ z + t & fe(x)), df) - @test formula(fems1.model[1]).rhs == @formula(y ~ 0 + fe(z) + fe(x)).rhs @test AFE.predictors(fems2.model[2])[1] == InterceptTerm{true}() end end From 5d8f415302a0c149ea1237572ea409d714db9545 Mon Sep 17 00:00:00 2001 From: yufongpeng <54415349+yufongpeng@users.noreply.github.com> Date: Thu, 29 Feb 2024 02:55:41 +0800 Subject: [PATCH 2/3] Bump to version 0.2.5 --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 086fe1f..824d1b9 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "AnovaFixedEffectModels" uuid = "e4965305-65d6-464d-9c03-ae3e5cffadab" authors = ["yufongpeng <54415349+yufongpeng@users.noreply.github.com> and contributors"] -version = "0.2.4" +version = "0.2.5" [deps] AnovaBase = "946dddda-6a23-4b48-8e70-8e60d9b8d680" From ef30afc00f8cbc47cc969f1dff405d0fc99d7813 Mon Sep 17 00:00:00 2001 From: yufongpeng <54415349+yufongpeng@users.noreply.github.com> Date: Thu, 29 Feb 2024 03:04:28 +0800 Subject: [PATCH 3/3] Add more tests --- test/runtests.jl | 2 ++ 1 file changed, 2 insertions(+) diff --git a/test/runtests.jl b/test/runtests.jl index c18a38a..5942670 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -56,8 +56,10 @@ isapprox(x::NTuple{N, Float64}, y::NTuple{N, Float64}, atol::NTuple{N, Float64} global aovfs2 = AFE.anova(fem0, fem1) global aovls = AnovaGLM.anova(lm0, lm1) @test !(@test_error test_show(aovf1)) + @test !(@test_error test_show(aovf1.anovamodel)) @test !(@test_error test_show(aovf2)) @test !(@test_error test_show(aovfs)) + @test !(@test_error test_show(aovfs.anovamodel)) @test nobs(aovf1) == nobs(aovl1) @test last(dof(aovf1)) == dof(aovl1)[end - 1] @test isapprox(first(deviance(aovf2)), first(deviance(aovl2)))