Rでゲノム・トランスクリプトーム解析

Transcription

1 版 第 3 部 :NGS 解析 ( 中 ~ 上級 ) ~ クラウド環境との連携 ロングリードデータの解析 ~ 東京大学 大学院農学生命科学研究科アグリバイオインフォマティクス教育研究プログラム門田幸二 ( かどたこうじ ) kadota@iu.a.u-tokyo.ac.jp Aug 03 06, NGS ハンズオン講習会

2 忘れたらココを思い出してね 利用プログラムの簡単な解説 DDBJ Pipeline (Nagasaki et al., DNA Res., 03) マッピングや de novo アセンブリを行ってくれるウェブツール ( クラウド解析環境 ) FaQCs (Lo and Chain, BMC Bioinformatics, 04) Quality Control 用プログラム クオリティフィルタリングやアダプター除去が主目的 FastQC Quality Control 用プログラム アダプターの混入など NGS データのクオリティチェックが主目的 HGAP (Chin et al., Nat Methods, 03) PacBio 用 de novo ゲノムアセンブラ.bax.h5 ファイルのみを入力として受け付ける KmerGenie (Chikhi and Medvedev, Bioinformatics, 04) de novo ゲノムアセンブリ時に用いる最適な k 値を算出してくれる 推定ゲノムサイズも返す Platanus (Kajitani et al., Genome Res., 04) de novo ゲノムアセンブラ 複数の k 値を利用するので KmerGenie とは無関係 SRA Toolkit sra 形式ファイル処理用のプログラム群 FASTQ に変換する fastq-dump 利用が主目的 Velvet (Zerbino and Birney, Genome Res., 008) de novo ゲノムアセンブラ 単一の k 値を利用するので KmerGenie と関係あり Aug 03 06, NGS ハンズオン講習会

3 おさらい 自習 解析データは 乳酸菌 (L. hokkaidonensis LOOC60 T ) ゲノム塩基配列 原著論文 Aug 03 06, NGS ハンズオン講習会 Tanizawa et al., BMC Genomics, 6: 40, 05 3

4 おさらい PacBio RS II データ (DRR04500) DRR04500 は登録内容に問題があったことが判明し消滅 4 セル分のデータ DRR に差し替えられている セルあたり約 5 万リード 4 セル分なので約 60 万リード Illumina MiSeq データ (DRR0450) paired-end ゲノムデータ リード長は forward 側と reverse 側共に 50 bp オリジナルは,97,30 リード 最初の 300,000 リードを解析 forward 側 (DRR0450sub_.fastq.gz) reverse 側 (DRR0450sub_.fastq.gz) FaQCs 実行結果 ( 第 6 回 W5-4) 300,000 リード 97,633 リード (W5-) forward 側 (QC..trimmed.fastq.gz) reverse 側 (QC..trimmed.fastq.gz) 自習 ファイルの場所 :~/Documents/DRR0450/result 乳酸菌ゲノム配列決定論文は 種類の NGS 機器から得られたデータを併用している 第 6 回は Illumina MiSeq データ (DRR0450) を そして第 7 回は PacBio データを取り扱っている 50 bp 50 bp Aug 03 06, NGS ハンズオン講習会 Tanizawa et al., BMC Genomics, 6: 40, 05 4

5 自習 第 6 回原稿 PDF の p45 複数の k 値で Velvet アセンブリを行い 主観で k=7 がいいと判断したところまでが の内容 今日は Velvet など単一の k 値を指定してアセンブリを実行する際に 客観的に最もよいと思われる k 値を出力してくれる KmerGenie のインストールの話からスタート KmerGenie は make コマンドを利用してインストールする 例目として また最適 k 値と同時に出力する ゲノムサイズ推定周辺を主目的として紹介 KmerGenie 出力結果の一部には で眺めた k-mer 出現頻度分布もある はゲノムサイズ推定の意義について Aug 03 06, NGS ハンズオン講習会 5

6 Contents 自習 W: ゲノムサイズ推定 (KmerGenie のインストールと利用 ) インストール single-end で実行 ( 結果の解説 ) paired-end で実行 ( 結果の解説 ) W: 配列長によるフィルタリング (Python プログラムの利用 ) DDBJ Pipeline(W3 から W7 まではほぼ省略 ) W8:k=3 での Velvet 実行結果の解析 W9:Platanus の実行 W0: 結果の解析 W:ACGT カウント ( 塩基ごとの出現頻度解析 ) ロングリード (PacBio) データと公共 DB W: PacBio 生データは bax.h5 形式 公共 DB は sra 形式と FASTQ 形式 W3: 公共 DB (DRA) の FASTQ ファイルを入力として FastQC W4:NCBI SRA (SRA) が提供する SRA Toolkit のインストール W5: 利用 (.sra.fastq への変換 ) W6:FastQC DDBJ Pipeline で HGAP を実行 W7:DDBJ Pipeline に解析したい生データファイル (.bax.h5) をアップロード W8:DDBJ Pipeline で HGAP を実行 W9:HGAP アセンブリ結果を眺める Aug 03 06, NGS ハンズオン講習会 6

7 ゲノムサイズ推定 自習 ゲノムサイズ推定 からスタート 3KmerGenie スライドを見るだけ 3 Aug 03 06, NGS ハンズオン講習会 7

8 自習 W-:KmerGenie 06 年 6 月 5 日現在の KmerGenie のバージョンは.706 だが 第 6 回ウェブ資料作成当時 (06 年 月 5 日 ) の KmerGenie (ver..698) をダウンロードできるのでそれで説明する で右クリックし ショートカットのコピーで wget でダウンロードするときに必要な URL 情報を取得できる 講習会では ver..698 の tar.gz ファイルをダウンロード済みですが サイズは約 400KB なので 時間をずらして休憩時間に最新版をダウンロードしてみてもよい Aug 03 06, NGSハンズオン講習会 Chikhi and Medvedev, Bioinformatics, 30: 3-37, 04 8

9 自習 W-: ダウンロードと解凍 ダウンロードと解凍 講習会ではver..698のtar.gzファイルをダウンロード済み 実際に行うのは の確認のみ Aug 03 06, NGS ハンズオン講習会 9

10 ゲノムサイズ推定 自習 もし.tar.gz ファイルを間違って消してしまったり その後の wget でしくじった場合は でリカバリしてください Aug 03 06, NGS ハンズオン講習会 0

11 W-: ダウンロードと解凍 解凍.tar.gz ファイルなので 解凍コマンドは W9-3 と同じ ここから手を動かす Aug 03 06, NGS ハンズオン講習会

12 W-: ダウンロードと解凍 解凍が無事終わると.tar.gz を除いた部分の名前のディレクトリが作成される Aug 03 06, NGS ハンズオン講習会

13 W-3:README 基本的には の場所で make と打てばいいが 今一度マニュアル (README) を less で眺める less の利用法は 第 3 回 W4-6 ( スペースキーを押してページスクロール q で抜ける ) Aug 03 06, NGS ハンズオン講習会 3

14 W-3:README 基本的な使用法は のような感じで.fastq.gz も読み込めることを確認 Aug 03 06, NGS ハンズオン講習会 4

15 W-3:README 前のスライドからスペースキーを 回押したあたりの画面 の赤枠部分に注目 デフォルトでは KmerGenie は k= までしか探索しない 探索の上限を 00 まで引き上げたい場合は make clean と打ったのち make k=00 と打て と書いている 3q と打って less を終了 3 Aug 03 06, NGS ハンズオン講習会 5

16 W-4:make clean まずは make clean 実行前後で特に変わった様子はないようだ Aug 03 06, NGS ハンズオン講習会 6

17 W-5:make k=00 引き続いて make k=00 約 0 秒 Aug 03 06, NGS ハンズオン講習会 7

18 W-5:make k=00 無事インストールできたようだ これを make が無事に通る などと呼ぶ Aug 03 06, NGS ハンズオン講習会 8

19 W-5:make k=00 06 年 月 5 日に作業している 個だったのが 個に増えている specialk というものが増えたもののようだ ls l 3 確かに specialk という実行ファイルができたて 他にも 4minia というディレクトリの内容もどこかが変更されたようだ 4 3 Aug 03 06, NGS ハンズオン講習会 9

20 W-6:README 利用法を見るべく もう一度 README ファイルを less する Aug 03 06, NGS ハンズオン講習会 0

21 W-6:README これが基本的な使い方 赤下線左端の./ を見て kmergenie のパスを通さなければ 絶対パスで書かないといけないのか と認識する そして 赤枠内の Optionally 以下の文章に目が留まる 実行プログラム ( この場合 kmergenie) のシンボリックリンクを /usr/local/bin に置けばいいことは第 4 回でも説明した 第 4 回 W9-5 で紹介した ln s の記述法は結構面倒だが make install と打てばいいらしいと解釈 Aug 03 06, NGS ハンズオン講習会

22 W-7:make install (q と打って less から抜けて )make install エラーが出て失敗していることがわかる この理由は 第 4 回 W9-5 でも説明しているが /usr/local/bin ディレクトリの所有者は root さんで root さん以外のヒトは書き込み権限が与えられていない ( つまり drwxr-xr-x となっており w がない ) そこに 3iu さんが書き込みを行おうとしたから Permission denied( 権限がない ) と文句を言われたというオチです 3 Aug 03 06, NGS ハンズオン講習会

23 W-7:make install sudo make install 赤下線のところでパスワードを打ち込むよう促されたら pass409 エラーが出ずに /usr/local/bin に kmergenie の実行ファイルを置けたようだ ( 無事パスを通せたようだ ) 記述内容を眺めることで make install の実体は /usr/local/bin/kmergenie があったならそれを rm で消して ln -s でパスを通す一連のコマンド実行のようだと学習する Aug 03 06, NGS ハンズオン講習会 3

24 W-7:make install パスを通す前は kmergenie のヘルプを表示させようとしても そんなコマンドはない と文句を言われるが 3 パスを通した後は文句を言われません./kmergenie と kmergenie の違いが分からないヒトは 第 4 回 W9 で復習 3 Aug 03 06, NGS ハンズオン講習会 4

25 W-7:kmergenie -h 4 kmergenie h の結果を眺める 乳酸菌 (haploid; 一倍体 ) の場合は気にしなくていいが ヒトなどの二倍体ゲノムの場合は --diploid オプションをつけないといけないようだ W0-6 で Velvet の結果として得たのは k=3 から 9 であった この範囲で best k-mer を探索したい場合は 3 l 3 -k 9 オプションを 4 の位置につけて実行すればいいようだ 3 Aug 03 06, NGS ハンズオン講習会 5

26 Contents W: ゲノムサイズ推定 (KmerGenie のインストールと利用 ) インストール single-end で実行 ( 結果の解説 ) paired-end で実行 ( 結果の解説 ) W: 配列長によるフィルタリング (Python プログラムの利用 ) DDBJ Pipeline(W3 から W7 まではほぼ省略 ) W8:k=3 での Velvet 実行結果の解析 W9:Platanus の実行 W0: 結果の解析 W:ACGT カウント ( 塩基ごとの出現頻度解析 ) ロングリード (PacBio) データと公共 DB W: PacBio 生データは bax.h5 形式 公共 DB は sra 形式と FASTQ 形式 W3: 公共 DB (DRA) の FASTQ ファイルを入力として FastQC W4:NCBI SRA (SRA) が提供する SRA Toolkit のインストール W5: 利用 (.sra.fastq への変換 ) W6:FastQC DDBJ Pipeline で HGAP を実行 W7:DDBJ Pipeline に解析したい生データファイル (.bax.h5) をアップロード W8:DDBJ Pipeline で HGAP を実行 W9:HGAP アセンブリ結果を眺める Aug 03 06, NGS ハンズオン講習会 6

27 W-8: 実行 一通りの実行方法がわかったので とりあえず W5-4 で作成した forward 側のファイル (QC..trimmed.fastq.gz) のみを入力として 3 から 9 の k-mer の探索範囲で kmergenie を実行 約 5 分 Aug 03 06, NGS ハンズオン講習会 7

28 W-8: 途中経過 リターンキーを押してから 確か 0 秒後くらいの状態 確かにオプションで指定した通りの範囲を探索してくれているようだ Aug 03 06, NGS ハンズオン講習会 8

29 W-8: 無事終了 終了後の状態 KmerGenie プログラムの主な目的は Velvet などのアセンブリプログラム実行時に最適だと思われる k-mer の推定 forward 側のみで実行した結果は k=87 がお勧めだという結果だった とりあえず ls Aug 03 06, NGS ハンズオン講習会 9

30 W-9: 確認 画面が一気に流れる KmerGenie は沢山のファイルをカレントディレクトリ上に吐き出すようだ ファイル群 histograms-k*.histo* は 個別の k 値の結果を見たいときに必要となるのだろう この html ファイルさえ見れば基本的に OK ゲスト OS(Bio-Linux) 上で firefox histograms_report.html とやってもよいが 画面が小さく見づらい ここでは 共有フォルダ経由でホスト OS の share フォルダにコピーして ホスト OS 上で html ファイルを眺める Aug 03 06, NGS ハンズオン講習会 30

31 W-9: 確認 ホスト OS 上で html ファイルを眺めた画面 アセンブリ時に用いるベストな k 値は 87 だという結果 が推定ゲノムサイズ,356,73 bp ( 約.4MB) という結果 Aug 03 06, NGS ハンズオン講習会 3

32 KmerGenie 結果解説 もう少し詳しく説明 横軸が調べた k-mer の k 値 指定した探索範囲は 3 から 9 なので妥当 縦軸の名称は Number of genomic k-mers genomic という形容詞がついていることからも想像できるが これが第 部初日 ( ) の最後のほうで示した シークエンスエラー由来のものを除いた k-mer の種類数 です Aug 03 06, NGS ハンズオン講習会 3

33 KmerGenie 結果解説 のあたりを眺めることで k=30-30 くらいの範囲でどの k-mer を用いてもゲノムサイズが.3-.4MB という結果になるのだろうと解釈する 但し これは genomic k-mers の種類数 なのでシークエンスエラー由来 k-mer のフィルタリングが実装されていない ( 甘い ) アセンブリプログラムだとそうはいかないだろう とも想像する Aug 03 06, NGS ハンズオン講習会 33

34 KmerGenie 結果解説 html ファイル自体は結構縦長である のあたりから オプションで指定した 3 から 9 の k 値の個別の解析結果が見られる 例えば が k=3 のときの k-mer 出現頻度分布 Aug 03 06, NGS ハンズオン講習会 34

35 KmerGenie 結果解説 KmerGenie は のシークエンスエラー由来 k-mer を除いた残りの k-mer の種類数を genomic と判定しているのだろう Aug 03 06, NGS ハンズオン講習会 35

36 KmerGenie 結果解説 ほぼ余談 実際には の付近のやたら沢山出現する k-mer もフィルタリングしている ( はず ) このあたりの多くは リピート配列由来 k-mer の領域だからです Aug 03 06, NGS ハンズオン講習会 36

37 Contents W: ゲノムサイズ推定 (KmerGenie のインストールと利用 ) インストール single-end で実行 ( 結果の解説 ) paired-end で実行 ( 結果の解説 ) W: 配列長によるフィルタリング (Python プログラムの利用 ) DDBJ Pipeline(W3 から W7 まではほぼ省略 ) W8:k=3 での Velvet 実行結果の解析 W9:Platanus の実行 W0: 結果の解析 W:ACGT カウント ( 塩基ごとの出現頻度解析 ) ロングリード (PacBio) データと公共 DB W: PacBio 生データは bax.h5 形式 公共 DB は sra 形式と FASTQ 形式 W3: 公共 DB (DRA) の FASTQ ファイルを入力として FastQC W4:NCBI SRA (SRA) が提供する SRA Toolkit のインストール W5: 利用 (.sra.fastq への変換 ) W6:FastQC DDBJ Pipeline で HGAP を実行 W7:DDBJ Pipeline に解析したい生データファイル (.bax.h5) をアップロード W8:DDBJ Pipeline で HGAP を実行 W9:HGAP アセンブリ結果を眺める Aug 03 06, NGS ハンズオン講習会 37

38 W-0:paired-end paired-end の場合に どう指定すればいいのか調べる ここでは 入力ファイルがあるディレクトリ上で KmerGenie の README ファイルを less で開く 打ち込み派のヒトは タブ補完をうまく利用すべし less の利用法は 第 3 回 W4-6 Aug 03 06, NGS ハンズオン講習会 38

39 W-0:paired-end 複数のファイルに分かれている場合はリストファイルを作成して入力として与えればよい 入力ファイルの順番 (order) やリードの向き (orientation) は気にしなくていいようだ q で less から抜ける Aug 03 06, NGS ハンズオン講習会 39

40 W-0:paired-end ls すると histograms* のファイルが沢山出て見づらいので削除する 実用上は削除候補ファイル群を ls やタブ補完でリストアップして 削除前に確認する ( 第 4 回の W-) Aug 03 06, NGS ハンズオン講習会 40

41 W-0:paired-end rm f で削除 * はワイルドカードの つ ( 第 4 回の W-) ワイルドカードは使いこなせるとかなり便利 ちなみに当面の目的は 3 paired-end のリストファイルを作成したい ここでは ワイルドカードを駆使して美しく作成する Tips を紹介 基本戦略は リストアップしたいもののみ ls し それをリダイレクトでファイルに書きだす です 実際には多少余分なものがあっても書きだした上で vi などのテキストエディタで余分なものを削除したりします 3 Aug 03 06, NGS ハンズオン講習会 4

42 W-0: ワイルドカード QC.* だと 赤枠で示す目的の pairedend ファイル以外のものが表示される QC.*.tri* でも目的外のファイルが表示される ( この場合は QC.*.gz でもうまくいくが )3 任意の数字 字を表す [0-9] を利用することで 目的の paired-end ファイルのみリストアップさせることができる アセンブリプログラム Platanus (ver...4; Kajitani et al., 04) のマニュアルを見るときに役立つ 3 Kajitani et al., Genome Res., 4: , 04 Aug 03 06, NGSハンズオン講習会 Chikhi and Medvedev, Bioinformatics, 30: 3-37, 04 4

43 W-0: ワイルドカード 3 でリストアップできることを確認したら 4 リダイレクト (>) を利用して任意のファイル名 ( ここでは list.txt) で保存すれば リストファイルの完成 この種のテクニックは頻用する 3 4 Aug 03 06, NGS ハンズオン講習会 43

44 W-:paired-end 実行 リストファイル list.txt があることを確認して kmergenie コマンドを実行 画面は数秒後の状態 うまくリストファイルを読み込めているようだ 約 0 分 Aug 03 06, NGS ハンズオン講習会 44

45 W-: 無事終了 終了後の状態 forward 側のみ (singleend) で実行した結果 (k=87) と違って paired-end のときは k=4 がお勧めとなっていることがわかる とりあえず ls Aug 03 06, NGS ハンズオン講習会 45

46 W-: 確認 W-9 と同様に画面が一気に流れる ここでは (single-end の結果ファイルが histograms_report.html として既に存在するので ) 共有フォルダ経由でホスト OS の share フォルダに histograms_report_paired.html という名前でコピーして ホスト OS 上で html ファイルを眺める Aug 03 06, NGS ハンズオン講習会 46

47 W-: 確認 single-end (k=87) と paired-end (k=4) で推奨 k 値は大きく異なるものの ゲノムサイズの推定値はほとんど変わらないことがわかる (single-end は bp paired-end は bp) Aug 03 06, NGS ハンズオン講習会 47

48 Contents W: ゲノムサイズ推定 (KmerGenie のインストールと利用 ) インストール single-end で実行 ( 結果の解説 ) paired-end で実行 ( 結果の解説 ) W: 配列長によるフィルタリング (Python プログラムの利用 ) DDBJ Pipeline(W3 から W7 まではほぼ省略 ) W8:k=3 での Velvet 実行結果の解析 W9:Platanus の実行 W0: 結果の解析 W:ACGT カウント ( 塩基ごとの出現頻度解析 ) ロングリード (PacBio) データと公共 DB W: PacBio 生データは bax.h5 形式 公共 DB は sra 形式と FASTQ 形式 W3: 公共 DB (DRA) の FASTQ ファイルを入力として FastQC W4:NCBI SRA (SRA) が提供する SRA Toolkit のインストール W5: 利用 (.sra.fastq への変換 ) W6:FastQC DDBJ Pipeline で HGAP を実行 W7:DDBJ Pipeline に解析したい生データファイル (.bax.h5) をアップロード W8:DDBJ Pipeline で HGAP を実行 W9:HGAP アセンブリ結果を眺める Aug 03 06, NGS ハンズオン講習会 48

49 第 6 回原稿 PDF の p45 の部分の話です 目的は アセンブリ結果からの短い配列の除外 Aug 03 06, NGS ハンズオン講習会 49

50 W-: ダウンロード multi-fasta 形式ファイルを入力として 指定した配列長未満の配列を除く Python プログラム (fastalengthfilter.py; 第 6 回筆頭著者作 ) をダウンロード 簡単な自作プログラムは ~/bin に置き そこにパスを通して使っている ( ヒトもいる ) そうすることで毎回 /usr/local/bin などにシンボリックリンクをはる手間が省ける Aug 03 06, NGS ハンズオン講習会 50

51 W-: ダウンロード ホームディレクトリに移動し bin ディレクトリがないことを確認して mkdir で作成 3bin ディレクトリに移動し 4fastaLengthFilter.py が存在する URL を指定して wget 5( 実行権限があればそうなるが ) ダウンロードしたファイルが緑色になっていないので 6 実行権限 ( 第 4 回 W3- ) が自分 ( この場合 iu) にもないことを認識する Aug 03 06, NGS ハンズオン講習会 5

52 W-:permission chmod で実行権限を変更 オプションで 755 とすると rwxr-xr-x になる 最初の 3 文字分は自分の権限に関する部分 読み込み (r) 書き込み (w) 実行 (x) の全てを許可するという意味で rwx となっている ここでは 自分以外には書き込み権限を与えない設定にしている 私はいつも 755 Aug 03 06, NGS ハンズオン講習会 5

53 W-3: パスを通す ~/bin にパスを通すための現状確認 ( 後でホームディレクトリの.zshrc ファイルをいじるので ) ホームディレクトリに移動 echo $PATH で現在通っているパスの確認 ( 第 4 回 W9-9) 確かに現状では ~/bin に相当する /home/iu/bin は存在しない 3Bio- Linux のデフォルトは zsh であることを一応確認 3 Aug 03 06, NGS ハンズオン講習会 53

54 W-3: パスを通す.zshrc ファイルの存在確認と 最後の 5 行分を表示 3 テキストエディタ gedit (vi や emacs でもよい ) を用いて 赤枠部分に :/home/iu/bin を追加して保存する gedit の使い方は第 4 回 W0-3 にもあり ~/Desktop/backup に.zshrc_org があります 消しちゃったヒトはそれで復旧してください 3 Aug 03 06, NGS ハンズオン講習会 54

55 W-3: パスを通す 赤枠部分に :/home/iu/bin を追加して保存し gedit を終了した後の状態 4 もう一度.zshrc ファイルの最後の 5 行分を表示 3 4 Aug 03 06, NGS ハンズオン講習会 55

56 W-3: パスを通す 確かに :/home/iu/bin が追加されている 3 4 Aug 03 06, NGS ハンズオン講習会 56

57 W-4: 確認 gedit.zshrc を実行したこのターミナルは まだターミナル起動時の環境設定のまま source コマンドで編集後の環境設定ファイル (.zshrc) を再読み込みしておく その後もう一度 echo $PATH で確認すると 3 確かに自分が追加した赤枠のものが見えている 基本的にこれでもう OK 3 Aug 03 06, NGS ハンズオン講習会 57

58 W-4: 確認 参考 以下はおまけです よく見ると の赤下線部分が重複しており 左側にあるものと全く同じ やらかしてしまったのではないかと思うが 重複しているだけで実際には問題ない だって 3 ここにも重複があるが 問題ないから 重複を防ぐには.zshrc 編集時に余分に追加されそうなところを削除しておけばよい つまり 3 Aug 03 06, NGS ハンズオン講習会 58

59 W-4: 確認 参考 例えば ~/hoge ディレクトリ以下のファイルにパスを通したい場合は の $PATH の部分のみ残して 赤枠の ~/hoge ディレクトリの絶対パスを追加すればよい Aug 03 06, NGS ハンズオン講習会 59

60 W-5: パスと改行コード where コマンドでパスが通っていることを念のため確認 file コマンドで Linux の改行コードになっているかどうか確認 ASCII text executable で終わっている場合は問題ないが ASCII text executable, with CRLF line terminators となっている場合は 第 4 回 W3-6 を参考にして nkf コマンドを用いて Linux の改行コードに変換しておく Aug 03 06, NGS ハンズオン講習会 60

61 W-6: ヘルプを表示 -h をつけて fastalengthfilter.py を実行 このプログラムが何をするものなのかについての簡単な description および 3 基本的な利用法 (Usage) が表示される 4 実は は偶然うまくいっただけ このプログラムに関しては -h オプションを受け付けてはおらず 入力ファイル ( 第 引数 ) と threshold( 第 引数 ) の つのオプションを同時入力しなかった場合に と 3 が表示される つまり 4 が正解 3 4 Aug 03 06, NGS ハンズオン講習会 6

62 W-7: 実行準備 fastalengthfilter.py を実行すべく W0-5 で作成した Velvet アセンブリ結果ファイル (contigs.fa) に近いディレクトリに移動し ls Aug 03 06, NGS ハンズオン講習会 6

63 W-7: 実行準備 W- で作成した KmerGenie 実行結果ファイル群 (histograms*) が一気に表示されて見づらいので 削除 とりあえず k= で実行した Velvet アセンブリ結果フォルダ hoge_ 中の contigs.fa を例題としてフィルタリングを行う Aug 03 06, NGS ハンズオン講習会 63

64 W-7: 実行準備 hoge_/contigs.fa の最初の 8 行分を表示 description 行は NODE_... みたいな感じになっていることがわかる wc 実行結果を表示 行数は 7,60 3 配列数は 3,76 で W0-6 と同じ 3 Aug 03 06, NGS ハンズオン講習会 64

65 W-8: 実行 hoge_/contigs.fa を入力として 300 bp 以上の配列のみ残した結果を contigs_300.fa というファイル名で出力 Aug 03 06, NGS ハンズオン講習会 65

66 W-8: 実行 行数が 7,60 行から 3,874 行に激減している! 配列数も 3,76 個から,937 個に激減 3 総塩基数も 5,78,04 bp から (83,550,937 =) 8,63 bp に激減していることがわかる 3 Aug 03 06, NGS ハンズオン講習会 66

67 W-8: 実行 4 の手順でも同じ結果になることを確認しており プログラムのバグではない 3 4 Aug 03 06, NGS ハンズオン講習会 67

68 W-8: 実行 行数が 3,874 行 配列数が,937 個に激減 3 description 行以外の行数も,937 行 という結果に違和感を覚えるかもしれない これは fastalengthfilter.py の出力ポリシーは 塩基配列部分が 配列 行だからである 実際 4 最後の 4 行を出力させると 配列分表示される 3 4 Aug 03 06, NGS ハンズオン講習会 68

69 W-9: 全て実行 参考 同じ作業を 他の k 値のアセンブリ結果に対しても実行 ないヒトはやったつもりでエアーハンズオン 入力ファイルがあるディレクトリをカレントディレクトリにしなかった理由は 複数のディレクトリ上にあるファイルを一気に処理できるようにしたかったからです Aug 03 06, NGS ハンズオン講習会 69

70 参考 W-9: 結果を確認 フィルタリング結果ファイルに対して 行数や ファイルサイズを表示 ないヒトはエアーハンズオン Aug 03 06, NGS ハンズオン講習会 70

71 W-9: 結果を確認 実際に行っているのは 赤枠内のコピペ 系統的に記述できている ( 違いはディレクトリ名部分のみ ) ことがわかる 赤枠部分のみからなるファイルは 一般的なシェルスクリプトといえる 私はこういうものを作って一気に系統的な解析を行います Aug 03 06, NGS ハンズオン講習会 7

72 W-0: これまでのまとめ 配列長 (< 300 bp) によるフィルタリング前は 最大で約 5.7MB (k= での Velvet アセンブリ結果 ) というゲノムサイズに達していたが フィルタリング後は最大でも約.6MB (k=5 の結果 ) となっていることがわかる Aug 03 06, NGS ハンズオン講習会 7

73 Contents W: ゲノムサイズ推定 (KmerGenie のインストールと利用 ) インストール single-end で実行 ( 結果の解説 ) paired-end で実行 ( 結果の解説 ) W: 配列長によるフィルタリング (Python プログラムの利用 ) DDBJ Pipeline(W3 から W7 まではほぼ省略 ) W8:k=3 での Velvet 実行結果の解析 W9:Platanus の実行 W0: 結果の解析 W:ACGT カウント ( 塩基ごとの出現頻度解析 ) ロングリード (PacBio) データと公共 DB W: PacBio 生データは bax.h5 形式 公共 DB は sra 形式と FASTQ 形式 W3: 公共 DB (DRA) の FASTQ ファイルを入力として FastQC W4:NCBI SRA (SRA) が提供する SRA Toolkit のインストール W5: 利用 (.sra.fastq への変換 ) W6:FastQC DDBJ Pipeline で HGAP を実行 W7:DDBJ Pipeline に解析したい生データファイル (.bax.h5) をアップロード W8:DDBJ Pipeline で HGAP を実行 W9:HGAP アセンブリ結果を眺める Aug 03 06, NGS ハンズオン講習会 73

74 DDBJ Pipeline DDBJ Pipeline で Velvet が利用可能である 実際に k=3 で de novo アセンブリを行う一連の手順を示し Bio-Linux 上の数値 (W0-5 や W-9) と同じ結果が得られて安心するところまでが W3 から W8 の内容 ダイジェストで紹介 Aug 03 06, NGS ハンズオン講習会 74

75 第 6 回原稿 PDF の p47 以降の話です DDBJ Pipeline は ウェブブラウザ経由で ( しかも無料で ) 遺伝研スパコン上で NGS 解析ができる大変ありがたいツールです Aug 03 06, NGS ハンズオン講習会 75

76 W3-:DDBJ pipeline しかるべき原著論文を適切に引用するのは エンドユーザの義務です! スライドを見るだけ Nagasaki et al., DNA Res., 0: , 03 Aug 03 06, NGS ハンズオン講習会 76

77 W3-:DDBJ pipeline DDBJ Pipeline 利用時には アカウントを作成しておかねばなりません その際に入力する 氏名 所属 利用目的 は 遺伝研スパコンのウェブサイト上で公開されるのでご注意くださいスライドを見るだけ Aug 03 06, NGS ハンズオン講習会 77

78 W3-:DDBJ pipeline つまり DDBJ Pipeline の新規アカウント作成時の登録情報 ( 氏名 所属 利用目的 ) は 遺伝研スパコンの 登録ユーザ情報の 3 のところで公開される 企業の方は 受託解析に使うことはできません 研究開発用途としてお使いください スライドを見るだけ 3 Aug 03 06, NGS ハンズオン講習会 78

79 W3-4: 手元のファイル おさらい Bio-Linux で行った Velvet の入力ファイルは の FaQCs 実行結果ファイル これを DDBJ Pipeline 上に FTP 経由でアップロードして Velvet を実行しようとしている スライドを見るだけ Aug 03 06, NGS ハンズオン講習会 79

80 W3-4: ファイル名に注意 連載第 回原稿の バイオインフォマティクス分野の常識 非常識 を再掲 非常識なファイル名で実行を試みて うまくいかないんですけど 的な質問は スライドを見るだけ Aug 03 06, NGS ハンズオン講習会 80

81 W5-:Select Tools デフォルトはマッピングになっているので de novo Assembly にチェックを入れて ページ下部に移動 Aug 03 06, NGS ハンズオン講習会 8

82 W5-:Select Tools Velvet にチェックをいれて NEXT Aug 03 06, NGS ハンズオン講習会 8

83 W5-3:Set QuerySet 解析したいデータにチェックをいれて Set as Pair-End 3NEXT 本来この画面は 複数のクエリファイルを使用する場合に それらを連結して つのクエリセットとしてジョブを実行するか 別々のクエリセットとして並列して実行するかを指定する画面である 今回はクエリファイル つを単独で使用するので この画面にはあまり意味がない 灰色は読まない 3 Aug 03 06, NGS ハンズオン講習会 83

84 W5-4:Set Ass.Options W7-3 の Bio-Linux 上での実行経験から DDBJ Pipeline では k=3 がデフォルトなのだろう W7-8 の velvetg 実行時は特にオプションを指定しなかったが で示されているようなオプションもあるのだと学ぶ 3 長さによる配列のフィルタリング (4 デフォルトは 00 bp) もやってくれるようだ 3 4 Aug 03 06, NGS ハンズオン講習会 84

85 W5-4:Set Ass.Options W-0( スライド 7) を眺め とりあえず k=3 で Bio-Linux 上で指定したオプションと同じにして (~4) 実行 5NEXT Aug 03 06, NGS ハンズオン講習会 85

86 W5-4:Set Ass.Options オプションの見栄えはプログラムごとに異なる これは Velvet の場合 Aug 03 06, NGS ハンズオン講習会 86

87 W5-5:Confirmation アセンブリが終了したら アカウント作成時に指定したアドレス宛にメールが送られる RUN 3OK 3 Aug 03 06, NGS ハンズオン講習会 87

88 W6-: 計算終了 DDBJ Pipeline から計算終了メールが届いたら DDBJ Pipeline に再度ログインして計算結果を眺める Aug 03 06, NGS ハンズオン講習会 88

89 W8-: フィルタリング前の結果 ページ下部に移動した後の状態 Download(66.5MB) をクリックすると velvet.zipという zip 圧縮ファイルをダウンロードできる 共有フォルダに保存してBio-Linux 上で眺める とオリジナル ( 第 6 回ウェブ資料 ) はなっているが 講習会では ~/Desktop/backup 上にある velvet.zipを取り扱うので 以降の内容は若干異なる Aug 03 06, NGS ハンズオン講習会 89

90 W8-:Bio-Linux で解凍 3 ~/Desktop/backup にある velvet.zip の存在確認 3unzip で解凍 4velvet フォルダ中の contigs.fa が Velvet の生の出力ファイル このあたりは手打ちでやってください 4 Aug 03 06, NGS ハンズオン講習会 90

91 W8-3: 行数とコンティグ数 velvet ディレクトリに移動 contigs.fa の 行数は 75,35 3 配列 ( コンティグ ) 数は 8,398 個 このあたりも手打ち 3 Aug 03 06, NGS ハンズオン講習会 9

92 W8-3: 行数とコンティグ数 velvet ディレクトリに移動 contigs.fa の 行数は 75,35 3 配列 ( コンティグ ) 数は 8,398 個 4 ウェブ上の数値と同じで安心 3 4 Aug 03 06, NGS ハンズオン講習会 9

93 W8-4: フィルタリング DDBJ Pipeline 実行結果ファイルを入力として 指定した配列長閾値未満の配列を取り除く Python プログラム fastalengthfilter.py を実行 300 bp 以上の配列は,449 個となり Bio-Linux 上で実行した結果と同じ (W-9) Aug 03 06, NGS ハンズオン講習会 93

94 Contents W: ゲノムサイズ推定 (KmerGenie のインストールと利用 ) インストール single-end で実行 ( 結果の解説 ) paired-end で実行 ( 結果の解説 ) W: 配列長によるフィルタリング (Python プログラムの利用 ) DDBJ Pipeline(W3 から W7 まではほぼ省略 ) W8:k=3 での Velvet 実行結果の解析 W9:Platanus の実行 W0: 結果の解析 W:ACGT カウント ( 塩基ごとの出現頻度解析 ) ロングリード (PacBio) データと公共 DB W: PacBio 生データは bax.h5 形式 公共 DB は sra 形式と FASTQ 形式 W3: 公共 DB (DRA) の FASTQ ファイルを入力として FastQC W4:NCBI SRA (SRA) が提供する SRA Toolkit のインストール W5: 利用 (.sra.fastq への変換 ) W6:FastQC DDBJ Pipeline で HGAP を実行 W7:DDBJ Pipeline に解析したい生データファイル (.bax.h5) をアップロード W8:DDBJ Pipeline で HGAP を実行 W9:HGAP アセンブリ結果を眺める Aug 03 06, NGS ハンズオン講習会 94

95 W9-:Select Tools DDBJ Pipeline では Platanus (ver...) も実行可能 チェックを入れて NEXT Aug 03 06, NGS ハンズオン講習会 95

96 W9-4:Set Ass. Options 乳酸菌ゲノム決定論文中では Platanus assembler ver. with the default settings と書かれている 赤枠の Step,, 3 の計 3 か所にオプションを指定する箇所があるが とりあえずここは空白として 細かい点はすっ飛ばして NEXT Aug 03 06, NGS ハンズオン講習会 96

97 W0-: 結果を眺める Platanus は 主に 3 ステップからなる (W9-4) が そのコマンドの実体がわかる 実行ログを眺めると k=3, 4, 5 などをいろいろ調べているのだろうと想像がつく 最後のほうが k=7, 0,,, 3 で終わっている Platanus はこれらの k 値の結果を全て取り込んだアセンブリ結果を返す multi-k genome assembler のカテゴリーに属する Velvet と違って k 値を指定するオプションが存在しないのは k 値の探索範囲も自動的に決められるから Aug 03 06, NGS ハンズオン講習会 97

98 W0-: 結果を眺める フィルタリング前の Platanus 実行結果ファイルは ここからダウンロードできる (de novo アセンブリの場合は ) おそらくどのダウンロードボタンを押しても同じものが得られると思うが 無難に最後のステップのところのものを選択 このあたりはプログラムごとに若干仕様が異なるため 統一的に見せるべくこのようになっているのだろう Aug 03 06, NGS ハンズオン講習会 98

99 W0-: 結果を眺める Platanus 実行結果の基本情報はここからみられる パッと見 フィルタリング前の段階ですでにコンティグ数も明らかに少なくよさそう 赤枠部分を拡大して Velvet の結果と比較 Aug 03 06, NGS ハンズオン講習会 99

100 W0-3:Velvet と比較 左表で調べた範囲の Velvet 実行結果 (W- 0) よりも 特に k 値を意識すらせずに実行した Platanus の結果のほうが明らかに優れていることがわかる 具体的には Velvet での最少配列数は 68 個 (k=7) であったが Platanus はフィルタリング前の状態で 7 個 そして 3Minimum contig size が 0 bp であることから 300 bp 未満でフィルタリングを行うと 間違いなく配列数が減ると予想する 3 Aug 03 06, NGS ハンズオン講習会 00

101 W0-4: ダウンロード Platanus 実行結果ファイル (platanusresult.zip) を 共有フォルダにダウンロード とオリジナル ( 第 6 回ウェブ資料 ) はなっているが 講習会では ~/Desktop/backup 上にある platanusresult.zip を取り扱うので 以降の内容は若干異なる Aug 03 06, NGS ハンズオン講習会 0

102 W0-5:Bio-Linux で解凍 ~/Desktop/backup にある platanusresult.zip の存在確認 このあたりは手打ちでやってください Aug 03 06, NGS ハンズオン講習会 0

103 W0-5:Bio-Linux で解凍 unzip コマンドで解凍 -q をつけて quiet で解凍 ls で確認 d オプションをつけることで ディレクトリの中身を表示させない Aug 03 06, NGS ハンズオン講習会 03

104 W0-5:Bio-Linux で解凍 platanusresult ディレクトリに移動し ls Platanus の最終結果ファイルは 3out_gapClosed.fa 3 Aug 03 06, NGS ハンズオン講習会 04

105 W0-6: 配列数確認 拡張子が.fa のファイルを表示 *.fa の配列数を一気に出力 3 確かに Platanus の最終結果ファイルである out_gapclosed.fa の配列数は 7 個であり 4DDBJ Pipeline のウェブサイト上の数値と一致する 3 4 Aug 03 06, NGS ハンズオン講習会 05

106 W0-6: 配列数確認 赤下線の配列数は のスライド 50 と同じです この 3 つは Platanus による 3 ステップからなる de novo アセンブリの 各ステップ実行後の出力ファイルでした Aug 03 06, NGS ハンズオン講習会 06

107 de novo アセンブリ 入力 :paired-end FASTQ ファイル おさらい :de novo アセンブリの手順は大まかに 3 つのステップからなる Step: Assembly contig contig contig3 contig4 contig5 Step: Scaffold NNNNN NNNN NNNNNN scaffold scaffold Step3: Gap close CA A T C GG G TA A NNNCA TNNC GGNNTA scaffold Aug 03 06, NGS ハンズオン講習会 scaffold 07

108 de novo アセンブリ 入力 :paired-end FASTQ ファイル 3 の 3 つは 各ステップ実行後のファイルだったことを思い出そう 配列数の変遷が妥当である という議論を にしました Step: Assembly out_contig.fa (349 配列 ) contig contig contig3 contig4 contig5 Step: Scaffold out_scaffold.fa (7 配列 ) NNNNN NNNN NNNNNN scaffold scaffold Step3: Gap close out_gapclosed.fa (7 配列 ) 3 CA A T C GG G TA A NNNCA TNNC GGNNTA scaffold Aug 03 06, NGS ハンズオン講習会 scaffold 08

109 参考 W0-7: フィルタリング out_gapclosed.fa を入力として 300 bp 未満の配列を除去するプログラムを実行 配列数は 5 個 3 原著論文中の結果 (53 配列 ) と似た個数が得られたことがわかる Aug 03 06, NGS ハンズオン講習会 3 Tanizawa et al., BMC Genomics, 6: 40, 05 09

110 W0-7: 総塩基数 参考 grep v > は > を含まない行を出力する (W8-; W0; W-8 など ) 300 bp 未満をフィルタリングする前は総塩基数が (,385, ) =,356,09 bp であった フィルタリングして配列数が 7 5 個に減ったあとの総塩基数は (,346,55 5) =,346,500 bp Aug 03 06, NGS ハンズオン講習会 0

111 Contents W: ゲノムサイズ推定 (KmerGenie のインストールと利用 ) インストール single-end で実行 ( 結果の解説 ) paired-end で実行 ( 結果の解説 ) W: 配列長によるフィルタリング (Python プログラムの利用 ) DDBJ Pipeline(W3 から W7 まではほぼ省略 ) W8:k=3 での Velvet 実行結果の解析 W9:Platanus の実行 W0: 結果の解析 W:ACGT カウント ( 塩基ごとの出現頻度解析 ) ロングリード (PacBio) データと公共 DB W: PacBio 生データは bax.h5 形式 公共 DB は sra 形式と FASTQ 形式 W3: 公共 DB (DRA) の FASTQ ファイルを入力として FastQC W4:NCBI SRA (SRA) が提供する SRA Toolkit のインストール W5: 利用 (.sra.fastq への変換 ) W6:FastQC DDBJ Pipeline で HGAP を実行 W7:DDBJ Pipeline に解析したい生データファイル (.bax.h5) をアップロード W8:DDBJ Pipeline で HGAP を実行 W9:HGAP アセンブリ結果を眺める Aug 03 06, NGS ハンズオン講習会

112 W-:ACGT カウント 参考 赤枠部分がA, C, G, T, Nの出現回数をカウントするRコード Step,, 3の配列長フィルタリング前の出力結果ファイルを入力としている のコピペ実行結果が次のスライド のスライド30-49で ホストOSの R 上で行った作業と同じことがLinux 上でもできることがわかってもらえればそれでよし スライド 0 まで省略 Aug 03 06, NGS ハンズオン講習会

113 参考 W-:ACGT カウント ACGT のカウントを調べたいファイルがあるディレクトリ上で R を起動 ( 第 5 回 W9-8) library(biostrings) と打って リターンキーを押す Aug 03 06, NGS ハンズオン講習会 3

114 参考 W-:ACGT カウント エラーなく読み込み完了 これで Biostrings パッケージが提供する ACGT の文字列をカウントする alphabetfrequency 関数などを利用可能な状態になった Aug 03 06, NGS ハンズオン講習会 4

115 参考 W-:ACGT カウント まずは PlatanusのStep 実行結果ファイル (out_contig.fa) を入力として A, C, G, T, Nの出現回数をカウントするRコード をコピペ実行 Aug 03 06, NGS ハンズオン講習会 5

116 参考 W-:ACGT カウント コピペ実行結果 の実行結果部分が A, C, G, T, N の塩基ごとの出現回数 N は つもないことがわかる Aug 03 06, NGS ハンズオン講習会 6

117 参考 W-:ACGT カウント の実行結果部分は 総塩基数を計算しているところ alphabetfrequency 関数実行によって得られた結果をhogeというものに格納している このhogeを入力として総和を計算する sum 関数を適用して全塩基数とみなしている Aug 03 06, NGS ハンズオン講習会 7

118 参考 W-:ACGT カウント Step 実行結果ファイル (out_scaffold.fa) を入力として実行した結果 Scaffolding によってコンティグ間が未知塩基 N で埋められたギャップで表現されるため N が存在 (49 個 ) するのは妥当 Aug 03 06, NGS ハンズオン講習会 8

119 W-:ACGT カウント 参考 Step3 実行結果ファイル (out_gapclosed.fa) を入力として実行した結果 Gap closingのおかげで この場合はNが0 個になっている 3 総塩基数,356,09 bpという数値は W0-7のwcコマンドから得られる結果と同じ 同じ目的を達成する上でも様々な手段がある 3 Aug 03 06, NGS ハンズオン講習会 9

120 W-:R の終了 参考 R の終了 ( 第 5 回 W9-8) Aug 03 06, NGS ハンズオン講習会 0

121 Contents W: ゲノムサイズ推定 (KmerGenie のインストールと利用 ) インストール single-end で実行 ( 結果の解説 ) paired-end で実行 ( 結果の解説 ) W: 配列長によるフィルタリング (Python プログラムの利用 ) DDBJ Pipeline(W3 から W7 まではほぼ省略 ) W8:k=3 での Velvet 実行結果の解析 W9:Platanus の実行 W0: 結果の解析 W:ACGT カウント ( 塩基ごとの出現頻度解析 ) ロングリード (PacBio) データと公共 DB W: PacBio 生データは bax.h5 形式 公共 DB は sra 形式と FASTQ 形式 W3: 公共 DB (DRA) の FASTQ ファイルを入力として FastQC W4:NCBI SRA (SRA) が提供する SRA Toolkit のインストール W5: 利用 (.sra.fastq への変換 ) W6:FastQC DDBJ Pipeline で HGAP を実行 W7:DDBJ Pipeline に解析したい生データファイル (.bax.h5) をアップロード W8:DDBJ Pipeline で HGAP を実行 W9:HGAP アセンブリ結果を眺める Aug 03 06, NGS ハンズオン講習会

122 おさらい PacBio RS II データ (DRR04500) DRR04500 は登録内容に問題があったことが判明し消滅 4 セル分のデータ DRR に差し替えられている セルあたり約 5 万リード 4 セル分なので約 60 万リード Illumina MiSeq データ (DRR0450) paired-end ゲノムデータ リード長は forward 側と reverse 側共に 50 bp オリジナルは,97,30 リード 最初の 300,000 リードを解析 forward 側 (DRR0450sub_.fastq.gz) reverse 側 (DRR0450sub_.fastq.gz) FaQCs 実行結果 ( 第 6 回 W5-4) 300,000 リード 97,633 リード (W5-) forward 側 (QC..trimmed.fastq.gz) reverse 側 (QC..trimmed.fastq.gz) ファイルの場所 :~/Documents/DRR0450/result 乳酸菌ゲノム配列決定論文は 種類の NGS 機器から得られたデータを併用している 第 6 回は Illumina MiSeq データ (DRR0450) を 連載第 7 回は PacBio データを取り扱っている Aug 03 06, NGS ハンズオン講習会 Tanizawa et al., BMC Genomics, 6: 40, 05

123 NGS 連載第 7 回は このあたりからスタート ロングリードの代表格である PacBio データとその周辺の話 オリジナルの第 7 回ウェブ資料 PDF と若干順番を入れ替えながら まずは事実の羅列から話を進めます Aug 03 06, NGS ハンズオン講習会 3

124 おさらい (NGS 連載第 3 回 ) 3 大公共 DB(NGS 分野 ) の提供ファイル形式 DDBJ SRA (DRA):sra 形式と FASTQ 形式 (bzip 圧縮 ) EMBL-EBI ENA (ENA):FASTQ 形式 (gzip 圧縮 ) NCBI SRA (SRA):sra 形式 大元は sra 形式なので FASTQ ファイルは sra ファイルから作成する (sra FASTQ) SRA Toolkit 中の fastq-dump で sra FASTQ の変換を行う FASTQ 提供サイト (DRA と ENA) は fastq-dump 実行時にクオリティフィルタリングなどを行っている このため FASTQ ファイル中のリード数は 元の sra ファイル中のリード数よりも若干少なくなる ( 第 3 回 W4) 実際 SRR6668 の sra ファイルは 35,073,834 リード FASTQ ファイルは 34,755,996 リードであった ( 第 3 回 W4-3) の第 3 回で行った議論は 基本的に Illumina データの話 PacBio データには通用しない 第 3 回 W0 から W 全部読んで説明 Aug 03 06, NGS ハンズオン講習会 4

125 W-7:bax.h5 ファイル まず PacBio の生データは sra でも FASTQ でもなく bax.h5 という形式のファイル ( 正確には PacBio RS II という機器が出力するファイル形式 ) サイズも巨大 具体的には セル分のみでも (747MB + 766MB + 90MB) =,44MB ( 約.4GB) に達する 3 これは DRR0543 の セル分のデータだが ファイル / セルではなく 3 つに分割されている点も最初は戸惑うポイント 3 Aug 03 06, NGS ハンズオン講習会 5

126 乳酸菌ゲノム配列決定論文は PacBio RS II で4セル分のデータを取得し de novoアセンブリを行った DRA 上では セルごとにDRR という計 4つのIDで登録されている その内のつであるDRR0543について解説しました DRR04500は登録内容に問題があったことが判明し消滅 4セル分のデータ DRR に差し替えられている おさらい PacBio RS II データ (DRR04500) セルあたり約 5 万リード 4 セル分なので約 60 万リード Illumina MiSeq データ (DRR0450) paired-end ゲノムデータ リード長は forward 側と reverse 側共に 50 bp オリジナルは,97,30 リード 最初の 300,000 リードを解析 forward 側 (DRR0450sub_.fastq.gz) reverse 側 (DRR0450sub_.fastq.gz) FaQCs 実行結果 ( 第 6 回 W5-4) 300,000 リード 97,633 リード (W5-) forward 側 (QC..trimmed.fastq.gz) reverse 側 (QC..trimmed.fastq.gz) ファイルの場所 :~/Documents/DRR0450/result Aug 03 06, NGS ハンズオン講習会 Tanizawa et al., BMC Genomics, 6: 40, 05 6

127 W-7:bax.h5 ファイル DDBJ Pipeline 上では PacBio 用のde novoアセンブリ用プログラムhgapを利用可能 しかし HGAPはbax.h5ファイルのみしか受け付けない HGAPを内包するPacBio 提供のSMRT Analysisシステムのインストールは超高難易度 ( 個人の感想です ) しかも数百 GBメモリ搭載マシンでないと動かせない (DDBJ Pipeline 上でDRR0543のみを入力としてHGAPを実行しても0GB 程度のメモリを要する ) 共同研究などで他人に頼むなど以外の場合 通常はDDBJ Pipeline 上でのHGAP 実行一択 Aug 03 06, NGS ハンズオン講習会 7

128 おさらい (NGS 連載第 3 回 ) 3 大公共 DB(NGS 分野 ) の提供ファイル形式 DDBJ SRA (DRA):sra 形式と FASTQ 形式 (bzip 圧縮 ) EMBL-EBI ENA (ENA):FASTQ 形式 (gzip 圧縮 ) NCBI SRA (SRA):sra 形式 大元は sra 形式なので FASTQ ファイルは sra ファイルから作成する (sra FASTQ) SRA Toolkit 中の fastq-dump で sra FASTQ の変換を行う FASTQ 提供サイト (DRA と ENA) は fastq-dump 実行時にクオリティフィルタリングなどを行っている このため FASTQ ファイル中のリード数は 元の sra ファイル中のリード数よりも若干少なくなる ( 第 3 回 W4) 実際 SRR6668 の sra ファイルは 35,073,834 リード FASTQ ファイルは 34,755,996 リードであった ( 第 3 回 W4-3) 公共 DB 上では bax.h5 ファイルは公開されていません sra や FASTQ ファイルが手元にあっても実質的には (DDBJ Pipeline で HGAP を実行できないので ) 無意味 しかも ( クオリティ値が PacBio よりも高い )Illumina データの場合は sra から FASTQ への変換時にリード数が若干減る程度だが Aug 03 06, NGS ハンズオン講習会 8

129 連載第 7 回 W3 3 大公共 DB(NGS 分野 ) の提供ファイル形式 DDBJ SRA (DRA):sra 形式とFASTQ 形式 (bzip 圧縮 ) EMBL-EBI ENA (ENA):FASTQ 形式 (gzip 圧縮 ) NCBI SRA (SRA):sra 形式 大元は sra 形式なので FASTQ ファイルは sra ファイルから作成する (sra FASTQ) SRA Toolkit 中の fastq-dump で sra FASTQ の変換を行う FASTQ 提供サイト (DRA と ENA) は fastq-dump 実行時にクオリティフィルタリングなどを行っている このため FASTQ ファイル中のリード数は 元の sra ファイル中のリード数よりも若干少なくなる ( 第 3 回 W4) 実際 SRR6668 の sra ファイルは 35,073,834 リード FASTQ ファイルは 34,755,996 リードであった ( 第 3 回 W4-3) PacBio データの場合は相当減る 具体的には 例えば DRA 上で同じ DRR0543 をダウンロードしても sra ファイルは 63,48 リード ( 実際には若干異なる?!) であるのに対し FASTQ ファイルだと 95 リードという結果になる これは減ったリード数ではなく 生き残ったリード数 Aug 03 06, NGS ハンズオン講習会 9

130 これはおそらく fastq-dump 実行時に指定しているオプションの閾値が Illuminaデータを合理的にフィルタリングするための条件のままで統一的にPacBioデータに対しても適用しているためであろう ( 推測の域 DDBJ SRA (DRA):sra 形式とFASTQ 形式 (bzip 圧縮を出ないが 論理的に考えればそのはず ) ) EMBL-EBI ENA (ENA):FASTQ 形式 (gzip 圧縮 ) 連載第 7 回 W3 3 大公共 DB(NGS 分野 ) の提供ファイル形式 NCBI SRA (SRA):sra 形式 大元は sra 形式なので FASTQ ファイルは sra ファイルから作成する (sra FASTQ) SRA Toolkit 中の fastq-dump で sra FASTQ の変換を行う FASTQ 提供サイト (DRA と ENA) は fastq-dump 実行時にクオリティフィルタリングなどを行っている このため FASTQ ファイル中のリード数は 元の sra ファイル中のリード数よりも若干少なくなる ( 第 3 回 W4) 実際 SRR6668 の sra ファイルは 35,073,834 リード FASTQ ファイルは 34,755,996 リードであった ( 第 3 回 W4-3) Aug 03 06, NGS ハンズオン講習会 30

131 DRA 上でダウンロードしたDRR0543の FASTQファイルが95リードになるのを FastQCで確認します ( 第 7 回 W3-) 次に ( 第 4 回 W4-5, W3-5, W4, W5でも紹介した ) apt-getを用いたプログラムインストールの応 DDBJ SRA (DRA):sra 形式とFASTQ 形式 (bzip 圧縮用編 ) ( 第 7 回 W4) として SRA Toolkitを解説 EMBL-EBI ENA (ENA):FASTQ 形式 (gzip 圧縮 ) この後の展開は 3 大公共 DB(NGS 分野 ) の提供ファイル形式 NCBI SRA (SRA):sra 形式 大元は sra 形式なので FASTQ ファイルは sra ファイルから作成する (sra FASTQ) SRA Toolkit 中のfastq-dumpで sra FASTQ の変換を行う FASTQ 提供サイト (DRAとENA) は fastq-dump 実行時にクオリティフィルタリングなどを行っている このため FASTQファイル中のリード数は 元のsraファイル中のリード数よりも若干少なくなる ( 第 3 回 W4) 実際 SRR6668のsraファイルは35,073,834リード FASTQファイルは34,755,996リードであった ( 第 3 回 W4-3) Aug 03 06, NGS ハンズオン講習会 3

132 SRA Toolkitは 第 7 回のオリジナルウェブ資料 (W5) ではFASTQファイルを生成して FastQCで確認する用途として利用しています しかし前述のように PacBioは入力がbax.h5 ファイルなので PacBioデータの性質を確認 DDBJ SRA (DRA):sra 形式とFASTQ 形式 (bzip 圧縮できてよかったね の域を出ません そして ) EMBL-EBI ENA (ENA):FASTQ 形式 (gzip 圧縮 ) Illuminaの場合は最初からFASTQファイルを NCBI SRA (SRA):sra 形式使えばよいので sraファイルを取り扱う意義は見出せません ( 第 3 回原稿 PDFのp36 左下 ) 念のため説明 3 大公共 DB(NGS 分野 ) の提供ファイル形式 大元は sra 形式なので FASTQ ファイルは sra ファイルから作成する (sra FASTQ) SRA Toolkit 中のfastq-dumpで sra FASTQ の変換を行う FASTQ 提供サイト (DRAとENA) は fastq-dump 実行時にクオリティフィルタリングなどを行っている このため FASTQファイル中のリード数は 元のsraファイル中のリード数よりも若干少なくなる ( 第 3 回 W4) 実際 SRR6668のsraファイルは35,073,834リード FASTQファイルは34,755,996リードであった ( 第 3 回 W4-3) Aug 03 06, NGS ハンズオン講習会 3

133 SRA Toolkitは あくまでもapt-getを用いたプログラムインストール ( 応用編 ; 第 7 回 W4) の例題という位置づけ apt-get( およびapt- cache) を使いこなして効率的にインストールを行うスキルがあれば de novo transcriptome DDBJ SRA (DRA):sra 形式とFASTQ 形式 (bzip 圧縮 assembly ) プログラムのつであるTrinityのイン EMBL-EBI ENA (ENA):FASTQ 形式 (gzip 圧縮 ) ストールも行えます ( ) NCBI SRA (SRA):sra 形式 念のため説明 3 大公共 DB(NGS 分野 ) の提供ファイル形式 大元は sra 形式なので FASTQ ファイルは sra ファイルから作成する (sra FASTQ) SRA Toolkit 中のfastq-dumpで sra FASTQ の変換を行う FASTQ 提供サイト (DRAとENA) は fastq-dump 実行時にクオリティフィルタリングなどを行っている このため FASTQファイル中のリード数は 元のsraファイル中のリード数よりも若干少なくなる ( 第 3 回 W4) 実際 SRR6668のsraファイルは35,073,834リード FASTQファイルは34,755,996リードであった ( 第 3 回 W4-3) Aug 03 06, NGS ハンズオン講習会 33

134 Contents W: ゲノムサイズ推定 (KmerGenie のインストールと利用 ) インストール single-end で実行 ( 結果の解説 ) paired-end で実行 ( 結果の解説 ) W: 配列長によるフィルタリング (Python プログラムの利用 ) DDBJ Pipeline(W3 から W7 まではほぼ省略 ) W8:k=3 での Velvet 実行結果の解析 W9:Platanus の実行 W0: 結果の解析 W:ACGT カウント ( 塩基ごとの出現頻度解析 ) ロングリード (PacBio) データと公共 DB W: PacBio 生データは bax.h5 形式 公共 DB は sra 形式と FASTQ 形式 W3: 公共 DB (DRA) の FASTQ ファイルを入力として FastQC W4:NCBI SRA (SRA) が提供する SRA Toolkit のインストール W5: 利用 (.sra.fastq への変換 ) W6:FastQC DDBJ Pipeline で HGAP を実行 W7:DDBJ Pipeline に解析したい生データファイル (.bax.h5) をアップロード W8:DDBJ Pipeline で HGAP を実行 W9:HGAP アセンブリ結果を眺める Aug 03 06, NGS ハンズオン講習会 34

135 W3-:FASTQ ダウンロード DRR0543 FASTQ 3bzip 圧縮 FASTQ ファイルをダウンロード 右クリックで ショートカットのコピー などで URL 情報を取得 ( 第 4 回 W9- や W8-) して wget してもよいし 共有フォルダ経由で Bio-Linux 上に置いてもよい とオリジナルのウェブ資料は書いてあるが 講習会用に既に ~/Desktop/backup にダウンロード済み それをコピーして利用 3 Aug 03 06, NGS ハンズオン講習会 35

136 W3-:FASTQ ダウンロード 今はこのあたりの話をしている W3- という情報を頼りに探す 3 赤枠をそのまま実行すれば wget は実行されないが 大量同時アクセスで DDBJ に迷惑をかけないように ~/Desktop/backup 上にあるファイルのコピーにしてください 3 Aug 03 06, NGS ハンズオン講習会 36

137 W3-:FASTQ ダウンロード 4 つまり ~/Documents/DRR0543 で wgetの部分はcpコマンドのほうのコピペにしてください ということ 4ディレクトリ作成は絶対に行い 指定された場所で作業を行うこと! 3 Aug 03 06, NGS ハンズオン講習会 37

138 W3-:FastQC FastQC ver は 第 4 回 W9- でインストールし fastqc というコマンドでパスを通している FastQC を実行し 共有フォルダ (~/Desktop/mac_share) に保存している 3 ファイルの確認 ヒトによっては他のものも見えるだろうが 最低限 4 が見えていれば問題ない 3 4 Aug 03 06, NGS ハンズオン講習会 38

139 W3-3: 結果を眺める 共有フォルダに保存することで 使いなれたホスト OS( この場合 Windows) 上で FastQC 実行結果ファイルを眺めることができる Aug 03 06, NGS ハンズオン講習会 39

140 W3-3: 結果を眺める FastQC 実行結果の解説は第 4 回 W8 と W7 および第 6 回 W4 にもあり 入力ファイル リード数は 95 3 配列長は bp の範囲であることがわかる 3 Aug 03 06, NGS ハンズオン講習会 40

141 W3-3: 結果を眺める Per base sequence quality この図の縦軸はクオリティスコア 赤線のスコア 0 を超えているかどうかが つの目安 Illumina HiSeq000 ( 第 4 回 W8) や MiSeq ( 第 6 回 W4) とは傾向が異なる 第 7 回 W6 では詳述しているが PacBio の公称スペックよりも全体的なクオリティが高い この理由は 現在眺めている 95 リードはクオリティフィルタリングで生き残ったリードに限定しているため Aug 03 06, NGS ハンズオン講習会 4

142 W3-3: 結果を眺める 横軸のリードポジションが 4000 bp あたりで黄色の縦棒がなくなっている この理由は 4000 bp 以上のリードが少数だからだと思われる それは の配列長分布 (Sequence Length Distribution) でも確認できる Aug 03 06, NGS ハンズオン講習会 4

143 W3-4: 配列長分布 配列長分布 (Sequence Length Distribution) このあたりで黄色の縦棒がなくなっているので おそらく 0 リードが黄色の縦棒の有無の閾値なのだろう Aug 03 06, NGS ハンズオン講習会 43

144 Contents W: ゲノムサイズ推定 (KmerGenie のインストールと利用 ) インストール single-end で実行 ( 結果の解説 ) paired-end で実行 ( 結果の解説 ) W: 配列長によるフィルタリング (Python プログラムの利用 ) DDBJ Pipeline(W3 から W7 まではほぼ省略 ) W8:k=3 での Velvet 実行結果の解析 W9:Platanus の実行 W0: 結果の解析 W:ACGT カウント ( 塩基ごとの出現頻度解析 ) ロングリード (PacBio) データと公共 DB W: PacBio 生データは bax.h5 形式 公共 DB は sra 形式と FASTQ 形式 W3: 公共 DB (DRA) の FASTQ ファイルを入力として FastQC W4:NCBI SRA (SRA) が提供する SRA Toolkit のインストール W5: 利用 (.sra.fastq への変換 ) W6:FastQC DDBJ Pipeline で HGAP を実行 W7:DDBJ Pipeline に解析したい生データファイル (.bax.h5) をアップロード W8:DDBJ Pipeline で HGAP を実行 W9:HGAP アセンブリ結果を眺める Aug 03 06, NGS ハンズオン講習会 44

145 NGS 連載第 7 回の 場所はこのあたり 次のスライドは 3 のリンク先ですが 基本的にコマンド打ち込み部分以外は リンク先にアクセスせずに スライドを眺めるだけでよい ネットワーク経由の作業なのでうまくいかない可能性もあるが その場合はスライドを眺めながらのエアーハンズオンに切り替え インストールに失敗しても その後に影響しないような構成にしてます 3 Aug 03 06, NGS ハンズオン講習会 45

146 W4-:SRA Toolkit 06 年 3 月 日現在の最新版は ver..5.7 プログラムは OS の種類ごとに用意されている Bio-Linux の実体は Ubuntu なので 3 ここ 話についてこれないヒトは 連載第 回を復習 4 のリンク先が次のスライド 4 3 Aug 03 06, NGS ハンズオン講習会 46

147 W4-:SRA Toolkit 目的の fastq-dump プログラムは よく使われるツール群 (Frequently Used Tools) の最初に位置する fastq-dump を利用したいがために SRA Toolkit をインストールするヒトがほとんどであろう 基本的には 3 を参考にインストールする クリック 3 Aug 03 06, NGS ハンズオン講習会 47

148 W4-:SRA Toolkit wget で tar.gz をダウンロードし 解凍するやり方が示されているが 折角なので第 4 回 W4-5, W3-5, W4, W5 で紹介した sudo apt-get install ソフトウェア名 で SRA Toolkit のインストールを行うやり方を伝授 Aug 03 06, NGS ハンズオン講習会 48

149 W4-:apt-cache 目的 : sudo apt-get install ソフトウェア名 で指定する名前を知りたい! やり方 : apt-cache n search キーワード で任意のキーワードを含むソフトウェア名をリストアップする ここでは SRA を含むソフトウェア名をリストアップ 欲しいソフトウェア名は sra-toolkit だということがわかる Aug 03 06, NGS ハンズオン講習会 49

150 W4-:apt-cache おまけ apt-cache n search SRA 実行結果として 小文字の sra を含むソフトウェア名もリストアップされたことから キーワード部分は 大文字でも小文字でもどちらでもいいのだろうと学習する また wc を追加することで ソフトウェア名が 6 個だったと認識 Aug 03 06, NGS ハンズオン講習会 50

151 W4-3:apt-get sudo apt-get install sra-toolkit root のパスワードを聞かれたら打ち込む ( 推奨手順通りだと pass409) Aug 03 06, NGS ハンズオン講習会 5

152 W4-3:apt-get Do you want to continue? と聞かれるので y イチイチ聞かれたくない場合は sudo apt-get y install sra-toolkit と -y オプションをつけておけばよい ( 第 4 回 W5-) Aug 03 06, NGS ハンズオン講習会 5

153 W4-3:apt-get インストール完了後の状態 06 年 06 月 9 日に気づいたこととして 赤下線部分とこの後のスライドを眺めればわかるが バージョンがかなり古い 06 年 7 月公開済みの原稿 PDF では apt-get の手順を推奨しているが 最新版 をインストールしたい場合は wget URL などプログラムのバージョンを自分で把握なり取り扱える手段でやりましょう SRA Toolkit はただの apt-get の練習台程度の位置づけでしかなかったが これもやったおかげで上記事実を知り こういうこともあるから気をつけねば という経験を積むことができた いろいろやるのは大事 Aug 03 06, NGS ハンズオン講習会 53

154 W4-4: 確認 インストール作業は ~/Documents/DRR0543 で行ったが sudo apt-get install ソフトウェア名 でやる場合は 基本的にどの作業ディレクトリ上でもよい インストール後は fastq-dump を使用可能 3 パスも既に通されている 3 Aug 03 06, NGS ハンズオン講習会 54

155 W4-5: パスが通っている sudo apt-get install sra-toolkit で無事インストール完了したあとは W4-4 で示したようにパスを通し終わった状態 それゆえ SRA Toolkit Installation and Configuration Guide 中の パスに関する注意書きは 気にしなくてもよい Aug 03 06, NGS ハンズオン講習会 55

156 Contents W: ゲノムサイズ推定 (KmerGenie のインストールと利用 ) インストール single-end で実行 ( 結果の解説 ) paired-end で実行 ( 結果の解説 ) W: 配列長によるフィルタリング (Python プログラムの利用 ) DDBJ Pipeline(W3 から W7 まではほぼ省略 ) W8:k=3 での Velvet 実行結果の解析 W9:Platanus の実行 W0: 結果の解析 W:ACGT カウント ( 塩基ごとの出現頻度解析 ) ロングリード (PacBio) データと公共 DB W: PacBio 生データは bax.h5 形式 公共 DB は sra 形式と FASTQ 形式 W3: 公共 DB (DRA) の FASTQ ファイルを入力として FastQC W4:NCBI SRA (SRA) が提供する SRA Toolkit のインストール W5: 利用 (.sra.fastq への変換 ) W6:FastQC DDBJ Pipeline で HGAP を実行 W7:DDBJ Pipeline に解析したい生データファイル (.bax.h5) をアップロード W8:DDBJ Pipeline で HGAP を実行 W9:HGAP アセンブリ結果を眺める Aug 03 06, NGS ハンズオン講習会 56

157 W5-6: デフォルトで実行 ここはダイジェスト版のみで眺めるだけ ( 詳細は第 7 回ウェブ資料 PDF) sra ファイルを入力として fastq-dump を実行すると fastq ファイルが得られる リード数情報は 3 を眺めるのでもわかる 4fastq ファイルの行数 (653,50) から辿れるリード数 ( 653,50/4 = 63,380 個 ) と同じ DRA 上で見られる数値 (63,48 リード ;W-5) と若干異なる理由は不明 3 4 Aug 03 06, NGS ハンズオン講習会 57

158 W6-3: 結果を眺める FastQC 実行結果 入力は さきほどの gzip 圧縮ファイル リード数は 63,380 3 配列長は bp の範囲であることがわかる 3 Aug 03 06, NGS ハンズオン講習会 58

159 W6-3: 結果を眺める Per base sequence quality この図の縦軸はクオリティスコア 赤線のスコア 0 を基準としてみると DRA から直接ダウンロードした 95 リードからなる FASTQ ファイルの FastQC 実行結果 (W3-3) と比べて 明らかに低くなっていることがわかる クオリティスコアは 平均で 0 弱程度 Aug 03 06, NGS ハンズオン講習会 59

160 W6-4:R で計算 塩基配列決定精度 ( エラー率 ) からクオリティスコアを R で計算する エラー率 3% のときは スコア 8.86 となるので 実際のスペック通りで安心 Rhoads and Au, Genomics Proteomics Bioinformatics, 3: 78-89, 05 Aug 03 06, NGS ハンズオン講習会 60

161 W6 は のあたりです FastQC 実行結果の html ファイルはここにあるので じっくり眺めたいヒトはどうぞ ~/Desktop/backup にもあります Aug 03 06, NGS ハンズオン講習会 6

162 Contents W: ゲノムサイズ推定 (KmerGenie のインストールと利用 ) インストール single-end で実行 ( 結果の解説 ) paired-end で実行 ( 結果の解説 ) W: 配列長によるフィルタリング (Python プログラムの利用 ) DDBJ Pipeline(W3 から W7 まではほぼ省略 ) W8:k=3 での Velvet 実行結果の解析 W9:Platanus の実行 W0: 結果の解析 W:ACGT カウント ( 塩基ごとの出現頻度解析 ) ロングリード (PacBio) データと公共 DB W: PacBio 生データは bax.h5 形式 公共 DB は sra 形式と FASTQ 形式 W3: 公共 DB (DRA) の FASTQ ファイルを入力として FastQC W4:NCBI SRA (SRA) が提供する SRA Toolkit のインストール W5: 利用 (.sra.fastq への変換 ) W6:FastQC DDBJ Pipeline で HGAP を実行 W7:DDBJ Pipeline に解析したい生データファイル (.bax.h5) をアップロード W8:DDBJ Pipeline で HGAP を実行 W9:HGAP アセンブリ結果を眺める Aug 03 06, NGS ハンズオン講習会 6

163 W-7:bax.h5 ファイル おさらい DDBJ Pipeline 上では PacBio 用の de novo アセンブリ用プログラム HGAP を利用可能 しかし HGAP は bax.h5 ファイルのみしか受け付けない Aug 03 06, NGS ハンズオン講習会 63

164 W7-:bax.h5 ファイル準備 DDBJ Pipeline にアップロードしたい bax.h5 ファイルを手元の PC にダウンロードしておく ここでは Desktop 上の share フォルダにダウンロード ファイルサイズは約.4GB に達するため それなりに時間はかかるだろう もちろん講習会ではやりません! 見るだけ! Aug 03 06, NGS ハンズオン講習会 64

165 W7-: アップロード WinSCP を用いて DDBJ Pipeline の実体である遺伝研スパコンにログインし デスクトップ share query フォルダに移動した状態 Aug 03 06, NGS ハンズオン講習会 65

166 W7-: アップロード アップロードしたい 3 つのファイルを query フォルダにドラッグ & ドロップ アップロードもそれなりに時間がかかる Aug 03 06, NGS ハンズオン講習会 66

167 W8-:HGAP 実行 de novo Assembly 3HGAP 4NEXT 3 4 Aug 03 06, NGS ハンズオン講習会 67

168 W8-:HGAP 実行 解析したい 3 つのファイルにチェックを入れて confirm Aug 03 06, NGS ハンズオン講習会 68

169 W8-: つのパラメータ と つの指定可能なパラメータがあります Aug 03 06, NGS ハンズオン講習会 69

170 W8-: パラメータの解説 つのパラメータに関する解説 Aug 03 06, NGS ハンズオン講習会 70

171 W8-3:HGAP 実行 ゲノムサイズは 乳酸菌の平均的なゲノムサイズである.5MB minimum length は Automatic Estimation( デフォルト ) を指定して 3NEXT PDF 原稿 p04 左下 3 Aug 03 06, NGS ハンズオン講習会 7

172 W8-3:HGAP 実行 RUN OK Aug 03 06, NGS ハンズオン講習会 7

173 W9-: 結果を眺める de novo Assembly Job ID 番号 (965) を頼りにすれば このページに辿り着ける 3 赤枠部分を見ると 4 コンティグ数は 4 つ このデータの正解は 3 つ ( chromosome and plasmids) 世間一般の評価通りのロングリード (PacBio) の長所がよく表れた結果 4 3 Aug 03 06, NGS ハンズオン講習会 73

174 W9-: 結果を眺める Total contig size は 乳酸菌の一般的なゲノムサイズと近く妥当 Maximum contig size のものが全体の 9 割以上を占めていることから これが乳酸菌の染色体ゲノムなのだろうと妄想する Aug 03 06, NGS ハンズオン講習会 74

175 W9-3: ダウンロード result.zip という zip 圧縮ファイルをダウンロードして その後の解析へと進んでいく Aug 03 06, NGS ハンズオン講習会 75

176 詳細は第 7 回原稿 PDF とウェブ資料をご覧ください この後の展開は W0:multi-FASTA ファイルの分割 プログラムによっては single-fasta のほうが取扱いやすい場合もある W:FASTQ ファイルの分割とクオリティスコア分布 FASTA ファイル用がうまく動かなくても 4 行で リードだということが分かっていれば Linux コマンドでどうにかなる という話 PacBio は FASTQ ファイルも出力する R で (single-)fastq ファイルを読み込んでクオリティスコア分布を描画し PacBio データは両末端のクオリティが低い傾向にあることを確認 環状化 ( ゲノム解読の finishing 作業の一部 ) アセンブリ結果として 最初と最後の末端部分が同じ配列の場合は 通常そのコンティグは環状と判断 それを確認するための基本的な考え方 手段 および環状化のノウハウを伝授 W:seqinr パッケージを用いて 仮想環状コンティグのドットプロットで感覚をつかむ W3: 重複配列の除去の感覚をつかむ ( これが環状化作業の実体 ) W4:dotter プログラムで実際の PacBio 出力結果に対して適用し 環状状態の概要を知る W5:Bio-Linux 上の blastn で 環状の同一コンティグ同士を DB 側と query 側にして実行 W6: 正確なアラインメントを眺め 切断箇所をクオリティスコア分布と合わせて判断する W7: 両端切断後のコンティグに対し クオリティスコア分布や dotter を再度実行して確認 W8:NCBI blast のやり方も示し 様々な解析手段を伝授 Aug 03 06, NGS ハンズオン講習会 76