macでインフォマティクス

macでインフォマティクス

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

冗長性がある配列データベースに対する正確なリードアライメントを行う KMA

 

 バイオインフォマティクスで最も古く、おそらく最も重要なツールは1つ以上の配列のアライメントである。アライメントは、ある配列が別の配列とどの程度類似しているかを知らせ、類似の配列パターンの存在量を定量化するために使用できる。見つかったパターンにアノテーションがある場合は、これをクエリシーケンスに送ることができる。したがって、アライメント方法を適切に選択することは、さらなる調査にとって重要である[ref.1]。

 過去10年間にわたり、ターゲット配列に対するrawリードの直接マッピングを可能にするいくつかのマッピング方法が開発されてきた。これまでの研究では、リードマッピングベースのアプローチの方が、最初にリードをアセンブルしてBLAST’edする従来の方法よりも優れていることが示されている。 BLASTは未処理のリードを直接マッピングするには遅すぎるので、このような方法が成功するかどうかはアセンブリの品質に大きく左右される。これはゲノムの反復領域および低深度領域内で問題となり、アセンブルされたコンティグ間にしばしばギャップを生じ、データの欠落を招く。今日、Bowtie 2 [ref.4]、BWA-MEM [ref.5]、およびMinimap 2 [ref.6]は、大きなリファレンスゲノムに対する高速なマッピングを可能にする一般的なマッピング方法である。これらの方法の全ては、配列のデータベース全体に対してマッピングするように適合させることができるが、最良の一致が複数ある場合にはランダム選択する問題に悩まされる。この同点問題を解決するために、前処理と後処理でこれらのランダム選択を解決するツールがいくつか開発されている。これらの例はSRST2 [ref.3]とMGmapper [ref.7]で、それぞれBowtie2とBWA-MEMのマッピングの前処理と後処理を行う。 Salmon [ref.8]などの他の方法は、伝統的な方法によって生成されたマッピングから相同テンプレートの配列決定レベルを推定するためにEMアルゴリズムwiki)を使用する。

 しかし、genomic epidemiologyで使用されている冗長なデータベースは、未加工のリードを直接マッピングすることに関して依然として課題を投げかけている。データベースが自然な進化のため新しいシーケンスで絶えず更新されるとき、結果は絶えず変化する状態にあり、クラスタリングを困難にする。データベースのこの機能は、リードが参照シーケンスの固有の部分をカバーするという保証がないため、rawリードの直接マッピングが面倒になる。バクテリアデータベースにマッピングする場合、冗長性が問題になるが、これらのデータベースはいくつかの遺伝子データベースほど冗長ではなく、ペアエンドの片方のペアは単一にマッピングされることが多いため、ペアエンドリードをより有効に利用できる。鋭敏な決定を下すことができるように、正確かつ迅速に冗長性に関する問題を解決するできる新しい方法論が必要である。この要求に応えて、新しいアラインメント方法、KMA(k-merアラインメント)が開発された。 KMAは、homologousなテンプレートを区別するために、新しい分類スキームConClaveを導入した。

 著者らは、冗長データベースに対してrawリードを直接マッピングすることを可能にする新しいアライメント方法KMA、およびスコアリングスキームConClaveを紹介する。 KMAはデータベース内での冗長性を可能にすることで既知のマッパーとは異なり、コンセンサス配列と結果の概要も生成する。 KMAは、現在のユーザーの要求と分析の課題に基づいて、できるだけ直感的でユーザーフレンドリーになるように作成された。 これを解決するために、KMAは5つの主なステップで動作する(trimming of reads、heuristic k-mer mapping、fine alignment、ConClave scoring and reference assembly)(論文図1を参照)。(以下略)

 

インストール

mac os10.14でテストした。

本体 Bitbucket

git clone https://bitbucket.org/genomicepidemiology/kma.git
cd kma && make
#パスの通ったディレクトリに移動
mv kma* /usr/local/bin/

> kma

$ kma

 Too few arguments handed

 Printing help message:

# KMA-1.1.7 mapps raw reads to a template database.

# Options are: Desc: Default: Requirements:

#

# -o Output file None REQUIRED

# -t_db Template DB None REQUIRED

# -i Input file name(s) STDIN

# -ipe Input paired end file name(s)

# -int Input interleaved file name(s)

# -k Kmersize DB defined

# -e evalue 0.05

# -ConClave ConClave version 1

# -mem_mode Use kmers to choose best

# template, and save memory False

# -ex_mode Searh kmers exhaustively False

# -ef Print additional features False

# -vcf Make vcf file, 2 to apply FT False/0

# -deCon Remove contamination False

# -dense Do not allow insertions

# in assembly False

# -ref_fsa Consensus sequnce will

# have "n" instead of gaps False

# -matrix Print assembly matrix False

# -a Print all best mappings False

# -mp Minimum phred score 30

# -5p Cut a constant number of

# nucleotides from the 5 prime. 0

# -Sparse Only count kmers False

# -Mt1 Map only to "num" template. 0 / False

# -ID Minimum ID 1.0%

# -ss Sparse sorting (q,c,d) q

# -pm Pairing method (p,u,f) u

# -fpm Fine Pairing method (p,u,f) u

# -apm Sets both pm and fpm u

# -shm Use shared DB made by kma_shm 0 (lvl)

# -swap Swap DB to disk 0 (lvl)

# -1t1 Skip HMM False

# -ck Count kmers instead of

# pseudo alignment False

# -ca Make circular alignments False

# -boot Bootstrap sequence False

# -bc Base calls should be

# significantly overrepresented. [True]

# -bc90 Base calls should be both

# significantly overrepresented,

# and have 90% agreement. False

# -bcNano Call bases at suspicious

# deletions, made for nanopore. False

# -and Both mrs and p_value thresholds

# has to reached to in order to

# report a template hit. or

# -mq Minimum mapping quality 0

# -mrs Minimum alignment score,

# normalized to alignment length 0.50

# -reward Score for match 1

# -penalty Penalty for mismatch -2

# -gapopen Penalty for gap opening -3

# -gapextend Penalty for gap extension -1

# -per Reward for pairing reads 7

# -cge Set CGE penalties and rewards False

# -t Number of threads 1

# -v Version

# -h Shows this help message

#

> kma_index

$ kma_index 

# Too few arguments handed.

# kma_index creates the databases needed to run KMA, from a list of fasta files given.

# Options are: Desc: Default:

#

# -i Input/query file name (STDIN: "--") None

# -o Output file Input/template file

# -batch Batch input file

# -deCon File with contamination (STDIN: "--") None/False

# -batchD Batch decon file

# -t_db Add to existing DB None/False

# -k Kmersize 16

# -k_t Kmersize for template identification 16

# -k_i Kmersize for indexing 16

# -ML Minimum length of templates kmersize (16)

# -CS Start Chain size 1 M

# -ME Mega DB False

# -NI Do not dump *.index.b False

# -Sparse Make Sparse DB ('-' for no prefix) None/False

# -ht Homology template 1.0

# -hq Homology query 1.0

# -and Both homolgy thresholds

# has to be reached or

# -v Version

# -h Shows this help message

#

> kma_shm

$ kma_shm 

# Too few arguments handed

# Printing help message:

# kma_shm sets up a shared database (sysV) for mapping with KMA.

# Options are: Desc: Default: Requirements:

#

# -t_db Template DB None REQUIRED

# -destroy Destroy shared DB False

# -shmLvl Level of shared memory 1

# -shm-h Explain shm levels

# -v Version

# -h Shows this help message

#

 

 

実行方法

1、index

kma_index -i templates.fsa.gz -o templates
  • -i     Input/query file name (STDIN: "--") None
  • -o    Output file Input/template file

出力

f:id:kazumaxneo:20190108162147j:plain

 

2、マッピング

ショートリードマッピング。1つのリードは1つのテンプレートのみマッピングされる。

#ペアエンド
kma -ipe pair*gz -t_db templates -o output -1t1 -t 20

#シングルエンドとペアエンド
kma -i singleEndReads.fq.gz -ipe pair*gz –o output -t_db templates -1t1 -t 20
  • -i     Input file name(s) STDIN 
  • -ipe     Input paired end file name(s)
  • -t_db   Template DB None REQUIRED
  • -1t1     Skip HMM False
  • -t     Number of threads

 

ロングリードマッピング

kma -i nanoporeReads.fq.gz –o output –t_db database/name -mp 20 -mrs 0.0 -bcNano -Mt1 2 -t 20
  • -mrs    Minimum alignment score
  • -bcNano   Call bases at suspicious
  • -Mt1    Map only to "num" template. 0 / False
  • -mp     Minimum phred score 30

 

出力はテンプレートとのアライメントファイル、アライメント結果から得られるコンセンサス配列などになります。詳細はBitbucketの方で確認して下さい。

引用

Rapid and precise alignment of raw reads against redundant databases with KMA
Clausen PTLC, Aarestrup FM, Lund O

BMC Bioinformatics. 2018 Aug 29;19(1):307