text:
@kernel function dup_kernel!(A, @Const(B))
    I = @index(Global)
    @inbounds A[I] = B[I] + B[I]
end

@kernel function split_kernel!(A, B, @Const(AB))
    I = @index(Global)
    @inbounds A[I] = AB[I].a
    @inbounds B[I] = AB[I].b
    nothing
end

@kernel function join_kernel!(AB, @Const(A), @Const(B))
    I = @index(Global)
    @inbounds AB[I] = (a=A[I], b=B[I])
end

@testset "StructArrays" begin
    TAs = ((Float64,  Array), (BigFloat, Array))
    if CUDA.has_cuda_gpu()
        TAs = (TAs..., (Float32, CuArray))
    end

    for (T, A) in TAs
        N = 4
        C = 3
        M = 5

        X = adapt(A, collect(reshape(oneunit(T):N*C*M, N, C, M)))

        tup = ntuple(i->view(X, :, i, :), C)
        b = StructArray{SMatrix{C, 1, T, C}}(tup)
        @test b isa StructArray
        @test Bennu.device(b) == Bennu.device(A)

        aos_b = collect(adapt(Array, b))
        @test aos_b isa Array
        @test Bennu.device(aos_b) == Bennu.device(Array)

        a = similar(b)
        aos_a = copy(aos_b)
        fill!(a, zero(eltype(a)))
        @test all(iszero.(a))
        @test a isa StructArray
        @test Bennu.device(a) == Bennu.device(A)

        a .= b .+ b
        @test isapprox(collect(adapt(Array, a)), aos_b .+ aos_b)

        ab = StructArray((a=a, b=b))
        @test ab isa StructArray
        @test Bennu.device(ab) == Bennu.device(A)
        @test all(ab.a .== a)
        @test all(ab.b .== b)

        ab2 = similar(ab, Int)
        @test eltype(ab2) == Int
        @test ab2 isa A
        @test Bennu.device(ab2) == Bennu.device(A)

        e = b .* 3
        @test isapprox(collect(adapt(Array, e)), aos_b .* 3 )
        @test e isa StructArray
        @test Bennu.device(e) == Bennu.device(A)

        c = norm.(b)
        @test isapprox(collect(adapt(Array, c)), norm.(aos_b))
        if T ≠ BigFloat
            @test c isa A
        end

        for s in ((), (1,), (1,3), (1,3,2))
            d = similar(b, SVector{2, Float32}, s)
            @test size(d) == s
            @test d isa StructArray
            @test Bennu.device(d) == Bennu.device(A)
        end

        a = similar(b)
        fill!(a, zero(eltype(a)))
        event = Event(Bennu.device(A))
        event = dup_kernel!(Bennu.device(A), 256)(a, b, ndrange=length(a),
                                                  dependencies = (event, ))
        wait(event)
        @test isapprox(collect(adapt(Array, a)), aos_b .* 2 )

        fill!(a, zero(eltype(a)))
        fill!(b, zero(eltype(b)))
        event = Event(Bennu.device(A))
        event = split_kernel!(Bennu.device(A), 256)(a, b, ab, ndrange=length(a),
                                                    dependencies = (event, ))
        wait(event)
        @test all(ab.a .== a)
        @test all(ab.b .== b)

        fill!(ab, (a=zero(eltype(a)), b=zero(eltype(b))))
        event = Event(Bennu.device(A))
        event = join_kernel!(Bennu.device(A), 256)(ab, a, b, ndrange=length(a),
                                                   dependencies = (event, ))
        wait(event)
        @test all(ab.a .== a)
        @test all(ab.b .== b)
    end
end


text:
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
#
# Description
#
#   Types and structures for the html backend.
#
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #

export HTMLDecoration, HTMLHighlighter, HTMLTableFormat

"""
    struct HTMLDecoration

Structure that defines parameters to decorate a table cell.

"""
@with_kw struct HTMLDecoration
    color::String           = ""
    background::String      = ""
    font_family::String     = ""
    font_weight::String     = ""
    font_decoration::String = ""
end

HTMLDecoration(color::String) = HTMLDecoration(color = color)

function Dict(d::HTMLDecoration)
    style = Dict{String,String}()

    !isempty(d.color)           && (style["color"]           = d.color)
    !isempty(d.background)      && (style["background"]      = d.background)
    !isempty(d.font_decoration) && (style["font-decoration"] = d.font_decoration)
    !isempty(d.font_family)     && (style["font-family"]     = d.font_family)
    !isempty(d.font_weight)     && (style["font-weight"]     = d.font_weight)

    return style
end

"""
    struct HTMLTableFormat

Format that will be used to print the HTML table. All parameters are strings
compatible with the corresponding HTML property.

# Fields

* `css`: CSS to be injected at the end of the `<style>` section.
* `header_colors`: A tuple with the foreground and background colors of the
                   header.
* `subheader_colors`: A tuple with the foreground and background colors of the
                      sub-headers
* `table_width`: Table width.

"""
@with_kw struct HTMLTableFormat
    css::String = """
    table, td, th {
        border-collapse: collapse;
        font-family: sans-serif;
    }

    td, th {
        border-bottom: 0;
        padding: 4px
    }

    tr:nth-child(odd) {
        background: #eee;
    }

    tr:nth-child(even) {
        background: #fff;
    }

    tr.header {
        background: navy;
        color: white;
        font-weight: bold;
    }

    tr.subheader, tr.subheaderLastRow {
        background: lightgray;
        color: black;
    }

    th.rowNumber, td.rowNumber {
        text-align: right;
    }
    """
    table_width::String = ""
end

#                                 Highlighters
# ==============================================================================

# Abstract type for all highlighters.
"""
    abstract type AbstractHTMLHighlighter

All HTLM highlighters must be a subtype of `AbstractHTMLHighlighter`. They API
dictates that they must implement two functions:

* `f`: Function with the signature `f(data,i,j)` in which should return `true`
       if the element `(i,j)` in `data` must be highlighter, or `false`
       otherwise.
* `fd`: Function with the signature `f(h,data,i,j)` in which `h` is the
        highlighter. This function must return the `HTMLDecoration` to be
        applied to the cell that must be highlighted.

"""
abstract type AbstractHTMLHighlighter end

"""
    struct HTMLHighlighter

Defines the default highlighter of a table when using the html backend.

# Fields

* `f`: Function with the signature `f(data,i,j)` in which should return `true`
       if the element `(i,j)` in `data` must be highlighter, or `false`
       otherwise.
* `fd`: Function with the signature `f(h,data,i,j)` in which `h` is the
        highlighter. This function must return the `HTMLDecoration` to be
        applied to the cell that must be highlighted.
* `decoration`: The `HTMLDecoration` to be applied to the highlighted cell if
                the default `fd` is used.

# Remarks

This structure can be constructed using two helpers:

    HTMLHighlighter(f::Function, decoration::HTMLDecoration)

    HTMLHighlighter(f::Function, fd::Function)

The first will apply a fixed decoration to the highlighted cell specified in
`decoration` whereas the second let the user select the desired decoration by
specifying the function `fd`.

"""
@with_kw struct HTMLHighlighter <: AbstractHTMLHighlighter
    # API
    f::Function
    fd::Function = (h,data,i,j)->h.decoration

    # Private
    decoration::HTMLDecoration = HTMLDecoration()
end

# Helper function to construct HTMLHighlighter.
HTMLHighlighter(f::Function, decoration::HTMLDecoration) =
    HTMLHighlighter(f = f, decoration = decoration)

HTMLHighlighter(f::Function, fd::Function) =
    HTMLHighlighter(f = f, fd = fd)

ext: jl


text:
# This file is a part of Julia. License is MIT: https://julialang.org/license

mutable struct SlotInfo
    defs::Vector{Int}
    uses::Vector{Int}
    any_newvar::Bool
end
SlotInfo() = SlotInfo(Int[], Int[], false)

function scan_entry!(result::Vector{SlotInfo}, idx::Int, @nospecialize(stmt))
    # NewVarNodes count as defs for the purpose
    # of liveness analysis (i.e. they kill use chains)
    if isa(stmt, NewvarNode)
        result[slot_id(stmt.slot)].any_newvar = true
        push!(result[slot_id(stmt.slot)].defs, idx)
        return
    elseif isexpr(stmt, :(=))
        if isa(stmt.args[1], SlotNumber)
            push!(result[slot_id(stmt.args[1])].defs, idx)
        end
        stmt = stmt.args[2]
    end
    if isa(stmt, Union{SlotNumber, TypedSlot})
        push!(result[slot_id(stmt)].uses, idx)
        return
    end
    for op in userefs(stmt)
        val = op[]
        if isa(val, Union{SlotNumber, TypedSlot})
            push!(result[slot_id(val)].uses, idx)
        end
    end
end

function lift_defuse(cfg::CFG, defuse)
    map(defuse) do slot
        SlotInfo(
            Int[block_for_inst(cfg, x) for x in slot.defs],
            Int[block_for_inst(cfg, x) for x in slot.uses],
            slot.any_newvar
        )
    end
end

@inline slot_id(s) = isa(s, SlotNumber) ? (s::SlotNumber).id : (s::TypedSlot).id
function scan_slot_def_use(nargs::Int, ci::CodeInfo, code::Vector{Any})
    nslots = length(ci.slotnames)
    result = SlotInfo[SlotInfo() for i = 1:nslots]
    # Set defs for arguments
    for var in result[1:(1+nargs)]
        push!(var.defs, 0)
    end
    for (idx, stmt) in Iterators.enumerate(code)
        scan_entry!(result, idx, stmt)
    end
    result
end

function renumber_ssa(stmt::SSAValue, ssanums::Vector{Any}, new_ssa::Bool=false, used_ssa::Union{Nothing, Vector{Int}}=nothing)
    id = stmt.id
    if id > length(ssanums)
        return stmt
    end
    val = ssanums[id]
    if isa(val, SSAValue) && used_ssa !== nothing
        used_ssa[val.id] += 1
    end
    return val
end

function renumber_ssa!(@nospecialize(stmt), ssanums::Vector{Any}, new_ssa::Bool=false, used_ssa::Union{Nothing, Vector{Int}}=nothing)
    isa(stmt, SSAValue) && return renumber_ssa(stmt, ssanums, new_ssa, used_ssa)
    return ssamap(val->renumber_ssa(val, ssanums, new_ssa, used_ssa), stmt)
end

function make_ssa!(ci::CodeInfo, code::Vector{Any}, idx, slot, @nospecialize(typ))
    (idx == 0) && return Argument(slot)
    stmt = code[idx]
    @assert isexpr(stmt, :(=))
    code[idx] = stmt.args[2]
    ci.ssavaluetypes[idx] = typ
    idx
end

struct UndefToken
end
const undef_token = UndefToken()

function new_to_regular(@nospecialize(stmt), new_offset::Int)
    if isa(stmt, NewSSAValue)
        return SSAValue(stmt.id + new_offset)
    end
    urs = userefs(stmt)
    for op in urs
        val = op[]
        if isa(val, NewSSAValue)
            op[] = SSAValue(val.id + new_offset)
        end
    end
    return urs[]
end

function fixup_slot!(ir::IRCode, ci::CodeInfo, idx::Int, slot::Int, @nospecialize(stmt::Union{SlotNumber, TypedSlot}), @nospecialize(ssa))
    # We don't really have the information here to get rid of these.
    # We'll do so later
    if ssa === undef_token
        insert_node!(ir, idx, Any, Expr(:throw_undef_if_not, ci.slotnames[slot], false))
        return undef_token
    end
    if !isa(ssa, Argument) && !(ssa === nothing) && ((ci.slotflags[slot] & SLOT_USEDUNDEF) != 0)
        insert_node!(ir, idx, Any, Expr(:undefcheck, ci.slotnames[slot], ssa))
    end
    if isa(stmt, SlotNumber)
        return ssa
    elseif isa(stmt, TypedSlot)
        return NewSSAValue(insert_node!(ir, idx, stmt.typ, PiNode(ssa, stmt.typ)).id - length(ir.stmts))
    end
end

function fixemup!(cond, rename, ir::IRCode, ci::CodeInfo, idx::Int, @nospecialize(stmt))
    if isa(stmt, Union{SlotNumber, TypedSlot}) && cond(stmt)
        return fixup_slot!(ir, ci, idx, slot_id(stmt), stmt, rename(stmt))
    end
    if isexpr(stmt, :(=))
        stmt.args[2] = fixemup!(cond, rename, ir, ci, idx, stmt.args[2])
        return stmt
    end
    if isa(stmt, PhiNode)
        for i = 1:length(stmt.edges)
            isassigned(stmt.values, i) || continue
            val = stmt.values[i]
            isa(val, Union{SlotNumber, TypedSlot}) || continue
            cond(val) || continue
            bb_idx = block_for_inst(ir.cfg, stmt.edges[i])
            from_bb_terminator = last(ir.cfg.blocks[bb_idx].stmts)
            stmt.values[i] = fixup_slot!(ir, ci, from_bb_terminator, slot_id(val), val, rename(val))
        end
        return stmt
    end
    if isexpr(stmt, :isdefined)
        val = stmt.args[1]
        if isa(val, Union{SlotNumber, TypedSlot})
            slot = slot_id(val)
            if (ci.slotflags[slot] & SLOT_USEDUNDEF) == 0
                return true
            else
                ssa = rename(val)
                if ssa === undef_token
                    return false
                elseif !isa(ssa, SSAValue) && !isa(ssa, NewSSAValue)
                    return true
                end
            end
            stmt.args[1] = ssa
        end
        return stmt
    end
    urs = userefs(stmt)
    for op in urs
        val = op[]
        if isa(val, Union{SlotNumber, TypedSlot}) && cond(val)
            x = fixup_slot!(ir, ci, idx, slot_id(val), val, rename(val))
            # We inserted an undef error node. Delete subsequent statement
            # to avoid confusing the optimizer
            if x === undef_token
                return nothing
            end
            op[] = x
        end
    end
    return urs[]
end

function fixup_uses!(ir::IRCode, ci::CodeInfo, code, uses::Vector{Int}, slot, @nospecialize(ssa))
    for use in uses
        code[use] = fixemup!(stmt->slot_id(stmt)==slot, stmt->ssa, ir, ci, use, code[use])
    end
end

function rename_uses!(ir::IRCode, ci::CodeInfo, idx::Int, @nospecialize(stmt), renames::Vector{Any})
    return fixemup!(stmt->true, stmt->renames[slot_id(stmt)], ir, ci, idx, stmt)
end

function strip_trailing_junk!(ci::CodeInfo, code::Vector{Any}, flags::Vector{UInt8})
    # Remove `nothing`s at the end, we don't handle them well
    # (we expect the last instruction to be a terminator)
    for i = length(code):-1:1
        if code[i] !== nothing
            resize!(code, i)
            resize!(ci.ssavaluetypes, i)
            resize!(ci.codelocs, i)
            resize!(flags, i)
            break
        end
    end
    # If the last instruction is not a terminator, add one. This can
    # happen for implicit return on dead branches.
    term = code[end]
    if !isa(term, GotoIfNot) && !isa(term, GotoNode) && !isa(term, ReturnNode)
        push!(code, ReturnNode())
        push!(ci.ssavaluetypes, Union{})
        push!(ci.codelocs, 0)
        push!(flags, 0x00)
    end
    nothing
end

struct DelayedTyp
    phi::NewSSAValue
end

# maybe use expr_type?
function typ_for_val(@nospecialize(x), ci::CodeInfo, spvals::SimpleVector, idx::Int, slottypes::Vector{Any})
    if isa(x, Expr)
        if x.head === :static_parameter
            return sparam_type(spvals[x.args[1]])
        elseif x.head === :boundscheck
            return Bool
        elseif x.head === :copyast
            return typ_for_val(x.args[1], ci, spvals, idx, slottypes)
        end
        return ci.ssavaluetypes[idx]
    end
    isa(x, GlobalRef) && return abstract_eval_global(x.mod, x.name)
    isa(x, SSAValue) && return ci.ssavaluetypes[x.id]
    isa(x, Argument) && return slottypes[x.n]
    isa(x, NewSSAValue) && return DelayedTyp(x)
    isa(x, QuoteNode) && return Const(x.value)
    isa(x, Union{Symbol, PiNode, PhiNode, SlotNumber, TypedSlot}) && error("unexpected val type")
    return Const(x)
end

struct BlockLiveness
    def_bbs::Vector{Int}
    live_in_bbs::Vector{Int}
end

# Run iterated dominance frontier
#
# The algorithm we have here essentially follows LLVM, which itself is a
# a cleaned up version of the linear-time algorithm described in
#
#  A Linear Time Algorithm for Placing phi-Nodes (by Sreedhar and Gao)
#
# The algorithm here, is quite straightforward. Suppose we have a CFG:
#
# A -> B -> D -> F
#  \-> C -------/
#
# and a corresponding dominator tree:
#
# A
# |- B - D
# |- C
# |- F
#
# Now, for every definition of our slot, we simply walk down the dominator
# tree and look for any edges that leave the sub-domtree rooted by our definition.
#
# E.g. in our example above, if we have a definition in `B`, we look at its successors,
#      which is only `D`, which is dominated by `B` and hence doesn't need a phi node.
#      Then we descend down the subtree rooted at `B` and end up in `D`. `D` has a successor
#      `F`, which is not part of the current subtree, (i.e. not dominated by `B`), so it
#      needs a phi node.
#
# Now, the key insight of that algorithm is that we have two defs, in blocks `A` and `B`,
# and `A` dominates `B`, then we do not need to recurse into `B`, because the set of
# potential backedges from a subtree rooted at `B` (to outside the subtree) is a strict
# subset of those backedges from a subtree rooted at `A` (out outside the subtree rooted
# at `A`). Note however that this does not work the other way. Thus, the algorithm
# needs to make sure that we always visit `B` before `A`.
function idf(cfg::CFG, liveness::BlockLiveness, domtree::DomTree)
    # This should be a priority queue, but TODO - sorted array for now
    defs = liveness.def_bbs
    pq = Tuple{Int, Int}[(defs[i], domtree.nodes[defs[i]].level) for i in 1:length(defs)]
    sort!(pq, by=x->x[2])
    phiblocks = Int[]
    # This bitset makes sure we only add a phi node to a given block once.
    processed = BitSet()
    # This bitset implements the `key insight` mentioned above. In particular, it prevents
    # us from visiting a subtree that we have already visited before.
    visited = BitSet()
    while !isempty(pq)
        # We pop from the end of the array - i.e. the element with the highest level.
        node, level = pop!(pq)
        worklist = Int[]
        push!(worklist, node)
        while !isempty(worklist)
            active = pop!(worklist)
            for succ in cfg.blocks[active].succs
                # Check whether the current root (`node`) dominates succ.
                # We are guaranteed that `node` dominates `active`, since
                # we've arrived at `active` by following dominator tree edges.
                # If `succ`'s level is less than or equal to that of `node`,
                # it cannot possibly be dominated by `node`. On the other hand,
                # since at this point we know that there is an edge from `node`'s
                # subtree to `succ`, we know that if succ's level is greater than
                # that of `node`, it must be dominated by `node`.
                succ_level = domtree.nodes[succ].level
                succ_level > level && continue
                # We don't dominate succ. We need to place a phinode,
                # unless liveness said otherwise.
                succ in processed && continue
                push!(processed, succ)
                if !(succ in liveness.live_in_bbs)
                    continue
                end
                push!(phiblocks, succ)
                # Basically: Consider the phi node we just added another
                # def of this value. N.B.: This needs to retain the invariant that it
                # is processed before any of its parents in the dom tree. This is guaranteed,
                # because succ_level <= level, which is the greatest level we have currently
                # processed. Thus, we have not yet processed any subtrees of level < succ_level.
                if !(succ in defs)
                    push!(pq, (succ, succ_level))
                    sort!(pq, by=x->x[2])
                end
            end
            # Recurse down the current subtree
            for child in domtree.nodes[active].children
                child in visited && continue
                push!(visited, child)
                push!(worklist, child)
            end
        end
    end
    phiblocks
end

function rename_incoming_edge(old_edge, old_to, result_order, bb_rename)
    new_edge_from = bb_rename[old_edge]
    if old_edge == old_to - 1
        # Could have been a crit edge break
        if new_edge_from < length(result_order) && result_order[new_edge_from + 1] == 0
            new_edge_from += 1
        end
    end
    new_edge_from
end

function rename_outgoing_edge(old_to, old_from, result_order, bb_rename)
    new_edge_to = bb_rename[old_to]
    if old_from == old_to - 1
        # Could have been a crit edge break
        if bb_rename[old_from] < length(result_order) && result_order[bb_rename[old_from]+1] == 0
            new_edge_to = bb_rename[old_from] + 1
        end
    end
    new_edge_to
end

function rename_phinode_edges(node, bb, result_order, bb_rename)
    new_values = Any[]
    new_edges = Any[]
    for (idx, edge) in pairs(node.edges)
        (edge == 0 || haskey(bb_rename, edge)) || continue
        new_edge_from = edge == 0 ? 0 : rename_incoming_edge(edge, bb, result_order, bb_rename)
        push!(new_edges, new_edge_from)
        if isassigned(node.values, idx)
            push!(new_values, node.values[idx])
        else
            resize!(new_values, length(new_values)+1)
        end
    end
    if length(new_edges) == 1
        return isassigned(new_values, 1) ? new_values[1] : PhiNode(Any[], Any[])
    else
        return PhiNode(new_edges, new_values)
    end
end

"""
    Sort the basic blocks in `ir` into domtree order (i.e. if bb`` is higher in
    the domtree than bb2, it will come first in the linear order). The resulting
    ir has the property that a linear traversal of basic blocks will also be a
    RPO traversal and in particular, any use of an SSA value must come after (by linear
    order) its definition.
"""
function domsort_ssa!(ir::IRCode, domtree::DomTree)
    # First compute the new order of basic blocks
    result_order = Int[]
    stack = Int[]
    node = 1
    ncritbreaks = 0
    nnewfallthroughs = 0
    while node !== -1
        push!(result_order, node)
        cs = domtree.nodes[node].children
        terminator = ir.stmts[last(ir.cfg.blocks[node].stmts)]
        iscondbr = isa(terminator, GotoIfNot)
        let old_node = node + 1
            if length(cs) >= 1
                # Adding the nodes in reverse sorted order attempts to retain
                # the original source order of the nodes as much as possible.
                # This is not required for correctness, but is easier on the humans
                if old_node in cs
                    # Schedule the fall through node first,
                    # so we can retain the fall through
                    append!(stack, reverse(sort(filter(x -> (x != old_node), cs))))
                    node = node + 1
                else
                    append!(stack, reverse(sort(cs)))
                    node = pop!(stack)
                end
            else
                if isempty(stack)
                    node = -1
                else
                    node = pop!(stack)
                end
            end
            if node != old_node && !isa(terminator, Union{GotoNode, ReturnNode})
                if isa(terminator, GotoIfNot)
                    # Need to break the critical edge
                    ncritbreaks += 1
                    push!(result_order, 0)
                else
                    nnewfallthroughs += 1
                end
            end
        end
    end
    bb_rename = IdDict{Int,Int}(i=>x for (x, i) in pairs(result_order) if i !== 0)
    new_bbs = Vector{BasicBlock}(undef, length(result_order))
    nstmts = 0
    for i in result_order
        if i !== 0
            nstmts += length(ir.cfg.blocks[i].stmts)
        end
    end
    result_stmts = Vector{Any}(undef, nstmts + ncritbreaks + nnewfallthroughs)
    result_types = Any[Any for i = 1:length(result_stmts)]
    result_ltable = fill(Int32(0), length(result_stmts))
    result_flags = fill(0x00, length(result_stmts))
    inst_rename = Vector{Any}(undef, length(ir.stmts))
    for i = 1:length(ir.new_nodes)
        push!(inst_rename, SSAValue(nstmts + i + ncritbreaks + nnewfallthroughs))
    end
    bb_start_off = 0
    crit_edge_breaks_fixup = Tuple{Int, Int}[]
    for (new_bb, bb) in pairs(result_order)
        if bb == 0
            @assert isa(result_stmts[bb_start_off+1], GotoNode)
            # N.B.: The .label has already been renamed when it was created.
            new_bbs[new_bb] = BasicBlock((bb_start_off+1):(bb_start_off+1), [new_bb-1], [result_stmts[bb_start_off+1].label])
            bb_start_off += 1
            continue
        end
        old_inst_range = ir.cfg.blocks[bb].stmts
        inst_range = (bb_start_off+1):(bb_start_off+length(old_inst_range))
        for (nidx, idx) in zip(inst_range, old_inst_range)
            inst_rename[idx] = SSAValue(nidx)
            stmt = ir.stmts[idx]
            if isa(stmt, PhiNode)
                result_stmts[nidx] = rename_phinode_edges(stmt, bb, result_order, bb_rename)
            else
                result_stmts[nidx] = stmt
            end
            result_types[nidx] = ir.types[idx]
            result_ltable[nidx] = ir.lines[idx]
            result_flags[nidx] = ir.flags[idx]
        end
        # Now fix up the terminator
        terminator = result_stmts[inst_range[end]]
        if isa(terminator, GotoNode)
            # Convert to implicit fall through
            if bb_rename[terminator.label] == new_bb + 1
                result_stmts[inst_range[end]] = nothing
            else
                result_stmts[inst_range[end]] = GotoNode(bb_rename[terminator.label])
            end
        elseif isa(terminator, GotoIfNot)
            # Check if we need to break the critical edge
            if bb_rename[bb + 1] != new_bb + 1
                @assert result_order[new_bb + 1] == 0
                # Add an explicit goto node in the next basic block (we accounted for this above)
                result_stmts[inst_range[end]+1] = GotoNode(bb_rename[bb+1])
            end
            result_stmts[inst_range[end]] = GotoIfNot(terminator.cond, bb_rename[terminator.dest])
        elseif !isa(terminator, ReturnNode)
            if isa(terminator, Expr) && terminator.head == :enter
                terminator.args[1] = bb_rename[terminator.args[1]]
            end
            if bb_rename[bb + 1] != new_bb + 1
                # Add an explicit goto node
                result_stmts[inst_range[end]+1] = GotoNode(bb_rename[bb+1])
                inst_range = first(inst_range):(last(inst_range)+1)
            end
        end
        bb_start_off += length(inst_range)
        local new_preds, new_succs
        let bb = bb, bb_rename = bb_rename, result_order = result_order
            new_preds = Int[rename_incoming_edge(i, bb, result_order, bb_rename) for i in ir.cfg.blocks[bb].preds if haskey(bb_rename, i)]
            new_succs = Int[rename_outgoing_edge(i, bb, result_order, bb_rename) for i in ir.cfg.blocks[bb].succs if haskey(bb_rename, i)]
        end
        new_bbs[new_bb] = BasicBlock(inst_range, new_preds, new_succs)
    end
    result_stmts = Any[renumber_ssa!(result_stmts[i], inst_rename, true) for i in 1:length(result_stmts)]
    cfg = CFG(new_bbs, Int[first(bb.stmts) for bb in new_bbs[2:end]])
    new_new_nodes = Vector{NewNode}(undef, length(ir.new_nodes))
    for i = 1:length(ir.new_nodes)
        entry = ir.new_nodes[i]
        new_new_nodes[i] = NewNode(inst_rename[entry.pos].id, entry.attach_after, entry.typ,
            renumber_ssa!(isa(entry.node, PhiNode) ?
                rename_phinode_edges(entry.node, block_for_inst(ir.cfg, entry.pos), result_order, bb_rename) : entry.node,
                inst_rename, true),
            entry.line)
    end
    new_ir = IRCode(ir, result_stmts, result_types, result_ltable, result_flags, cfg, new_new_nodes)
    return new_ir
end

function compute_live_ins(cfg::CFG, defuse)
    # We remove from `uses` any block where all uses are dominated
    # by a def. This prevents insertion of dead phi nodes at the top
    # of such a block if that block happens to be in a loop
    ordered = Tuple{Int, Int, Bool}[(x, block_for_inst(cfg, x), true) for x in defuse.uses]
    for x in defuse.defs
        push!(ordered, (x, block_for_inst(cfg, x), false))
    end
    ordered = sort(ordered, by=x->x[1])
    bb_defs = Int[]
    bb_uses = Int[]
    last_bb = last_def_bb = 0
    for (_, bb, is_use) in ordered
        if bb != last_bb && is_use
            push!(bb_uses, bb)
        end
        last_bb = bb
        if last_def_bb != bb && !is_use
            push!(bb_defs, bb)
            last_def_bb = bb
        end
    end
    # To obtain live ins from bb_uses, recursively add predecessors
    extra_liveins = BitSet()
    worklist = Int[]
    for bb in bb_uses
        append!(worklist, filter(p->p != 0 && !(p in bb_defs), cfg.blocks[bb].preds))
    end
    while !isempty(worklist)
        elem = pop!(worklist)
        (elem in bb_uses || elem in extra_liveins) && continue
        push!(extra_liveins, elem)
        append!(worklist, filter(p->p != 0 && !(p in bb_defs), cfg.blocks[elem].preds))
    end
    append!(bb_uses, extra_liveins)
    BlockLiveness(bb_defs, bb_uses)
end

function recompute_type(node::Union{PhiNode, PhiCNode}, ci::CodeInfo, ir::IRCode, spvals::SimpleVector, slottypes::Vector{Any})
    new_typ = Union{}
    for i = 1:length(node.values)
        if isa(node, PhiNode) && !isassigned(node.values, i)
            if !isa(new_typ, MaybeUndef)
                new_typ = MaybeUndef(new_typ)
            end
            continue
        end
        typ = typ_for_val(node.values[i], ci, spvals, -1, slottypes)
        was_maybe_undef = false
        if isa(typ, MaybeUndef)
            typ = typ.typ
            was_maybe_undef = true
        end
        @assert !isa(typ, MaybeUndef)
        while isa(typ, DelayedTyp)
            typ = types(ir)[typ.phi::NewSSAValue]
        end
        new_typ = tmerge(new_typ, was_maybe_undef ? MaybeUndef(typ) : typ)
    end
    return new_typ
end

function construct_ssa!(ci::CodeInfo, code::Vector{Any}, ir::IRCode, domtree::DomTree, defuse, nargs::Int, spvals::SimpleVector,
                        slottypes::Vector{Any})
    cfg = ir.cfg
    left = Int[]
    catch_entry_blocks = Tuple{Int, Int}[]
    for (idx, stmt) in pairs(code)
        if isexpr(stmt, :enter)
            push!(catch_entry_blocks, (block_for_inst(cfg, idx), block_for_inst(cfg, stmt.args[1])))
        end
    end

    exc_handlers = IdDict{Int, Int}()
    # Record the correct exception handler for all cricitcal sections
    for (enter_block, exc) in catch_entry_blocks
        exc_handlers[enter_block+1] = exc
        # TODO: Cut off here if the terminator is a leave corresponding to this enter
        for block in dominated(domtree, enter_block+1)
            exc_handlers[block] = exc
        end
    end

    phi_slots = Vector{Int}[Vector{Int}() for _ = 1:length(ir.cfg.blocks)]
    phi_nodes = Vector{Pair{NewSSAValue,PhiNode}}[Vector{Pair{NewSSAValue,PhiNode}}() for _ = 1:length(cfg.blocks)]
    phi_ssas = SSAValue[]
    phicnodes = IdDict{Int, Vector{Tuple{SlotNumber, NewSSAValue, PhiCNode}}}()
    for (_, exc) in catch_entry_blocks
        phicnodes[exc] = Vector{Tuple{SlotNumber, NewSSAValue, PhiCNode}}()
    end
    @timeit "idf" for (idx, slot) in Iterators.enumerate(defuse)
        # No uses => no need for phi nodes
        isempty(slot.uses) && continue
        # TODO: Restore this optimization
        if false # length(slot.defs) == 1 && slot.any_newvar
            if slot.defs[] == 0
                typ = slottypes[idx]
                ssaval = Argument(idx)
                fixup_uses!(ir, ci, code, slot.uses, idx, ssaval)
            elseif isa(code[slot.defs[]], NewvarNode)
                typ = MaybeUndef(Union{})
                ssaval = nothing
                for use in slot.uses[]
                    insert_node!(ir, use, Union{}, Expr(:throw_undef_if_not, ci.slotnames[idx], false))
                end
                fixup_uses!(ir, ci, code, slot.uses, idx, nothing)
            else
                val = code[slot.defs[]].args[2]
                typ = typ_for_val(val, ci, spvals, slot.defs[], slottypes)
                ssaval = SSAValue(make_ssa!(ci, code, slot.defs[], idx, typ))
                fixup_uses!(ir, ci, code, slot.uses, idx, ssaval)
            end
            continue
        end
        @timeit "liveness" (live = compute_live_ins(cfg, slot))
        for li in live.live_in_bbs
            cidx = findfirst(x->x[2] == li, catch_entry_blocks)
            if cidx !== nothing
                # The slot is live-in into this block. We need to
                # Create a PhiC node in the catch entry block and
                # an upsilon node in the corresponding enter block
                node = PhiCNode(Any[])
                phic_ssa = NewSSAValue(insert_node!(ir, first_insert_for_bb(code, cfg, li), Union{}, node).id - length(ir.stmts))
                push!(phicnodes[li], (SlotNumber(idx), phic_ssa, node))
                # Inform IDF that we now have a def in the catch block
                if !(li in live.def_bbs)
                    push!(live.def_bbs, li)
                end
            end
        end
        phiblocks = idf(cfg, live, domtree)
        for block in phiblocks
            push!(phi_slots[block], idx)
            node = PhiNode()
            ssa = NewSSAValue(insert_node!(ir, first_insert_for_bb(code, cfg, block), Union{}, node).id - length(ir.stmts))
            push!(phi_nodes[block], ssa=>node)
        end
        push!(left, idx)
    end
    # Perform SSA renaming
    initial_incoming_vals = Any[
        if 0 in defuse[x].defs
            Argument(x)
        elseif !defuse[x].any_newvar
            undef_token
        else
            SSAValue(-1)
        end for x in 1:length(ci.slotnames)
    ]
    worklist = Tuple{Int, Int, Vector{Any}}[(1, 0, initial_incoming_vals)]
    visited = BitSet()
    type_refine_phi = BitSet()
    @timeit "SSA Rename" while !isempty(worklist)
        (item, pred, incoming_vals) = pop!(worklist)
        # Rename existing phi nodes first, because their uses occur on the edge
        # TODO: This isn't necessary if inlining stops replacing arguments by slots.
        for idx in cfg.blocks[item].stmts
            stmt = code[idx]
            if isexpr(stmt, :(=))
                stmt = stmt.args[2]
            end
            isa(stmt, PhiNode) || continue
            for (edgeidx, edge) in pairs(stmt.edges)
                from_bb = edge == 0 ? 0 : block_for_inst(cfg, edge)
                from_bb == pred || continue
                isassigned(stmt.values, edgeidx) || break
                stmt.values[edgeidx] = rename_uses!(ir, ci, edge, stmt.values[edgeidx], incoming_vals)
                break
            end
        end
        # Insert phi nodes if necessary
        for (idx, slot) in Iterators.enumerate(phi_slots[item])
            ssaval, node = phi_nodes[item][idx]
            incoming_val = incoming_vals[slot]
            if incoming_val == SSAValue(-1)
                # Optimistically omit this path.
                # Liveness analysis would probably have prevented us from inserting this phi node
                continue
            end
            push!(node.edges, pred)
            if incoming_val == undef_token
                resize!(node.values, length(node.values)+1)
            else
                push!(node.values, incoming_val)
            end
            # TODO: Remove the next line, it shouldn't be necessary
            push!(type_refine_phi, ssaval.id)
            if isa(incoming_val, NewSSAValue)
                push!(type_refine_phi, ssaval.id)
            end
            typ = incoming_val == undef_token ? MaybeUndef(Union{}) : typ_for_val(incoming_val, ci, spvals, -1, slottypes)
            old_entry = ir.new_nodes[ssaval.id]
            if isa(typ, DelayedTyp)
                push!(type_refine_phi, ssaval.id)
            end
            new_typ = isa(typ, DelayedTyp) ? Union{} : tmerge(old_entry.typ, typ)
            ir.new_nodes[ssaval.id] = NewNode(old_entry.pos, old_entry.attach_after, new_typ, node, old_entry.line)
            incoming_vals[slot] = ssaval
        end
        (item in visited) && continue
        # Record phi_C nodes if necessary
        if haskey(phicnodes, item)
            for (slot, ssa, _) in phicnodes[item]
                incoming_vals[slot_id(slot)] = ssa
            end
        end
        # Record initial upsilon nodes if necessary
        eidx = findfirst(x->x[1] == item, catch_entry_blocks)
        if eidx !== nothing
            for (slot, _, node) in phicnodes[catch_entry_blocks[eidx][2]]
                ivalundef = incoming_vals[slot_id(slot)] === undef_token
                unode = ivalundef ? UpsilonNode() : UpsilonNode(incoming_vals[slot_id(slot)])
                typ = ivalundef ? MaybeUndef(Union{}) : slottypes[slot_id(slot)]
                push!(node.values,
                    NewSSAValue(insert_node!(ir, first_insert_for_bb(code, cfg, item),
                                 typ, unode, true).id - length(ir.stmts)))
            end
        end
        push!(visited, item)
        for idx in cfg.blocks[item].stmts
            stmt = code[idx]
            (isa(stmt, PhiNode) || (isexpr(stmt, :(=)) && isa(stmt.args[2], PhiNode))) && continue
            if isa(stmt, NewvarNode)
                incoming_vals[slot_id(stmt.slot)] = undef_token
                code[idx] = nothing
            else
                stmt = rename_uses!(ir, ci, idx, stmt, incoming_vals)
                if stmt === nothing && isa(code[idx], Union{ReturnNode, GotoIfNot}) && idx == last(cfg.blocks[item].stmts)
                    # preserve the CFG
                    stmt = ReturnNode()
                end
                code[idx] = stmt
                # Record a store
                if isexpr(stmt, :(=)) && isa(stmt.args[1], SlotNumber)
                    id = slot_id(stmt.args[1])
                    val = stmt.args[2]
                    typ = typ_for_val(val, ci, spvals, idx, slottypes)
                    # Having undef_token appear on the RHS is possible if we're on a dead branch.
                    # Do something reasonable here, by marking the LHS as undef as well.
                    if val !== undef_token
                        incoming_vals[id] = SSAValue(make_ssa!(ci, code, idx, id, typ))
                    else
                        code[idx] = nothing
                        incoming_vals[id] = undef_token
                    end
                    if haskey(exc_handlers, item)
                        exc = exc_handlers[item]
                        cidx = findfirst(x->slot_id(x[1]) == id, phicnodes[exc])
                        if cidx !== nothing
                            node = UpsilonNode(incoming_vals[id])
                            if incoming_vals[id] === undef_token
                                node = UpsilonNode()
                                typ = MaybeUndef(Union{})
                            end
                            push!(phicnodes[exc][cidx][3].values,
                                NewSSAValue(insert_node!(ir, idx, typ, node, true).id - length(ir.stmts)))
                        end
                    end
                end
            end
        end
        for succ in cfg.blocks[item].succs
            push!(worklist, (succ, item, copy(incoming_vals)))
        end
    end
    # Delete any instruction in unreachable blocks (except for terminators)
    for bb in setdiff(BitSet(1:length(cfg.blocks)), visited)
        for idx in cfg.blocks[bb].stmts
            if isa(code[idx], Union{GotoNode, GotoIfNot, ReturnNode})
                code[idx] = ReturnNode()
            else
                code[idx] = nothing
            end
        end
    end
    # Convert into IRCode form
    new_code = ir.stmts
    ssavalmap = Any[SSAValue(-1) for _ in 1:(length(ci.ssavaluetypes)+1)]
    result_types = Any[Any for _ in 1:length(new_code)]
    # Detect statement positions for assignments and construct array
    for (bb, idx) in bbidxiter(ir)
        stmt = code[idx]
        # Convert GotoNode/GotoIfNot/PhiNode to BB addressing
        if isa(stmt, GotoNode)
            new_code[idx] = GotoNode(block_for_inst(cfg, stmt.label))
        elseif isa(stmt, GotoIfNot)
            new_dest = block_for_inst(cfg, stmt.dest)
            if new_dest == bb+1
                # Drop this node - it's a noop
                new_code[idx] = stmt.cond
            else
                new_code[idx] = GotoIfNot(stmt.cond, new_dest)
            end
        elseif isexpr(stmt, :enter)
            new_code[idx] = Expr(:enter, block_for_inst(cfg, stmt.args[1]))
        elseif isexpr(stmt, :leave) || isexpr(stmt, :(=)) || isexpr(stmt, :return) ||
            isexpr(stmt, :meta) || isa(stmt, NewvarNode)
            new_code[idx] = stmt
        else
            ssavalmap[idx] = SSAValue(idx)
            result_types[idx] = ci.ssavaluetypes[idx]
            if isa(stmt, PhiNode)
                edges = Any[edge == 0 ? 0 : block_for_inst(cfg, edge) for edge in stmt.edges]
                new_code[idx] = PhiNode(edges, stmt.values)
            else
                new_code[idx] = stmt
            end
        end
    end
    for (_, nodes) in phicnodes
        for (_, ssa, node) in nodes
            new_typ = Union{}
            # TODO: This could just be the ones that depend on other phis
            push!(type_refine_phi, ssa.id)
            new_idx = ssa.id
            node = ir.new_nodes[new_idx]
            for i = 1:length(node.node.values)
                orig_typ = typ = typ_for_val(node.node.values[i], ci, spvals, -1, slottypes)
                @assert !isa(typ, MaybeUndef)
                while isa(typ, DelayedTyp)
                    typ = types(ir)[typ.phi::NewSSAValue]
                end
                new_typ = tmerge(new_typ, typ)
            end
            ir.new_nodes[new_idx] = NewNode(node.pos, node.attach_after, new_typ, node.node, node.line)
        end
    end
    # This is a bit awkward, because it basically duplicates what type
    # inference does. Ideally, we'd just use this representation earlier
    # to make sure phi nodes have accurate types
    changed = true
    while changed
        changed = false
        for new_idx in type_refine_phi
            node = ir.new_nodes[new_idx]
            new_typ = recompute_type(node.node, ci, ir, spvals, slottypes)
            if !(node.typ ⊑ new_typ) || !(new_typ ⊑ node.typ)
                ir.new_nodes[new_idx] = NewNode(node.pos, node.attach_after, new_typ, node.node, node.line)
                changed = true
            end
        end
    end
    result_types = Any[isa(result_types[i], DelayedTyp) ? types(ir)[result_types[i].phi::NewSSAValue] : result_types[i] for i in 1:length(result_types)]
    new_nodes = NewNode[let node = ir.new_nodes[i]
            typ = isa(node.typ, DelayedTyp) ? types(ir)[node.typ.phi::NewSSAValue] : node.typ
            NewNode(node.pos, node.attach_after, typ, node.node, node.line)
        end for i in 1:length(ir.new_nodes)]
    # Renumber SSA values
    new_code = Any[new_to_regular(renumber_ssa!(new_code[i], ssavalmap), length(ir.stmts)) for i in 1:length(new_code)]
    new_nodes = NewNode[let node = new_nodes[i]
            NewNode(node.pos, node.attach_after, node.typ,
            new_to_regular(renumber_ssa!(node.node, ssavalmap), length(ir.stmts)),
            node.line)
        end for i in 1:length(new_nodes)]
    ir = IRCode(ir, new_code, result_types, ir.lines, ir.flags, ir.cfg, new_nodes)
    @timeit "domsort" ir = domsort_ssa!(ir, domtree)
    return ir
end

max_stars_repo_name: asbisen/julia

max_issues_repo_name: asbisen/julia

max_forks_repo_name: asbisen/julia


text:
"""
Original publication:
A. Ghasemi, H. R. Rabiee, M. Fadaee, M. T. Manzuri,
and M. H. Rohban. Active learning from positive and
unlabeled data. In 2011 IEEE 11th International
Conference on Data Mining Workshops, pages
244–250, Dec 2011.
"""
struct ExpectedMaximumEntropyPQs <: DataBasedPQs
    p_x::Array{Float64, 1}
    bw_method::Union{String, U} where U <: Real
    ExpectedMaximumEntropyPQs(x::Array{T, 2}, bw_method="scott") where T <: Real = new(multi_kde(x, bw_method)(x), bw_method)
end

function qs_score(qs::ExpectedMaximumEntropyPQs, x::Array{T, 2}, labels::Dict{Symbol, Array{Int, 1}})::Array{Float64, 1} where T <: Real
    haskey(labels, :Lin) || throw(MissingLabelTypeException(:Lin))
    p_x_inlier = multi_kde(x[:, labels[:Lin]], qs.bw_method)(x)
    a = p_x_inlier ./ qs.p_x
    valid_score(a) = (-a^2 * log(a) + a + (a - 1)^2 * log(1 - a)) / (2 * a)
    # fallback for densities with values > 1
    return [(1 - s) <= 0 ? -Inf : valid_score(s) for s in a]
end


text:
### A Pluto.jl notebook ###
# v0.16.0

using Markdown
using InteractiveUtils

# ╔═╡ 7a3de02b-c423-496e-bf92-9981d71e5eab
begin
    using Arrow
	using Chain
	using DataFrames
	using MixedModels
	using StatsBase
	using StatsModels
	using StatsModels: ContrastsCoding
end

# ╔═╡ 09a33b05-568e-427a-bd66-812d271d1791
md"""
# Contrast Coding of Visual Attention Effects
## Date: 2021-05-29 (revised: 2021-09-16)
"""

# ╔═╡ 28c37807-b83f-4c3e-848d-a7401b4eda29
md"""
## Example data
We take the `KWDYZ` dataset ([Kliegl et al., 2011; Frontiers](https://doi.org/10.3389/fpsyg.2010.00238)). This is an experiment looking at three effects of visual cueing under four different cue-target relations (CTRs). Two horizontal rectangles are displayed above and below a central fixation point or they displayed in vertical orientation to the left and right of the fixation point.  Subjects react to the onset of a small visual target occuring at one of the four ends of the two rectangles. The target is cued validly on 70% of trials by a brief flash of the corner of the rectangle at which it appears; it is cued invalidly at the three other locations 10% of the trials each.

We specify three contrasts for the four-level factor CTR that are derived from spatial, object-based, and attractor-like features of attention. They map onto sequential differences between appropriately ordered factor levels. Interestingly, a different theoretical perspective, derived from feature overlap, leads to a different set of contrasts. Can the results refute one of the theoretical perspectives?

We also have a dataset from a replication and extension of this study (Kliegl, Kuschela, & Laubrock, 2015). Both data sets are available in [R-package RePsychLing](https://github.com/dmbates/RePsychLing/tree/master/data/) (Baayen et al., 2014).

## Preprocessing
"""

# ╔═╡ ccf701bb-cf16-4339-bd93-64ffe87435bb
begin
	dat1 = DataFrame(Arrow.Table("./data/kwdyz11.arrow"))
	dat1 = select(dat1, :subj => :Subj, :tar => :CTR, :rt);

	# Descriptive statistics
	cellmeans = combine(groupby(dat1, [:CTR]),
                             :rt => mean, :rt  => std, :rt => length,
                             :rt => (x -> std(x)/sqrt(length(x))) => :rt_semean)
end

# ╔═╡ cb63bd4b-5f0e-49d1-a257-2b83c1687975
md"""## SeqDiffCoding
This contrast corresponds to `MASS::contr.sdif()` in R.
"""

# ╔═╡ 8ffba70f-fc97-47a2-b973-3c0f2047eb9b
begin
	cntr1 = Dict(
    :CTR  => SeqDiffCoding(levels=["val", "sod", "dos", "dod"]),
    :Subj => Grouping()
	);

	formula = @formula  rt ~ 1 + CTR + (1 + CTR | Subj);
	m1 = fit(MixedModel, formula, dat1, contrasts=cntr1)
end

# ╔═╡ 4e9359a3-2533-4ce8-b2fd-26411e84f59c
md"""
## HypothesisCoding
This contrast corresponds to `MASS::contr.sdif()` in R.
A general solution (not inverse of last contrast)
"""

# ╔═╡ b00e953a-b279-4ec3-8744-5518d1242d8e
begin
	cntr1b = Dict(
    :CTR => HypothesisCoding([-1  1  0  0
                               0 -1  1  0
                               0  0  1 -1],
            levels=["val", "sod",  "dos", "dod"],
            labels=["spt", "obj", "grv"])
	);

	m1b = fit(MixedModel, formula, dat1, contrasts=cntr1b)
end

# ╔═╡ 7a3749b3-d2e9-4b53-9e57-975ee428b416
md"""

Controlling the ordering of levels for contrasts:

1.  kwarg `levels=` to order the levels; the first is set as the baseline.
2.  kwarg `base=` to fix the baseline level.

The assignment of random factors such as `Subj` to `Grouping()` is only necessary when the sample size is very large and leads to an out-of-memory error; it is included only in the first example for reference.

## DummyCoding
Thi contrast corresponds to `contr.treatment()` in R
"""

# ╔═╡ 6393d21d-a3a6-4463-8981-052440ba7887
begin
	cntr2 = Dict(:CTR => DummyCoding(base= "val"));

	m2 = fit(MixedModel, formula, dat1, contrasts=cntr2)
end

# ╔═╡ d087afe8-9d2b-489e-9160-206f22486505
md"""
This contrast has the disadvantage that the intercept returns the mean of the level specified as `base`, default is the first level, not the GM.

## YchycaeitCoding

The contrasts returned by `DummyCoding` may be what you want. Can't we have them, but also the GM rather than the mean of the base level?  Yes, we can!  I call this "You can have your cake and it eat, too"-Coding (YchycaeitCoding).
"""

# ╔═╡ bb97e862-24ad-4ce5-9011-429a0e69e5ab
begin
	cntr2b = Dict(
    :CTR => HypothesisCoding([-1  1  0  0
                              -1  0  1  0
                              -1  0  0  1],
            levels=["val", "sod",  "dos", "dod"])
	);

	m2b = fit(MixedModel, formula, dat1, contrasts=cntr2b)
end

# ╔═╡ 71528268-41fe-4de5-8725-1c65c4839c65
md"""
Just relevel the factor or move the column with -1s for a different base.

## EffectsCoding - corresponds to `contr.sum()` in R
"""

# ╔═╡ e2372cb1-03d4-4404-bfef-640ad5653792
begin
	cntr3 = Dict(:CTR => EffectsCoding(base= "dod"));

	m3 = fit(MixedModel, formula, dat1, contrasts=cntr3)
end

# ╔═╡ 142c86d3-e340-494c-8b52-8f2013993a19
md"## HelmertCoding"

# ╔═╡ 0bb71abe-433c-44e3-bff4-eb26e63bee2d
begin
	cntr4 = Dict(:CTR => HelmertCoding());

	fit(MixedModel, formula, dat1, contrasts=cntr4)
end

# ╔═╡ 0cf43493-af05-4744-a154-d7209251ae8f
md"## Reverse HelmertCoding"

# ╔═╡ 2d778d8e-949c-4a7b-a774-770c6bc4bce6
begin
	cntr4b = Dict(:CTR => HelmertCoding(levels=reverse(levels(dat1.CTR))));

	fit(MixedModel, formula, dat1, contrasts=cntr4b)
end

# ╔═╡ a71103ca-a680-45ac-ba02-61a5660cc157
md"""
Helmert contrasts are othogonal.

## AnovaCoding - Anova contrasts are orthogonal.

### A(2) x B(2)

An A(2) x B(2) design can be recast as an F(4) design with the levels (A1-B1, A1-B2, A2-B1, A2-B2). The following contrast specifiction returns estimates for the main effect of A, the main effect of B, and the interaction of A and B. In a figure With A on the x-axis and the levels of B shown as two lines, the interaction tests the null hypothesis that the two lines are parallel. A positive coefficient implies overadditivity (diverging lines toward the right) and a negative coefficient underadditivity (converging lines).
"""

# ╔═╡ 00bedbae-0bc4-4e4e-85d5-91ffc8ea77ee
begin
	cntr5 = Dict(
    :CTR => HypothesisCoding([-1  -1 +1  +1          # A
                              -1  +1 -1  +1          # B
                              +1  -1 -1  +1],        # A x B
            levels=["val", "sod",  "dos", "dod"],
            labels=["A", "B", "AxB"])
	);
	m5 = fit(MixedModel, formula, dat1, contrasts=cntr5)
end

# ╔═╡ f504edce-35de-4938-a3c6-08d4fc8f2e66
md"""
It is also helpful to see the corresponding layout of the four means for the interaction of A and B (i.e., the third contrast)

```
        B1     B2
   A1   +1     -1
   A2   -1     +1
```

Thus, interaction tests whether the difference between main diagonal and minor diagonal is different from zero.
=#

### A(2) x B(2) x C(2)

Going beyond the four level factor; it is also helpful to see the corresponding layout of the eight means for the interaction of A and B and C.

```
          C1              C2
      B1     B2        B1     B2
 A1   +1     -1   A1   -1     +1
 A2   -1     +1   A2   +1     -1
```

### A(2) x B(2) x C(3)

TO BE DONE
"""

# ╔═╡ a12c0a44-f02e-479e-94f9-c5d09c3b25ae
md"""
## NestedCoding

An A(2) x B(2) design can be recast as an F(4) design with the levels (A1-B1, A1-B2, A2-B1, A2-B2).  The following contrast specifiction returns an estimate for the main effect of A and the effects of B nested in the two levels of A. In a figure With A on the x-axis and the levels of B shown as two lines, the second contrast tests whether A1-B1 is different from A1-B2 and the third contrast tests whether A2-B1 is different from A2-B2.
"""

# ╔═╡ ab8b2f3a-c8d0-4f74-943b-9426d1c65deb
begin
	cntr8 = Dict(
    :CTR => HypothesisCoding([-1  -1 +1  +1
                              -1  +1  0   0
                               0   0 +1  -1],
            levels=["val", "sod",  "dos", "dod"],
            labels=["do_so", "spt", "grv"])
	);
	m8 = fit(MixedModel, formula, dat1, contrasts=cntr8)
end

# ╔═╡ aef896f9-83d6-45f7-ac28-ba7d1ab7064b
md"""
The three contrasts for one main effect and two nested contrasts are orthogonal. There is no test of the interaction (parallelism).

## Other orthogonal contrasts

For factors with more than four levels there are many options for specifying orthogonal contrasts as long as one proceeds in a top-down strictly hiearchical fashion.

Suppose you have a factor with seven levels and let's ignore shifting colummns. In this case, you have six options for the first contrast, that is 6 vs. 1, 5 vs.2 , 4 vs. 3, 3 vs. 4, 2 vs. 5, and 1 vs. 6 levels.  Then, you specify orthogonal contrasts for partitions with more than 2 elements and so on. That is, you don't specify a contrast that crosses an earlier partition line.

In the following example, after an initial 4 vs 3 partitioning of levels, we specify `AnovaCoding` for the left and `HelmertCoding` for the right partition.
"""

# ╔═╡ 1a9a46c9-42fb-4e27-87af-e9278547eae6
begin
	cntr9 = Dict(
    :CTR => HypothesisCoding(
    [-1/4 -1/4 -1/4 -1/4 +1/3 +1/3 +1/3
     -1/2 -1/2 +1/2 +1/2   0    0    0
     -1/2 +1/2 -1/2 +1/2   0    0    0
     +1/2 -1/2 -1/2 +1/2   0    0    0
       0    0    0    0   -1   +1    0
       0    0    0    0  -1/2 -1/2   1
     ],
    levels=["A1", "A2",  "A3", "A4", "A5", "A6", "A7"],
    labels=["c567.1234", "B", "C", "BxC", "c6.5", "c6.56"])
	)
end

# ╔═╡ 93fcc1b6-98a8-41b5-a84d-46e5bd2df7b9
md"""
There are two rules that hold for all orthogonal contrasts:

1. The weights within rows sum to zero.
2. For all pairs of rows, the sum of the products of weights in the same columns sums to zero.

# Appendix: Summary (Dave Kleinschmidt)

[StatsModels](https://juliastats.org/StatsModels.jl/latest/contrasts/)

StatsModels.jl provides a few commonly used contrast coding schemes,
some less-commonly used schemes, and structs that allow you to manually
specify your own, custom schemes.

## Standard contrasts

The most commonly used contrasts are `DummyCoding` and `EffectsCoding`
(which are similar to `contr.treatment()` and `contr.sum()` in R,
respectively).

## "Exotic" contrasts (rk: well ...)

We also provide `HelmertCoding` and `SeqDiffCoding` (corresponding to
base R's `contr.helmert()` and `MASS::contr.sdif()`).

## Manual contrasts

**ContrastsCoding()**

There are two ways to manually specify contrasts. First, you can specify
them **directly** via `ContrastsCoding`. If you do, it's good practice
to specify the levels corresponding to the rows of the matrix, although
they can be omitted in which case they'll be inferred from the data.

**HypothesisCoding()**

A better way to specify manual contrasts is via `HypothesisCoding`, where each
row of the matrix corresponds to the weights given to the cell means of the
levels corresponding to each column (see [Schad et
al. 2020](https://doi.org/10.1016/j.jml.2019.104038) for more information).
"""

# ╔═╡ 00000000-0000-0000-0000-000000000001
PLUTO_PROJECT_TOML_CONTENTS = """
[deps]
Arrow = "69666777-d1a9-59fb-9406-91d4454c9d45"
Chain = "8be319e6-bccf-4806-a6f7-6fae938471bc"
DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
MixedModels = "ff71e718-51f3-5ec2-a782-8ffcbfa3c316"
StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
StatsModels = "3eaba693-59b7-5ba5-a881-562e759f1c8d"

[compat]
Arrow = "~1.6.2"
Chain = "~0.4.8"
DataFrames = "~1.2.2"
MixedModels = "~4.1.1"
StatsBase = "~0.33.10"
StatsModels = "~0.6.25"
"""

# ╔═╡ 00000000-0000-0000-0000-000000000002
PLUTO_MANIFEST_TOML_CONTENTS = """
# This file is machine-generated - editing it directly is not advised

[[ArgTools]]
uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f"

[[Arrow]]
deps = ["ArrowTypes", "BitIntegers", "CodecLz4", "CodecZstd", "DataAPI", "Dates", "Mmap", "PooledArrays", "SentinelArrays", "Tables", "TimeZones", "UUIDs"]
git-tree-sha1 = "b00e6eaba895683867728e73af78a00218f0db10"
uuid = "69666777-d1a9-59fb-9406-91d4454c9d45"
version = "1.6.2"

[[ArrowTypes]]
deps = ["UUIDs"]
git-tree-sha1 = "a0633b6d6efabf3f76dacd6eb1b3ec6c42ab0552"
uuid = "31f734f8-188a-4ce0-8406-c8a06bd891cd"
version = "1.2.1"

[[Artifacts]]
uuid = "56f22d72-fd6d-98f1-02f0-08ddc0907c33"

[[Base64]]
uuid = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f"

[[BenchmarkTools]]
deps = ["JSON", "Logging", "Printf", "Statistics", "UUIDs"]
git-tree-sha1 = "42ac5e523869a84eac9669eaceed9e4aa0e1587b"
uuid = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf"
version = "1.1.4"

[[BitIntegers]]
deps = ["Random"]
git-tree-sha1 = "f50b5a99aa6ff9db7bf51255b5c21c8bc871ad54"
uuid = "c3b6d118-76ef-56ca-8cc7-ebb389d030a1"
version = "0.2.5"

[[Bzip2_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "19a35467a82e236ff51bc17a3a44b69ef35185a2"
uuid = "6e34b625-4abd-537c-b88f-471c36dfa7a0"
version = "1.0.8+0"

[[Chain]]
git-tree-sha1 = "cac464e71767e8a04ceee82a889ca56502795705"
uuid = "8be319e6-bccf-4806-a6f7-6fae938471bc"
version = "0.4.8"

[[ChainRulesCore]]
deps = ["Compat", "LinearAlgebra", "SparseArrays"]
git-tree-sha1 = "bdc0937269321858ab2a4f288486cb258b9a0af7"
uuid = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4"
version = "1.3.0"

[[CodecBzip2]]
deps = ["Bzip2_jll", "Libdl", "TranscodingStreams"]
git-tree-sha1 = "2e62a725210ce3c3c2e1a3080190e7ca491f18d7"
uuid = "523fee87-0ab8-5b00-afb7-3ecf72e48cfd"
version = "0.7.2"

[[CodecLz4]]
deps = ["Lz4_jll", "TranscodingStreams"]
git-tree-sha1 = "59fe0cb37784288d6b9f1baebddbf75457395d40"
uuid = "5ba52731-8f18-5e0d-9241-30f10d1ec561"
version = "0.4.0"

[[CodecZlib]]
deps = ["TranscodingStreams", "Zlib_jll"]
git-tree-sha1 = "ded953804d019afa9a3f98981d99b33e3db7b6da"
uuid = "944b1d66-785c-5afd-91f1-9de20f533193"
version = "0.7.0"

[[CodecZstd]]
deps = ["TranscodingStreams", "Zstd_jll"]
git-tree-sha1 = "d19cd9ae79ef31774151637492291d75194fc5fa"
uuid = "6b39b394-51ab-5f42-8807-6242bab2b4c2"
version = "0.7.0"

[[Compat]]
deps = ["Base64", "Dates", "DelimitedFiles", "Distributed", "InteractiveUtils", "LibGit2", "Libdl", "LinearAlgebra", "Markdown", "Mmap", "Pkg", "Printf", "REPL", "Random", "SHA", "Serialization", "SharedArrays", "Sockets", "SparseArrays", "Statistics", "Test", "UUIDs", "Unicode"]
git-tree-sha1 = "727e463cfebd0c7b999bbf3e9e7e16f254b94193"
uuid = "34da2185-b29b-5c13-b0c7-acf172513d20"
version = "3.34.0"

[[CompilerSupportLibraries_jll]]
deps = ["Artifacts", "Libdl"]
uuid = "e66e0078-7015-5450-92f7-15fbd957f2ae"

[[Crayons]]
git-tree-sha1 = "3f71217b538d7aaee0b69ab47d9b7724ca8afa0d"
uuid = "a8cc5b0e-0ffa-5ad4-8c14-923d3ee1735f"
version = "4.0.4"

[[DataAPI]]
git-tree-sha1 = "ee400abb2298bd13bfc3df1c412ed228061a2385"
uuid = "9a962f9c-6df0-11e9-0e5d-c546b8b5ee8a"
version = "1.7.0"

[[DataFrames]]
deps = ["Compat", "DataAPI", "Future", "InvertedIndices", "IteratorInterfaceExtensions", "LinearAlgebra", "Markdown", "Missings", "PooledArrays", "PrettyTables", "Printf", "REPL", "Reexport", "SortingAlgorithms", "Statistics", "TableTraits", "Tables", "Unicode"]
git-tree-sha1 = "d785f42445b63fc86caa08bb9a9351008be9b765"
uuid = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
version = "1.2.2"

[[DataStructures]]
deps = ["Compat", "InteractiveUtils", "OrderedCollections"]
git-tree-sha1 = "7d9d316f04214f7efdbb6398d545446e246eff02"
uuid = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8"
version = "0.18.10"

[[DataValueInterfaces]]
git-tree-sha1 = "bfc1187b79289637fa0ef6d4436ebdfe6905cbd6"
uuid = "e2d170a0-9d28-54be-80f0-106bbe20a464"
version = "1.0.0"

[[Dates]]
deps = ["Printf"]
uuid = "ade2ca70-3891-5945-98fb-dc099432e06a"

[[DelimitedFiles]]
deps = ["Mmap"]
uuid = "8bb1440f-4735-579b-a4ab-409b98df4dab"

[[Distributed]]
deps = ["Random", "Serialization", "Sockets"]
uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b"

[[Distributions]]
deps = ["ChainRulesCore", "FillArrays", "LinearAlgebra", "PDMats", "Printf", "QuadGK", "Random", "SparseArrays", "SpecialFunctions", "Statistics", "StatsBase", "StatsFuns"]
git-tree-sha1 = "c2dbc7e0495c3f956e4615b78d03c7aa10091d0c"
uuid = "31c24e10-a181-5473-b8eb-7969acd0382f"
version = "0.25.15"

[[DocStringExtensions]]
deps = ["LibGit2"]
git-tree-sha1 = "a32185f5428d3986f47c2ab78b1f216d5e6cc96f"
uuid = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae"
version = "0.8.5"

[[Downloads]]
deps = ["ArgTools", "LibCURL", "NetworkOptions"]
uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6"

[[ExprTools]]
git-tree-sha1 = "b7e3d17636b348f005f11040025ae8c6f645fe92"
uuid = "e2ba6199-217a-4e67-a87a-7c52f15ade04"
version = "0.1.6"

[[FillArrays]]
deps = ["LinearAlgebra", "Random", "SparseArrays", "Statistics"]
git-tree-sha1 = "7c365bdef6380b29cfc5caaf99688cd7489f9b87"
uuid = "1a297f60-69ca-5386-bcde-b61e274b549b"
version = "0.12.2"

[[Formatting]]
deps = ["Printf"]
git-tree-sha1 = "8339d61043228fdd3eb658d86c926cb282ae72a8"
uuid = "59287772-0a20-5a39-b81b-1366585eb4c0"
version = "0.4.2"

[[Future]]
deps = ["Random"]
uuid = "9fa8497b-333b-5362-9e8d-4d0656e87820"

[[GLM]]
deps = ["Distributions", "LinearAlgebra", "Printf", "Reexport", "SparseArrays", "SpecialFunctions", "Statistics", "StatsBase", "StatsFuns", "StatsModels"]
git-tree-sha1 = "f564ce4af5e79bb88ff1f4488e64363487674278"
uuid = "38e38edf-8417-5370-95a0-9cbb8c7f171a"
version = "1.5.1"

[[HTTP]]
deps = ["Base64", "Dates", "IniFile", "Logging", "MbedTLS", "NetworkOptions", "Sockets", "URIs"]
git-tree-sha1 = "60ed5f1643927479f845b0135bb369b031b541fa"
uuid = "cd3eb016-35fb-5094-929b-558a96fad6f3"
version = "0.9.14"

[[IniFile]]
deps = ["Test"]
git-tree-sha1 = "098e4d2c533924c921f9f9847274f2ad89e018b8"
uuid = "83e8ac13-25f8-5344-8a64-a9f2b223428f"
version = "0.5.0"

[[InteractiveUtils]]
deps = ["Markdown"]
uuid = "b77e0a4c-d291-57a0-90e8-8db25a27a240"

[[InvertedIndices]]
deps = ["Test"]
git-tree-sha1 = "15732c475062348b0165684ffe28e85ea8396afc"
uuid = "41ab1584-1d38-5bbf-9106-f11c6c58b48f"
version = "1.0.0"

[[IrrationalConstants]]
git-tree-sha1 = "f76424439413893a832026ca355fe273e93bce94"
uuid = "92d709cd-6900-40b7-9082-c6be49f344b6"
version = "0.1.0"

[[IteratorInterfaceExtensions]]
git-tree-sha1 = "a3f24677c21f5bbe9d2a714f95dcd58337fb2856"
uuid = "82899510-4779-5014-852e-03e436cf321d"
version = "1.0.0"

[[JLLWrappers]]
deps = ["Preferences"]
git-tree-sha1 = "642a199af8b68253517b80bd3bfd17eb4e84df6e"
uuid = "692b3bcd-3c85-4b1f-b108-f13ce0eb3210"
version = "1.3.0"

[[JSON]]
deps = ["Dates", "Mmap", "Parsers", "Unicode"]
git-tree-sha1 = "8076680b162ada2a031f707ac7b4953e30667a37"
uuid = "682c06a0-de6a-54ab-a142-c8b1cf79cde6"
version = "0.21.2"

[[JSON3]]
deps = ["Dates", "Mmap", "Parsers", "StructTypes", "UUIDs"]
git-tree-sha1 = "b3e5984da3c6c95bcf6931760387ff2e64f508f3"
uuid = "0f8b85d8-7281-11e9-16c2-39a750bddbf1"
version = "1.9.1"

[[JSONSchema]]
deps = ["HTTP", "JSON", "ZipFile"]
git-tree-sha1 = "b84ab8139afde82c7c65ba2b792fe12e01dd7307"
uuid = "7d188eb4-7ad8-530c-ae41-71a32a6d4692"
version = "0.3.3"

[[LazyArtifacts]]
deps = ["Artifacts", "Pkg"]
uuid = "4af54fe1-eca0-43a8-85a7-787d91b784e3"

[[LibCURL]]
deps = ["LibCURL_jll", "MozillaCACerts_jll"]
uuid = "b27032c2-a3e7-50c8-80cd-2d36dbcbfd21"

[[LibCURL_jll]]
deps = ["Artifacts", "LibSSH2_jll", "Libdl", "MbedTLS_jll", "Zlib_jll", "nghttp2_jll"]
uuid = "deac9b47-8bc7-5906-a0fe-35ac56dc84c0"

[[LibGit2]]
deps = ["Base64", "NetworkOptions", "Printf", "SHA"]
uuid = "76f85450-5226-5b5a-8eaa-529ad045b433"

[[LibSSH2_jll]]
deps = ["Artifacts", "Libdl", "MbedTLS_jll"]
uuid = "29816b5a-b9ab-546f-933c-edad1886dfa8"

[[Libdl]]
uuid = "8f399da3-3557-5675-b5ff-fb832c97cbdb"

[[LinearAlgebra]]
deps = ["Libdl"]
uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"

[[LogExpFunctions]]
deps = ["DocStringExtensions", "IrrationalConstants", "LinearAlgebra"]
git-tree-sha1 = "3d682c07e6dd250ed082f883dc88aee7996bf2cc"
uuid = "2ab3a3ac-af41-5b50-aa03-7779005ae688"
version = "0.3.0"

[[Logging]]
uuid = "56ddb016-857b-54e1-b83d-db4d58db5568"

[[Lz4_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "5d494bc6e85c4c9b626ee0cab05daa4085486ab1"
uuid = "5ced341a-0733-55b8-9ab6-a4889d929147"
version = "1.9.3+0"

[[Markdown]]
deps = ["Base64"]
uuid = "d6f4376e-aef5-505a-96c1-9c027394607a"

[[MathOptInterface]]
deps = ["BenchmarkTools", "CodecBzip2", "CodecZlib", "JSON", "JSONSchema", "LinearAlgebra", "MutableArithmetics", "OrderedCollections", "SparseArrays", "Test", "Unicode"]
git-tree-sha1 = "575644e3c05b258250bb599e57cf73bbf1062901"
uuid = "b8f27783-ece8-5eb3-8dc8-9495eed66fee"
version = "0.9.22"

[[MathProgBase]]
deps = ["LinearAlgebra", "SparseArrays"]
git-tree-sha1 = "9abbe463a1e9fc507f12a69e7f29346c2cdc472c"
uuid = "fdba3010-5040-5b88-9595-932c9decdf73"
version = "0.7.8"

[[MbedTLS]]
deps = ["Dates", "MbedTLS_jll", "Random", "Sockets"]
git-tree-sha1 = "1c38e51c3d08ef2278062ebceade0e46cefc96fe"
uuid = "739be429-bea8-5141-9913-cc70e7f3736d"
version = "1.0.3"

[[MbedTLS_jll]]
deps = ["Artifacts", "Libdl"]
uuid = "c8ffd9c3-330d-5841-b78e-0817d7145fa1"

[[Missings]]
deps = ["DataAPI"]
git-tree-sha1 = "2ca267b08821e86c5ef4376cffed98a46c2cb205"
uuid = "e1d29d7a-bbdc-5cf2-9ac0-f12de2c33e28"
version = "1.0.1"

[[MixedModels]]
deps = ["Arrow", "DataAPI", "Distributions", "GLM", "JSON3", "LazyArtifacts", "LinearAlgebra", "Markdown", "NLopt", "PooledArrays", "ProgressMeter", "Random", "SparseArrays", "StaticArrays", "Statistics", "StatsBase", "StatsFuns", "StatsModels", "StructTypes", "Tables"]
git-tree-sha1 = "f318e42a48ec0a856292bafeec6b07aed3f6d600"
uuid = "ff71e718-51f3-5ec2-a782-8ffcbfa3c316"
version = "4.1.1"

[[Mmap]]
uuid = "a63ad114-7e13-5084-954f-fe012c677804"

[[Mocking]]
deps = ["ExprTools"]
git-tree-sha1 = "748f6e1e4de814b101911e64cc12d83a6af66782"
uuid = "78c3b35d-d492-501b-9361-3d52fe80e533"
version = "0.7.2"

[[MozillaCACerts_jll]]
uuid = "14a3606d-f60d-562e-9121-12d972cd8159"

[[MutableArithmetics]]
deps = ["LinearAlgebra", "SparseArrays", "Test"]
git-tree-sha1 = "3927848ccebcc165952dc0d9ac9aa274a87bfe01"
uuid = "d8a4904e-b15c-11e9-3269-09a3773c0cb0"
version = "0.2.20"

[[NLopt]]
deps = ["MathOptInterface", "MathProgBase", "NLopt_jll"]
git-tree-sha1 = "d80cb3327d1aeef0f59eacf225e000f86e4eee0a"
uuid = "76087f3c-5699-56af-9a33-bf431cd00edd"
version = "0.6.3"

[[NLopt_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "2b597c46900f5f811bec31f0dcc88b45744a2a09"
uuid = "079eb43e-fd8e-5478-9966-2cf3e3edb778"
version = "2.7.0+0"

[[NetworkOptions]]
uuid = "ca575930-c2e3-43a9-ace4-1e988b2c1908"

[[OpenSpecFun_jll]]
deps = ["Artifacts", "CompilerSupportLibraries_jll", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "13652491f6856acfd2db29360e1bbcd4565d04f1"
uuid = "efe28fd5-8261-553b-a9e1-b2916fc3738e"
version = "0.5.5+0"

[[OrderedCollections]]
git-tree-sha1 = "85f8e6578bf1f9ee0d11e7bb1b1456435479d47c"
uuid = "bac558e1-5e72-5ebc-8fee-abe8a469f55d"
version = "1.4.1"

[[PDMats]]
deps = ["LinearAlgebra", "SparseArrays", "SuiteSparse"]
git-tree-sha1 = "4dd403333bcf0909341cfe57ec115152f937d7d8"
uuid = "90014a1f-27ba-587c-ab20-58faa44d9150"
version = "0.11.1"

[[Parsers]]
deps = ["Dates"]
git-tree-sha1 = "438d35d2d95ae2c5e8780b330592b6de8494e779"
uuid = "69de0a69-1ddd-5017-9359-2bf0b02dc9f0"
version = "2.0.3"

[[Pkg]]
deps = ["Artifacts", "Dates", "Downloads", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "REPL", "Random", "SHA", "Serialization", "TOML", "Tar", "UUIDs", "p7zip_jll"]
uuid = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f"

[[PooledArrays]]
deps = ["DataAPI", "Future"]
git-tree-sha1 = "cde4ce9d6f33219465b55162811d8de8139c0414"
uuid = "2dfb63ee-cc39-5dd5-95bd-886bf059d720"
version = "1.2.1"

[[Preferences]]
deps = ["TOML"]
git-tree-sha1 = "00cfd92944ca9c760982747e9a1d0d5d86ab1e5a"
uuid = "21216c6a-2e73-6563-6e65-726566657250"
version = "1.2.2"

[[PrettyTables]]
deps = ["Crayons", "Formatting", "Markdown", "Reexport", "Tables"]
git-tree-sha1 = "0d1245a357cc61c8cd61934c07447aa569ff22e6"
uuid = "08abe8d2-0d0c-5749-adfa-8a2ac140af0d"
version = "1.1.0"

[[Printf]]
deps = ["Unicode"]
uuid = "de0858da-6303-5e67-8744-51eddeeeb8d7"

[[ProgressMeter]]
deps = ["Distributed", "Printf"]
git-tree-sha1 = "afadeba63d90ff223a6a48d2009434ecee2ec9e8"
uuid = "92933f4c-e287-5a05-a399-4b506db050ca"
version = "1.7.1"

[[QuadGK]]
deps = ["DataStructures", "LinearAlgebra"]
git-tree-sha1 = "12fbe86da16df6679be7521dfb39fbc861e1dc7b"
uuid = "1fd47b50-473d-5c70-9696-f719f8f3bcdc"
version = "2.4.1"

[[REPL]]
deps = ["InteractiveUtils", "Markdown", "Sockets", "Unicode"]
uuid = "3fa0cd96-eef1-5676-8a61-b3b8758bbffb"

[[Random]]
deps = ["Serialization"]
uuid = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"

[[RecipesBase]]
git-tree-sha1 = "44a75aa7a527910ee3d1751d1f0e4148698add9e"
uuid = "3cdcf5f2-1ef4-517c-9805-6587b60abb01"
version = "1.1.2"

[[Reexport]]
git-tree-sha1 = "45e428421666073eab6f2da5c9d310d99bb12f9b"
uuid = "189a3867-3050-52da-a836-e630ba90ab69"
version = "1.2.2"

[[Rmath]]
deps = ["Random", "Rmath_jll"]
git-tree-sha1 = "bf3188feca147ce108c76ad82c2792c57abe7b1f"
uuid = "79098fc4-a85e-5d69-aa6a-4863f24498fa"
version = "0.7.0"

[[Rmath_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "68db32dff12bb6127bac73c209881191bf0efbb7"
uuid = "f50d1b31-88e8-58de-be2c-1cc44531875f"
version = "0.3.0+0"

[[SHA]]
uuid = "ea8e919c-243c-51af-8825-aaa63cd721ce"

[[SentinelArrays]]
deps = ["Dates", "Random"]
git-tree-sha1 = "54f37736d8934a12a200edea2f9206b03bdf3159"
uuid = "91c51154-3ec4-41a3-a24f-3f23e20d615c"
version = "1.3.7"

[[Serialization]]
uuid = "9e88b42a-f829-5b0c-bbe9-9e923198166b"

[[SharedArrays]]
deps = ["Distributed", "Mmap", "Random", "Serialization"]
uuid = "1a1011a3-84de-559e-8e89-a11a2f7dc383"

[[ShiftedArrays]]
git-tree-sha1 = "22395afdcf37d6709a5a0766cc4a5ca52cb85ea0"
uuid = "1277b4bf-5013-50f5-be3d-901d8477a67a"
version = "1.0.0"

[[Sockets]]
uuid = "6462fe0b-24de-5631-8697-dd941f90decc"

[[SortingAlgorithms]]
deps = ["DataStructures"]
git-tree-sha1 = "b3363d7460f7d098ca0912c69b082f75625d7508"
uuid = "a2af1166-a08f-5f64-846c-94a0d3cef48c"
version = "1.0.1"

[[SparseArrays]]
deps = ["LinearAlgebra", "Random"]
uuid = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"

[[SpecialFunctions]]
deps = ["ChainRulesCore", "LogExpFunctions", "OpenSpecFun_jll"]
git-tree-sha1 = "a322a9493e49c5f3a10b50df3aedaf1cdb3244b7"
uuid = "276daf66-3868-5448-9aa4-cd146d93841b"
version = "1.6.1"

[[StaticArrays]]
deps = ["LinearAlgebra", "Random", "Statistics"]
git-tree-sha1 = "3240808c6d463ac46f1c1cd7638375cd22abbccb"
uuid = "90137ffa-7385-5640-81b9-e52037218182"
version = "1.2.12"

[[Statistics]]
deps = ["LinearAlgebra", "SparseArrays"]
uuid = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"

[[StatsAPI]]
git-tree-sha1 = "1958272568dc176a1d881acb797beb909c785510"
uuid = "82ae8749-77ed-4fe6-ae5f-f523153014b0"
version = "1.0.0"

[[StatsBase]]
deps = ["DataAPI", "DataStructures", "LinearAlgebra", "Missings", "Printf", "Random", "SortingAlgorithms", "SparseArrays", "Statistics", "StatsAPI"]
git-tree-sha1 = "8cbbc098554648c84f79a463c9ff0fd277144b6c"
uuid = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
version = "0.33.10"

[[StatsFuns]]
deps = ["ChainRulesCore", "IrrationalConstants", "LogExpFunctions", "Reexport", "Rmath", "SpecialFunctions"]
git-tree-sha1 = "46d7ccc7104860c38b11966dd1f72ff042f382e4"
uuid = "4c63d2b9-4356-54db-8cca-17b64c39e42c"
version = "0.9.10"

[[StatsModels]]
deps = ["DataAPI", "DataStructures", "LinearAlgebra", "Printf", "ShiftedArrays", "SparseArrays", "StatsBase", "StatsFuns", "Tables"]
git-tree-sha1 = "3fa15c1f8be168e76d59097f66970adc86bfeb95"
uuid = "3eaba693-59b7-5ba5-a881-562e759f1c8d"
version = "0.6.25"

[[StructTypes]]
deps = ["Dates", "UUIDs"]
git-tree-sha1 = "8445bf99a36d703a09c601f9a57e2f83000ef2ae"
uuid = "856f2bd8-1eba-4b0a-8007-ebc267875bd4"
version = "1.7.3"

[[SuiteSparse]]
deps = ["Libdl", "LinearAlgebra", "Serialization", "SparseArrays"]
uuid = "4607b0f0-06f3-5cda-b6b1-a6196a1729e9"

[[TOML]]
deps = ["Dates"]
uuid = "fa267f1f-6049-4f14-aa54-33bafae1ed76"

[[TableTraits]]
deps = ["IteratorInterfaceExtensions"]
git-tree-sha1 = "c06b2f539df1c6efa794486abfb6ed2022561a39"
uuid = "3783bdb8-4a98-5b6b-af9a-565f29a5fe9c"
version = "1.0.1"

[[Tables]]
deps = ["DataAPI", "DataValueInterfaces", "IteratorInterfaceExtensions", "LinearAlgebra", "TableTraits", "Test"]
git-tree-sha1 = "d0c690d37c73aeb5ca063056283fde5585a41710"
uuid = "bd369af6-aec1-5ad0-b16a-f7cc5008161c"
version = "1.5.0"

[[Tar]]
deps = ["ArgTools", "SHA"]
uuid = "a4e569a6-e804-4fa4-b0f3-eef7a1d5b13e"

[[Test]]
deps = ["InteractiveUtils", "Logging", "Random", "Serialization"]
uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40"

[[TimeZones]]
deps = ["Dates", "Future", "LazyArtifacts", "Mocking", "Pkg", "Printf", "RecipesBase", "Serialization", "Unicode"]
git-tree-sha1 = "6c9040665b2da00d30143261aea22c7427aada1c"
uuid = "f269a46b-ccf7-5d73-abea-4c690281aa53"
version = "1.5.7"

[[TranscodingStreams]]
deps = ["Random", "Test"]
git-tree-sha1 = "216b95ea110b5972db65aa90f88d8d89dcb8851c"
uuid = "3bb67fe8-82b1-5028-8e26-92a6c54297fa"
version = "0.9.6"

[[URIs]]
git-tree-sha1 = "97bbe755a53fe859669cd907f2d96aee8d2c1355"
uuid = "5c2747f8-b7ea-4ff2-ba2e-563bfd36b1d4"
version = "1.3.0"

[[UUIDs]]
deps = ["Random", "SHA"]
uuid = "cf7118a7-6976-5b1a-9a39-7adc72f591a4"

[[Unicode]]
uuid = "4ec0a83e-493e-50e2-b9ac-8f72acf5a8f5"

[[ZipFile]]
deps = ["Libdl", "Printf", "Zlib_jll"]
git-tree-sha1 = "c3a5637e27e914a7a445b8d0ad063d701931e9f7"
uuid = "a5390f91-8eb1-5f08-bee0-b1d1ffed6cea"
version = "0.9.3"

[[Zlib_jll]]
deps = ["Libdl"]
uuid = "83775a58-1f1d-513f-b197-d71354ab007a"

[[Zstd_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "cc4bf3fdde8b7e3e9fa0351bdeedba1cf3b7f6e6"
uuid = "3161d3a3-bdf6-5164-811a-617609db77b4"
version = "1.5.0+0"

[[nghttp2_jll]]
deps = ["Artifacts", "Libdl"]
uuid = "8e850ede-7688-5339-a07c-302acd2aaf8d"

[[p7zip_jll]]
deps = ["Artifacts", "Libdl"]
uuid = "3f19e933-33d8-53b3-aaab-bd5110c3b7a0"
"""

# ╔═╡ Cell order:
# ╟─09a33b05-568e-427a-bd66-812d271d1791
# ╠═7a3de02b-c423-496e-bf92-9981d71e5eab
# ╟─28c37807-b83f-4c3e-848d-a7401b4eda29
# ╠═ccf701bb-cf16-4339-bd93-64ffe87435bb
# ╟─cb63bd4b-5f0e-49d1-a257-2b83c1687975
# ╠═8ffba70f-fc97-47a2-b973-3c0f2047eb9b
# ╟─4e9359a3-2533-4ce8-b2fd-26411e84f59c
# ╠═b00e953a-b279-4ec3-8744-5518d1242d8e
# ╟─7a3749b3-d2e9-4b53-9e57-975ee428b416
# ╠═6393d21d-a3a6-4463-8981-052440ba7887
# ╟─d087afe8-9d2b-489e-9160-206f22486505
# ╠═bb97e862-24ad-4ce5-9011-429a0e69e5ab
# ╟─71528268-41fe-4de5-8725-1c65c4839c65
# ╠═e2372cb1-03d4-4404-bfef-640ad5653792
# ╟─142c86d3-e340-494c-8b52-8f2013993a19
# ╠═0bb71abe-433c-44e3-bff4-eb26e63bee2d
# ╟─0cf43493-af05-4744-a154-d7209251ae8f
# ╠═2d778d8e-949c-4a7b-a774-770c6bc4bce6
# ╟─a71103ca-a680-45ac-ba02-61a5660cc157
# ╠═00bedbae-0bc4-4e4e-85d5-91ffc8ea77ee
# ╟─f504edce-35de-4938-a3c6-08d4fc8f2e66
# ╟─a12c0a44-f02e-479e-94f9-c5d09c3b25ae
# ╠═ab8b2f3a-c8d0-4f74-943b-9426d1c65deb
# ╟─aef896f9-83d6-45f7-ac28-ba7d1ab7064b
# ╠═1a9a46c9-42fb-4e27-87af-e9278547eae6
# ╟─93fcc1b6-98a8-41b5-a84d-46e5bd2df7b9
# ╟─00000000-0000-0000-0000-000000000001
# ╟─00000000-0000-0000-0000-000000000002