파트 3: 코호트에 대한 조인트 호출¶
AI 지원 번역 - 자세히 알아보기 및 개선 제안
이전에는 각 샘플의 데이터를 독립적으로 처리하는 샘플별 변이 호출 파이프라인을 구축했습니다. 이제 파트 1(다중 샘플 사용 사례)에서 설명한 대로 조인트 변이 호출을 구현하도록 파이프라인을 확장하겠습니다.
이 섹션을 시작하는 방법
이 과정의 이 섹션은 파트 1: 방법 개요, 파트 2: 샘플별 변이 호출을 완료하고 작동하는 genomics.nf 파이프라인이 있다고 가정합니다.
파트 2를 완료하지 않았거나 이 파트를 새로 시작하려면 파트 2 해결책을 시작점으로 사용할 수 있습니다.
nf4-science/genomics/ 디렉토리 내에서 다음 명령을 실행하세요:
cp solutions/part2/genomics-2.nf genomics.nf
cp solutions/part2/nextflow.config .
cp solutions/part2/modules/* modules/
이렇게 하면 완전한 샘플별 변이 호출 워크플로우를 얻을 수 있습니다. 다음 명령을 실행하여 성공적으로 실행되는지 테스트할 수 있습니다:
과제¶
이 과정의 이번 파트에서는 다음을 수행하도록 워크플로우를 확장할 것입니다:
- Samtools를 사용하여 각 BAM 입력 파일에 대한 인덱스 파일 생성
- 각 BAM 입력 파일에 대해 GATK HaplotypeCaller를 실행하여 샘플별 게놈 변이 호출의 GVCF 생성
- 모든 GVCF를 수집하고 GenomicsDB 데이터 저장소로 결합
- 결합된 GVCF 데이터 저장소에 대해 조인트 유전자형 분석을 실행하여 코호트 수준의 VCF 생성
이는 파트 1: 방법 개요(다중 샘플 사용 사례를 다루는 두 번째 섹션)에서 설명한 방법을 구현하며, 파트 2: 샘플별 변이 호출에서 생성한 워크플로우를 직접 기반으로 합니다.
학습 계획¶
이를 두 단계로 나누었습니다:
- 샘플별 변이 호출 단계를 수정하여 GVCF를 생성합니다. 프로세스 명령과 출력을 업데이트하는 것을 다룹니다.
- 샘플별 GVCF를 결합하고 유전자형을 분석하는 조인트 유전자형 분석 단계를 추가합니다.
collect()연산자, 명령줄 구성을 위한 Groovy 클로저, 그리고 다중 명령 프로세스를 소개합니다.
이는 파트 1: 방법 개요의 두 번째 섹션에서 컨테이너에서 수동으로 실행했던 단계를 자동화합니다.
팁
올바른 작업 디렉토리에 있는지 확인하세요:
cd /workspaces/training/nf4-science/genomics
1. 샘플별 변이 호출 단계를 수정하여 GVCF 생성¶
파트 2의 파이프라인은 VCF 파일을 생성하지만, 조인트 호출에는 GVCF 파일이 필요합니다. GVCF 변이 호출 모드를 활성화하고 출력 파일 확장자를 업데이트해야 합니다.
파트 1의 GVCF 변이 호출 명령을 상기해 보세요:
gatk HaplotypeCaller \
-R /data/ref/ref.fasta \
-I /data/bam/reads_mother.bam \
-O reads_mother.g.vcf \
-L /data/ref/intervals.bed \
-ERC GVCF
파트 2에서 래핑한 기본 HaplotypeCaller 명령과 비교하면, 차이점은 -ERC GVCF 매개변수와 .g.vcf 출력 확장자입니다.
1.1. HaplotypeCaller가 GVCF를 내보내도록 지시하고 출력 확장자 업데이트¶
modules/gatk_haplotypecaller.nf 모듈 파일을 열어 두 가지를 변경하세요:
- GATK HaplotypeCaller 명령에
-ERC GVCF매개변수를 추가합니다; - GATK 규칙에 따라 출력 파일 경로를 해당하는
.g.vcf확장자를 사용하도록 업데이트합니다.
-ERC GVCF를 추가할 때 이전 줄 끝에 백슬래시(\)를 추가해야 합니다.
또한 새 파일 확장자와 일치하도록 출력 블록을 업데이트해야 합니다.
명령 출력을 .vcf에서 .g.vcf로 변경했으므로, 프로세스 output: 블록도 동일한 변경 사항을 반영해야 합니다.
1.2. 프로세스 출력 블록에서 출력 파일 확장자 업데이트¶
또한 새 GVCF 출력을 반영하도록 워크플로우의 게시 및 출력 구성을 업데이트해야 합니다.
1.3. 새 GVCF 출력에 대한 게시 대상 업데이트¶
이제 VCF 대신 GVCF를 생성하므로, 더 설명적인 이름을 사용하도록 워크플로우의 publish: 섹션을 업데이트해야 합니다.
또한 명확성을 위해 GVCF 파일을 자체 하위 디렉토리로 구성하겠습니다.
이제 일치하도록 출력 블록을 업데이트하세요.
1.4. 새 디렉토리 구조에 맞게 출력 블록 업데이트¶
또한 GVCF 파일을 gvcf 하위 디렉토리에 넣도록 output 블록을 업데이트해야 합니다.
모듈, 게시 대상, 출력 블록이 모두 업데이트되었으므로 변경 사항을 테스트할 수 있습니다.
1.5. 파이프라인 실행¶
변경 사항이 작동하는지 확인하기 위해 워크플로우를 실행하세요.
명령 출력
Nextflow 출력은 이전과 동일하게 보이지만, .g.vcf 파일과 해당 인덱스 파일이 이제 하위 디렉토리로 구성되어 있습니다.
디렉토리 내용 (심볼릭 링크 축약됨)
results/
├── gvcf/
│ ├── reads_father.bam.g.vcf -> */27/0d7eb9*/reads_father.bam.g.vcf
│ ├── reads_father.bam.g.vcf.idx -> */27/0d7eb9*/reads_father.bam.g.vcf.idx
│ ├── reads_mother.bam.g.vcf -> */e4/4ed55e*/reads_mother.bam.g.vcf
│ ├── reads_mother.bam.g.vcf.idx -> */e4/4ed55e*/reads_mother.bam.g.vcf.idx
│ ├── reads_son.bam.g.vcf -> */08/e95962*/reads_son.bam.g.vcf
│ └── reads_son.bam.g.vcf.idx -> */08/e95962*/reads_son.bam.g.vcf.idx
└── indexed_bam/
├── reads_father.bam -> */9a/c7a873*/reads_father.bam
├── reads_father.bam.bai -> */9a/c7a873*/reads_father.bam.bai
├── reads_mother.bam -> */f1/8d8486*/reads_mother.bam
├── reads_mother.bam.bai -> */f1/8d8486*/reads_mother.bam.bai
├── reads_son.bam -> */cc/fbc705*/reads_son.bam
└── reads_son.bam.bai -> */cc/fbc705*/reads_son.bam.bai
GVCF 파일 중 하나를 열어 스크롤하면 GATK HaplotypeCaller가 요청한 대로 GVCF 파일을 생성했는지 확인할 수 있습니다.
핵심 정리¶
도구 명령의 출력 파일명을 변경하면, 프로세스 output: 블록과 게시/출력 구성을 일치하도록 업데이트해야 합니다.
다음 단계¶
채널의 내용을 수집하고 단일 입력으로 다음 프로세스에 전달하는 방법을 학습합니다.
2. 조인트 유전자형 분석 단계 추가¶
이제 샘플별 GVCF를 수집하고, GenomicsDB 데이터 저장소로 결합한 다음, 조인트 유전자형 분석을 실행하여 코호트 수준의 VCF를 생성해야 합니다.
파트 1에서 다룬 것처럼, 이것은 두 도구 작업입니다: GenomicsDBImport가 GVCF를 결합하고, GenotypeGVCFs가 최종 변이 호출을 생성합니다.
두 도구를 GATK_JOINTGENOTYPING이라는 단일 프로세스로 래핑하겠습니다.
파트 1의 두 명령을 상기해 보세요:
gatk GenomicsDBImport \
-V reads_mother.g.vcf \
-V reads_father.g.vcf \
-V reads_son.g.vcf \
-L /data/ref/intervals.bed \
--genomicsdb-workspace-path family_trio_gdb
첫 번째 명령은 샘플별 GVCF와 인터벌 파일을 받아 GenomicsDB 데이터 저장소를 생성합니다.
두 번째 명령은 해당 데이터 저장소와 참조 게놈을 받아 최종 코호트 수준의 VCF를 생성합니다.
컨테이너 URI는 HaplotypeCaller와 동일합니다: community.wave.seqera.io/library/gatk4:4.5.0.0--730ee8817e436867.
2.1. 입력 설정¶
조인트 유전자형 분석 프로세스에는 아직 없는 두 종류의 입력이 필요합니다: 임의의 코호트 이름과 모든 샘플의 수집된 GVCF 출력을 함께 묶은 것입니다.
2.1.1. cohort_name 매개변수 추가¶
코호트에 대한 임의의 이름을 제공해야 합니다.
교육 시리즈의 후반부에서 이러한 종류의 작업에 샘플 메타데이터를 사용하는 방법을 배우겠지만, 지금은 params를 사용하여 CLI 매개변수를 선언합니다.
이 매개변수를 최종 출력 파일 이름 지정에 사용하겠습니다.
2.1.2. 테스트 프로필에 cohort_name의 기본값 추가¶
또한 테스트 프로필에 cohort_name 매개변수의 기본값을 추가합니다:
| nextflow.config | |
|---|---|
다음으로, 샘플별 출력을 수집하여 함께 처리할 수 있도록 해야 합니다.
2.1.3. 샘플 전체에서 HaplotypeCaller 출력 수집¶
GATK_HAPLOTYPECALLER의 출력 채널을 새 프로세스에 직접 연결하면, Nextflow는 각 샘플 GVCF에 대해 프로세스를 개별적으로 호출합니다.
세 개의 GVCF(및 해당 인덱스 파일)를 모두 묶어서 Nextflow가 단일 프로세스 호출에 모두 함께 전달하도록 하려고 합니다.
collect() 채널 연산자를 사용하여 이를 수행할 수 있습니다.
GATK_HAPLOTYPECALLER 호출 직후 workflow 본문에 다음 줄을 추가하세요:
이를 분석하면:
.out속성을 사용하여GATK_HAPLOTYPECALLER의 출력 채널을 가져옵니다.- 섹션 1에서
emit:을 사용하여 출력의 이름을 지정했기 때문에,.vcf로 GVCF를 선택하고.idx로 인덱스 파일을 선택할 수 있습니다. 이름이 지정된 출력이 없으면.out[0]과.out[1]을 사용해야 합니다. collect()연산자는 모든 파일을 단일 요소로 묶으므로,all_gvcfs_ch에는 세 개의 GVCF가 모두 함께 포함되고,all_idxs_ch에는 세 개의 인덱스 파일이 모두 함께 포함됩니다.
Nextflow가 실행을 위해 모든 입력 파일을 함께 스테이징하므로 인덱스 파일이 GVCF와 함께 존재하게 되어, GVCF와 해당 인덱스 파일을 (튜플로 함께 유지하는 것과 반대로) 별도로 수집할 수 있습니다.
팁
view() 연산자를 사용하여 채널 연산자를 적용하기 전후에 채널의 내용을 검사할 수 있습니다.
2.2. 조인트 유전자형 분석 프로세스 작성 및 워크플로우에서 호출¶
파트 2에서 사용한 것과 동일한 패턴을 따라, 모듈 파일에 프로세스 정의를 작성하고, 워크플로우로 가져온 다음, 방금 준비한 입력에 대해 호출하겠습니다.
2.2.1. 각 GVCF에 -V 인수를 제공하는 문자열 구성¶
프로세스 정의를 채우기 시작하기 전에 해결해야 할 한 가지가 있습니다.
GenomicsDBImport 명령은 각 GVCF 파일에 대해 별도의 -V 인수를 기대합니다:
gatk GenomicsDBImport \
-V reads_mother.bam.g.vcf \
-V reads_father.bam.g.vcf \
-V reads_son.bam.g.vcf \
...
-V ${all_gvcfs_ch}로 작성하면, Nextflow는 단순히 파일명을 연결하고 명령의 해당 부분은 다음과 같이 보일 것입니다:
하지만 문자열이 다음과 같이 보여야 합니다:
중요한 것은, 수집된 채널에 있는 파일이 무엇이든 이 문자열을 동적으로 구성해야 한다는 것입니다. Nextflow(Groovy를 통해)는 이를 수행하는 간결한 방법을 제공합니다:
이를 분석하면:
all_gvcfs.collect { gvcf -> "-V ${gvcf}" }는 각 파일 경로를 반복하고-V를 앞에 추가하여["-V A.g.vcf", "-V B.g.vcf", "-V C.g.vcf"]를 생성합니다..join(' ')은 공백으로 연결합니다:"-V A.g.vcf -V B.g.vcf -V C.g.vcf".- 결과는 로컬 변수
gvcfs_line(def로 정의됨)에 할당되며, 명령 템플릿에 보간할 수 있습니다.
이 줄은 명령 템플릿 전에 프로세스의 script: 블록 내부에 들어갑니다.
script:와 명령 템플릿의 여는 """ 사이에 임의의 Groovy 코드를 배치할 수 있습니다.
그러면 프로세스의 script: 블록에서 전체 문자열을 gvcfs_line으로 참조할 수 있습니다.
2.2.2. 조인트 유전자형 분석 프로세스의 모듈 작성¶
다음으로, 전체 프로세스 작성을 다룰 수 있습니다.
modules/gatk_jointgenotyping.nf를 열고 프로세스 정의의 개요를 살펴보세요.
위에 제공된 정보를 사용하여 프로세스 정의를 작성한 다음, 아래 "후" 탭의 해결책과 비교하여 작업을 확인하세요.
여기서 주목할 만한 몇 가지가 있습니다.
이전과 마찬가지로, 명령이 직접 참조하지 않더라도 여러 입력이 나열되어 있습니다: all_idxs, ref_index, ref_dict.
이들을 나열하면 Nextflow가 명령에 나타나는 파일과 함께 작업 디렉토리에 이러한 파일을 스테이징하도록 보장하며, GATK는 명명 규칙에 따라 이를 찾을 것으로 예상합니다.
gvcfs_line 변수는 위에서 설명한 Groovy 클로저를 사용하여 GenomicsDBImport에 대한 -V 인수를 구성합니다.
이 프로세스는 터미널에서와 마찬가지로 두 명령을 순차적으로 실행합니다.
GenomicsDBImport는 샘플별 GVCF를 데이터 저장소로 결합한 다음, GenotypeGVCFs는 해당 데이터 저장소를 읽고 최종 코호트 수준의 VCF를 생성합니다.
GenomicsDB 데이터 저장소(${cohort_name}_gdb)는 프로세스 내에서만 사용되는 중간 산출물입니다; 출력 블록에는 나타나지 않습니다.
이를 완료하면 프로세스를 사용할 준비가 됩니다. 워크플로우에서 사용하려면 모듈을 가져오고 프로세스 호출을 추가해야 합니다.
2.2.3. 모듈 가져오기¶
기존 import 문 아래에 genomics.nf에 import 문을 추가하세요:
이제 프로세스를 워크플로우 범위에서 사용할 수 있습니다.
2.2.4. 프로세스 호출 추가¶
collect() 줄 다음에 워크플로우 본문에 GATK_JOINTGENOTYPING 호출을 추가하세요:
| genomics.nf | |
|---|---|
이제 프로세스가 완전히 연결되었습니다. 다음으로 출력이 게시되는 방식을 구성합니다.
2.3. 출력 처리 구성¶
조인트 VCF 출력을 게시해야 합니다. 조인트 유전자형 분석 결과에 대한 게시 대상과 출력 블록 항목을 추가하세요.
2.3.1. 조인트 VCF에 대한 게시 대상 추가¶
워크플로우의 publish: 섹션에 조인트 VCF와 해당 인덱스를 추가하세요:
이제 일치하도록 출력 블록을 업데이트하세요.
2.3.2. 조인트 VCF에 대한 출력 블록 항목 추가¶
조인트 VCF 파일에 대한 항목을 추가하세요. 이것이 최종 출력이므로 결과 디렉토리의 루트에 배치하겠습니다.
프로세스, 게시 대상, 출력 블록이 모두 준비되었으므로 완전한 워크플로우를 테스트할 수 있습니다.
2.4. 워크플로우 실행¶
모든 것이 작동하는지 확인하기 위해 워크플로우를 실행하세요.
명령 출력
처음 두 단계는 이전 실행에서 캐시되었으며, 새로운 GATK_JOINTGENOTYPING 단계는 세 샘플 모두의 수집된 입력에 대해 한 번 실행됩니다.
최종 출력 파일인 family_trio.joint.vcf(및 해당 인덱스)는 결과 디렉토리에 있습니다.
디렉토리 내용 (심볼릭 링크 축약됨)
results/
├── family_trio.joint.vcf -> */a6/7cc8ed*/family_trio.joint.vcf
├── family_trio.joint.vcf.idx -> */a6/7cc8ed*/family_trio.joint.vcf.idx
├── gvcf/
│ ├── reads_father.bam.g.vcf -> */27/0d7eb9*/reads_father.bam.g.vcf
│ ├── reads_father.bam.g.vcf.idx -> */27/0d7eb9*/reads_father.bam.g.vcf.idx
│ ├── reads_mother.bam.g.vcf -> */e4/4ed55e*/reads_mother.bam.g.vcf
│ ├── reads_mother.bam.g.vcf.idx -> */e4/4ed55e*/reads_mother.bam.g.vcf.idx
│ ├── reads_son.bam.g.vcf -> */08/e95962*/reads_son.bam.g.vcf
│ └── reads_son.bam.g.vcf.idx -> */08/e95962*/reads_son.bam.g.vcf.idx
└── indexed_bam/
├── reads_father.bam -> */9a/c7a873*/reads_father.bam
├── reads_father.bam.bai -> */9a/c7a873*/reads_father.bam.bai
├── reads_mother.bam -> */f1/8d8486*/reads_mother.bam
├── reads_mother.bam.bai -> */f1/8d8486*/reads_mother.bam.bai
├── reads_son.bam -> */cc/fbc705*/reads_son.bam
└── reads_son.bam.bai -> */cc/fbc705*/reads_son.bam.bai
조인트 VCF 파일을 열면 워크플로우가 예상된 변이 호출을 생성했는지 확인할 수 있습니다.
이제 자동화되고 완전히 재현 가능한 조인트 변이 호출 워크플로우가 완성되었습니다!
참고
제공된 데이터 파일은 염색체 20의 아주 작은 부분만 다룬다는 점을 명심하세요. 실제 변이 호출 세트의 크기는 수백만 개의 변이로 계산됩니다. 그래서 교육 목적으로는 아주 작은 데이터 하위 집합만 사용합니다!
핵심 정리¶
채널에서 출력을 수집하고 단일 입력으로 다른 프로세스에 묶는 방법을 알게 되었습니다. 또한 Groovy 클로저를 사용하여 명령줄을 구성하는 방법과 단일 프로세스에서 여러 명령을 실행하는 방법도 알게 되었습니다.
다음 단계¶
스스로를 크게 칭찬하세요! Nextflow for Genomics 과정을 완료했습니다.
학습한 내용을 검토하고 다음 단계를 알아보려면 최종 과정 요약으로 이동하세요.