Skip to content
Closed
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 3 additions & 0 deletions .github/workflows/nextflow-test.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -5,6 +5,9 @@ on:
pull_request:
branches:
- 'main'
push:
branches:
- 'v2'
workflow_dispatch:

env:
Expand Down
5 changes: 4 additions & 1 deletion .gitignore
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,7 @@ runs/*
results*
work
out-*
out/*

## nf-test
.nf-test/
Expand Down Expand Up @@ -34,4 +35,6 @@ notebooks/*
tmp

## vscode
.vscode/*
.vscode/*

.e2e_test_tmp/
Empty file added assets/NO_FILE
Empty file.
78 changes: 78 additions & 0 deletions bin/bulk_to_export.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,78 @@
#!/usr/bin/env python3
"""
Bridge: BULK_TO_EXPORT

Builds a per-cell export table (the format REPERTOIRE / MASTER_SUMMARY consume) from
pseudobulk clonotype data, so those reports can run in the VDJ-only route where there is
no GEX Seurat object.

The repertoire report derives per-clone cell counts by COUNTING ROWS grouped by
(sample, clone_id) — i.e. it treats each row as one cell. So each clonotype with
duplicate_count = N is expanded into N rows (one synthetic "cell" each). Without a GEX
object there is no cell-type annotation, so `annot` is a constant ("Unannotated").

Input : concatenated canonical clonotype table (junction_aa, v_call, j_call,
duplicate_count, ..., sample) — e.g. ANNOTATE_FROM_CONCAT's concat_cdr3_sorted.
Optional: samplesheet CSV mapping sample -> patient_id (else patient_id = sample).

Output: export_cells.tsv with columns REPERTOIRE/MASTER_SUMMARY resolve:
cell_id, sample, patient_id, clone_id, clone_size, annot, has_tcr, paired_tcr
"""

import sys
import argparse
import pandas as pd


def main():
ap = argparse.ArgumentParser(description="BULK_TO_EXPORT bridge")
ap.add_argument("concat_cdr3", help="concatenated canonical clonotype TSV (has 'sample')")
ap.add_argument("--samplesheet", default=None, help="CSV with sample,patient_id (optional)")
ap.add_argument("--out", default="export_cells.tsv")
args = ap.parse_args()

df = pd.read_csv(args.concat_cdr3, sep="\t", low_memory=False)

for col in ("junction_aa", "duplicate_count", "sample"):
if col not in df.columns:
sys.exit(f"[BULK_TO_EXPORT] input missing required column '{col}'. Have: {list(df.columns)}")

# Clone identity = CDR3b + V gene (falls back to junction_aa alone if v_call absent)
if "v_call" in df.columns:
df["clone_id"] = df["junction_aa"].astype(str) + "_" + df["v_call"].astype(str)
else:
df["clone_id"] = df["junction_aa"].astype(str)

df["clone_size"] = pd.to_numeric(df["duplicate_count"], errors="coerce").fillna(0).astype(int)
df = df[df["clone_size"] > 0].copy()
if df.empty:
sys.exit("[BULK_TO_EXPORT] no clonotypes with positive counts.")

# patient_id: from samplesheet if given, else = sample
if args.samplesheet:
ss = pd.read_csv(args.samplesheet)
pcol = next((c for c in ("patient_id", "patient") if c in ss.columns), None)
if pcol and "sample" in ss.columns:
pmap = dict(zip(ss["sample"].astype(str), ss[pcol].astype(str)))
df["patient_id"] = df["sample"].astype(str).map(pmap).fillna(df["sample"].astype(str))
else:
df["patient_id"] = df["sample"].astype(str)
else:
df["patient_id"] = df["sample"].astype(str)

# Expand each clonotype into `clone_size` per-cell rows (one row = one cell)
expanded = df.loc[df.index.repeat(df["clone_size"])].reset_index(drop=True)
expanded["annot"] = "Unannotated" # no GEX → no cell-type label
expanded["has_tcr"] = "TRUE"
expanded["paired_tcr"] = "FALSE"
expanded["cell_id"] = ["cell_%d" % i for i in range(len(expanded))]

out_cols = ["cell_id", "sample", "patient_id", "clone_id", "clone_size",
"annot", "has_tcr", "paired_tcr"]
expanded[out_cols].to_csv(args.out, sep="\t", index=False)
print(f"[BULK_TO_EXPORT] {df['clone_id'].nunique()} clonotypes across "
f"{df['sample'].nunique()} sample(s) -> {len(expanded)} cell rows in {args.out}")


if __name__ == "__main__":
main()
91 changes: 91 additions & 0 deletions bin/cluster_to_sc.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,91 @@
#!/usr/bin/env python3
"""
Bridge 2: CLUSTER_TO_SC (optional)
Maps TCRtoolkit bulk cluster assignments (CDR3b-level) onto single cells
using CDR3b as the join key, producing per-cell TSVs compatible with
SCRATCH-TCR's Consensus Clustering module inputs.

Usage:
cluster_to_sc.py <export_cells.tsv> <giana_clusters> <gliph_clusters>

Pass 'NO_FILE' for any cluster file that is not available.

Outputs (written only when input is not NO_FILE):
giana_export_cells.tsv — per-cell with giana_cluster column
gliph2_export_cells.tsv — per-cell with gliph2_cluster column
"""

import sys
import pandas as pd


NO_FILE_SENTINEL = 'NO_FILE'


def extract_beta_cdr3(ctaa):
"""Extract beta chain CDR3 from CTaa field (format: 'alpha;beta')."""
if pd.isna(ctaa):
return None
parts = str(ctaa).split(';')
return parts[1].strip().upper() if len(parts) > 1 else None


def map_clusters(cells: pd.DataFrame, cluster_file: str, cluster_col_name: str,
out_file: str) -> None:
"""
Join CDR3b-level cluster assignments from a TCRtoolkit output file onto
per-cell data and write the result.

The cluster file is expected to have:
- First column: CDR3b amino acid sequence (may have a comment header)
- 'cluster' column (or second column if 'cluster' is absent)
"""
try:
cl = pd.read_csv(cluster_file, sep='\t', comment='#')
except Exception as e:
print(f"[Bridge 2] Warning: could not read {cluster_file}: {e}")
return

cdr3_col = cl.columns[0]
cluster_col = 'cluster' if 'cluster' in cl.columns else cl.columns[1]

cl_map = (
cl[[cdr3_col, cluster_col]]
.rename(columns={cdr3_col: 'CDR3b', cluster_col: cluster_col_name})
.drop_duplicates('CDR3b')
)
cl_map['CDR3b'] = cl_map['CDR3b'].str.upper()

merged = cells.merge(cl_map, on='CDR3b', how='left')
merged.to_csv(out_file, sep='\t', index=False)
assigned = merged[cluster_col_name].notna().sum()
print(f"[Bridge 2] {out_file}: {assigned}/{len(merged)} cells assigned a {cluster_col_name}")


def main():
if len(sys.argv) < 4:
print("Usage: cluster_to_sc.py <export_cells.tsv> <giana_clusters> <gliph_clusters>")
sys.exit(1)

export_cells_file = sys.argv[1]
giana_file = sys.argv[2]
gliph_file = sys.argv[3]

cells = pd.read_csv(export_cells_file, sep='\t', low_memory=False)

if 'CTaa' in cells.columns:
cells['CDR3b'] = cells['CTaa'].apply(extract_beta_cdr3)
elif 'junction_aa' in cells.columns:
cells['CDR3b'] = cells['junction_aa'].str.upper()
else:
raise ValueError("export_cells.tsv must contain 'CTaa' or 'junction_aa'")

if giana_file and giana_file != NO_FILE_SENTINEL:
map_clusters(cells, giana_file, 'giana_cluster', 'giana_export_cells.tsv')

if gliph_file and gliph_file != NO_FILE_SENTINEL:
map_clusters(cells, gliph_file, 'gliph2_cluster', 'gliph2_export_cells.tsv')


if __name__ == '__main__':
main()
Loading
Loading