P2 validation: oracle test and large-scale raw X/Y check
- tests/oracle/test_p2_raw_xy.py: IDAT-assembled raw X/Y and controls equal the arrays in oracle GTCs - validate_p1_large.py: phase-P2 comparison against vendor GTCs with diagnostics (channel swap, bead-count disagreements) - cloud trial: also runs `simurg ingest` on all samples (time/memory), disk 250 GB
This commit is contained in:
1 parent
4c0ea5e04f
commit
7b83155d56
4 files changed
+141
-5
No files matched your search
@@ -16,7 +16,7 @@ REGION=${REGION:-us-east-1}
|
||||
TYPE=${TYPE:-m7g.4xlarge}
|
||||
MAX_PRICE=${MAX_PRICE:-0.30} # USD/hour spot cap (on-demand is 0.6528)
|
||||
MAX_MINUTES=${MAX_MINUTES:-360} # hard cap: 6 hours
|
||||
DISK_GB=${DISK_GB:-200}
|
||||
DISK_GB=${DISK_GB:-250}
|
||||
BUDGET_USD=${BUDGET_USD:-10}
|
||||
TAG=simurg-trial
|
||||
ACCOUNT=$(aws sts get-caller-identity --query Account --output text)
|
||||
|
||||
@@ -33,7 +33,7 @@ echo "=== start $(date -u +%FT%TZ) on $(uname -m), $(nproc) cpus ==="
|
||||
|
||||
export DEBIAN_FRONTEND=noninteractive
|
||||
apt-get update -q
|
||||
apt-get install -y -q pigz git python3-venv curl unzip
|
||||
apt-get install -y -q pigz git python3-venv curl unzip time
|
||||
snap install aws-cli --classic
|
||||
S3="/snap/bin/aws s3"
|
||||
( while sleep 300; do upload_log; done ) &
|
||||
@@ -60,4 +60,14 @@ echo "=== validation $(date -u +%T) ==="
|
||||
--egt $D/HumanOmni2.5-4v1-Multi_B.egt --csv $D/HumanOmni2.5-4v1_B.csv --out /work/summary.json \
|
||||
--workers "$(nproc)" || echo "validation reported differences (see summary.json)"
|
||||
$S3 cp /work/summary.json "s3://$BUCKET/$RUN/summary.json" --only-show-errors
|
||||
|
||||
echo "=== simurg ingest at scale $(date -u +%T) ==="
|
||||
/usr/bin/time -v .venv/bin/simurg ingest --bpm $D/HumanOmni2.5-4v1-Multi_B.bpm --idats $D/data \
|
||||
--out /work/ingest --threads "$(nproc)" --quiet 2> /work/ingest_time.txt || echo "ingest exit code $?"
|
||||
tail -25 /work/ingest_time.txt
|
||||
for f in ingest_summary.tsv run_manifest.json; do
|
||||
$S3 cp "/work/ingest/$f" "s3://$BUCKET/$RUN/ingest_$f" --only-show-errors || true
|
||||
done
|
||||
$S3 cp /work/ingest_time.txt "s3://$BUCKET/$RUN/ingest_time.txt" --only-show-errors || true
|
||||
echo "ingest output size: $(du -sh /work/ingest | cut -f1)"
|
||||
echo "=== done $(date -u +%T) ==="
|
||||
@@ -103,6 +103,59 @@ def check_idat_pair(args: tuple[str, str]) -> dict[str, Any]:
|
||||
return res
|
||||
|
||||
|
||||
_P2: dict[str, Any] = {}
|
||||
|
||||
|
||||
def _p2_init(manifest: Any) -> None:
|
||||
_P2.clear()
|
||||
_P2["m"] = manifest
|
||||
|
||||
|
||||
def check_raw_xy(args: tuple[str, str, str]) -> dict[str, Any]:
|
||||
"""P2: raw X/Y + controls assembled from IDATs vs the arrays stored in the GTC (TOC 1000/1001/500/501)."""
|
||||
from simurg.formats.gtc import read_gtc
|
||||
from simurg.formats.idat import read_idat
|
||||
from simurg.signal.raw_xy import LocusMapCache, assemble_raw_xy
|
||||
|
||||
gtc, grn, red = args
|
||||
res: dict[str, Any] = {"gtc": gtc, "ok": True}
|
||||
try:
|
||||
m = _P2["m"]
|
||||
cache = _P2.setdefault("cache", LocusMapCache(m))
|
||||
ref = read_gtc(gtc, check_write_complete=False)
|
||||
raw = assemble_raw_xy(read_idat(grn), read_idat(red), m, cache=cache)
|
||||
dx = raw.raw_x != ref.raw_x
|
||||
dy = raw.raw_y != ref.raw_y
|
||||
diff = dx | dy
|
||||
res.update(
|
||||
n_diff=int(diff.sum()),
|
||||
n_missing=raw.n_missing,
|
||||
vendor_zero=int(((ref.raw_x == 0) & (ref.raw_y == 0)).sum()),
|
||||
# plan §20 U1: would swapping channels explain the differences?
|
||||
swapped=int(((raw.raw_x == ref.raw_y) & (raw.raw_y == ref.raw_x) & diff).sum()),
|
||||
ours_missing_vendor_nonzero=int((raw.missing & ((ref.raw_x != 0) | (ref.raw_y != 0))).sum()),
|
||||
vendor_zero_ours_present=int((~raw.missing & (ref.raw_x == 0) & (ref.raw_y == 0)).sum()),
|
||||
controls_equal=bool(
|
||||
np.array_equal(raw.controls_x, ref.controls_x) and np.array_equal(raw.controls_y, ref.controls_y)
|
||||
),
|
||||
)
|
||||
if res["n_diff"]:
|
||||
idx = np.flatnonzero(diff)[:5].tolist()
|
||||
res["examples"] = [
|
||||
{
|
||||
"locus": m.names[i],
|
||||
"assay": int(m.assay_types[i]),
|
||||
"ours": [int(raw.raw_x[i]), int(raw.raw_y[i])],
|
||||
"vendor": [int(ref.raw_x[i]), int(ref.raw_y[i])],
|
||||
}
|
||||
for i in idx
|
||||
]
|
||||
res["ok"] = res["n_diff"] == 0 and res["controls_equal"]
|
||||
except Exception as e:
|
||||
res.update(ok=False, error=f"{type(e).__name__}: {e}")
|
||||
return res
|
||||
|
||||
|
||||
def main() -> int:
|
||||
ap = argparse.ArgumentParser()
|
||||
ap.add_argument("--root", required=True, type=Path)
|
||||
@@ -112,6 +165,7 @@ def main() -> int:
|
||||
ap.add_argument("--out", required=True, type=Path)
|
||||
ap.add_argument("--workers", type=int, default=os.cpu_count() or 1)
|
||||
ap.add_argument("--limit", type=int, default=0)
|
||||
ap.add_argument("--no-p2", action="store_true", help="skip the phase-P2 raw X/Y check")
|
||||
a = ap.parse_args()
|
||||
|
||||
from simurg.formats.bpm import read_bpm
|
||||
@@ -164,6 +218,28 @@ def main() -> int:
|
||||
"normalization_lookups": m.normalization_lookups().tolist() == b.normalization_lookups,
|
||||
"n_bins": len(np.unique(m.norm_bin_ids())),
|
||||
}
|
||||
if not a.no_p2:
|
||||
import multiprocessing as mp
|
||||
import re
|
||||
|
||||
key_re = re.compile(r"(\d+_r\d{2}c\d{2})", re.IGNORECASE)
|
||||
triples = []
|
||||
for g in gtcs:
|
||||
mk = key_re.search(Path(g).name)
|
||||
k = mk.group(1).lower() if mk else Path(g).stem.lower()
|
||||
if k in grn and k in red:
|
||||
triples.append((g, str(grn[k]), str(red[k])))
|
||||
summary["p2_matched_gtc_idat"] = len(triples)
|
||||
with ProcessPoolExecutor(
|
||||
a.workers, mp_context=mp.get_context("spawn"), initializer=_p2_init, initargs=(m.lite(),)
|
||||
) as ex:
|
||||
p2 = list(ex.map(check_raw_xy, triples, chunksize=4))
|
||||
summary["p2_ok"] = sum(r["ok"] for r in p2)
|
||||
summary["p2_total_diff_loci"] = sum(r.get("n_diff", 0) for r in p2)
|
||||
summary["p2_swapped_loci"] = sum(r.get("swapped", 0) for r in p2)
|
||||
summary["p2_ours_missing_vendor_nonzero"] = sum(r.get("ours_missing_vendor_nonzero", 0) for r in p2)
|
||||
summary["p2_vendor_zero_ours_present"] = sum(r.get("vendor_zero_ours_present", 0) for r in p2)
|
||||
summary["p2_failures"] = [r for r in p2 if not r["ok"]][:30]
|
||||
e = read_egt(a.egt)
|
||||
with open(a.egt, "rb") as fh:
|
||||
be = module.ClusterFile.read_cluster_file(fh)
|
||||
@@ -200,7 +276,10 @@ def main() -> int:
|
||||
summary["gtc_transforms_vs_bins"] = [len(g0.normalization_transforms), summary["bpm"]["n_bins"]]
|
||||
|
||||
summary["seconds"] = round(time.time() - t0, 1)
|
||||
passed = (
|
||||
p2_passed = a.no_p2 or (
|
||||
summary.get("p2_ok") == summary.get("p2_matched_gtc_idat", -1) and summary.get("p2_matched_gtc_idat", 0) > 0
|
||||
)
|
||||
passed = p2_passed and (
|
||||
summary["gtc_ok"] == len(gtcs)
|
||||
and summary["idat_ok"] == len(pairs)
|
||||
and all(v for k, v in summary["bpm"].items() if isinstance(v, bool))
|
||||
@@ -213,8 +292,19 @@ def main() -> int:
|
||||
print(
|
||||
json.dumps(
|
||||
{
|
||||
k: summary[k]
|
||||
for k in ("n_gtc", "gtc_ok", "gtc_roundtrip_identical", "n_idat_pairs", "idat_ok", "passed", "seconds")
|
||||
k: summary.get(k)
|
||||
for k in (
|
||||
"n_gtc",
|
||||
"gtc_ok",
|
||||
"gtc_roundtrip_identical",
|
||||
"n_idat_pairs",
|
||||
"idat_ok",
|
||||
"p2_matched_gtc_idat",
|
||||
"p2_ok",
|
||||
"p2_total_diff_loci",
|
||||
"passed",
|
||||
"seconds",
|
||||
)
|
||||
}
|
||||
)
|
||||
)
|
||||
|
||||
Reference in new issue
Block a user