macでインフォマティクス

macでインフォマティクス

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

アセンブリ過程でロングリードをフィルタリングする fpa

 

以前、ロングリードのアセンブリ前処理ツール yacrdを紹介した。

今回はアセンブリ過程でフィルタリングして出力を調節するfpaを紹介する。

 

以下のフィルタリングが行える (Githubより)。

  • internal match
  • containment
  • dovetails
  • self matching
  • read name match against regex
  • length of overlap
  • length of read in overlap

 

Pierre MarijonさんのBlog (More details and some experiment are present)

https://blog.pierre.marijon.fr/binary-mapping-format/

 

インストール

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

ビルド依存

  • Rust >= 1.32
  • libgz
  • libbzip2
  • liblzma

その他 (option)

  • minimap2
  • miniasm

本体 Github

#bioconda (link) (*1)
conda install -c bioconda -y fpa

> fpa -h

# fpa -h

fpa 0.5 Krabby

Pierre Marijon <pierre.marijon@inria.fr>

fpa take long read mapping information and filter them

 

USAGE:

    fpa [OPTIONS] [SUBCOMMAND]

 

FLAGS:

    -h, --help       Prints help information

    -V, --version    Prints version information

 

OPTIONS:

    -z, --compression-out <compression-out>

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

            lzma, no]

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

    -i, --input <input>                                    Path to input file, use '-' for stdin [default: -]

        --internal-threshold <internal-match-threshold>

            A match is internal match if overhang length > match length * internal threshold this option set internal

            match [default: 0.8]

    -o, --output <output>                                  Path to output file, use '-' for stdout [default: -]

 

SUBCOMMANDS:

    drop      fpa drop mapping match this constraints

    gfa       fpa generate a overlap graph in gfa1 format with mapping passing filter

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

    index     fpa generate a index of mapping passing filter

    keep      fpa keep only mapping match this constraints

    rename    fpa rename reads with name you chose or with incremental counter

> fpa keep -h

# fpa keep -h

fpa-keep

fpa keep only mapping match this constraints

 

USAGE:

    fpa keep [FLAGS] [OPTIONS] [SUBCOMMAND]

 

FLAGS:

    -c, --containment      Keep only containment mapping

    -d, --dovetail         Keep only dovetail mapping

    -h, --help             Prints help information

    -i, --internalmatch    Keep only internal mapping

    -m, --same-name        Keep only mapping where reads have same name

    -V, --version          Prints version information

 

OPTIONS:

    -l, --length-lower <length_lower>                      Keep only mapping with length lower than value

    -L, --length-upper <length_upper>                      Keep only mapping with length upper than value

    -n, --name-match <name_match>                          Keep only mapping where one reads match with regex

    -s, --sequence-length-lower <sequence_length_lower>

            Keep only mapping where one reads have length lower than value

 

    -S, --sequence-length-upper <sequence_length_upper>

            Keep only mapping where one reads have length upper than value

fpa drop -h

# fpa drop -h

fpa-drop 

fpa drop mapping match this constraints

 

USAGE:

    fpa drop [FLAGS] [OPTIONS] [SUBCOMMAND]

 

FLAGS:

    -c, --containment      Drop containment mapping

    -d, --dovetail         Drop dovetail mapping

    -h, --help             Prints help information

    -i, --internalmatch    Drop internal mapping

    -m, --same-name        Drop mapping where reads have same name

    -V, --version          Prints version information

 

OPTIONS:

    -l, --length-lower <length_lower>                      Drop mapping with length lower than value

    -L, --length-upper <length_upper>                      Drop mapping with length upper than value

    -n, --name-match <name_match>                          Drop mapping where one reads match with regex

    -s, --sequence-length-lower <sequence_length_lower>    Drop mapping where one reads have length lower than value

    -S, --sequence-length-upper <sequence_length_upper>    Drop mapping where one reads have length upper than value

fpa gfa -h

# fpa gfa -h

fpa-gfa 

fpa generate a overlap graph in gfa1 format with mapping passing filter

 

USAGE:

    fpa gfa [FLAGS] [OPTIONS] [SUBCOMMAND]

 

FLAGS:

    -c, --containment      Keep containment overlap

    -h, --help             Prints help information

    -i, --internalmatch    Keep internal match overlap

    -V, --version          Prints version information

 

OPTIONS:

    -o, --output <output>    Write mapping passing filter in gfa1 graph format in path passed as parameter

fpa index -h

# fpa index -h

fpa-index 

fpa generate a index of mapping passing filter

 

USAGE:

    fpa index [OPTIONS] [SUBCOMMAND]

 

FLAGS:

    -h, --help       Prints help information

    -V, --version    Prints version information

 

OPTIONS:

    -f, --filename <filename>    Write index of mapping passing filter in path passed as parameter

    -t, --type <type>            Type of index, only reference read when it's query, target or both of them [default:

                                 both]  [possible values: query, target, both]

> fpa rename -h

# fpa rename -h

fpa-rename 

fpa rename reads with name you chose or with incremental counter

 

USAGE:

    fpa rename [OPTIONS] [SUBCOMMAND]

 

FLAGS:

    -h, --help       Prints help information

    -V, --version    Prints version information

 

OPTIONS:

    -i, --input <input>      Rename reads with value in path passed as parameter

    -o, --output <output>    Write rename table in path passed as parameter

 

実行方法

1、 keep - 制約に合致するマッピングを出力

2、 drop - 制約に合致するマッピングを除去

minimap/miniasmのロングリードアセンブリの途中にfpaのフィルタリングを加えるように設計されている。そのため、fpaはpafファイルを入力とし、pafファイルを出力する。miniasmによるRaw assemblyの第一ステップはminimapまたはminimap2によるリード同士のオーバーラップ(all-vs-all read self-mappings)なので、この第一ステップの出力(minimapの出力はdefaultではPAF)を受け取り、fpaでフィルタリングしてminiasmでアセンブリ(コンカテネートしてunitigをGFA出力(エラーレートはinputのリードと同じ))を行う。

例えばONTのrawアセンブリを以下の流れで実行するとする。
minimap2 -t 18 -x ava-ont nanopore.fq.gz nanopore.fq.gz | gzip -1 > reads.paf.gz
miniasm -f nanopore.fq.gz reads.paf.gz > raw_assebly.gfa


#この流れにfpaを組み込む。"fpa drop -l <NUM>"で指定以下のリードを捨てて出力。
#1000bp以上のリードだけpaf出力。
minimap2 -t 18 -x ava-ont nanopore.fq.gz nanopore.fq.gz | fpa drop -l 1000 | gzip - > filtered.paf.gz
#それからminiasmを実行する

#fpa dropなら逆の条件でフィルタリングできる。fpa keep -l <NUM>で指定以下のリードのみキープして出力
#10000bp以下のリードだけpaf出力。
minimap2 -t 18 -x ava-ont nanopore.fq.gz nanopore.fq.gz | fpa keep -l 10000 | gzip - > filtered.paf.gz

#つまりkeepとdropは逆の挙動をする。1000bp以上のリードだけキープして出力するなら"keep -L <NUM>"か"drop -l <NUM>"。keepを使うなら
minimap2 -t 18 -x ava-ont nanopore.fq.gz nanopore.fq.gz | fpa keep -L 1000 | gzip - > filtered.paf.gz

 

サイズセレクト

 1000bp以上10000-bp以下のリードのみにフィルタリング。

#keepを使うなら
minimap2 -t 18 -x ava-ont nanopore.fq.gz nanopore.fq.gz | fpa keep -L 1000 -l 10000| gzip - > filtered.paf.gz

#dropを使うなら
minimap2 -t 18 -x ava-ont nanopore.fq.gz nanopore.fq.gz | fpa drop -l 1000 -L 10000| gzip - > filtered.paf.gz

 

dovetail(蟻継ぎ)のオーバーラップを持つリードのみにフィルタリング。

minimap2 long_read.fasta long_read.fasta | fpa keep -d | gzip - > filtered.paf.gz
  • -d     Keep only dovetail mapping
       

インターナルのオーバーラップを持つリードのみにフィルタリング。

minimap2 long_read.fasta long_read.fasta | fpa keep -i | gzip - > filtered.paf.gz
  • -i     Keep only internal mapping 

 

汚染配列を除去。

minimap2 long_read.fasta long_read.fasta | fpa keep -c | gzip - > filtered.paf.gz
  • -c   Keep only containment mapping

 

汚染配列除去とサイズセレクト。

minimap2 long_read.fasta long_read.fasta | fpa keep -c keep -L 1000 | gzip - > filtered.paf.gz
  • -c   Keep only containment mapping

 

 

他の例はGithubを参考にしてください。 

 

fpa gfaコマンドも組み込めば、アセンブリグラフのGFA(version1)で出力できる。bandage等を使えばグラフを可視化できる。

minimap2 long_read.fasta long_read.fasta | fpa keep -d gfa -o dovetail.gfa > dovetail.paf

あくまでアセンブリ前提に設計されているので、GFAの出力(-o)以外にリダイレクトでPAFも出力する。fpa gfa -oで出力したGFAに、配列情報は含まれないので注意する。

 

 

fpaはロングリードのアセンブリを前提に設計されているので、アセンブリ前にロングリードを前処理したいならyacrdを使う。

引用

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

 


 

*1

2019 7/4現在、macosだと0.4までしか入らない。v0.4だとkeepなどの一部コマンドが存在しない。 keepなどを使いたければソースからビルドする。