Open/filled data points based whether inside or outside exon

Author: micknudsen

src/gcnvplot/core.py

@@ -373,7 +373,10 @@ def build_background(args: argparse.Namespace) -> int:
 
 
 def rows_for_plot(
-    args: argparse.Namespace, region: Interval, background_summary: BackgroundSummary
+    args: argparse.Namespace,
+    region: Interval,
+    background_summary: BackgroundSummary,
+    transcript: TranscriptAnnotation | None = None,
 ) -> list[dict[str, object]]:
     """Join one sample with the background for plotting."""
     counts = parse_read_counts(args.read_counts)
@@ -413,6 +416,9 @@ def rows_for_plot(
             "bg_upper_norm": bg_upper,
             "background_n": int(bg["N"]),
         }
+        row["overlaps_exon"] = transcript is None or any(
+            overlaps(interval[1], interval[2], exon.start, exon.end) for exon in transcript.exons
+        )
         row["signal"] = log2_ratio(sample_value, bg_median, 0.01)
         row["signal_lower"] = log2_ratio(bg_lower, bg_median, 0.01)
         row["signal_upper"] = log2_ratio(bg_upper, bg_median, 0.01)
@@ -505,7 +511,9 @@ def y_for(value: float) -> float:
         ".grid { stroke: #d8dee4; stroke-width: 1; }",
         ".ribbon { fill: #c8d2df; opacity: 0.85; }",
         ".sample { fill: none; stroke: #127c78; stroke-width: 2.2; }",
-        ".point { fill: #127c78; stroke: white; stroke-width: 1.2; }",
+        ".point { stroke-width: 1.6; }",
+        ".point-filled { fill: #127c78; stroke: white; }",
+        ".point-open { fill: white; stroke: #127c78; }",
         ".baseline { stroke: #5b6472; stroke-width: 2; stroke-dasharray: 6 5; }",
         ".highlight-band { fill: #f59e0b; opacity: 0.28; stroke: #b45309; stroke-width: 1.1; }",
         ".panel { fill: #ffffff; stroke: #d8dee4; stroke-width: 1.2; }",
@@ -564,8 +572,9 @@ def y_for(value: float) -> float:
             f"SIGNAL={fmt(float(row['signal']))}\n"
             f"BACKGROUND_N={row['background_n']}"
         )
+        point_class = "point point-filled" if bool(row["overlaps_exon"]) else "point point-open"
         elements.append(
-            f'<circle class="point" cx="{x:.2f}" cy="{y:.2f}" r="4.8">'
+            f'<circle class="{point_class}" cx="{x:.2f}" cy="{y:.2f}" r="4.8">'
             f"<title>{html.escape(title_text)}</title></circle>"
         )
 
@@ -770,12 +779,10 @@ def add_arrow(x_pos: float) -> None:
     elements[3] = f'<rect class="panel" x="{panel_x}" y="{panel_y}" width="{panel_width}" height="{panel_h}"/>'
     elements.extend(
         [
-            f'<rect x="{panel_x + 16}" y="{legend_y - 13}" width="15" height="15" fill="#c8d2df" opacity="0.85"/>',
-            f'<text class="panel-text" x="{panel_x + 40}" y="{legend_y}">background</text>',
-            f'<line x1="{panel_x + 16}" y1="{legend_y + 23}" x2="{panel_x + 34}" y2="{legend_y + 23}" stroke="#5b6472" stroke-width="2" stroke-dasharray="6 5"/>',
-            f'<text class="panel-text" x="{panel_x + 40}" y="{legend_y + 27}">expected baseline</text>',
-            f'<circle class="point" cx="{panel_x + 24}" cy="{legend_y + 48}" r="5"/>',
-            f'<text class="panel-text" x="{panel_x + 40}" y="{legend_y + 52}">sample signal</text>',
+            f'<circle class="point point-filled" cx="{panel_x + 24}" cy="{legend_y}" r="5"/>',
+            f'<text class="panel-text" x="{panel_x + 40}" y="{legend_y + 4}">Overlaps exon</text>',
+            f'<circle class="point point-open" cx="{panel_x + 24}" cy="{legend_y + 24}" r="5"/>',
+            f'<text class="panel-text" x="{panel_x + 40}" y="{legend_y + 28}">Outside exon</text>',
             "</svg>",
         ]
     )
@@ -793,7 +800,7 @@ def plot_sample(args: argparse.Namespace) -> int:
         region = (transcript.contig, transcript.start, transcript.end)
 
     background_summary = parse_background(args.background)
-    rows = rows_for_plot(args, region, background_summary)
+    rows = rows_for_plot(args, region, background_summary, transcript=transcript)
     if not rows:
         raise SystemExit("No intervals with background statistics overlap the selected region.")
 

tests/test_cli.py

@@ -76,6 +76,16 @@ def svg_highlight_bands(svg: str) -> list[tuple[float, float, float, float]]:
     return bands
 
 
+def svg_point_class_by_interval(svg: str) -> dict[str, str]:
+    """Return point classes keyed by interval title."""
+    classes: dict[str, str] = {}
+    pattern = r'<circle class="([^"]*point[^"]*)"[^>]*><title>([^<]+)</title></circle>'
+    for point_class, title_text in re.findall(pattern, svg):
+        interval = title_text.split("&#10;", 1)[0]
+        classes[interval] = point_class
+    return classes
+
+
 def test_main_version(capsys: pytest.CaptureFixture[str]) -> None:
     with pytest.raises(SystemExit) as exc_info:
         cli.main(["--version"])
@@ -341,8 +351,9 @@ def test_plot_transcript_writes_exon_track_and_gene_name(
     assert "Exons" in svg
     assert "Highlight" in svg
     assert "log2(sample/background)" in svg
+    assert "Overlaps exon" in svg
+    assert "Outside exon" in svg
     assert ">2</text>" in svg
-    assert ">background</text>" in svg
     assert svg.count('class="panel-separator"') == 2
     assert svg.index("Sample") < svg.index("Gene") < svg.index("Transcript") < svg.index("Region") < svg.index("Highlight") < svg.index("Exons")
     highlight_bands = svg_highlight_bands(svg)
@@ -358,6 +369,78 @@ def test_plot_transcript_writes_exon_track_and_gene_name(
     )
 
 
+def test_plot_transcript_uses_open_points_for_non_exonic_intervals(
+    tmp_path: Path, capsys: pytest.CaptureFixture[str]
+) -> None:
+    background_path = tmp_path / "background.tsv"
+    background_path.write_text(
+        "\n".join(
+            [
+                "# normalization=median-of-ratios",
+                "# baseline=median-positive-count",
+                "# samples=2",
+                "# lower_percentile=5",
+                "# upper_percentile=95",
+                "CONTIG\tSTART\tEND\tBASELINE_MEDIAN\tN\tBG_NORM_MEAN\tBG_NORM_MEDIAN\tBG_NORM_SD\tBG_NORM_P5\tBG_NORM_P95",
+                "chr1\t100\t120\t10\t2\t10\t10\t0\t10\t10",
+                "chr1\t121\t179\t20\t2\t20\t20\t0\t20\t20",
+                "chr1\t180\t220\t30\t2\t30\t30\t0\t30\t30",
+            ]
+        )
+        + "\n",
+        encoding="utf-8",
+    )
+    sample_path = tmp_path / "sample.tsv"
+    write_read_counts(
+        sample_path,
+        [
+            ("chr1", 100, 120, 12),
+            ("chr1", 121, 179, 18),
+            ("chr1", 180, 220, 33),
+        ],
+    )
+    gtf_path = tmp_path / "transcript.gtf.gz"
+    write_gtf_gz(
+        gtf_path,
+        [
+            'chr1\tsource\ttranscript\t100\t220\t.\t+\t.\tgene_id "GENE1"; gene_name "MYGENE"; transcript_id "TX1";',
+            'chr1\tsource\texon\t100\t120\t.\t+\t.\tgene_id "GENE1"; gene_name "MYGENE"; transcript_id "TX1"; exon_number "1";',
+            'chr1\tsource\texon\t180\t220\t.\t+\t.\tgene_id "GENE1"; gene_name "MYGENE"; transcript_id "TX1"; exon_number "2";',
+        ],
+    )
+    output_path = tmp_path / "plot.svg"
+
+    assert (
+        cli.main(
+            [
+                "plot",
+                "--read-counts",
+                str(sample_path),
+                "--background",
+                str(background_path),
+                "--transcript",
+                "TX1",
+                "--gtf",
+                str(gtf_path),
+                "--output",
+                str(output_path),
+            ]
+        )
+        == 0
+    )
+
+    svg = output_path.read_text(encoding="utf-8")
+    point_classes = svg_point_class_by_interval(svg)
+    assert point_classes["chr1:100-120"] == "point point-filled"
+    assert point_classes["chr1:121-179"] == "point point-open"
+    assert point_classes["chr1:180-220"] == "point point-filled"
+
+    assert capsys.readouterr().out == (
+        "Plotted intervals: 3\n"
+        f"Wrote: {output_path}\n"
+    )
+
+
 def test_plot_reverse_transcript_arrows_do_not_touch_exons(
     tmp_path: Path, capsys: pytest.CaptureFixture[str]
 ) -> None: