ゲノム工学実習 - 荒川和晴, 慶應義塾大学先端生命科学研究所 (Kazuharu Arakawa, Institute for Advanced Biosciences, Keio University)

ゲノム工学実習大学院科目

【開講日程】 2018年度秋学期特定期間集中【担当教員】荒川　和晴

【前提科目（推奨）】B6161:基礎分子生物学４ 
【前提科目（推奨）】B6160:基礎分子生物学３
【前提科目（推奨）】B6159:基礎分子生物学２
【前提科目（推奨）】B6158:基礎分子生物学１
【前提科目（推奨）】34190:基礎分子生物学４
【前提科目（推奨）】34180:基礎分子生物学３
【前提科目（推奨）】34170:基礎分子生物学２
【前提科目（推奨）】34160:基礎分子生物学１
【前提科目（推奨）】B3215:生命科学実験の基礎
【前提科目（推奨）】C2038:遺伝子解析実習
【前提科目（推奨）】34130:遺伝子解析実習

【開講場所】 TTCK 【授業形態】講義、実習
【履修条件】TTCK生のみ履修可
【連絡先】 gaou@sfc.keio.ac.jp

注意: 配布資料などはSFC-SFSの授業ページで公開します。

DNA解析技術の飛躍的な向上により、微生物程度のゲノム解析はもはや「誰でも」「どこでも」可能なレベルにまで簡単になってきている。特に、携帯型ナノポアシークエンサーの登場は初期投資をほぼ必要とせずに長鎖DNAの解析を安価に可能とした。このような現状を踏まえれば、微生物程度のゲノムであれば遺伝子単位ではなくもはやゲノム単位でDNAを解析することが第一選択肢となる時代が到来していることを意味する。そこで、本実習では任意の微生物から長鎖DNAを抽出・精製し、ナノポアシークエンサーにてDNAを読み取り、それをバイオインフォマティクスによりアセンブル・アノテーションし、解析可能なゲノム情報にして、さらにそれをGenome Reportsの形にして国際誌に投稿するまでの全過程を学ぶ。

授業シラバス

主題と目標／授業の手法など

前半の過程では実際にナノポアシークエンサーに適した長鎖DNAをシーケンスする実験を実習として行い、後半ではシーケンスされたDNAをコンピュータを用いて解析する。よって、実験・バイオインフォマティクス双方の過程を学ぶが、知識としては少なくとも実験の経験があれば構わない。

教材・参考文献

参考文献：

荒川和晴(企画), "どこでも　誰でも　より長く　ナノポアシークエンサーが研究の常識を変える！", 実験医学 2018年1月号 Vol.36 No.1

提出課題・試験・成績評価の方法など

実験ノート及び最終レポートをもって評価する

履修上の注意

実験経験のあるTTCK生のみ履修可。

授業計画

第1回イントロダクション

ナノポアシーケンスと、ゲノム解析の流れについて講義します。

第2回長鎖DNA抽出 1

ナノポアシーケンス用長鎖DNAを抽出します。

第3回長鎖DNA抽出 2

ナノポアシーケンス用長鎖DNAを抽出します。

第4回長鎖DNA抽出 3

ナノポアシーケンス用長鎖DNAを抽出します。

第5回長鎖DNA抽出 4

ナノポアシーケンス用長鎖DNAを抽出します。

第6回長鎖DNA QC 1

長鎖DNAの品質をパルスフィールド電気泳動を用いて検証します。

第7回長鎖DNA QC 2

長鎖DNAの品質をパルスフィールド電気泳動を用いて検証します。

第8回長鎖DNA QC 3

長鎖DNAの品質をパルスフィールド電気泳動を用いて検証します。

第9回ナノポアライブラリ作製

ナノポアシーケンス用ライブラリを作成します。

第10回ナノポアシーケンシング

ライブラリをシーケンスにかけます。

第11回ゲノムアセンブリー

得られたゲノムをアセンブルします。

第12回エラー補正

Nanopolishを用いてエラー補正します。

第13回ゲノムアノテーション

D-FASTを用いてゲノムをアノテーションします。

第14回 Genome Report執筆

これまでに得られたデータをGenome Reportの形にまとめます。

その他

毎回実験ノートをまとめ、次回の準備をする

ゲノムのアセンブリー

環境設定

ソフトウェアのインストールにはcondaが便利です。https://www.anaconda.com/products/distribution

$ wget https://repo.anaconda.com/archive/Anaconda3-2022.10-Linux-x86_64.sh

$ bash Anaconda3-2022.10-Linux-x86_64.sh 

$ printenv SHELL を実行した結果が/bin/bashではない場合、$ bashを起動。

$ conda activate base 

$ conda install -c bioconda bbmap

bbmapのところをcanu, nanoplot, などとソフトウェアの名前に変更してインストール。

サーバはkingを使えると良いですが、Torqueの使い方がわからない場合にはsmith6/7を使っても良いです（授業以外での利用はご遠慮ください）。

ファイル

冨田研ファイルサーバの

/home/gaou/gew/

の該当する年度のフォルダの下にそれぞれのバーコードに相当するナノポア配列(fastq)、illuminaフォルダ以下にバーコードに対応するIllumina配列ファイル(fastq)があります。guppyソフトウェアでベースコール後、バーコードのdemultiplexを行なっています。各自自分のホームにファイルをコピーして以降の解析を実施してください。

リードのフィルタリング

まず、現状ではリードが多すぎるので、だいたいx50~x100になるようにリードを調整します。この時、長いリードは残したいので、

awk 'BEGIN {OFS = "\n"} {header = $0 ; getline seq ; getline qheader ; getline qseq ; if (length(seq) >= 10000 ) {print header, seq, qheader, qseq}}' < input.fq > filtered-10000.fastq

のように（ここのinput.fqを入力ファイル, あとは10000を任意の数字に）すると、任意の長さ以上の配列だけを取得できます。あるいは、BBMapのreformat.shを使って、

/home/gaou/kumamushi/software/bbmap/reformat.sh in=BC01.fq out=BC01-filter10k.fq minlength=10000 qin=33

のようにします。

以下のコマンド(BBMap: https://jgi.doe.gov/data-and-tools/bbtools/bb-tools-user-guide/bbmap-guide/) で統計値を算出して、カバレッジを合わせます。

/home/gaou/kumamushi/software/bbmap/stats.sh filtered-10000.fastq

基本的な配列の統計を可視化するには、NanoPlot (https://github.com/wdecoster/NanoPlot) が便利です。

[gaou@smith7 NanoPlot]$ NanoPlot -t 24 --fastq BC01.fq -p BC01

Canuでのアセンブリー

次にcanu ( https://github.com/marbl/canu/releases )をインストールします。Canuは最新版だと場合によってはまだバイナリが配布されていないので、自分でコンパイルが必要です。その場合、一つ前のバージョンだとLinux版のバイナリが配布されています。指示にしたがってコンパイルするか、バイナリの場合は解凍してください。コンパイルの場合、Javaのバージョンに依存性があるので、あまり古いOSだとコンパイルできない可能性がある点ご注意ください。kingではワーニングたくさんでますが問題なくコンパイルできました。

その後canuでアセンブリーを行います。メモリやCPUの関係上kingでやるのが良いかと思いますが、（king利用についてはこちらを参照。https://www.bioinfo.ttck.keio.ac.jp/wordpress/?page_id=2383)

qsub -I -l nodes=1:ppn=32

でkingのノードにログインし、

 /home/gaou/kumamushi/software/canu-2.2/bin/canu -nanopore BC01.fastq -d BC01 -p BC01 -fast useGrid=false genomeSize=4m maxThreads=8

のように打って(BC01.fastqは上でフィルタリングしたfastqのファイル名に、BC01のところは自分のバーコード、あるいは任意の名前に変えてくださいね）アセンブリーを実行してください。-fastオプションは精度を犠牲にして実行速度を上げるオプションで、今回くらいのカバレッジ（x100)があればつけても問題なくアセンブリーできると思いますが、時間がかかっても構わないなら外してください。genomeSizeオプション（頭にハイフンをつけない点に注意）は予想ゲノムサイズより少し大きめを設定してください。バクテリアの場合大抵4mで良いと思います。maxThreadsは、使用するサーバに合わせて設定してください。kingの場合、32に設定してください。

大きな問題がなければこれで数時間でアセンブリーが終了します。

/home/gaou/kumamushi/software/bbmap/stats.sh BC01/BC01.contigs.fasta

と打って、ちゃんとアセンブリーが終了したか確認します。

末端処理

まず環状化ができているかを確認します。

grep ">" BC01.contigs.fasta

で各contigのFASTAヘッダを見て、右端のsuggestCircularがyesになっているか確認します。長さが10kbpに満たないsuggestCircularがnoのcontigはゴミの可能性が高いですので、多くの場合破棄して構いません。

suggestCircularがyesのものは、染色体かプラスミドの可能性が高いです。長さが1Mbpに満たないものはプラスミドの可能性が高いので、適当に数kbp分をコピーしてNCBI BLASTに投げてみましょう。プラスミド配列にヒットするようでしたら高い確率でプラスミドと言えます。

ここでsuggestCircularがyesになっていても、環状であることを確認しているだけで環状化されているわけではありません。最初の50文字程度で自身の配列内を検索して、末端部分に該当する場所を見つけてください。基本的に、その後に続く配列が先と部分と末端の該当部分で一致するはずです。一致を確認したら、末端の一致部分を削除します。

アセンブリークオリティの検証

アセンブリーのクオリティ確認はN50だけではだめで、ちゃんとゲノムとして全遺伝情報がカバーできているか、を確認する必要があります。このためにはCEGMAやBUSCOといった手法(http://kazumaxneo.hatenablog.com/entry/2017/07/19/145640) を用います。これらのソフトウェアはインストールが面倒なのですが、理研が開発しているgVolanteというウェブサーバが非常に簡単に使えるようにしてくれています。https://gvolante.riken.jpBUSCOはv.1がバクテリアに対応しています。

バクテリア用にはCheckMというソフトウェアの方がより詳細にcompletenessを計算できます。CheckMはDFAST Quality Controlツール ( https://dfast.nig.ac.jp/dqc/submit/ )で簡単に検証できるので、Genome ReportにはBUSCOよりもCheckMの値を載せた方が良いです。

ローカルで CheckMを利用する際にはこちらを参考にしてください( https://kazumaxneo.hatenablog.com/entry/2017/09/22/012544 )。

ナノポアの配列のみでアセンブルした場合、大抵の場合BUSCOスコアが目安となる90%程度を著しく下回ります。これは、多くの場合indel系のエラーが解決されないことに起因します。

エラーコレクション

エラー補正には、ナノポアリードを使う場合と、Illuminaリードを使う場合の二通りがあります。ただし、Illuminaリードがを用いた方が確実に良い結果が得られます。今回この授業ではIlluminaシーケンスを用意しましたので、これを使ってエラーコレクションを行いましょう。新規で読んでいてIlluminaリードがない場合にはnanoporeオンリーのエラーコレクションを行います（非推奨）。

0. IlluminaデータのフィルタリングIlluminaのデータは正確ですが、とはいえ末端の配列などはクオリティが高くない場合があります。クオリティが悪い配列でエラーコレクションしようとすると余計にエラーを導入することになりかねないので、事前にクオリティでフィルタリングをした方が良いです。この時にはfastpというソフトウェアを使うとほぼ全自動でやってくれます。(https://kazumaxneo.hatenablog.com/entry/2018/05/21/111947)

 ~/kumamushi/software/fastp -i BC01_S1_L001_R1_001.fastq -I BC01_S1_L001_R2_001.fastq -o fastpR1.fq -O fastpR2.fq -g -c -w 8

配列が多すぎる場合にはだいたい100xくらいの量にサブサンプリングした方がいい場合があります。300bp pairedの配列で、ゲノムサイズが4Mの場合、100xだと600000ペアの配列があれば十分です。

/home/gaou/kumamushi/software/bbmap/reformat.sh in1=fastpR1.fq in2=fastpR2.fq out1=subsampledR1.fq out2=subsampledR2.fq samplereadstarget=600000 qin=auto qout=auto

1. Illuminaでの補正Pilonを使います。https://github.com/broadinstitute/pilon/wiki基本的にはBWAでmapping後pilonをかけるだけです。

まずはmapping。

/home/gaou/kumamushi/software/bwa-0.7.11/bwa index BC01.contigs.fasta
/home/gaou/kumamushi/software/bwa-0.7.11/bwa mem -t 8 BC01.contigs.fasta BC01_S1_merged_R1.fq  | /home/gaou/miniconda3/envs/hypo/bin/samtools view -@ 4 -b -o aln.bam - 
/home/gaou/miniconda3/envs/hypo/bin/samtools sort -T sort.tmp -o aln.sorted.bam -@ 4 aln.bam
/home/gaou/miniconda3/envs/hypo/bin/samtools index aln.sorted.bam

次にpilon。

java -Xms8g -jar /home/gaou/kumamushi/software/pilon-1.23.jar --genome BC01.contigs.fasta --bam aln.sorted.bam --threads 4 --output pilon1

pilon後に再度gVolanteでBUSCOスコアを算出すると、大幅に向上していることが確認できます。ただし、一回では不十分なことが多いので、ここで作成したエラーコレクション後のファイル(pilon1.fasta)に再度Illuminaリードをマッピングし、pilonをかけ直し、BUSCOスコアが向上する限りこれを繰り返します。

IlluminaデータをSRAからダウンロードする場合、SRA形式からFASTQ形式に変換する必要があります。

/home/gaou/kumamushi/software/sratoolkit.2.8.2-1-centos_linux64/bin/fastq-dump SRR390728 --split-files

ペアドエンドの場合には--split-filesオプションをつけてください。

2. Nanoporeでの補正nanopolishを使います。https://github.com/jts/nanopolishnanopolish indexでfastqとfast5を対応付け、fastqをbwaでリファレンスにマッピングし、実行します。マニュアルだとparallelを用いた方法が書いてありますが、以下のように実行するとシングルで実行できます。が、非常に時間がかかるので、parallelをインストールできる場合parallelを使って16並列 x 4スレッド、あるいは8並列 x 8スレッドくらいでやったほうがいいかもしれません。

nanopolish variants --consensus -r BC01.dedup.fastq -b BC01/reads.sorted.bam -g BC01/BC01.contigs.fasta -o BC01.nanopolished.fasta -t 64 -q dcm,dam -w tig00000001:1-3640229

マニュアルに書いていませんが、この時、-q dcm,damオプションをつけることは重要です。これはDNAのメチル化を考慮に入れたベースコールをするオプションで、これが入るとエラーコレクションの精度が大分向上します。（バクテリアゲノムは特にメチレーションが多いため）

kingだと例によってwarningは出ますが問題なくコンパイルできます。
fastqエントリがduplicateだと言われる場合： /home/gaou/bin/fastq-dedup.pl BC01.fastq > BC01.dedup.fastq
nanopolish index時は -fオプションでsequencing_summaryの場所を必ず指定すること（でないと異様に遅い）
nanopolishの最適カバレッジはx100~200です。これ以上ある場合にはx200以下までダウンサンプリングしてください。
nanopolish indexとbwa memはいずれもかなり時間がかかるので、同時にかけると良いでしょう。
nanopolish は異様に時間かかるので、先にRaconで様子を見てみるのもアリです。http://kazumaxneo.hatenablog.com/entry/2018/03/22/013006

Ishikawa S, Huang M, Tomita A, Kurihara Y, Watanabe R, Iwai H, Arakawa K*, "Complete Genome Sequences of Four Bacteria Isolated from the Gut of a Spiny Ant (Polyrhachis lamellidens)", Microbiol Resour Announc, 2022, 11:e00333-22. (Publisher)

Takeda T, Fukumitsu N, Yuzawa S, Arakawa K*, "Complete Genome Sequence of Streptomyces albus Strain G153", Microbiol Resour Announc, 2022, 11:e00332-22. (Publisher)

2020年

Takahashi H, Yang J, Yamamoto H, Fukuda S, Arakawa K*, "Complete Genome Sequence of Adlercreutzia equolifaciens subsp. celatus DSM18785", Microbiol Resour Announc, 2021, 10:e00354-21. (Publisher)

Warashina, T, Yamamura S, Suzuki H, Amachi S, Arakawa K, "Complete Genome Sequence of Geobacter sp. Strain SVR, an Antimonate-reducing Bacterium Isolated from Antimony-rich Mine Soil", Microbiol Resour Announc, 2021, 10:e00142-21. (Publisher)

2019年

Takeyama N, Huang M, Sato K, Galipon J, Arakawa K*, "Complete Genome Sequence of Halomonas hydrothermalis Strain Slthf2, a Halophilic Bacterium Isolated from a Deep-Sea Hydrothermal-Vent Environment", Microbiol Resour Announc, 2020, 9:e00294-20. (Publisher)

Takahashi Y, Takahashi H, Galipon J, Arakawa K*, "Complete Genome Sequence of Halomonas meridiana Strain Slthf1, Isolated from a Deep-Sea Thermal Vent", Microbiol Resour Announc, 2020, 9:e00292-20. (Publisher)

Seo K, Tanaka K, Fukuda S, Arakawa K*, "Complete Genome Sequences of Two Cutibacterium acnes Strains Isolated from an Orthopedic Surgical Site", Microbiol Resour Announc, 2020, 9:e00290-20. (Publisher)

Kurihara Y, Kawai S, Sakai A, Galipon J, Arakawa K*, "Complete Genome Sequence of Halomonas meridiana Strain Eplume2, Isolated from a Hydrothermal Plume in the Northeast Pacific Ocean", Microbiol Resour Announc, 2020, 9:e00330-20. (Publisher)

Inoue H, Shibata S, Ii K, Inoue J, Fukuda S, Arakawa K, "Complete Genome Sequence of Bifidobacterium longum Strain Jih1, Isolated from Human Feces", Microbiol Resour Announc, 2020, 9:e00319-20. (Publisher)

Nishimura K, Ikarashi M, Yasuda Y, Sato M, Cano Guerrero M, Galipon J, Arakawa K, "Complete Genome Sequence of Sphingomonas paucimobilis Strain Kira, Isolated from Human Neuroblastoma SH-SY5Y Cell Cultures Supplemented with Retinoic Acid.", Microbiol Resour Announc, 2021, 10(6):e01156-20. (PubMed)

2018年

Tsurumaki M, Deno S, Galipon J, Arakawa K*, "Complete Genome Sequence of Halophilic Deep-Sea Bacterium Halomonas axialensis Strain Althf1", Microbiol Resour Announc, 2019, 8:e00839-19. https://mra.asm.org/content/8/31/e00839-19

Evans-Yamamoto D, Takeuchi N, Masuda T, Murai Y, Onuma Y, Mori H, Masuyama N, Ishiguro S, Yachie N, Arakawa K*, "Complete genome sequence of Psychrobacter sp. strain KH172YL61, isolated from deep-sea sediments in the Nankai Trough, Japan", Microbiol Resour Announc, 2019, 8:e00326-19. https://mra.asm.org/content/8/16/e00326-19

Nagata S, Ii KM, Tsukimi T, Miura MC, Galipon J, Arakawa K*, "Complete genome sequence of Halomonas olivaria, a moderately halophilic bacterium isolated from olive processing effluents, obtained by nanopore sequencing", Microbiol Resour Announc, 2019, 8:e00144-19. https://mra.asm.org/content/8/18/e00144-19

Saito M, Nishigata A, Galipon J, Arakawa K*, "Complete Genome Sequence of Halomonas sulfidaeris Strain Esulfide1 Isolated from a Metal Sulfide Rock at a Depth of 2,200 Meters, Obtained Using Nanopore Sequencing", Microbiol Resour Announc, 2019, 8(23):e00327-19. https://mra.asm.org/content/8/23/e00327-19

Murai Y, Masuda T, Onuma Y, Evans-Yamamoto D, Takeuchi N, Mori H, Masuyama N, Ishiguro S, Yachie N, Arakawa K*, "Complete Genome Sequence of Bacillus sp. Strain KH172YL63, Isolated from Deep-Sea Sediment", Microbiol Resour Announc, 2020, 9:e00291-20. (Publisher)

ゲノム工学実習 大学院科目

メニュー

研究関連リンク

授業

ソフトウェア

データベース

所属学会

リンク

連絡先

ゲノム工学実習大学院科目