{"spec_id":"heatmap-clustered","library":"makie","language":"julia","code":"# anyplot.ai\n# heatmap-clustered: Clustered Heatmap\n# Library: makie 0.21.9 | Julia 1.11.9\n# Quality: 79/100 | Created: 2026-09-05\n\nusing CairoMakie\nusing Colors\nusing ColorSchemes\nusing Random\n\nRandom.seed!(42)\n\n# --- Theme tokens -----------------------------------------------------------\nconst THEME    = get(ENV, \"ANYPLOT_THEME\", \"light\")\nconst PAGE_BG  = THEME == \"light\" ? colorant\"#FAF8F1\" : colorant\"#1A1A17\"\nconst INK      = THEME == \"light\" ? colorant\"#1A1A17\" : colorant\"#F0EFE8\"\nconst INK_SOFT = THEME == \"light\" ? colorant\"#4A4A44\" : colorant\"#B8B7B0\"\nconst _MIDPOINT = THEME == \"light\" ? colorant\"#FAF8F1\" : colorant\"#1A1A17\"\nconst IMPRINT_DIV = cgrad([colorant\"#AE3030\", _MIDPOINT, colorant\"#4467A3\"])\nconst IMPRINT_PALETTE = [colorant\"#009E73\", colorant\"#C475FD\", colorant\"#4467A3\"]\n\n# --- Hierarchical clustering (Ward's method, Lance-Williams update) --------\n\nfunction pairwise_sqdist(X::AbstractMatrix{<:Real})\n    n = size(X, 1)\n    D = zeros(Float64, n, n)\n    for i in 1:n, j in (i + 1):n\n        d = sum((X[i, :] .- X[j, :]) .^ 2)\n        D[i, j] = d\n        D[j, i] = d\n    end\n    return D\nend\n\nfunction ward_linkage(D::Matrix{Float64})\n    n = size(D, 1)\n    dist = fill(Inf, 2n - 1, 2n - 1)\n    dist[1:n, 1:n] .= D\n    sizes = ones(Int, 2n - 1)\n    active = collect(1:n)\n    merge_a = zeros(Int, n - 1)\n    merge_b = zeros(Int, n - 1)\n    merge_h = zeros(Float64, n - 1)\n\n    for step in 1:(n - 1)\n        best = Inf\n        bi, bj = active[1], active[2]\n        for ai in 1:length(active), aj in (ai + 1):length(active)\n            i, j = active[ai], active[aj]\n            d = dist[min(i, j), max(i, j)]\n            if d < best\n                best = d\n                bi, bj = i, j\n            end\n        end\n\n        new_id = n + step\n        merge_a[step] = bi\n        merge_b[step] = bj\n        merge_h[step] = best\n\n        ni, nj = sizes[bi], sizes[bj]\n        for a in active\n            if a == bi || a == bj\n                continue\n            end\n            nm = sizes[a]\n            d_im = dist[min(bi, a), max(bi, a)]\n            d_jm = dist[min(bj, a), max(bj, a)]\n            dist[min(a, new_id), max(a, new_id)] = ((ni + nm) * d_im + (nj + nm) * d_jm - nm * best) / (ni + nj + nm)\n        end\n\n        sizes[new_id] = ni + nj\n        filter!(x -> x != bi && x != bj, active)\n        push!(active, new_id)\n    end\n\n    return merge_a, merge_b, merge_h\nend\n\nfunction collect_leaves!(order, xpos, merge_a, merge_b, n, node_id)\n    if node_id <= n\n        push!(order, node_id)\n        xpos[node_id] = Float64(length(order))\n    else\n        step = node_id - n\n        collect_leaves!(order, xpos, merge_a, merge_b, n, merge_a[step])\n        collect_leaves!(order, xpos, merge_a, merge_b, n, merge_b[step])\n        xpos[node_id] = (xpos[merge_a[step]] + xpos[merge_b[step]]) / 2\n    end\n    return nothing\nend\n\nfunction dendrogram_segments(merge_a, merge_b, merge_h, xpos, n)\n    pos_pts = Float64[]\n    height_pts = Float64[]\n    for step in 1:length(merge_h)\n        a, b, h = merge_a[step], merge_b[step], merge_h[step]\n        ha = a <= n ? 0.0 : merge_h[a - n]\n        hb = b <= n ? 0.0 : merge_h[b - n]\n        pa, pb = xpos[a], xpos[b]\n        append!(pos_pts, (pa, pa, NaN, pb, pb, NaN, pa, pb, NaN))\n        append!(height_pts, (ha, h, NaN, hb, h, NaN, h, h, NaN))\n    end\n    return pos_pts, height_pts\nend\n\n# --- Data --------------------------------------------------------------------\n# Gene expression across control + two treatment conditions, four replicates\n# each. Genes fall into four latent co-expression modules so clustering\n# recovers a block structure once rows and columns are reordered.\nn_genes = 16\nn_samples = 12\n\nconditions = repeat([\"Ctrl\", \"TreatA\", \"TreatB\"], inner = 4)\nreplicate = repeat(1:4, outer = 3)\ncolumn_labels = [conditions[i] * \"-\" * string(replicate[i]) for i in 1:n_samples]\nrow_labels = [\"Gene \" * lpad(string(g), 2, \"0\") for g in 1:n_genes]\n\ngene_module = repeat(1:4, inner = 4)\nmodule_pattern = Dict(\n    1 => Dict(\"Ctrl\" => 2.0, \"TreatA\" => -1.5, \"TreatB\" => -1.0),\n    2 => Dict(\"Ctrl\" => -1.5, \"TreatA\" => 2.0, \"TreatB\" => 0.5),\n    3 => Dict(\"Ctrl\" => -0.5, \"TreatA\" => -0.5, \"TreatB\" => 2.2),\n    4 => Dict(\"Ctrl\" => 0.2, \"TreatA\" => -0.2, \"TreatB\" => 0.1),\n)\n\nexpression = zeros(Float64, n_genes, n_samples)\nfor g in 1:n_genes, s in 1:n_samples\n    expression[g, s] = module_pattern[gene_module[g]][conditions[s]] + randn() * 0.6\nend\n\nrow_mean = [sum(expression[g, :]) / n_samples for g in 1:n_genes]\nrow_std = [sqrt(sum((expression[g, :] .- row_mean[g]) .^ 2) / n_samples) for g in 1:n_genes]\nzscore = (expression .- row_mean) ./ row_std\nzmax = maximum(abs.(zscore))\n\n# --- Clustering ---------------------------------------------------------------\nrow_merge_a, row_merge_b, row_merge_h = ward_linkage(pairwise_sqdist(zscore))\ncol_merge_a, col_merge_b, col_merge_h = ward_linkage(pairwise_sqdist(Matrix(zscore')))\n\nrow_order = Int[]\nrow_xpos = Dict{Int,Float64}()\ncollect_leaves!(row_order, row_xpos, row_merge_a, row_merge_b, n_genes, 2 * n_genes - 1)\n\ncol_order = Int[]\ncol_xpos = Dict{Int,Float64}()\ncollect_leaves!(col_order, col_xpos, col_merge_a, col_merge_b, n_samples, 2 * n_samples - 1)\n\nzscore_reordered = zscore[row_order, col_order]\nrow_labels_reordered = row_labels[row_order]\ncol_labels_reordered = column_labels[col_order]\n\n# Ctrl/TreatA/TreatB group-color annotation strip, reordered alongside columns.\ncondition_colors = Dict(\n    \"Ctrl\" => IMPRINT_PALETTE[1], \"TreatA\" => IMPRINT_PALETTE[2], \"TreatB\" => IMPRINT_PALETTE[3],\n)\ncol_group_img = reshape([condition_colors[conditions[i]] for i in col_order], n_samples, 1)\n\nrow_pos, row_height = dendrogram_segments(row_merge_a, row_merge_b, row_merge_h, row_xpos, n_genes)\ncol_pos, col_height = dendrogram_segments(col_merge_a, col_merge_b, col_merge_h, col_xpos, n_samples)\nrow_max_height = maximum(row_merge_h)\ncol_max_height = maximum(col_merge_h)\n\n# --- Plot ----------------------------------------------------------------------\nfig = Figure(\n    size            = (1200, 1200),\n    fontsize        = 14,\n    backgroundcolor = PAGE_BG,\n)\n\nLabel(\n    fig[1, 1:3], \"heatmap-clustered · julia · makie · anyplot.ai\";\n    fontsize = 20, color = INK, font = :bold,\n)\n\ncol_dendro_ax = Axis(\n    fig[2, 2];\n    backgroundcolor = PAGE_BG,\n    limits          = (0.5, n_samples + 0.5, 0.0, col_max_height * 1.05),\n)\nlines!(col_dendro_ax, col_pos, col_height; color = INK_SOFT, linewidth = 1.6)\nhidedecorations!(col_dendro_ax)\nhidespines!(col_dendro_ax)\n\nLegend(\n    fig[2, 3],\n    [PolyElement(color = c) for c in IMPRINT_PALETTE],\n    [\"Ctrl\", \"TreatA\", \"TreatB\"],\n    \"Group\";\n    labelcolor = INK_SOFT, titlecolor = INK, framevisible = false,\n    labelsize = 12, titlesize = 12, patchsize = (12, 12),\n)\n\ncol_group_ax = Axis(\n    fig[3, 2];\n    backgroundcolor = PAGE_BG,\n    limits          = (0.5, n_samples + 0.5, 0.0, 1.0),\n)\nimage!(col_group_ax, 0.5 .. (n_samples + 0.5), 0.0 .. 1.0, col_group_img; interpolate = false)\nhidedecorations!(col_group_ax)\nhidespines!(col_group_ax)\n\nrow_dendro_ax = Axis(\n    fig[4, 1];\n    backgroundcolor = PAGE_BG,\n    limits          = (0.0, row_max_height * 1.05, 0.5, n_genes + 0.5),\n    xreversed       = true,\n    yreversed       = true,\n)\nlines!(row_dendro_ax, row_height, row_pos; color = INK_SOFT, linewidth = 1.6)\nhidedecorations!(row_dendro_ax)\nhidespines!(row_dendro_ax)\n\nheat_ax = Axis(\n    fig[4, 2];\n    backgroundcolor      = PAGE_BG,\n    limits               = (0.5, n_samples + 0.5, 0.5, n_genes + 0.5),\n    yreversed            = true,\n    xticks               = (1:n_samples, col_labels_reordered),\n    yticks               = (1:n_genes, row_labels_reordered),\n    xticklabelrotation   = pi / 4,\n    xticklabelcolor      = INK_SOFT,\n    yticklabelcolor      = INK_SOFT,\n    xticklabelsize       = 12,\n    yticklabelsize       = 12,\n    yaxisposition        = :right,\n    xgridvisible         = false,\n    ygridvisible         = false,\n    topspinevisible      = false,\n    rightspinevisible    = false,\n    leftspinevisible     = false,\n    bottomspinevisible   = false,\n)\nhm = heatmap!(\n    heat_ax, 1:n_samples, 1:n_genes, permutedims(zscore_reordered);\n    colormap = IMPRINT_DIV, colorrange = (-zmax, zmax),\n)\n\n# Thin theme-adaptive cell grid + outer frame so near-zero cells (which sit at\n# the midpoint color, equal to the page background) stay visually separated\n# from each other and from the canvas in both themes.\ngrid_color = RGBAf(INK.r, INK.g, INK.b, 0.15)\nv_segs = Point2f[]\nfor gx in 0.5:1.0:(n_samples + 0.5)\n    push!(v_segs, Point2f(gx, 0.5), Point2f(gx, n_genes + 0.5))\nend\nh_segs = Point2f[]\nfor gy in 0.5:1.0:(n_genes + 0.5)\n    push!(h_segs, Point2f(0.5, gy), Point2f(n_samples + 0.5, gy))\nend\nlinesegments!(heat_ax, v_segs; color = grid_color, linewidth = 0.75)\nlinesegments!(heat_ax, h_segs; color = grid_color, linewidth = 0.75)\n\nColorbar(\n    fig[4, 3], hm;\n    label = \"Expression (row z-score)\", labelcolor = INK,\n    ticklabelcolor = INK_SOFT, width = 18,\n)\n\ncolsize!(fig.layout, 1, Relative(0.14))\ncolsize!(fig.layout, 2, Relative(0.62))\ncolsize!(fig.layout, 3, Relative(0.24))\nrowsize!(fig.layout, 2, Relative(0.13))\nrowsize!(fig.layout, 3, Relative(0.05))\nrowsize!(fig.layout, 4, Relative(0.70))\n\n# --- Save -----------------------------------------------------------------------\nsave(\"plot-$(THEME).png\", fig; px_per_unit = 2)\n"}