16. 짧은 조각에서는 E-value 가 아무도 거르지 못한다 · 정교수의 mRNA백신 디자인 첫걸음 I

앞 강의 마지막 두 줄로 돌아옵니다.

E-value > 10⁻⁴ 그리고 bit score < 100

이 강에서 할 일은 하나입니다 — 이 기준을 논문이 고른 스무 개 서열에 실제로 대어 보는 것입니다.

E-value 가 무엇인지 다시

BLAST 의 E-value 는 이렇게 생겼습니다.

E = m × n ÷ 2^비트점수

기호	무엇
m	찾는 조각의 길이 — 여기서는 9 · 15 · 16
n	뒤지는 데이터베이스의 크기 (아미노산 총수)
비트점수	일치가 얼마나 좋은가 — 높을수록 좋은 일치

여기서 m 이 분자에 있습니다. 조각이 짧으면 E 도 작아질 것 같지만, 훨씬 크게 움직이는 것은 분모입니다. 조각이 짧으면 비트점수가 클 수 없고, 그러면 분모가 작아져 E 가 커집니다.

100% 일치인데 전부 통과합니다. 체가 작동하지 않습니다.

Python 브라우저에서 실행

# 9잔기가 100% 똑같을 때 E-value 가 얼마나 되는지 계산한다.
# BLOSUM62 자기 자신 점수(대각선)와 Karlin-Altschul 표준 상수를 쓴다.
import math

B62 = {"A": 4, "R": 5, "N": 6, "D": 6, "C": 9, "Q": 5, "E": 5, "G": 6,
       "H": 8, "I": 4, "L": 4, "K": 5, "M": 5, "F": 6, "P": 7, "S": 4,
       "T": 5, "W": 11, "Y": 7, "V": 4}
lam = 0.318          # BLOSUM62 무갭 lambda
K = 0.134            # BLOSUM62 무갭 K
n = 9e7              # 사람 단백질체 + 장내세균 79종 아미노산 총수 (어림)

def evalue(pep):
    raw = sum(B62[c] for c in pep)                 # 100% 일치일 때의 원점수
    bits = (lam * raw - math.log(K)) / math.log(2)
    return raw, bits, len(pep) * n / 2 ** bits

peps = ["KTYSNITIA", "MEAAYTASL", "LPYFHNINY", "VQMTFVISV"]
print("9잔기 CTL 에피토프가 사람 단백질과 100% 똑같다고 하자")
print()
print("%-12s %5s %7s %12s %s" % ("조각", "원점수", "비트", "E-value", "판정"))
for p in peps:
    raw, bits, e = evalue(p)
    verdict = "통과" if e > 1e-4 else "탈락"
    print("%-12s %5d %7.1f %12.3g %s" % (p, raw, bits, e, verdict))
print()
print("기준은 'E-value > 1e-4 면 통과(안 닮음)' 이다.")
print("100% 똑같은데도 전부 통과한다.")

얼마나 모자라나

E-value 가 10⁻⁴ 밑으로 내려가려면 비트점수가 43.7 을 넘어야 합니다. 9잔기가 100% 똑같을 때 나오는 비트점수는 22~28 입니다.

필요한 비트점수

43.7 이상

9잔기의 천장

약 23~28 — 전부 트립토판이어도 48 이 한계

천장이 기준보다 낮습니다. 그러니까 어떤 9잔기 조각도 이 체에 걸릴 수 없습니다.

스무 개 전부가 두 기준을 통과합니다 — 100% 일치라고 가정해도요.

Python 브라우저에서 실행

# 스무 개 전부에 같은 계산을 해 본다. 15잔기·16잔기도 포함해서.
import math
B62 = {"A": 4, "R": 5, "N": 6, "D": 6, "C": 9, "Q": 5, "E": 5, "G": 6,
       "H": 8, "I": 4, "L": 4, "K": 5, "M": 5, "F": 6, "P": 7, "S": 4,
       "T": 5, "W": 11, "Y": 7, "V": 4}
lam, K, n = 0.318, 0.134, 9e7

ctl = ["KTYSNITIA", "TYSNITIAM", "LPYFHNINY", "VQMTFVISV",
       "VSIPTNFSF", "IPTNFSFGL", "MEAAYTASL", "LASFAAIPF"]
htl = ["TTARIKKIYPAFVLG", "TIVISPTTTARIKKI", "YFHNINYYSVIPRSP",
       "ISPTTTARIKKIYPA", "WTAGLASFAAIPFAQ", "LYFFHVGYQPTALIN",
       "RYVAMDVKFENITTK"]
lbl = ["FPLSMASYLRPGSTGP", "SMTERVQMTFVISVTF", "DVESGVYSVSSFEAKS",
       "PGSTGPTAEFNYRQDF", "SDGQSASMTERVQMTF"]

passed = 0
total = 0
maxbits = 0.0
for group, peps in (("CTL", ctl), ("HTL", htl), ("LBL", lbl)):
    for p in peps:
        raw = sum(B62[c] for c in p)
        bits = (lam * raw - math.log(K)) / math.log(2)
        e = len(p) * n / 2 ** bits
        total += 1
        maxbits = max(maxbits, bits)
        if e > 1e-4:
            passed += 1

print("스무 개가 사람 단백질과 100% 똑같다고 가정했을 때")
print("  E-value > 1e-4 기준을 통과한 것:", passed, "/", total)
print("  가장 높은 비트점수:", round(maxbits, 1))
print("  bit score < 100 기준을 통과한 것:", total, "/", total)
print()
print("두 기준 모두 스무 개 전부를 통과시킨다.")
print("(이 계산은 우리가 한 것이다 — 논문에 없다. 인용하지 말 것)")

조각 길이와 데이터베이스 크기를 움직여 보세요. 어느 길이부터 체가 작동하기 시작하는지 보입니다.

불러오는 중…

그래서 무엇이 문제인가

논문이 거짓말을 한 것이 아닙니다. 기준을 적었고, 도구를 돌렸고, 결과가 '통과' 였습니다. 전부 사실입니다.

문제는 이 단계가 아무것도 걸러 내지 못했는데도 걸러 낸 것처럼 보인다는 데 있습니다. 독자는 '자가면역 위험을 확인했다' 고 읽습니다.

통과했다는 말과 걸러졌다는 말은 다르다

짧은 펩타이드의 자기유사성을 제대로 보려면 E-value 가 아니라 정확 일치 탐색이나 MHC 제시까지 본 비교가 필요합니다. 그런 도구가 따로 있습니다.

🤔 확인 문제

9잔기 조각에서 E-value 기준이 작동하지 않는 까닭은?

🤔 확인 문제

이 발견을 가장 정확히 표현한 것은?

🤔 확인 문제

짧은 펩타이드의 자기유사성을 보는 더 나은 방법은?