[BioInformatics] NGS Data Analysis Workflow - 시퀀싱된 데이터를 파이썬으로 다뤄보자

NGS(Next Generation Sequencing)가 뭐지 ?
NGS Workflow
시퀀싱 완료된 데이터 형태
FASTQ에 대해서
FastQC & Trimming
Alignment 정렬

NGS(Next Generation Sequencing)가 뭐지 ?

차세대 염기서열 분석(NGS; Next Generation Sequencing)이란 기존의 고전적인 염기서열 분석 방식인 Sanger Sequencing 이후 개발된 기술로, 한 번의 실험에서 수백만 개 이상의 DNA 서열을 병렬로 빠르게 읽어낼 수 있는 고속, 고효율 기술입니다.

기존의 방법은 직접 염기서열분석법이고, 분석하고자 하는 부위를 PCR 증폭해야 해서 목표 범위가 넓어질 경우 많은 시간과 노력, 비용이 소모되어 낮은 효율성을 가지고 있습니다. 하지만 NGS는 DNA 가닥을 각각 하나씩 분석함으로써 기존의 직접 염기서열분석법에 비해 매우 빠르고 저렴하게 염기서열이 가능하다는 장점이 있습니다.

NGS는 3가지의 종류가 있습니다.

1. 패널 검사(Targeted Sequencing, Panels)

원하는 유전자 부위(일부 DNA 또는 RNA)를 선별해서 분석합니다.

다중 유전성 패널 검사의 경우, 여러 유전자가 원인이 되는 다양한 질환에서 진단율이 높아 비용 대비 효율성이 높습니다.

유전성 질환, 혈액 종양, 암 질환 등에서 널리 활용되고 있다고 합니다.

2. 엑솜 검사, 전장 엑솜 염기서열분석 (Whole Exome sequencing, WES)

전체 유전체 중 2만여 개에 달하는 유전자의 단백질을 코딩하는 엑손(exon) 부위를 분석하는 것으로 전체 유전체의 1% 정도입니다. 질환과 연관성이 알려진 대부분의 변이는 엑손 부위에 발생하기 때문에 소요 시간 및 분석 비용 대비 효과적인 방법입니다. 패널 검사의 경우 목표한 부위에 포함된 유전자 변이만 검출이 가능한데 반해 엑솜 검사에서는 예측하지 못했던 유전자 변이까지도 발견이 가능합니다.

3. 지놈 검사, 전장 유전체 염기서열분석 (Whole Genome Sequencing, WGS)

한 종의 유전 정보를 저장하는 유전체를 전부 분석합니다.

엑솜 검사에서는 엑손과 매우 인접한 부위의 인트론(intron)에 존재하는 변이까지만 검출이 가능해서 엑손과 멀리 떨어진 부위의 인트론과 조절 부위(regulatory region)에 있는 변이까지 알아보기 위해서는 지놈 검사가 필요해요.

엑솜 검사와 지놈 검사에서는 예상했던 유전자 이외에 다른 유전자에서 예기치 못하게 돌연변이를 발견할 수 있습니다. 그렇기에 환자가 가진 증상만으로 진단하기 어렵거나 원인 유전자가 다양하여 일반적인 패널 검사로 검사가 충분하지 못할 때 유용하게 사용됩니다.

NGS Workflow

시퀀싱된 FASTQ 데이터를 처리하는 프로세스는 아래와 같습니다.

크게 보면 Raw한 FASTQ 데이터를 받아서 품질 검사를 진행합니다(FastQC).

품질이 좋지 않다면 Trimming 이라는 품질 개선 단계를 거치고 STAR 기법을 사용한 Alignment 단계를 따릅니다.

여기까지 마치면 파이썬을 통해서 발현량 등을 분석할 수 있는 데이터가 나오게 됩니다.

시퀀싱 완료된 데이터 형태

- Illumina 장비는 일반적으로 FASTQ 형식으로 데이터 제공
- PacBio 장비는 일반적으로 FASTQ 또는 BAM 파일 제공
- Oxford Nanopore는 전기신호 기반 FAST5 파일을 제공한 후 이를 FASTQ로 변환 가능

FASTQ에 대해서

NCBI에서 제공하는 SRA-tool은 데이터를 쉽게 다운로드하고 관리하는 도구입니다.
FASTQ 형식으로 데이터를 제공하며 데이터 품질 관리 및 데이터 포맷 전환을 제공하는 기본적이고 핵심적인 도구입니다.
아래의 명령어를 통해 필요한 툴을 설치합니다.

# Conda로 설치 (추천)
conda install -c bioconda sra-tools

# Ubuntu
sudo apt-get install sra-toolkit

상적으로 설치가 끝나면 sra-tools를 사용할 수 있습니다.

기본적인 fastq-dump 사용 명령어는 아래와 같습니다.

fastq-dump [옵션] SRR_ID

옵션은 여러가지가 있는데, 가장 많이 사용하는 옵션들에 대해서 살펴보면 다음과 같습니다.

`--split-files` : Paired-end 데이터일 때 Forward와 Reverse read 분리
ex. `fastq-dump --split-files SRR1234567`
`--gzip` : 압축된 gz 형식으로 FASTQ 저장 (용량 절약)
ex. `fastq-dump --gzip SRR1234567`

`fastq-dump --split-files --gzip SRR1234567`

위의 명령어는 SRR1234567 이라는 일반적인 데이터를 다운로드 받을 수 있는 명령어 입니다. 원하는 질병의 SRR은 NCBI SRA 홈페이지(https://www.ncbi.nlm.nih.gov/sra)에서 찾을 수 있습니다.

데이터를 다운받으면 저렇게 두 개의 데이터가 디렉토리에 나타납니다.

이제 데이터를 분석하기 전에 FastQC & Trimming를 통해 다운로드 받은 FASTQ 파일의 품질을 확인하는 것이 좋습니다.

# FastQC 설치 (conda 환경에서)
conda install -c bioconda fastqc

# 품질 확인 진행
fastqc SRR390728.fastq

FastQC & Trimming

FastQC를 실행하면 저렇게 .html 파일과 zip 파일이 생성되는데, 중요한 파일은 .html 파일입니다.

생성된 HTML 파일을 열어보면 아래와 같은 화면이 나오게 됩니다.

FASTQC 리포트를 간단하게 해석하기 위해서 중요한 항목에 대해서만 알아보겠습니다.

보고서에서 중요한 항목:
- Per Base Sequence Quality
    - 각 염기 위치 별 품질(Phred Score)을 확인
    - 품질이 20 이상이면 양호, 20 미만 구간이 많다면 추가 trimming 필요
- Per Sequence Quality Scores
    - 전반적인 데이터 품질 분포 확인 가능
- Sequence Length Distribution
    - 서열 길이의 분포를 통해 trimming이 필요한지 판단 가능
- Overrepresented sequences
    - 중복된 서열이 많은 경우 오염 가능성 의심
- 전반적으로 빨간색 ❌나 노란색 경고⚠️ 항목이 많다면, 품질 저하 가능성이 있으므로 Trimming 필수 수행 권장합니다.

왼쪽의 Summary를 보면 'Per base sequence quality'가 ❌ 표시되어 있으므로 trimming을 해주는게 좋겠네요.

(BIO) min@minyoungxi:~/main/Bio/DNA practice$ sudo apt install trim-galore

Trimming은 FASTQ 파일에서 품질이 낮거나 신뢰할 수 없는 부분을 잘라내는 작업을 말합니다. 즉, NGS로 얻은 FASTQ 데이터의 불필요하거나 품질이 낮은 염기서열을 제거하는 과정입니다.

Trimming 과정에서는 Quality-based trimming(품질 기반 트리밍), Adapter removal(어댑터 제거)와 같은 일들을 수행합니다. 전자의 경우 염기서열 끝부분에서 품질이 낮은 영역을 제거하는 작업이고, 후자의 경우 DNA 조각을 시퀀싱 할 때 사용하는 어댑터가 데이터에 섞이는 경우 어댑터 서열을 제거합니다.

Trimming을 위해서 `trim_galore` 을 설치하고 아래와 같이 명령어를 입력하여 Trimming을 진행합니다.

Trimming이 수행되면 아래와 같은 결과가 나옵니다.

1. Total reads processed
- 처리된 총 리드 수 (reads 수) : 7,178,576개의 리드를 처리했습니다.
2. Reads with adapters
- 어댑터 서열이 포함된 리드 수 : 3,245,087개 (전체의 45.2%)의 리드에 어댑터 서열이 있어서 이를 제거했습니다.
- 이 수치가 높으면 실험 과정에서 사용한 어댑터가 데이터에 많이 포함되었다는 뜻이며, 트리밍이 매우 잘 된 것입니다.
3 . Reads written (passing filters)
- 트리밍 후 최종적으로 남겨진 리드 수 : 7,178,576개 (100%) 모두 남겨졌습니다.
- 모든 리드가 최소 품질 기준을 통과하여 삭제된 리드는 없었습니다.

Alignment 정렬

정렬(Alignment) 이란, FASTQ 파일에서 읽은 짧은 서열(Read)을 기존에 알려진 표준 유전체(Reference Genome)에 정확한 위치로 매핑(mapping)하는 과정입니다.

정렬을 통해 다음의 분석을 수행할 수 있습니다.
- 유전자 발현량 분석 (어떤 유전자에서 얼마나 발현되었는지)
- 변이(variants) 발견 및 분석
- 특정 유전자의 발현 패턴 분석 및 질병 예측 분석 가능

1. STAR 설치

conda install -c bioconda star

2. Reference Genome 준비
먼저 표준 유전체(reference genome)를 준비하고, STAR가 사용할 수 있도록 색인(index)을 생성해야 합니다.

아래는 인간 유전체(GRCh38) 기준의 예시

Reference Genome 다운로드

wget http://ftp.ensembl.org/pub/release-110/fasta/homo_sapiens/dna/Homo_sapiens.GRCh38.dna.primary_assembly.fa.gz
gunzip Homo_sapiens.GRCh38.dna.primary_assembly.fa.gz

GTF 파일 (유전자 위치 정보) 다운로드

wget http://ftp.ensembl.org/pub/release-110/gtf/homo_sapiens/Homo_sapiens.GRCh38.110.gtf.gz
gunzip Homo_sapiens.GRCh38.110.gtf.gz

Reference Genome

Reference Genome(기준 유전체)이란 한 생물종의 전체 유전체 정보를 대표적으로 정리한 표준적인 유전체 서열을 의미합니다.

쉽게 말해, 한 종(예: 인간, 쥐, 효모 등)을 대표할 수 있는 표준적인 DNA 서열을 하나의 완전한 세트로 정리한 것입니다.
- 예를 들어, 인간의 대표적인 레퍼런스 유전체는 GRCh38 (Genome Reference Consortium Human Build 38) 입니다.

Reference Genome은 유전체 분석의 기준 역할을 합니다.

특히 NGS 분석에서는 짧게 읽힌 서열(read)을 어디에 매핑(위치)할지 결정하는 기준이 됩니다.

즉, 짧은 서열(Reads, FASTQ 데이터) → Reference genome에 정렬 → 특정 유전자 또는 변이 분석 의 흐름으로 실행됩니다.

FASTQ 데이터 (짧은 서열 Reads)
    ↓ Alignment (STAR 등)
Reference Genome (FASTA) → Reads의 위치를 결정하는 기준
    ↓ 발현량 정량화 (featureCounts, Salmon 등)
Reference Genome annotation (GTF) → 유전자별로 reads 개수를 계산하는 기준

3. STAR Index 만들기
유전체에 대한 색인(index)을 만들어줍니다. (최초 1회만 수행)

mkdir star_index

STAR --runThreadN 8 \
     --runMode genomeGenerate \
     --genomeDir star_index \
     --genomeFastaFiles Homo_sapiens.GRCh38.dna.primary_assembly.fa \
     --sjdbGTFfile Homo_sapiens.GRCh38.110.gtf \
     --sjdbOverhang 100

--runThreadN: 사용할 CPU 코어 수

--sjdbOverhang: 리드 길이 -1 (대체로 100 정도 권장)

4. 정렬(Alignment) 수행

이제 트리밍한 파일을 정렬합니다.
예시 파일명: `SRR390728_1_trimmed.fq`

mkdir aligned

STAR --runThreadN 8 \
     --genomeDir star_index \
     --readFilesIn SRR390728_1_trimmed.fq \
     --outFileNamePrefix aligned/SRR390728_ \
     --outSAMtype BAM SortedByCoordinate

`--runThreadN 8`: CPU를 8개 사용할 때 (필요에 따라 숫자 조절 가능)
`--genomeDir star_index`: 방금 생성한 index가 저장된 폴더
`--readFilesIn`: 입력 FASTQ (pair로 존재한다면 둘 다 넣어줌)
`--outFileNamePrefix`: 결과 파일 저장 위치와 이름 지정
`--outSAMtype BAM SortedByCoordinate`: BAM 형태로 정렬되어 저장

5. Alignment 결과 확인하기

aligned/
    ├── SRR390728_Aligned.sortedByCoord.out.bam  # 최종 결과 파일
    ├── SRR390728_Log.final.out                  # 정렬 통계 정보
    ├── SRR390728_Log.out                        # STAR 수행 로그

- bam: Binary Alignment Map의 약자로, 정렬된 결과 파일이며 이후 발현량 정량화의 입력 데이터로 사용됩니다.
- Log.final.out:
    STAR의 정렬 정확도와 성능을 요약한 통계정보.
    이 파일을 열어서 alignment의 품질을 확인할 수 있습니다.

출력된 결과 해석

samtools flagstat SRR390728_Aligned.sortedByCoord.out.bam

495008 + 0 in total (QC-passed reads + QC-failed reads)
412848 + 0 primary
82160 + 0 secondary
0 + 0 supplementary
0 + 0 duplicates
0 + 0 primary duplicates
82160 + 0 secondary

지금 결과에서 주목할 수치는 다음과 같습니다:

- 전체 read 수: 495,008개

- 정확히 정렬된 reads (properly paired): 94.41%
    → 일반적으로 이 수치가 80% 이상이면 양호한 데이터입니다.
    → 94.41%는 매우 우수한 수치입니다.

- Singleton 비율이 약 5.59%로 적당한 수준입니다.
    - 일반적으로 10% 미만이면 문제가 없습니다.
결론적으로, 현재 데이터는 정렬 품질이 아주 좋은 상태입니다.

---

다음 포스팅에서는 전처리 완료된 'SRR390728_Aligned.sortedByCoord.out.bam' 데이터를 파이썬으로 불러와서 분석하는 방법을 알아보겠습니다.

저작자표시 비영리 변경금지

내 블로그 - 관리자 홈 전환	`Q` `Q`
새 글 쓰기	`W` `W`

글 수정 (권한 있는 경우)	`E` `E`
댓글 영역으로 이동	`C` `C`

이 페이지의 URL 복사	`S` `S`
맨 위로 이동	`T` `T`
티스토리 홈 이동	`H` `H`
단축키 안내	`Shift` + `/` `⇧` + `/`

[BioInformatics] NGS Data Analysis Workflow - 시퀀싱된 데이터를 파이썬으로 다뤄보자

NGS(Next Generation Sequencing)가 뭐지 ?

NGS Workflow

시퀀싱 완료된 데이터 형태

FASTQ에 대해서

FastQC & Trimming

Alignment 정렬

Reference Genome

티스토리툴바

단축키

내 블로그

블로그 게시글

모든 영역