macでインフォマティクス

macでインフォマティクス

NGS関連のインフォマティクス情報についてまとめています。

植物のRNA seqデータからvirus配列を検出する Kodoja

 

 ウイルス感染は、食物と燃料のために栽培される作物で特に重要な問題である。ウイルスは収量と品質の大きな損失を引き起こし、その結果、ウイルスは重要な経済的悪影響を及ぼす[ref.1]。英国では、ポテトウイルスYは年間3,000〜4,000万ポンドのジャガイモの収穫損失[ref.2]を引き起こし、アジアではイネに感染するウイルス(ライスグラッシースタントウイルスなど)は年間1億2,000万ドルの作物損失を引き起こす[ref.3]。これらの例は、高速で正確なウイルス検出方法の必要性を強調している。ウイルス感染の症状には黄変や発育阻害が含まれる場合があるが、多くの場合、症状は存在しないか、他の要因によって隠されている。場合によっては、植物ウイルスが相乗的に相互作用して、新しいまたはより深刻な病気の症状を引き起こす[ref.4]。 1つの例は、ラズベリーの果物がもろくなる複雑な病気で、これは2つのウイルス、 Raspberry bushy dwarf virusおよびRaspberry潜伏ウイルス[ref.5]の存在によって引き起こされる。
 作物が新しい地理的位置で栽培され、農業慣行が強化されると、新しいウイルスが確立され、既存のウイルスが宿主範囲を広げるリスクが高まる。したがって、植物ウイルスの診断は、将来の食料安全保障の観点から重要性が高まっている分野である。
ウイルスを検出するための標準的な分子技術には、逆転写(RT)PCRに基づく方法が含まれる。ただし、このような手法では既知のウイルス検出のみ可能になる。つまり、各テストは1つのウイルスまたはごく少数の関連ウイルスに固有である[ref.6]。さらに、ウイルスゲノムが進化し、テストが時間の経過とともに無効になる可能性があり、病気の診断が制限される。このような制限は、最近、複数のウイルスを仮説なしで同時に検出するための次世代シーケンス(NGS)メソッドの使用により克服された[ref.7]。植物ウイルスの大部分はRNAを遺伝物質として持ち、DNAゲノムを持つウイルスはRNA転写産物を生成する。さらに、真核生物のsiRNAはRNA干渉を介して抗ウイルス免疫を誘導し、この過程でウイルス由来のsiRNAが宿主に濃縮される[ref.8]。したがって、RNAおよびsmall RNA(sRNA)シーケンスの両方が、植物のウイルス検出に効果的な方法である。ただし、これは2つの重要な要素に依存している:(a)堅牢なRNA抽出およびエンリッチメントプロトコル、および(b)ウイルス同定のための高速で堅牢なバイオインフォマティクスツール。
 RNA抽出およびエンリッチメントプロトコルの範囲、およびバイオインフォマティクスのワークフローは、以前にヒトの臨床サンプル用に開発された(レビューについては[ref.9]を参照)。最近、このような研究により、48時間以内にウイルス性疾患の診断と実用的な臨床管理が行われた[ref.10]。この臨床研究で使用されたワークフローは、ウイルス診断ツールに必要な2つの主要要素で構成されていた。(a)宿主ヌクレオチド配列の特定と除去、(b)ウイルス配列の特定。ただし、ヒトゲノムには十分なアノテーションが付けられており(ホスト配列の簡単な削除が可能)、またヒトウイルスデータは配列データベースでより一般的のため、臨床サンプルでのウイルス検出は植物よりも簡単である。それに比べて、多くの作物のゲノムは不完全であるかアノテーションが不十分であり、植物ウイルス配列はデータベースでは十分に表現されていない。

 最近、NGSデータからウイルス検出に利用可能なバイオインフォマティクスツールとワークフローをレビューした[ref.9]。大部分はヒトのNGSデータに最適化されており、ウイルスの識別にはシーケンスIDを超えるものはほとんどなく、インストールや使用に多くの重要な計算知識が必要であると結論付けた。 Taxonmer [ref.11]とVirusDetect [ref.12]の2つのツールはWebサーバーとして利用でき、植物[ref.2]からのRNAシーケンシングデータの分析の可能性を提供している。ただし、レビューでは、公開された3つのツールが植物データでテストされたが、植物のウイルス検出に焦点を当てたプロジェクトでは使用されていないという事実が強調された。代わりに、このアプローチは一般に分析中により高い柔軟性を提供するため、ワークフローの外部でスタンドアロンマッピングおよびアセンブリアルゴリズムを使用した。
 ウイルス識別ワークフローは、(a)低品質のリードおよびアダプターシーケンスのトリミングを含む、rawデータファイルの品質管理対策の実施、(b)ホストシーケンスの識別、および(c)ウイルスシーケンスの識別が可能である必要がある。既知のウイルスの識別は、既存のウイルス配列のデータベースにマッピングすることで実行できるが、新しい株や新しいウイルスの識別には、ワークフローを超えた専門知識と追加の分析が必要である。
 公開されているウイルス検出ワークフローの多くは、アセンブリマッピングアルゴリズムを使用してウイルスシーケンスを識別する。ただし、アセンブリマッピングの両方が非常に計算集約的である可能性があるため、ワークフローは
大規模なデータセットでは実行時間が長くなる可能性がある。また、アセンブリおよびマッピング方法では、未アセンブリのリードが識別されないままになる。 RNA-seqデータセットでウイルスリードを識別するもう1つの方法は、k-merプロファイリングを使用することである。これはTaxonomer [ref.11]で正常に実装されている。 RNAおよびDNAシーケンスは長さkの複数の部分文字列に分割できる。このようにして、シーケンスをk-merプロファイルで表すことができ、これらのプロファイルをtaxaに割り当てして比較できる。 k-merプロファイルは、バイオインフォマティクスの類似検索の範囲で使用されている。メタゲノミクスでは配列間のアラインメントのない類似性分析が可能である[ref.13]。分類プロファイリングでは、ビニング法はk-merプロファイルを使用して配列をクラスター化し、ドラフトゲノム回復を可能にする。このような方法は、リファレンスゲノムを使用せずに、集団内のウイルスハプロタイプの同定にもうまく適用されている。
 ここで紹介するKodojaワークフローは、NGSデータセットを使用する幅広い研究者に適用できる独自の機能セットを組み合わせている。目的は、アセンブリマッピングの方法を超えたワークフローを開発することであった。これは、植物データセット用に特に最適化され、非生物情報学者がアクセスできる。 Kodojaは、混合された(植物とウイルス、細菌、真菌の両方のヌクレオチドを含む)RNA-seqデータからウイルスシーケンスを識別することができるワークフローである。 Kodojaは、(a)植物NGSデータに固有であり、(b)ヌクレオチドレベルでk-merプロファイリングを使用し、既存のツールKraken [ref.16] ]およびKaiju [ref.17]を使い、(c)はBioconda [ref.18]を介してローカルインストールして使用でき、(d)Galaxy Webベースの分析環境[ref.19]内のツールとしても使用できる。

 

f:id:kazumaxneo:20191207200849p:plain

Kodoja workflow.  論文より転載。

 

インストール

ubuntu18.04のPython 3.7環境でテストした(docker使用、ホストOS macos10.14)。

依存

You can use Python 2.7 or Python 3, specifically Kodoja has been tested on Python 3.6.

  • FastQC v0.11.5
  • Trimmomatic v0.36
  • Kraken v1.0
  • Kaiju v1.5.0

Python packages:

  • numpy v1.9
  • biopython v1.67
  • pandas v0.14
  • ncbi-genome-download v0.2.6

本体 Github

#bioconda (link)
conda install -c bioconda kodoja

kodoja_search.py  -h

# kodoja_search.py  -h

usage: kodoja_search.py [-h] [--version] -o OUTPUT_DIR -d1 KRAKEN_DB -d2

                        KAIJU_DB -r1 READ1 [-r2 READ2] [-f DATA_FORMAT]

                        [-t THREADS] [-s HOST_SUBSET] [-m TRIM_MINLEN]

                        [-a TRIM_ADAPT] [-q KRAKEN_QUICK] [-p]

                        [-c KAIJU_SCORE] [-l KAIJU_MINLEN] [-i KAIJU_MISMATCH]

 

Kodoja Search is a tool intended to identify viral sequences

in a FASTQ/FASTA sequencing run by matching them against both Kraken and

Kaiju databases.

 

optional arguments:

  -h, --help            show this help message and exit

  --version             show program's version number and exit

  -o OUTPUT_DIR, --output_dir OUTPUT_DIR

                        Output directory path, required

  -d1 KRAKEN_DB, --kraken_db KRAKEN_DB

                        Kraken database path, required

  -d2 KAIJU_DB, --kaiju_db KAIJU_DB

                        Kaiju database path, required

  -r1 READ1, --read1 READ1

                        Read 1 file path, required

  -r2 READ2, --read2 READ2

                        Read 2 file path

  -f DATA_FORMAT, --data_format DATA_FORMAT

                        Sequence data format (default fastq)

  -t THREADS, --threads THREADS

                        Number of threads (default 1)

  -s HOST_SUBSET, --host_subset HOST_SUBSET

                        Subset sequences with this tax id from results

  -m TRIM_MINLEN, --trim_minlen TRIM_MINLEN

                        Trimmomatic minimum length

  -a TRIM_ADAPT, --trim_adapt TRIM_ADAPT

                        Illumina adapter sequence file

  -q KRAKEN_QUICK, --kraken_quick KRAKEN_QUICK

                        Number of minium hits by Kraken

  -p, --kraken_preload  Kraken preload database

  -c KAIJU_SCORE, --kaiju_score KAIJU_SCORE

                        Kaju alignment score

  -l KAIJU_MINLEN, --kaiju_minlen KAIJU_MINLEN

                        Kaju minimum length

  -i KAIJU_MISMATCH, --kaiju_mismatch KAIJU_MISMATCH

                        Kaju allowed mismatches

 

The main output of ``kodoja_search.py`` is a file called ``virus_table.txt``

in the specified output directory. This is a plain text tab-separated table,

the columns are as follows:

 

1. Species name,

2. Species NCBI taxonomy identifier (TaxID),

3. Number of reads assigned by *either* Kraken or Kaiju to this species,

4. Number of Reads assigned by *both* Kraken and Kaiju to this species,

5. Genus name,

6. Number of reads assigned by *either* Kraken or Kaiju to this genus,

7. Number of reads assigned by *both* Kraken and Kaiju to this genus.

 

The output directory includes additional files, including ``kodoja_VRL.txt``

(a table listing the read identifiers used) which is intended mainly as

input to the ``kodoja_retrieve.py`` script.

 

See also https://github.com/abaizan/kodoja/wiki/Kodoja-Manual

kodoja_build.py -h

# kodoja_build.py -h

usage: kodoja_build.py [-h] [--version] -o OUTPUT_DIR [-t THREADS]

                       [-p HOST_TAXID] [-d DOWNLOAD_PARALLEL] [-n]

                       [-e [EXTRA_FILES [EXTRA_FILES ...]]]

                       [-x [EXTRA_TAXIDS [EXTRA_TAXIDS ...]]] [-v]

                       [-b KRAKEN_TAX] [-k KRAKEN_KMER] [-m KRAKEN_MINIMIZER]

                       [-a DB_TAG]

 

Kodoja Build creates a database for use with Kodoja Search.

 

optional arguments:

  -h, --help            show this help message and exit

  --version             show program's version number and exit

  -o OUTPUT_DIR, --output_dir OUTPUT_DIR

                        Output directory path, required

  -t THREADS, --threads THREADS

                        Number of threads, default=1

  -p HOST_TAXID, --host_taxid HOST_TAXID

                        Host tax ID

  -d DOWNLOAD_PARALLEL, --download_parallel DOWNLOAD_PARALLEL

                        Parallel genome download, default=4

  -n, --no_download     Genomes have already been downloaded

  -e [EXTRA_FILES [EXTRA_FILES ...]], --extra_files [EXTRA_FILES [EXTRA_FILES ...]]

                        List of extra files added to "extra" dir

  -x [EXTRA_TAXIDS [EXTRA_TAXIDS ...]], --extra_taxids [EXTRA_TAXIDS [EXTRA_TAXIDS ...]]

                        List of taxID of extra files

  -v, --all_viruses     Build databases with all viruses (not plant specific)

  -b KRAKEN_TAX, --kraken_tax KRAKEN_TAX

                        Path to taxonomy directory

  -k KRAKEN_KMER, --kraken_kmer KRAKEN_KMER

                        Kraken kmer size, default=31

  -m KRAKEN_MINIMIZER, --kraken_minimizer KRAKEN_MINIMIZER

                        Kraken minimizer size, default=15

  -a DB_TAG, --db_tag DB_TAG

                        Suffix for databases

 

See also https://github.com/abaizan/kodoja/wiki/Kodoja-Manual

kodoja_retrieve.py -h

# kodoja_retrieve.py -h

usage: kodoja_retrieve.py [-h] [--version] -o FILE_DIR -r1 READ1 [-r2 READ2]

                          [-f USER_FORMAT] [-t TAXID] [-g] [-s]

 

Kodoja Retrieve is used with the output of Kodoja Search to

give subsets of your input sequencing reads matching viruses.

 

optional arguments:

  -h, --help            show this help message and exit

  --version             show program's version number and exit

  -o FILE_DIR, --file_dir FILE_DIR

                        Path to directory of kodoja_search results, required

  -r1 READ1, --read1 READ1

                        Read 1 file path, required

  -r2 READ2, --read2 READ2

                        Read 2 file path, default: False

  -f USER_FORMAT, --user_format USER_FORMAT

                        Sequence data format, default: fastq

  -t TAXID, --taxID TAXID

                        Virus tax ID for subsetting, default: All viral

                        sequences

  -g, --genus           Include sequences classified at genus

  -s, --stringent       Only subset sequences identified by both tools

 

The main output of ``kodoja_search.py`` is a file called ``virus_table.txt``

(a table summarising the potential viruses found), but the specified output

directory will also contain ``kodoja_VRL.txt`` (a table listing the read

identifiers). This second file is used as input to ``kodoja_retrieve.py``

along with the original sequencing read files.

 

A sub-directory ``subreads/`` will be created in the output directory,

which will include either FASTA or FASTQ files named as follows:

 

* ``subset_files/virus_all_sequences1.fasta`` FASTA output

* ``subset_files/virus_all_sequences1.fastq`` FASTQ output

 

And, for paired end datasets,

 

* ``subset_files/virus_all_sequences2.fasta`` FASTA output

* ``subset_files/virus_all_sequences2.fastq`` FASTQ output

 

However, if the ``-t 12345`` option is used rather than ``virus_all_...``

the files will be named ``virus_12345_...`` instead.

 

See also https://github.com/abaizan/kodoja/wiki/Kodoja-Manual

 

 

データベースの準備

ここではオーサーらがビルドしたデータベース(krakenのDB)を使用する。 

mkdir kodojaDB_v1.0
cd kodojaDB_v1.0
wget https://zenodo.org/record/1406071/files/kodojaDB_v1.0.tar.gz
tar -xvf kodojaDB_v1.0.tar.gz

解凍してできたkaijuDB/とkrakenDB/を使う。

 

実行方法

kodoja_search.py --kraken_db kodojaDB_v1.0/krakenDB --kaiju_db kodojaDB_v1.0/kaijuDB \
-r1 pair_1.fq -r2 pair_2.fq -t 8 -o outdir
  • -r1    path to the single-end or first paired-end file (required)
  • -r2    path to second paired-end file (default=False)
  • -f      specify the file-type for file1 ("fasta" or "fastq" - default='fastq')
  • -o     path to the results folder (required)
  • -t      number of threads on cluster (default=1)
  • --kraken_db   path to kraken database (required)
  • --kaiju_db       path to kaiju database, nodes.dmp and names.dmp files (required)

 ランが終わるとvirus_table.txtに結果がまとめられる。

 

引用

Kodoja: A workflow for virus detection in plants using k-mer analysis of RNA-sequencing data

Amanda Baizan-Edge, Peter Cock, Stuart MacFarlane, Wendy McGavin, Lesley Torrance, Susan Jones

J Gen Virol. 2019 Mar;100(3):533-542