2026-03-11 ⋯ Bioinformatics Engineering #3 변이 Annotation
Annotation이란? VCF 파일에는 "chr1:12345에서 A가 G로 바뀌었다"는 정보만 있습니다. Annotation은 이 변이에 생물학적 의미를 부여하는 작업입니다. 어떤 유전자에 위치하는지, 엑손인지 인트론인지, 아미노산이 바뀌는지, 그 변화가 단백질 기능에 해로운지, 다른 환자에서도 발견된 변이인지 등의 정보를 추가합니다. Annotation 도구 선택하기 가장 널리 사용되는 도구는 VEP(Variant Effect Predictor), ANNOVAR, SnpEff입니다. VEP는 Ensembl에서 개발하여 유전자 모델이 풍부하고 업데이트가 빠릅니다. ANNOVAR는 다양한 데이터베이스 통합이 쉽고 속도가 빠릅니다. SnpEff는 설치가 간단하고 메모리 사용이 적습니다. 국립암센터처럼 대규모 데이터를 다루는 기관에서는 VEP가 표준으로 많이 사용됩니다. 캐시 데이터 설치하기 Annotation 도구는 유전자 정보, 전사체 정보, 단백질 서열 등을 담은 캐시 데이터가 필요합니다. VEP의 경우 homo_sapiens_vep_110_GRCh38.tar.gz 같은 파일을 다운로드하여 지정된 경로에 설치합니다. 캐시에는 Ensembl/RefSeq 유전자 모델, 규제 요소 정보, 단백질 도메인 정보 등이 포함됩니다. 온라인 모드로도 실행 가능하지만, 대량 데이터 처리 시 캐시를 사용해야 속도가 확보됩니다. 기본 Annotation 실행하기 VEP를 예로 들면, vep -i input.vcf -o output.vcf --cache --dir_cache 캐시경로 --fasta 참조유전체.fa --vcf 형식으로 실행합니다. 도구가 각 변이의 좌표를 유전자 모델과 대조하여, 해당 변이가 어떤 전사체(transcript)의 어떤 위치에 해당하는지 계산합니다. 출력 파일의 INFO 필드에 CSQ(Consequence) 정보가 추가됩니다. Consequence 이해하기 Annotation 결과에서 가장 중요한 것은 변이의 consequence(결과/영향)입니다. missense_variant는 아미노산이 다른 아미노산으로 바뀌는 것, nonsense_variant(stop_gained)는 조기 종결 코돈이 생기는 것, frameshift_variant는 삽입/결실로 읽기 틀이 바뀌는 것, synonymous_variant는 아미노산 변화 없이 코돈만 바뀌는 것입니다. 일반적으로 frameshift, nonsense, splice_site 변이가 단백질 기능에 가장 큰 영향을 줍니다. 추가 Annotation 플러그인 적용하기 기본 annotation 외에 플러그인으로 추가 정보를 붙일 수 있습니다. SIFT와 PolyPhen은 missense 변이가 단백질 기능에 해로운지 예측합니다. CADD는 변이의 유해성을 종합 점수로 제공합니다. dbNSFP는 수십 개의 예측 알고리즘 결과를 한 번에 추가합니다. VEP 실행 시 --plugin CADD,파일경로 형식으로 지정하면 해당 정보가 출력에 포함됩니다. 공개 데이터베이스 정보 추가하기 임상적 해석을 위해 공개 데이터베이스 정보를 추가합니다.
ClinVar는 변이의 임상적 의미(pathogenic, likely_pathogenic, uncertain_significance, likely_benign, benign)를 제공합니다. 특정 변이가 유전 질환과 연관되어 있는지 알 수 있습니다.
COSMIC(Catalogue of Somatic Mutations in Cancer)은 해당 변이가 다른 암 환자에서 몇 번 발견되었는지, 어떤 암종에서 흔한지 정보를 제공합니다. 빈번하게 발견되는 변이는 driver mutation일 가능성이 높습니다.
gnomAD(Genome Aggregation Database)는 일반 인구에서의 변이 빈도를 제공합니다. 인구 집단의 1% 이상에서 발견되는 변이는 보통 질병 원인이 아닌 다형성(polymorphism)으로 간주합니다. 결과 파싱 및 정리하기 Annotation 출력은 복잡한 형식으로 되어 있습니다. VEP의 경우 CSQ 필드에 파이프(|)로 구분된 수십 개의 값이 들어 있습니다. 이를 분석에 사용하려면 파싱하여 별도의 컬럼으로 분리해야 합니다. 또한 한 변이가 여러 전사체에 영향을 줄 수 있으므로, canonical transcript를 기준으로 대표 annotation을 선택하거나, 가장 심각한 consequence를 가진 것을 선택하는 규칙을 정해야 합니다. 품질 검증하기 최종 annotation 결과를 검증합니다. 알려진 병원성 변이(예: BRCA1의 특정 변이)가 제대로 ClinVar pathogenic으로 annotation되었는지 확인합니다. 유전자명(Hugo Symbol)이 최신 명명법을 따르는지, 좌표와 REF/ALT가 원본과 일치하는지 점검합니다. 이 검증을 통해 파이프라인의 신뢰성을 확보합니다.각각 이런식으로 전체 흐름 요약도 적어줘 정리 VCF 파일에서 시작하여, VEP나 ANNOVAR 같은 도구와 해당 참조유전체 버전의 캐시 데이터를 준비합니다. 기본 annotation을 실행하여 각 변이의 유전자 위치, 전사체 영향, consequence(missense, frameshift 등)를 파악하고, SIFT/PolyPhen/CADD 같은 플러그인으로 기능적 영향을 예측합니다. ClinVar로 임상적 의미를, COSMIC으로 암 연관성을, gnomAD로 인구 집단 빈도를 추가합니다. 복잡한 출력을 파싱하여 분석 가능한 형태로 정리하고, 알려진 변이로 검증하여 파이프라인의 신뢰성을 확보합니다. 이렇게 annotation된 데이터는 변이의 생물학적·임상적 의미를 해석하고, 치료 타깃 발굴이나 환자 층화에 활용됩니다.
2026-03-11 ⋯ Bioinformatics Engineering #2 Lift-over
참조유전체 Lift-over 유전체 데이터는 참조유전체(reference genome)를 기준으로 좌표가 기록됩니다. 그런데 참조유전체는 계속 업데이트됩니다. GRCh37(hg19)에서 chr1:12345 위치에 있던 변이가 GRCh38(hg38)에서는 chr1:12456일 수 있습니다. 참조유전체 버전이 개선되면서 일부 서열이 추가되거나 수정되기 때문입니다. Lift-over는 한 참조유전체 버전의 좌표를 다른 버전의 좌표로 변환하는 작업입니다. 왜 Lift-over가 필요한지 이해하기 다기관 연구에서 A병원은 GRCh37로, B병원은 GRCh38로 데이터를 생성했다면, 두 데이터를 통합하려면 좌표 기준을 맞춰야 합니다. 또한 ClinVar, COSMIC 같은 공개 데이터베이스가 GRCh38 기준이라면, GRCh37 데이터에 annotation을 추가하기 전에 좌표를 변환해야 합니다. 좌표가 맞지 않으면 동일한 변이를 서로 다른 변이로 인식하거나, annotation이 잘못 매핑됩니다. Chain 파일 이해하기 Lift-over의 핵심은 chain 파일입니다. Chain 파일은 두 참조유전체 버전 간의 서열 대응 관계를 기록한 매핑 파일입니다. 예를 들어 "GRCh37의 chr1:10000-20000 구간은 GRCh38의 chr1:10500-20500에 대응한다"는 정보가 담겨 있습니다. UCSC Genome Browser에서 hg19ToHg38.over.chain.gz 같은 형식으로 제공됩니다. 반대 방향(GRCh38→GRCh37)이 필요하면 hg38ToHg19.over.chain.gz를 사용합니다. 입력 파일 준비하기 Lift-over할 파일을 준비합니다. VCF 파일이라면 염색체, 위치, REF, ALT 정보가 포함되어 있습니다. BED 파일이라면 염색체, 시작 위치, 끝 위치가 탭으로 구분되어 있습니다. 중요한 것은 입력 파일의 현재 참조유전체 버전을 정확히 아는 것입니다. VCF 헤더의 라인이나 라인에서 확인할 수 있습니다. Lift-over 도구 실행하기 가장 많이 사용되는 도구는 UCSC liftOver와 CrossMap입니다. UCSC liftOver는 BED 형식 입력에 최적화되어 있고, CrossMap은 VCF, BAM, BED 등 다양한 형식을 지원합니다. CrossMap을 예로 들면, CrossMap.py vcf chain파일 입력.vcf 타겟참조유전체.fa 출력.vcf 형식으로 실행합니다. 도구가 chain 파일을 참조하여 각 좌표를 새 버전으로 변환합니다. 변환 실패 케이스 처리하기 모든 좌표가 성공적으로 변환되지는 않습니다. 일부 영역은 새 참조유전체에서 삭제되었거나 크게 재배열되어 매핑이 불가능합니다. 이런 변이들은 "unmapped" 파일로 분리됩니다. 연구 목적에 따라 이 변이들을 제외하거나, 수동으로 검토하거나, 원본 좌표를 별도 컬럼에 보존할 수 있습니다. 일반적으로 전체 변이의 1-5% 정도가 매핑 실패합니다. 변환 결과 검증하기 Lift-over 후에는 결과를 검증해야 합니다. 몇 개의 변이를 무작위로 선택하여 UCSC Genome Browser나 Ensembl에서 양쪽 버전의 좌표를 직접 확인합니다. 또한 변환 전후의 REF 염기가 해당 참조유전체와 일치하는지 확인합니다. GRCh38 좌표인데 REF 염기가 GRCh38 서열과 다르다면 변환 오류입니다. 정리 원본 파일의 현재 참조유전체 버전을 확인하고, 목표 버전으로의 변환이 필요한지 판단합니다. UCSC에서 적절한 chain 파일을 다운로드하고, CrossMap이나 liftOver 도구로 좌표 변환을 실행합니다. 변환 실패한 변이들은 unmapped 파일로 분리하여 별도 검토하고, 성공한 변이들은 새 참조유전체 서열과 REF 염기가 일치하는지 검증합니다. 이렇게 좌표가 통일된 데이터는 다기관 데이터 통합이나 공개 데이터베이스 annotation에 활용할 수 있게 됩니다.
2026-03-11 ⋯ Bioinformatics Engineering #1 VCF to MAF 변환
VCF to MAF 변환 VCF(Variant Call Format)는 개별 샘플에서 발견된 유전체 변이를 기록하는 파일 형식입니다. MAF(Mutation Annotation Format)는 여러 샘플의 변이를 통합하고, 각 변이에 유전자명, 단백질 변화, 임상적 의미 등의 annotation을 추가한 형식입니다. 쉽게 말해 VCF가 "어디에 변이가 있다"를 기록한다면, MAF는 "그 변이가 어떤 유전자에 있고 어떤 영향을 주는가"까지 설명합니다. VCF 파일 VCF 파일을 열면 헤더와 데이터 라인으로 구성되어 있습니다. 헤더는 시작하며 파일의 메타정보를 담고 있고, 데이터 라인은 각 변이의 염색체 위치(CHROM), 좌표(POS), 참조 염기(REF), 대체 염기(ALT), 품질 점수(QUAL) 등을 탭으로 구분하여 기록합니다. 예를 들어 chr1의 12345번 위치에서 A가 G로 바뀌었다면, 그 정보가 한 줄에 담깁니다. 참조유전체 버전 확인하기 VCF 파일이 어떤 참조유전체를 기준으로 만들어졌는지 확인해야 합니다. 현재 주로 사용되는 버전은 GRCh37(hg19)과 GRCh38(hg38)입니다. 만약 VCF가 GRCh37 기준인데 annotation 데이터베이스가 GRCh38을 사용한다면, 먼저 좌표를 변환(Lift-over)해야 합니다. 이 과정에서 UCSC의 liftOver 도구나 CrossMap을 사용합니다. Annotation 수행하기 변이의 생물학적 의미를 알기 위해 annotation을 수행합니다. 가장 많이 사용되는 도구는 VEP(Variant Effect Predictor)입니다. VEP는 각 변이가 어떤 유전자의 어떤 부위(exon, intron, UTR 등)에 위치하는지, 아미노산 변화가 발생하는지, 그 변화가 단백질 기능에 어떤 영향을 주는지 예측합니다. VEP를 실행하면 원본 VCF에 annotation 정보가 추가된 파일이 생성됩니다. vcf2maf 도구로 변환하기 vcf2maf는 VEP로 annotation된 VCF를 MAF 형식으로 변환하는 표준 도구입니다. 이 도구는 VCF의 각 변이를 MAF의 필수 컬럼(Hugo_Symbol, Chromosome, Start_Position, End_Position, Variant_Classification, Variant_Type, Reference_Allele, Tumor_Seq_Allele2 등)에 맞게 재구성합니다. 실행 시 참조유전체 FASTA 파일과 VEP 캐시 경로를 지정해야 합니다. 체세포 vs 생식세포 변이 분류하기 암 연구에서는 체세포 변이(somatic mutation)와 생식세포 변이(germline variant)를 구분하는 것이 중요합니다. 체세포 변이는 암세포에서만 발생한 후천적 변이이고, 생식세포 변이는 환자가 태어날 때부터 가진 변이입니다. 일반적으로 종양 샘플과 정상 샘플을 비교하여, 종양에만 있는 변이를 체세포로, 둘 다에 있는 변이를 생식세포로 분류합니다. MAF 파일의 Mutation_Status 컬럼에 이 정보가 기록됩니다. Public DB annotation 추가하기 MAF 파일에 임상적 의미를 더하기 위해 공개 데이터베이스 정보를 추가합니다. ClinVar는 변이의 병원성(pathogenic, benign 등)을, COSMIC은 해당 변이가 다른 암 환자에서 얼마나 자주 발견되었는지를, gnomAD는 일반 인구에서의 빈도를 제공합니다. 이 정보들이 추가되면 "이 변이가 암 발생에 얼마나 중요한가"를 판단할 수 있게 됩니다. 품질 검증 및 정제하기 최종 MAF 파일의 품질을 검증합니다. 필수 컬럼에 결측값이 없는지, 좌표가 올바른 범위 내에 있는지, 변이 분류(Variant_Classification)가 표준 용어를 사용하는지 확인합니다. 또한 분석 목적에 따라 특정 변이 유형(예: synonymous 변이 제외, 특정 VAF 이하 제외)을 필터링하여 최종 분석용 MAF를 생성합니다. 정리 원본 VCF에서 시작하여, 참조유전체 버전을 맞추고, VEP로 유전자/단백질 영향을 annotation하고, vcf2maf로 형식을 변환하고, 체세포/생식세포를 분류하고, 공개 DB 정보를 추가하고, 품질 검증을 거쳐 최종 MAF가 완성됩니다. 이 MAF 파일은 여러 환자의 변이를 통합 분석하거나, AI 모델 학습에 사용하거나, 임상적 의사결정을 지원하는 데 활용됩니다.
2026-01-26 ⋯ RNA-seq 전처리 #3 Read counting
read counting alignment output file인 accepted_hits.bam 파일을 사용해서 gene-level count 테이블을 생성하기. HTSeq-count 도구를 이용했고 run_htseq_count.sh 스크립트로 작업. 실행하면 샘플에 대해 Count/count_150-2.txt와 같은 gene-level count 파일이 생성된다. HTSeq-count 오류 트러블슈팅 accepted_hits.bam 파일을 HTSeq-count에 입력했는데, 아래와 같은 경고/에러가 발생하며 count 생성이 중단됨 HTSeq-count는 paired-end 모드에서 read1/read2가 입력 스트림에서 연속된 줄로 들어오는 것을 전제로 mate를 매칭함. accepted_hits.bam은 일반적으로 좌표 기준 정렬(coordinate-sorted) 상태이기 때문에, mate 쌍이 파일 내에서 떨어져 존재할 수 있고(서로 다른 위치로 정렬되면 더 흔함), 그 결과 HTSeq가 “인접 라인에서 mate를 못 찾았다”고 판단하며 오류가 발생한다. 입력 BAM을 read name 기준으로 정렬(name-sorted) 한 뒤, HTSeq-count에 전달하도록 변경함 즉, mate끼리 항상 연속된 상태로 들어가도록 보장해서 pair_alignments 에러를 방지함. 정리 - Alignment 결과(accepted_hits.bam)를 입력으로 gene-level read counting 수행
- hg38.refGene.gtf annotation 기준으로 HTSeq-count 적용
- header 포함을 위해 samtools view -H + cat 파이프라인 사용
- BAM을 SAM으로 스트리밍 변환하면서 paired-end 처리 중 mate가 인접하지 않은(coordinate-sorted) BAM 입력으로 인해 pair_alignments 오류가 발생함
- samtools sort -n을 이용해 name-sorted BAM으로 변환 후 HTSeq-count를 수행하여 mate adjacency를 보장하고 오류를 해결.
2026-01-25 ⋯ RNA-seq 전처리 #2 Alignment
Alignment 작업 hg38 참조 유전체를 기준으로 DHT 약물 처리된 18개 샘플의 paired-end RNA-seq 데이터 정렬 작업
정렬 도구로는 TopHat 2.0.13을 사용하였고 대규모 데이터 처리 시간을 고려하여 40개 CPU 스레드를 활용하고 백그라운드에서 쉘과 분리된 상태(disown)로 실행함.
Alignment 실행 로그는 파일로 리다이렉션하여 추후 확인 가능하도록 설정했고 파일 길이가 평균 100M으로 소요시간이 매우 커서 장시간 작업을 위해 터미널 종료 후에도 백그라운드 실행 유지하도록 & + disown으로 작업. Topaht 오류 트러블슈팅 TopHat2는 구형 코드로써 python2만을 지원하고 python3를 사용하는 서버 환경에서는 올바르게 실행되지 않는 문제가 발생. 설치된 TopHat 디렉토리 구조를 확인한 결과, 실제 실행 파일은 아래 경로의 tophat.py임을 확인하였다. 해당 파일을 직접 호출할 수 있도록 .bashrc에 tophatpy라는 alias를 추가했고 alias 적용 후에는 tophat 대신 tophatpy 명령어를 통해 정상적으로 alignment 수행이 가능했다. 정리 hg38 기준으로 DHT 처리 RNA-seq 샘플 18개를 TopHat2로 정렬 작업. 대규모 연산을 위해 멀티스레드(40 threads) 및 백그라운드 실행 구성했다. python2 의존성으로 인한 실행 오류를 alias 설정으로 해결했다.
2026-01-25 ⋯ RNA-seq 전처리 #1 Quality control
quality control DHT 약물 처리된 18개 샘플의 paired-end RNA-seq 데이터를 quality control 수행. FastQC 도구를 사용했고 run_fastqc.sh 스크립트 만들어서 작업. 품질 레포트 확인 결과 모든 샘플에 공통적으로 read 5’ 말단 1-10bp 구간의 염기 조성이 균일하지 않음. 참고 논문에 따라 1-10bp에 atgc가 서열 내용이 무작위로 분포하는 부분을 Trimming 작업. trim_fastq.sh 스크립트 만들어서 작업함. 정리 18개 샘플의 paired-end RNA-seq 데이터에 대해서 Fastqc를 사용해서 품질 평가를 수행. 품질 레포트 결과에 따라 5′ 말단 10bp trimming을 수행하여 alignment 입력 RNA-seq의 일관성을 확보했다.
2025-07-28 ⋯ MAFFT 작업 #2 MAFFT 실행
1. Objective Influenza의 Reference squence는 길이가 fix되어있지만, - 각 sequence는 삽입/탈락 mutation이 일어남에 따라 모두 길이가 같지 않다. - 이 길이를 맞춰주는 padding을 하기 위해 MAFFT를 이용해 정렬(Multiple Sequence Alignment)한다. 2. MAFFT 실행 bash script
2025-07-28 ⋯ MAFFT 작업 #1 Fasta 파일 전처리
1. Load package 2. Objective Influenza type A의 H1N1 strain의 fasta 파일을 확인해보면? `>`로 시작하는 행에 해당 시퀀스의 메타데이터가 있고 다음 `>`로 시작하는 행 이전까지 해당 시퀀스 정보가 있다. - `>`로 시작하는 행을 |로 분리했을때 제일 마지막값에 유전자 정보가 있다. Raw 데이터는 아래와 같이 구성돼있는데 유전자가 섞여있는 fasta 파일을 분석하기 편하게 유전자별로 분리해줄려고한다. 3. Process fasta file 8개 gene segment별로 fasta 파일을 분리함.