19강 / 20강

실습 종합 — 코드랩 0~3단계를 줄 단위로

<code>rosie_mct_mrna_vaccine_codelab.ipynb</code> · 3시간 워크숍의 앞 절반

워크숍에서 여는 노트북은 Google Colab 에서 런타임 → 모두 실행 으로 그대로 돕니다. 여기서는 앞 네 단계를 줄 단위로 읽습니다.

단계	하는 일	이 강의
0	도구 설치 — bwa · minimap2 · samtools · bcftools	12강
1	Ensembl 에서 개 KIT 유전자좌 20 kb 받기	6·7강
2	합성 FASTQ 만들기 + 드라이버 심기	8·9·10강
3	정렬 → 정렬(sort) → 중복표시 → 커버리지	11~17강

0단계 — 왜 도구를 '찾는' 코드가 필요한가

def _resolve_tools():

    micromamba = os.path.expanduser("~/.ghbio-micromamba/envs/ghbio-gen/bin")

    for n in ("bwa", "minimap2", "samtools", "bcftools"):

        p = shutil.which(n)

Colab 에서는 관리자 권한이 있으므로 apt install 이 됩니다. 학교 서버나 연구실 워크스테이션에서는 대개 안 됩니다. 그래서 이미 설치된 것을 찾는 길을 함께 둡니다.

이것이 재현성의 첫 단추입니다 — "내 컴퓨터에서는 되는데" 를 없애는 일입니다.

2단계 — '심어 두고 다시 찾기'

simulate_reads(ref, out_r1, out_r2, coverage, read_len,

               driver_pos0, driver_alt, vaf, err=0.002, insert=350, seed=1)

정상(혈액)은 야생형만, 종양은 드라이버를 VAF 0.847 로 가진 리드를 만듭니다.

이 방식의 값어치는 정답을 알고 있다는 것입니다. 파이프라인이 그것을 되찾아 내면 파이프라인이 옳게 돈다는 뜻입니다. 실제 데이터에서는 절대 알 수 없는 확인입니다.

기대값을 먼저 적어 두고 돌리는 것. 워크숍에서 학생에게 가장 먼저 가르치는 습관입니다.

Python 브라우저에서 실행

# 노트별 매개변수가 결과를 어떻게 정하나 — 심기 전에 계산해 본다
L, read_len, insert = 20_000, 100, 350

for cov in (20, 30, 60):
    n_pairs = max(1, int(cov * L / (2 * read_len)))
    print(f"커버리지 {cov:>3}x  →  짝 {n_pairs:>6,}개  "
          f"= 리드 {n_pairs*2:>7,}개  ≈ FASTQ {n_pairs*2*(read_len*2+60)/1e6:>5.1f} MB")

print()
vaf, cov = 0.847, 60
alt = cov * vaf
print(f"VAF {vaf} · 깊이 {cov}x  →  드라이버 자리에서 변이 리드 약 {alt:.0f}개")
print(f"                              참조 리드 약 {cov-alt:.0f}개")
print()
print("이 숫자를 미리 알고 4단계로 갑니다. 나오는 값이 이것과 다르면 어딘가 잘못된 것입니다.")

3단계 — 한 줄씩

명령	무엇	확인할 것
`samtools faidx kit_locus.fa`	참조 색인 (7강)	`.fai` 가 생겼는가
`bwa index kit_locus.fa`	정렬용 FM-index (11강)	`.bwt .sa .amb .ann .pac` 다섯 개
`bwa mem -R '@RG…'`	정렬 (12강)	읽기 그룹이 붙었는가 — `samtools view -H`
`samtools sort -o …`	좌표 정렬 (13강)	헤더에 `SO:coordinate`
`samtools index`	BAM 색인 (13강)	`.bai` 가 생겼는가
`samtools markdup`	중복 표시 (16강)	중복률이 얼마인가
`samtools depth`	커버리지 (17강)	종양 : 정상 ≈ 2 : 1 인가

이 일곱 줄이 이 권의 요약입니다. BAM 하나를 믿을 수 있게 만드는 것이 I 권이 한 일입니다.

Python 브라우저에서 실행

# 3단계가 끝났을 때 무엇을 확인해야 하나 — 점검표를 코드로
checks = [
    ("참조 색인 .fai 가 있다",             True,  "없으면 아무 도구도 시작 못 한다"),
    ("BAM 헤더에 @RG SM: 이 있다",         True,  "없으면 II 권의 뺄셈이 시료를 못 가른다"),
    ("BAM 헤더에 SO:coordinate 가 있다",   True,  "없으면 index 가 안 걸린다"),
    (".bai 색인이 있다",                    True,  "없으면 구간 조회가 파일 전체 훑기가 된다"),
    ("종양 평균 깊이 / 정상 평균 깊이 ≈ 2", True,  "다르면 시료가 바뀌었을 수 있다"),
    ("중복률이 30% 미만이다",               True,  "넘으면 라이브러리 복잡도가 모자란 것"),
    ("드라이버 자리 깊이가 20x 이상",       True,  "얕으면 4단계에서 못 찾는다"),
]
ok = 0
for name, passed, why in checks:
    mark = "✔" if passed else "✘"
    print(f" {mark}  {name}")
    print(f"      {why}")
    ok += passed
print()
print(f"{ok} / {len(checks)} 통과 — 전부 통과해야 4단계로 갑니다.")
print("하나라도 어긋나면, 뒤에서 나오는 결과는 전부 의심스럽습니다.")

워크숍 진행 — 앞 90분

시간	무엇	참고 강
0~15분	로지 이야기 · 비만세포종이 무엇인지	1~3강
15~25분	Colab 열기 · 0단계 실행 (도구 설치)	12강
25~40분	1단계 — 개 참조 유전체를 실제로 받아 본다	6·7강
40~60분	2단계 — FASTQ 를 만들고 열어 본다	8·9강
60~85분	3단계 — 정렬 · 중복 · 커버리지	11~17강
85~90분	점검표 일곱 줄 확인	이 강

뒤 90분은 II 권이 맡습니다.

🤔 확인 문제

2단계에서 변이를 '심어 두고' 4단계에서 '다시 찾는' 방식의 값어치는?

🤔 확인 문제

3단계가 끝난 뒤 반드시 확인해야 할 것으로 가장 중요한 것은?

해설 SM: 이 없으면 다음 권의 뺄셈이 종양과 정상을 구별하지 못합니다. 그리고 정렬을 다시 돌리는 것은 실제 WGS 에서 몇 시간짜리 일입니다. 다른 항목들(.fai·.bai·SO:coordinate)은 없으면 바로 오류가 나므로 오히려 안전합니다.