Nanopore Long-read RNA-seq 분석 따라하기

🧬 RNA-seq 실전 연습｜Nanopore Long-read RNA-seq 분석 따라하기 (2025 최신 논문 기반)

“논문 속 RNA-seq 코드를 직접 구현해보고 싶다면?”
오늘은 2025 Nature Methods의 롱리드 RNA-seq 벤치마크 논문을 바탕으로
실제 분석 파이프라인을 Python 코드로 재현해봅니다.

(실제 데이터 대신 예시 FASTQ로 구현 가능하며, 초보자도 실행할 수 있는 버전입니다.)

---

📦 1️⃣ 환경 준비

RNA-seq 분석의 첫걸음은 환경 세팅부터!
논문에서도 Python + Bioconda 조합을 권장합니다.

# Conda 환경 생성
conda create -n rna_longread python=3.10
conda activate rna_longread

# 필수 패키지 설치
conda install -c bioconda minimap2 samtools gffread
pip install pandas numpy matplotlib

 

---

🧬 2️⃣ FASTQ → 정렬(Alignment)

논문에서는 minimap2를 사용하여 Nanopore 리드를
reference transcriptome에 정렬했습니다.

# reference genome과 annotation 파일
REF=reference_genome.fa
READS=nanopore_reads.fastq

# minimap2로 정렬 수행
minimap2 -ax splice -uf -k14 $REF $READS | samtools sort -o aligned.bam
samtools index aligned.bam

 

🧩 -uf -k14 옵션은 Nanopore direct RNA-seq용으로 권장된 파라미터입니다.
(Chen et al., Nature Methods, 2025 기준)

---

🔍 3️⃣ 전사체(Transcript) 어셈블리

정렬 결과를 바탕으로, gffread로 전사체를 조립합니다.

 

# GTF 생성
gffread aligned.bam -T -o transcripts.gtf

📌 논문에서는 Isoform-level 분석을 위해
“StringTie2 / TALON / FLAIR” 중 하나를 비교 사용했습니다.
여기선 간단히 gffread로 시연합니다.

---

📊 4️⃣ Isoform 수준 정량화 (Python Practice)

정렬된 BAM 파일에서 isoform별 리드 수를 계산하는
간단한 Python 연습 코드입니다.

 

import pysam
import pandas as pd

bam_file = "aligned.bam"
gtf_file = "transcripts.gtf"

bam = pysam.AlignmentFile(bam_file, "rb")

isoform_counts = {}
for read in bam:
    if read.reference_name:
        isoform_counts[read.reference_name] = isoform_counts.get(read.reference_name, 0) + 1

df = pd.DataFrame(list(isoform_counts.items()), columns=["Transcript", "ReadCount"])
df = df.sort_values("ReadCount", ascending=False)

# 상위 10개 isoform 출력
print(df.head(10))
df.to_csv("isoform_counts.csv", index=False)

💡 실제 논문에서는 isoform-level TPM 계산까지 포함했지만,
위 코드는 연습용으로 리드 카운트를 세는 기본 구조입니다.

---

📈 5️⃣ 품질(QC) 시각화

논문에서는 리드 길이 분포·품질 점수를 함께 분석했습니다.
아래는 연습용 Python 시각화 코드입니다.

 

import matplotlib.pyplot as plt
import numpy as np
import gzip

fastq = "nanopore_reads.fastq"
lengths = []

with open(fastq, "r") as f:
    for i, line in enumerate(f):
        if i % 4 == 1:  # 시퀀스 라인
            lengths.append(len(line.strip()))

plt.hist(lengths, bins=50)
plt.title("Nanopore Read Length Distribution")
plt.xlabel("Read length (bp)")
plt.ylabel("Count")
plt.show()

📊 평균 리드 길이와 분포를 확인하면,
샘플의 품질과 스플라이싱 분석 가능성을 빠르게 파악할 수 있습니다.

---

🧠 6️⃣ 논문과의 비교 포인트

항목논문 수행 내용연습 코드 버전

정렬기	minimap2 (splice mode)	동일
정량화	TALON/FLAIR 기반 TPM	단순 read count
QC	Read length / Error rate / Isoform diversity	Length histogram
분석 목표	Isoform diversity benchmark	기본 파이프라인 학습용

---

💬 정리하며

🧩 “논문 속 RNA-seq 코드를 그대로 따라 하려면?”
핵심은 정렬 → 조립 → 정량화 → QC 시각화 네 단계입니다.
오늘의 코드만으로도
RNA-seq 데이터가 어떻게 구조화·해석되는지 체감할 수 있을 거예요.

---

✨ 마무리 Tip

• 실제 데이터는 ENA의 공개 Nanopore RNA-seq 샘플을 이용하세요.
  (검색 키워드: Nanopore direct RNA-seq human sample)
• 더 정확한 isoform 정량화에는 TALON, FLAIR, or StringTie2를 사용하세요.
• Python 분석과 함께 R의 DESeq2 / edgeR를 결합하면 통계적 유의성 검증까지 가능합니다.

---

📌 추천 관련 포스팅

• 🔗 2025 RNA-seq 기술 총정리｜롱리드·공간 전사체학 비교
• 🔗 AI 기반 RNA-seq 자동 분석 툴 3가지 리뷰