Gene-to-distance-to-plaque correlation — microglia only

Microglia object provided by Kai Saito with Akhil HALO annotations

Published

April 27, 2026

By restricting the analysis to microglia only, we remove the composition confound. If a gene still correlates with distance within a single cell type, the association more likely reflects a genuine distance-dependent transcriptional programme rather than the fact that microglia are (or are not) concentrated near plaques.

Setup and data loading

Code 
source("R/setup.R")
library(ggrepel)
library(patchwork)
library(ggVennDiagram)
library(gridExtra)

data <- load_slides_and_distances("microglia")
slide1_mg <- data$slide1; slide2_mg <- data$slide2
cell_distances_1 <- data$distances_1; cell_distances_2 <- data$distances_2
cell_distances_all <- data$distances_all

message(sprintf(
  "Microglia: slide 1 = %d cells, slide 2 = %d cells (total = %d)",
  ncol(slide1_mg), ncol(slide2_mg), ncol(slide1_mg) + ncol(slide2_mg)
))
Code 
# Add proximity labels (50 um threshold) — needed for treatment DE later
add_proximity <- function(seurat_obj, cell_distances) {
  prox <- cell_distances |>
    mutate(
      prox_caa   = if_else(dist_caa < 50,   "proximal", "distal"),
      prox_paren = if_else(dist_paren < 50,  "proximal", "distal"),
      prox_any   = if_else(pmin(dist_caa, dist_paren) < 50, "proximal", "distal")
    ) |>
    column_to_rownames("cell") |>
    select(dist_caa, dist_paren, prox_caa, prox_paren, prox_any)

  AddMetaData(seurat_obj, prox)
}

slide1_mg <- add_proximity(slide1_mg, cell_distances_1)
slide2_mg <- add_proximity(slide2_mg, cell_distances_2)

cell_distances_all <- cell_distances_all |>
  mutate(
    prox_caa   = if_else(dist_caa < 50,   "proximal", "distal"),
    prox_paren = if_else(dist_paren < 50,  "proximal", "distal"),
    prox_any   = if_else(dist_any < 50,    "proximal", "distal")
  )

Cell counts per group

Code 
# sample_ID and Treatment.Group are already in cell_distances_all (joined in load_slides_and_distances)
cell_distances_all |>
  group_by(sample_ID, Treatment.Group, prox_any) |>
  summarise(n = n(), .groups = "drop") |>
  pivot_wider(names_from = prox_any, values_from = n, values_fill = 0) |>
  mutate(total = proximal + distal) |>
  kbl(caption = "Microglia per sample — proximity to any plaque (50 um threshold)") |>
  kable_styling("striped", full_width = FALSE)
Microglia per sample — proximity to any plaque (50 um threshold)
sample_ID Treatment.Group distal proximal total
KK4_464 Adu 1488 2788 4276
KK4_465 IgG 840 2565 3405
KK4_492 Adu 2020 3057 5077
KK4_496 IgG 1052 2682 3734
KK4_502 Adu 1059 2551 3610
KK4_504 IgG 925 1976 2901

Correlation function and plot

For each gene we compute a Spearman rank correlation (rho) between two vectors that share one entry per cell:

  1. The cell’s raw transcript count for that gene.
  2. The cell’s Euclidean distance (in um) to its nearest plaque of the specified type.

Spearman correlation operates on ranks rather than raw values, so it captures any monotonic relationship (not just linear) and is robust to the extreme zero-inflation typical of spatial transcriptomics count data. We set exact = FALSE because the large number of tied zeros makes an exact permutation p-value computationally infeasible.

A gene that is highly expressed near plaques (short distance) and lowly expressed far away (long distance) will have a negative rho (high expression co-occurs with low distance). Conversely, a gene expressed preferentially far from plaques will have a positive rho.

To make the results more intuitive,so that “plaque-enriched” genes appear at the top of a ranked list rather than the bottom, we multiply rho by -1 to obtain rho_flip. After flipping, positive values mean “enriched near plaques” and negative values mean “enriched far from plaques”.

We apply the Benjamini-Hochberg (BH) procedure to control the false discovery rate, producing an adjusted p-value (padj). A gene is called significant at padj < 0.05.

Final classification. Each gene is labelled as:

  • enriched_near — padj < 0.05 and rho_flip > 0 (upregulated close to plaques)
  • enriched_far — padj < 0.05 and rho_flip < 0 (upregulated far from plaques)
  • ns — not significant

compute_gene_distance_correlation() is centralised in R/setup.R (sourced in setup chunk).

Code 
plot_gene_distance_correlation <- function(
    results,
    title          = "Gene-distance correlation",
    n_top          = 5,
    genes_to_label = NULL
) {
  top_pos <- results |> filter(sig == "enriched_near") |> slice_min(rank, n = n_top)
  top_neg <- results |> filter(sig == "enriched_far")  |> slice_max(rank, n = n_top)
  auto_label <- bind_rows(top_pos, top_neg)

  if (!is.null(genes_to_label)) {
    user_genes <- results |> filter(gene %in% genes_to_label)
    auto_label <- bind_rows(auto_label, user_genes) |> distinct(gene, .keep_all = TRUE)
  }

  colour_map <- c(
    "enriched_near" = "#E07B39",
    "enriched_far"  = "#7B68AE",
    "ns"            = "grey70"
  )

  ggplot(results, aes(x = rank, y = rho_flip)) +
    geom_point(aes(colour = sig), size = 2, alpha = 0.5) +
    geom_hline(yintercept = 0, linetype = "dashed", colour = "grey40") +
    geom_text_repel(
      data          = auto_label,
      aes(label = gene, colour = sig),
      size          = 5,
      fontface      = "italic",
      box.padding   = 0.4,
      max.overlaps  = 20,
      show.legend   = FALSE,
      seed = 23
    ) +
    scale_colour_manual(
      values = colour_map,
      labels = c(
        "enriched_near" = "Enriched near plaque (padj < 0.05)",
        "enriched_far"  = "Enriched far from plaque (padj < 0.05)",
        "ns"            = "Not significant"
      )
    ) +
    labs(
      title  = title,
      x      = "Gene rank",
      y      = "Spearman rho (flipped)",
      colour = NULL
    ) +
    theme_minimal(base_size = 14) +
    theme(
      legend.position = "bottom",
      panel.grid.minor = element_blank(),
      panel.grid.major.x = element_blank())
}

Run analysis — microglia only

Code 
seurat_list <- list(slide1_mg, slide2_mg)

cor_caa   <- compute_gene_distance_correlation(seurat_list, cell_distances_all, distance_col = "dist_caa")
cor_paren <- compute_gene_distance_correlation(seurat_list, cell_distances_all, distance_col = "dist_paren")
cor_any   <- compute_gene_distance_correlation(seurat_list, cell_distances_all, distance_col = "dist_any")

# qs2::qs_save(cor_caa,   "spatial/mg_cor_caa.qs2")
# qs2::qs_save(cor_paren, "spatial/mg_cor_paren.qs2")
# qs2::qs_save(cor_any,   "spatial/mg_cor_any.qs2")
Code 
seurat_list <- list(slide1_mg, slide2_mg)

dist_adu <- cell_distances_all |> filter(Treatment.Group == "Adu")
dist_igg <- cell_distances_all |> filter(Treatment.Group == "IgG")

for (trt in c("adu", "igg")) {
  dist_trt <- if (trt == "adu") dist_adu else dist_igg
  for (dcol in c("dist_caa", "dist_paren", "dist_any")) {
    res    <- compute_gene_distance_correlation(seurat_list, dist_trt, distance_col = dcol)
    suffix <- sub("dist_", "", dcol)
    qs2::qs_save(res, paste0("spatial/mg_cor_", suffix, "_", trt, ".qs2"))
  }
}
Code 
cor_caa   <- qs2::qs_read("spatial/mg_cor_caa.qs2")
cor_paren <- qs2::qs_read("spatial/mg_cor_paren.qs2")
cor_any   <- qs2::qs_read("spatial/mg_cor_any.qs2")

Rank plots

The microglia-specific gene-distance correlation, split by plaque type reveals that microglia activate fundamentally different transcriptional programmes near CAA vs parenchymal plaques.

Parenchymal plaque (right panel)

Near plaque (orange): Classic DAM (Disease-Associated Microglia) signature:

  • Cst7 (rho ~0.45) — the canonical DAM marker, cystatin F, involved in lysosomal processing
  • Ctsb, Cd63 — lysosomal/endosomal genes: microglia near parenchymal plaques are actively phagocytosing amyloid
  • Lyz2 — lysozyme, antimicrobial/inflammatory
  • Lpl — lipoprotein lipase, lipid metabolism rewiring in DAM

Far from plaque (purple): Homeostatic microglia signature:

  • P2ry12, Tmem119 — the two defining homeostatic microglia markers. Their strong negative rho confirms that microglia lose their resting identity as they approach parenchymal plaques.
  • Maf, Ptgs1, Lpcat2 — transcription factor and lipid metabolism genes associated with the surveilling/quiescent state

CAA plaque (left panel)

The effect sizes are much smaller (~0.08 vs ~0.45) and the gene signature is completely different.

Near CAA (orange): This is a lipid/cholesterol metabolism programme, not a phagocytic one. CAA sits in vessel walls, and microglia nearby appear to be responding to vascular lipid dysregulation rather than engulfing amyloid.

  • Glul — glutamine synthetase, typically an astrocyte gene but also involved in glutamate-glutamine cycling at the vasculature
  • Abhd17a — depalmitoylase, involved in protein trafficking
  • Srebf2, Hmgcs1 — cholesterol biosynthesis pathway (SREBP2 is the master regulator, HMGCS1 is its direct target)
  • Acsl6 — fatty acid activation

Far from CAA (purple):

  • Lgals3bp, Ftl1 — galectin-binding and ferritin (iron storage)
  • Abca1, Abcg1 — cholesterol efflux transporters (ABC family).
Code 
p_caa   <- plot_gene_distance_correlation(cor_caa,   title = "CAA plaque")
p_paren <- plot_gene_distance_correlation(cor_paren, title = "Parenchymal plaque")

(p_caa | p_paren) +
  plot_annotation( title = "Plaque proximity gene-distance correlation")+
  plot_layout(guides = "collect") &
  theme(legend.position = "bottom")

For each gene, Spearman’s rank correlation coefficient (ρ) was calculated between per-cell gene expression levels (raw count data) and the minimum distance to the nearest amyloid plaque. To facilitate interpretation, correlation coefficients were sign-inverted such that positive values indicate higher gene expression in cells located closer to plaques, whereas negative values indicate higher expression at increasing distances from plaques. Genes were ranked by the inverted correlation coefficient. Statistical significance was assessed using two-sided Spearman correlation tests followed by Benjamini–Hochberg false discovery rate correction. Genes were classified as plaque-enriched, plaque-distal, or not significant based on adjusted P values. The five most positive and five most negative correlations are labeled
Code 
near_caa   <- cor_caa   |> filter(padj < 0.05, rho_flip > 0) |> pull(gene)
near_paren <- cor_paren |> filter(padj < 0.05, rho_flip > 0) |> pull(gene)

ggVennDiagram(
  list(CAA = near_caa, Parenchymal = near_paren),
  label_alpha = 0, label = "count", shape_id = "201") +
  labs(title = "Enriched-near genes (padj < 0.05)") +
  theme(legend.position = "none")+
  scale_fill_distiller(palette = "RdBu")

Overlap of genes significantly enriched near plaques (padj < 0.05, positive rho_flip) between CAA and parenchymal analyses.
Code 
cat("Common enriched-near genes in both CAA and parenchymal analyses:\n",
  paste(intersect(near_caa, near_paren), collapse = ", "))
Common enriched-near genes in both CAA and parenchymal analyses:
 H2-Eb1, Bhlhe40, Gapdh, Arsb, Vegfa, Cd74, Pkm, Slc2a1, H2-Ab1, Ldlr, Hkdc1, Cd14
Code 
plot_gene_distance_correlation(cor_any, title = "Microglia — distance to nearest plaque (any type)")

Top genes

Code 
show_top <- function(results, label) {
  top <- bind_rows(
    results |> filter(sig == "enriched_near") |> slice_min(rank, n = 10),
    results |> filter(sig == "enriched_far")  |> slice_max(rank, n = 10)
  ) |>
    select(gene, rho, rho_flip, pval, padj, rank, sig) |>
    mutate(across(c(rho, rho_flip), \(x) round(x, 4)),
           across(c(pval, padj), \(x) ifelse(x < 1e-300, "< 1e-300", formatC(x, format = "e", digits = 2))))

  top |>
    kbl(caption = label) |>
    kable_styling("striped", full_width = FALSE)
}

show_top(cor_caa,   "Top genes — microglia, distance to CAA")
Top genes — microglia, distance to CAA
gene rho rho_flip pval padj rank sig
Glul -0.0743 0.0743 1.63e-29 2.79e-27 1 enriched_near
Abhd17a -0.0715 0.0715 1.76e-27 2.02e-25 2 enriched_near
Hmgcs1 -0.0711 0.0711 3.58e-27 3.07e-25 3 enriched_near
Srebf2 -0.0644 0.0644 1.43e-22 6.63e-21 4 enriched_near
Acsl6 -0.0643 0.0643 1.55e-22 6.63e-21 5 enriched_near
Ids -0.0630 0.0630 1.07e-21 3.66e-20 6 enriched_near
Reln -0.0601 0.0601 7.76e-20 2.42e-18 7 enriched_near
H2-Eb1 -0.0553 0.0553 4.63e-17 1.32e-15 8 enriched_near
Gstm1 -0.0551 0.0551 5.89e-17 1.55e-15 9 enriched_near
Pla2g7 -0.0512 0.0512 7.42e-15 1.82e-13 10 enriched_near
Lgals3bp 0.0744 -0.0744 1.26e-29 2.79e-27 343 enriched_far
Ftl1 0.0667 -0.0667 4.26e-24 2.92e-22 342 enriched_far
Abca1 0.0646 -0.0646 1.11e-22 6.37e-21 341 enriched_far
Abcg1 0.0637 -0.0637 4.29e-22 1.63e-20 340 enriched_far
Tmsb4x 0.0472 -0.0472 8.35e-13 1.68e-11 339 enriched_far
C1qa 0.0470 -0.0470 9.49e-13 1.76e-11 338 enriched_far
C1qc 0.0437 -0.0437 3.42e-11 4.69e-10 337 enriched_far
Lyz2 0.0391 -0.0391 3.09e-09 3.54e-08 336 enriched_far
Ctsd 0.0352 -0.0352 9.23e-08 8.62e-07 335 enriched_far
Apod 0.0330 -0.0330 5.36e-07 4.27e-06 334 enriched_far
Code 
show_top(cor_paren, "Top genes — microglia, distance to parenchymal")
Top genes — microglia, distance to parenchymal
gene rho rho_flip pval padj rank sig
Cst7 -0.4642 0.4642 < 1e-300 < 1e-300 1 enriched_near
Ctsb -0.3495 0.3495 < 1e-300 < 1e-300 2 enriched_near
Cd63 -0.3415 0.3415 < 1e-300 < 1e-300 3 enriched_near
Lyz2 -0.3349 0.3349 < 1e-300 < 1e-300 4 enriched_near
Lpl -0.3321 0.3321 < 1e-300 < 1e-300 5 enriched_near
Ftl1 -0.3207 0.3207 < 1e-300 < 1e-300 6 enriched_near
Myo1e -0.3159 0.3159 < 1e-300 < 1e-300 7 enriched_near
Gpnmb -0.2763 0.2763 < 1e-300 < 1e-300 8 enriched_near
Gusb -0.2680 0.2680 < 1e-300 < 1e-300 9 enriched_near
Hexa -0.2651 0.2651 < 1e-300 < 1e-300 10 enriched_near
P2ry12 0.3402 -0.3402 < 1e-300 < 1e-300 343 enriched_far
Tmem119 0.3169 -0.3169 < 1e-300 < 1e-300 342 enriched_far
Lpcat2 0.2798 -0.2798 < 1e-300 < 1e-300 341 enriched_far
Ptgs1 0.2662 -0.2662 < 1e-300 < 1e-300 340 enriched_far
Maf 0.2373 -0.2373 7.85e-292 1.50e-290 339 enriched_far
Bin1 0.2231 -0.2231 2.01e-257 3.45e-256 338 enriched_far
Ids 0.2097 -0.2097 6.44e-227 9.60e-226 337 enriched_far
Cx3cr1 0.1999 -0.1999 6.85e-206 9.40e-205 336 enriched_far
Glul 0.1994 -0.1994 5.67e-205 7.49e-204 335 enriched_far
Selplg 0.1944 -0.1944 8.91e-195 1.09e-193 334 enriched_far
Code 
show_top(cor_any,   "Top genes — microglia, distance to any plaque")
Top genes — microglia, distance to any plaque
gene rho rho_flip pval padj rank sig
Cst7 -0.4719 0.4719 < 1e-300 < 1e-300 1 enriched_near
Ctsb -0.3533 0.3533 < 1e-300 < 1e-300 2 enriched_near
Cd63 -0.3437 0.3437 < 1e-300 < 1e-300 3 enriched_near
Lpl -0.3381 0.3381 < 1e-300 < 1e-300 4 enriched_near
Lyz2 -0.3376 0.3376 < 1e-300 < 1e-300 5 enriched_near
Myo1e -0.3190 0.3190 < 1e-300 < 1e-300 6 enriched_near
Ftl1 -0.3165 0.3165 < 1e-300 < 1e-300 7 enriched_near
Gpnmb -0.2752 0.2752 < 1e-300 < 1e-300 8 enriched_near
Gusb -0.2695 0.2695 < 1e-300 < 1e-300 9 enriched_near
Hexa -0.2660 0.2660 < 1e-300 < 1e-300 10 enriched_near
P2ry12 0.3417 -0.3417 < 1e-300 < 1e-300 343 enriched_far
Tmem119 0.3276 -0.3276 < 1e-300 < 1e-300 342 enriched_far
Lpcat2 0.2828 -0.2828 < 1e-300 < 1e-300 341 enriched_far
Ptgs1 0.2707 -0.2707 < 1e-300 < 1e-300 340 enriched_far
Maf 0.2318 -0.2318 3.22e-278 6.14e-277 339 enriched_far
Bin1 0.2286 -0.2286 2.42e-270 4.37e-269 338 enriched_far
Ids 0.2026 -0.2026 1.28e-211 1.83e-210 337 enriched_far
Cx3cr1 0.1994 -0.1994 6.89e-205 9.45e-204 336 enriched_far
Tgfbr1 0.1942 -0.1942 2.15e-194 2.74e-193 335 enriched_far
Selplg 0.1928 -0.1928 1.39e-191 1.65e-190 334 enriched_far

Treatment-split correlations

Same analysis run separately in Aducanumab and IgG microglia.

Code 
cor_caa_adu   <- qs2::qs_read("spatial/mg_cor_caa_adu.qs2")
cor_caa_igg   <- qs2::qs_read("spatial/mg_cor_caa_igg.qs2")
cor_paren_adu <- qs2::qs_read("spatial/mg_cor_paren_adu.qs2")
cor_paren_igg <- qs2::qs_read("spatial/mg_cor_paren_igg.qs2")
cor_any_adu   <- qs2::qs_read("spatial/mg_cor_any_adu.qs2")
cor_any_igg   <- qs2::qs_read("spatial/mg_cor_any_igg.qs2")

Rank plots by treatment

Code 
plot_gene_distance_correlation(cor_caa_adu, title = "CAA — Aducanumab") |
plot_gene_distance_correlation(cor_caa_igg, title = "CAA — IgG")

Code 
plot_gene_distance_correlation(cor_paren_adu, title = "Parenchymal — Aducanumab") |
plot_gene_distance_correlation(cor_paren_igg, title = "Parenchymal — IgG")

Code 
plot_gene_distance_correlation(cor_any_adu, title = "Any plaque — Aducanumab") |
plot_gene_distance_correlation(cor_any_igg, title = "Any plaque — IgG")

Adu vs IgG scatter

Code 
plot_treatment_scatter <- function(cor_adu, cor_igg, title) {
  df <- inner_join(
    cor_adu |> select(gene, rho_flip_adu = rho_flip, sig_adu = sig),
    cor_igg |> select(gene, rho_flip_igg = rho_flip, sig_igg = sig),
    by = "gene"
  ) |>
    mutate(
      concordance = case_when(
        sig_adu == "enriched_near" & sig_igg == "enriched_near" ~ "both near",
        sig_adu == "enriched_far"  & sig_igg == "enriched_far"  ~ "both far",
        sig_adu != "ns"            & sig_igg != "ns"            ~ "discordant",
        sig_adu != "ns"            | sig_igg != "ns"            ~ "one sig",
        TRUE                                                     ~ "ns"
      )
    )

  colour_map <- c(
    "both near"  = "#E07B39",
    "both far"   = "#7B68AE",
    "discordant" = "#C0392B",
    "one sig"    = "#5DADE2",
    "ns"         = "grey80"
  )

  label_df <- df |> filter(concordance %in% c("both near", "both far", "discordant"))

  ggplot(df, aes(x = rho_flip_igg, y = rho_flip_adu)) +
    geom_point(aes(colour = concordance), size = 1.8, alpha = 0.6) +
    geom_abline(slope = 1, intercept = 0, linetype = "dashed", colour = "grey40") +
    geom_hline(yintercept = 0, colour = "grey70", linewidth = 0.3) +
    geom_vline(xintercept = 0, colour = "grey70", linewidth = 0.3) +
    geom_text_repel(
      data         = label_df,
      aes(label = gene, colour = concordance),
      size         = 3.5,
      fontface     = "italic",
      max.overlaps = 20,
      show.legend  = FALSE
    ) +
    scale_colour_manual(values = colour_map) +
    labs(
      title  = title,
      x      = "rho_flip — IgG",
      y      = "rho_flip — Aducanumab",
      colour = NULL
    ) +
    theme_minimal(base_size = 13) +
    theme(
      legend.position  = "bottom",
      panel.grid.minor = element_blank(),
      panel.grid.major = element_line(colour = "grey92")
    )
}

plot_treatment_scatter(cor_caa_adu,   cor_caa_igg,   "CAA — Adu vs IgG") |
plot_treatment_scatter(cor_paren_adu, cor_paren_igg, "Parenchymal — Adu vs IgG") |
plot_treatment_scatter(cor_any_adu,   cor_any_igg,   "Any plaque — Adu vs IgG")

Candidate gene heatmap

The rank plots above identified two biologically distinct distance-dependent programmes. To compare their effect sizes and significance side by side across both plaque types, we visualise the rho_flip values for the named genes in a single heatmap. Candidates were drawn directly from the labelled genes in the rank plots:

  • DAM / phagocytic programme (near parenchymal plaque): Cst7, Ctsb, Cd63, Lyz2, Lpl — the top-ranked enriched-near genes for parenchymal distance.
  • Vascular-lipid programme (near CAA): Srebf2, Hmgcs1, Acsl6, Glul, Abhd17a — the top-ranked enriched-near genes for CAA distance.
  • Homeostatic programme (far from parenchymal plaque): P2ry12, Tmem119, Maf, Ptgs1, Lpcat2 — the top-ranked enriched-far genes for parenchymal distance.
  • Cholesterol efflux / iron (far from CAA): Abca1, Abcg1, Lgals3bp, Ftl1 — the top-ranked enriched-far genes for CAA distance.

The heatmap is intended to show whether genes that are strongly associated with one plaque type also show a weaker or opposing signal for the other.

Code 
mg_candidate_genes <- c(
  "Cst7", "Ctsb", "Cd63", "Lyz2", "Lpl",
  "Srebf2", "Hmgcs1", "Acsl6", "Glul", "Abhd17a",
  "P2ry12", "Tmem119", "Maf", "Ptgs1", "Lpcat2",
  "Abca1", "Abcg1", "Lgals3bp", "Ftl1"
)

add_sig_and_gene_factor <- function(df) {
  df |>
    filter(gene %in% mg_candidate_genes) |>
    mutate(
      sig_label = case_when(
        padj < 0.001 ~ "***",
        padj < 0.01  ~ "**",
        padj < 0.05  ~ "*",
        TRUE         ~ ""
      ),
      gene = factor(gene, levels = rev(mg_candidate_genes))
    )
}

build_heatmap <- function(df, x_col, title) {
  ggplot(df, aes(x = .data[[x_col]], y = gene, fill = rho_flip)) +
    geom_tile(colour = "white", linewidth = 0.5) +
    geom_text(aes(label = sig_label), size = 4, vjust = 0.75) +
    scale_fill_gradient2(
      low  = "#7B68AE", mid = "white", high = "#E07B39",
      midpoint = 0,
      name = "rho_flip\n(+ = near plaque)"
    ) +
    labs(title = title, x = NULL, y = NULL) +
    theme_minimal(base_size = 13) +
    theme(axis.text.x = element_text(angle = 30, hjust = 1))
}

overall_df <- bind_rows(
  cor_caa   |> mutate(plaque = "CAA"),
  cor_paren |> mutate(plaque = "Parenchymal")
) |>
  add_sig_and_gene_factor() |>
  mutate(plaque = factor(plaque, levels = c("CAA", "Parenchymal")))

by_trt_df <- bind_rows(
  cor_caa_adu   |> mutate(column = "CAA-Adu"),
  cor_caa_igg   |> mutate(column = "CAA-IgG"),
  cor_paren_adu |> mutate(column = "Paren-Adu"),
  cor_paren_igg |> mutate(column = "Paren-IgG")
) |>
  add_sig_and_gene_factor() |>
  mutate(column = factor(column, levels = c("CAA-Adu", "CAA-IgG",
                                            "Paren-Adu", "Paren-IgG")))

p_overall <- build_heatmap(overall_df, "plaque", "Overall (Adu + IgG)")
p_by_trt  <- build_heatmap(by_trt_df,  "column", "Split by treatment")

(p_overall / p_by_trt) +
  plot_annotation(
    title    = "Microglia candidate gene-distance correlations",
    subtitle = "* padj < 0.05  ** padj < 0.01  *** padj < 0.001"
  ) +
  plot_layout(guides = "collect", heights = c(1, 1)) &
  theme(legend.position = "right")

Spearman rho_flip for microglia candidate genes. Top: pooled across treatments, columns = plaque type. Bottom: split by treatment, columns = plaque type x treatment. Positive (orange) = enriched near plaque; negative (purple) = enriched far from plaque. Stars indicate significance (* padj < 0.05, ** padj < 0.01, *** padj < 0.001) computed per column.