macでインフォマティクス

macでインフォマティクス

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

メタゲノムのビニングを行う COCACOLA

 

 アセンブリはコンティグを生成するが、それ以上の分類学的なプロファイリングや機能解析のためには、OTUに分類することが重要である。このOTUクラスタリングはビニングとも呼ばれる。しかしコンティグの正確なビニングは、ゲノム中のリピート配列、シークエンシングエラー、同じ種内の株間レベル変動などによるキメラアセンブリなどの理由のために依然として困難である。現在利用可能なビニング手法は、分類手法とクラスタリング手法に大別される。分類アプローチ(Classification approaches )はtaxonomy依存であり、すなわち、コンティグまたはリードから有意なtaxonへの割り当てに参照するデータベースが必要である。分類は、配列同一性による相同性(最も一般的な分類学的祖先に割り当てるMEGANなど)、またはオリゴヌクレオチド組成パターンおよび分類学的cladeのようなゲノムシグネチャーに基づく分類(PhyloPythiaやkrakenなど)がある。クラスタリングアプローチは参照データベースまたは分類学情報は必要ない。クラスタリングアプローチは、GC含量、テトラマー組成(論文より Albertsen et al、2013; Chatterji et al、2008; Yang et al、2010)またはInterpolated Markov Models(Kelley and Salzberg、2010)からのコンティグカバレッジプロファイル(Baran and Halperin、2012; Wu and Ye、2011)を利用する。

 最近ではコンティグのカバレッジプロファイルを用いてビニングする方法が開発されている。この考え方は、2つのコンティグが同じゲノムに由来する場合、複数のメタゲノムサンプルにわたるそれらのカバレッジプロファイルが高度に相関していることである。この方法は、コンティグプロファイルをテトラマー頻度と統合することによってさらに改善され得る。GroopM(Imelfort et al。、2014)は、視覚化されインタラクティブなパイプライン持ち、ユーザーは専門家の下でビンをマージして分割できる。他方、熟練者の介入がなければ、GroopMの自動ビニング結果はCONCOCTほど満足できるものではない(Alneberg et al、2014)。 CONCOCT(Alneberg et al、2014)は、コンティグをビンにクラスタリングするためにガウス混合モデル(GMM)を利用する。また、CONCOCTは、変分ベイズモデルの選択によって最適なOTU数を自動的に決定する。 MetaBAT(Kang et al、2015)は、対のコンティグについての積分距離を計算し、次いで、modifiled K-medoids algorithmによって反復的にコンティグをクラスター化する。 MaxBin(Wu et al。、2015)は同じゲノム間およびゲノム間の距離の分布を比較している。

 COCACOLAはシーケンスの構成、カバレッジ、Co-alignment、および複数のサンプルにわたるペアエンドリンケージを組み込んでいる。デフォルトでは、COCACOLAは、複数のサンプルにわたって配列合成とカバレッジを使用してビニングを実行する。CONCOCT、GroopM、MaxBin、MetaBATと比較されており、Precision、recallでよりよい成績を出し、また、スケーラブルかつ高速であると示されている。

 

 

インストール

 ubuntu16.0.4にて、conda createでpython2.7.15環境を作ってテストした。

本体 Github

ここではPythonバージョンをインストールする。

#conda createで仮想環境を作る
conda create -n cocacola_env -y python=2.7.14
source activate cocacola_env

#他の依存
conda install -y numpy scipy pandas scikit-learn cvxopt

python版COCACOLAはdropboxからダウンロードする。

link => here

docker imageもある。

docker pull bilalarxd/cocacola

python COCACOLA-python/cocacola.py

# python COCACOLA-python/cocacola.py 

usage: cocacola.py [-h] [--contig_file CONTIG_FILE]

                   [--abundance_profiles ABUNDANCE_PROFILES]

                   [--composition_profiles COMPOSITION_PROFILES]

                   [--edge_list EDGE_LIST] [--output OUTPUT]

                   [--clusters CLUSTERS]

cocacola.py: error: Data is missing, add file(s) using --contig_file <contig_file> and/or --abundance_profiles <abund_profiles> and/or --composition_profiles <comp_profiles>

(cocacola_env) root@c875fa54ab81:/data# python COCACOLA-python/cocacola.py -h

usage: cocacola.py [-h] [--contig_file CONTIG_FILE]

                   [--abundance_profiles ABUNDANCE_PROFILES]

                   [--composition_profiles COMPOSITION_PROFILES]

                   [--edge_list EDGE_LIST] [--output OUTPUT]

                   [--clusters CLUSTERS]

 

optional arguments:

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

  --contig_file CONTIG_FILE

                        The contigs file.

  --abundance_profiles ABUNDANCE_PROFILES

                        The abundance profiles, containing a table where each

                        row correspond to a contig, and each column correspond

                        to a sample. All values are separated with tabs.

  --composition_profiles COMPOSITION_PROFILES

                        The composition profiles, containing a table where

                        each row correspond to a contig, and each column

                        correspond to the kmer composition of particular kmer.

                        All values are separated with comma.

  --edge_list EDGE_LIST

                        The edges encoding either the co-alignment or the

                        pair-end linkage information, one row for one edge in

                        the format: contig_name_A contig_name_B weight. The

                        edge is undirected.

  --output OUTPUT       The output file, storing the binning result. If not

                        specified, the result is displayed directly on the

                        console.

  --clusters CLUSTERS   specify the number of clusters. If not specified, the

                        cluster number is estimated by single-copy genes.

 

 

テストラン

contigファイル、abundanceファイルとk-merプロファイルファイルを指定する。

cd COCACOLA-python/
python cocacola.py \
--contig_file data/SpeciesMock/input/SpeciesMock_Contigs_cutup_10K_nodup_filter_1K.fasta \
--abundance_profiles data/SpeciesMock/input/cov_inputtableR.tsv \
--composition_profiles data/SpeciesMock/input/kmer_4.csv \
--output data/SpeciesMock/result.csv
  • --abundance_profiles    ABUNDANCE_PROFILES:  The abundance profiles, containing a table where each row correspond to a contig, and each column correspond to a sample. All values are separated with tabs.

  • --composition_profiles    The composition profiles, containing a table where each row correspond to a contig, and each column correspond to the kmer composition of particular kmer. All values are separated with comma.

 

引用

COCACOLA: binning metagenomic contigs using sequence COmposition, read CoverAge, CO-alignment and paired-end read LinkAge

Yang Young Lu Ting Chen Jed A. Fuhrman Fengzhu Sun

Bioinformatics, Volume 33, Issue 6, 15 March 2017, Pages 791–798

 

 

メタゲノムシーケンシングリードをアセンブリしてvirusゲノム配列を探す自動パイプライン virMine

 

 真核生物および原核生物とは対照的に、ウイルスゲノムはごく一部のみがシーケンシングされ特徴付けられている。ウイルスのメタゲノム研究は、地球上でのウイルスの多様性についての理解を深めるうえで極めて重要である。海水(Breitbart et al、2002; Yooseph et al、2007; Hurwitz&Sullivan、2013; Brum et al、2015; Coutinho et al、2017; Zeigler Allen et al、2017;レビューBrum&Sullivan、2015)、土壌(Fierer et al、2007; Zablocki et al、2014; Adriaenssens et al、2017;レビューPratama&Van Elsas、2018)、淡水 (López-Bueno et al., 2009; López-Bueno et al., 2015; Roux et al., 2012; see review Bruder et al., 2016),、およびヒト微生物叢 (e.g., Reyes et al., 2010; Minot et al., 2011; Minot et al., 2013; Pride et al., 2012; Hannigan et al., 2015; Santiago-Rodriguez et al., 2015; Miller-Ensminger et al., 2018; see review Abeles & Pride, 2014) などの多数の生息地が調査されている。 最近、ヒトの微生物叢のウイルスメンバー(レビュー、Barr、2017; Keen&Dantas、2018年参照)および海洋環境(レビュー、Breitbart et al、2018年参照)が、かつて考えられていたよりも極めて重要な役割を果たすことを明らかにした。探索された環境にかかわらず、産生された圧倒的多数のウイルス配列は、特徴付けられたウイルス種に対して配列相同性を示さない。よく研究されている海洋ウイルス群集でさえ、予測されるコード領域の60%以上が完全に新規である(Coutinho et al、2017)。

 メタゲノミクスは、遺伝子マーカー(例えば、16S rRNA遺伝子)および未培養真核生物および原核生物種のゲノムを同定するのに有益であるが(Hug et al、2016)、viromesの調査は特有の課題に直面している(Bruder et al、2016; Rose et al 、2016)。第一に、細胞生物とは異なり、ウイルスに普遍的に保存された遺伝子はない。ウイルスは高度の遺伝的多様性にまたがり、本質的にモザイクである(Hatfull、2008年)。第二に、精製ビリオンをシーケンシングするときでさえ、シーケンシングデータはしばしば非ウイルス(ホスト)DNAを含む。これは、ウイルスゲノムDNAがしばしばサンプル中の宿主細胞または他の生物よりも桁違いに少ないという事実によってさらに複雑になる。ウイルスメタゲノミクスの実験手順の開発(例えば、ConceiçãoNeto et al、2015; Hayes et al、2017; Lewandowska et al、2017)に加えて、いくつかのバイオインフォマティクスソリューションが、ウイルス配列の検出を助けるために作成されている。混合コミュニティ(例えばRoux et al、2015; Hatzopoulos、Watkins&Putonti、2016; Yamashita、Sekizuka&Kuroda、2016; Ren et al、2017; Amgarten et al、2018;総説Hurwitz et al、2018参照) ; Nooij et al、2018)。第三に、現存のウイルスデータ貯蔵所はウイルス種を十分含んでいない。したがって、細菌メタゲノム解析のためのものなどの、配列相同性を同定することに依存するツール(レビューNayfach&Pollard、2016を参照)は、viromes研究における用途が限られている。

 単一または少数のウイルス種を含有する試料からのウイルスゲノムの同定は、非ウイルス配列の大きなバックグラウンドが存在する場合でさえも比較的簡単である。そのような例は、臨床サンプルの潜在的なウイルス性病原体の検索であろう。 VIP(Li et al、2016b)、VirAmp(Wan et al、2015)、およびVirFind(Ho&Tzanetakis、2014)を含むソフトウェアツールは、そのような場合のために特別に設計された。しかしながら、それらは既知のウイルス分類群の単離に限定されている。複雑なウイルスコミュニティは、さらに大きな課題を抱えている。通常、2つのアプローチのうちの1つが取られる。第1のアプローチは、配列属性、例えばそれらのヌクレオチド使用プロファイル(Ren et al、2017)および/またはコンティグカバレッジに基づいてメタゲノムデータセットからコンティグを同定する (see reviews Sharon & Banfield, 2013; Garza & Dutilh, 2015; Sangwan, Xia & Gilbert, 2016)。第2のより頻繁に追求されている方法は、既知のウイルス配列、例えば、Phage Eco-Locator(Aziz et al、2011)、VIROME(Wommack et al、2012)、MetaVir(Roux et al、 2014)、VirSorter(Roux et al、2015)、MetaPhinder(Jurtz et al、2016)、VirusSeeker(Zhao et al、2017)、およびFastViromeExplorer(Tithi et al、2018)である。ツールMARVELは、ゲノムの特徴(遺伝子密度およびストランドシフト)および配列相同性に基づいてテールファージ配列を予測しながら、2つのアプローチを統合している(Amgarten et al、2018)。どのようなアプローチを取ったとしても、手作業によるキュレーションと検査はプロセスの中で重要なステップとなる。

 (一部略)

 ここで本著者らはメタゲノムデータセット内のウイルスゲノムの同定のためにvirMineを提示する。 virMineは発見プロセスを自動化する。生のシーケンスからクオリティチェックを経てアセンブリアノテーションを行う。 virMineには、公に利用可能なさまざまなツールとユーザー定義の基準が組み込まれている。地球上でのウイルスの多様性についての限られた知識に基づいてウイルスの「シグネチャ」を検索する以前のバイオインフォマティクスツールとは対照的に、virMineは細胞生物に利用可能な豊富な配列データを利用する。したがって、ウイルス(バクテリオファージおよび真核ウイルス)の発見は、ウイルスではないことがわかっているものを除外するプロセスを通じて行われる。次に、「非ウイルス性」ではない配列(すなわち、推定ウイルス配列)をウイルス配列のデータベースと比較する。この比較は、既知のウイルス配列と下流の分析のための新規ウイルスを表し得るものと類似の推定ウイルス配列を区別する。このツールのベータ版を用いて、尿中メタゲノムデータセットからウイルス配列を単離した(Garretto et al、2018)。ここでは、4つのケーススタディを使用してこのツールの有用性を説明する:合成データセット、gut マイクロバイオーム、尿中viromes、および淡水viromes、既知のウイルスの新型ならびに新規のウイルスゲノムの同定。

 このパイプラインは、PythonとBioPythonライブラリを使用して、既存のツールと新しいアルゴリズムを統合している(Cock et al、2009)。論文図1はvirMineが採用しているプロセスを示している。ツールの重要な側面はその柔軟性である。virMineはモジュラー設計されていて、ユーザーが個々に機能性にアクセスするか、または完全なパイプラインを実行することを可能にする。このリリースにはいくつかの方法が組み込まれている(論文表1)が、新しいツールを簡単に追加することができる。さらに、ターゲットを絞った分析を容易にするために、プログラミングの専門知識がなくてもユーザーはフィルタリングオプションとカスタマイズオプションを利用できる。

 

f:id:kazumaxneo:20190706220609j:plain

Figure 1: Overview of virMine pipeline. 論文より転載

 

 

 インストール

dockerイメージをビルドしてテストした(ホストOS macos10.12)。

本体 Github

#dockerイメージをbuildして使う。
git clone https://github.com/thatzopoulos/virMine.git
cd virMine/
docker build -t virmine .

#dockerを立ち上げてから

>python2.7 /virMine.py -h

# python2.7 /virMine.py -h

usage: virMine.py [assembly options] [filter options] [databases] -o output_path

 

optional arguments:

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

  -a {spades,metaspades,megahit,all3}, --assembler {spades,metaspades,megahit,all3}

                        Assembly method. all3 will test all three and pick the

                        best based on N50 scores (Paired-end reads only).

  -A <filename>, --assembled_contigs <filename>

                        Contigs file. Start analysis from an existing

                        assembly.

  -s <filename>, --single_reads <filename>

                        Single reads.

  -p <filename> <filename>, --paired_end_reads <filename> <filename>

                        Paired-end reads. List both read files.

  -o <directory>, --output_path <directory>

                        Directory to store all the resulting files (required)

  -t <int>, --num_threads <int>

                        Number of processors to use.

  -b {blastx,blastp}, --blast_option {blastx,blastp}

                        BLAST method to classify as viral or non-viral. Note:

                        blastp is faster.

  --version             show program's version number and exit

 

filter options:

  -m <int>, --min_contig_size <int>

                        Minimum contig size.

  -M <int>, --max_contig_size <int>

                        Maximum contig size.

  -c <int>, --min_coverage <int>

                        Minimum coverage.

  -cov <float>, --min_SPAdes_cov <float>

                        Minimum SPAdes cov value.

  -g <filename>, --genes_of_interest <filename>

                        PROTEIN sequences of interest. FASTA format. Contig

                        must contain homology to genes of interest.

 

database options:

  -v <filename>, --viral <filename>

                        Viral PROTEIN sequence collection. FASTA format.

                        (required)

  -nv <filename>, --non_viral <filename>

                        Non-viral PROTEIN sequence collection. FASTA format.

                        (required)

 

 

データベースの準備

ここではdockerイメージをランして中で行っている

#docker内で行う
docker run -itv $PWD:/data -w data virmine

#libraryがないので先に導入
pip install ftputil

#データベースダウンロードとビルド(カレントのパスは/dataとして)*1
git clone https://github.com/thatzopoulos/virMine.git
python2.7 virMine/virmine_make_dbs.py

no_phage_bact.fastaとall_viral_faa.fastaができる。 

 

 

実行方法

データベース準備ステップですでにdocker runしているとして進める。

python2.7 virMine.py -a spades \
-p pair_1.fastq pair_2.fastq \
-v all_viral_faa.fasta -nv no_phage_bact.fasta \
-o outputFolder
  • -v    Viral PROTEIN sequence collection. FASTA format (required)
  • -nv   Non-viral PROTEIN sequence collection. FASTA format (required)
  • -p     Paired-end reads. List both read files
  • -o     Directory to store all the resulting files (required)
  • -a     assembler {spades,metaspades,megahit,all3 

データが巨大だったり複雑だとspadesではヘビーになっていく。必要に応じてmegahitに切り替える。

 

 

引用

virMine: automated detection of viral sequences from complex metagenomic samples
Garretto A, Hatzopoulos T, Putonti C

PeerJ. 2019 Apr 10;7:e6695

 

関連

 

 


 

 

*1ビルドにはかなりの時間がかかるので、データベースををビルドしてから抜け、commitしてローカルに保管しておくと2回目以降使いやすい。

 

似た名前のツールがありますが別のものです。ご注意ください。

メタゲノムのファージ配列分析webサーバー VirMiner - macでインフォマティクス

 

 

タイプ株を中心にバクテリアの表現型情報をまとめたデータベース BacDive

2020 10/16  タイトル変更、誤字修正

2021 10/7 画像一部更新

 

 原核生物は、研究開発との関連性が高い多種多様な表現型形質を発現する。バクテリアメタデータホットスポットとしてよく利用できるのは、最初の(一次)文献で報告された種の説明と、生物資源センター(BRC)によって管理されているデータベースである(ref.1)。それでも、両方の情報源にアクセスすることは困難であり、そして特定の表現型データの検索はしばしば面倒である。全ゲノム配列の数が急速に増えており(2018年7月5日現在までに、10 957ゲノムがNCBIゲノムデータベースに登録されている)、細菌表現型の予測および分析にますます使用されている。 Microbial trait analyzer(ref.2)のようなツールは、ゲノム配列から表現型の属性を自動的に予測し、シーケンシングされた生物の潜在的な表現型についての迅速な概観を可能にする。それにもかかわらず、これらのアプローチは、モデル生物の機能的配列情報に対する相同性の枠組みに依存しており、それはしばしば調査対象の生物とは遠く離れている。したがって、実験室の表現型データによる結果の検証は不可欠になる。さらに、包括的でよく構造化された表現型データを使用してアルゴリズムの精度を向上させ、それによって微生物の生物多様性に対するより良い予測と新しい洞察を可能にすることができる。

 BacDive - the Bacterial Diversity Metadatabase(https://bacdive.dsmz.de)は、広範囲のバクテリアおよび古細菌株にわたってメタデータを集め構造化し、メタデータを提供している。主な情報源は、BRCの内部データベース(例:DSMZ、CIP)、未発表の研究者のコレクション(例:Reichenbach collection of Myxobacteria、Wink compendium of Actinobacteria (ref.3))および一次文献の種の記述(例:International Journal of Systematic and Evolutionary Microbiology)である 。これらの情報源にアクセスし、関連するデータを整理して組み合わせることで、BacDiveは分類学、形態学、生理学、培養、地理的起源、原核生物の応用および分子データに関する包括的なデータセットを提供する。全てのデータは株に関連しており、それぞれの参考文献にリンクされている。

 2012年4月の発表以来(ref.3)、BacDiveは実質的にさらに発展している。データフィールドの数が増え、検索機能が改善され、Restful Webサービスが実装された(ref.4)。最近の開発は、一次文献から手作業でアノテーション付けしたデータおよびBRCの内部ファイルおよび個人的なコレクションから動員したデータを用いた株データセットの修正に焦点を当てていた(ref.1)。データ内容およびデータフィールドの増加、ならびにそれに付随するグラフィカルユーザインターフェース(GUI)および検索機能への適応について、以下に説明する。

 現在のBacDiveリリース(2018年7月)は、2637属および13 569種に分布する63 669の菌株を対象としている。 BacDiveに含まれる12 715株は、その種のタイプ株を表す。 データベースはProkaryotic Nomenclature Up-to-Date (PNU)データベースにリストされている14 395の正当に公表されている型株の88%をカバーしている。 過去数年にわたり、既存の株データセットの充実は、新しい株に関する情報の収集よりも優先されていた。 したがって、2015年にBacDiveがカバーした株の総数(53 978)と比較すると、新しい株の増加は18%だが、1株あたりのデータ量はすべてのセクションで大幅に増加した(論文表1)。 セクションの形態と生理学(485%)、培養と成長条件(391%)、そして単離、サンプリングと環境情報(247%)は特に増加している。

新しい種類のデータ

 Analytical Profile Index(API)テスト、脂肪酸プロファイル、抗生物質感受性テストなどの表現型テストのデータを統合することで、新しいデータタイプがBacDiveに導入された。(一部略)

 The API® strip testは、微生物を迅速に識別するための十分に確立された生理学的試験システムである。製造業者bioMérieuxには22のテストストリップを提供しており、各テストストリップは最大50のシングルテストで構成されている。(一部略)現在、BacDiveは、5237細菌株のデータを含む8977API®テストデータを提供している。これは世界最大のAPI®テスト結果のコレクションであり、多様な細菌株に関する表現型情報を提供している。

 最近、脂肪酸プロファイルと抗生物質感受性のデータもBacDiveに統合された。どちらのデータタイプも分類法、したがってBacDiveの大規模なユーザーグループにとって非常に重要である。最初に、395の脂肪酸プロファイルおよび475の抗生物質感受性試験の結果が集められた。このコレクションは継続的に修正されている。

文献データ

 最も包括的だが同時にメタデータのための最も多様な情報源は一次文献である。メタデータはさまざまなジャーナルに散在している。標準化と一意の識別子がなく、文学からの動員は困難である。細菌および古細菌メタデータについては、the International Journal of Systematic and Evolutionary Microbiology (IJSEM) が特に重要であり、最大80%の新種の説明がこのジャーナルに掲載されている。 IJSEMの種の記述は表現型メタデータに富んでいるだけでなく、かなり保存された構造も持っている。これにより、ジャーナルはメタデータ注釈の主要な情報源になる。 IJSEMの種の記述に基づいて、最大152の異なるデータフィールドに手動でデータをアノテーションするシステムを設定した。これらのデータは、コントロールされた語彙と割り当てられたグローバル識別子を使用して変換された。最初に、523種類の株に関するデータが収集された(ref.1)。それ以来、本著者らはすべての新しく発表された種の説明からのデータに継続的にアノテーションをつけることによって進めた。これまでのところ、1387株のデータには手動でアノテーションが付けられていた。 2017年に、Barberán et al.,(ref.5)は5129株および最大30のデータフィールドについてのIJSEMからのアノテーション付きデータを含むデータセットを発表した。慎重な調整後、この供給源からの4081株のデータがBacDiveに組み込まれた。現在のリリースでは、BacDiveは合計5468株についてのIJSEMの種の説明からのアノテーション付きデータを提供している(https://bacdive.dsmz.de/advsearch?site=advsearch&searchparams%5B84%5D%5Bsearchterm%5D=IJSEM&advsearch=sea

(以下略) 

 

チュートリアル

https://bacdive.dsmz.de/tutorials

BacDive Quickstart tutorial (音が鳴ります)

 

他にもいくつか動画が準備されています。

https://bacdive.dsmz.de/tutorials

 

wiki

https://en.wikipedia.org/wiki/BacDive

 


使い方

https://bacdive.dsmz.de/ にアクセスする。 

f:id:kazumaxneo:20211007223715p:plain

(画像更新)

 

タイプ株を中心に様々な細菌の情報がまとめられている。

f:id:kazumaxneo:20211007224150p:plain

(グラフの要素をクリックするとその情報ソースにジャンプする)

 

その場で素早く検索できる。f:id:kazumaxneo:20190705230031p:plain

 

 

Acinetobacter baylyi A7のページを開いた。

f:id:kazumaxneo:20190705230130p:plain

 

strain A7はプロテオバクテリア門のガンマプロテオバクテリア綱、アシネトバクター属のAcinetobacter baylyiという菌であることが分かる。すべての情報には左端に引用元文献のリンクがついている。

f:id:kazumaxneo:20190705230354p:plain

 

資化できる炭素/窒素源、耐性など様々なテスト(API®テスト)の結果がまとめられている。

f:id:kazumaxneo:20190705230929p:plain

API®テストについて

f:id:kazumaxneo:20190705231102p:plain

 

培養条件もまとめられる。

f:id:kazumaxneo:20190705231405p:plain

 

単離場所

f:id:kazumaxneo:20190705231424p:plain

 

ページ下には引用文献のリンクがついている。中にはもしかして間違っていたり、結果が揺れ動くようなメタ情報(+/-)もあるかもしれない。そのような時に引用文献のリンクがあるととても役に立つ。

 

他の機能について

API® test finder

メニューからAPI® test finderを選択

f:id:kazumaxneo:20190705231819p:plain

 

各テスト項目の結果から菌を絞り込める。

f:id:kazumaxneo:20190705231950p:plain

空白部分をクリックして+か-か選んでいく。

 

テスト内容が不明な場合、テストの文字部分をクリックすることで解説ページにジャンプできる。

f:id:kazumaxneo:20190705232218p:plain

 

TAXplorer

taxonomyの枠組みに従ってバクテリアアーキアを検索できる。

f:id:kazumaxneo:20190705235254p:plain

 

Advanced search

詳細な条件を指定した検索

f:id:kazumaxneo:20190705235533p:plain

 

Isolation sources

f:id:kazumaxneo:20190705235614p:plain

 

左上のFilterボタンをクリックすると、菌を単離した環境から絞り込める。

f:id:kazumaxneo:20190705235937p:plain

 

Enviromental => Aquatic => Freshwaterを選択して、Submit Filterボタンをクリック。

f:id:kazumaxneo:20190706001220p:plain

該当する菌が表示されている。

 

テーブル上のde-/select all displayed strainsボタンをクリックして、該当する菌の行を全選択した。

f:id:kazumaxneo:20190706001406p:plain

 

この状態でテーブル上のadd selectionをクリックすれば、該当する菌のメタデータを含むCSVをダウンロードできる。

f:id:kazumaxneo:20190706001634p:plain

exportボタンをクリックすればCSV形式でダウンロードされる。

 

Mapをクリックすると国ごとの単離された菌体数を視覚化できる。

f:id:kazumaxneo:20190706000505p:plain

単離された菌体数は色の濃さで表される。

 

単離環境の内訳はkronaで視覚化されている(右上のopen chartをクリック)。

f:id:kazumaxneo:20190706000843p:plain

 

菌のwikiでもよく引用ソースになっています。情報の精度が高いため、メタゲノムやコンタネーションの分析にも利用できると思います。

引用

BacDive in 2019: bacterial phenotypic data for High-throughput biodiversity analysis
Lorenz Christian Reimer Anna Vetcininova Joaquim Sardà Carbasse Carola Söhngen Dorothea Gleim Christian Ebeling Jörg Overmann
Nucleic Acids Research, Volume 47, Issue D1, 08 January 2019, Pages D631–D636

 
BacDive--The Bacterial Diversity Metadatabase in 2016
Söhngen C, Podstawka A, Bunk B, Gleim D, Vetcininova A, Reimer LC, Ebeling C, Pendarovski C, Overmann J

Nucleic Acids Res. 2016 Jan 4;44(D1):D581-5


BacDive--the Bacterial Diversity Metadatabase
Söhngen C, Bunk B, Podstawka A, Gleim D, Overmann J

Nucleic Acids Res. 2014 Jan;42(Database issue):D592-9

 

関連


アセンブリのグラフを分析する Asgan

 

 Asgan - [As] sembly [G] raphs [An] alyzer - は、アセンブリグラフを分析するためのツールである。 このツールはGFA形式の2つのアセンブリグラフを入力として受け取り、そのグラフの最小セットの相同配列(シンテニーパス)を見つけ、見つかったパスに基づいてさまざまな統計を計算する。

 

詳細はこちら

 


 

 インストール

ubuntu16.0.4のminiconda3-4.3.30環境でテストした(docker使用、ホストOS macos10.12)。*1

依存

  • networkx、decorator等
  • graphviz (graph visualization)
conda install -y decorator networkx 
conda install -y graphviz

> dot -?

$ dot -?

Usage: dot [-Vv?] [-(GNE)name=val] [-(KTlso)<val>] <dot files>

(additional options for neato)    [-x] [-n<v>]

(additional options for fdp)      [-L(gO)] [-L(nUCT)<val>]

(additional options for memtest)  [-m<v>]

(additional options for config)  [-cv]

 

 -V          - Print version and exit

 -v          - Enable verbose mode 

 -Gname=val  - Set graph attribute 'name' to 'val'

 -Nname=val  - Set node attribute 'name' to 'val'

 -Ename=val  - Set edge attribute 'name' to 'val'

 -Tv         - Set output format to 'v'

 -Kv         - Set layout engine to 'v' (overrides default based on command name)

 -lv         - Use external library 'v'

 -ofile      - Write output to 'file'

 -O          - Automatically generate an output filename based on the input filename with a .'format' appended. (Causes all -ofile options to be ignored.) 

 -P          - Internally generate a graph of the current plugins. 

 -q[l]       - Set level of message suppression (=1)

 -s[v]       - Scale input by 'v' (=72)

 -y          - Invert y coordinate in output

 

 -n[v]       - No layout mode 'v' (=1)

 -x          - Reduce graph

 

 -Lg         - Don't use grid

 -LO         - Use old attractive force

 -Ln<i>      - Set number of iterations to i

 -LU<i>      - Set unscaled factor to i

 -LC<v>      - Set overlap expansion factor to v

 -LT[*]<v>   - Set temperature (temperature factor) to v

 

 -m          - Memory test (Observe no growth with top. Kill when done.)

 -m[v]       - Memory test - v iterations.

 

 -c          - Configure plugins (Writes $prefix/lib/graphviz/config 

               with available plugin information.  Needs write privilege.)

 -?          - Print usage and exit

 本体 Github

git clone --recurse-submodules https://github.com/epolevikov/Asgan
make -C Asgan/lib/minimap2

 > python asgan.py -h

$ python asgan.py -h

usage: asgan.py [-h] [--input-query INPUT_QUERY] [--input-target INPUT_TARGET]

                [--out-dir OUT_DIR]

 

optional arguments:

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

  --input-query INPUT_QUERY

  --input-target INPUT_TARGET

  --out-dir OUT_DIR

 

 

テストラン

2つのアセンブリのGFAファイルを指定する。

cd Asgan
python asgan.py \
--input-query=test/NCTC9016-Flye.gfa \
--input-target=test/NCTC9016-Canu.gfa \
--out-dir=output

テストデータはCanuとflyeのアセンブリ結果のGFAを指定している。graphvizに対応したgzファイル(.dot)が出力される。

 

stats

f:id:kazumaxneo:20190704235628p:plain

synteny_paths.txt

f:id:kazumaxneo:20190704235715p:plain

 

視覚化

cd output/
dot -Tpng -O adjacency_graph_query.gv
dot -Tpng -O adjacency_graph_target.gv

 

 

adjacency_graph_target.gv.png

下の図は、フォワードのパスが実線、リバースのパスが破線で表されている。ターゲット(ここではcanu)のアセンブリグラフパスは、フォワードが-2, -1, -4, +3、リバースが-3, +4, +1, +2のパスになっている。

f:id:kazumaxneo:20190622175529p:plain

 

adjacency_graph_query.gv.png

クエリとターゲットのパス番号は対応している。番号がついていない黒い線はターゲットのcanuには対応しないパスで、リピート解像されていない。

f:id:kazumaxneo:20190622175542p:plain

 

アセンブルを真面目に考えたい人にとっては便利なツールですね。まだ開発中のようですが、早めに紹介しました。Documentationが出てきたら詳細を追記します。

引用

GitHub - epolevikov/Asgan: A tool for analysis of assembly graphs

 

https://zenodo.org/record/3198701#.XQ3v8C_ANts

 

参考

https://dev.classmethod.jp/tool/graphviz-beginner/

 

ubuntu graphviz-dev パッケージ

Ubuntu – Package Search Results -- graphviz

 

*1

仮想環境でテストしましたが、導入時トラブりました。いきなり本番環境にいれないほうがいいと思います。

メタゲノムcontigのビニングとアノテーションwebサーバー BusyBee Web

 

 メタゲノムシーケンシング、すなわち微生物混合群集から無差別に抽出されたDNAの全ゲノムシーケンシングは、分類学的組成および環境マイクロバイオームの機能的可能性を研究するために首尾よく使用されてきた(ref.1-4)。従来の単離培養工程の独立性は、費用および時間の削減、ならびにこれまで人工実験室条件下での培養の試みに抵抗してきた微生物を特徴付けることを可能にするので、しばしば利点と考えられる(ref.5,6)。メタゲノムシーケンシングは主に基礎研究に使用されてきたが、臨床現場でのその可能性は最近実証された(ref.7,8)。さらに、第3世代シーケンシング技術、例えばPacific Biosciences(PacBio)またはOxford Nanopore Technologies(ONT)の論文が発表されており、微生物混合群集の長期にわたる研究に基づいた研究が可能になっている(ref.9–11)。

 計算を用いたメタゲノムシーケンシングデータからの個々の生物(またはclosely relatedな生物集団)のレベルゲノム配列の回収は「ビニング」と呼ばれる。現在のビニング手法の本体は、(i)基準依存型手法と(ii)基準非依存型手法とに大別することができる。リファレンス依存のビニング手法は、通常、非常に短い実行時間と、高い感度および精度によって特徴付けられる(ref.12〜16)。しかしながら、これらのアプローチは、設計上、データベースに存在する参考文献の一部であるかまたはそれにclosely relatedな生物に由来する配列に対して最もよく機能し、そしてこれまで特徴付けられていない微生物に由来するゲノム配列でチャレンジングになる。対照的に、リファレンスに依存しないビニング手法は、入力データのみから配列クラスター構造を推測するため(ref.17–20)、主に配列構成に基づいている。事前知識に依存せず、複数のサンプルにわたるabundanceの共変動に頼る手法が出現している(ref.21〜24)。それらのリファレンス非依存性のために、これらのアプローチは現在のリファレンスゲノムデータベースでは限定的にしかない環境の分析に特に有用であり、しばしば「未分類」の配列の分解を可能にする。しかしながら、リファレンスに依存しないビニングは、かなりの量のCPU時間、特定の閾値を超えるシーケンス長、例えば1000 bp以上、理想的には独立したサンプルを必要とすることが多い。 さまざまなビニングWebサーバーが存在するが、これらは主にリファレンス依存のアプローチ(ref.15、25〜27)に基づいているか、専用のコンピューティングリソースおよび/またはユーザートレーニング(28、29)が必要となる事前の計算が必要である。

 この論文では、メタゲノムシーケンシングデータセットのbootstrapped supervised binning(BSB)を実装するWebアプリケーションであるBusyBee Webサーバーを紹介することによって、現在利用可能なリファレンスに依存しないビニングツールを拡張する。本ビニングアプローチは、リファレンスデータベースに頼るのではなく、入力からトレーニングデータを「ブートストラップ」することによって、教師なしと教師付きの機械学習アプローチを組み合わせたものである。 BusyBee Webは、入力として単一のFASTAフォーマットのファイルのみを必要とし、シーケンスを母集団分解ビンに自動デコンボリューションする。 BSBでは、教師なしアプローチを使用して、クラスタがシーケンスのサブセットに対してde novoで定義される(ref.30–32)。このステップに続いて、応答/従属変数としてクラスターラベルを使用して、ランダムフォレストベースの分類器を学習する(教師部分)。ビニングをさらに加速するために、監視されていない部分(圧縮)の間にデータ点がそれらの最近傍(代表者)の代表として働くようにランダムにサンプリングされる任意選択の「圧縮」ステップが実施される。代表者とその仲間は、その後、それぞれの代表者のde novoクラスターラベル(展開)と組み合わせて、監視対象部分で使用される。したがって、トレーニングセットのサイズは、ランダムにサンプリングされた代表的なデータポイントを使用するだけの場合と比較して大きくなる。最終的に、すべてのシーケンス(デフォルトでは500 bp以上)には、ブートストラップ訓練された分類器を使用してラベルが割り当てられ、それによって最終的なビンのセットが定義される。クラスタリング/ビニング結果の検査のために、データ固有の推測された構造の2D散布図がユーザに提示される。これを補足するために、ビン品質、すなわち完成度、汚染度、および歪みの不均一性の推定値が計算され視覚化される。さらに、配列はKrakenを用いて分類学的にアノテーション付けられ、抗生物質耐性遺伝子の機能的アノテーションが行われる。すべてのビニングおよびアノテーションステップはWebサーバーによってユーザーに対して透過的に自動的に実行されるため、専用のコンピューティングリソースや特別なユーザートレーニングは必要ない。さらに、カスタムのシーケンス毎のアノテーションをユーザがアップロードすることができる。特定の関心のあるシーケンスをハイライトするために、そしてBusyBee Webは生成された結果をダウンロードするためのオプションを提供する。

 

Tutorial

https://ccb-microbe.cs.uni-saarland.de/busybee/tutorial/

Github

 

使い方

https://ccb-microbe.cs.uni-saarland.de/busybee/ にアクセスする。

f:id:kazumaxneo:20190703230408p:plain

Submit new jobをクリック。

f:id:kazumaxneo:20190704002159p:plain

 

 FASTA配列 (contig / long read) をアップロードする(最大100MB)。

f:id:kazumaxneo:20190704002438p:plain

バグ防止のため、defaultではシーケンス名が短く変更される(ヘッダーが全て20文字以下ならOFFにしてO.K)。

 

Binningだけでなくアノテーションも実行するならEnabledに切り替える。

f:id:kazumaxneo:20190704002532p:plain

Taxonomic annotationにはkrakenを使ったtaxnomic assignmentが実行される。Functional annotationにはResfamsを使った抗生物質耐性遺伝子アノテーション

 

必要に応じてパラメータを変更する。

f:id:kazumaxneo:20190704002535p:plain

submitボタンをクリックしてジョブを開始する。サンプルの複雑さとサーバーの混雑度に応じて、結果が出るまで数分〜1時間ほどかかる。

 

結果

ビンング結果は二次元の散布図として視覚化される。

複雑性が低いデータ。

f:id:kazumaxneo:20190703231137p:plain

複雑性が低いサンプルであれば、ラン時に cluster point threshold をdefaultの2000-bpから500-bpなどに下げることで回収率を上げることができる。

 

 

散布図をドラッグして囲むと拡大できる。

f:id:kazumaxneo:20190704005731p:plain

右上のボタンからも拡大可能。右上の弧を描いた矢印ボタンをクリックすると拡大がリセットされる。

 

右のメニューから表示にするBinを選択できる。

f:id:kazumaxneo:20190704005732p:plain

 

Bin1とBin5を非表示にした。

f:id:kazumaxneo:20190704010323p:plain

 

Taxonに切り替えれば、プロットは分類階級によって色分けされる。

f:id:kazumaxneo:20190704005705p:plain

 

Functional annotationをONにして解析していれば、右のメニューから抗生物質耐性遺伝子のプロットを表示できる。

f:id:kazumaxneo:20190704011728p:plain

バンコマイシン耐性関連遺伝子を視覚化した。抗生物質耐性遺伝子のプロットサイズは図の下にあるPoint size annotation pointsバーから変更可能。

 

複雑性が高いデータ。右側の交差してしまっているビンを非表示にした。

f:id:kazumaxneo:20190704011012p:plain

 

配列数とビンクオリティ

f:id:kazumaxneo:20190704011241p:plain

 

各ビン/クラスターの分類組成

f:id:kazumaxneo:20190704011245p:plain


 ビニング結果は下のボタンからダウンロードできる。

f:id:kazumaxneo:20190704011455p:plain

引用
BusyBee Web: metagenomic data analysis by bootstrapped supervised binning and annotation
Cedric C. Laczny, Christina Kiefer, Valentina Galata, Tobias Fehlmann, Christina Backes, Andreas Keller

Nucleic Acids Res. 2017 Jul 3; 45(Web Server issue): W171–W179

 

アセンブリの前処理としてロングリードのキメラ領域(低オーバーラップ領域)を除く yacrd

2019 コマンドの誤り修正

2020 3/30  バージョンによるコマンドの違いを記載

2020 3/31 version0.6.0のコマンドを一番下に追記

2020 4/23 論文追記

 

 第三世代DNAシーケンシング法(PacBio、オックスフォードナノポア)は、リファレンスゲノムの構築(デノボアセンブリ)のための重要な技術となりつつある。この種のデータに対する新しいバイオインフォマティクス手法が急速に登場している。
 一部のロングリードアセンブラは、アセンブリ前にリードに対してエラー訂正を実行する。訂正は、第3世代リードの高いエラー率を減らし、アセンブリを扱いやすくするのに役立つが、時間とメモリを消費するステップでもある。最近のアセンブラ(例:Li(2016); Ruan and Li(2019)など)は、未訂正の未加工のリードを直接アセンブルする方法を見つけた。したがって、ここでは未訂正アセンブリのみに焦点を当てる。この設定では、アセンブリの品質は、キメラリードと非常に誤った領域(Myers、2015)の影響を受ける。
 DASCRUBBERプログラム(Myers、2017)は、リードの「スクラビング」の概念を導入した。これは、他の方法で塩基を修正することを試みることなく、リード中の問題のある領域を迅速に除去する。その考えは、リードをスクラブすることは訂正よりも軽量の操作であり、それゆえ高性能で訂正のないゲノムアセンブラに適しているということである。
 DASCRUBBERは、リードのall-against-allマッピングを実行し、リードごとにパイルアップを作成する。次に、マッピング品質を分析して、推定上高いエラー率の領域を決定する。これは、パイルアップ内の他のリードからの同等で高品質の領域に置き換えられる。 MiniScrub(LaPierre et al、2018)は、オーバーラップ検出に使用されたアンカーの位置を記録するために、Minimap2(Li、2017)の修正版を使用する別のスクラビングツールである。 MiniScrubはリードごとにアンカー位置をイメージに変換する。次に、畳み込みニューラルネットワークが低品質のリード領域を検出して削除する。
 リードスクラビングのさらに上流にあるもう1つの問題は、リード間のオーバーラップの計算である。オーバーラップの保存はディスクを大量に消費し、著者らの知る限りでは、その潜在的に高いディスクスペースを最適化する試みはこれまでになかった。
 本稿では、ロングリードのアセンブラの初期段階を一緒に最適化する2つのツールを紹介する。 1つは、高速で効果的なリードのスクラビング用のyacrd(for Yet Another Chimeric Read Detector)である。もう1つは、リード間で検出される重複をフィルター処理するfpa(フィルターペアアライメント)である。

 DASCRUBBERやMiniScrubと同様に、yacrdはリードの低品質領域は他のリードでは十分にサポートされていないという仮定に基づいている。そのような領域を検出するために、yacrdはMinimap2を使用してall-against-allリードマッピングを実行してから、各リードのカバレッジを計算する。 DASCRUBBERおよびMiniScrubとは対照的に、yacrdはMinimap2によって与えられるおおよその位置マッピング情報のみを使用する。これは時間のかかるアライメントステップを回避する。これはベースレベルのアライメントを持たないことを犠牲にしているが、これはスクラビングを実行するのに十分であることが判明している。カバレッジが特定のしきい値(デフォルトでは4に設定されている)を下回った場所でリードが分割され、カバレッジの低い領域が完全に削除される。

 

f:id:kazumaxneo:20190703004847p:plain
Githubより転載

 

2022/03/02

 

2020 3/3追記

https://twitter.com/pierre_marijon/status/1234783942986944512

f:id:kazumaxneo:20200303215507p:plain

 

インストール

ubuntu16.0.4のminicona3.4.0.5環境でテストした。

依存

  • Rust in stable channel
  • libgz
  • libbzip2
  • liblzma
  • minimap2 (リード同士のマッピングに必要)
conda install -c bioconda -y minimap2

本体 Github 

#bioconda (link) 0.6とはコマンドが異なる。0.5.1を導入してテストした。
conda install -c bioconda -y yacrd==0.5.1

#cargo
cargo install yacrd

yacrd -h

$ yacrd -h

yacrd 0.5.1 Omanyte

Pierre Marijon <pierre.marijon@inria.fr>

Yet Another Chimeric Read Detector

 

USAGE:

    yacrd [SUBCOMMAND]

 

FLAGS:

    -h, --help       Prints help information

    -V, --version    Prints version information

 

SUBCOMMANDS:

    chimeric     In chimeric mode yacrd detect chimera if coverage gap are in middle of read

    help         Prints this message or the help of the given subcommand(s)

    scrubbing    In scrubbing mode yacrd remove all part of read not covered

> yacrd chimeric -h

$ yacrd chimeric -h

yacrd-chimeric 

In chimeric mode yacrd detect chimera if coverage gap are in middle of read

 

USAGE:

    yacrd chimeric [FLAGS] [OPTIONS]

 

FLAGS:

    -j, --json       Yacrd report are write in json format

    -h, --help       Prints help information

    -V, --version    Prints version information

 

OPTIONS:

    -i, --input <input>...

            Mapping input file in PAF or MHAP format (with .paf or .mhap extension), use - for read standard input (no

            compression allowed, paf format by default) [default: -]

    -o, --output <output>

            Path where yacrd report are writen, use - for write in standard output same compression as input or use

            --compression-out [default: -]

    -f, --filter <filter>...

            Create a new file {original_path}_fileterd.{original_extension} with only not chimeric records, format

            support fasta|fastq|mhap|paf

    -e, --extract <extract>...

            Create a new file {original_path}_extracted.{original_extension} with only chimeric records, format support

            fasta|fastq|mhap|paf

    -s, --split <split>...

            Create a new file {original_path}_splited.{original_extension} where chimeric records are split, format

            support fasta|fastq

    -F, --format <format>                                  Force the format used [possible values: paf, mhap]

    -c, --chimeric-threshold <chimeric-threshold>

            Overlap depth threshold below which a gap should be created [default: 0]

 

    -n, --not-covered-threshold <not-covered-threshold>

            Coverage depth threshold above which a read are marked as not covered [default: 0.80]

 

        --filtered-suffix <filtered-suffix>

            Change the suffix of file generate by filter option [default: _filtered]

 

        --extracted-suffix <extracted-suffix>

            Change the suffix of file generate by extract option [default: _extracted]

 

        --splited-suffix <splited-suffix>

            Change the suffix of file generate by split option [default: _splited]

 

    -C, --compression-out <compression-out>

            Output compression format, the input compression format is chosen by default [possible values: gzip, bzip2,

            lzma, no]

> yacrd scrubbing -h

$ yacrd scrubbing -h

yacrd-scrubbing 

In scrubbing mode yacrd remove all part of read not covered

 

USAGE:

    yacrd scrubbing [FLAGS] [OPTIONS] --mapping <mapping> --report <report> --scrubbed <scrubbed> --sequence <sequence>

 

FLAGS:

    -j, --json       Yacrd report are write in json format

    -h, --help       Prints help information

    -V, --version    Prints version information

 

OPTIONS:

    -m, --mapping <mapping>

            Path to mapping file in PAF or MHAP format (with .paf or .mhap extension, paf format by default)

 

    -s, --sequence <sequence>

            Path to sequence you want scrubbed, format support fasta|fastq

 

    -r, --report <report>                                  Path where yacrd report are writen

    -S, --scrubbed <scrubbed>                              Path where scrubbed read are write [default: -]

    -c, --chimeric-threshold <chimeric-threshold>

            Overlap depth threshold below which a gap should be created [default: 0]

 

    -n, --not-covered-threshold <not-covered-threshold>

            Coverage depth threshold above which a read are marked as not covered [default: 0.80]

 

    -M, --mapping-format <format>                          Force the format used [possible values: paf, mhap]

 

 

実行方法

1、chimeric   カバレッジギャップがリードの中央にある場合にキメラ検出

 

キメラ検出。必要に応じて一括セッティングオプション"-x"を使用する (e.g., "-x ava-ont"))。12スレッド使用。

minimap2 -t 12 reads.fq.gz reads.fq.gz | yacrd chimeric -o reads.yacrd

 

キメラフィルタリングされたリードを出力

minimap2 -t 12 reads.fq.gz reads.fq.gz > mapping.paf
yacrd chimeric -i mapping.paf -f reads.fasta > reads.yacrd # produce reads_fileterd.fasta
  • -f, --filter <filter>    Create a new file {original_path}_fileterd.{original_extension} with only not chimeric records, format support fasta|fastq|mhap|paf

 

キメラリードのみ出力

minimap2 -t 12 reads.fq.gz reads.fq.gz > mapping.paf
yacrd chimeric -i mapping.paf -e reads.fasta > reads.yacrd # produce reads_extracted.fasta
  • -e, --extract <extract>   Create a new file {original_path}_extracted {original_extension} with only chimeric records, format support fasta|fastq|mhap|paf

 

キメラ領域をsplitして除き、全リード出力

minimap2 -t 12 reads.fq.gz reads.fq.gz > mapping.paf
yacrd chimeric -i mapping.paf -s reads.fasta > reads.yacrd # produce reads_splited.fasta
  • -s, --split <split>    Create a new file {original_path}_splited.{original_extension} where chimeric records are split, format support fasta|fastq 

 

2、scrubbing    -  chimericとは違い、カバレッジギャップがリードのどの領域にあってもキメラ検出

ONTのキメラ検出。12スレッド使用。

minimap2 -x ava-ont -g 500 -t 12 reads.fq.gz reads.fq.gz > overlap.paf
yacrd scrubbing -c 4 -n 0.4 -m overlap.paf -s reads.fasta -S reads_scrubbed.fasta -r scrubbed_report.yacrd

minimap2

  • -g    stop chain enlongation if there are no minimizers in INT-bp [5000]

yacrd

  • -c     Overlap depth threshold below which a gap should be created [default: 0]
  • -n     Coverage depth threshold above which a read are marked as not covered [default: 0.80]
  • -m    Path to mapping file in PAF or MHAP format (with .paf or .mhap extension, paf format by default)
  • -s     Path to sequence you want scrubbed, format support fasta|fastq
  • -S     Path where scrubbed read are write [default: -]
  • -r      Path where yacrd report are writen

 

PacbioのP6-C4のキメラ検出。12スレッド使用。

minimap2 -x ava-pb -g 800 -t 12 reads.fq.gz reads.fq.gz > overlap.paf
yacrd scrubbing -c 4 -n 0.4 -m overlap.paf -s reads.fasta -S reads_scrubbed.fasta -r scrubbed_report.yacrd 

 

PacbioのSequelのキメラ検出。12スレッド使用。

minimap2 -x ava-pb -g 800 -t 12 reads.fq.gz reads.fq.gz > overlap.paf
yacrd scrubbing -c 4 -n 0.4 -m overlap.paf -s reads.fasta -S reads_scrubbed.fasta -r scrubbed_report.yacrd 

 

Preprintでは、キメラリードを除くことでヒトゲノムと線虫のRedBean (旧 wtdbg2)とminiasmを使ったde novo assemblyのパフォーマンスが大きく伸びることが示されています。

fpaは別に紹介します。 

 

 

2020 03/30

v0.6からコマンドが大きく変わったので、使い方を簡単に書いておきます。

#ONT-readsのall versus all overlap
minimap2 -x ava-ont -g 500 -t 16 ONT.fq.gz reads.fq.gz > overlap.paf

#キメラがないか分析
yacrd -i overlap.paf -o reads.yacrd

#キメラやオーバーラップなしリードの表記があれば、以下のコマンドでフィルタリング
#その前にリードをfasta変換しておく。
seqkit fq2fa ONT.fq.gz > ONT.fa

#yacrdのサブコマンド;filter実行
yacrd -i overlap.paf -o reads.yacrd filter -i ONT.fa -o reads.filter.fasta

引用

yacrd and fpa: upstream tools for long-read genome assembly

Pierre Marijon, Rayan Chikhi, Jean-Stéphane Varré 

bioRxiv preprint first posted online Jun. 18, 2019

 

yacrd and fpa: upstream tools for long-read genome assembly
Pierre Marijon, Rayan Chikhi, Jean-Stéphane Varré
Bioinformatics, 2020 

 

関連


メタゲノムデータセットをタンパク質レベルでアセンブリし、ホモログサーチを行う GRASP2

 

 メタゲノミクスは、特定の微生物群集のゲノム含有量を研究するための培養に依存しないアプローチである。典型的なメタゲノミクス研究では、環境サンプルから微生物のDNAが抽出され、次世代シークエンシング(NGS)技術を使用してシークエンシングされる。中程度または複雑度の高い微生物群集から得られたメタゲノミクスデータの分析は、構成微生物の分類学的およびゲノム的多様性が高いこと、およびしばしば不十分なおよび/または不均一なシーケンシングカバレッジのために困難である。そのようなデータセットde novoアセンブリ(シーケンスリードから個々の微生物ゲノムを再構築することを目的とする)は、しばしば不完全で断片的であり、下流の機能的および分類学的分析を妨げる。たとえば、MetaHIT(Metagenomics of the Human Intestinal Tract)プロジェクトデータアセンブリでは、コンティグN50はわずか2.2kbであり、リードの53%以上は未アセンブルのままだった[ref.1]。

 アセンブリに依存しない分析方法は、利用可能なデータベースに対してそれらを検索することによって個々のリードに直接アノテーションを付ける。これらのデータベースは、completeシーケンシングされたゲノム、タンパク質およびタンパク質ドメイン、ならびにアノテーション付き分類法を有するマーカー遺伝子を含み得る。データベースに対する有意なヒットは、データベース内のリードと配列の間の相同性を示唆し、個々のリードの機能と分類法を推論し、続いてコミュニティ全体の機能を予測することを可能にする。ただし、このアプローチは既存のデータベースの完全性に大きく依存している。実際には、そのようなデータベースは、ほとんどの環境における微生物の多様性が十分に特徴付けられていないかまたはシーケンシングされていない、単純でよく研究されたコミュニティを除いて、めったに利用できない。この場合、ほとんどのデータベース検索では、中程度または遠隔相同性の検出が行われる。これは、closely relatredの検出と比較してより困難である。そのようなホモログ検索は、限られた量の機能的および分類学的シグナルしか含まないリード長の短さによってさらに混乱する。

 本著者らは、失われた機能的および分類学的シグナルを再構築するためにリードを重ね合わせそれらを長くすることでこの問題に取り組むことを試みた。リードのオーバーラップ検出は、デノボアセンブリで採用されているものと同様である。しかし、ほとんどのアセンブラのように過度のエラー修正や積極的なグラフフィルタリング、プルーニングは行わなかったため、豊富なゲノムのほとんどの多型を保持することができる。さらに、(FragGeneScan [ref.2]を使用して)ヌクレオチドリードから翻訳されたショートペプチド配列に対してリードオーバーラップ検出を行った。このアプローチは、アミノ酸空間における同義置換の崩壊により、タンパク質配列をより効果的に再構築することが示されている[ref.3]。そのような直感に基づいて、著者らは以前に概念的な重なりグラフにおいてクエリ/リファレンスタンパク質のホモログに対応するパスウェイを見つけることを目的としたGRAPS(ガイド付き参照ベースの短いペプチドのアセンブリ)と呼ばれる同時アラインメントおよびアセンブリアルゴリズムを開発した。。この意味で、GRASPはリファレンスのホモログタンパク質配列を再構築する遺伝子中心アセンブリツールとして使用することができる。それはまた、アセンブリされたコンティグをテンプレートとして使用することによって、個々の相同なリードを動員するためのホモログ検索ツールとしても使用され得る。ベンチマークデータは、GRASPが既存のホモログ検索ツールの再現率を精度を落とすことなく約40%から約60%に向上させたことを示した[ref.4]。 RAPSearch [ref.5]やDIAMOND [ref.6]のような最近開発されたホモログ検索ツールは、計算効率の向上に焦点を当てており、再現率のパフォーマンスは低い。したがって、効率的な同族体検索ツールは多数利用可能だが、GRASPは依然として高い再現率と全体的なパフォーマンスを誇るツールである。

 しかし、GRASPで行われるde novoアセンブリのため、BLASTのような伝統的なリードベースのホモログ検索ツールよりはるかに遅くなる。さらに、GRASPはまた、アセンブリ工程を補助するために使用される接尾辞配列データ構造を記憶するために大量のメモリを必要とする。したがって、GRASPは、人間の口腔環境から生成されたデータなど、比較的単純で小さなメタゲノムデータセットにのみ適用されていた。 GRASPの実用的な実用性を拡張するために、本著者らはGRASPxと呼ばれる新しいアライメントとアセンブリの同時アルゴリズムを開発した[ref.7]。 GRASPxは、GRASPの本来のパフォーマンスを犠牲にすることなく、GRASPの速度を最大30倍向上させる。 GRASPxもGRASPと同量のメモリを使用する。パフォーマンスが大幅に向上したにもかかわらず、GRASPxは他のホモログ検索ツールよりも多くの時間とスペースを必要としていた。

 この論文では、GRASPxの計算効率とスペース効率をさらに向上させるために、GRASP2と呼ばれる、完全に再設計された同時アライメントとアセンブリアルゴリズムを紹介する。 GRASPx、BLAST、PSI-BLAST、およびFASTMを含むホモログ検索アルゴリズムのセットに対してGRASP2をベンチマークした。ベンチマークの結果は、GRASP2がGRASPxよりも8倍高速であり、インデックス作成フェーズとアセンブリ/検索フェーズの両方で8倍メモリ使用量が少ないことを示した。 GRASP2はGRASPxと同じ高い性能を持っている。そしてGRASP2は他のホモログ検索ツールよりもかなり優れている。これらの結果は、本発明者らの新しいアルゴリズムが、同時アラインメントおよびアセンブリアルゴリズムの実行時間および所要スペースを効果的に減らすことができることを示唆している。結果として得られるソフトウェアGRASP2は、その高い性能と著しく改善された計算効率のために大きな応用の可能性を持っている。

 最初にGRASPxのアルゴリズムフレームワークを要約して、ここでさらに改善された制限を特定する。 GRASPxは、2つの主な段階、すなわち効率的なオーバーラップ検出を容易にするためにリードセット全体から接尾辞配列を構築するためのindexing段階、およびリファレンスタンパク質のホモログに対応するパスウェイについて概念的重なりグラフを検索するアセンブリ/検索段階を含む。(以下略)

 

 

インストール

ubuntu16.04でテストした(docker使用、ホストOS macos10.12)。

依存

  • GCC version 2.8.0+
  • Boost library version 1.46.0+
  • cmake
  • make
#手っ取り早く全部入れる
apt update && apt install libboost-all-dev cmake make gcc

#Fraggenescanなどシーケンシングリードからproteinを予測する方法も必要
#bioconda (link)
conda install -c bioconda -y fraggenescan

GRASP2本体 SourceForge 

ターボールをダウンロードして解凍、ビルドする。

tar -zxvf GRASP2_V0.02.tar
cd GRASP2_release/
cmake CMakeLists.txt
make -j 16
cd bin/

> ./grasp2-assemble --help

# ./grasp2-assemble --help

Usage: grasp-assemble [query] [peptide_db] [contig_out]

 

List of options:

  --help                 print the help message

  --query arg            query sequences (in FASTA format)

  --db arg               short-peptide reads (in FASTA format)

  --contig_out arg       assembled contigs output (in FASTA format)

  --alignment_out arg    alignment output (between query and contigs, in 

                         BLAST-style format)

  --index arg (=index)   working directory where the indexing files should be 

                         found

  --gap_open arg (=-10)  gap open penalty for sequence alignment

  --gap_extend arg (=-1) gap extension penalty for sequence alignment

  --band_size arg (=40)  band size for sequence alignment

  --e_value arg (=10)    e-value cutoff (Karlin-Altschul statistics)

  --orphan arg           align orphan reads (true/false; default false)

  --num_threads arg (=1) maximum number of threads to be used

  --verbose arg          print intermediate information (true/false; default 

                         true)

 

> ./grasp2-build --help

# ./grasp2-build --help

Usage: grasp-build [peptide_db (FASTA)]

 

List of options:

  --help                     print the help message

  --db_file arg              short-peptide reads (in FASTA format)

  --index arg (=index)       working directory for indexing file dump

  --extension_len arg (=10)  minimum overlap length for path extension

  --neighbor_score arg (=11) neighbor score for 3-mer seed matches

  --num_threads arg (=1)     maximum number of threads to be used

  --verbose arg              print intermediate information (default true)

 

GRASPxも入れる。

同様にbuildする。

cd GRASPx_src/
#Makefileの2行目をカレントパスに修正してからmake
make -j 16

 

実行方法

1、FragGeneScanを使いgenomeやショートリードからアミノ酸配列を得る。ここではアセンブリ配列contig.fastaを指定。

run_FragGeneScan.pl -genome=contigs.fasta -out=output -complete=0 -train=illumina_5 -thread=20
  • -genome=        sequence file name including the full path
  • -complete=    1 if the sequence file has complete genomic sequences. 0 if the sequence file has short sequence reads
  • train=     file name that contains model parameters. [complete] for complete genomic sequences or short sequence reads without sequencing error.  [sanger_5] for Sanger sequencing reads with about 0.5% error rate. [sanger_10] for Sanger sequencing reads with about 1% error rate. [454_10] for 454 pyrosequencing reads with about 1% error rate. [454_30] for 454 pyrosequencing reads with about 3% error rate. [illumina_5] for Illumina sequencing reads with about 0.5% error rate. [illumina_10] for Illumina sequencing reads with about 1% error rate.

output.faやgffファイル等が出力される。

 

2、Indexing 

bin/grasp2-build Short-peptides.fa

 

3、アセンブリ

bin/grasp2-assemble [Query_proteins] [Short-peptide_reads] [Contig_output]

 

4、リードをアセンブリしたcontigにマッピング

graps-map [Short-peptide_reads] [Contig_output] [Map_output]

テスト中。

 

引用

GRASP2: fast and memory-efficient gene-centric assembly and homolog search for metagenomic sequencing data

Cuncong Zhong, Youngik Yang, Shibu Yooseph
BMC Bioinformatics 2019 20 (Suppl 11) :276