From 3094a4dd2cf36a8fa6fc610ad7daf32ad85965c6 Mon Sep 17 00:00:00 2001 From: Sergio-ote Date: Wed, 9 Sep 2026 14:58:04 +0200 Subject: [PATCH 1/2] fix: relocate legends in UMAP and spatial plots to avoid overlapping issues --- docs/rules.md | 6 +- workflow/scripts/annotate_cells.py | 45 +++- workflow/scripts/explore_genes_common.py | 79 ++++--- workflow/scripts/explore_genes_integrated.py | 38 +++- workflow/scripts/explore_genes_sample.py | 41 +++- workflow/scripts/leiden_analysis.py | 35 ++- workflow/scripts/plotting_legends.py | 219 +++++++++++++++++++ workflow/scripts/qc_sweep.py | 36 ++- workflow/scripts/sample_report.py | 28 ++- workflow/scripts/spatial_niches.py | 42 +++- workflow/scripts/subcluster.py | 25 ++- 11 files changed, 497 insertions(+), 97 deletions(-) create mode 100644 workflow/scripts/plotting_legends.py diff --git a/docs/rules.md b/docs/rules.md index 309277f..4c159a8 100644 --- a/docs/rules.md +++ b/docs/rules.md @@ -189,7 +189,7 @@ Generates integrated (Harmony) and unintegrated results, with UMAPs coloured by #### `sample_report` Assembles a compact, multi-page PDF with graphical outputs per sample, ideal to share with internal and external collaborators. -Each page contains the most important data-descriptive information: cluster / annotation / region / niche UMAPs and matching spatial maps, bar plots for absolute and relative composition, and a dot plot with the top-10 expressed markers per cell type. Only annotation sources containing real labels are shown; in runs without a manual TSV annotation, the external or ingest annotation also drives cell-type composition. +Each page contains the most important data-descriptive information: cluster / annotation / region / niche UMAPs and matching spatial maps, bar plots for absolute and relative composition, and a dot plot with the top-10 expressed markers per cell type. Only annotation sources containing real labels are shown; in runs without a manual TSV annotation, the external or ingest annotation also drives cell-type composition. Categorical legends use reserved bands below their plots, with adaptive rows and wrapped long labels, and always remain on the same PDF page as the corresponding panels. --- @@ -202,12 +202,12 @@ Integrated exploration scales are standardized to ease interpretation of the res #### `explore_genes_integrated` Exploration of the integrated samples. -Loads the Harmony-integrated object once and computes AUCell scores for the input gene or gene set. It then writes shared expression ranges and produces integrated UMAP/dotplot PNGs per gene or gene set. +Loads the Harmony-integrated object once and computes AUCell scores for the input gene or gene set. It then writes shared expression ranges and produces integrated UMAP/dotplot PNGs per gene or gene set. Multi-panel categorical legends are kept in dedicated bands within each composite. #### `explore_genes_sample` Exploration of individual samples. -Using the shared expression ranges from the integrated exploration, produces per-sample spatial composites (expression + cell type + region) and stacked dotplot composites for each query entry. +Using the shared expression ranges from the integrated exploration, produces per-sample spatial composites (expression + cell type + region) and stacked dotplot composites for each query entry. Multi-panel categorical legends are kept in dedicated bands within each composite. --- diff --git a/workflow/scripts/annotate_cells.py b/workflow/scripts/annotate_cells.py index 3a394b0..05d10bd 100644 --- a/workflow/scripts/annotate_cells.py +++ b/workflow/scripts/annotate_cells.py @@ -37,6 +37,11 @@ record_active_annotation_columns, ) from external_metadata import optional_input_path, read_cell_metadata +from plotting_legends import ( + category_labels, + move_legend_to_axis, + panel_figure, +) # ── Logging ────────────────────────────────────────────────────────────────── @@ -433,12 +438,21 @@ def generate_annotation_plots(adata, annot_key, label, plots_dir, sample_id, annot_cols_present = active_annotations n_panels = len(annot_cols_present) if n_panels >= 2: - fig, axes = plt.subplots(1, n_panels, figsize=(8 * n_panels, 6), - gridspec_kw={"wspace": 0.5}) - for ax, col in zip(axes, annot_cols_present): + labels_by_panel = [category_labels(adata, col) for col in annot_cols_present] + fig, plot_axes, legend_axes = panel_figure( + n_panels, + labels_by_panel, + panel_width=8, + plot_height=6, + wspace=0.18, + ) + axes = plot_axes[0] + for ax, legend_ax, col in zip(axes, legend_axes, annot_cols_present): + title = col.replace("cell_type_", "") sc.pl.umap(adata, color=col, size=UMAP_POINT_SIZE, frameon=False, - title=col.replace("cell_type_", ""), ax=ax, show=False, - legend_fontsize=6, na_in_legend=False) + title=title, ax=ax, show=False, + legend_fontsize=12, na_in_legend=False) + move_legend_to_axis(ax, legend_ax, title=title) plt.savefig(os.path.join(plots_dir, f"UMAP_comparison_{sample_id}.png"), dpi=DPI, bbox_inches="tight") plt.close() @@ -446,13 +460,24 @@ def generate_annotation_plots(adata, annot_key, label, plots_dir, sample_id, # Side-by-side SPATIAL comparison (only if >1 annotation method) if n_panels >= 2: try: - fig, axes = plt.subplots(1, n_panels, figsize=(8 * n_panels, 6), - gridspec_kw={"wspace": 0.5}) - for ax, col in zip(axes, annot_cols_present): + labels_by_panel = [ + category_labels(adata, col) for col in annot_cols_present + ] + fig, plot_axes, legend_axes = panel_figure( + n_panels, + labels_by_panel, + panel_width=8, + plot_height=6, + wspace=0.18, + ) + axes = plot_axes[0] + for ax, legend_ax, col in zip(axes, legend_axes, annot_cols_present): + title = col.replace("cell_type_", "") sc.pl.spatial(adata, color=col, spot_size=SPATIAL_POINT_SIZE, frameon=False, - title=col.replace("cell_type_", ""), + title=title, library_id=library_id, ax=ax, show=False, - legend_fontsize=6, na_in_legend=False) + legend_fontsize=12, na_in_legend=False) + move_legend_to_axis(ax, legend_ax, title=title) plt.savefig(os.path.join(plots_dir, f"spatial_comparison_{sample_id}.png"), dpi=DPI, bbox_inches="tight") plt.close() diff --git a/workflow/scripts/explore_genes_common.py b/workflow/scripts/explore_genes_common.py index d44f0fb..b0c4787 100644 --- a/workflow/scripts/explore_genes_common.py +++ b/workflow/scripts/explore_genes_common.py @@ -12,6 +12,7 @@ import matplotlib.pyplot as plt from matplotlib.patches import Patch from PIL import Image +from plotting_legends import legend_layout log = logging.getLogger(__name__) @@ -75,24 +76,6 @@ def make_score_adata(adata, score_col): ) -def _compact_legend(ax, title="", per_col=20, fontsize=5): - """Recast an axis legend into multiple narrow columns so a high niche count - does not blow up the panel (~per_col entries per column).""" - leg = ax.get_legend() - if leg is None: - return - handles = (leg.legend_handles if hasattr(leg, "legend_handles") - else getattr(leg, "legendHandles", [])) - labels = [t.get_text() for t in leg.get_texts()] - if not handles: - return - ncol = max(1, (len(labels) + per_col - 1) // per_col) - ax.legend(handles, labels, loc="center left", bbox_to_anchor=(1.0, 0.5), - ncol=ncol, fontsize=fontsize, frameon=False, title=title, - title_fontsize=fontsize + 1, handletextpad=0.3, - columnspacing=0.6, labelspacing=0.25, borderaxespad=0.2) - - def annotate_ct_region_dotplot(annot_key, annotation_colors, region_colors): """Add cell type and region colour bars BETWEEN the row labels and the dot grid, with colour legends below. tick_params(pad) creates space @@ -181,8 +164,12 @@ def annotate_ct_region_dotplot(annot_key, annotation_colors, region_colors): def create_annotation_legend(annot_key, annotation_colors, region_colors, out_path, dpi): - """Create a standalone legend image for cell type + region colour bars. - Legends are stacked vertically so they never overlap.""" + """Create a temporary legend band composited beneath its dotplots. + + Cell-type and region legends use dedicated rows with adaptive height, so long + labels cannot overlap. The caller embeds this image in the same output PNG and + deletes the temporary file; it is never a separate report page. + """ try: from matplotlib.patches import Patch @@ -213,29 +200,41 @@ def create_annotation_legend(annot_key, annotation_colors, region_colors, if not unique_cts and not unique_regs: return - # Two axes stacked vertically — one legend per row, no overlap - n_rows = (1 if unique_cts else 0) + (1 if unique_regs else 0) - fig_leg, axes_leg = plt.subplots(n_rows, 1, figsize=(12, 1.5 * n_rows), - squeeze=False) - row = 0 + groups = [] if unique_cts: - ax = axes_leg[row, 0] - ax.set_axis_off() - ct_patches = [Patch(facecolor=c, label=n) for n, c in unique_cts.items()] - ax.legend(handles=ct_patches, title="Cell type", - loc="center", fontsize=12, title_fontsize=13, - frameon=True, edgecolor="lightgray", - ncol=max(1, len(unique_cts) // 4 + 1)) - row += 1 - + groups.append(("Cell type", unique_cts)) if unique_regs: - ax = axes_leg[row, 0] + groups.append(("Region", unique_regs)) + + layouts = [legend_layout(list(entries), preferred_rows=4) + for _, entries in groups] + heights = [height for _, _, height in layouts] + fig_leg = plt.figure(figsize=(12, sum(heights))) + grid = fig_leg.add_gridspec(len(groups), 1, height_ratios=heights, hspace=0.15) + for row, ((title, entries), (labels, ncol, _)) in enumerate( + zip(groups, layouts) + ): + ax = fig_leg.add_subplot(grid[row, 0]) ax.set_axis_off() - reg_patches = [Patch(facecolor=c, label=n) for n, c in unique_regs.items()] - ax.legend(handles=reg_patches, title="Region", - loc="center", fontsize=12, title_fontsize=13, - frameon=True, edgecolor="lightgray", - ncol=max(1, len(unique_regs) // 3 + 1)) + patches = [ + Patch(facecolor=color, label=label) + for label, color in entries.items() + ] + ax.legend( + handles=patches, + labels=labels, + title=title, + loc="upper center", + fontsize=12, + title_fontsize=13, + frameon=True, + edgecolor="lightgray", + ncol=ncol, + handletextpad=0.35, + columnspacing=0.8, + labelspacing=0.35, + borderaxespad=0.0, + ) fig_leg.savefig(out_path, dpi=dpi, bbox_inches="tight") plt.close(fig_leg) diff --git a/workflow/scripts/explore_genes_integrated.py b/workflow/scripts/explore_genes_integrated.py index 84876cf..d87d30a 100644 --- a/workflow/scripts/explore_genes_integrated.py +++ b/workflow/scripts/explore_genes_integrated.py @@ -28,7 +28,12 @@ from explore_genes_common import ( read_tsv_to_dict, classify_entries, apply_annotation_palette, apply_region_palette, make_score_adata, annotate_ct_region_dotplot, create_annotation_legend, - composite_vertical, _compact_legend, + composite_vertical, +) +from plotting_legends import ( + category_labels, + move_legend_to_axis, + panel_figure, ) # Large cell type × region dotplots can exceed PIL's default pixel limit @@ -165,7 +170,8 @@ def generate_dotplots(adata, var_names, annot_key, has_regions, out_dir, plt.savefig(tmp, dpi=dpi, bbox_inches="tight") tmp_files.append(tmp) - # Standalone annotation legend (before plt.close so gcf() works) + # Temporary annotation-legend band (before plt.close so gcf() works). + # It is composited into the same output PNG below and then deleted. if annotation_colors or region_colors: tmp_leg = tempfile.NamedTemporaryFile(suffix=".png", delete=False).name create_annotation_legend(annot_key, @@ -202,8 +208,20 @@ def generate_umap_composite(adata, color_col, annot_key, has_regions, out_path, multi-column legend for high niche counts.""" has_niche = bool(niche_col and niche_col in adata.obs.columns) n_panels = 2 + (1 if has_regions else 0) + (1 if has_niche else 0) - fig, axes = plt.subplots(1, n_panels, figsize=(8 * n_panels, 7), - gridspec_kw={"wspace": 0.5}) + panel_keys = [None, annot_key] + if has_regions: + panel_keys.append("region_annotation") + if has_niche: + panel_keys.append(niche_col) + labels_by_panel = [category_labels(adata, key) for key in panel_keys] + fig, plot_axes, legend_axes = panel_figure( + n_panels, + labels_by_panel, + panel_width=8, + plot_height=6, + wspace=0.18, + ) + axes = plot_axes[0] try: sc.pl.umap(adata, color=color_col, size=UMAP_POINT_SIZE, frameon=False, @@ -214,9 +232,10 @@ def generate_umap_composite(adata, color_col, annot_key, has_regions, out_path, try: sc.pl.umap(adata, color=annot_key, size=UMAP_POINT_SIZE, frameon=False, - title="Cell types", legend_fontsize=6, + title="Cell types", legend_fontsize=12, na_in_legend=False, ax=axes[1], show=False) + move_legend_to_axis(axes[1], legend_axes[1], title="Cell types") except Exception as e: log.warning(" UMAP celltype failed: %s", e) @@ -224,9 +243,10 @@ def generate_umap_composite(adata, color_col, annot_key, has_regions, out_path, if has_regions: try: sc.pl.umap(adata, color="region_annotation", size=UMAP_POINT_SIZE, frameon=False, - title="Regions", legend_fontsize=6, + title="Regions", legend_fontsize=12, na_in_legend=False, ax=axes[idx], show=False) + move_legend_to_axis(axes[idx], legend_axes[idx], title="Regions") except Exception as e: log.warning(" UMAP region failed: %s", e) idx += 1 @@ -240,9 +260,11 @@ def generate_umap_composite(adata, color_col, annot_key, has_regions, out_path, adata.uns[f"{niche_col}_colors"] = [ niche_palette.get(str(c), "#cccccc") for c in cats] sc.pl.umap(adata, color=niche_col, size=UMAP_POINT_SIZE, frameon=False, - title="Spatial niches", legend_fontsize=5, + title="Spatial niches", legend_fontsize=12, na_in_legend=False, ax=axes[idx], show=False) - _compact_legend(axes[idx], title="Niche") + move_legend_to_axis( + axes[idx], legend_axes[idx], title="Spatial niches", fontsize=12 + ) except Exception as e: log.warning(" UMAP niche failed: %s", e) idx += 1 diff --git a/workflow/scripts/explore_genes_sample.py b/workflow/scripts/explore_genes_sample.py index f6de536..8b31700 100644 --- a/workflow/scripts/explore_genes_sample.py +++ b/workflow/scripts/explore_genes_sample.py @@ -33,7 +33,12 @@ from explore_genes_common import ( read_tsv_to_dict, classify_entries, apply_annotation_palette, apply_region_palette, make_score_adata, annotate_ct_region_dotplot, create_annotation_legend, - composite_vertical, _compact_legend, + composite_vertical, +) +from plotting_legends import ( + category_labels, + move_legend_to_axis, + panel_figure, ) # Large cell type × region dotplots can exceed PIL's default pixel limit @@ -72,8 +77,21 @@ def generate_spatial_composite(adata, color_col, annot_key, has_regions, multi-column legend so a high niche count does not distort the composite.""" has_niche = bool(niche_col and niche_col in adata.obs.columns) n_panels = 2 + (1 if has_regions else 0) + (1 if has_niche else 0) - fig, axes = plt.subplots(1, n_panels, figsize=(8 * n_panels, 7), - gridspec_kw={"wspace": 0.5}) + panel_keys = [None, annot_key] + if has_regions: + panel_keys.append("region_annotation") + if has_niche: + panel_keys.append(niche_col) + labels_by_panel = [category_labels(adata, key) for key in panel_keys] + fig, plot_axes, legend_axes = panel_figure( + n_panels, + labels_by_panel, + panel_width=8, + plot_height=6, + wspace=0.18, + top_margin=0.8, + ) + axes = plot_axes[0] fig.suptitle(sample_id, fontsize=16, fontweight="bold", y=1.02) # (1) Expression / AUCell @@ -90,8 +108,9 @@ def generate_spatial_composite(adata, color_col, annot_key, has_regions, try: sc.pl.spatial(adata, color=annot_key, spot_size=SPATIAL_POINT_SIZE, frameon=False, title="Cell types", library_id=library_id, - legend_fontsize=6, na_in_legend=False, + legend_fontsize=12, na_in_legend=False, ax=axes[1], show=False) + move_legend_to_axis(axes[1], legend_axes[1], title="Cell types") except Exception as e: log.warning(" Spatial celltype failed: %s", e) axes[1].set_title("Cell types (failed)") @@ -103,8 +122,9 @@ def generate_spatial_composite(adata, color_col, annot_key, has_regions, sc.pl.spatial(adata, color="region_annotation", spot_size=SPATIAL_POINT_SIZE, frameon=False, title="Regions", library_id=library_id, - legend_fontsize=6, na_in_legend=False, + legend_fontsize=12, na_in_legend=False, ax=axes[idx], show=False) + move_legend_to_axis(axes[idx], legend_axes[idx], title="Regions") except Exception as e: log.warning(" Spatial region failed: %s", e) axes[idx].set_title("Regions (failed)") @@ -122,9 +142,11 @@ def generate_spatial_composite(adata, color_col, annot_key, has_regions, niche_palette.get(str(c), "#cccccc") for c in cats] sc.pl.spatial(adata, color=niche_col, spot_size=SPATIAL_POINT_SIZE, frameon=False, title="Spatial niches", library_id=library_id, - legend_fontsize=5, na_in_legend=False, + legend_fontsize=12, na_in_legend=False, ax=axes[idx], show=False) - _compact_legend(axes[idx], title="Niche") + move_legend_to_axis( + axes[idx], legend_axes[idx], title="Spatial niches", fontsize=12 + ) except Exception as e: log.warning(" Spatial niche failed: %s", e) axes[idx].set_title("Spatial niches (failed)") @@ -197,8 +219,9 @@ def generate_dotplot_composite(adata, var_names, annot_key, has_regions, plt.savefig(tmp, dpi=dpi, bbox_inches="tight") tmp_files.append(tmp) - # Standalone annotation legend (before plt.close so gcf() can - # still read the dotplot's ytick labels) + # Temporary annotation-legend band (before plt.close so gcf() can + # still read the dotplot's ytick labels). It is composited into the + # same output PNG below and then deleted. if annotation_colors or region_colors: tmp_leg = tempfile.NamedTemporaryFile(suffix=".png", delete=False).name create_annotation_legend(annot_key, diff --git a/workflow/scripts/leiden_analysis.py b/workflow/scripts/leiden_analysis.py index b745123..eee2ee5 100644 --- a/workflow/scripts/leiden_analysis.py +++ b/workflow/scripts/leiden_analysis.py @@ -23,6 +23,17 @@ from sklearn.metrics import silhouette_samples, silhouette_score from sklearn.preprocessing import StandardScaler +try: + _here = os.path.dirname(os.path.abspath(__file__)) +except NameError: # very old Snakemake + _here = os.getcwd() +sys.path.insert(0, _here) +from plotting_legends import ( + category_labels, + grid_figure, + move_legend_to_axis, +) + log_handlers = [logging.StreamHandler(sys.stderr)] if hasattr(snakemake, "log"): @@ -69,10 +80,21 @@ def read_tsv_to_dict(tsv_path): def split_umap(adata, split_by, ncol=None, nrow=None, **kwargs): categories = adata.obs[split_by].cat.categories - ncol = ncol or len(categories) + ncol = ncol or min(4, len(categories)) nrow = nrow or int(np.ceil(len(categories) / ncol)) - fig, axs = plt.subplots(nrow, ncol, figsize=(5 * ncol, 4 * nrow)) - axs = np.atleast_1d(axs).flatten() + color_key = kwargs.get("color") + if isinstance(color_key, (list, tuple)): + color_key = color_key[0] if color_key else None + labels = category_labels(adata, color_key) + fig, axes, legend_ax = grid_figure( + nrow, + ncol, + labels, + cell_width=5, + cell_height=4, + max_legend_columns=min(8, 2 * ncol), + ) + axs = axes.ravel() for i, cat in enumerate(categories): sc.pl.umap( adata[adata.obs[split_by] == cat], @@ -80,7 +102,12 @@ def split_umap(adata, split_by, ncol=None, nrow=None, **kwargs): ) for j in range(i + 1, len(axs)): axs[j].set_visible(False) - plt.tight_layout() + move_legend_to_axis( + axs[:len(categories)], + legend_ax, + title=str(color_key).replace("_", " ").title() if color_key else None, + max_columns=min(8, 2 * ncol), + ) try: diff --git a/workflow/scripts/plotting_legends.py b/workflow/scripts/plotting_legends.py new file mode 100644 index 0000000..a00650d --- /dev/null +++ b/workflow/scripts/plotting_legends.py @@ -0,0 +1,219 @@ +"""Helpers for keeping categorical legends inside multi-panel figures.""" + +from __future__ import annotations + +import math +import textwrap + +import matplotlib.pyplot as plt +import numpy as np + + +DEFAULT_LEGEND_FONTSIZE = 12 + + +def category_labels(adata, obs_key): + """Return displayed category labels in their plotting order.""" + if not obs_key or obs_key not in adata.obs.columns: + return [] + values = adata.obs[obs_key] + categories = getattr(values.dtype, "categories", None) + if categories is not None: + labels = categories + else: + labels = values.dropna().unique() + return [str(label) for label in labels] + + +def _wrap_label(label, width=24): + """Wrap only the displayed legend text; underlying category names are unchanged.""" + return textwrap.fill( + str(label), + width=width, + break_long_words=True, + break_on_hyphens=False, + ) + + +def legend_layout(labels, wrap_width=24, preferred_rows=5, max_columns=4, + fontsize=DEFAULT_LEGEND_FONTSIZE): + """Return wrapped labels, column count, and estimated height in inches.""" + wrapped = [_wrap_label(label, wrap_width) for label in labels] + if not wrapped: + return wrapped, 1, 0.0 + + ncol = min(max_columns, max(1, math.ceil(len(wrapped) / preferred_rows))) + rows = math.ceil(len(wrapped) / ncol) + line_counts = [label.count("\n") + 1 for label in wrapped] + + # Matplotlib fills legends column-first. Estimate each displayed row using the + # tallest wrapped entry that can land on it, then reserve a title/margin band. + row_lines = [] + for row in range(rows): + row_lines.append( + max( + (line_counts[col * rows + row] + for col in range(ncol) + if col * rows + row < len(line_counts)), + default=1, + ) + ) + # Reserve physical space in proportion to the rendered font size. The + # multiplier includes Matplotlib's default line spacing plus row padding. + line_height = 0.30 * (fontsize / 10) + title_height = 0.55 * (fontsize / 10) + height = title_height + sum(line_height * lines for lines in row_lines) + return wrapped, ncol, max(0.8, height) + + +def legend_height(labels_by_panel, minimum=0.8): + """Height needed by the largest panel legend.""" + heights = [legend_layout(labels)[2] for labels in labels_by_panel] + return max([minimum, *heights]) + + +def panel_figure(n_panels, labels_by_panel, plot_rows=1, panel_width=8.0, + plot_height=6.0, wspace=0.18, hspace=0.25, + top_margin=0.0): + """Create plot axes plus a dedicated legend row in the same figure. + + The returned legend axes are part of the figure itself. Consequently, a PDF + writer always keeps each legend on the same page as its plots. + """ + legend_h = legend_height(labels_by_panel) + fig = plt.figure( + figsize=(panel_width * n_panels, + plot_height * plot_rows + legend_h + top_margin) + ) + grid = fig.add_gridspec( + plot_rows + 1, + n_panels, + height_ratios=[plot_height] * plot_rows + [legend_h], + wspace=wspace, + hspace=hspace, + ) + plot_axes = np.asarray( + [[fig.add_subplot(grid[row, col]) for col in range(n_panels)] + for row in range(plot_rows)] + ) + legend_axes = [fig.add_subplot(grid[plot_rows, col]) for col in range(n_panels)] + for ax in legend_axes: + ax.set_axis_off() + return fig, plot_axes, legend_axes + + +def grid_figure(nrows, ncols, legend_labels, cell_width=5.0, cell_height=4.0, + wspace=0.18, hspace=0.25, top_margin=0.0, + max_legend_columns=8): + """Create a plot grid with one shared legend band on the same figure.""" + _, _, legend_h = legend_layout( + legend_labels, + max_columns=max_legend_columns, + ) + legend_h = max(0.8, legend_h) + fig = plt.figure( + figsize=(cell_width * ncols, + cell_height * nrows + legend_h + top_margin) + ) + grid = fig.add_gridspec( + nrows + 1, + ncols, + height_ratios=[cell_height] * nrows + [legend_h], + wspace=wspace, + hspace=hspace, + ) + axes = np.asarray( + [[fig.add_subplot(grid[row, col]) for col in range(ncols)] + for row in range(nrows)] + ) + legend_ax = fig.add_subplot(grid[nrows, :]) + legend_ax.set_axis_off() + return fig, axes, legend_ax + + +def draw_legend_on_axis( + legend_ax, + handles, + labels, + title=None, + fontsize=DEFAULT_LEGEND_FONTSIZE, + wrap_width=24, + preferred_rows=5, + max_columns=4, +): + """Draw handles in a reserved legend axis using adaptive rows and wrapping.""" + legend_ax.set_axis_off() + if not handles: + return + wrapped, ncol, _ = legend_layout( + labels, + wrap_width=wrap_width, + preferred_rows=preferred_rows, + max_columns=max_columns, + fontsize=fontsize, + ) + legend_ax.legend( + handles, + wrapped, + title=title, + loc="upper center", + bbox_to_anchor=(0.5, 1.0), + ncol=ncol, + fontsize=fontsize, + title_fontsize=fontsize + 1, + frameon=False, + handletextpad=0.35, + columnspacing=0.8, + labelspacing=0.35, + borderaxespad=0.0, + ) + + +def move_legend_to_axis( + source_axes, + legend_ax, + title=None, + fontsize=DEFAULT_LEGEND_FONTSIZE, + wrap_width=24, + preferred_rows=5, + max_columns=4, +): + """Move unique entries from Scanpy legends into a dedicated in-figure axis. + + Any duplicate legends on the supplied source axes are removed. ``source_axes`` + can therefore contain matching UMAP and spatial panels while only one legend + is rendered below them. + """ + if not isinstance(source_axes, (list, tuple, np.ndarray)): + source_axes = [source_axes] + + handles = [] + labels = [] + seen_labels = set() + for ax in source_axes: + legend = ax.get_legend() + if legend is None: + continue + axis_handles = ( + list(legend.legend_handles) + if hasattr(legend, "legend_handles") + else list(getattr(legend, "legendHandles", [])) + ) + axis_labels = [text.get_text() for text in legend.get_texts()] + for handle, label in zip(axis_handles, axis_labels): + if label not in seen_labels: + handles.append(handle) + labels.append(label) + seen_labels.add(label) + legend.remove() + + draw_legend_on_axis( + legend_ax, + handles, + labels, + title=title, + fontsize=fontsize, + wrap_width=wrap_width, + preferred_rows=preferred_rows, + max_columns=max_columns, + ) diff --git a/workflow/scripts/qc_sweep.py b/workflow/scripts/qc_sweep.py index b394a12..414cd3a 100644 --- a/workflow/scripts/qc_sweep.py +++ b/workflow/scripts/qc_sweep.py @@ -21,6 +21,13 @@ import pandas as pd import scanpy as sc +try: + _here = os.path.dirname(os.path.abspath(__file__)) +except NameError: # very old Snakemake + _here = os.getcwd() +sys.path.insert(0, _here) +from plotting_legends import grid_figure, move_legend_to_axis + # ── Logging ────────────────────────────────────────────────────────────────── log_handlers = [logging.StreamHandler(sys.stderr)] if hasattr(snakemake, "log"): @@ -202,9 +209,20 @@ def make_landing(label, vals, op, ths, out_path): nrows = 2 if ingest_vals is not None else 1 ncols = 1 + len(ths) - fig, axes = plt.subplots(nrows, ncols, figsize=(5.6 * ncols, 5.0 * nrows), - squeeze=False, - gridspec_kw={"wspace": 0.6, "hspace": 0.35}) + legend_labels = list(ingest_palette) + if ths: + legend_labels.append("removed") + fig, axes, legend_ax = grid_figure( + nrows, + ncols, + legend_labels, + cell_width=5.6, + cell_height=5.0, + wspace=0.18, + hspace=0.3, + top_margin=0.6, + max_legend_columns=min(8, 2 * ncols), + ) def draw(ax, color, title, **extra): try: @@ -221,7 +239,7 @@ def draw(ax, color, title, **extra): adata.uns["_qc_colors"] = [ingest_palette[c] for c in adata.obs["_qc"].cat.categories] draw(axes[r][0], "_qc", "all cells (ingest)", - legend_fontsize=5, na_in_legend=False) + legend_fontsize=12, na_in_legend=False) for j, t in enumerate(ths): mask = (vals < t) if op == "below" else (vals > t) title = f"{label} {t}\n{int(mask.sum()):,} removed" @@ -234,7 +252,7 @@ def draw(ax, color, title, **extra): adata.obs["_qc"] = pd.Categorical(tmp, categories=cats) adata.uns["_qc_colors"] = [ingest_palette[c] for c in cats] draw(axes[r][j + 1], "_qc", title, na_color="#e8e8e8", - na_in_legend=False, legend_fontsize=5) + na_in_legend=False, legend_fontsize=12) r += 1 draw(axes[r][0], None, "all cells") @@ -250,6 +268,14 @@ def draw(ax, color, title, **extra): adata.uns["_qc_colors"] = ["#e41a1c"] draw(axes[r][j + 1], "_qc", title, na_color="#e8e8e8", na_in_legend=False) + move_legend_to_axis( + axes.ravel(), + legend_ax, + title="Cell status / ingest cell type", + fontsize=12, + max_columns=min(8, 2 * ncols), + ) + if "_qc" in adata.obs: del adata.obs["_qc"] adata.uns.pop("_qc_colors", None) diff --git a/workflow/scripts/sample_report.py b/workflow/scripts/sample_report.py index d8ec6ff..6f1a129 100644 --- a/workflow/scripts/sample_report.py +++ b/workflow/scripts/sample_report.py @@ -33,6 +33,11 @@ from composition_barplots import (draw_stacked_composition, find_niche_column, build_niche_palette) from annotation_utils import active_annotation_columns, has_meaningful_annotation +from plotting_legends import ( + category_labels, + move_legend_to_axis, + panel_figure, +) # ── Logging ────────────────────────────────────────────────────────────────── @@ -128,24 +133,37 @@ def _build_page1(adata, sample_id, keys, library_id, has_regions=False): panels.append((niche, "Spatial niche")) panels = [(k, t) for k, t in panels if k] n = len(panels) - fig, axes = plt.subplots(2, n, figsize=(8 * n, 14), - gridspec_kw={"wspace": 0.75, "hspace": 0.35}, - squeeze=False) + labels_by_panel = [category_labels(adata, key) for key, _ in panels] + fig, axes, legend_axes = panel_figure( + n, + labels_by_panel, + plot_rows=2, + panel_width=8, + plot_height=6, + wspace=0.18, + hspace=0.28, + top_margin=1.0, + ) fig.suptitle(f"Sample: {sample_id}", fontsize=18, fontweight="bold", y=0.99) for j, (key, title) in enumerate(panels): try: sc.pl.umap(adata, color=key, size=UMAP_POINT_SIZE, frameon=False, title=title, - legend_fontsize=6, na_in_legend=False, + legend_fontsize=12, na_in_legend=False, ax=axes[0, j], show=False) except Exception as e: log.warning(" UMAP %s failed: %s", title, e); axes[0, j].set_title(f"{title} (failed)") try: sc.pl.spatial(adata, color=key, spot_size=SPATIAL_POINT_SIZE, frameon=False, title=title, library_id=library_id, - legend_fontsize=6, na_in_legend=False, + legend_fontsize=12, na_in_legend=False, ax=axes[1, j], show=False) except Exception as e: log.warning(" Spatial %s failed: %s", title, e); axes[1, j].set_title(f"{title} (failed)") + move_legend_to_axis( + [axes[0, j], axes[1, j]], + legend_axes[j], + title=title, + ) return fig diff --git a/workflow/scripts/spatial_niches.py b/workflow/scripts/spatial_niches.py index 430e015..2266471 100644 --- a/workflow/scripts/spatial_niches.py +++ b/workflow/scripts/spatial_niches.py @@ -35,6 +35,7 @@ # Shared composition-barplot helpers (scripts/ is on sys.path for script: rules) sys.path.insert(0, os.path.dirname(os.path.abspath(__file__))) from composition_barplots import composition_pair, build_niche_palette # noqa: E402 +from plotting_legends import draw_legend_on_axis, grid_figure # noqa: E402 # ── Logging ────────────────────────────────────────────────────────────────── @@ -159,8 +160,14 @@ def _spatial_scatter(adata, sample_key, color_key, out_path, dpi, color_map): n = len(samples) ncol = min(4, n) nrow = int(np.ceil(n / ncol)) - fig, axes = plt.subplots(nrow, ncol, figsize=(4 * ncol, 4 * nrow), - squeeze=False) + fig, axes, legend_ax = grid_figure( + nrow, + ncol, + [str(c) for c in cats], + cell_width=4, + cell_height=4, + max_legend_columns=min(8, 2 * ncol), + ) for ax in axes.ravel(): ax.set_axis_off() for i, s in enumerate(samples): @@ -182,9 +189,14 @@ def _spatial_scatter(adata, sample_key, color_key, out_path, dpi, color_map): ax.set_xticks([]); ax.set_yticks([]) handles = [plt.Line2D([0], [0], marker="o", linestyle="", markersize=5, color=colour[c], label=str(c)) for c in cats] - fig.legend(handles=handles, loc="center left", bbox_to_anchor=(1.0, 0.5), - frameon=False, title="Spatial niche", fontsize=7, title_fontsize=8) - fig.tight_layout() + draw_legend_on_axis( + legend_ax, + handles, + [str(c) for c in cats], + title="Spatial niche", + fontsize=12, + max_columns=min(8, 2 * ncol), + ) fig.savefig(out_path, dpi=dpi, bbox_inches="tight") plt.close(fig) @@ -198,7 +210,15 @@ def _spatial_scatter_highlight(adata, sample_key, color_key, target, out_path, d n = len(samples) ncol = min(4, n) nrow = int(np.ceil(n / ncol)) - fig, axes = plt.subplots(nrow, ncol, figsize=(4 * ncol, 4 * nrow), squeeze=False) + highlight_labels = [f"niche {target}", "other"] + fig, axes, legend_ax = grid_figure( + nrow, + ncol, + highlight_labels, + cell_width=4, + cell_height=4, + max_legend_columns=min(8, 2 * ncol), + ) for ax in axes.ravel(): ax.set_axis_off() for i, s in enumerate(samples): @@ -226,9 +246,13 @@ def _spatial_scatter_highlight(adata, sample_key, color_key, target, out_path, d color=hi_color, label=f"niche {target}"), plt.Line2D([0], [0], marker="o", linestyle="", markersize=5, color=bg_color, label="other")] - fig.legend(handles=handles, loc="center left", bbox_to_anchor=(1.0, 0.5), - frameon=False, fontsize=7) - fig.tight_layout() + draw_legend_on_axis( + legend_ax, + handles, + highlight_labels, + fontsize=12, + max_columns=min(8, 2 * ncol), + ) fig.savefig(out_path, dpi=dpi, bbox_inches="tight") plt.close(fig) diff --git a/workflow/scripts/subcluster.py b/workflow/scripts/subcluster.py index ea01cd4..0c3a32f 100644 --- a/workflow/scripts/subcluster.py +++ b/workflow/scripts/subcluster.py @@ -35,6 +35,7 @@ _here = os.getcwd() sys.path.insert(0, _here) from composition_barplots import save_stacked_composition, find_niche_column, build_niche_palette +from plotting_legends import category_labels, grid_figure, move_legend_to_axis # ── Logging ────────────────────────────────────────────────────────────────── @@ -66,10 +67,21 @@ def _resolution_range(rmin, rmax, step): def split_umap(adata, split_by, ncol=None, nrow=None, **kwargs): categories = adata.obs[split_by].cat.categories - ncol = ncol or len(categories) + ncol = ncol or min(4, len(categories)) nrow = nrow or int(np.ceil(len(categories) / ncol)) - fig, axs = plt.subplots(nrow, ncol, figsize=(5 * ncol, 4 * nrow)) - axs = np.atleast_1d(axs).flatten() + color_key = kwargs.get("color") + if isinstance(color_key, (list, tuple)): + color_key = color_key[0] if color_key else None + labels = category_labels(adata, color_key) + fig, axes, legend_ax = grid_figure( + nrow, + ncol, + labels, + cell_width=5, + cell_height=4, + max_legend_columns=min(8, 2 * ncol), + ) + axs = axes.ravel() for i, cat in enumerate(categories): sc.pl.umap( adata[adata.obs[split_by] == cat], @@ -77,7 +89,12 @@ def split_umap(adata, split_by, ncol=None, nrow=None, **kwargs): ) for j in range(i + 1, len(axs)): axs[j].set_visible(False) - plt.tight_layout() + move_legend_to_axis( + axs[:len(categories)], + legend_ax, + title=str(color_key).replace("_", " ").title() if color_key else None, + max_columns=min(8, 2 * ncol), + ) def recalculate_qc_metrics(adata): From f61c6c7e692f46c101646a929d8224cb4d25d2c5 Mon Sep 17 00:00:00 2001 From: Sergio-ote Date: Wed, 9 Sep 2026 15:43:38 +0200 Subject: [PATCH 2/2] fix: clustree and silhouette display --- workflow/scripts/leiden_analysis.py | 74 ++++++++++++++++++++++++++--- 1 file changed, 67 insertions(+), 7 deletions(-) diff --git a/workflow/scripts/leiden_analysis.py b/workflow/scripts/leiden_analysis.py index eee2ee5..4ca9332 100644 --- a/workflow/scripts/leiden_analysis.py +++ b/workflow/scripts/leiden_analysis.py @@ -16,10 +16,12 @@ import matplotlib matplotlib.use("Agg") import matplotlib.pyplot as plt +from matplotlib.colors import is_color_like, to_hex import numpy as np import pandas as pd from pyclustree import clustree import scanpy as sc +from scanpy.plotting import palettes as scanpy_palettes from sklearn.metrics import silhouette_samples, silhouette_score from sklearn.preprocessing import StandardScaler @@ -67,6 +69,43 @@ SPATIAL_POINT_SIZE = float(getattr(snakemake.params, "spatial_point_size", 20)) +def _resolution_range(rmin, rmax, step): + """Return the same inclusive resolution ladder used during preprocessing.""" + n = round((rmax - rmin) / step) + return [round(rmin + i * step, 10) for i in range(n + 1)] + + +def _ensure_scanpy_palette(adata, obs_key): + """Return category-to-colour mapping used by subsequent Scanpy plots.""" + categories = list(adata.obs[obs_key].cat.categories) + color_key = f"{obs_key}_colors" + stored = list(adata.uns.get(color_key, [])) + valid = ( + len(stored) == len(categories) + and all(is_color_like(color) for color in stored) + ) + + if valid: + colors = [to_hex(color, keep_alpha=True) for color in stored] + else: + n_categories = len(categories) + prop_cycle = plt.rcParams["axes.prop_cycle"].by_key()["color"] + if len(prop_cycle) >= n_categories: + colors = prop_cycle[:n_categories] + elif n_categories <= 20: + colors = scanpy_palettes.default_20[:n_categories] + elif n_categories <= 28: + colors = scanpy_palettes.default_28[:n_categories] + elif n_categories <= len(scanpy_palettes.default_102): + colors = scanpy_palettes.default_102[:n_categories] + else: + colors = ["grey"] * n_categories + colors = [to_hex(color) for color in colors] + + adata.uns[color_key] = colors + return {str(category): color for category, color in zip(categories, colors)} + + def read_tsv_to_dict(tsv_path): with open(tsv_path) as fh: reader = csv.DictReader(fh, delimiter="\t") @@ -154,6 +193,7 @@ def _set_palette(adata, obs_key, colors_cfg): adata.uns[f"{obs_key}_colors"] = [cd.get(str(c), "#cccccc") for c in cats] _set_palette(adata, "leiden", ANNOTATION_COLORS) + leiden_palette = _ensure_scanpy_palette(adata, "leiden") # Apply custom region palette if configured if REGION_COLORS and "region_annotation" in adata.obs.columns: @@ -172,20 +212,32 @@ def _set_palette(adata, obs_key, colors_cfg): # ── Clustering evaluation (silhouette) ─────────────────────────────── log.info("Silhouette analysis …") embedding_scaled = StandardScaler().fit_transform(adata.obsm["X_pca"]) - labels = adata.obs["leiden"].astype(int) + labels = adata.obs["leiden"].astype(str).to_numpy() sil_values = silhouette_samples(embedding_scaled, labels) sil_avg = silhouette_score(embedding_scaled, labels) plt.figure(figsize=(12, 6)) - unique_labels = np.unique(labels) - colours = plt.cm.tab10(np.linspace(0, 1, len(unique_labels))) + unique_labels = [str(label) for label in adata.obs["leiden"].cat.categories] x_pos = 0 - for label, colour in zip(unique_labels, colours): + for label in unique_labels: vals = np.sort(sil_values[labels == label]) + if not len(vals): + continue xs = np.arange(len(vals)) + x_pos - plt.bar(xs, vals, color=colour, label=f"Cluster {label}", width=1, alpha=0.7) + plt.bar( + xs, + vals, + color=leiden_palette[label], + label=f"Cluster {label}", + width=1, + ) x_pos = xs[-1] + 5 - plt.axhline(y=sil_avg, color="red", linestyle="--", label="Average Silhouette") + plt.axhline( + y=sil_avg, + color="black", + linestyle="--", + label="Average Silhouette", + ) plt.ylabel("Silhouette Score") plt.xlabel("Cluster and Cell Index") plt.title(f"Silhouette (resolution {resolution}) – avg {sil_avg:.3f}") @@ -197,14 +249,22 @@ def _set_palette(adata, obs_key, colors_cfg): plt.close() # Clustree (uses all pre-computed leiden columns) - resolutions = np.arange(RES_SCAN_MIN, RES_SCAN_MAX, RES_SCAN_STEP).round(1) + resolutions = _resolution_range(RES_SCAN_MIN, RES_SCAN_MAX, RES_SCAN_STEP) leiden_keys = [f"leiden_{str(r).replace('.', '_')}" for r in resolutions] available_keys = [k for k in leiden_keys if k in adata.obs.columns] + missing_keys = [k for k in leiden_keys if k not in adata.obs.columns] + if missing_keys: + log.warning("Clustree is missing screened columns: %s", missing_keys) if len(available_keys) >= 2: fig = clustree(adata, available_keys, edge_weight_threshold=0.00, show_fraction=True) fig.savefig(os.path.join(clust_eval_dir, f"clustree_{sample_id}.png"), dpi=DPI, bbox_inches="tight") plt.close() + else: + log.warning( + "Clustree requires at least two available screened resolutions; found %s", + available_keys, + ) # ── QC plots ───────────────────────────────────────────────────────── log.info("QC plots …")