16강 / 20강

DIGGIT ④ — Regulator score, 그리고 논문과의 정면 대조

프로토콜 Phase I-3 결과 · 우리는 어디까지 재현했나

마지막 조립입니다. 후보 유전자 하나마다 이렇게 셉니다.

Regulator score

핵심 전사인자 하나하나에 대해, 아래 세 조건을 모두 만족하면 1점

  ① MINDy  — 그 TF 의 표적 관계를 유의하게 변조한다
  ② aQTL   — 조절자의 복제수가 그 TF 의 활성과 연동된다
  ③ 독립성 — 조절자가 그 TF 자신과 강하게 상관되어 있지 않다
              (상관되어 있으면 '변조'가 아니라 '동반 조절'일 수 있다)

Regulator score = 1점을 받은 핵심 TF 의 개수

③ 이 왜 필요한지 짚어 둡시다. 조절자 M 이 전사인자 T 와 강하게 상관되어 있다면, M 으로 표본을 가르는 것은 사실상 T 로 표본을 가르는 것과 같습니다. 그러면 T—표적 관계가 달라지는 것이 당연하고, 변조의 증거가 되지 못합니다.

고등학생논문과의 정면 대조

지금까지 계산한 것을 논문 결과와 견준다.

확인	내용
같은 후보가 나오나	—
순위가 비슷한가	—
다르면 왜 다른가	가장 중요

우리는 62위로 나왔다

논문의 핵심 후보가 우리 계산에서는 62위 다.

재현 시도의 목적은 같은 답을 얻는 것이 아니라 차이의 원인을 아는 것이다

전체 DIGGIT 은 브라우저에 무거워서, 결과는 위젯에 담아 두었습니다.

Python 브라우저에서 실행

import gene

# 전체 DIGGIT 은 브라우저에서 돌리기에 무겁습니다.
# 여기서는 계산의 '규모'만 확인하고, 결과는 아래 위젯에서 봅니다.
W = gene.pivotal_tfs()
piv = W[(W.act_p < 0.05) & (W.expr_p > 0.05)]
fc  = gene.fcnv()
cand = fc[(fc.p < 1e-6) & (fc.rho > 0)]

print("DIGGIT 의 계산 규모")
print(f"  핵심 전사인자      : {len(piv):>6} 개")
print(f"  F-CNV 통과 후보    : {len(cand):>6} 개")
print(f"  검사할 (후보,TF) 쌍: {len(piv)*len(cand):>6} 개")
print(f"  각 쌍마다 표적 100개 × 두 무리 = 상관계수 {len(piv)*len(cand)*200:,} 번")
print()
print("이 계산은 브라우저에서 몇 분 걸립니다.")
print("결과는 아래 위젯에 미리 담아 두었으니 그것으로 확인하세요.")

우리가 얻은 순위표

아래 위젯에서 우리 결과와 논문 Table S2A 를 나란히 볼 수 있습니다. 먼저 결론부터 정직하게 적어 둡니다.

재현된 것

· 논문이 발표한 마스터 조절자 목록이
우리 순위표 위쪽에 유의하게 몰려 있다
(Mann-Whitney p = 0.0015)
· DDX54 는 후보 768개 중 62위 — 상위 8%
· DDX54 의 발현이 면역 사막에서 유의하게 높다
(log2FC +0.25, p = 0.025) — 논문 Fig 2D 와 일치

재현되지 않은 것

· DDX54 가 1위로 나오지는 않았다
· 우리 상위권은 세포주기·복제 관련 유전자가 차지
· 논문 상위 10개 중 일부(PCNA, RAP1B, CDK4)는
우리 순위에서 한참 아래

이것을 실패라고 부르지 않습니다. 부분 재현이라고 부르고, 왜 그만큼인지를 따지는 것이 재현 연구의 본체입니다.

차이가 생긴 이유 — 여섯 갈래

① 표본을 가른 기준이 다르다
논문은 H&E 병리 이미지 딥러닝(Saltz 등)의 brisk/non-brisk, 우리는 발현 기반 림프구 침윤 점수. 다른 환자 집합을 비교한 셈입니다. 이 갈래만은 추측하지 않고 바로 아래에서 숫자로 잽니다 — 4강에 실어 둔 그 TIL 지도로요.

② 전사인자 목록이 다르다
AnimalTFDB 3.0 vs Human TFs v1.01. CHAF1A 처럼 한쪽에만 있는 유전자가 있습니다.

③ 유전자 우주가 다르다
논문은 전체 전사체, 우리는 1,807개 조각. 조각 밖의 조절자는 후보에조차 오르지 못합니다 (논문 9위 RHOF 가 그렇습니다).

④ 부트스트랩을 생략했다
ARACNE 의 100회 재추출을 돌리지 않아 레귤론이 논문보다 덜 안정적입니다.

⑤ MI 추정 방식이 다르다
원 ARACNE 는 적응적 분할 커널 추정, 우리는 순위-가우시안 근사.

⑥ DIGGIT 구현이 다르다
원 DIGGIT 은 Bioconductor 패키지로 아주 많은 세부 설정이 있습니다. 우리는 논문에 적힌 설명(F-CNV → MINDy → aQTL → 개수 세기)을 읽고 다시 구현했습니다.

① 을 숫자로 재 본다 — 우리 선과 논문 선은 얼마나 다른가

여섯 갈래 중 첫째만은 직접 잴 수 있습니다. 4강에서 TCIA 의 TIL 지도 542장을 실어 두었으니까요.

논문의 brisk / non-brisk 는 병리의사가 지도를 눈으로 보고 붙인 라벨이라 우리에게 없습니다. 대신 지도에서 우리가 직접 센 TIL 비율에 5% 선을 그어 그 판정을 흉내 냅니다. 폐선암 458장에서 이 5% 선은 논문의 병리 라벨을 81% 맞히고 AUC 는 0.896 입니다 — 완벽하지는 않지만, 병리의사의 눈을 대신 세우기에는 충분한 자입니다.

셋 중 하나가 어긋납니다. 논문과 우리는 79명 중 26명을 반대쪽 무리에 넣고 시작한 셈입니다. 이 표본으로 만든 시그니처 → 네트워크 → 활성 → 순위가 논문과 똑같이 나오는 편이 오히려 이상하죠.

Python 브라우저에서 실행

import pandas as pd, numpy as np, gene
from scipy import stats

# 4강에서 연 그 지도들 - 환자당 한 장
til = pd.read_csv("TCGA_LUAD_TILmap.tsv", sep="\t")
til = til.sort_values("tissue_patches").drop_duplicates("participant", keep="last")
m = gene.meta().copy()
m["participant"] = [s[:12] for s in m.index]
j = m.merge(til, on="participant")

# H&E 지도로 병리 판정을 흉내낸다 - TIL 비율 5% 를 경계로 brisk / non-brisk
j["HE"] = np.where(j.til_pct >= 5, "brisk-like", "non-brisk-like")

sub = j[j.group.isin(["desert", "inflamed"])]
ct = pd.crosstab(sub.group, sub.HE)
print("우리가 발현으로 가른 무리(행) vs H&E 지도로 가른 무리(열)")
print(ct.to_string())
print()

agree = ct.loc["desert", "non-brisk-like"] + ct.loc["inflamed", "brisk-like"]
print(f"두 기준이 같은 쪽으로 부른 표본 : {agree} / {ct.values.sum()}"
      f"  ({100*agree/ct.values.sum():.0f}%)")
print(f"우리가 '사막'이라 부른 {ct.loc['desert'].sum()}명 중 "
      f"{ct.loc['desert','brisk-like']}명은 H&E 로 보면 림프구가 많다")
print(f"우리가 '염증'이라 부른 {ct.loc['inflamed'].sum()}명 중 "
      f"{ct.loc['inflamed','non-brisk-like']}명은 H&E 로 보면 거의 비어 있다")
print()
print("독립성 검정 : Fisher p =", f"{stats.fisher_exact(ct.values)[1]:.4f}",
      "- 무관하지는 않지만, 같은 자가 아니다")

고등학생왜 62위인가 — 이유들

이유	강
축소된 네트워크	7강 — 부트스트랩이 적음
단순한 MI 추정	6강
암종을 섞음	3강
자기 발현이 안 변함	12~13강
표본 규모	—

네 번째가 흥미롭다

DDX54 는 자기 발현이 안 변한다. 그러니 발현에 기대는 단계에서 불리하다.

그런데 그것이 이 방법으로 찾은 이유이기도 하다

발현 분석으로는 아예 못 찾는다. 활성으로 보면 보인다 — 62위라도.

같은 성질이 어떤 방법에서는 약점이고 어떤 방법에서는 발견의 이유가 된다

왼쪽 환자는 우리가 면역 사막이라 부른 사람인데 조직의 절반이 붉습니다. 오른쪽 환자는 우리가 면역 염증이라 부른 사람인데 붉은 조각이 하나도 없습니다. 논문이라면 두 사람의 자리가 서로 바뀌었을 겁니다.

Python 브라우저에서 실행

import gene, pandas as pd
import matplotlib.pyplot as plt

til = pd.read_csv("TCGA_LUAD_TILmap.tsv", sep="\t")
til = til.sort_values("tissue_patches").drop_duplicates("participant", keep="last")
m = gene.meta().copy()
m["participant"] = [s[:12] for s in m.index]
j = m.merge(til, on="participant")

# 두 기준이 정면으로 어긋난 표본 두 명을 눈으로 본다
fig, axes = plt.subplots(1, 2, figsize=(9.5, 4.2))
for ax, (pid, tag) in zip(axes, [("TCGA-78-7536", "our group: desert"),
                                 ("TCGA-05-4410", "our group: inflamed")]):
    a = (plt.imread(f"TILmap_{pid}.png")[:, :, :3] * 255).round().astype(int)
    r, g, b = a[:, :, 0], a[:, :, 1], a[:, :, 2]
    n_til = int(((r == 255) & (g == 0) & (b == 0)).sum())
    n_tis = n_til + int(((r == 0) & (g == 0) & (b == 255)).sum())
    pct   = 100 * n_til / n_tis
    row   = j[j.participant == pid].iloc[0]
    print(f"{pid} - {tag} (LI = {row.LI:+.2f}) -> H&E TIL {pct:.1f}%")
    ax.imshow(a)
    ax.axis("off")
    ax.set_title(f"{pid}\n{tag}\nH&E TIL {pct:.1f}%", fontsize=10)
fig.tight_layout()

그래도 발현 점수를 버리지 않는 이유

H&E 지도가 더 옳은 자라서 논문이 그것을 쓴 것은 아닙니다. 둘은 서로 다른 것을 봅니다.

H&E 지도 — 슬라이드 한 장에서 림프구가 어디에 있는지까지 봅니다. 대신 그 한 장이 종양 전체를 대표하는지는 아무도 보증하지 않습니다 (위 오른쪽 환자의 지도는 조직 조각이 5,618개뿐입니다 — 아주 작게 잘린 슬라이드죠).

발현 점수 — RNA 를 뽑은 조각 전체의 평균이라 표본 대표성은 낫지만, 위치 정보가 없습니다. 종양 안에 들어와 있는 림프구와 가장자리에 멈춰 선 림프구를 구별하지 못합니다.

논문의 결론이 "안으로 못 들어온다" 는 공간적 주장이었으니 H&E 쪽을 쓴 것이 맞습니다. 우리는 그 자를 통째로 쓸 수 없어 대신 쓸 수 있는 자를 골랐고, 그 대가가 이 26명입니다.

🎮 순위표 대조

우리 결과와 논문 Table S2A 를 나란히 놓고, 논문의 마스터 조절자들이 우리 순위표 어디에 있는지 표시합니다. DDX54 도 눈에 띄게 표시됩니다.

순위는 못 맞혔지만, 논문이 든 유전자들이 실제로 면역 사막에서 높다는 것은 확인됩니다.

Python 브라우저에서 실행

import gene

# 논문 주장 하나는 우리 자료로도 독립적으로 확인됩니다 —
# "DDX54 는 마스터 조절자 후보 중 유일하게 두 무리 사이 발현 차이가 유의했다" (Fig 2D)
sig = gene.signature()
mr  = gene.paper_mr()
desert = mr[mr.phenotype == "desert"].head(10)

print("논문 상위 10개 마스터 조절자의 — 우리 자료에서의 발현 차이 (사막 − 염증)")
print(f"{'유전자':<10}{'논문점수':>8}{'log2FC':>10}{'p':>10}")
for _, row in desert.iterrows():
    g = row.gene
    if g in sig.index:
        r = sig.loc[g]
        star = " *" if r.p < 0.05 else ""
        print(f"{g:<10}{row.score:>8}{r.log2FC:>+10.3f}{r.p:>10.4f}{star}")
    else:
        print(f"{g:<10}{row.score:>8}{'조각에 없음':>20}")
print()
print("* = p < 0.05 — DDX54 를 포함해 여러 개가 우리 자료에서도 유의하게 높습니다.")

순위는 못 맞혔습니다. 그러면 생물학도 틀렸을까요?

여기서 멈추면 "재현 실패" 라고 쓰고 끝납니다. 하지만 우리가 실제로 못 맞힌 것은 순위 하나입니다. 논문의 진짜 주장은 순위가 아니라 이것이었어요 —

DDX54 가 높은 종양은 발암·면역회피 경로가 켜져 있다

이건 DIGGIT 없이도 물어볼 수 있는 질문입니다. 표본을 DDX54 발현으로 가르고, 그 두 무리 사이에서 경로가 움직였는지만 보면 되니까요. 논문 저자들도 정확히 그렇게 했습니다 — SI Fig. S4 입니다.

그리고 중요한 점: 저자들은 그 확인을 TCGA 가 아닌 다른 자료로 했습니다. 같은 자료로 같은 결론을 다시 내면 그건 확인이 아니라 메아리니까요.

코호트	무엇	표본	단위
GSE81089	스웨덴 우메오 NSCLC 코호트	종양 197	FPKM
GSE171415	Sherlock-Lung (NCI) 폐암	종양 35	TPM
TCGA-LUAD	우리가 3강부터 써 온 조각	종양 512	RSEM

둘 다 GEO 에서 그대로 받아, 우리가 쓰는 1,807 유전자 조각으로 줄여 실어 두었습니다. 줄였기 때문에 치를 대가가 있는데, 그건 곧 보게 됩니다.

고등학생62위를 어떻게 읽어야 하나

읽기	평가
재현 실패	너무 극단적
재현 성공	과장
부분 재현	정확

62위는 상위 몇 %인가

전사인자가 1,500개면 상위 4% 다.

무작위였다면 기대 순위가 750위 다.

그러니 신호는 있다

다만 1위는 아니다.

순위를 백분율로 바꿔 보면 그것이 우연 수준인지 알 수 있다

서로 다른 나라 · 다른 시퀀싱 · 다른 정규화. 그래도 DDX54 는 환자마다 다릅니다.

Python 브라우저에서 실행

# 외부 코호트 두 개를 열어 본다 - TCGA 와 완전히 다른 환자들이다
import gene
import pandas as pd

for 이름 in ["GSE81089", "GSE171415"]:
    M = gene.ext(이름)
    d = M.loc["DDX54"]
    print(f"{이름:<12} {M.shape[0]:>5} 유전자 x {M.shape[1]:>4} 종양   "
          f"DDX54 log2 발현 {d.min():.2f} ~ {d.max():.2f} (중앙값 {d.median():.2f})")

T = gene.expr()
d = T.loc["DDX54"]
print(f"{'TCGA-LUAD':<12} {T.shape[0]:>5} 유전자 x {T.shape[1]:>4} 종양   "
      f"DDX54 log2 발현 {d.min():.2f} ~ {d.max():.2f} (중앙값 {d.median():.2f})")
print()
print("세 코호트 모두 DDX54 발현이 환자마다 크게 다릅니다 - 가를 수 있다는 뜻입니다.")

이제 각 코호트에서 DDX54 상위 3분의 1 과 하위 3분의 1 을 갈라, 그 사이에서 경로 여섯 개가 움직였는지 봅니다. 농축 검정은 새로 짜지 않습니다 — 11강에서 만든 순위 농축을 그대로 씁니다.

gene.pathway_table(M)
   → DDX54 상위/하위로 가르고 · t 통계량으로 줄 세우고 · 여섯 경로를 순열검정

NES 가 +면 DDX54 가 높은 종양에서 그 경로가 켜져 있다는 뜻입니다.

Python 브라우저에서 실행

# 논문 SI Fig. S4 를 우리 손으로 - 외부 코호트 197명
import gene
import pandas as pd
pd.set_option("display.width", 100)

표 = gene.pathway_table(gene.ext("GSE81089"))
print("GSE81089 (스웨덴 NSCLC 197) - DDX54 상위 66 vs 하위 66")
print(표.round(3))

표본 512명에서는 신호가 더 크고, 35명에서는 대부분 사라집니다.

Python 브라우저에서 실행

# 우리 TCGA 조각에서도 같은가 - 그리고 표본이 적으면 어떻게 되나
import gene
import pandas as pd
pd.set_option("display.width", 100)

for 이름, M in [("TCGA-LUAD 512", gene.expr()),
                ("GSE171415 35", gene.ext("GSE171415"))]:
    print(f"=== {이름} ===")
    print(gene.pathway_table(M)[["NES", "p"]].round(3))
    print()

세 코호트를 나란히 놓으면 이렇습니다.

경로	GSE81089 (197)	TCGA-LUAD (512)	GSE171415 (35)	논문
MYC_TARGETS	+2.78 ✔	+4.08 ✔	−1.30 —	증가
EMBRYONIC_STEM	+2.85 ✔	+4.21 ✔	+1.38 —	증가
TNFA_NFKB	−1.80	−1.46	−1.77	증가
EMT	−1.78	−2.22	−1.15	증가
IL6_JAK_STAT3	+0.84	−0.98	−0.80	증가
WNT_BETA_CATENIN	검정할 수 없었다 — 아래 참조			증가

둘은 표본이 넉넉한 두 코호트에서 크고 일관되게 재현되고, 셋은 오히려 반대 방향입니다. 이것이 우리가 실제로 얻은 결과이고, 그대로 적는 것이 재현 연구입니다.

고등학생TIL 지도와의 대조

4강에서 실어 둔 TCIA TIL 지도 542장 을 여기서 쓴다.

확인	내용
DDX54 활성이 높으면	TIL 이 적나
독립적 자료	발현이 아니라 조직 사진

이것이 독립적 검증이다

발현으로 후보를 찾고 조직 사진 으로 확인한다.

서로 다른 종류의 자료 이니 같은 편향이 반복될 이유가 없다.

그리고 18강의 태아 흉내와 연결된다

인간 자료에서의 근거가 여기 에 있다.

다른 종류의 자료가 같은 답을 주면 그 답이 훨씬 믿을 만해진다

🤔 확인 문제

논문의 순위(1위)를 재현하지 못했는데도 외부 코호트에서 Myc·배아줄기 시그니처가 재현된 것은 무엇을 뜻하나요?

🤔 확인 문제

우리 결과에서 논문의 마스터 조절자 목록이 상위에 유의하게 몰려 있다는 것은 무엇을 뜻하나요?

🤔 확인 문제

논문 9위 RHOF 가 우리 순위표에 아예 없는 이유는?

🤔 확인 문제

우리 DIGGIT 상위권을 세포주기 유전자들이 차지하는 현상의 원인은?