テロメアは、真核生物の染色体末端に見られる繰り返し配列で、多くの細胞で分裂時にその平均長が短くなることから、「体内時計」と考えられている。テロメアの長さの異常は、老化との関連に加え、複数の癌やテロメア短小症候群との関連や、様々な疾患の危険因子としてよく知られている。テロメアの長さを測定する方法の多くは、全染色体の平均的な長さを報告するが、特定の染色体群の異常が特定の疾患のバイオマーカーとなることが知られている。テロメアは繰り返しが多いため、ショートリードのシーケンシング手法ではこの分解能での特性解析は困難であり、ロングリードのシーケンシング手法でもなお難しい。
本著者らは、ロングリードシーケンスデータから染色体特異的なテロメア長を報告する手法であるTelogatorを発表する。Telogatorが、様々なリード長やエラー率のシミュレーションデータから、染色体特異的なテロメア長を検出することができることを実証する。次に、Telogatorを10個の生殖細胞サンプルに適用し、平均テロメア長を報告する際にショートリード法との高い相関性を示した。さらに、一般的なサブテロメアリアレンジメントを調べ、これらのハプロタイプを持つサンプルでテロメア/サブテロメアの境界をanchorするために必要な最小リード長を特定した。
インストール
ubuntu18.04でコマンドの使用方法のみ確認した。
依存
他に必要なツール
- pbmm2
- samtools
git clone https://github.com/zstephens/telogator.git
$ python3 telogator/grab_reads_from_subtel_aln.py -h
usage: grab_reads_from_subtel_aln.py [-h] --bam input.bam --fa input.fa --out output.fa --bed subtel.bed [--samtools samtools] [--readtype CCS / CLR / SRA] [--newlines yes / no]
grab_subreads_from_t2t-and-subtel_aln.py
optional arguments:
-h, --help show this help message and exit
--bam input.bam * Input BAM (default: None)
--fa input.fa * Raw subreads (.fa or .fa.gz) (default: None)
--out output.fa * Output reads (.fa or .fa.gz) (default: None)
--bed subtel.bed * Subtel regions (default: None)
--samtools samtools /path/to/samtools (default: samtools)
--readtype CCS / CLR / SRA
Read name format (default: CLR)
--newlines yes / no Read data split across multiple lines? (default: yes)
> python3 telogator/telogator.py -h
usage: telogator.py [-h] -i in.sam / in.p / - -o output/ [-k default_kmers.tsv] [-l 5000] [-t 0.9] [-p 0.5] [-pq 0.25] [--sa largest] [--sm mapq] [--job my_job n_jobs] [--plot]
[--debug]
Telogator v1.0
optional arguments:
-h, --help show this help message and exit
-i in.sam / in.p / - * Long reads aligned to subtel ref (- for stdin) (default: None)
-o output/ * Path to output directory (default: None)
-k default_kmers.tsv Telomere kmers to search for (default: )
-l 5000 Min read length (default: 5000)
-t 0.9 Max frac of read that can be tel (default: 0.9)
-p 0.5 Telomere signal threshold (0-1) (default: 0.5)
-pq 0.25 Max minor p/q fraction in tel region (default: 0.25)
--sa largest Subtel/tel anchoring strategy (default: largest)
--sm mapq Repeated matches trimming strategy (default: mapq)
--job my_job n_jobs Job splitting for parallelization (default: (0, 0))
--plot Create read plots (default: False)
--debug Print results for each read as its processed (default: False)
実行方法
1、Telomere-to-telomere consortiumのレポジトリ(link)からT2Tアセンブリをダウンロードする。ここでは最新のv1.1をダウンロードした。それにサブテロメア領域のアセンブリを追加してインデックスを付ける。
git clone https://github.com/zstephens/telogator.git
gzip -dc telogator/resources/stong_subtels.fa.gz > stong_subtels.fa
cat chm13.draft_v1.1.fasta stong_subtels.fa > t2t-and-subtel.fa
samtools faidx t2t-and-subtel.fa
2、PacBio CLRのリードを作成した全ゲノムリファレンスにアラインする。レポジトリではPacBio のpbmm2アライナー(BLASRの後継)を使っている。
pbmm2 align t2t-and-subtel.fa clr-reads.fa.gz aln.bam --preset SUBREAD --sort
3、サブテロメア領域のリードを抽出する。
python3 telogator/grab_reads_from_subtel_aln.py \
--bam aln.bam \
--fa clr-reads.fa.gz \
--out subtel-reads.fa.gz \
--bed resources/subtel_regions.bed
テストした時は、このステップで回収されたリードがゼロになってしまった。
4、抽出したリードを、サブテロメア領域の専用リファレンス(レポジトリに含まれる)にアラインする。
gzip -dc telogator/resources/t2t-telogator-ref.fa.gz > t2t-telogator-ref.fa
pbmm2 align t2t-telogator-ref.fa subtel-reads.fa.gz subtel_aln.bam --preset SUBREAD --sort
5、専用リファレンスにアラインされたbamを使ってTelogatorを実行する。
samtools view subtel_aln.bam | python3 telogator/telogator.py -i - -o telogator_out/
6、プロットを作成し、レポートを出力する。
python3 merge_jobs.py -i telogator_out/
T2Tアセンブリ(CHM13hTERTヒト細胞株(リンク))のpacbioデータが公開されているので、これを使うことができます(リンク)。ONTのロングリードを使う場合の注意点はレポジトリに書かれています。興味がある方は確認して下さい。
引用
Telogator: a method for reporting chromosome-specific telomere lengths from long reads
Zachary Stephens, Alejandro Ferrer, Lisa Boardman, Ravishankar K Iyer, Jean-Pierre A Kocher
Bioinformatics, Published: 12 January 2022
参考
関連