Enhance the tbl property of a parametricbootstrap result

src/bootstrap.jl

@@ -253,7 +253,7 @@
     a breaking change.
 """
 function allpars(bsamp::MixedModelFitCollection{T}) where {T}
-    fits, λ, fcnames = bsamp.fits, bsamp.λ, bsamp.fcnames
+    (; fits, λ, fcnames) = bsamp
     npars = 2 + length(first(fits).β) + sum(map(k -> (k * (k + 1)) >> 1, size.(bsamp.λ, 2)))
     nresrow = length(fits) * npars
     cols = (
@@ -392,12 +392,12 @@
 end
 
 """
+    setθ!(bsamp::MixedModelFitCollection, θ::AbstractVector)
     setθ!(bsamp::MixedModelFitCollection, i::Integer)
 
 Install the values of the i'th θ value of `bsamp.fits` in `bsamp.λ`
 """
-function setθ!(bsamp::MixedModelFitCollection, i::Integer)
-    θ = bsamp.fits[i].θ
+function setθ!(bsamp::MixedModelFitCollection{T}, θ::AbstractVector{T}) where {T}
     offset = 0
     for (λ, inds) in zip(bsamp.λ, bsamp.inds)
         λdat = _getdata(λ)
@@ -410,6 +410,10 @@
     return bsamp
 end
 
+function setθ!(bsamp::MixedModelFitCollection, i::Integer)
+    return setθ!(bsamp, bsamp.θ[i])
+end
+
 _getdata(x::Diagonal) = x
 _getdata(x::LowerTriangular) = x.data
 
@@ -521,6 +525,7 @@
 
 """
     tidyσs(bsamp::MixedModelFitCollection)
+
 Return a tidy (row)table with the estimates of the variance components (on the standard deviation scale) spread into columns
 of `iter`, `group`, `column` and `σ`.
 """
@@ -544,40 +549,114 @@
     return result
 end
 
-function pbstrtbl(bsamp::MixedModelFitCollection{T}) where {T}
-    (; fits, λ, inds) = bsamp
-    row1 = first(fits)
-    cnms = [:obj, :σ]
-    pos = Dict{Symbol,UnitRange{Int}}(:obj => 1:1, :σ => 2:2)
-    βsz = length(row1.β)
-    append!(cnms, _generatesyms('β', βsz))
-    lastpos = 2 + βsz
-    pos[:β] = 3:lastpos
-    σsz = sum(m -> size(m, 1), bsamp.λ)
-    append!(cnms, _generatesyms('σ', σsz))
-    pos[:σs] = (lastpos + 1):(lastpos + σsz)
-    lastpos += σsz
-    θsz = length(row1.θ)
-    append!(cnms, _generatesyms('θ', θsz))
-    pos[:θ] = (lastpos + 1):(lastpos + θsz)
-    tblrowtyp = NamedTuple{(cnms...,),NTuple{length(cnms),T}}
-    val = sizehint!(tblrowtyp[], length(bsamp.fits))
-    v = Vector{T}(undef, length(cnms))
-    for (i, r) in enumerate(bsamp.fits)
-        v[1] = r.objective
-        v[2] = coalesce(r.σ, one(T))
-        copyto!(view(v, pos[:β]), r.β)
-        fill!(view(v, pos[:σs]), zero(T))
-        copyto!(view(v, pos[:θ]), r.θ)
-        setθ!(bsamp, i)
-        vpos = first(pos[:σs])
-        for l in λ
-            for λr in eachrow(l)
-                v[vpos] = r.σ * norm(λr)
-                vpos += 1
+_nρ(d::Diagonal) = 0
+_nρ(t::LowerTriangular) = kchoose2(size(t.data, 1))
+
+function σρnms(λ)
+    σsyms = _generatesyms('σ', sum(first ∘ size, λ))
+    ρsyms = _generatesyms('ρ', sum(_nρ, λ))
+    val = Symbol[]
+    for l in λ
+        for _ in axes(l, 1)
+            push!(val, popfirst!(σsyms))
+        end
+        for _ in 1:_nρ(l)
+            push!(val, popfirst!(ρsyms))
+        end
+    end
+    return val
+end
+
+"""
+    _offsets(bsamp::MixedModelBootstrap)
+
+Return a Tuple{4,Int} of offsets for β, σρ, and θ values in the table columns
+
+The initial `2` is for the `:obj` and `:σ` columns. The last element is the total number of columns.
+"""
+function _offsets(bsamp::MixedModelBootstrap)
+    (; fits, λ) = bsamp
+    (; β, θ) = first(fits)
+    σρ = σρnms(λ)
+    lβ = length(β)
+    lθ = length(θ)
+    syms = [:obj, :σ]
+    append!(syms, _generatesyms('β', lβ))
+    append!(syms, σρnms(λ))
+    append!(syms, _generatesyms('θ', lθ))
+    return Tuple(syms), cumsum((2, lβ, length(σρ), lθ))
+end
+
+function σρ!(v::AbstractVector, off::Integer, d::Diagonal, σ)
+    diag = d.diag
+    diag *= σ
+    return copyto!(v, off + 1, d.diag)
+end
+
+function σρ!(v::AbstractVector{T}, off::Integer, t::LowerTriangular{T}, σ::T) where {T}
+    dat = t.data
+    for i in axes(dat, 1)
+        ssqr = zero(T)
+        for j in 1:i
+            ssqr += abs2(dat[i, j])
+        end
+        len = sqrt(ssqr)
+        v[off += 1] = σ * len
+        if len > 0
+            for j in 1:i
+               dat[i, j] /= len
+            end
+        end
+    end
+    for i in axes(dat, 1)
+        for j in 1:(i - 1)
+            s = zero(T)
+            for k in 1:i
+                s += dat[i, k] * dat[j, k]
             end
+            v[off += 1] = s
         end
-        push!(val, tblrowtyp(v))
+    end
+    return v
+end
+
+function _allpars!(
+    v::AbstractVector{T},
+    bsamp::MixedModelBootstrap{T},
+    i::Integer,
+    offsets::NTuple{4,Int},
+) where {T}
+    fiti = bsamp.fits[i]
+    setθ!(bsamp, i)
+    λ = bsamp.λ
+    v[1] = fiti.objective
+    v[2] = σ = coalesce(fiti.σ, one(T))
+    off = 2
+    for b in fiti.β
+        v[off += 1] = b
+    end
+#    copyto!(v, 3, fiti.β)
+    off = offsets[2]
+    for lam in λ
+        σρ!(v, off, lam, σ)
+        off += size(lam, 1) + _nρ(lam)
+    end
+    off = offsets[3]
+    for t in fiti.θ
+        v[off += 1] = t
+    end
+#    copyto!(v, offsets[3] + 1, fiti.θ)
+    return v
+end
+
+function pbstrtbl(bsamp::MixedModelFitCollection{T}) where {T}
+    fits = bsamp.fits
+    syms, offsets = _offsets(bsamp)
+    nsym = length(syms)
+    val = NamedTuple{syms,NTuple{nsym,T}}[]
+    v = sizehint!(Vector{T}(undef, nsym), size(fits, 1))
+    for i in axes(fits, 1)
+        push!(val, NamedTuple{syms,NTuple{nsym,T}}(_allpars!(v, bsamp, i, offsets)))
     end
     return val
 end

src/profile/utilities.jl

@@ -24,10 +24,12 @@ Utility to generate a vector of Symbols of the form :<tag><index> from a tag and
 
 The indices are left-padded with zeros to allow lexicographic sorting.
 """
-function _generatesyms(tag::Char, len::Integer)
+function _generatesyms(tag::AbstractString, len::Integer)
     return Symbol.(string.(tag, lpad.(Base.OneTo(len), ndigits(len), '0')))
 end
 
+_generatesyms(tag::Char, len::Integer) = _generatesyms(string(tag), len)
+
 function TableColumns(m::LinearMixedModel{T}) where {T}
     nmvec = [:ζ]
     positions = Dict(:ζ => 1:1)