NGS解析(RNAseq)
坊農 秀雅(ぼうのう ひでまさ)
大学共同利用機関法人 情報・システム研究機構
ライフサイエンス統合データベースセンター
(DBCLS)
今日のレシピ
»NGS は Now Generation
Sequencer!
»NGSデータ解析の概要!
»具体的な解析例→ハンズオン!
»おまけ
4
NGS(Next Generation
Sequencer)は既に Now
Generation Sequencer
DNA塩基配列解読の超高速化
• かつてはSanger法
• 最近は「次世代シーケンサー(NGS)」
‒Illumina: Sequence By Synthesis
• http://www.youtube.com/watch?v=womKfikWlxM !
‒Life Technologies(現 ThermoFisher Scientific)
•ヌクレオチドがDNA鎖に取り込まれる過程でポリメラー
ゼによって放出される水素イオンを検出
• http://www.youtube.com/watch?v=MxkYa9XCvBQ !
‒PacBio: 一分子・リアルタイム(SMRT®)検出
• http://www.youtube.com/watch?v=NHCJ8PtYCFc
6
MiSeq → NextSeq
• Illumina社のデスクトップ次世代シーケンサ
NextSeq 500
• ヒト全ゲノムリシーケンス(3Gb ゲノム)
‒coverage: 30x
‒リード長: 2 x 150塩基
‒時間: 29h
‒データ出力: 100-120Gb/run
‒(右図はMiSeq)
7
NGSからの生データ
• FASTQフォーマットのファイル
‒4行/readが基本単位
‒Illumina NextSeq 500
‒3億リードx4行
‒=12億行
SRR001356.1 2023DAAXX:5:1:123:563 length=33
TGTCGGTCCAGCTCGGCCTTGGGCTCCGTTTTC
+SRR001356.1 2023DAAXX:5:1:123:563 length=33
-IIIIIIII8IIIIIIIIIII6IIIIIIIII9I
@SRR001356.2 2023DAAXX:5:1:123:476 length=33
TCTGAACCCGACTCCCTTTCGATCGGCCGCGGG
+SRR001356.2 2023DAAXX:5:1:123:476 length=33
IIIIIIIIIIIIIIIIIIIIIGIIIIIII-III
@SRR001356.3 2023DAAXX:5:1:121:746 length=33
GTGGCAGCGTTTTTGGGCCCGCCGCTTGCCGTT
+SRR001356.3 2023DAAXX:5:1:121:746 length=33
IIIII&IIIIIIIIIIIIIIIIHI1IIIIIIII
• ファイルサイズも4Gbyte/file超
‒FAT32フォーマットでは扱えない
• いわゆる「開く」ことが不可能
8
NGSでできること
• RNA転写量測定!
–RNAseq(transcriptome sequencing)!
• DNA結合タンパク質の結合配列の解析!
–ChIPseq(ヒストンや転写因子)!
• ChIPはChromatin immunoprecipitationの略!
• 多型解析!
–Exome(exon限定), Re-sequence!
• その他、塩基配列解読が伴うさまざまな応用
9
NGSによるデータ解析の概要
トランスクリプトーム解析
(RNAseq)
RNAseqとは？
• 「次世代シーケンサを利用して、サンプル中
の RNA の中身に関する情報を得るために cDNA をシーケンシングする方法」!
–http://en.wikipedia.org/wiki/RNA-Seqより勝手に翻訳!
• Whole transcriptome shutgun sequencing (WTSS) や!
• Transcriptome sequencingとも
11
SRR001356.1 2023DAAXX:5:1:123:563 length=33!
TGTCGGTCCAGCTCGGCCTTGGGCTCCGTTTTC!
+SRR001356.1 2023DAAXX:5:1:123:563 length=33!
-IIIIIIII8IIIIIIIIIII6IIIIIIIII9I!
@SRR001356.2 2023DAAXX:5:1:123:476 length=33!
TCTGAACCCGACTCCCTTTCGATCGGCCGCGGG!
+SRR001356.2 2023DAAXX:5:1:123:476 length=33!
IIIIIIIIIIIIIIIIIIIIIGIIIIIII-III!
@SRR001356.3 2023DAAXX:5:1:121:746 length=33!
GTGGCAGCGTTTTTGGGCCCGCCGCTTGCCGTT!
+SRR001356.3 2023DAAXX:5:1:121:746 length=33!
IIIII&IIIIIIIIIIIIIIIIHI1IIIIIIII
FASTQ
RNAseq データ
解析の流れ1
ゲノム
1.tophat
(bowtie)
2.cufflinks
.fa
ゲノムに対する多重配列アラインメント
.gtf
ゲノムアノ
テーション
予測転写単位ごとの
（推定）発現量情報
.bam
遺伝子アノ
テーション
3.cummeRbund
12
種々のデータフォーマット
ファイルフォーマット ファイル拡張子
1
2
3
4
5
6
7
FASTA
FASTQ
SRA/SRA-lite
SAM/BAM
GTF(GFF)
BED
VCF
.fa .fasta
.fq .fastq
.sra .lite.sra
.sam .bam
.gtf .gff
.bed
.vcf
13
1. FASTA
• FASTAというプログラムで使われる配列データ形式!
–プレーンテキスト。ファイル拡張子: .fa .fasta など!
• 1行目に“>”で始まる1行のヘッダ行!
• 2行目以降に実際のシーケンス文字列
>gi|5524211|gb|AAD44166.1| cytochrome b [Elephas maximus maximus]!
LCLYTHIGRNIYYGSYLYSETWNTGIMLLLITMATAFMGYVLPWGQMSFWGATVITNLFSAIPYIGTNLV!
EWIWGGFSVDKATLNRFFAFHFILPFTMVALAGVHLTFLHETGSNNPLGLTSDSDKIPFHPYYTIKDFLG!
LLILILLLLLLALLSPDMLGDPDNHMPADPLNTPLHIKPEWYFLFAYAILRSVPNKLGGVLALFLSIVIL!
GLMPFLHTSKHRSMMLRPLSQALFWTLTMDLLTLTWIGSQPVEYPYTIIGQMASILYFSIILAFLPIAGX!
IENY
参考: http://ja.wikipedia.org/wiki/FASTA
14
2. FASTQ
• NGSデータの配列データ形式のデファクトスタンダード!
–プレーンテキスト。ファイル拡張子: .fq .fastq など!
• 1行目に“@”で始まる1行のヘッダ行!
• 2行目に実際の塩基配列!
• 3行目に”+”!
• 4行目に2行目に記述した配列のクオリティ値
@SEQ_ID!
GATTTGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTTT!
+!
!''*((((***+))%%%++)(%%%%).1***-+*''))**55CCF>>>>>>CCCCCCC65
参考: http://ja.wikipedia.org/wiki/Fastq
15
3. SRA, SRA-lite
• FASTQ形式の代わりに使われている、NGS
配列データ配布フォーマット!
–配列拡張子: .sra .lite.sra !
• SRA-toolkitを使ってFASTQを生成できる!
–http://www.ncbi.nlm.nih.gov/Traces/sra/?view=software
fastq-dump -A SRR233129 SRR233129.lite.sra
16
4. SAM/BAM
• ゲノムマッピングしたときに生成されるアラ
インメントのフォーマット!
–リファレンスゲノム配列に対するアラインメント!
• SAMはプレーンテキスト(ASCII)形式なのに
対して、BAMはバイナリ(binary)形式!
1:497:R:-272+13M17D24M!
19:20389:F:275+18M2D19M!
19:20389:F:275+18M2D19M!
9:21597+10M2I25M:R:-209!
113!
99!
147!
83!
1!
1!
1!
1!
497! 37!
17644!0!
17919!0!
21678!0!
37M! 15!
37M! =!
18M2D19M!
8M2I27M!
100338662! 0!
CGGGTCTGACCTGAGGAGAACTGTGCTCCGCCTTCAG! 0;==-==9;>>>>>=>>>>>>>>>>>=
17919!314! TATGACTGCTAATAATACCTACACATGTTAGAACCAT!
>>>>>>>>>>>>>>>>>>>><<>>><<>>4:
=!
17644!-314! GTAGTACCAACTGTAAGTCCTTATCTTCATACTTTGT!
;44999;499<8<8<<<8<<><<<<><
=!
21469!-244! CACCACATCACATATACCAAGCCTGGCTGTGTCTTCT!
<;9<<5><<<<><<<>><<><>><9>
参考: http://genome.sph.umich.edu/wiki/SAM
17
5. GTF(GFF)
• General Transfer Format. GFF(General
Feature Format)のversion2!
• ゲノムアノテーションのフォーマット!
–例: ゲノム上のどこに遺伝子があるか
X!
X!
X!
X!
X!
X!
Ensembl!
Ensembl!
Ensembl!
Ensembl!
Ensembl!
Ensembl!
Repeat!2419108! 2419128! 42! .! .! hid=trf; hstart=1; hend=21!
Repeat!2419108! 2419410! 2502! -! .! hid=AluSx; hstart=1; hend=303!
Repeat!2419108! 2419128! 0! .! .! hid=dust; hstart=2419108; hend=2419128!
Pred.trans.!2416676! 2418760! 450.19!-! 2! genscan=GENSCAN00000019335!
Variation! 2413425! 2413425! .! +! .! !
Variation! 2413805! 2413805! .! +! .
参考: http://asia.ensembl.org/info/website/upload/gff.html
18
6. BED
• ゲノムアノテーションのフォーマット!
–例: ゲノム上のどこに遺伝子があるか
track name=pairedReads description="Clone Paired Reads" useScore=1!
chr22 1000 5000 cloneA 960 + 1000 5000 0 2 567,488, 0,3512!
chr22 2000 6000 cloneB 900 - 2000 6000 0 2 433,399, 0,3601
参考: http://genome.ucsc.edu/FAQ/FAQformat.html#format1
19
7. VCF
• Variant Call Format!
• 配列の多型を記述するフォーマット
##fileformat=VCFv4.0!
##fileDate=20110705!
##reference=1000GenomesPilot-NCBI37!
##phasing=partial!
##INFO=<ID=NS,Number=1,Type=Integer,Description="Number of Samples With Data">!
##INFO=<ID=DP,Number=1,Type=Integer,Description="Total Depth">!
##INFO=<ID=AF,Number=.,Type=Float,Description="Allele Frequency">!
##INFO=<ID=AA,Number=1,Type=String,Description="Ancestral Allele">!
##INFO=<ID=DB,Number=0,Type=Flag,Description="dbSNP membership, build 129">!
##INFO=<ID=H2,Number=0,Type=Flag,Description="HapMap2 membership">!
##FILTER=<ID=q10,Description="Quality below 10">!
##FILTER=<ID=s50,Description="Less than 50% of samples have data">!
##FORMAT=<ID=GQ,Number=1,Type=Integer,Description="Genotype Quality">!
##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">!
##FORMAT=<ID=DP,Number=1,Type=Integer,Description="Read Depth">!
##FORMAT=<ID=HQ,Number=2,Type=Integer,Description="Haplotype Quality">!
#CHROM POS ID
REF ALT QUAL FILTER INFO
FORMAT
Sample1
Sample2
Sample3!
2
4370 rs6057 G A
29 .
NS=2;DP=13;AF=0.5;DB;H2
GT:GQ:DP:HQ 0|0:48:1:52,51 1|0:48:8:51,51 1/1:43:5:.,.!
2
7330 .
T A
3 q10 NS=5;DP=12;AF=0.017
GT:GQ:DP:HQ 0|0:46:3:58,50 0|1:3:5:65,3 0/0:41:3!
2
110696 rs6055 A G,T 67 PASS NS=2;DP=10;AF=0.333,0.667;AA=T;DB GT:GQ:DP:HQ 1|2:21:6:23,27 2|1:2:0:18,2 2/2:35:4!
2
130237 .
T .
47 .
NS=2;DP=16;AA=T
GT:GQ:DP:HQ 0|0:54:7:56,60 0|0:48:4:56,51 0/0:61:2!
2
134567 microsat1 GTCT G,GTACT 50 PASS NS=2;DP=9;AA=G
GT:GQ:DP 0/1:35:4
0/2:17:2
1/1:40:3
参考: http://en.wikipedia.org/wiki/Variant_Call_Format
20
NGSに比べてマイクロアレイ
• 2000年前後に使われ初めて、ある程度（技術と
して）枯れてきた!
–参考: 公共データベースの登録数の推移!
–遺伝子発現バンク(GEO)目次 http://lifesciencedb.jp/geo/!
• 本もそれなりに出ている
21
マイクロアレイ解析の流れ
遺伝子アノ
テーション
oligoprobeに対応する
遺伝子ごとの発現量
Genespring
22
SRR001356.1 2023DAAXX:5:1:123:563 length=33!
TGTCGGTCCAGCTCGGCCTTGGGCTCCGTTTTC!
+SRR001356.1 2023DAAXX:5:1:123:563 length=33!
-IIIIIIII8IIIIIIIIIII6IIIIIIIII9I!
@SRR001356.2 2023DAAXX:5:1:123:476 length=33!
TCTGAACCCGACTCCCTTTCGATCGGCCGCGGG!
+SRR001356.2 2023DAAXX:5:1:123:476 length=33!
IIIIIIIIIIIIIIIIIIIIIGIIIIIII-III!
@SRR001356.3 2023DAAXX:5:1:121:746 length=33!
GTGGCAGCGTTTTTGGGCCCGCCGCTTGCCGTT!
+SRR001356.3 2023DAAXX:5:1:121:746 length=33!
IIIII&IIIIIIIIIIIIIIIIHI1IIIIIIII
FASTQ
RNAseq データ
解析の流れ
ゲノム
1.tophat
(bowtie)
2.cufflinks
.fa
ゲノムに対する多重配列アラインメント
.gtf
ゲノムアノ
テーション
予測転写単位ごとの
（推定）発現量情報
.bam
遺伝子アノ
テーション
3.cummeRbund
23
マイクロアレイのデータ形式の実際
• タブ区切りテキスト!
–数万(=スポットの数)行!
• (古い)Excelでも「開ける」!
–Excel2003の行数制限内!
• コマンドライン操作なしで中身が直接見れる
24
データ解析に必要なもの
解析ソフト
マイクロアレイ NGS(RNAseq)
+++
+++
遺伝子
+++
アノテーション
+++
ゲノム
アノテーション
++
ゲノム配列
-
++
コマンドライ
ン操作
+
+++
計算機パワー
+
+++
25
マイクロアレイとの違い: RPKM
• Reads Per Kilobase per Million mapped reads!
• ノーマライズした遺伝子発現量!
–100万リード数マップされたとき、転写産物を
1000塩基長としたときのマップされたリード数!
• FPKMもほぼ同じ!
–Fragments Per Kilobase of exon per Million
mapped fragments!
–!
• Reference: Nat Methods, 5(7):621-628.
26
具体的な解析例
ハンズオン!
次世代シークエンス解析スタンダード!
NGSのポテンシャルを活かしきるWET&DRY
• 通称、「愛本」!
• https://www.yodosha.co.jp/jikkenigaku/book/9784758101912/
28
The cat way
• 理化学研究所の二階堂愛さんのブログ!
–http://cat.hackingisbelieving.org/lecture/!
tophat
! -p 8 -r 100 -o output_dir/iPS_01 bowtie2_indexes/mm9 iPS_01_1.fastq
!
cuffdiff
-p 24 ensembl_gene.gtf !
-L iPS_01,iPS_02,hESC_01,hESC_02,Fibroblast_01,Fibroblast_02!
-o! results iPS_01.bam,iPS_2.bam hESC_1.bam,hESC_2.bam
Fibroblast_01.bam,Fibroblast_02.bam!
!
• オープンソースソフトウェア!
–Tuxedo suite!
• bowtie,tophat,cufflinks!
–R + Bioconductor
29
SRR001356.1 2023DAAXX:5:1:123:563 length=33!
TGTCGGTCCAGCTCGGCCTTGGGCTCCGTTTTC!
+SRR001356.1 2023DAAXX:5:1:123:563 length=33!
-IIIIIIII8IIIIIIIIIII6IIIIIIIII9I!
@SRR001356.2 2023DAAXX:5:1:123:476 length=33!
TCTGAACCCGACTCCCTTTCGATCGGCCGCGGG!
+SRR001356.2 2023DAAXX:5:1:123:476 length=33!
IIIIIIIIIIIIIIIIIIIIIGIIIIIII-III!
@SRR001356.3 2023DAAXX:5:1:121:746 length=33!
GTGGCAGCGTTTTTGGGCCCGCCGCTTGCCGTT!
+SRR001356.3 2023DAAXX:5:1:121:746 length=33!
IIIII&IIIIIIIIIIIIIIIIHI1IIIIIIII
FASTQ
RNAseq データ
解析の流れ1
ゲノム
1.tophat
(bowtie)
2.cufflinks
.fa
ゲノムに対する多重配列アラインメント
.gtf
ゲノムアノ
テーション
予測転写単位ごとの
（推定）発現量情報
.bam
遺伝子アノ
テーション
3.cummeRbund
30
0. ゲノム＆アノテーション
ファイル取得
• Illumina iGenomes
31
Illumina iGenomes
• 超でかいですが、とりあえず必要なのは(human
UCSC hg19の場合)!
• (インデックス済みゲノム配列)!
–Homo_sapiens/UCSC/hg19/Sequence/Bowtie2Index/ 以下!
–これをコピーしておきましょう!
% cp
-r Bowtie2Index/ 141025/
!
• (ゲノムアノテーション)!
–Homo_sapiens/UCSC/hg19/Annotation/Genes/genes.gtf !
–これもコピーしておきましょう
% cp genes.gtf 141025/
32
1. fastqファイル取得
• すでに.fastqファイルを持っている場合はpass!
• .sra 形式のファイルは fastq-dump で変換!
–fastq-dumpはhomebrewのsratoolkitに含まれています!
% ! brew install sratoolkit -v
–インストール失敗した人はtapしていないからかも!
% brew tap homebrew/science -v
!
% ! fastq-dump *.sra
–分オーダーで時間がかかります!
• .fastqなファイルをゲット
33
2. tophat(bowtie)!
でゲノムにマッピング
• tophatをインストール!
% brew install tophat -v
!
• bowtie2, samtoolsとか同時に色々入ります!
• 実行!!
%! tophat -p 8 -o tophat/CNTL Bowtie2Index/
genome GSE50491/SRR960409.fastq
!
• tophatからbowtieを呼び出してゲノム配列にマッピ
ングを行います
35
tophat2
• 出力ディレクトリがないと怒れました。作りましょう!
% ! mkdir tophat
• 再度実行!
% ! tophat -p 8 -o tophat/BRD4 Bowtie2Index/
genome
GSE50491/SRR960410.fastq
!
• オプション!
• p: 使用CPU数!
• r: ペアエンド間の距離の指定(ペアエンドシークエンス
の場合に指定)!
• 1fastqあたり小1時間かかりました(MacProで-p 12)
36
3. cufflinksで発現量を定量
• cufflinksをインストール!
% brew
install cufflinks -v
!
• 実行!!
% cufflinks -p 8 -o cufflinks_CNTL tophat/
!
CNTL/accepted_hits.bam!
%! cufflinks -p 8 -o cufflinks_BRD4 tophat/
BRD4/accepted_hits.bam
!
• 1bamあたり1.5時間程かかりました(MacProで-p 12)
38
cuffmerge
• cuffmergeでassemblyをmerge!
% cuffmerge -p8 —s hg19.fa assemblies.txt
!
• cuffcompare!
% ! cuffcompare -p8 —s hg19.fa -r genes.gtf
merged_asm/merged.gtf
!
• こちらを使うのは新規transcriptやslice variant
の発見を行う場合!
• そうでない場合は次のcuffdiffを
39
cuffdiff
% cuffdiff -p 8 -L CNTL,BRD4 -o cuffdiff
genes.gtf tophat/CNTL/accepted_hits.bam
tophat/BRD4/accepted_hits.bam
• homebrewで入れたバイナリではエラー…!
• 仕方なく、バイナリ配布のものを使用!
• http://cufflinks.cbcb.umd.edu/
% ~/Downloads/cufflinks-2.2.1.OSX_x86_64/
cuffdiff -p 8 -L CNTL,BRD4 -o cuffdiff
genes.gtf tophat/CNTL/accepted_hits.bam
tophat/BRD4/accepted_hits.bam
40
4. cummeRbundでRの世界に
% R
• cummeRbundをインストール!
> !source("http://bioconductor.org/biocLite.R")!
> biocLite("cummeRbund")
!
• 可能なグラフ in R Graphical Manual!
• http://rgm3.lab.nig.ac.jp/RGM/R_image_list?
package=cummeRbund&init=true!
> browseVignettes("cummeRbund")
42
cummeRbund(続)
• manualを参照!
>
>
>
>
library(“cummeRbund")!
!setwd("~/Desktop/141025")!
!cuff.dir <- "cuffdiff"!
cuff <- readCufflinks(dir=cuff.dir)
!
• 発現の密度分布プロット!
> dens <- csDensity(genes(cuff))!
> !dens
• CSV形式でFPKM値を出力
> gene.matrix <- fpkmMatrix(genes(cuff))!
> write.csv(gene.matrix, file="fpkm.csv")
43
4. おまけ
統計解析環境R
• Rを使ったトランスクリプトーム解析!
–(Rで)マイクロアレイデータ解析!
• http://www.iu.a.u-tokyo.ac.jp/~kadota/r.html!
–(Rで)塩基配列解析!
• http://www.iu.a.u-tokyo.ac.jp/~kadota/r_seq.html!
!
• トランスクリプトーム解析 by 門田幸二 from 共立出版!
–http://www.kyoritsu-pub.co.jp/bookdetail/9784320123700
45
有償のソフトウェアの利用
•CLC Genomics workbench!
•Agilent!
–Avadis NGS!
–GeneSpring!
•TIBCO Spotfire
46
SRR001356.1 2023DAAXX:5:1:123:563 length=33!
TGTCGGTCCAGCTCGGCCTTGGGCTCCGTTTTC!
+SRR001356.1 2023DAAXX:5:1:123:563 length=33!
-IIIIIIII8IIIIIIIIIII6IIIIIIIII9I!
@SRR001356.2 2023DAAXX:5:1:123:476 length=33!
TCTGAACCCGACTCCCTTTCGATCGGCCGCGGG!
+SRR001356.2 2023DAAXX:5:1:123:476 length=33!
IIIIIIIIIIIIIIIIIIIIIGIIIIIII-III!
@SRR001356.3 2023DAAXX:5:1:121:746 length=33!
GTGGCAGCGTTTTTGGGCCCGCCGCTTGCCGTT!
+SRR001356.3 2023DAAXX:5:1:121:746 length=33!
IIIII&IIIIIIIIIIIIIIIIHI1IIIIIIII
FASTQ
RNAseq データ
解析の流れ1
ゲノム
1.tophat
(bowtie)
2.cufflinks
.fa
ゲノムに対する多重配列アラインメント
.gtf
ゲノムアノ
テーション
予測転写単位ごとの
（推定）発現量情報
.bam
遺伝子アノ
テーション
3.cummeRbund
47
統合TVに動画チュートリアルが
• CLC Genomics Workbench でショートリー
ドのマッピングを行う!
–http://togotv.dbcls.jp/20110628.html
48
RNA-seq by Avadis NGS
• http://togotv.dbcls.jp/20111124.html
49
ChIP-seq by Avadis NGS
• http://togotv.dbcls.jp/20120626.html
50
GeneSpring
• https://www.youtube.com/user/GeneSpringTV
51
Spotfireによるcuffdiff出力の可視化
% cuffdiff -p 8 Caenorhabditis_elegans.WBcel215.69.gtf
-L N2,UV -o cuffdiff SRR454084.bam SRR454085.bam
52
© Copyright 2024 ExpyDoc
