macでインフォマティクス

macでインフォマティクス

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

nanocorrectでロングリードを修正する

nanocorrectはナノポアリードをpolishする方法論。速度が遅いのが欠点らしく、後継としてnaonpolishが発表されている(リンク)。 

  

インストール

依存

全てbrewで導入できる。

Github

https://github.com/jts/nanocorrect

 

ラン

最初にDALIGNERのデータベースを作る。

nanocorrect-overlap INPUT=reads.fasta NAME=nc

 ncは自分が決めた名前である。

全リードをpolishする。

python nanocorrect.py nc all > corrected.fasta

論文では99.5%リファレンスのE.coliの塩基と合致するconitgが得られたと記載されている。 また複数回して精度を上げるのもありらしい。

 

 

引用

A complete bacterial genome assembled de novo using only nanopore sequencing data

Nicholas J Loman, Joshua Quick & Jared T Simpson

Nature Methods 12, 733–735 (2015)

 

SlnCで高速に変異入りリードをシュミレートする

 SlnCは最も多い変異であるSNV、indel、CNVをシミュレートできるNGSのリードシミュレーションツール。マルチコアに対応しており、ARTのようなツールと比較して高速にカバレッジのディープなデータセットを発生させることができる。

 

 

ダウンロード

依存

brewで導入できる。

SourceForge (Binaryのみ)

https://sourceforge.net/p/sincsimulator/code/ci/master/tree/

git clone https://git.code.sf.net/p/sincsimulator/code sincsimulator-code

注: mac OSXのようなDarwin系列OSでは動作しないので、cent OSに導入した。binaryにパスを通しておく。

 

ラン

 Step 1: Quality profile generation

genProfile -R <read tag(1 for R1, 2 for R2)> -l <read length> <input.txt>順番で記載する。

genProfile -R 1 -l 100 input.txt

 

 

Step 2: Simulation of SNPs, INDELs, CNVs

SInC_simulate ref.fasta

Example:

./SInC_simulate -S 0.002 -I 0.0001 -p 2 -l 1000 -u 150000 -t 2

 

  

Step 3: Read generation

SInC_readGen [options] <in.ref.fa> <read_1_profile.txt> <read_2_prof.txt>  

Example:

./SInC_readGen -C 5 -T 1 -R 100 chr22_allele_1.fa 100_bp_read1_profile.txt 100_bp_read2_profile.txt

./SInC_readGen -C 5 -T 1 -R 100 chr22_allele_2.fa 100_bp_read1_profile.txt 100_bp_read2_profile.txt 

 

 

 

 

 

引用

SInC: an accurate and fast error-model based simulator for SNPs, Indels and CNVs coupled with a read generator for short-read sequence data.

Pattnaik S, Gupta S, Rao AA, Panda B1.

BMC Bioinformatics. 2014 Feb 5;15:40. doi: 10.1186/1471-2105-15-40.

 

SInC: an accurate and fast error-model based simulator for SNPs, Indels and CNVs coupled with a read generator for short-read sequence data. - PubMed - NCBI

 

 

シュードゲノムのシミューレーター Simulome

 

Simulomeは2017年に発表されたbacteria向けの遺伝子のシミュレートツールである。gene情報を与えることで、標準では一部の遺伝子に限定してシミュレートする。具体的には、遺伝子の長さの分布を調べ、その平均と標準偏差から遺伝子のサンプリングをお行い、サンプリングされた遺伝子のみ出力してpseudo-genomeを作成する。このとき、SNPsやindelのフラグをつけることで、指向的な変異を導入することが可能になっている。また、単純にゲノム全体を(変異を入れて)シミュレートすることも可能である。

 

 wiki

https://github.com/price0416/Simulome/wiki

マニュアル

ダウンロードしたフォルダにマニュアルPDFが含まれている。

ダウンロード

依存

  • Python 2.7.2
  • Biopython 1.6.1+
  • BLAST 2.3.0+ 

Github

GitHub - price0416/Simulome: Simulome, genome and mutant variant simulator.

 

ラン

sample_dataディレクトリのデータを使いテストランを行う。500遺伝子サンプリングしてpusedo-genomeを作成する。1遺伝子につき10のSNPsを導入する。

python2.7 simulome.py --genome=sample_data/ecoli_genome.fasta --anno=sample_data/ecoli_anno.gtf --output=ecoli_simulation -- num_genes=500 --snp=TRUE --num_snp=10
  • --genome= File representing genome. FASTA nucleotide format.
  • --anno= File containing genome annotation information in GTF/GFF3 format.
  • --output= Creates a folder named with the supplied prefix containing output files.
  • --num_snp= The number of SNPs to simulate per gene. This argument is required for SNP run mode.

全ゲノムのシミュレートは--num_snpフラグは除き、--whole_genome TRUEをつける。

 

 

500遺伝子サンプリングしてpusedo-genomeを作成する。1遺伝子につき、1つの100-bp挿入を導入する。

python2.7 simulome.py --genome=sample_data/ecoli_genome.fasta --anno=sample_data/ecoli_anno.gtf --output=ecoli_simulation -- num_genes=500 --snp=TRUE --num_snp=10 --indel 3 --ins_len=100
  • --indel File representing genome. FASTA nucleotide format.

 Possible values are:

  1 = Insertions only.

  2 = Deletions only.

  3 = Both insertions and deletions.

  • --ins_len= Length of insertion events. Required for insertion mode.

全ゲノムのシミュレートは

 

Synonymous mutationとnonsynonymous mutationsの指定してアミノ酸変化を起こす塩基置換の割合を指定したり、duplication eventを導入したりできます。かなりのオプションがあるので、詳細はPDFマニュアルを確認してください。

 

  

引用

Simulome: a genome sequence and variant simulator.

Adam Price and Cynthia Gibas

Bioinformatics, 33(12), 2017, 1876–1878 

 

 

MAFFTでマルチプルアライメントを行う

MAFFTはマルチプルアライメントツール。t-coffeeやclustal omegaより高速に動作するとされ、数百のrRNA配列に対してマルチプルアライメントを実行する例が記載されて入る。

 

公式ページ

クイックマニュアル

Manpage of MAFFT

Tips

https://mafft.cbrc.jp/alignment/software/tips0.html

 

インストール

macOSX版ダウンロード

MAFFT for Mac OS X - a multiple sequence alignment program

pkgファイルを叩き、指示に従ってインストールする。brewで導入することもできる。

 

ラン

公式サイトから59のrRNAの配列のFASTAをダウンロードできる(リンク)。 これを使ってみる。

mafft --auto ex1.txt > output

 

またはユーザーと対話式でランもできる。

mafft

 

user$ mafft

 

---------------------------------------------------------------------

 

   MAFFT v7.305b (2016/Aug/16)

 

        Copyright (c) 2016 Kazutaka Katoh

        MBE 30:772-780 (2013), NAR 30:3059-3066 (2002)

        http://mafft.cbrc.jp/alignment/software/

---------------------------------------------------------------------

 

 

Input file? (fasta format)

ex1.txt

Output file?

@ out

Output format?

  1. Clustal format / Sorted

  2. Clustal format / Input order

  3. Fasta format   / Sorted

  4. Fasta format   / Input order

  5. Phylip format  / Sorted

  6. Phylip format  / Input order

@ 1

OK. arguments = --reorder

 

 

Strategy?

  1. --auto

  2. FFT-NS-1 (fast)

  3. FFT-NS-2 (default)

  4. G-INS-i  (accurate)

  5. L-INS-i  (accurate)

  6. E-INS-i  (accurate)

@ 1

OK. arguments = --auto --reorder

 

 

Additional arguments? (--ep # --op # --kappa # etc)

@

command=

"/usr/local/bin/mafft"  --auto --clustalout --reorder "ex1.txt" > "out"

OK?

@ [Y] y

 

 

web版もいくつかの機関で構築されています。ターミナルでランできなければこちらを使ってみてください。

MAFFT - a multiple sequence alignment program

 

 

*Sierra環境でconfig errorが出たが、"unset MAFFT_BINARIES" を打つことで修復された。

 

 

引用

MAFFT: a novel method for rapid multiple sequence alignment based on fast Fourier transform.

Katoh, Misawa, Kuma and Miyata

Nucleic Acids Res. 30:3059-3066, 2002

 

MAFFT version 5: improvement in accuracy of multiple sequence alignment.

Katoh, Kuma, Toh and Miyata

Nucleic Acids Res. 33:511-518, 2005

 

PartTree: an algorithm to build an approximate tree from a large number of unaligned sequences.

Katoh and Toh

Bioinformatics 23:372-374, 2007

 

MAFFTを使ってマルチプルアラインメントを行う 統合TV(togotv)|生命科学系DB・ツール使い倒し系チャンネル

トランスポゾン検出ツール6 Tangram

 Tangramはトランスポゾンの検出に特化した構造変化検出ツール。SV検出で用いられるread-pairとsplit-readのアルゴリズムを使い高感度にトランスポゾンを検出する。1000ゲノムでもmobile element検出ツールとして用いられた。トランスポゾン検出ツールは様々報告されているが、Tangramはsplit read情報とread pair情報を両方使い、1塩基の精度でトランスポゾン挿入位置を検出することが可能である。

 

 

 

インストール

github

https://github.com/jiantao/Tangram

 

依存

  • g++ 4.2.0 and above
  • zlib
  • pthread lib

 

 

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

git clone git://github.com/jiantao/Tangram.git 
cd src
make

 bin/にパスを通しておく。

 

 

ラン

以下のフローでトランスポゾンを検出する。

 

f:id:kazumaxneo:20171017213212j:plain 

Gitのマニュアルより

 

MOSAIKでアライメントしたza tag付きのbamが必要。なければ、bamにZA -tagを付加する作業を最初に行う。fastaのindexがなければ、ここで自動作成される。

step1 ZA tagの付加。

tangram_bam -i input.bam -r ref.fa -o ZA-tagged.bam

 <@/&:MQ1:MQ2:SP_REF:NUM_MAP:CIGAR:MD>のようなTagがbamに付加される。詳細はGitのマニュアルを参考にしてください。

 

step2 bamのスキャン。

tangram_scan -in ZA-tagged.bam -dir output

 

 

step3 index。

tangram_index -ref -sp output/input.bam -out output/

 

step4 トランスポゾンの検出。

tangram_detect

 

step5 フィルタリング。

tangram_filter

 

作成中。

 

 

引用

Tangram: a comprehensive toolbox for mobile element insertion detection

Jiantao Wu, Wan-Ping Lee, Alistair Ward, Jerilyn A Walker, Miriam K Konkel, Mark A Batzer and Gabor T Marth

BMC Genomics201415:795

TagDust2によるアダンプタートリミング

 TgaDust2は、アダプター、バーコード、単純リピートなどの不要な情報を見つけて除去するツール。2009年にTagDDustが発表され、その後2015年にTagDust2が発表された。

 

公式サイト

TagDust

 

 

インストール

brewで導入できる。

brew install TagDust

brewではTagDust2が導入されている。

user$ tagdust

 

Tagdust 2.33, Copyright (C) 2015 Timo Lassmann <timolassmann@gmail.com>

 

Usage:   tagdust [options] <file>  -o <output prefix> 

 

Options:

-Q                      FLT       confidence threshold [20].    

-l                      STR       log file directory name.      

-start                  INT       start of search area [0].     

-end                    INT       end of search area [length of sequence].

-format                 STR       format of input sequence file.

-minlen                INT       minimal accepted read length [16].

-ref                    STR       reference fasta file to be compared against[].

-fe                     INT       number of errors allowed when comparing to reference[2].

-dust                   INT       remove low complexity sequences. [100].

-e                      FLT       expected sequencer error rate [0.05].

-o                      STR       output file name.             

-a                      STR       output file for artifacts [NA].

-t                      INT       number of threads [8].        

-show_finger_seq         NA       print fingerprint as sequence (default is as base 4 number).

-h/help                  NA       print help.                   

-v/version               NA       print version number.         

-1                      STR       type of the first HMM building block.

-2                      STR       type of the second HMM building block.

-...                    STR       type of the . . . HMM building block.

 

 

 

ラン

tagdust adapter.fasta input.fastq -e 0.05 -fe 2 -o output.fastq -a artifactual.fastq -t 12
  • -dust remove low complexity sequences. [100].
  • -e expected sequencer error rate [0.05].
  • -o output file name.
  • -a output file for artifacts [NA].
  • -t number of threads [8].
  • -fe number of errors allowed when comparing to reference[2].
  • -minlen minimal accepted read length [16].

 

 Scytheのオーサーらは、Scytheでアダプター除去し、Scytheで除けなかった5'側のアダプターなどをTagDustで除去する2段構えのアダプター除去プロセスを提案しています(その場合のワークフローは、Scythe-> TagDust2 -> クオリティトリミング-> QC、の流れにすべきとしています)。Scytheは以前紹介しています。


 

 

 

引用

TagDust--a program to eliminate artifacts from next generation sequencing data.

Lassmann T1, Hayashizaki Y, Daub CO.

Bioinformatics. 2009 Nov 1;25(21):2839-40. doi: 10.1093/bioinformatics/btp527. Epub 2009 Sep 7.

 

TagDust2: a generic method to extract reads from sequencing data.

Lassmann T1,2.

BMC Bioinformatics. 2015 Jan 28;16:24. doi: 10.1186/s12859-015-0454-y.

 

http://catway.jp/bioinformatics/qc/removeseq.html

 

embossのseqretでFASTAを修復する

FASTAをいじっていると、何らかの拍子に構造がおかしくなってsamtoolsのindexでsegmentation errorを起こすことがある。途中に空行ができていたり、特殊文字が入っていたり、何らかの理由があるわけだが、embossのseqretを使うと簡単に修復することができる。seqretは入力ファイルを分析し、パースして標準的なNCBIFASTA形式で出力することに使われるコマンドである。

 

公式サイト

http://emboss.sourceforge.net/apps/release/6.6/emboss/apps/seqret.html

 

インストール

embossbrewで導入できる。

brew install emboss

 

ラン

seqret 

入力のFASTAと出力のFASTA名を順番に入力する。

 user$ seqret

Read and write (return) sequences

Input (gapped) sequence(s): input.fasta 

output sequence(s) [chr.fasta]:out.fa

 

これだけでFASTAを修復できる。

 

 

引用

EMBOSS: The European Molecular Biology Open Software Suite

Rice P1, Longden I, Bleasby A.

Trends Genet. 2000 Jun;16(6):276-7.

 

http://seqanswers.com/forums/showthread.php?t=2352