#!/usr/bin/env python3 """ AdipoGlowIVIS - 오프라인 CLI 시연 (표준 라이브러리만 사용) rodent IVIS 광학영상 지방 depot ROI 정량 도구. 이 main.py 는 numpy/pandas/scipy 가 설치되지 않아도 동작합니다. ( app.py 는 streamlit 기반 GUI 로 별도. ) 사용법: python3 main.py --demo # data/ 의 샘플로 전체 분석 파이프라인 출력 python3 main.py --help # 도움말 본 도구는 연구용·참고용이며, 정량 결과는 사용자 검증이 필요합니다. 샘플 데이터는 합성(synthetic) 데이터입니다. """ import sys import os import csv import json import math import argparse HERE = os.path.dirname(os.path.abspath(__file__)) DATA_DIR = os.path.join(HERE, "data") DEFAULT_CSV = os.path.join(DATA_DIR, "sample_roi.csv") DEFAULT_REF = os.path.join(DATA_DIR, "reporter_reference.json") DISCLAIMER = ("본 도구는 연구용·참고용이며, 정량 결과는 사용자 검증이 필요합니다. " "샘플 데이터는 합성(synthetic) 데이터입니다.") # normalize 기준값: 노출 30초, binning 8 을 reference 로 삼아 # 이미지 간 비교 가능한 단위로 환산한다. REF_EXPOSURE = 30.0 REF_BINNING = 8.0 # --------------------------------------------------------------------------- # 1. ingest + 노출/binning normalize + QC # --------------------------------------------------------------------------- def load_csv(path): rows = [] with open(path, newline="") as f: reader = csv.DictReader(f) for r in reader: rec = { "animal_id": r["animal_id"].strip(), "group": r["group"].strip(), "timepoint": r["timepoint"].strip(), "depot": r["depot"].strip(), "total_flux_ps": _f(r["total_flux_ps"]), "avg_radiance": _f(r["avg_radiance"]), "exposure_sec": _f(r["exposure_sec"]), "binning": _f(r["binning"]), "f_stop": _f(r["f_stop"]), "substrate_post_min": _f(r["substrate_post_min"]), "is_background": int(float(r["is_background"])), } rows.append(rec) return rows def _f(x): if x is None or str(x).strip() == "": return None return float(x) def normalize_signal(value, exposure_sec, binning): """노출시간·binning 으로 normalize. 신호는 노출시간에 비례, binning 면적(픽셀 합)에 비례한다고 가정. normalized = value * (REF_EXPOSURE/exposure) * (REF_BINNING/binning) -> 노출 길수록 down-scale, binning 클수록 down-scale. """ if exposure_sec in (None, 0) or binning in (None, 0): return None return value * (REF_EXPOSURE / exposure_sec) * (REF_BINNING / binning) def add_normalized(rows): for r in rows: r["norm_total_flux"] = normalize_signal( r["total_flux_ps"], r["exposure_sec"], r["binning"]) r["norm_avg_radiance"] = normalize_signal( r["avg_radiance"], r["exposure_sec"], r["binning"]) return rows def qc_flags(rows): """QC: saturation 추정, 노출 불일치, 기질 시점 미기록.""" flags = [] # saturation: avg_radiance 가 비현실적으로 큰지 (시연 임계 1.0e6) SAT = 1.0e6 exposures = sorted(set(r["exposure_sec"] for r in rows if r["exposure_sec"])) for r in rows: rid = "%s/%s/%s" % (r["animal_id"], r["timepoint"], r["depot"]) if r["avg_radiance"] and r["avg_radiance"] >= SAT: flags.append((rid, "SATURATION 의심 (avg_radiance >= %.0e)" % SAT)) if r["substrate_post_min"] is None and not r["is_background"]: flags.append((rid, "기질 주입 후 경과시간 미기록")) if len(exposures) > 1: flags.append(("DATASET", "노출시간 불일치 감지: %s (normalize 적용됨)" % ", ".join("%gs" % e for e in exposures))) return flags # --------------------------------------------------------------------------- # 2. 배경 차감 + depot 비교 # --------------------------------------------------------------------------- def background_map(rows): """(animal, timepoint) -> 배경 normalized avg_radiance / total_flux.""" bg = {} for r in rows: if r["is_background"]: key = (r["animal_id"], r["timepoint"]) bg[key] = { "rad": r["norm_avg_radiance"], "flux": r["norm_total_flux"], } return bg def subtract_background(rows): bg = background_map(rows) out = [] for r in rows: if r["is_background"]: continue key = (r["animal_id"], r["timepoint"]) b = bg.get(key, {"rad": 0.0, "flux": 0.0}) rr = dict(r) rr["corr_avg_radiance"] = (r["norm_avg_radiance"] or 0.0) - (b["rad"] or 0.0) rr["corr_total_flux"] = (r["norm_total_flux"] or 0.0) - (b["flux"] or 0.0) out.append(rr) return out def depot_table(corr_rows): """depot별 신호 테이블 (total flux vs avg radiance 병기).""" agg = {} for r in corr_rows: d = r["depot"] agg.setdefault(d, {"flux": [], "rad": []}) agg[d]["flux"].append(r["corr_total_flux"]) agg[d]["rad"].append(r["corr_avg_radiance"]) table = [] for d, v in agg.items(): table.append({ "depot": d, "n": len(v["flux"]), "mean_corr_total_flux": mean(v["flux"]), "mean_corr_avg_radiance": mean(v["rad"]), }) order = {"iBAT": 0, "iWAT": 1, "eWAT": 2} table.sort(key=lambda x: order.get(x["depot"], 99)) return table # --------------------------------------------------------------------------- # 3. 종단 within-animal fold + 기질 동역학 # --------------------------------------------------------------------------- def within_animal_fold(corr_rows, baseline_tp="baseline"): """동일 개체 baseline 대비 fold-change (paired). depot별.""" # index: (animal, depot, timepoint) -> corr_avg_radiance idx = {} timepoints = set() for r in corr_rows: idx[(r["animal_id"], r["depot"], r["timepoint"])] = r["corr_avg_radiance"] timepoints.add(r["timepoint"]) folds = [] for (animal, depot, tp), val in idx.items(): if tp == baseline_tp: continue base = idx.get((animal, depot, baseline_tp)) if base in (None, 0): folds.append({"animal_id": animal, "depot": depot, "timepoint": tp, "fold": None, "flag": "baseline 누락"}) else: folds.append({"animal_id": animal, "depot": depot, "timepoint": tp, "fold": val / base, "flag": ""}) folds.sort(key=lambda x: (x["animal_id"], x["depot"])) return folds def substrate_kinetics(rows): """기질(luciferin) 주입 후 경과시간 기록 + peak 시점 정렬 정보.""" times = [r["substrate_post_min"] for r in rows if r["substrate_post_min"] is not None and not r["is_background"]] if not times: return {"recorded": False} return { "recorded": True, "min_post_min": min(times), "max_post_min": max(times), "mean_post_min": mean(times), "note": "peak 정렬: 측정 시점이 다르면 peak(통상 luciferin 10-15분) 기준 보정 권장", } # --------------------------------------------------------------------------- # 4. reporter / model reference + QC flag (load) # --------------------------------------------------------------------------- def load_reference(path): with open(path) as f: return json.load(f) def classify_pharmacology(fold_rows): """약력 분류: BAT 활성형 vs 지방염증 감소형 vs depot 선택형. depot별 평균 fold 로 단순 규칙 분류 (시연용).""" by_depot = {} for fr in fold_rows: if fr["fold"] is None: continue by_depot.setdefault(fr["depot"], []).append(fr["fold"]) means = {d: mean(v) for d, v in by_depot.items()} ibat = means.get("iBAT", 1.0) iwat = means.get("iWAT", 1.0) ewat = means.get("eWAT", 1.0) if ibat >= 1.5 and ibat > ewat * 1.2: label = "BAT 활성형 (iBAT 신호 우세 증가)" elif ewat <= 0.7: label = "지방염증 감소형 (eWAT 신호 감소)" elif ibat >= 1.2 and (iwat >= 1.2 or ewat <= 0.9): label = "depot 선택형 (depot별 차등 반응)" else: label = "변화 미미 / 미분류" return {"means": means, "label": label} # --------------------------------------------------------------------------- # 5. cohort 통계 (표준 라이브러리 t-test 구현) # --------------------------------------------------------------------------- def mean(xs): xs = [x for x in xs if x is not None] return sum(xs) / len(xs) if xs else 0.0 def std(xs, ddof=1): xs = [x for x in xs if x is not None] n = len(xs) if n <= ddof: return 0.0 m = mean(xs) return math.sqrt(sum((x - m) ** 2 for x in xs) / (n - ddof)) def paired_ttest(a, b): """paired t-test, 표준 라이브러리만. (a,b 동일 길이) t 분포 p값은 정규근사로 근사 보고 (scipy 없음).""" diffs = [x - y for x, y in zip(a, b)] n = len(diffs) if n < 2: return {"t": None, "df": n - 1, "p_approx": None, "n": n} md = mean(diffs) sd = std(diffs, ddof=1) if sd == 0: return {"t": float("inf"), "df": n - 1, "p_approx": 0.0, "n": n} t = md / (sd / math.sqrt(n)) # 정규근사 양측 p (df 큰 경우 근사). 작은 n 에서는 보수적 참고치. p = 2 * (1 - _norm_cdf(abs(t))) return {"t": t, "df": n - 1, "p_approx": p, "n": n, "mean_diff": md} def _norm_cdf(z): return 0.5 * (1 + math.erf(z / math.sqrt(2))) def cohort_stats(corr_rows, depot="iBAT", timepoint="week2"): """군 × depot × 시점 요약 + 군간(unpaired 느낌) 비교용 그룹별 평균. 여기서는 군 내 baseline vs timepoint paired t-test 를 depot별로 수행.""" # group -> animal -> {baseline, tp} out = {} by_group = {} for r in corr_rows: if r["depot"] != depot: continue by_group.setdefault(r["group"], {}).setdefault(r["animal_id"], {}) by_group[r["group"]][r["animal_id"]][r["timepoint"]] = r["corr_avg_radiance"] for g, animals in by_group.items(): base, tp = [], [] for a, tps in animals.items(): if "baseline" in tps and timepoint in tps: base.append(tps["baseline"]) tp.append(tps[timepoint]) stat = paired_ttest(tp, base) out[g] = { "n": len(base), "mean_baseline": mean(base), "mean_%s" % timepoint: mean(tp), "paired_t": stat["t"], "p_approx": stat["p_approx"], } return out # --------------------------------------------------------------------------- # 리포트 (국문/영문) # --------------------------------------------------------------------------- def build_report(depot_tbl, fold_rows, pharm, cohort): ko = [] ko.append("[방법] rodent IVIS depot ROI total flux·average radiance 를 " "노출(%gs)·binning(%g) 기준으로 normalize 하고, 동일 개체·시점의 " "배경 ROI 를 차감하였다. depot별 신호 테이블, 동일 개체 baseline 대비 " "fold-change, 군별 paired t-test(정규근사) 를 산출하였다." % (REF_EXPOSURE, REF_BINNING)) ko.append("[결과] depot별 평균 fold: " + ", ".join("%s=%.2f" % (k, v) for k, v in pharm["means"].items()) + ". 약력 분류: %s." % pharm["label"]) en = [] en.append("[Methods] Rodent IVIS depot ROI total flux and average radiance " "were normalized to %gs exposure / binning %g, and per-animal/timepoint " "background ROI was subtracted. Depot signal tables, within-animal " "baseline fold-change, and per-group paired t-tests (normal approx.) " "were computed." % (REF_EXPOSURE, REF_BINNING)) en.append("[Results] Mean fold per depot: " + ", ".join("%s=%.2f" % (k, v) for k, v in pharm["means"].items()) + ". Pharmacology class: %s." % pharm["label"]) return "\n".join(ko), "\n".join(en) # --------------------------------------------------------------------------- # CLI demo # --------------------------------------------------------------------------- def run_demo(csv_path=DEFAULT_CSV, ref_path=DEFAULT_REF): line = "=" * 64 print(line) print("AdipoGlowIVIS — rodent IVIS 지방 depot 광학영상 정량 (CLI demo)") print(line) print(DISCLAIMER) print() rows = load_csv(csv_path) rows = add_normalized(rows) print("[1] Ingest + 노출/binning normalize") print(" 총 %d ROI rows 로드 (배경 포함). reference 노출=%gs, binning=%g" % (len(rows), REF_EXPOSURE, REF_BINNING)) sample = [r for r in rows if r["depot"] == "iBAT"][:2] for r in sample: print(" %s/%s iBAT: raw avg_rad=%.3e exp=%gs -> norm=%.3e" % (r["animal_id"], r["timepoint"], r["avg_radiance"], r["exposure_sec"], r["norm_avg_radiance"])) print() print("[1b] QC flags") for rid, msg in qc_flags(rows): print(" - %s: %s" % (rid, msg)) print() corr = subtract_background(rows) print("[2] 배경 차감 + depot 비교 (background-corrected, normalized)") print(" %-6s %4s %18s %18s" % ("depot", "n", "mean_total_flux", "mean_avg_rad")) for t in depot_table(corr): print(" %-6s %4d %18.3e %18.3e" % (t["depot"], t["n"], t["mean_corr_total_flux"], t["mean_corr_avg_radiance"])) print(" (total flux=전체 신호, avg radiance=밀도 → ROI 크기 의존성 분리)") print() folds = within_animal_fold(corr) print("[3] 종단 within-animal fold-change (week2 / baseline, paired, avg_radiance)") by_gd = {} for fr in folds: if fr["fold"] is None: print(" %s %s: FLAG=%s" % (fr["animal_id"], fr["depot"], fr["flag"])) continue by_gd.setdefault(fr["depot"], []).append(fr["fold"]) for depot, vals in by_gd.items(): print(" %-6s mean fold=%.2f (n=%d, range %.2f-%.2f)" % (depot, mean(vals), len(vals), min(vals), max(vals))) sk = substrate_kinetics(rows) if sk["recorded"]: print(" 기질 동역학: 주입 후 %g~%g분 (평균 %.1f분). %s" % (sk["min_post_min"], sk["max_post_min"], sk["mean_post_min"], sk["note"])) print() ref = load_reference(ref_path) print("[4] reporter/model reference 내장 (로드 확인)") for name, info in ref["reporters"].items(): print(" - %s: primary=%s, substrate=%s, peak %s분" % (name, info["primary_depot"], info["substrate"], info["typical_peak_min"])) print() print("[5] cohort 통계 (iBAT, week2 vs baseline, 군별 paired t-test)") cohort = cohort_stats(corr, depot="iBAT", timepoint="week2") for g, s in cohort.items(): p = s["p_approx"] pstr = "%.4f" % p if p is not None else "NA" print(" [%s] n=%d baseline=%.3e week2=%.3e paired_t=%s p~%s" % (g, s["n"], s["mean_baseline"], s["mean_week2"], ("%.2f" % s["paired_t"]) if s["paired_t"] is not None else "NA", pstr)) pharm = classify_pharmacology(folds) print() ko, en = build_report(depot_table(corr), folds, pharm, cohort) print("[리포트 / Report]") print(" " + ko.replace("\n", "\n ")) print() print(" " + en.replace("\n", "\n ")) print() print(line) print("demo 완료. (p값은 scipy 미사용 정규근사 — 정식 분석시 scipy paired t-test 권장)") print(line) def main(argv=None): parser = argparse.ArgumentParser( prog="main.py", description="AdipoGlowIVIS — rodent IVIS 지방 depot 광학영상 정량 (오프라인 CLI, 표준 라이브러리만).", epilog=DISCLAIMER, ) parser.add_argument("--demo", action="store_true", help="data/ 의 샘플 데이터로 전체 분석 파이프라인을 stdout 출력") parser.add_argument("--csv", default=DEFAULT_CSV, help="ROI CSV 경로 (기본: data/sample_roi.csv)") parser.add_argument("--ref", default=DEFAULT_REF, help="reporter reference JSON 경로 (기본: data/reporter_reference.json)") args = parser.parse_args(argv) if args.demo: run_demo(args.csv, args.ref) else: parser.print_help() print("\n힌트: python3 main.py --demo 로 시연을 실행하세요.") if __name__ == "__main__": main()