macでインフォマティクス

macでインフォマティクス

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

Minhashをメタゲノム解析へ応用する CMash

 

 Minhashは、2つの集合の類似性をJaccard指数(集合の和に対する交点の大きさの比として定義される)の観点から推定する確率的な手法である。この手法は、対象となる集合の大きさが似ている場合に最も優れた性能を発揮し、集合の大きさが大きく異なる場合には性能が大幅に低下することを実証した。本論文では、大きさの異なる集合のJaccard indexを推定するのに適した、「containment minhash approach」と呼ばれる新しい効率的な手法を紹介する。これは、高速なメンバーシップクエリのための別の確率的手法(特にBloomフィルタ)を利用することで実現している。containment Minhash法の推定誤差の確率を導出し、古典的なMinhash法を大幅に改善することを示す。また、時間的・空間的な複雑さの点でも大きな改善が見られる。応用例として、この手法を用いてメタゲノムデータセット中の生物の有無を検出し、非常に小さく、存在量の少ない微生物の存在を検出できることを示す。

 

インストール

Github

mamba install -c bioconda cmash -y

> MakeDNADatabase.py  -h

usage: MakeDNADatabase.py [-h] [-p PRIME] [-t THREADS] [-n NUM_HASHES]

                          [-k K_SIZE] [-i]

                          in_file out_file

 

This script creates training/reference sketches for each FASTA/Q file listed

in the input file.

 

positional arguments:

  in_file               Input file: file containing (absolute) file names of

                        training genomes.

  out_file              Output training database/reference file (in HDF5

                        format)

 

optional arguments:

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

  -p PRIME, --prime PRIME

                        Prime (for modding hashes) (default: 9999999999971)

  -t THREADS, --threads THREADS

                        Number of threads to use (default: 128)

  -n NUM_HASHES, --num_hashes NUM_HASHES

                        Number of hashes to use. (default: 500)

  -k K_SIZE, --k_size K_SIZE

                        K-mer size (default: 21)

  -i, --intersect_nodegraph

                        Optional flag to export Nodegraph file (bloom filter)

                        containing all k-mers in the training database. Saved

                        in same location as out_file. This is to be used with

                        QueryDNADatabase.py (default: False)

~                                                                                                                                                                                             

~                                                                                                                                                                                             

~

> MakeNodeGraph.py -h

usage: MakeNodeGraph.py [-h] [-fp FP_RATE] [-i INTERSECT_NODEGRAPH]

                        [-k K_SIZE] [-t THREADS]

                        in_file out_dir

 

This script will create node graph for a given k-mer size and query file (can

be used as input to QueryDNADatabase.py)

 

positional arguments:

  in_file               Input file: FASTQ/A file (can be gzipped).

  out_dir               Output directory

 

optional arguments:

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

  -fp FP_RATE, --fp_rate FP_RATE

                        False positive rate. (default: 0.0001)

  -i INTERSECT_NODEGRAPH, --intersect_nodegraph INTERSECT_NODEGRAPH

                        Location of Node Graph. Will only insert query k-mers

                        in bloom filter if they appear anywhere in the

                        training database. Note that the Jaccard estimates

                        will now be J(query intersect union_i training_i,

                        training_i) instead of J(query, training_i), but will

                        use significantly less space (unfortunately will also

                        disable threading). (default: None)

  -k K_SIZE, --k_size K_SIZE

                        K-mer size (default: 21)

  -t THREADS, --threads THREADS

                        Number of threads to use (default: 128)

>  QueryDNADatabase.py -h

usage: QueryDNADatabase.py [-h] [-t THREADS] [-f] [-fp FP_RATE]

                           [-ct CONTAINMENT_THRESHOLD] [-c CONFIDENCE]

                           [-ng NODE_GRAPH] [-b] [-i]

                           in_file training_data out_csv

 

This script creates a CSV file of similarity indicies between the input file

and each of the sketches in the training/reference file.

 

positional arguments:

  in_file               Input file: FASTQ/A file (can be gzipped).

  training_data         Training/reference data (HDF5 file created by

                        MakeTrainingDatabase.py)

  out_csv               Output CSV file

 

optional arguments:

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

  -t THREADS, --threads THREADS

                        Number of threads to use (default: 128)

  -f, --force           Force creation of new NodeGraph. (default: False)

  -fp FP_RATE, --fp_rate FP_RATE

                        False positive rate. (default: 0.0001)

  -ct CONTAINMENT_THRESHOLD, --containment_threshold CONTAINMENT_THRESHOLD

                        Only return results with containment index above this

                        value (default: 0.02)

  -c CONFIDENCE, --confidence CONFIDENCE

                        Desired probability that all results were returned

                        with containment index above threshold [-ct] (default:

                        0.95)

  -ng NODE_GRAPH, --node_graph NODE_GRAPH

                        NodeGraph/bloom filter location. Used if it exists; if

                        not, one will be created and put in the same directory

                        as the specified output CSV file. (default: None)

  -b, --base_name       Flag to indicate that only the base names (not the

                        full path) should be saved in the output CSV file

                        (default: False)

  -i, --intersect_nodegraph

                        Option to only insert query k-mers in bloom filter if

                        they appear anywhere in the training database. Note

                        that the Jaccard estimates will now be J(query

                        intersect union_i training_i, training_i) instead of

                        J(query, training_i), but will use significantly less

                        space. (default: False)

 

 

実行方法

リファレンス/トレーニングデータベースを作成し、ブルームフィルターを形成し、データベースに問い合わせるという流れで使う。

1、fasta/fastqのフルパスを指定した1行ずつ指定したファイルを指定してデータベースを作成。

MakeDNADatabase.py list.txt TrainingDatabase.h5

 

2、(任意)ブルームフィルタを作成

MakeNodeGraph.py TrainingDatabase.h5 outdir

 

3、クエリファイルMetagenome.faとデータベースとのJaccard indexの推定値を得る

QueryDNADatabase.py Metagenome.fa TrainingDatabase.h5 output.csv

 

引用

Improving MinHash via the containment index with applications to metagenomic analysis
David Koslickia, Hooman Zabetib

https://doi.org/10.1016/j.amc.2019.02.018

 

関連


 

包括的な遺伝子セットのエンリッチメント解析ウェブサーバー Enrichr

 

 エンリッチメント解析は、ゲノムワイド実験で得られた遺伝子セットを解析するための一般的な手法である。ここでは、Enrichrと呼ばれるこの分野のツールの1つを大幅に更新した。Enrichrには、現在、解析やダウンロードが可能な多様な遺伝子セットライブラリの大規模なコレクションが含まれている。Enrichrは現在、102の遺伝子セットライブラリから180 184のアノテーションされた遺伝子セットを含んでいる。Enrichrには、ファジーセットの送信機能、BEDファイルのアップロード、アプリケーション・プログラミング・インターフェースの改善、クラスターグラムとしての結果の視覚化などの新機能が追加されている。Enrichrは、キュレーションされた遺伝子セットの包括的なリソースであり、さらなる生物学的発見のために生物学的知識を蓄積する検索エンジンでもある。Enrichrは次のサイトで自由に利用できる: http://amp.pharm.mssm.edu/Enrichr.

 

help

https://maayanlab.cloud/Enrichr/help#basics

 

webサービス

https://maayanlab.cloud/Enrichr/にアクセスする。ヒトとマウス(M. musculus)に対応している。

f:id:kazumaxneo:20211212005959p:plain

D. melanogasterC. elegansS. cerevisiaeD. rerioにも対応している。

https://maayanlab.cloud/modEnrichr/

f:id:kazumaxneo:20211213085937p:plain

 

Entrez遺伝子シンボルのセットをテキストボックスにペーストする。各遺伝子シンボルは1行ずつ入力する。ここではExampleのSTAT3データを指定した(変動遺伝子TOP100 )。

f:id:kazumaxneo:20211213090640p:plain

 

 

 

大規模なデータセットの発現行列と比較して共発現している遺伝子のエンリッチメント解析を行う。ここではARCHS4(ヒトとマウスの公開RNA-seqデータの大部分を均一な条件で解析して公開している)を指定。

f:id:kazumaxneo:20211213090743p:plain

submitをクリック。

出力例

結果ページでは、分析結果が充実したカテゴリーに分けて表示される。左端はTranscriptionのタブになっている。

f:id:kazumaxneo:20211213204154p:plain

結果ページにある各エンリッチメント・タームは、その遺伝子セット・ライブラリによって構成されている。これらのライブラリは、発表された研究や主要な生物学的および生物医学的オンラインデータベースなどの多くの情報源から構築されている。Enrichrのために作成され、Enrichrでしか利用できないものもある。

 

Transcriptionのタブでは、転写の調節に関わる遺伝子セットライブラリのエンリッチメント解析がリストアップされている。

f:id:kazumaxneo:20211213204658p:plain

 

クリックすると展開される。ちなみに、モバイル端末では左右のスワイプで表示内容を切り替えることができる。

Bar Graph

f:id:kazumaxneo:20211213204954p:plain

棒グラフの長さは、その特定の遺伝子セットや用語の重要性を表している。色の明るさも同じで、色が明るいほどその用語の重要性が高いことを示している。初期設定では、すべてのフィギュアの配色が赤になっているが、歯車アイコンをクリックするとカラーピッカーが表示され、配色を変更することができる(棒グラフ、Clustergramのグリッド、ネットワークに対応)。

 

Table

f:id:kazumaxneo:20211213205013p:plain

Enrichrでは、エンリッチメント結果を報告するために4つのスコア:p-value、q-value、rank (Z-score)、およびcombined scoreを実装しているが(詳細はbackgroundについて)、データテーブルではこれらのデータの一部を見ることができる。特定の用語を検索して結果をフィルタリングすることもできる。

 

Clustergram

f:id:kazumaxneo:20211213204942p:plain

グリッドの各マスは用語を表し、他の用語との遺伝子セットの類似性に基づいて配置されている。マス目が明るいほどその用語の重要度が高いことを意味している。。グリッドをクリックすると、隣接する用語との相関スコアに基づいてグリッドを着色し、有意な用語を白いドットで表す別の表示に切り替えることができる。

 

Enrichrでは、結果を他の人と簡単に共有することができるようになっている。それには、左のメニューボックスの共有アイコンをクリックする。これにより、分析結果への一時的なダイレクトリンクが提供され、他の人と共有することができる。

f:id:kazumaxneo:20211213211454p:plain

 

    Appyter

f:id:kazumaxneo:20211213205139p:plain

 

Enrichr Appyterを使用すると、エンリッチメント解析の結果をすぐに公開できるようなビジュアライゼーションを作成できる。

Appyter使用例

f:id:kazumaxneo:20211213205337p:plain

f:id:kazumaxneo:20211213205434p:plain

f:id:kazumaxneo:20211213205401p:plain

f:id:kazumaxneo:20211213205415p:plain

(Appytersカタログは、生物医学研究者がデータを様々な方法で分析し、視覚化することを可能にするソフトウェアプラットフォーム)

 

左から2番目のタブはPathwayとなっている。KEGGPanther、Wikipathway、Reactomeなどの結果が表示される。

f:id:kazumaxneo:20211213205713p:plain

KEGG Human

f:id:kazumaxneo:20211213205731p:plain

 

左から3番目のタブはOntologies。GOなどがある。

f:id:kazumaxneo:20211213205803p:plain

4番目のタブはDiseases。ClinVarやCOVID-19などがある。

f:id:kazumaxneo:20211213205948p:plain

 

5番目のタブはcell type。

f:id:kazumaxneo:20211213210537p:plain

 

6番目のタブはmisc (多種類)。

f:id:kazumaxneo:20211213210636p:plain

 

7番目のタブはlegacy。KEGG2016などがある。

f:id:kazumaxneo:20211213210805p:plain

 

右端のタブはCrowd(雑多、混ぜ物)。

f:id:kazumaxneo:20211213210953p:plain

 

Enrichrでは、ユーザーアカウントを登録してログインすることもできる。アカウントを持つことで、結果をEnrichrサーバーに永久保存することができ、再アップロードの必要なくセットに戻ることができる。使用するには、ログイン後、セットを保存または共有ボタンをクリックする。これにより、自動的にアカウントと関連付けられる(helpより)。

引用
Enrichr: a comprehensive gene set enrichment analysis web server 2016 update
Maxim V. Kuleshov, Matthew R. Jones, Andrew D. Rouillard, Nicolas F. Fernandez, Qiaonan Duan, Zichen Wang, Simon Koplev, Sherry L. Jenkins, Kathleen M. Jagodnik, Alexander Lachmann, Michael G. McDermott, Caroline D. Monteiro, Gregory W. Gundersen, Avi Ma'ayan

Nucleic Acids Res. 2016 Jul 8; 44(Web Server issue): W90–W97

 

Enrichr: interactive and collaborative HTML5 gene list enrichment analysis tool
Edward Y Chen, Christopher M Tan, Yan Kou, Qiaonan Duan, Zichen Wang, Gabriela Vaz Meirelles, Neil R Clark, Avi Ma'ayan

BMC Bioinformatics. 2013 Apr 15;14:128

シーケンスアラインメントやHMMER3のHMMプロファイルをlogoで視覚化する skylign

 

 ロゴは、分子生物学において、配列の保存パターンをコンパクトなグラフで表現するためによく用いられる。ロゴは、配列アラインメントや隠れマルコフモデルに含まれる情報を、各位置に文字のスタックを描くことで表現する。スタックの高さはその位置の保存度に対応し、スタック内の各文字の高さはその位置でのその文字の頻度に依存する。

 本著者らはSkylignと呼ばれる新しいツールとウェブサーバーを発表する。このツールは、配列アラインメントとHMMプロファイルの両方のロゴを作成するための統一されたフレームワークを提供する。静的な画像ファイルに加えて、SkylignはWebページに含めるための新しいインタラクティブなロゴプロットを作成する。これらのインタラクティブなロゴは、スクロール、ズーム、および基本的な値の検査が可能になっている。Skylignは、冗長な配列をダウンウェイトし、観測されたカウントをinformed priors(事前承認)と組み合わせることにより、配列アライメントにおけるサンプリングバイアスを回避することができる。また、ギャップパラメータの表現を簡素化し、オプションとして、位置の保存性の代替計算に基づいて文字の高さをスケーリングすることができる。

 Skylignは、ウェブサイト、RESTfulインターフェイスを持つスクリプト可能なウェブサービス、およびダウンロード可能なソフトウェアパッケージとして提供されている。Skylignのインタラクティブなロゴは、数行のHTMLマークアップで簡単にWebページに組み込むことができる。Skylignは、http://skylign.orgで利用できる。

 

Help

http://skylign.org/help/

 

webサービス

http://skylign.org/にアクセスする。

f:id:kazumaxneo:20211211223359p:plain

 

HMMモデルを指定する。ここではKEGGのKofamKOALAのKEGG IDマッピング時に作成されたHMMモデルを使った。

f:id:kazumaxneo:20211211223638p:plain

 

出力例

f:id:kazumaxneo:20211211223751p:plain

 

 

クリックするとその位置のアミノ酸の推定出現確率が表示される。

f:id:kazumaxneo:20211211224055p:plain

HMMを提供した場合、出現確率はHMMから直接抽出される。アラインメントファイルを提供した場合、Probabilityはスタート時に選べる3つの計算方法のいずれかにより推定される。

 

酵素タンパク質であれば、活性部位のアミノ酸残基などに注目して見ていく。ロゴがほぼ見えない部位は、保存されていないアミノ酸残基(ループ構造を形成する部位など)を意味する。結果はPNGSVGなどのフォーマットでダウンロードできる。

引用

Skylign: a tool for creating informative, interactive logos representing sequence alignments and profile hidden Markov models
Travis J Wheeler, Jody Clements & Robert D Finn 
BMC Bioinformatics volume 15, Article number: 7 (2014)

 

ショートリードとロングリードによりトランスクリプトームアセンブリの構造回復とアバンダンス推定を行う StringTieの新しいバージョン

 

 トランスクリプトームのアセンブリには、short-read RNA sequencingとlong-read RNA sequencingのそれぞれに長所と短所がある。ショートリードは精度が高い反面、複数のエクソンにまたがることができない。Long-read技術は、完全な長さの転写産物を捉えることができるが、エラー率が高く、スプライスサイトを誤認することが多く、スループットが低いため、定量化が困難である。ここでは、ハイブリッドリードアセンブリを実行するStringTieの新リリースを紹介する。ロングリードとショートリードの両方の強みを活用することで、StringTieによるハイブリッドリードアセンブリは、ロングリードのみ、またはショートリードのみのアセンブリよりも精度が高く、いくつかのデータセットでは、ロングリードデータのみのアセンブリよりも大幅に高い精度を得ながら、正しくアセンブリされた転写産物の数を2倍以上に増やすことができる。ここでは、シロイヌナズナ、ムスキュラス、ヒトのシミュレーションデータおよび実データを用いて、精度の向上を実証した。また、ハイブリッド・リード・アセンブリは、アセンブリ前にロングリードを修正するよりも精度が高く、かつ大幅に高速であることを示す。StringTieは、オープンソースソフトウェアとして、https://github.com/gpertea/stringtie から自由に入手できる。

 

HP

https://ccb.jhu.edu/software/stringtie/index.shtml

 

 

インストール

Github

git clone https://github.com/gpertea/stringtie
cd stringtie
make release

> ./stringtie 

stringtie <in.bam ..> [-G <guide_gff>] [-l <prefix>] [-o <out.gtf>] [-p <cpus>]

 [-v] [-a <min_anchor_len>] [-m <min_len>] [-j <min_anchor_cov>] [-f <min_iso>]

 [-c <min_bundle_cov>] [-g <bdist>] [-u] [-L] [-e] [--viral] [-E <err_margin>]

 [--ptf <f_tab>] [-x <seqid,..>] [-A <gene_abund.out>] [-h] {-B|-b <dir_path>}

 [--mix] [--conservative] [--rf] [--fr]

Assemble RNA-Seq alignments into potential transcripts.

Options:

 --version : print just the version at stdout and exit

 --conservative : conservative transcript assembly, same as -t -c 1.5 -f 0.05

 --mix : both short and long read data alignments are provided

        (long read alignments must be the 2nd BAM/CRAM input file)

 --rf : assume stranded library fr-firststrand

 --fr : assume stranded library fr-secondstrand

 -G reference annotation to use for guiding the assembly process (GTF/GFF)

 --ptf : load point-features from a given 4 column feature file <f_tab>

 -o output path/file name for the assembled transcripts GTF (default: stdout)

 -l name prefix for output transcripts (default: STRG)

 -f minimum isoform fraction (default: 0.01)

 -L long reads processing; also enforces -s 1.5 -g 0 (default:false)

 -R if long reads are provided, just clean and collapse the reads but

    do not assemble

 -m minimum assembled transcript length (default: 200)

 -a minimum anchor length for junctions (default: 10)

 -j minimum junction coverage (default: 1)

 -t disable trimming of predicted transcripts based on coverage

    (default: coverage trimming is enabled)

 -c minimum reads per bp coverage to consider for multi-exon transcript

    (default: 1)

 -s minimum reads per bp coverage to consider for single-exon transcript

    (default: 4.75)

 -v verbose (log bundle processing details)

 -g maximum gap allowed between read mappings (default: 50)

 -M fraction of bundle allowed to be covered by multi-hit reads (default:1)

 -p number of threads (CPUs) to use (default: 1)

 -A gene abundance estimation output file

 -E define window around possibly erroneous splice sites from long reads to

    look out for correct splice sites (default: 25)

 -B enable output of Ballgown table files which will be created in the

    same directory as the output GTF (requires -G, -o recommended)

 -b enable output of Ballgown table files but these files will be 

    created under the directory path given as <dir_path>

 -e only estimate the abundance of given reference transcripts (requires -G)

 --viral : only relevant for long reads from viral data where splice sites

    do not follow consensus (default:false)

 -x do not assemble any transcripts on the given reference sequence(s)

 -u no multi-mapping correction (default: correction enabled)

 -h print this usage message and exit

 --ref/--cram-ref reference genome FASTA file for CRAM input

 

Transcript merge usage mode: 

  stringtie --merge [Options] { gtf_list | strg1.gtf ...}

With this option StringTie will assemble transcripts from multiple

input files generating a unified non-redundant set of isoforms. In this mode

the following options are availabl

  -G <guide_gff>   reference annotation to include in the merging (GTF/GFF3)

  -o <out_gtf>     output file name for the merged transcripts GTF

                    (default: stdout)

  -m <min_len>     minimum input transcript length to include in the merge

                    (default: 50)

  -c <min_cov>     minimum input transcript coverage to include in the merge

                    (default: 0)

  -F <min_fpkm>    minimum input transcript FPKM to include in the merge

                    (default: 1.0)

  -T <min_tpm>     minimum input transcript TPM to include in the merge

                    (default: 1.0)

  -f <min_iso>     minimum isoform fraction (default: 0.01)

  -g <gap_len>     gap between transcripts to merge together (default: 250)

  -i               keep merged transcripts with retained introns; by default

                   these are not kept unless there is strong evidence for them

  -l <label>       name prefix for output transcripts (default: MSTRG)

 

Error: no input file provided!

 

 

 

実行方法

ショートリードRNAseqとロングリードRNAseqのbamファイルを指定する。

stringtie --mix -o mix_reads.out.gtf mix_short.bam mix_long.bam

 

bamファイルは座標ソートされている必要があります。Githubで確認して下さい。

引用

Improved Transcriptome Assembly Using a Hybrid of Long and Short Reads with StringTie
Alaina Shumate, Brandon Wong, Geo Pertea,  Mihaela Pertea

bioRxiv, Posted December 10, 2021

 

関連


生命科学実験の検索・提案のためのウェブアプリケーション LEXAS

2021 12/11 誤字修正

 

 細胞生物学では,研究者は関連する論文を読み,記述されている実験や結果を検討することでウェットな実験を計画する。今日、研究者は実験を計画するために長い時間をかけて文献を調査している。

 実験計画を加速するために、本著者らはLEXAS (Life-Science EXperiment seArch and Suggestion)というウェブアプリケーションを開発した。LEXASは、生物医学実験の記述をキュレーションし、次に行うことができる遺伝子に関する実験を提案する。LEXASを開発するために、まずPubMed Central(PMC)に保存されている生物医学論文のフルテキストから実験の記述を検索した。これらの検索された実験と生物医学のナリッジベースやデータベースを用いて、次の実験を提案する機械学習モデルを学習した。このモデルは、合理的な遺伝子だけでなく、関心のある遺伝子といくつかの重要な特徴を共有する新規の遺伝子も次の実験のターゲットとして提案することができる。

 LEXASはhttps://lexas.f.u-tokyo.ac.jp/にて、ユーザーに検索と提案の2つのインターフェースを提供している。検索インターフェースでは、ユーザーは実験の説明の包括的なリストを見つけることができ、提案インターフェースでは、ユーザーは可能な実験方法とともに分析可能な遺伝子のリストを見つけることができる。ソースコードhttps://github.com/lexas-f-utokyo/lexas にある。

 

Documentation

https://lexas.f.u-tokyo.ac.jp/documentation.html

 

webサービス

LEXAS - Lifescience experiment search and suggestion -にアクセスする。

f:id:kazumaxneo:20211210094856p:plain

サーチをクリック。

 

用語や遺伝子名を入力する。例のTP53を入力した。

f:id:kazumaxneo:20211210095045p:plain

 

クリックすると類似する候補が出る。

f:id:kazumaxneo:20211210095142p:plain

カテゴリを選ぶ。ここではBioinformaticsを選ぶ。

f:id:kazumaxneo:20211210095325p:plain

Goをクリック。

 

検索結果。PMCID(PMCにアーカイブされている文献の識別子)とヒットした単語を含む文が表示される。初期は信頼度スコア順にソートされている。

f:id:kazumaxneo:20211210095416p:plain



PMCIDをクリックするとPMCの該当する文献にジャンプする。

f:id:kazumaxneo:20211210095726p:plain

 

LEXASのもう1つの機能であるSuggestも同じようにして使用できる。キーワードを入力する(ヒットがないこともある)。

f:id:kazumaxneo:20211210095939p:plain

右側のウィンドウでは、提案に使用する2つの機械学習モデルを選択する。1つは信頼性の高いモデルで、BioGRIDなどの7つの主要な生物医学データベースと、Gene Ontologyなどの4つの知識ベースを使用している。このモデルは、公表された知識体系に沿った信頼性の高い提案を求める人に適している。対照的に、もう一つの探索的モデルは、生物医学データベースのみを使用して構築されている。このモデルは、遺伝子の比較的客観的な特徴に基づいて、新規の接続を探している人に適している(論文より)。

 

出力例

f:id:kazumaxneo:20211210100026p:plain

Suggestでは、ある遺伝子の実験後に分析可能な遺伝子のリストを見つけることができる。

 

引用

LEXAS: a web application for life science experiment search and suggestion
Kei K Ito, Yoshimasa Tsuruoka, Daiju Kitagawa

bioRxiv, Posted December 07, 2021

 

関連


オックスフォードナノポアシークエンシングリードのトリミングツール ProwlerTrimmer

 

 トリミングおよびフィルタリングツールは、配列アラインメントの精度を高め、結果の信頼性を向上させるため、DNAシーケンス解析において有用である。オックスフォード・ナノポア・テクノロジー(ONT)のトリミングおよびフィルタリングツールは、現在のところ初歩的なもので、一般的にはリード全体の平均品質に基づいてリードをフィルタリングしているだけである。これでは、高品質な配列の領域を含むリードが廃棄されてしまう。ここでは、ショートリードデータのトリミングに使用されるアルゴリズムにヒントを得て、ウィンドウベースのアプローチを使用するトリマー「Prowler」を提案する。重要な事は、オプションでトリミングされたセクションをNsに置き換えることで、位相とリードの長さの情報を保持することである。
 Prowlerを哺乳類とバクテリアのデータセットに適用し、アラインメントとアセンブリそれぞれに対する効果を評価した。Nanofiltでフィルタリングしたデータと比較して、Prowlerでトリミングしたデータのアラインメントはエラー率が低く、マッピングされたリードの数も多かった。Prowlerでトリミングされたデータのアセンブリは、Nanofiltでフィルタリングされたデータよりもエラー率が低かったが、アセンブリの連続性が少し(some)犠牲になっていた。

 

 

インストール

Github

git clone https://github.com/ProwlerForNanopore/ProwlerTrimmer.git
cd ProwlerTrimmer/

python TrimmerLarge.py -h

$ python TrimmerLarge.py -h

Usage: TrimmerLarge.py [options]

 

Options:

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

  -f FILE, --file=FILE  the name of the file you want to trim, wihtout the

                        folderpath

  -i DIRECTORY, --infolder=DIRECTORY

                        the folderpath where your file to be trimmed is

                        located (default = cwd)

  -o DIRECTORY, --outfolder=DIRECTORY

                        the folderpath where your want to save the trimmed

                        file (default = cwd)

  -w INT, --windowsize=INT

                        change the size of the trimming window (default=

                        100bp)

  -l INT, --minlen=INT  change the minimum acceptable numer of bases in a read

                        (default=100)

  -m [S/D], --trimmode=[S/D]

                        select trimming algorithm: S for static  or D for

                        dynamic (default=)

  -q INT, --qscore=INT  select the phred quality score trimming threshold

                        (default=7)

  -d INT, --datamax=INT

                        select a maximum data subsample in MB (default=0,

                        entire file))

  -r [.fasta/.fastq], --outformat=[.fasta/.fastq]

                        select output format of trimmed file (fastq or fasta)

                        (default=.fastq)

  -c [L/T/LT], --clip=[L/T/LT]

                        select L to clip leading Ns, T to trim trialing Ns and

                        LT to trim both (default=LT)

  -g [U0/F1/F2...], --fragments=[U0/F1/F2...]

                        select fragmentation mode (default=U0)

 

 

実行方法

fastqファイルを指定するか、fastqファイルを含むディレクトリを指定する。

python TrimmerLarge.py -i ONT.fq -w 1000 -l 1000 -o outprefix
  • -f     the name of the file you want to trim, wihtout the folderpath
  • -r [.fasta | .fastq]     select output format of trimmed file (fastq or fasta)   (default=.fastq)

  •  -i     the folderpath where your file to be trimmed is located (default = cwd)

  • -w    change the size of the trimming window (default= 100bp)

  • -l    change the minimum acceptable numer of bases in a read (default=100) 

     

 

引用

Prowler: A novel trimming algorithm for Oxford Nanopore sequence data
Simon Lee,  Loan T Nguyen,  Ben J Hayes,  Elizabeth Ross
Bioinformatics, Published: 02 September 2021

 

関連


ゲノムアセンブリと遺伝地図を統合するツール Chromonomer

 

 新しいリファレンスゲノムの配列決定とコンピュータによるアセンブリのペースは加速している。しかし、DNAシーケンシング技術やアセンブルソフトウェアツールは進化し続けているが、反復配列などのゲノムの生物学的特徴や、シーケンシングライブラリの調製に伴う分子アーティファクトは、断片的なアセンブルやキメラ状のアセンブルを引き起こす可能性がある。このような欠陥を放置しておくと、ゲノムの構造や機能の理解が進まないばかりか、最悪の場合、研究を大きく誤らせることにもなりかねない。幸いなことに、遺伝地図(特にRADseqから得られるマーカー密度の高い地図)や、近縁種から得られる保存されたオルソログ遺伝子の順序など、追加の独立した情報を統合することで、リンクされていない無秩序な断片をつなぎ合わせたり、間違って結合されたリファレンスゲノムを再構築したりすることができる。これらのプロセスを自動化するためのツールセットを紹介する。このツールは、アセンブリや遺伝地図へのあらゆる変更を追跡し、ウェブベースのグラフィカルな視覚化を用いてユーザーがこれらの変更を精査できるようになっている。Chromonomerは、ユーザーが定義したリファレンスゲノム、遺伝子マーカーのマップ、およびオプションとして保存されたシンテニー情報を用いて、染色体モデルの改良されたリファレンスゲノム、すなわち「chromonome」を構築する。Chromonomerの性能を、特性や品質レベルが異なるゲノムアセンブリや遺伝地図で実証する。

 

manual

http://catchenlab.life.illinois.edu/chromonomer/manual/#intro

マニュアルより

Chromonomerは、ゲノムアセンブリと遺伝地図を統合するために設計されたプログラムです。Chromonomerは、ローカルアセンブリの順序と比較して、遺伝地図上で順序が狂っているマーカーを識別して削除したり、遺伝地図に従って正しく組み立てられていないscaffoldsを識別して、そのscaffoldsを分割したりすることに努力します。

 

インストール

HPからダウンロードしてビルドする。

http://catchenlab.life.illinois.edu/chromonomer/

f:id:kazumaxneo:20210725002113p:plain

 

cd chromonomer_x.xx
./configure --prefix=/usr/local/bin
make -j
sudo make install

> ./chromonomer

$ ./chromonomer 

You must specify a path to the genetic map linkage group definition file.

chromonomer 1.13

chromonomer --map map_path --alns bam_path --agp agp_path --out_path out_path [--verbose] [--depth path] [--gtf path --orth_gtf path --orthologs path] [-h]

--map <path>: TSV file containing the genetic map linkage group definitions.

--alns <path>: SAM or BAM formatted alignments linking map markers to genome contigs.

--agp <path>: AGP formatted file defining scaffold layout including any assembly gaps.

--fasta <path>: optional, supply a scaffold-based FASTA that will be translated according to the integrated assembly.

 

Virtual Breakpoint Options:

--depth <path>: supply a Samtools-formatted, TSV depth-of-coverage file for use in inferring scaffold breakpoints.

--depth_stdevs <num>: number of standard deviations from the mean to call a coverage window a breakpoint (default 5 stdevs).

--depth_win_size <num>: size of sliding window to use to call mean depth of coverage (default 10Kb).

 

Ordering by conserved synteny options:

--gtf <path>: supply a GTF file describing the location of genes in the scaffold-level assembly.

--gff <path>: supply a GFF file describing the location of genes in the scaffold-level assembly (supply either GTF or GFF, not both).

--orth_gtf <path>: supply a GTF file describing the location of orthologous genes in a related assembly.

--orth_gff <path>: supply a GFF file describing the location of orthologous genes in a related assembly (supply either GTF or GFF, not both).

--orthologs <path>: supply a two-column, TSV file describing the orthologous relationships between pairs of genes contained in the two GTF files.

 

Genome Correction Options:

--rescaffold: Re-scaffold contiguous sequence based on the genetic map (will break contigs).

 

Output Options:

--out_path: path to write chromonomer output.

--verbose: turn on detailed console output for each scaffold.

--scaffold_prefix <string>: text string to use as a common naming prefix when creating new scaffolds (default 'CHRR').

--description <string>: supply a description of this chromonomer execution that will be written into the output log.

--join_gap_size <num>: number of Ns to insert in between scaffold joins (default: 100).

 

Operational Options:

--min_markers <num>: minimum number of markers required to anchor a scaffold on a particular

map node. Will not remove scaffolds only present on one node (default 2).

--allpaths_agp: The AllPaths-LG assembler produces a non-standard AGP file. This flag will change

the AGP parser to accommodate it.

--disable_splitting: do not split scaffolds. When markers indicate a conflict in where to place a scffold, do

not split the scaffold, instead place it in the integration where the majority of markers are.

 

Filtering Options:

--scaffold_wl <path>: only process scaffolds in this file.

 

HTML/PHP Options:

--data_version <string>: used to create versioned PHP files so multiple runs can be

executed and compared in the web interface.

--html_prefix <string>: any text that should be prepended onto links in the PHP output

to accommodate where they are located in the web server document root.

h,--help: display this help message.

 

 

ApacheウェブサーバでChromonomerウェブインタフェースを有効にするため、以下の内容のconfigファイル;chromonomer.confを作成。

<Directory "/usr/local/share/chromonomer/php">
Order deny,allow
Deny from all
Allow from all
Require all granted
</Directory>

Alias /chromonomer "/usr/local/share/chromonomer/php"

 

このconfigファイルを/etc/apache2/conf-available/に置く(apacheが入ってないなら"sudo apt install apache2")。さらにapache2/conf-enabled/にシンボリックリンクを張ってapacheを再起動。

sudo cp chromonomer.conf /etc/apache2/conf-available/
sudo ln -s /etc/apache2/conf-available/chromonomer.conf /etc/apache2/conf-enabled/
sudo apachectl restart

 

 

実行手順

http://catchenlab.life.illinois.edu/chromonomer/manual/#exec

1、新しく作成したリファレンスゲノムに対してマーカーをアラインメントする。

アラインメントを行うには、各マーカーの配列を含むFASTAファイルを用意する。各マーカーのIDは、遺伝地図を記述するために用意したマーカーのIDと一致していなければならない。Stacks(HP)を使用してマーカーを生成した場合、UNIXコマンドを使用してStacksカタログから各マーカーのコンセンサス配列をエクスポートできる。マーカー遺伝子のFASTAファイルが用意できたら、BWAなどのアライナーを使い、ゲノムにマッピングしてBAMファイルを作る。

 

2、遺伝地図を記述したタブ区切りのリスト、つまり地図中のマーカーとその連鎖グループおよびcentiMorganの位置を記したファイルが必要になる。Chromonomerは、Stacks( http://catchenlab.life.illinois.edu/stacks )を使用してRADデータから構築された遺伝地図を扱うように設計されているが、マーカーがゲノムにアラインメントされている遺伝地図であれば、どのようなものでも動作する。

 

3、Chromonomerの出力を格納するディレクトリを作成する。例えば、20150603。Webインターフェイスの下にも同じ名前のディレクトリを作成する。

 

4、Chromonomerを実行し、適切な入出力パスと--data_versionフラグで作成したディレクトリを指定する。例えば以下のようになる。

chromonomer -p ~/research/20150603_linkage_map.tsv \
    -o ~/research/20150603/ -s ~/research/markers.sam \
    -a ~/research/final.assembly.agp --data_version 20150603

 

引用

Chromonomer: A Tool Set for Repairing and Enhancing Assembled Genomes Through Integration of Genetic Maps and Conserved Synteny
Julian Catchen, Angel Amores, Susan Bassham

G3 (Bethesda). 2020 Nov 5;10(11):4115-4128

 

Chromonomer: a tool set for repairing and enhancing assembled genomes through integration of genetic maps and conserved synteny
Julian Catchen, Angel Amores, Susan Bassham

bioRxiv, Posted February 05, 2020

 

関連