macでインフォマティクス

macでインフォマティクス

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

contigやシーケンシングリードのリファレンスへのアラインメントを複数の方法で視覚化する Alvis

2019 6/10 誤字修正

2019 6/21 リンク追加

 

 2セットの配列間のアラインメントを見つけることは、バイオインフォマティクスにおける基本的な作業である。ロングリードの解析、アセンブリ結果の評価、またはターゲットキャプチャープロトコルの評価では、リファレンスゲノムまたは遺伝子セットに対するアライメントが必要になる。アラインメントを計算するためのさまざまなツールが存在し、これらのツールはさまざまな出力形式を生成する。一般的なフォーマットのほとんどは、直観的な人間の理解を伝えるのではなく、コンピューターによる解析を容易にするために設計された大きなタブ区切りリストで構成されている。それでも、染色体全体のコンティグのレイアウト調査、染色体や遺伝子全体のリードカバレッジデプスの描写など、多くの種類の分析でアライメントの視覚的表現が役立つ。そのような分析はまた、ゲノムの切断部分の人工的結合によって形成されたキメラの存在を明らかにすることができる。
 今日まで、アラインメントの可視化に利用可能な多くのツールがあり、それらはユーザーが直面しているタスクに応じてそれらの特定の機能が異なる。 1つの可能性はensembl [ref.1]やIGV [ref.2]のようなゲノムブラウザを使うことである。ゲノムブラウザは、単なる配列アラインメントの視覚化以上のものを実行し、ユーザーがゲノム配列とアノテーションデータを閲覧、検索、分析することを可能にする。ゲノムブラウザはユーザにグラフィカルインタフェースを提供し、その多くはWebブラウザを通じてオンラインで利用できる。これは、コマンドラインに慣れていないユーザに役立つ。しかしながら、ゲノムブラウザで可能な広範囲のタスクはそれらが提供するインタフェースにおいていくらかの複雑さを必要とし、それらが生成する画像は通常publication品質のものではない。
ツールCircos [ref.3]は、publication用、特に単一ゲノムに関連する多くのデータセットを表示するため(例:[ref.4])、およびゲノム間の構造変化のための高品質画像を作成するための一般的なツールである。ただし、これらの図は表面的には魅力的に見えるが、使い過ぎると明瞭さに欠ける可能性があり、通常は個々のコンティグの配置を視覚化するための縮尺が間違っている。
 MUMmerパッケージ[ref.5]に含まれているmummerplotツールは、2セットのシーケンス間のアラインメントのドットプロットを作成することができ、リアレンジメントの容易な評価を可能にする。これらのプロットはシンプルで、明瞭で読みやすく、そしてプログラムそれ自体は使いやすいシンプルなコマンドラインインターフェースを持っている。しかしながら、それらはユーザがMUMmerパッケージに含まれるアラインメントソフトウェアを使用することを必要とし、そして多くのクエリーまたは標的配列がある場合には読みにくくなる可能性がある。
 もう1つの選択肢は、BLASTのWebインターフェースによって提供されます[ref.6]。これにより、ユーザーはBLASTを使用してNCBI ntデータベースに対して照会されるシーケンスを入力できる。ベストヒットへのアラインメントが図に表示される(図1aと同様)。しかしながら、この方法は、多数の配列を分析したり、一度にいくつかのアラインメントを比較するには不適切である。
 ここでは、さまざまな入力ファイルの種類から4種類のpublication品質の図を生成できる、ロングリードとアセンブリの配置を視覚化するためのツールであるAlvisを紹介する。 Alvisは、ユーザーによるアラインメントの柔軟なフィルタリングを可能にし、潜在的なキメラリードまたはコンティグを自動的に強調表示する。さらに、ユーザは、潜在的なキメラ配列のみを見て、そのようなすべての配列のリストおよびおおよその結合位置を含むテキストファイルを入手することを選択してもよい。このオプションをロングリードシーケンスデータに使用すると、アセンブリの正確性と連続性が向上することを示す。

 AlvisはJavaで書かれており、Java Runtime Environmentを備えたどのプラットフォームでも実行できる。

 Alvisは次のフォーマットの入力を受け付ける:BLASTテーブル形式、BWA [ref.7]および他のアライナからのSAMファイル、minimap2 [ref.8]からのPAFファイル、MUMmerの.coordsおよび.tilingファイル[ref.5]、BLAT [ref.9]からのPSLファイル。さらに、このソフトウェアは他のフォーマットを簡単に実装できるように設計されている。ダイアグラムはSVGまたはTeXフォーマットで出力できる。後者の場合、PDFファイルを生成するためにLaTeXコンパイラが必要である。
 ユーザーは、ノイズを除去するためにアラインメントをフィルタリングすることも選べる。これは、デフォルトではクエリ長の1%未満であるアライメントを破棄することによって達成される(この値はユーザーによって変更できる)。ユーザはまた、キメラリードから生じた可能性があるアラインメントのみを表示する選択もできる。(一部略)各アライメントの長さは、クエリ配列の長さの少なくとも10%でなければならない。これらのデフォルト値はユーザーが調整できる。 Alvisは、結合のおおよその位置とともに、潜在的にキメラリードまたは連続のIDを含む別のプレーンテキストファイルを出力することもできる。
 Alvisはディストリビューションに含まれているJarファイルを使ってコマンドラインから実行できる。たとえば、minimap2ファイルからカバレッジマップを作成するには、次のようにする。
Java -jar Alvis.jar -inputfmt paf -outputfmt tex -type coverageMap -coverageType long -in alignments.paf -outdir output -out outprefix

  

manual

Welcome to Alvis’s documentation! — Alvis documentation

 

インストール

ubuntu16.04のpython2.7.14環境でテストした(docker使用、ホストOS macos10.14)。

依存

本体 Github

git clone https://github.com/SR-Martin/alvis.git
cd alvis/dist/

java -jar Alvis.jar -h

$ java -jar Alvis.jar -h

Alvis

 

Alvis - ALignment VISualisation tool

Syntax: Java -jar .../Alvis.jar [options]

 

Options:

    -type alignment|contigalignment|coveragemap|genomecoverage  Type of diagram to draw.

    -inputfmt psl|coords|tiling|sam|paf|blast  Format of alignment file to use.

    -outputfmt tex|svg  Format for output.

    -in <filename>  Filename for alignment file.

    -outdir <directory>  Output directory - where the diagrams go.

    -out <prefix>  Prefix for output file.

    -tsizes <int>  The sizes of target contigs. To be used when inputfmt is sam and values are not provided  in the header of the SAM file.

    -blastfmt <format string>  The format of blast file. To be used when inputfmt is blast and format of  blast file is user-defined.

    -alignmentQueryName <query sequence>  Name of query sequence for detailed contig alignment diagram.

    -alignmentTargetName <target sequence>  Name of target sequence for detailed contig alignment diagram.

    -coverageType square|long  Type of coverage map to draw for coveragemap diagram.

    -binSize <int>  Size of bin for genomecoverage diagrams.

    -filter  Filter small alignments.

    -minAlignmentPC  Minimum size of alignments, as a percent of the query length, to keep when filtering.

    -chimeras  Only draw potentially chimeric contigs in contig alignment diagram.

    -printChimeras  Output a file containing coordinate information for chimeric contigs.

 

 

テストラン

入ファイルと出力ディレクトリのほか、入力と出力のフォーマットを指定して実行する。

cd alvis/
mkdir output

java -jar dist/Alvis.jar -type alignment \
-inputfmt paf -outputfmt svg -outdir output -out prefix \
-in tutorial_data/tutorial_data/alignments_sample.paf

-type <name>では、4種類の視覚化方式のいずれかを指定する(詳細)。

  • alignment
  • contigAlignment
  • coveragemap
  • genomecoverage

入力と出力に指定できるフォーマットは以下の通り(documentより)。出力はsvgとtaxのみになる。

f:id:kazumaxneo:20190609001051p:plain

 サブコマンドによっても利用できるフォーマットはわずかに異なるものがある。

 

 

実行方法

4つ順番に説明する。いずれの方式でも、contigやリードは前もってリファレンスにアラインメントしてcoordsやpafを得ておく必要がある。

 

1、alignment

alignmentはリファレンス(ゲノム、ターゲット遺伝子等)とターゲット配列(contigsやロングリード)のアラインメントに使う。contigやlong readは前もってminimap2、mummer、blast等でリファレンスにアラインメントし、PAF、coords、blast table出力しておく(SAMはcoveragemapとgenomecoverageのみ対応)。

java -jar dist/Alvis.jar -type alignment \
-inputfmt paf -outputfmt svg -outdir output -out prefix \
-in contig_alignment.paf

上の赤いバーがリファレンス、青のboxがcontigになる。右端にcontig名が並んでいる。

f:id:kazumaxneo:20190609014206p:plain

実際に試してみればわかるが、alignmentコマンドでゲノムとロングリードのアラインメントを行ってしまうと、リードが多いため大量にsvgファイルが出力される*1。この方式は(closely relatedなfinishした)ゲノムとアセンブルして得た長いcontigとのアラインメントに向いているように思う。

 

"-filter"や"-chimera"をつけることで、出力は調整できます(下の例参照)。

 

 

2、contigAlignment

alignmentとは異なり、contig中心のアラインメント表示になる。

java -jar dist/Alvis.jar -type contigAlignment \
-inputfmt paf -outputfmt svg -outdir output -out prefix \
-in contig_alignment.paf

下の図のように、各contigが上から並び、リファレンスゲノムのどのクロモソームやplasidとアラインメントしたか色でわかるようになっている。色は図の一番下に表示されたレジェンドの色に対応している。さらに色はグラデーションになっており、アラインメントされたクロモソームの大まかな位置を色の濃さで判断できる。

f:id:kazumaxneo:20190609124145p:plain

複数ポジションにマップされていると、box内部で複数バーに分けられ、それぞれの色がアサインされる。

 

SVG出力だと、boxをクリックすることでアライメントポジション等の情報を表示できる。

f:id:kazumaxneo:20190609130538p:plain

 

キメラコンティグと判定されると右端に"C"がつく。

f:id:kazumaxneo:20190609125339p:plain

 

オプション "-alignmentQueryName"でcontig名、"-alignmentTargetName"でリファレンスのchromosome名を指定すると、指定したペア配列間だけの詳細アラインメントが出力される。

java -jar dist/Alvis.jar -type contigAlignment \
-inputfmt paf -outputfmt svg -outdir output -out node1_vs_chr \
-alignmentQueryName NODE_1_length_352258_cov_15.530624 \
-alignmentTargetName chr1 \
-in aln.paf

f:id:kazumaxneo:20190609130903p:plain

 


3、coveragemap

Coverage Mapはリファレンスの(リード)マッピングデプスを表すために使用する。各ターゲットコンティグについて、各ピクセルがターゲットコンティグ内の単一位置の範囲を表すヒートマップ画像が生成される。デフォルトでは、行の上から下、左から右に正方形のイメージが生成される。

java -jar dist/Alvis.jar -type coveragemap \
-inputfmt paf -outputfmt svg -outdir output -out node1_vs_chr \
-in aln.paf

f:id:kazumaxneo:20190609195949j:plain

 

4、genomecoverage

java -jar dist/Alvis.jar -type genomecoverage \
-inputfmt paf -outputfmt svg -outdir output -out node1_vs_chr \
-in aln.paf

出力はsvgのみとなる。 

テスト時はエラーが出た。 

 

A.thalianaのナノポアのシーケンシングリードをflyeでアセンブリし、リファレンスにマッピング

#flyeのコマンドは省略

#paf作成
minimap2 -t 16 -x asm5 Arabidopsis_thaliana.TAIR10.dna_sm.toplevel.fa \
flye_assembly/scaffolds.fasta -t 20 > aln.paf

#視覚化1。リファレンスの1%以上の長さのアライメントだけ出力する。
mkdir Arabi
java -jar Alvis.jar -type alignment -inputfmt paf -outputfmt svg \
-outdir Arabi/ -out alignment -in aln.paf \
-filter -minAlignmentPC 1

出力ディレクトリにSVGができる。

f:id:kazumaxneo:20190609190941j:plain

chromosomeスケールの長いアライメントに限定しているので、ファイル数は少ない。

 

svgファイルの1つをブラウザで開いてみた。scaffold_44もscaffold_40もどちらもギャップのある状態でchr1とアラインメントしているため、boxが複数の行に別れて表示されている。

f:id:kazumaxneo:20190609191000j:plain

contigAlignmentコマンドでも視覚化する。

#視覚化2
mkdir Arabi2
java -jar Alvis.jar -type contigAlignment -inputfmt paf -outputfmt svg \
-outdir Arabi2/ -out contigAlignment -in aln.paf \
-filter -minAlignmentPC 1

下の図の一番上がscaffold_44。

f:id:kazumaxneo:20190609195638j:plain

 

ACT(まとめ)でもscaffold_44とchr1のアラインメントを可視化してみる。

python bwast/bwast.py At_chr1.fasta scaffolds_44.fasta -a 

上がchr1、下がscaffolds44。

f:id:kazumaxneo:20190609192627j:plain

ACT(blastn)の方もよく見るとgapがあるが、gapの範囲がscaffoldsと比べてずっと短く目視では分かりにくい。Alvisの方がgapが強調表示されていて、アラインメントのブレイクポイントをつかみやすくなっている。

例2

キメラアラインメントの検出

"-chimera"オプション使用することで、キメラであると考えられるアライメントのみを表示できる。 ここでキメラ配列は、異なるchr、または同じchrの異なる遺伝子座のアラインメントによってクエリの配列の90%以上がカバーされており、また、それぞれのアラインメントがクエリ配列の10%以上をカバーしているものと定義される。値は"-minChimeraCoveragePC"と"-minChimeraAlignmentPC"で変更可能。(*2)

#long readをリファレンスにmappingしてmap作成しておく

mkdir chimera
java -jar Alvis.jar -type alignment -inputfmt paf -outputfmt svg \
-outdir chimera/ -out alignment -in aln.paf \
-chimera -printChimeras

 キメラリードのfastqをスプリットしてfastq出力する。

python dist/chimera_splitter.py -i scaffolds.fasta -c chimeras.txt -o out

Alvis.jarの"-printChimeras"がうまく動作せずchimeras.txtが出力されなかった。修正でき次第追記します。ペーパー中では疑わしいリードを分割し、再アセンブルすることで連続性を改善できることを指摘しています。

 

引用

Alvis: a tool for contig and read ALignment VISualisation and chimera detection

Samuel Martin, Richard Mark Leggett

bioRxiv preprint first posted online Jun. 6, 2019

 


 関連1

アラインメント


関連2

視覚化


 

 

*1

LaTeXを使いtexをpdf変換しても同じ。

 

*2

ペーパー中にあるように、環状DNAの終わりと先頭をまたぐアラインメントもキメラと判定されるので、環状リファレンスとわかっている場合、ラン後にそれらのキメラアラインメントをフィルタイングする必要がある。