macでインフォマティクス

macでインフォマティクス

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

ロングリードを使ってcontigのscaffoldingを行う LINKS

2015 GigaScience Nanopore long read assembly Pacbio Linked reads

2019 6/1 help追加、インストール追記

 

　ロングリードシークエンシング技術はここ数年で急速に成熟し、ゲノムアセンブリへのロングリードの利点は欠かせないものになった[論文より　ref.1]。最近、複数グループがエラーの多いロングリードから完全なバクテリアゲノムへのデノボアセンブリが可能であることを示している[ref.2-4]。マイクロ流体力学、エレクトロニクス、ナノポア技術の急速な進歩により、ポータブルなロングリードシーケンシング技術が私たちの目の前にある[ref.5]。 （一部略）

　本発明者らは、E.coliのONTデータセットとして公開されたデータ[ref.6 pubmed6]を使い、R7ケミストリーデータ（Full 2Dと呼ばれる）の特性を評価した。著者らは、これらのリードが平均して77.1±10.6％のsequence identityを有することを観察した（大腸菌K-12 MG1655に対して50％以上のsequence identityを有する1714のリード11Mbpにおいて）。間違いなくこの全体的なクオリティは低いにもかかわらず、完成したゲノムと比較して、リードにおいて正確なk塩基の連続したストレッチが頻繁に存在する。これらのストレッチは特異性を持つには十分な長さだが、エラーフリーになるにはきわめて短い（k = 15、論文図1）。著者らはこの特性を利用して、LINKSを開発した。LINKSは、ONTリードからペアのkmersを抽出し、それらを用いてコンティグ対を連結する（図2）。

　(以下略)

 

f:id:kazumaxneo:20181016070732p:plain

 LINKS algorithm. 論文より転載。

 

HP

http://www.bcgsc.ca/platform/bioinfo/software/links

 

インストール

ubuntu18.04のminiconda3-4.0.5でテストした。

#anacondaを使っているならcondaで導入できる
conda install -c bioconda -y links

#/usr/bin/timeを使うので入れておく
apt install time

> LINKS --help

usage（Githubより）

Usage: ./LINKS [v1.8.6]
-f  sequences to scaffold (Multi-FASTA format with each sequence on a single line, required)
-s  file-of-filenames, full path to long sequence reads or MPET pairs [see below] (Multi-FASTA/fastq format, required)
-m  MPET reads (default -m 1 = yes, default = no, optional
	! DO NOT SET IF NOT USING MPET. WHEN SET, LINKS WILL EXPECT A SPECIAL FORMAT UNDER -s
	! Paired MPET reads in their original outward orientation <- -> must be separated by ":"
	  >template_name
	  ACGACACTATGCATAAGCAGACGAGCAGCGACGCAGCACG:ATATATAGCGCACGACGCAGCACAGCAGCAGACGAC
-d  distance between k-mer pairs (ie. target distances to re-scaffold on. default -d 4000, optional)
	Multiple distances are separated by comma. eg. -d 500,1000,2000,3000
-k  k-mer value (default -k 15, optional)
-t  step of sliding window when extracting k-mer pairs from long reads
(default -t 2, optional)
	Multiple steps are separated by comma. eg. -t 10,5
-o  offset position for extracting k-mer pairs (default -o 0, optional)
-e  error (%) allowed on -d distance   e.g. -e 0.1  == distance +/- 10%
(default -e 0.1, optional)
-l  minimum number of links (k-mer pairs) to compute scaffold (default -l 5, optional) 
-a  maximum link ratio between two best contig pairs (default -a 0.3, optional)
	 *higher values lead to least accurate scaffolding*
-z  minimum contig length to consider for scaffolding (default -z 500, optional)
-b  base name for your output files (optional)
-r  Bloom filter input file for sequences supplied in -s (optional, if none provided will output to .bloom)
	 NOTE: BLOOM FILTER MUST BE DERIVED FROM THE SAME FILE SUPPLIED IN -f WITH SAME -k VALUE
	 IF YOU DO NOT SUPPLY A BLOOM FILTER, ONE WILL BE CREATED (.bloom)
-p  Bloom filter false positive rate (default -p 0.001, optional; increase to prevent memory allocation errors)
-x  Turn off Bloom filter functionality (-x 1 = yes, default = no, optional)
-v  Runs in verbose mode (-v 1 = yes, default = no, optional)


Notes:
-s K12_F2D.fof specifies a file-of-filenames (text file) listing: K12_full2dONT_longread.fa (see ./test folder)
   -f and -s : sequences must be on a SINGLE line with no linebreaks
   eg.  
   >LONGREAD-1
   AATACAATAGACGCACA...ATGAACGCAGACTTACAG
   >LONGREAD-2
   TGTGCTCTCTGTAATGTTC...ATACAGAACACGCAGCCAAGCGA

-x When turned on (-x 1), LINKS will run with a behaviour equivalent to v1.3 (no Bloom filters).  
This may be useful for large genome assembly drafts and when long reads are extremely high quality.

> LINKS

$ LINKS

Usage: /Users/user/.pyenv/versions/anaconda2-4.2.0/bin/LINKS [v1.5.2]

-f  sequences to scaffold (Multi-FASTA format, required)

-s  file-of-filenames, full path to long sequence reads (Multi-FASTA/fastq format, required)

-d  distance between k-mer pairs (eg. target distance to re-scaffold on. default -d 4000, optional)

-k  k-mer value (default -k 15, optional)

-t  step of sliding window when extracting k-mer pairs from long reads (default -t 2, optional)

-o  offset position for extracting k-mer pairs (default -o 0, optional)

-e  error (%) allowed on -d distance   e.g. -e 0.1  == distance +/- 10% (default -e 0.1, optional)

-l  minimum number of links (k-mer pairs) to compute scaffold (default -l 5, optional)

-a  maximum link ratio between two best contig pairs (default -a 0.3, optional)

    *higher values lead to least accurate scaffolding*

-b  base name for your output files (optional)

-r  Bloom filter input file for sequences supplied in -s (optional, if none provided will output to .bloom)

    NOTE: BLOOM FILTER MUST BE DERIVED FROM THE SAME FILE SUPPLIED IN -f WITH SAME -k VALUE

    IF YOU DO NOT SUPPLY A BLOOM FILTER, ONE WILL BE CREATED (.bloom)

-p  Bloom filter false positive rate (default -p 0.0001, optional; increase to prevent memory allocation errors)

-x  Turn off Bloom filter functionality (-x 1 = yes, default = no, optional)

-v  Runs in verbose mode (-v 1 = yes, default = no, optional)

 

テストラン

git clone https://github.com/bcgsc/LINKS.git
cd LINKS/links_v1.8.6/test/
./runall.sh

シーケンスデータをダウンロードして実行します。かなり時間がかかります。またメモリを130GBほど使います。

 

実行方法

LINKS -f contig.fa -s K12_F2D.fof -b links_output
  •  -f   sequences to scaffold (Multi-FASTA format with each sequence on a single line, required)
  • -s   file-of-filenames, full path to long sequence reads or MPET pairs [see below] (Multi-FASTA/fastq format, required)
  • -b   base name for your output files (optional)

  

引用

LINKS: Scalable, alignment-free scaffolding of draft genomes with long reads

Warren RL, Yang C, Vandervalk BP, Behsaz B, Lagman A, Jones SJ, Birol I

Gigascience. 2015 Aug 4;4:35