PowerPoint プレゼンテーション

Size: px

Start display at page:

Download "PowerPoint プレゼンテーション"

ひろみやまのかみしゃ
5 years ago
Views:

1 V1 次世代シークエンサ実習 II

2 本講義の内容 Reseq 解析 RNA-seq 解析公開データ取得クオリティコントロールマッピング変異検出公開データ取得クオリティコントロールマッピング発現定量 FPKM を算出します 2

3 R N A - s e q とはメッセンジャー RNA(mRNA) をキャプチャして次世代シーケンサーでシーケンシングする手法リファレンスがある生物種の場合 : 既知遺伝子にマッピングするリファレンスにマッピングして遺伝子構造を同定する本日の内容リファレンスがない生物種の場合 : アセンブリングして転写物構造を予測しそれに対してマッピングする近いゲノムのリファレンスにマッピングする Copyright Amelieff Corporation. All Rights Reserved. 3

4 R N A - s e q 解析 : パイプラインスプライスジャンクションを考慮したマッピングを行います 4

5 R N A - s e q 解析 : パイプライン今日は一部のコマンドを実行しますヒートマップ GO 解析パスウェイ解析など様々 5

6 R N A - s e q 解析でできること発現量の定量比較新規転写物新規スプライシングバリアントの探索 RNA-seq がマイクロアレイと比較して優れている点新規転写物や融合遺伝子が検出可 SNV small Indel も検出可プローブの設計を必要としない ( 非モデル生物にも対応可 ) Copyright Amelieff Corporation. All Rights Reserved. 6

7 R N A - s e q 解析 : データ酵母のゲノムのリファレンス取得リファレンス配列の fasta のみではなくマッピングソフトのインデックスファイルや遺伝子情報ファイルも一緒に圧縮されて公開しています 7

8 R N A - s e q 解析 : データ TRY! 酵母のゲノムのリファレンス配列を確認 $ cd /home/admin1409/genome/saccharomyces_cerevis iae/ncbi/build3.1 $ ll Sequence 8

9 R N A - s e q 解析 : データ TRY! 酵母のゲノムのリファレンス配列を確認 $ ll Sequence/Bowtie2Index 今回はこちらのインデックスを使用します 9

R N A - s e q 解析 : データシーケンスデータ取得 http://trace.ddbj.nig.

10 R N A - s e q 解析 : データシーケンスデータ取得 DDBJ の Sequence Read Archive Search 10

11 R N A - s e q 解析 : データシーケンスデータ取得 Accession に SRR と入力 Search 11

12 R N A - s e q 解析 : データシーケンスデータ取得今回は 1 サンプルで実行しますが発現比較する場合には複数サンプル必要で replicate も多いほうが良いですここからダウンロード実験の詳細 Navigation エリアの Experiment SRX をクリック 12

13 R N A - s e q 解析 : データシーケンスデータ取得他にもシーケンサのプラットフォームやリード長などの情報も記載されています転写産物 13

14 R N A - s e q 解析 : データシーケンスデータ取得 ( 実行済み ) ダウンロードします $ wget ftp://ftp.ddbj.nig.ac.jp/ddbj_database/dra/fastq/sra055/sr A055683/SRX157933/SRR518891_1.fastq.bz2 $ wget ftp://ftp.ddbj.nig.ac.jp/ddbj_database/dra/fastq/sra055/sr A055683/SRX157933/SRR518891_2.fastq.bz2 14

15 R N A - s e q 解析 : データシーケンスデータ取得 ( 実行済み ) 解凍して先頭 1000 リードを抽出します $ bunzip2 SRR518891_1.fastq.bz2 $ bunzip2 SRR518891_2.fastq.bz2 $ head SRR518891_1.fastq > 1K_SRR518891_1.fastq $ head SRR518891_2.fastq > 1K_SRR518891_2.fastq 15

16 R N A - s e q 解析 : データ TRY! シーケンスデータを確認 $ cd /home/admin1409/amelieff/ngs $ ll 16

17 TRY! R N A - s e q 解析 : クオリティコントロールシーケンスデータのクオリティを確認 FastQC を実行します $ mkdir rnaseq $ fastqc -o rnaseq -f fastq 1K_SRR518891_1.fastq 1K_SRR518891_2.fastq fastqc_report.html をウェブブラウザで開きます $ firefox rnaseq/1k_srr518891_1_fastqc/fastqc_report.html $ firefox rnaseq/1k_srr518891_2_fastqc/fastqc_report.html 以降の解析は片側のリードのみ使用します 17

18 TRY! R N A - s e q 解析 : クオリティコントロール PolyA/T tail が存在最初の 1 塩基 18

19 応用 ) P o l y A /T t a i l の混入 RNA-seq(mRNA) では 3' 末端に PolyA/T tail がついている転写物をシーケンシングするため 19

20 TRY! R N A - s e q 解析 : クオリティコントロール最初の 1 塩基を削除 fastx_trimmerの使用方法を確認する $ fastx_trimmer -h 20

21 TRY! R N A - s e q 解析 : クオリティコントロール最初の 1 塩基を削除 1 塩基目を削除する $ cd rnaseq $ fastx_trimmer -f 2 -i../1k_srr518891_1.fastq -o 1K_SRR518891_1_s.fastq -Q33 21

22 TRY! R N A - s e q 解析 : クオリティコントロール PolyA/T tail を除去 3' 端にAを5 連続以上含むリード数がどのくらいあるか調べる $ grep "AAAAA$" 1K_SRR518891_1_s.fastq wc -l 39 fastx_clipperの使用方法を確認する $ fastx_clipper h 22

23 TRY! R N A - s e q 解析 : クオリティコントロール PolyA/T tail を除去 PolyA/T tail を除去する $ fastx_clipper -a AAAAA -i 1K_SRR518891_1_s.fastq -o 1K_SRR518891_1_s_notail.fastq -Q 33 Prinseq など各シーケンスの PolyA/T の数に合わせて除去するソフトもあります 23

24 TRY! R N A - s e q 解析 : クオリティコントロールクオリティの低いリードを除外 3 末端からクオリティ 20 未満の塩基をトリミングし長さが 30 塩基未満になったリードを破棄するさらに 80% 以上の塩基がクオリティー 20 以上のリードのみを抽出する $ fastq_quality_trimmer -t 20 -l 30 -Q 33 -i 1K_SRR518891_1_s_notail.fastq fastq_quality_filter -q 20 -p 80 -Q 33 -o 1K_SRR518891_1_clean.fastq 24

25 TRY! R N A - s e q 解析 : クオリティコントロールクリーニング結果の確認 FastQCを実行します $ fastqc -f fastq 1K_SRR518891_1_clean.fastq fastqc_report.html をウェブブラウザで開きます $ firefox 1K_SRR518891_1_clean_fastqc/fastqc_report.html クリーニング前後のリード数とクオリティの変化を確認してください 25

26 TRY! R N A - s e q 解析 : クオリティコントロールクリーニング前クリーニング後サンプルや調整方法シーケンサの特徴にあわせてクリーニング項目や条件を工夫しています 26

27 応用 ) アダプタ配列を除外するソフト 1 cutadapt 2 FastX-Toolkit(fastxclipper) 3 tagcleaner 指定した配列はどのソフトでも除ける fastx_clipper は部分配列もかなり除けたがリード数も 1/10 以下に減るためアダプタ以外の配列も除いている可能性がある tagcleaner は一度に 1 アダプタ配列しか指定できない $ cutadapt -b TCTCGTATGCCGTCTTC -b CTACAGTCCGACGA -m 10 -n 2 $FILE.fastq 1> $FILE_cutadapt.fastq オプション m n O e これより短くなったものは破棄同リードへのアダプタ出現回数マッチ領域の最少長エラー塩基数/ マッチ領域長の最大 27

28 R N A - s e q 解析 : マッピング TRY! TopHat の使い方を確認 $ tophat スプライシングを考慮してマッピングするため既知の遺伝子情報を使用することもできます 28

29 R N A - s e q 解析 : マッピング TRY! マッピング今回はリード数が少ないためマッピング基準を緩めています $ tophat -o 1K_SRR g 3 -N 10 --read-edit-dist 10 --read-gap-length 10 /home/admin1409/genome/saccharomyces_cerevisiae/ncbi /build3.1/sequence/bowtie2index/genome 1K_SRR518891_1_clean.fastq BAMとインデックス BEDなどが作成されます $ ls 1K_SRR N/--read-mismatches --read-edit-dist --read-gap-length Final read alignments having more than these many mismatches are discarded. Final read alignments having more than these many edit distance are discarded. Final read alignments having more than these many total length of gaps are discarded. 29

30 R N A - s e q 解析 : マッピング TRY! マッピング率 $ less 1K_SRR518891/align_summary.txt マッピング率 30

ポイント )T o p H a t のアルゴリズム 1. リードをペアエンドでリファレンスにマッピングする 2.

マッピング結果をもとに転写構造をアセンブリングする http://en.wikipedia.

32 R N A - s e q 解析 : マッピング TRY! マッピング結果の可視化 $ samtools index accepted_hits.bam $ igv.sh accepted_hits.bam を IGV で表示してください 32

33 R N A - s e q 解析 : マッピング TRY! マッピング結果の可視化 Position の例 : III:177, ,500 33

34 ポイント ) 遺伝子の発現量遺伝子上にマップされたリード数長い遺伝子ほどマップされるリードは多くなる ( 遺伝子間のバイアス ) サンプル量の多いランほどマップされるリードは多くなる ( ラン間のバイアス ) これらのバイアスを補正してから発現量を比較する必要があります発現量としてよく使われる指標 RPKM(Reads Per Kilobase per Million mapped reads) FPKM(Fragments Per Kilobase of exon per Million mapped fragments) どちらも発現量をエクソン長と全マッピング数で補正した値 FPKM = raw counts 1,000,000 all reads 1,000 gene length 34

35 R N A - s e q 解析 : 発現定量 TRY! Cufflinks の使い方を確認 $ cufflinks アセンブルのガイドとして既知の遺伝子情報を使用することもできます 35

--min-frags-per-transfrag 2 -o 1K_SRR518891

36 R N A - s e q 解析 : 発現定量 TRY! 発現量を計算今回はリード数が少ないため検出基準を緩めています $ cufflinks --min-frags-per-transfrag 2 -o 1K_SRR K_SRR518891/accepted_hits.bam $ ll h 1K_SRR fpkm_tracking ファイルが作成されます --min-frags-per-transfrag minimum number of fragments needed for new transfrags 36

37 R N A - s e q 解析 : 発現定量 TRY! 発現量を計算 $ less 1K_SRR518891/genes.fpkm_tracking 4 列目が Gene ID 10 列目が FPKM です 37

応用 ) サンプル間比較サンプルごとに発現量を計算したあとサンプルごとに発現している遺伝子が違うため比較の基準とする遺伝子リストを作成します $ cuffmerge -o COMPARE -g Genes.gtf -s genome.fa transcript.gtf.txt Group1/S1/transcript.gtf Group1/S2/transcript.

38 応用 ) サンプル間比較サンプルごとに発現量を計算したあとサンプルごとに発現している遺伝子が違うため比較の基準とする遺伝子リストを作成します $ cuffmerge -o COMPARE -g Genes.gtf -s genome.fa transcript.gtf.txt Group1/S1/transcript.gtf Group1/S2/transcript.gtf Group2/S3/transcript.gtf Group2/S4/transcript.gtf transcript.gtf.txt 各サンプルの cufflinks 結果を羅列したファイル発現量を比較します $ cuffdiff -o COMPARE -L Group1, Group2 Genes.gtf Group1/S1/accepted_hits.bam, Group1/S2/accepted_hits.bam Group2/S3/accepted_hits.bam, Group2/S4/accepted_hits.bam 38

39 本講義の内容 Reseq 解析 RNA-seq 解析公開データ取得クオリティコントロールマッピング変異検出公開データ取得クオリティコントロールマッピング発現定量 39

40 ご清聴ありがとうございました 40

PowerPoint プレゼンテーション

PowerPoint プレゼンテーション V1 次世代シークエンサ実習 II 本講義にあたって代表的な解析の流れを紹介します論文でよく使用されているツールを使用しますコマンドを沢山実行しますスペルミスが心配な方はコマンド例がありますのでコピーして実行してください /home/admin1409/amelieff/ngs/reseq_command.txt マークのコマンドは実行してください実行が遅れてもあせらずに応用や課題の間に追い付いてください