support rank deficiency and reduce dependency on Makie recipes (#86)

* start work on rank deficiency

* initial pass at conversion to recipeless implementation

* version and compat bump

* Indexable

* coefplot

* ridgeplot

* ridge2d

* YAS

* CI update

* handle rank deficiency

* add tests

* YAS

* fix a scope warning

* eliminate progress message in precompilation

.github/workflows/ci.yml

@@ -22,12 +22,12 @@ jobs:
     runs-on: ${{ matrix.os }}
     strategy:
       matrix:
-        version: [1.8]
+        version: [1.9]
         arch: [x64]
         os: [ubuntu-22.04] # macos-10.15, windows-2019
     steps:
       - name: Checkout
-        uses: actions/checkout@v2
+        uses: actions/checkout@v4
       - name: Julia Setup
         uses: julia-actions/setup-julia@v1
         with:
@@ -44,7 +44,7 @@ jobs:
         env:
           GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }}
       - name: Percy Upload
-        if: ${{ matrix.version == '1.8' }}
+        if: ${{ matrix.version == '1.9' }}
         run: |
           ls ./test/output/ # useful for debugging
           npx @percy/cli upload ./test/output
@@ -53,4 +53,6 @@ jobs:
       - name: Coverage Process
         uses: julia-actions/julia-processcoverage@v1
       - name: Coverage Upload
-        uses: codecov/codecov-action@v3
+        uses: codecov/codecov-action@v4
+        with:
+          token: ${{ secrets.CODECOV_TOKEN }}

.github/workflows/documenter.yml

@@ -25,7 +25,7 @@ jobs:
     # Run on push's or non-draft PRs
     if: (github.event_name == 'push') || (github.event.pull_request.draft == false)
     steps:
-      - uses: actions/checkout@v2
+      - uses: actions/checkout@v4
       - uses: julia-actions/setup-julia@v1
         with:
           version: 1

.github/workflows/style.yml

@@ -18,7 +18,7 @@ jobs:
     runs-on: ubuntu-latest
     strategy:
       matrix:
-        julia-version: [1.8]
+        julia-version: [1.9]
     steps:
       - uses: julia-actions/setup-julia@latest
         with:

Project.toml

@@ -1,7 +1,7 @@
 name = "MixedModelsMakie"
 uuid = "b12ae82c-6730-437f-aff9-d2c38332a376"
 authors = ["Phillip Alday <[email protected]>", "Douglas Bates <[email protected]>", "contributors"]
-version = "0.3.28"
+version = "0.4.0"
 
 [deps]
 BSplineKit = "093aae92-e908-43d7-9660-e50ee39d5a0a"
@@ -19,11 +19,11 @@ StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
 [compat]
 BSplineKit = "0.15, 0.16, 0.17"
 DataFrames = "1"
-Distributions = "0.21, 0.22, 0.23, 0.24, 0.25"
+Distributions = "0.25"
 KernelDensity = "0.6.3"
 Makie = "0.20"
 MixedModels = "4.14"
 PrecompileTools = "1"
 SpecialFunctions = "1, 2"
-StatsBase = "0.31, 0.32, 0.33, 0.34"
-julia = "1.8"
+StatsBase = "0.33, 0.34"
+julia = "1.9"

src/MixedModelsMakie.jl

@@ -36,6 +36,10 @@ export RanefInfo,
        zetaplot,
        zetaplot!
 
+# from https://github.com/MakieOrg/Makie.jl/issues/2992
+const Indexable = Union{Makie.Figure,Makie.GridLayout,Makie.GridPosition,
+                        Makie.GridSubposition}
+
 include("utilities.jl")
 include("shrinkage.jl")
 include("caterpillar.jl")
@@ -48,7 +52,7 @@ include("recipes.jl")
 
 @setup_workload begin
     model = fit(MixedModel, @formula(reaction ~ 1 + days + (1 + days | subj)),
-                MixedModels.dataset(:sleepstudy))
+                MixedModels.dataset(:sleepstudy); progress=false)
     @compile_workload begin
         caterpillar(model)
         coefplot(model)

src/coefplot.jl

@@ -1,65 +1,59 @@
 """
-    coefplot(x::MixedModel;
-             conf_level=0.95, vline_at_zero=true, show_intercept=true, attributes...)
-    coefplot(x::MixedModelBootstrap;
-             conf_level=0.95, vline_at_zero=true, show_intercept=true, attributes...)
+    coefplot(x::Union{MixedModel,MixedModelBootstrap}; kwargs...)::Figure
+    coefplot!(fig::$(Indexable), x::Union{MixedModel,MixedModelBootstrap}; 
+              kwargs...)
+    coefplot!(ax::Axis, Union{MixedModel,MixedModelBootstrap};
+              conf_level=0.95, vline_at_zero=true, show_intercept=true, attributes...)
 
 Create a coefficient plot of the fixed-effects and associated confidence intervals.
 
-!!! note
-    This functionality is implemented using [Makie recipes](https://makie.juliaplots.org/v0.15.0/recipes.html)
-    and thus there are also additional auto-generated methods for `coefplot` and `coefplot!` that may be useful
-    when constructing more complex figures.
+Inestimable coefficients (coefficients removed by pivoting in the rank deficient case) are excluded.
+
+`attributes` are passed onto `scatter!` and `errorbars!`.
+
+The mutating methods return the original object.
+
+!!! note 
+    Inestimable coefficients (coefficients removed by pivoting in the rank deficient case) 
+    are excluded.
 """
-function coefplot(x::Union{MixedModel,MixedModelBootstrap};
-                  conf_level=0.95,
-                  vline_at_zero=true,
-                  show_intercept=true,
-                  attributes...)
+function coefplot(x::Union{MixedModel,MixedModelBootstrap}; show_intercept=true, kwargs...)
     # need to guarantee a min height of 150
     fig = Figure(; size=(640, max(150, 75 * _npreds(x; show_intercept))))
-    ax = Axis(fig[1, 1])
-    pl = coefplot!(ax, x; conf_level, vline_at_zero, show_intercept, attributes...)
-    return Makie.FigureAxisPlot(fig, ax, pl)
+    coefplot!(fig, x; kwargs...)
+    return fig
 end
 
-@recipe(CoefPlot, x) do scene
-    return Attributes(; conf_level=0.95, vline_at_zero=true, show_intercept=true)
+"""$(@doc coefplot)"""
+function coefplot!(fig::Indexable, x::Union{MixedModel,MixedModelBootstrap}; kwargs...)
+    ax = Axis(fig[1, 1])
+    coefplot!(ax, x; kwargs...)
+    return fig
 end
 
-function Makie.plot!(ax::Axis, P::Type{<:CoefPlot}, allattrs::Makie.Attributes, x)
-    plot = Makie.plot!(ax.scene, P, allattrs, x)
+"""$(@doc coefplot)"""
+function coefplot!(ax::Axis, x::Union{MixedModel,MixedModelBootstrap};
+                   conf_level=0.95,
+                   vline_at_zero=true,
+                   show_intercept=true,
+                   attributes...)
+    ci = confint_table(x, conf_level; show_intercept)
+    y = nrow(ci):-1:1
+    xvals = ci.estimate
+    xlabel = @sprintf "Estimate and %g%% confidence interval" conf_level * 100
+
+    attributes = merge((; xlabel), attributes)
+    ax.xlabel = attributes.xlabel
+
+    scatter!(ax, xvals, y; attributes...)
+    errorbars!(ax, xvals, y, xvals .- ci.lower, ci.upper .- xvals;
+               direction=:x, attributes...)
+    vline_at_zero && vlines!(ax, 0; color=(:black, 0.75), linestyle=:dash)
 
-    if haskey(allattrs, :title)
-        ax.title = allattrs.title[]
-    end
-    if haskey(allattrs, :xlabel)
-        ax.xlabel = allattrs.xlabel[]
-    else
-        ax.xlabel = @sprintf "Estimate and %g%% confidence interval" (allattrs.conf_level[] *
-                                                                      100)
-    end
-    if haskey(allattrs, :ylabel)
-        ax.ylabel = allattrs.ylabel[]
-    end
     reset_limits!(ax)
-    show_intercept = allattrs.show_intercept[]
     cn = _coefnames(x; show_intercept)
     nticks = _npreds(x; show_intercept)
     ax.yticks = (nticks:-1:1, cn)
     ylims!(ax, 0, nticks + 1)
-    allattrs.vline_at_zero[] && vlines!(ax, 0; color=(:black, 0.75), linestyle=:dash)
-    return plot
-end
-
-function Makie.plot!(plot::CoefPlot{<:Tuple{Union{MixedModel,MixedModelBootstrap}}})
-    model_or_boot, conf_level = plot[1][], plot.conf_level[]
-    ci = confint_table(model_or_boot, conf_level)
-    plot.show_intercept[] || filter!(:coefname => !=("(Intercept)"), ci)
-    y = nrow(ci):-1:1
-    xvals = ci.estimate
-    scatter!(plot, xvals, y)
-    errorbars!(plot, xvals, y, xvals .- ci.lower, ci.upper .- xvals; direction=:x)
-
-    return plot
+    return ax
 end

src/ridge.jl

@@ -1,15 +1,9 @@
-
-@recipe(RidgePlot, x) do scene
-    return Attributes(;
-                      conf_level=0.95,
-                      vline_at_zero=true,
-                      show_intercept=true)
-end
-
 """
-    ridgeplot(x::MixedModelBootstrap;
-              conf_level=0.95, vline_at_zero=true, show_intercept=true,
-              attributes...)
+    ridgeplot(x::Union{MixedModel,MixedModelBootstrap}; kwargs...)::Figure
+    ridgeplot!(fig::$(Indexable), x::Union{MixedModel,MixedModelBootstrap}; 
+              kwargs...)
+    ridgeplot!(ax::Axis, Union{MixedModel,MixedModelBootstrap};
+              conf_level=0.95, vline_at_zero=true, show_intercept=true, attributes...)
 
 Create a ridge plot for the bootstrap samples of the fixed effects.
 
@@ -18,84 +12,89 @@ Densities are normalized so that the maximum density is always 1.
 The highest density interval corresponding to `conf_level` is marked with a bar at the bottom of each density.
 Setting `conf_level=missing` removes the markings for the highest density interval.
 
-!!! note
-    This functionality is implemented using [Makie recipes](https://makie.juliaplots.org/v0.15.0/recipes.html)
-    and thus there are also additional auto-generated methods for `ridgeplot` and `ridgeplot!` that may be useful
-    when constructing more complex figures.
+`attributes` are passed onto [`coefplot`](@ref), `band!` and `lines!`.
+
+The mutating methods return the original object.
+
+!!! note 
+    Inestimable coefficients (coefficients removed by pivoting in the rank deficient case) 
+    are excluded.
 """
-function ridgeplot(x::MixedModelBootstrap;
-                   conf_level=0.95,
-                   vline_at_zero=true,
-                   show_intercept=true,
-                   attributes...)
+function ridgeplot(x::MixedModelBootstrap; show_intercept=true, kwargs...)
     # need to guarantee a min height of 200
     fig = Figure(; size=(640, max(200, 100 * _npreds(x; show_intercept))))
-    ax = Axis(fig[1, 1])
-    if !ismissing(conf_level)
-        pl = coefplot!(ax, x; conf_level, vline_at_zero, show_intercept, color=:black,
-                       attributes...)
-    end
-    pl = ridgeplot!(ax, x; vline_at_zero, conf_level, show_intercept, attributes...)
-    return Makie.FigureAxisPlot(fig, ax, pl)
+    return ridgeplot!(fig, x; show_intercept, kwargs...)
 end
 
-function Makie.plot!(ax::Axis, P::Type{<:RidgePlot}, allattrs::Makie.Attributes, x)
-    plot = Makie.plot!(ax.scene, P, allattrs, x)
+"""$(@doc ridgeplot)"""
+function ridgeplot!(fig::Indexable, x::MixedModelBootstrap; kwargs...)
+    ax = Axis(fig[1, 1])
+    ridgeplot!(ax, x; kwargs...)
+    return fig
+end
 
-    if haskey(allattrs, :title)
-        ax.title = allattrs.title[]
-    end
-    if haskey(allattrs, :xlabel)
-        ax.xlabel = allattrs.xlabel[]
+"""
+    _color(s::Symbol)
+    _color(p::Pair)
+    
+Extract the color part out of either a color name or a `(color, alpha)` pair.
+"""
+_color(s) = s
+_color(p::Pair) = first(p)
+
+"""$(@doc ridgeplot)"""
+function ridgeplot!(ax::Axis, x::MixedModelBootstrap;
+                    conf_level=0.95,
+                    vline_at_zero=true,
+                    show_intercept=true,
+                    attributes...)
+    xlabel = if !ismissing(conf_level)
+        @sprintf "Normalized bootstrap density and %g%% confidence interval" (conf_level *
+                                                                              100)
     else
-        lab = if !ismissing(allattrs.conf_level[])
-            @sprintf "Normalized bootstrap density and %g%% confidence interval" (allattrs.conf_level[] *
-                                                                                  100)
-        else
-            "Normalized bootstrap density"
-        end
-        ax.xlabel = lab
-    end
-    if haskey(allattrs, :ylabel)
-        ax.ylabel = allattrs.ylabel[]
+        "Normalized bootstrap density"
     end
-    reset_limits!(ax)
-    show_intercept = allattrs.show_intercept[]
 
-    # check conf_level so that we don't double print
-    # if coefplot took care of it for us
-    if ismissing(allattrs.conf_level[])
-        cn = _coefnames(x; show_intercept)
-        nticks = _npreds(x; show_intercept)
-        ax.yticks = (nticks:-1:1, cn)
-        ylims!(ax, 0, nticks + 1)
-        allattrs.vline_at_zero[] && vlines!(ax, 0; color=(:black, 0.75), linestyle=:dash)
+    if !ismissing(conf_level)
+        coefplot!(ax, x; conf_level, vline_at_zero, show_intercept, color=:black,
+                  attributes...)
     end
-    return plot
-end
 
-function Makie.plot!(plot::RidgePlot{<:Tuple{MixedModelBootstrap}})
-    boot, conf_level = plot[1][], plot.conf_level[]
-    df = transform!(DataFrame(boot.β), :coefname => ByRow(string) => :coefname)
-    plot.show_intercept[] || filter!(:coefname => !=("(Intercept)"), df)
+    attributes = merge((; xlabel, color=:black), attributes)
+    band_attributes = merge(attributes, (; color=(_color(attributes.color), 0.3)))
+
+    ax.xlabel = attributes.xlabel
+
+    df = transform!(DataFrame(x.β), :coefname => ByRow(string) => :coefname)
+    filter!(:coefname => in(_coefnames(x; show_intercept)), df)
     group = :coefname
     densvar = :β
     gdf = groupby(df, group)
 
-    y = length(gdf):-1:1
-
     dens = combine(gdf, densvar => kde => :kde)
 
     for (offset, row) in enumerate(reverse(eachrow(dens)))
         # multiply by 0.95 so that the ridges don't overlap
         dd = 0.95 * row.kde.density ./ maximum(row.kde.density)
-        lower = Observable(Point2f.(row.kde.x, offset))
-        upper = Observable(Point2f.(row.kde.x, dd .+ offset))
-        band!(plot, lower, upper; color=(:black, 0.3))
-        lines!(plot, upper; color=(:black, 1.0))
+        lower = Point2f.(row.kde.x, offset)
+        upper = Point2f.(row.kde.x, dd .+ offset)
+        band!(ax, lower, upper; band_attributes...)
+        lines!(ax, upper; attributes...)
+    end
+
+    # check conf_level so that we don't double print
+    # if coefplot took care of it for us
+    if ismissing(conf_level)
+        cn = _coefnames(x; show_intercept)
+        nticks = _npreds(x; show_intercept)
+        ax.yticks = (nticks:-1:1, cn)
+        ylims!(ax, 0, nticks + 1)
+        vline_at_zero && vlines!(ax, 0; color=(:black, 0.75), linestyle=:dash)
     end
 
-    return plot
+    reset_limits!(ax)
+
+    return ax
 end
 
 # """