macでインフォマティクス

macでインフォマティクス

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

完全に自動化された 16S･18Sメタアンプリコン解析パイプライン AmpWrap

2025 amplicon sequence Bioinformatics Advances ASV (amplicon sequence variant) Nanopore long read snakemake automated pipeline

 

 次世代シーケンシング革命は、群集構成を探索するための効率的かつ費用対効果の高い方法としてメタバーコーディングの確立を推進した。原核生物の16S rRNA遺伝子などの分類マーカー遺伝子のアンプリコンシーケンシングは、ハイスループット分類プロファイリングのための効率的な方法を提供する。ロングリード技術の登場により、16S rRNA遺伝子の一部領域ではなく全体のシーケンシングが可能になり、種レベルの解像度を達成する可能性が高まった。このような実験は費用対効果が高く拡張性も高いものの、統合されたユーザーフレンドリーな分析ワークフローの欠如が依然として大きなボトルネックとなっている。現在のパイプラインでは、複雑な依存関係を持つ複数のツールを使用する必要があり、パラメータの最適化は手動で行われることが多く、再現性と全体的な効率性が制限されている。これらの制限に対処するため、イルミナとナノポアの両方のアンプリコンを解析するために設計された自動化されたワンラインワークフローであるAmpWrapを開発した。このワークフローは、ユーザーの手間を最小限に抑え、ノイズを低減しながら最大限のリード数と情報を保持するようにトリミングパラメータを自動的に最適化する。

 

インストール

Github

git clone https://github.com/LDoni/AmpWrap.git
cd AmpWrap/ampwrap/
mamba env create -f ampwrap.yml
conda activate ampwrap
chmod +x setup.sh
./setup.sh
bash setup.sh

#docker
docker pull ghcr.io/ldoni/ampwrap:v1.1.1

>  ampwrap -h

Usage: ampwrap <workflow> [args]

    workflow: 'short' or 'long' (required)

    args:     arguments for the selected workflow

 

To see workflow-specific help:

    ampwrap short --help

    ampwrap long --help

 

>  ampwrap short --help

 

   █████╗ ███╗   ███╗██████╗                        ██╗

  ██╔══██╗████╗ ████║██╔══██╗                       ╚██╗

  ███████║██╔████╔██║██████╔╝█████╗█████╗█████╗█████╗╚██╗

  ██╔══██║██║╚██╔╝██║██╔═══╝ ╚════╝╚════╝╚════╝╚════╝██╔╝

  ██║  ██║██║ ╚═╝ ██║██║                            ██╔╝

  ╚═╝  ╚═╝╚═╝     ╚═╝╚═╝                            ╚═╝

    ██╗                 ██╗    ██╗██████╗  █████╗ ██████╗

   ██╔╝                 ██║    ██║██╔══██╗██╔══██╗██╔══██╗

  ██╔╝█████╗█████╗█████╗██║ █╗ ██║██████╔╝███████║██████╔╝

  ╚██╗╚════╝╚════╝╚════╝██║███╗██║██╔══██╗██╔══██║██╔═══╝

   ╚██╗                 ╚███╔███╔╝██║  ██║██║  ██║██║

    ╚═╝                  ╚══╝╚══╝ ╚═╝  ╚═╝╚═╝  ╚═╝╚═╝

 

 

Welcome to ampwrap_short

 

usage: AmpWrap_short [-h] [-i INPUT_DIRECTORY [INPUT_DIRECTORY ...]] [-a FORWARD_PRIMER] [-A REVERSE_PRIMER] [-o OUTPUT_DIRECTORY] [-l AMPLEN] [-d DB_CHOICE] [-q] [--db-info] [--cutadapt_trim_by_length]

                     [--cutadapt_disable_discard_untrimmed] [--dada2_params DADA2_PARAMS] [--dry-run-tune]

 

ampwrap_short: process raw sequencing files and assign taxonomy.

 

optional arguments:

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

  -i INPUT_DIRECTORY [INPUT_DIRECTORY ...], --input_directory INPUT_DIRECTORY [INPUT_DIRECTORY ...]

                        Input directory with raw fq/fq.gz/fastq/fastq.gz files, add multiple space-separated paths for multiple runs e.g. -i path/to/run1 path/to/run2 path/to/run3

  -a FORWARD_PRIMER, --forward_primer FORWARD_PRIMER

                        Forward primer sequence

  -A REVERSE_PRIMER, --reverse_primer REVERSE_PRIMER

                        Reverse primer sequence

  -o OUTPUT_DIRECTORY, --output_directory OUTPUT_DIRECTORY

                        Output directory (default: ./ampwrap_short)

  -l AMPLEN, --amplen AMPLEN

                        Length of the amplified sequence target (not including primers-->e.g., R926 - F515 = 411 - pF19 - pR20 = 372)

  -d DB_CHOICE, --db_choice DB_CHOICE

                        Database for taxonomy assignment. Choices: decipher_silva138 decipher_gtdb226 decipher_rdp18 dada2_silva_genus138 dada2_RDP_genus19 dada2_GG2_genus09 dada2_RefSeq_RDPv16 dada2_GTDB_r202

                        dada2_pr2_5.1.1 dada2_18S_silva132

  -q, --quiet           Minimal standard output

  --db-info             Show available databases and exit

  --cutadapt_trim_by_length

                        Cutadapt trim primers by length

  --cutadapt_disable_discard_untrimmed

                        Deactivate cutadapt option --discard-untrimmed

  --dada2_params DADA2_PARAMS

                        Semicolon-separated list of DADA2 parameters to override Figaro (e.g., 'truncLen=c(240,200);maxEE=c(2,2)'

  --dry-run-tune        Run only Figaro tuning and show first 5 results, then clean up

 

> ampwrap long --help

 

 █████╗ ███╗   ███╗██████╗                                                                  ██╗

██╔══██╗████╗ ████║██╔══██╗                                                                 ╚██╗

███████║██╔████╔██║██████╔╝█████╗█████╗█████╗█████╗█████╗█████╗█████╗█████╗█████╗█████╗█████╗╚██╗

██╔══██║██║╚██╔╝██║██╔═══╝ ╚════╝╚════╝╚════╝╚════╝╚════╝╚════╝╚════╝╚════╝╚════╝╚════╝╚════╝██╔╝

██║  ██║██║ ╚═╝ ██║██║                                                                      ██╔╝

╚═╝  ╚═╝╚═╝     ╚═╝╚═╝                                                                      ╚═╝

 

  ██╗                                                           ██╗    ██╗██████╗  █████╗ ██████╗

 ██╔╝                                                           ██║    ██║██╔══██╗██╔══██╗██╔══██╗

██╔╝█████╗█████╗█████╗█████╗█████╗█████╗█████╗█████╗█████╗█████╗██║ █╗ ██║██████╔╝███████║██████╔╝

╚██╗╚════╝╚════╝╚════╝╚════╝╚════╝╚════╝╚════╝╚════╝╚════╝╚════╝██║███╗██║██╔══██╗██╔══██║██╔═══╝

 ╚██╗                                                           ╚███╔███╔╝██║  ██║██║  ██║██║

  ╚═╝                                                            ╚══╝╚══╝ ╚═╝  ╚═╝╚═╝  ╚═╝╚═╝

 

 

 

 

 

Welcome to ampwrap_long

usage: AmpWrap_long [-h] [-i INPUT_DIRECTORY] [-o OUTPUT_DIRECTORY] [-tr {porechop,porechop_abi}] [-d {rdp,silva,emu}] [-q] [-t THREADS] [--nl-min-len NL_MIN_LEN] [--nl-max-len NL_MAX_LEN] [--nl-min-qual NL_MIN_QUAL]

                    [--emu-min-ab EMU_MIN_AB] [--cutadapt-forward CUTADAPT_FORWARD] [--cutadapt-reverse CUTADAPT_REVERSE] [--db-info] [-c] [-kc]

 

Process raw sequencing files for downstream analysis.

 

optional arguments:

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

  -i INPUT_DIRECTORY, --input-directory INPUT_DIRECTORY

                        Specify the input directory containing raw fq/fq.gz/fastq/fastq.gz files

  -o OUTPUT_DIRECTORY, --output-directory OUTPUT_DIRECTORY

                        Specify the output directory (default: ./ampwrap_long_results)

  -tr {porechop,porechop_abi}, --trimming {porechop,porechop_abi}

                        Specify the method to remove adapters, skip if unspecified

  -d {rdp,silva,emu}, --db-choice {rdp,silva,emu}

                        Specify the database choice for taxonomy assignment for EMU: silva, rdp or emu_db (default: emu)

  -q, --quiet           Minimal standard output

  -t THREADS, --threads THREADS

                        Number of threads

  --nl-min-len NL_MIN_LEN

                        Filter amplicons on a minimum read length. Default=1200

  --nl-max-len NL_MAX_LEN

                        Filter amplicons on a maximum read length. Default=1800

  --nl-min-qual NL_MIN_QUAL

                        Filter on a minimum average read quality score. Default=10

  --emu-min-ab EMU_MIN_AB

                        Filter species with relative abundance below threshold; .01 = 1 percent. Default=0.0001

  --cutadapt-forward CUTADAPT_FORWARD

                        Forward primer sequence to trim it with cutadapt, no trimming if unspecified

  --cutadapt-reverse CUTADAPT_REVERSE

                        Reverse primer sequence to trim it with cutadapt, no trimming if unspecified

  --db-info             Print available databases with version and checksum info

  -c, --chimera         Specify chimera flag to remove chimeras, skip if unspecified

  -kc, --keep-counts    Include estimated read counts in EMU output

usage: AmpWrap_long [-h] [-i INPUT_DIRECTORY] [-o OUTPUT_DIRECTORY] [-tr {porechop,porechop_abi}] [-d {rdp,silva,emu}] [-q] [-t THREADS] [--nl-min-len NL_MIN_LEN] [--nl-max-len NL_MAX_LEN] [--nl-min-qual NL_MIN_QUAL]

                    [--emu-min-ab EMU_MIN_AB] [--cutadapt-forward CUTADAPT_FORWARD] [--cutadapt-reverse CUTADAPT_REVERSE] [--db-info] [-c] [-kc]

 

 

テストラン

cd AmpWrap/test/

#16S V4V5
bash test_short.sh

#18S V4
bash test_real_short_data18S.sh

#Nanopore long read
bash test_long.sh

単純にfastqのダウンロード後､以下のコマンドが実行されている(16S V4V5)

ampwrap short -i short -a GTGYCAGCMGCCGCGGTAA -A CCGYCAATTYMTTTRAGTTT -l 372  -o short_test_output

 

テストの出力 #16S V4V5

results/

ASVs_taxonomy.tsv: 各ASVの系統分類結果（from Kingdom to Species）

ASVs.fa: 各ASVの塩基配列データ

ASVs_counts.tsv: サンプルごとのAmplicon Sequence Variantカウントデータ。(どのサンプルにどの菌が何リード存在するか)

phyloseq_object.rds:  phyloseq パッケージ読み込み用オブジェクトファイル

read-count-tracking.tsv: フィルタリングやDADA2工程で残ったリードの推移

final_report.txt: 解析全体のサマリーレポート

 

run_1/には､fastqディレクトリ1個目のQCや中間出力のfastqなどが保存される｡

multiQC report (raw qc)

 

実行方法

ショートリードは以下のフローで解析される｡

  • FastQCによる初期品質管理とMultiQCによるQCレポート生成
  • Cutadaptを用いたシーケンシングリードからのプライマー除去
  • FastQCによるCutadapt後の品質管理とMultiQCによるQCレポート生成
  • FIGAROを用いたDADA2の最適なトリミングパラメータの決定
  • DADA2によるアンプリコン配列バリアント推定

fastqファイルを配置したディレクトリを指定してランする｡イルミナの標準的なファイル命名則: <group>_<sample>_S##_L###_R[12]_001.(fastq|fq)[.gz]

に対応していれば自動で認識する｡

準備ができたら実行する｡ショートリードにはampwrap shortコマンドを使う｡"-i"でfastqを配置したディレクトリ､そして"-a"と"-F"でアンプリコン配列増幅に使用したプライマーを指定する｡レポジトリには典型的なV3ーV4領域のPrimer例があるので､これを指定している｡(論文になっているデータは論文を確認すること､自分の外注データなら､納品時の書類を確認すること)

#実行､シングルディレクトリ指定
ampwrap short -i fastq_dir -a CCTACGGGNGGCWGCAG -A GACTACHVGGGTATCTAATCC -l 411

#複数プロジェクトのディレクトリ指定(スペース区切り指定)
ampwrap short -i input_directory1 input_directory2 -a CCTACGGGNGGCWGCAG -A GACTACHVGGGTATCTAATCC -l amplicon_length --cutadapt_trim_by_length
  •  -i    Input directory with raw fq/fq.gz/fastq/fastq.gz files, add multiple space-separated paths for multiple runs e.g. -i path/to/run1 path/to/run2 path/to/run3 
  • -l    Length of the amplified sequence target (not including primers -->e.g., R926 - F515 = 411 - pF19 - pR20 = 372)
  • -a    Forward primer sequence
  • -A   Reverse primer sequence
  • -o    Output directory (default: ./ampwrap_short) 
  • -d    Database for taxonomy assignment. Choices: decipher_silva138 decipher_gtdb226 decipher_rdp18 dada2_silva_genus138 dada2_RDP_genus19 dada2_GG2_genus09 dada2_RefSeq_RDPv16 dada2_GTDB_r202
                            dada2_pr2_5.1.1 dada2_18S_silva132

複数実行モードでは、"--cutadapt_trim_by_length" を使用してプライマーが固定長トリミングされるようにすることが推奨されている｡

 

ナノポアロングリードは以下のフローで解析される｡

  • FastQCによる初期品質管理とMultiQCによるQCレポート生成
  • PorechopまたはPorechop_ABIによるアダプター除去
  • Cutadaptによるプライマー除去
  • NanoFiltによる長さと品質によるリードのフィルタリング
  • FastQCによるフィルタリング後の品質管理とMultiQCによるQCレポート生成
  • yacrdによるキメラの検出とフィルタリング（オプション）
  • EMUによる分類と存在量推定

ロングリードにはampwrap longコマンドを使う｡

ampwrap long -i input_directory -o output_directory

 

レポジトリより

  • 16S rRNA遺伝子の V4 や V3-V4 領域は、すべての菌種で全く同じ長さではなく、数〜数十塩基の個体差がある｡もし -lの数値を短めに設定してしまうとFIGAROはリードを切りすぎてしまう｡数値は期待される最も長いサイズを使用するのがベスト(増幅長 － Primerx2長)

コメント

リード長が揃っていないとエラーになりやすいので注意してください｡

引用

AmpWrap: a one-line fully automated amplicon metabarcoding 16S and 18S rRNA gene analysis Open Access

Lapo Doni , Alessia Marotta , Luigi Vezzulli , Emanuele Bosi

Bioinformatics Advances, Volume 5, Issue 1, 2025

 

関連