Skip to content

Commit dc7d8f0

Browse files
authored
Cohorts: run_cohort — N subjects per run (#8)
* feat(cohort): CohortSubjectResult + CohortResult core accessors * feat(cohort): combined to_dataframe/to_csv/to_fasta with union columns + unique headers * feat(cohort): subject-id + counts resolution helpers * feat(cohort): Experiment.run_cohort happy path (clone-per-subject loop) * test(cohort): auto-ids, counts override/zero, mutual-exclusion, cartridge-hash, mid-loop failure safety * test(cohort): mixed novel/plain subjects + per-subject refdata/validate + e2e * fix(cohort): CohortResult.records returns independent dict copies; test allele-weights exclusion * docs(cohort): run_cohort guide section + remove from limitations * fix(cohort): apply with_metadata per subject (+collision guard); reject mapping/bytes counts; empty-export schema
1 parent 8b0c237 commit dc7d8f0

4 files changed

Lines changed: 676 additions & 2 deletions

File tree

site_docs/guides/genotype.md

Lines changed: 56 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -409,6 +409,62 @@ the allele frequencies are cartridge- or uniform-sourced — supplying an
409409
`model_id` / `source` / `version`, the plane's `model_checksum`, per-segment gene
410410
counts, the novel-allele count, and `source_field` (`"DataConfig.genotype_priors"`).
411411

412+
## Cohorts
413+
414+
A single genotype models one subject. To simulate a **cohort** — many subjects,
415+
each with their own genotype — use `run_cohort`. It runs the single-subject path
416+
once per subject and collects the results:
417+
418+
```python
419+
import GenAIRR as ga
420+
import GenAIRR.data as gdata
421+
from GenAIRR.genotype import Genotype
422+
423+
cfg = gdata.HUMAN_IGH_OGRDB
424+
# one sampled genotype per donor
425+
donors = [Genotype.sample(cfg, seed=s, subject_id=f"donor_{s}") for s in range(5)]
426+
427+
cohort = ga.Experiment.on(cfg).recombine().run_cohort(
428+
donors, n_per_subject=200, seed=0, expose_provenance=True)
429+
430+
cohort.subject_ids # ['donor_0', ..., 'donor_4']
431+
len(cohort) # 1000 total records
432+
cohort.result_for("donor_2") # that donor's SimulationResult
433+
cohort.to_csv("cohort.csv") # combined, subject-tagged, unique sequence_id
434+
```
435+
436+
`run_cohort` returns a `CohortResult`:
437+
438+
- `.subject_ids` / `.genotypes` / `.results` — per-subject, in input order.
439+
- `.result_for(sid)` / `.refdata_for(sid)` — one subject's `SimulationResult` and
440+
the reference it ran against (preserved per subject, so `validate_records` and
441+
novel-allele truth calls stay correct even when subjects differ).
442+
- `.records` — a fresh, subject-tagged, `sequence_id`-namespaced concatenation
443+
(each `sequence_id` is `"{subject_id}_{...}"`, so combined AIRR/FASTA export
444+
never collides).
445+
- `.to_dataframe()` / `.to_csv()` / `.to_fasta()` — combined export over a stable
446+
union of columns.
447+
448+
**Record counts.** `n_per_subject` applies to all subjects; pass `counts` (a
449+
parallel list, same length as the genotypes) to vary per-subject repertoire
450+
sizes. A count of `0` is allowed — that subject appears in the cohort with zero
451+
records.
452+
453+
```python
454+
cohort = ga.Experiment.on(cfg).recombine().run_cohort(
455+
donors, counts=[500, 200, 0, 1000, 300], seed=0)
456+
```
457+
458+
**Subject IDs.** Taken from each genotype's `subject_id`; if none are set they are
459+
auto-assigned `subject_0..N-1`. Mixed (some set, some not) or duplicate IDs raise.
460+
461+
**Determinism.** Each subject gets an independent sub-seed derived from `seed`, so
462+
a cohort is fully reproducible and subjects are independent.
463+
464+
`run_cohort` is mutually exclusive with `with_genotype`, `restrict_alleles`, and
465+
`recombine(*_allele_weights=...)` (the genotype owns allele expression), and — in
466+
this release — is not combined with `receptor_revision` or clonal forks.
467+
412468
## Novel / private alleles
413469

414470
Individuals carry germline alleles that aren't in any reference — *private* or
@@ -616,8 +672,6 @@ deletion calls.
616672

617673
The genotype foundation is deliberately scoped. Deferred to later work:
618674

619-
- **Cohorts** — many subjects, each with their own genotype, in one run
620-
(`with_genotype` is single-subject; `result.genotypes` is a one-element list).
621675
- **External loaders** — importing genotypes from VDJbase / TIgGER / IgDiscover /
622676
partis output.
623677
- **Same-haplotype receptor revision**`receptor_revision` with a genotype is

src/GenAIRR/cohort.py

Lines changed: 157 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,157 @@
1+
"""Cohort orchestration results — N subjects each with their own diploid
2+
genotype, produced by :meth:`GenAIRR.Experiment.run_cohort`.
3+
4+
``CohortResult`` is built from per-subject ``CohortSubjectResult`` entries and
5+
exposes per-subject access (``result_for`` / ``refdata_for``) plus combined
6+
export. It explicitly stores each subject's refdata because ``SimulationResult``
7+
does not preserve it (needed for ``validate_records`` and novel-allele subjects).
8+
See ``.private/specs/2026-06-17-genotype-cohorts-design.md``.
9+
"""
10+
from __future__ import annotations
11+
12+
from dataclasses import dataclass
13+
from typing import Any, Dict, List
14+
15+
16+
def _resolve_subject_ids(raw_ids: List[Any]) -> List[str]:
17+
"""Resolve per-subject IDs: all-None -> subject_0..N-1; all-present ->
18+
require unique after str() normalization; mixed -> raise."""
19+
none_count = sum(1 for i in raw_ids if i is None)
20+
n = len(raw_ids)
21+
if none_count == n:
22+
return [f"subject_{i}" for i in range(n)]
23+
if none_count != 0:
24+
raise ValueError(
25+
"run_cohort: some genotypes have a subject_id and others don't; set "
26+
"subject_id on all genotypes or none")
27+
ids = [str(i) for i in raw_ids]
28+
if len(ids) != len(set(ids)):
29+
raise ValueError("run_cohort: duplicate subject_id after normalization")
30+
return ids
31+
32+
33+
def _check_count(c, where: str) -> int:
34+
if isinstance(c, bool) or not isinstance(c, int) or c < 0:
35+
raise ValueError(f"{where} must be an int >= 0, got {c!r}")
36+
return c
37+
38+
39+
def _resolve_counts(n_genotypes: int, n_per_subject, counts) -> List[int]:
40+
"""Resolve per-subject record counts. ``counts`` (a parallel sequence)
41+
overrides ``n_per_subject`` when supplied; entries are validated as
42+
non-bool ints >= 0."""
43+
_check_count(n_per_subject, "n_per_subject")
44+
if counts is None:
45+
return [int(n_per_subject)] * n_genotypes
46+
from collections.abc import Mapping
47+
if isinstance(counts, (str, bytes, bytearray, Mapping)):
48+
raise ValueError(
49+
"run_cohort: counts must be a parallel list/tuple of ints, not a "
50+
f"{type(counts).__name__}")
51+
counts = list(counts)
52+
if len(counts) != n_genotypes:
53+
raise ValueError(
54+
f"run_cohort: counts length {len(counts)} != number of genotypes "
55+
f"{n_genotypes}")
56+
return [_check_count(c, f"counts[{i}]") for i, c in enumerate(counts)]
57+
58+
59+
@dataclass(frozen=True)
60+
class CohortSubjectResult:
61+
"""One subject's slice of a cohort run."""
62+
63+
subject_id: str
64+
genotype: Any # GenAIRR.genotype.Genotype (snapshot)
65+
result: Any # GenAIRR.result.SimulationResult
66+
refdata: Any # the refdata this subject ran against (base or effective)
67+
seed: int # derived per-subject sub-seed
68+
count: int # records requested for this subject
69+
70+
71+
class CohortResult:
72+
"""Combined result of a :meth:`Experiment.run_cohort` run.
73+
74+
``subjects`` is the single source of truth; ``subject_ids`` / ``genotypes`` /
75+
``results`` are derived so parallel lists can't drift."""
76+
77+
def __init__(self, subjects: List[CohortSubjectResult]):
78+
self._subjects = list(subjects)
79+
self._by_id = {s.subject_id: s for s in self._subjects}
80+
81+
@property
82+
def subjects(self) -> List[CohortSubjectResult]:
83+
return list(self._subjects)
84+
85+
@property
86+
def subject_ids(self) -> List[str]:
87+
return [s.subject_id for s in self._subjects]
88+
89+
@property
90+
def genotypes(self) -> List[Any]:
91+
return [s.genotype for s in self._subjects]
92+
93+
@property
94+
def results(self) -> List[Any]:
95+
return [s.result for s in self._subjects]
96+
97+
def result_for(self, subject_id: str):
98+
return self._by_id[subject_id].result
99+
100+
def refdata_for(self, subject_id: str):
101+
return self._by_id[subject_id].refdata
102+
103+
@property
104+
def records(self) -> List[Dict[str, Any]]:
105+
"""A fresh concatenated list of every subject's record dicts (each
106+
already ``subject_id``-tagged and ``sequence_id``-namespaced). Each dict
107+
is a shallow copy, so neither the returned list nor mutating a record in
108+
it affects the per-subject results — edit those via ``result_for``."""
109+
out: List[Dict[str, Any]] = []
110+
for s in self._subjects:
111+
out.extend(dict(rec) for rec in s.result.records)
112+
return out
113+
114+
def __len__(self) -> int:
115+
return sum(len(s.result) for s in self._subjects)
116+
117+
def __repr__(self) -> str:
118+
return f"<CohortResult subjects={len(self._subjects)} records={len(self)}>"
119+
120+
# ── combined export ─────────────────────────────────────────────
121+
def to_dataframe(self, *, airr_strict: bool = False):
122+
"""Combined DataFrame over every subject's records. Column set is the
123+
stable union across subjects (pandas fills missing keys with NaN)."""
124+
import pandas as pd
125+
126+
from .result import _DEFAULT_COLUMN_ORDER
127+
128+
records = self.records
129+
if not records:
130+
# Empty cohort: preserve the default AIRR schema plus cohort-owned
131+
# columns so an empty export still has a usable header (parity with
132+
# SimulationResult.to_dataframe on an empty result).
133+
cols = list(_DEFAULT_COLUMN_ORDER)
134+
for extra in ("subject_id", "haplotype"):
135+
if extra not in cols:
136+
cols.append(extra)
137+
return pd.DataFrame(columns=cols)
138+
if airr_strict:
139+
from .result import _to_airr_strict
140+
records = [_to_airr_strict(r) for r in records]
141+
return pd.DataFrame(records)
142+
143+
def to_csv(self, path: str, *, airr_strict: bool = False) -> None:
144+
"""Write the combined records as CSV (union columns guaranteed by the
145+
DataFrame)."""
146+
self.to_dataframe(airr_strict=airr_strict).to_csv(path, index=False)
147+
148+
def to_fasta(self, path: str) -> None:
149+
"""Write one FASTA record per concatenated record. Header is exactly
150+
``>{sequence_id}`` (the namespaced id); body is the record's
151+
``sequence``. Unlike ``SimulationResult.to_fasta`` (whose headers come
152+
from the enumerate index), this keeps cohort headers globally unique."""
153+
with open(path, "w", encoding="utf-8") as fh:
154+
for rec in self.records:
155+
sid = rec.get("sequence_id", "")
156+
seq = rec.get("sequence", "")
157+
fh.write(f">{sid}\n{seq}\n")

src/GenAIRR/experiment.py

Lines changed: 129 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -2948,6 +2948,135 @@ def _build_simulator(
29482948
allow_curatable_refdata=allow_curatable_refdata,
29492949
)
29502950

2951+
def run_cohort(
2952+
self,
2953+
genotypes,
2954+
*,
2955+
n_per_subject: int = 1,
2956+
seed: int = 0,
2957+
counts=None,
2958+
strict: bool = False,
2959+
expose_provenance: bool = False,
2960+
validate_records: bool = False,
2961+
allow_curatable_refdata: Optional[bool] = None,
2962+
) -> "CohortResult":
2963+
"""Run a cohort: N subjects, each with its own diploid genotype, in one
2964+
call. A Python loop around the single-subject genotype path — each
2965+
subject is compiled and run independently, records are tagged with
2966+
``subject_id`` and given a namespaced ``sequence_id``, and the per-subject
2967+
``SimulationResult`` (with its own refdata) is collected into a
2968+
:class:`~GenAIRR.cohort.CohortResult`.
2969+
2970+
``n_per_subject`` applies to every subject; ``counts`` (a parallel
2971+
sequence, same length as ``genotypes``) overrides it per subject and may
2972+
contain ``0`` (that subject appears with zero records). Subject IDs are
2973+
taken from each genotype, or auto-assigned ``subject_0..N-1`` when all are
2974+
unset; mixed/duplicate IDs raise. Per-subject sub-seeds are derived
2975+
deterministically from ``seed``.
2976+
2977+
Mutually exclusive with :meth:`with_genotype`, :meth:`restrict_alleles`,
2978+
and ``recombine(*_allele_weights=...)`` (the genotype owns allele
2979+
expression). Not supported with :meth:`receptor_revision` or clonal forks
2980+
in this release.
2981+
"""
2982+
import copy as _copy
2983+
import random as _random
2984+
2985+
from .cohort import (
2986+
CohortResult,
2987+
CohortSubjectResult,
2988+
_resolve_counts,
2989+
_resolve_subject_ids,
2990+
)
2991+
from .genotype import Genotype
2992+
from .result import SimulationResult
2993+
2994+
gts = list(genotypes)
2995+
if not gts:
2996+
raise ValueError("run_cohort: genotypes must be a non-empty sequence")
2997+
for g in gts:
2998+
if not isinstance(g, Genotype):
2999+
raise TypeError(
3000+
f"run_cohort: every element must be a Genotype, got "
3001+
f"{type(g).__name__}")
3002+
3003+
# Mutual exclusions — mirror with_genotype / compile.
3004+
if self._genotype is not None:
3005+
raise ValueError(
3006+
"run_cohort() and with_genotype() are mutually exclusive")
3007+
if any(v is not None for v in self._locks.values()):
3008+
raise ValueError(
3009+
"run_cohort() and restrict_alleles() are mutually exclusive")
3010+
if self._user_allele_weights_set:
3011+
raise ValueError(
3012+
"run_cohort() and recombine(*_allele_weights=...) are mutually "
3013+
"exclusive: the genotype owns allele expression")
3014+
if any(isinstance(s, _ReceptorRevisionStep) for s in self._steps):
3015+
raise ValueError(
3016+
"run_cohort() is not supported with receptor_revision() in this "
3017+
"release (the revision pass is not haplotype-aware)")
3018+
if self._has_clonal_fork():
3019+
raise ValueError(
3020+
"run_cohort() is not supported together with expand_clones() / "
3021+
"clonal_lineage() / clonal_repertoire() in this release")
3022+
# with_metadata() is applied per subject below, but it must not overwrite
3023+
# cohort-owned columns.
3024+
_cohort_owned = {"subject_id", "sequence_id", "haplotype"}
3025+
_md_collision = _cohort_owned & set(self._metadata)
3026+
if _md_collision:
3027+
raise ValueError(
3028+
f"run_cohort: with_metadata keys {sorted(_md_collision)} are "
3029+
f"cohort-owned (reserved); rename them")
3030+
3031+
# Cartridge-hash check (same as with_genotype).
3032+
live_hash = self._refdata.content_hash()
3033+
for g in gts:
3034+
if g._source_hash != live_hash:
3035+
raise ValueError(
3036+
"run_cohort: a genotype was built against a different cartridge "
3037+
f"(content hash {g._source_hash!r} != experiment {live_hash!r})")
3038+
3039+
resolved_counts = _resolve_counts(len(gts), n_per_subject, counts)
3040+
subject_ids = _resolve_subject_ids([g.subject_id for g in gts])
3041+
base_rng = _random.Random(seed)
3042+
3043+
subjects = []
3044+
for g, sid, count in zip(gts, subject_ids, resolved_counts):
3045+
sub_seed = base_rng.getrandbits(63)
3046+
snap = g._snapshot()
3047+
snap.subject_id = sid
3048+
# Clone the uncompiled experiment; never mutate self.
3049+
exp_i = _copy.copy(self)
3050+
exp_i._genotype = snap
3051+
compiled = exp_i.compile(allow_curatable_refdata=allow_curatable_refdata)
3052+
refdata_i = compiled.refdata
3053+
if count == 0:
3054+
# run() rejects n < 1; build an empty result and stamp the
3055+
# genotype manually (run_records would otherwise do it).
3056+
res = SimulationResult.from_outcomes(
3057+
[], refdata_i, expose_provenance=expose_provenance)
3058+
res._genotypes = [snap]
3059+
else:
3060+
res = compiled.run_records(
3061+
n=count, seed=sub_seed, strict=strict,
3062+
expose_provenance=expose_provenance,
3063+
validate_records=validate_records)
3064+
# Apply with_metadata() per subject (parity with run_records), then
3065+
# namespace sequence_id so combined export never collides. Metadata
3066+
# is stamped first; the cohort-owned sequence_id rewrite wins (and a
3067+
# collision was already rejected up front).
3068+
if self._metadata:
3069+
for rec in res.records:
3070+
for key, value in self._metadata.items():
3071+
rec[key] = value
3072+
for rec in res.records:
3073+
rec["sequence_id"] = f"{sid}_{rec.get('sequence_id', '')}"
3074+
subjects.append(CohortSubjectResult(
3075+
subject_id=sid, genotype=snap, result=res, refdata=refdata_i,
3076+
seed=sub_seed, count=count))
3077+
3078+
return CohortResult(subjects)
3079+
29513080
def run_records(
29523081
self,
29533082
*,

0 commit comments

Comments
 (0)