macでインフォマティクス

macでインフォマティクス

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

シングルセルのアセンブルツール HyDA

 大部分の微生物は一般的な培地では培養不能で、シングルセルシーケンスは微生物の洞察を得るための重要な方法となっている。シングルセルシーケンスには、全ゲノム増幅によってDNAをランダム増幅させる必要があるが、シーケンスバイアスが非常に大きいという問題点がある。

Hybrid DeNovo Assembler(HyDA)はシングルセルのアセンブラ。ディープ・シングルセル・シーケンスを少し浅い数個のシングルセル・シーケンスに置き換え、コアセンブリーすることで、シーケンシング・コストを大幅に増加させることなく、サンプル損失のリスクを回避してハイクオリティなアセンブルを行うことができるとされ、論文中では、複数細胞から控えめにシーケンスしたfastqをmixし、SPAdesおよびIDBA-UDと比較してより長いコンティグを作ることを示している。入力データセットのシェアされた配列と排他配列を検出して分離することによって、キメラアセンブリを回避する機能を持つ。

 

 公式ページ(シングルセルシーケンスデータリンクは下の方)

Synergistic single cell sequencing

 

インストール

cent OSに導入した。

本体 SourceForge

https://sourceforge.net/projects/hyda/?source=typ_redirect

ソースコードをダウンロードし、解凍してビルドする。

 

tar xvf squeezambler-2.0.3-hyda-1.3.1.tar
cd squeezambler-2.0.3-hyda-1.3.1/
make clean; make #k-merのmaxは64、もっと増やすならMAXK=Nで指定する
cd bin/

 > ./import -h

$ ./import -h

Usage: import [options]

       Output (FASTA by default) file name is identified by -O option.

       Input files are identified by -q and -f options.

       Example: import -f=lib1-1.fasta,lib1-2.fasta -q=lib2-1.fastq,lib2-2.fastq -f=lib3.fasta -O=reads.fasta

       Above, lib1-1.fasta and lib1-2.fasta are paired reads in library 0, lib2-1.fastq and lib2-2.fastq are paired reads in library 1, and lib3.fasta are single reads in library 2.

 

 

-V,--version      prints the version

-h,--help      shows this help

-O,--output      =output file

-q,--fastqlibrary      =comma separated list, without space, of FASTQ file names in the library

-f,--fastalibrary      =comma separated list, without space, of FASTA file names in the library

-t,--truncate      =bad quality flag, truncate trailing nucleotides with the bad quality value in FASTQ files

-l,--minlen      =minimum length of a read (default 20)

-n,--maxn      =maximum number of N's in the relevant part of a read (default 4)

  --fastqout      fastq output

 

> ./assemble-unitig -h

$ ./assemble-unitig 

Error: inputs not specified, try 'assemble-unitig -h' for help.

[uesaka@cyano bin]$ ./assemble-unitig -h

Usage: assemble-unitig [options] reads-file

       Output files prefix is identified by -O option.

 

 

-V,--version      prints the version

-h,--help      shows this help

-O,--output      =output files prefix (mandatory)

-C,--config      =coloring config file, if absent all libraries are colored 0

-k,--kmer      =k as in kmer (mandatory)

-p,--procs      =number of threads per processing node (default 1)

-q,--fastq      input is in FASTQ format

-m,--machines      =number of machines (only needed when distributing among multiple machines)

-s,--slice      =0-base slice number (only needed when distributing among multiple machines)

-S,--single      single strand

  --exportkmers      export kmers with their multiplicities

 

 

assemble-finish -h

$ ./assemble-finish -h

Usage: assemble-finish [options] reads-file

       Input unitigs file is identified by -U option.

       Output files prefix is identified by -O option.

 

 

-V,--version      prints the version

-h,--help      shows this help

-U,--unitig      =input unitigs file ([collective] output of assemble-unitig [for multiple slices]) (mandatory)

-O,--output      =output files prefix, may be the same as assemble-unitig output prefix (mandatory)

-C,--config      =pairing config file

-c,--cov      =minimum average coverage of a condensed node, as a multiple of the assembly coverage mean (default 1.0)

-l,--lowcovlen      =maximum length of a low average coverage condensed node to be removed (default 1000nt)

-p,--procs      =number of threads

-q,--fastq      input is in FASTQ format

-S,--single      single strand

-a,--exportasm      export assembly graph for downstream applications, e.g. align-longreads

  --mincontiglen      =minimum length of an output contig (default k)

  --exportgraph      =max contig size, export those nodes of the condensed graph that are shorter than this value, in DOT format

  --exportlowcov      export removed low coverage condensed nodes

  --save      =file to save condensed graph into

  --load      =file to load condensed graph from

  --separatecolors      output contigs of each color into a separate file

  --nolowcovunitigs      do not even initially load those unitigs that have coverage less than 2

  --mincontigcov      =minimum coverage of an output contig [applicable only when --separatecolors used] (default 1.0)

  --keepcolor0      keeps all nodes in color 0 no matter what their coverages are. That is the cut-off for color 0 is 0.

bin/にパスを通しておく(echo export PATH=\$PATH:`pwd`\ >> ~/.bash_profile && source ~/.bash_profile)。

 

 

ラン

アセンブル前に fastqやfastaを HyDAの FASTA/FASTQ フォーマットに変換する必要がある。それにはimportコマンドを使う。公式からダウンロードできるテストデータをランしてみる。/Smith/PairedEnd/K04の4組のペアードエンドfastqを使う。

import -q=2_ACTTGA_L002_R1_001.fastq,2_ACTTGA_L002_R2_001.fastq \
-q=2_ACTTGA_L002_R1_002.fastq,2_ACTTGA_L002_R2_002.fastq \
-q=2_ACTTGA_L002_R1_003.fastq,2_ACTTGA_L002_R2_003.fastq \
-q=2_ACTTGA_L002_R1_004.fastq,2_ACTTGA_L002_R2_004.fastq \
-O=output.fasta
  •  -q --fastqlibrary =comma separated list, without space, of FASTQ file names in the library
  • -O --output =output file

output.fastaが出力される。 

 

k=11でアセンブルする。

assemble-unitig -k=11 -p=40 -O=out output.fasta
  • -p =number of threads per processing node (default 1)
  • -k =k as in kmer (mandatory)
  • -O =output files prefix (mandatory)

logと out.unitigsができる。

 

Finishing。-Uで先ほどの出力を指定する。

assemble-finish -p=40 -c=0.5 -U=out.unitigs -O=final output.fasta
  • -p =number of threads
  • -U =input unitigs file ([collective] output of assemble-unitig
  • -O =output files prefix, may be the same as assemble-unitig output prefix (mandatory)
  • -c =minimum average coverage of a condensed node, as a multiple of the assembly coverage mean (default 1.0)

final.extcontigsができる。

 

 

 

引用

Efficient Synergistic Single-Cell Genome Assembly

Narjes S. Movahedi, Mallory Embree, Harish Nagarajan, Karsten Zengler, and Hamidreza Chitsaz

Front Bioeng Biotechnol. 2016; 4: 42. Published online 2016 May 23.

 

シングルセルシーケンスのカバレッジバイアスを見積もる Preseq

 単一細胞レベルで変異を調べるには、単一細胞のシーケンス決定技術が必要になる。このシングルセルシーケンスの技術は、腫瘍細胞のシーケンスや未培養の細菌集団の細胞の多様性を調べるような研究にも用いられてきた。また、着床前遺伝子診断などに利用して、重篤な変異のない細胞を選抜することに利用できる可能性がある。ただし、単一細胞に存在するのDNAの分子量は非常に小さく、ヒトでは1細胞に100ピコグラム以下のDNAしかない。シーケンスにはナノオーダーのDNAが必要なため、シングルセルのシーケンスには一般に全ゲノムのランダム増幅が必要になる。その方法として非PCRのMultiple displacement amplification (MDA)で使用されるDNAポリメラーゼφ29の利用がよく知られているが、φ29のプライミング効率および伸長速度がヌクレオチド含量に依存し、不均一な増幅をもたらすという知見も報告されている(Preseq論文)。 増幅時の強い増幅バイアスにより、シーケンスのカバレッジには強いバイアスがかかり、カバレッジの大きな変動によるリードがまったくない領域、対立遺伝子の片方しかシーケンスされなかったことによる、多型の消失などが起こる。

 Preseqはduplicate readsでないリード(= distinct reads)の割合を計算することで、シングルセルシーケンスデータが比較的良好なのか、duplicate readsがかなりの割合を占めるのか見積もるツール。Preseqにかけて、増幅バイアスがひどいのかどうか分析することで、今後さらにシーケンスしても意味がないのか、全領域で一定のシーケンスをできそうなのか最低限のコストで推定したり、良好なライブラリ作製法を比較検討するための基礎データにすることができる。

 

公式サイト

Preseq | The Smith Lab

マニュアル

http://smithlabresearch.org/wp-content/uploads/manual.pdf

 

インストール

cent OSに導入した。

Github

https://github.com/smithlabcode/preseq

リリース(リンク)からpreseq_v2.0.2.tar.bz2をダウンロードして解凍して中に入りビルドする。

 make all

> ./preseq

$ ./preseq

preseq:  a program for analyzing library complexity

Version: 2.0.0

 

Usage:   preseq <command> [OPTIONS]

 

<command>: c_curve    generate complexity curve for a library

           lc_extrap  predict the yield for future experiments

           gc_extrap  predict genome coverage low input

                      sequencing experiments

           bound_pop  lower bound on population size

 

 

 

ラン

bedは使用前にソートしておく必要がある(chromosome, start position, end position, and strand)。sortを以下のオプション付きでかける。

sort -k 1,1 -k 2,2n -k 3,3n -k 6,6 input.bed > input.sort.bed

 

予想される歩留まりを計算する。

preseq c_curve -o output.txt input.sort.bed

 

 

90%信頼区間に設定して、ペアードエンドのbamファイルを指定してランする。ステップサイズは5万(リード)とする。

preseq lc_extrap -o yield.txt -B input.sort.bam -P -s 50000 -c 0.90

出力は、total readとdistinct read、その90%信頼区間の上限値、下限値となる。

f:id:kazumaxneo:20180215193316j:plain

 

外挿してdintinct readsを予測している。例えば100万リードシーケンスすればexpected distinct reads (duplicationでないリード) の数は994642と期待できる。duplicationが多い質の低いシングルセルシーケンスだと、total readsの数を増やしてもdistinct readsの数は飽和して伸びないという結果が出るはずである。詳細は公式マニュアル(後半のグラフ)を確認してください。SRAの公共データを数パーセントダウンサンプリングして良好なデータかどうか推測しています。

 

 

引用

Modeling genome coverage in single-cell sequencing

Timothy Daley and Andrew D. Smith.

Bioinformatics. 2014 Nov 15; 30(22): 3159–3165.

 

NGSのリード情報を使いスキャッホールドのギャップを埋める FGAP

 

FGAPはドラフトゲノムのギャップを埋めるためのツール。BLASTを使用して、複数のコンティグをドラフトゲノムアセンブリに対して整列させ、ギャップを埋めるために最良のシーケンスを検出する。ヒトchr14では、ギャップの数を35%減少させたと述べられている。

 

インストール

cent OSに導入した。

依存

  • blast+
  • MCR

MATLAB Runtimeは本体とともにSorceForgeからダウンロードできる。ダウンロードしたら解凍して、インストールディレクトリを指定して./installMCR.shを実行する。homeのMCR/にインストールするなら

cd MCR_LINUX64b/
./installMCR.sh /home/uesaka/MCR/

 

 

本体 SourceForge

https://sourceforge.net/projects/fgap/files/?source=navbar

 

 

ラン

上で導入したMATLAB Runtimeのインストールディレクトリを指定してランする。

 

テストデータを実行する。454のリードでアセンブルしたcontigを伸ばすために、454のリードと、illumina hiseqのリードを指定している。

 ./run_fgap.sh MCR/v717/ -d sample_data/DRAFT_ecoli_hiseq454.fasta -a "sample_data/DATASET_ecoli_hiseq.fasta,sample_data/DATASET_ecoli_454.fasta" -b blast/

1分程度でテストランは終わる。中間ファイルも含めて複数ファイルが作られる。

> ls -alth

$ ls -alth

total 355M

drwxr-xr-x   7 uesaka user   4.0K Feb 14 11:03 .

-rw-r--r--   1 uesaka user   2.2K Feb 14 11:03 output_fgap.stats

-rw-r--r--   1 uesaka user   4.5M Feb 14 11:03 output_fgap_4.fasta

-rw-r--r--   1 uesaka user   4.5M Feb 14 11:03 output_fgap.final.fasta

-rw-r--r--   1 uesaka user   8.5K Feb 14 11:03 output_fgap_4.log

-rw-r--r--   1 uesaka user   4.5M Feb 14 11:03 output_fgap_3.fasta

-rw-r--r--   1 uesaka user    38K Feb 14 11:03 output_fgap_3.log

-rw-r--r--   1 uesaka user   4.5M Feb 14 11:03 output_fgap_2.fasta

-rw-r--r--   1 uesaka user    87K Feb 14 11:03 output_fgap_2.log

-rw-r--r--   1 uesaka user   4.5M Feb 14 11:03 output_fgap_1.fasta

-rw-r--r--   1 uesaka user   221K Feb 14 11:03 output_fgap_1.log

drwx------ 207 uesaka uesaka  36K Feb 14 11:03 ..

 

statsファイル。

$ cat output_fgap.stats

-------------------- GENERAL STATS --------------------

 

Closed gaps (N): 97

 

Before FGAP: 

 Gaps: 123

 Sequences: 73

 Length: 4554392bp 

 GC: 50.7047%

 N50: 172167

 Min: 314

 Max: 414074

 Ns: 2621

 

After FGAP: 

 Gaps: 26

 Sequences: 73

 Length: 4555558bp 

 GC: 50.7279%

 N50: 172155

 Min: 314

 Max: 414012

 Ns: 654

 

Inserted: 3512bp

Removed : 2346bp

 

Closed gaps by each dataset:

 sample_data/DATASET_ecoli_hiseq.fasta: 28 gaps

 sample_data/DATASET_ecoli_454.fasta: 69 gaps

 

Closed gaps by type:

 Negative gap: 0 gaps

 Zero gap: 0 gaps

 Positive gap: 97 gaps

 

 finalのFASTAは、N50はほぼ変化なしだが、gapが123から73に減っている。

 

 

引用

FGAP: an automated gap closing tool.

Piro VC, Faoro H, Weiss VA, Steffens MB, Pedrosa FO, Souza EM, Raittz RT1.

BMC Res Notes. 2014 Jun 18;7:371.

 

bamの分析に使うバイオインフォマティクスのツールキット goleft

 

goleftはMIT licence下で提供されているバイオイオンフォのツールキット。GO言語で構築されている。

 

インストール

Github

https://github.com/brentp/goleft

リリース(リンク)からosx向けバイナリーをダウンロードできる。パスの通ったディレクトリに移動しておく。

chmod u+x goleft-osx.dms
mv goleft-osx.dms /usr/local/bin/goleft
goleft #動作確認

$ goleft

goleft Version: 0.1.17

 

covstats   : coverage stats across bams by sampling

depth      : parallelize calls to samtools in user-defined windows

depthwed   : matricize output from depth to n-sites * n-samples

indexcov   : quick coverage estimate using only the bam index

indexsplit : create regions of even coverage across bams/crams

samplename : report samplename(s) from a bam's SM tag

goleftはBiocondaでも提供されている(リンク)。

 

 

ラン

 1、covstats   bamをサンプリングしてカバレッジとインサートサイズをレポートする。

$ goleft covstats -h

coverage insert_mean insert_sd insert_5th insert_95th template_mean template_sd pct_unmapped pct_bad_reads pct_duplicate pct_proper_pair read_length bam sample

Usage: goleft [--n N] [--regions REGIONS] BAMS [BAMS ...]

 

Positional arguments:

  BAMS                   bam for which to estimate coverage

 

Options:

  --n N, -n N            number of reads to sample for length [default: 1000000]

  --regions REGIONS, -r REGIONS

                         optional bed file to specify target regions

  --help, -h             display this help and exit

 

bamを指定してランする。bedファイルで特定の領域や染色体だけ解析もできる。

goleft covstats pair.bam

カバレッジ、インサートサイズとそのSDなどが出力される。

f:id:kazumaxneo:20180213214311j:plain

cutとcolumnに渡して整形表示。

 

 

 2、depth  一定のウィンドウサイズでbamのカバレッジを計算する。

$ goleft depth -h

Usage: goleft [--windowsize WINDOWSIZE] [--maxmeandepth MAXMEANDEPTH] [--ordered] [--q Q] [--chrom CHROM] [--mincov MINCOV] [--stats] [--reference REFERENCE] [--processes PROCESSES] [--bed BED] --prefix PREFIX BAM

 

Positional arguments:

  BAM                    bam for which to calculate depth

 

Options:

  --windowsize WINDOWSIZE, -w WINDOWSIZE

                         window size in which to calculate high-depth regions [default: 250]

  --maxmeandepth MAXMEANDEPTH, -m MAXMEANDEPTH

                         windows with depth > than this are high-depth. The default reports the depth of all regions.

  --ordered, -o          force output to be in same order as input even with -p.

  --q Q, -Q Q            mapping quality cutoff [default: 1]

  --chrom CHROM, -c CHROM

                         optional chromosome to limit analysis

  --mincov MINCOV        minimum depth considered callable [default: 4]

  --stats, -s            report sequence stats [GC CpG masked] for each window

  --reference REFERENCE, -r REFERENCE

                         path to reference fasta

  --processes PROCESSES, -p PROCESSES

                         number of processors to parallelize.

  --bed BED, -b BED      optional file of positions or regions to restrict depth calculations.

  --prefix PREFIX        prefix for output files depth.bed and callable.bed

  --help, -h             display this help and exit

 

goleft depth --reference ref.fa --prefix output input.bam 

prefixで指定したファイルに指定ウィンドウサイズのカバレッジが出力される。

 

 

 3、indexcov   bam.baiからのカバレッジの超高速な推定。30サンプルを30秒で解析可能。

$ goleft indexcov -h

Usage: goleft --directory DIRECTORY [--includegl] [--excludepatt EXCLUDEPATT] [--sex SEX] [--chrom CHROM] [--fai FAI] BAM [BAM ...]

 

Positional arguments:

  BAM                    bam(s) or crais for which to estimate coverage

 

Options:

  --directory DIRECTORY, -d DIRECTORY

                         directory for output files

  --includegl, -e        plot GL chromosomes like: GL000201.1 which are not plotted by default

  --excludepatt EXCLUDEPATT [default: ^chrEBV$|^NC|_random$|Un_|^HLA\-|_alt$|hap\d$]

  --sex SEX, -X SEX      comma delimited names of the sex chromosome(s) used to infer sex. Set to '' if no sex chromosomes are present. [default: X,Y]

  --chrom CHROM, -c CHROM

                         optional chromosome to extract depth. default is entire genome.

  --fai FAI, -f FAI      fasta index file. Required when crais are used.

  --help, -h             display this help and exit

 

複数のbamをディレクトリに準備してランする。

goleft indexcov --directory output/ bam_data/*.bam

関係ないファイルがあるとエラーになったので、bamとbaiだけ集めたほうがよいかもしれない。htmlで結果は出力される。

f:id:kazumaxneo:20180213220113j:plain

上では2サンプルだけです。 公式のほうがわかりやすい結果を載せています(リンク)。

 

他にも、bamをほぼ等しい量のデータでN個の領域に分割するindexsplitや(計算を並列化するために使う)、SM tagからsamplenameをレポートするsamplenameがある。

 

引用

https://github.com/brentp/goleft

 

 

ペアエンドリードを使いミスアセンブリを検出する misFinder

 ミスアセンブルを検出するツールにはQuest、GAGEなどがあるが、これらのツールはミスアセンブルとリファンレスの違いを区別せず全て報告する。そのため、ミスアセンブルのみを工夫するには通常さらなる工夫が必要になる。 一方、CGALやALE、REAPRはDe novoのアプローチでミスアセンブルを検出する。すなわち、リードをアセンブルした配列にアライメントし、カバレッジ情報、インサートサイズの変動などの情報を使い、異常なアセンブルがないか解析する方法となる。ただしこの方法は、リアルデータのカバレッジのばらつきやscaffoldsのギャップなどでfalse positiveが出る傾向にある。

 misFinderは、リファレンスベースの手法とDe novoの手法をハイブリッドに使うミスアセンブルの検出ツール。複数の情報を適用して、バイアスの少ない検出ができるとされる。

 

 

インストール

Github

https://github.com/hitbio/misFinder

git clone https://github.com/hitbio/misFinder.git 
cd misFinder/
chmod a+x autogen.sh
./autogen.sh
cd bin/
./mf -h

$ ./mf -h

misFinder: v0.4.05.05

Released : Oct 17, 2015

 

Usage: mf <command> [option]

    metrics     Compute the assembly metrics

    misass      Compute mis-assemblies

    all         Do all the above in turn

 

PROGRAM OPTIONS:

  1) metrics -- compute the metrics:

    -conf <FILE>       Configuration file. It is required.

    -m <INT>           The minimal query length. Default is 100.

    -pt <FLOAT>        The minimal identity percentage for matched queries and 

                       matched segments. Default is 0.95.

    -t <INT>           The number of threads for the alignment between queries 

                       and subjects. Default is the number of CPU cores.

    -o <STR>

    -out <STR>         Output directory for the output files. Default is "./"

    -h

    -help              Show help information.

  2) misass -- compute mis-assemblies:

    -conf <FILE>       Configuration file. It is required.

    -i <INT>           Minimal indel size. Default is 5 bp.

    -t <INT>           The number of threads for reads alignment. Default is

                       the number of CPU cores.

    -sc <INT>          Single-cell paired-end data flag. Default is 0.

                       0: standard genomic DNA prepared from culture;

                       1: single-cell data.

    -o <STR>

    -out <STR>         Output directory for the output files. Default is "./"

    -h

    -help              Show help information.

 

パスを通しておく。

 

ラン

confファイルに参照のリファレンスFASTAアセンブルして得たFASTA、paired-endのfastqのパスを記載する。

> cat conf_example

f:id:kazumaxneo:20180213195334j:plain

 

confファイルを指定してランする。

mf all -conf config_file -o output

出力はミスアアセンブルのレポートと、それが修復されたFASTAファイルなどになる。

 

 

引用

misFinder: identify mis-assemblies in an unbiased manner using reference and paired-end reads.

Zhu X, Leung HC, Wang R, Chin FY, Yiu SM, Quan G, Li Y、 Zhang R, Jiang Q, Liu B, Dong Y, Zhou G, Wang Y.

BMC Bioinformatics. 2015 Nov 16;16:386.

 

 

 

 

 

 

 

 

リファンレンスガイドのトランスクリプトのアセンブル TransComb

 

TransCombは、junction graphに基づいて開発されたゲノムガイドのアセンブルツール。ペアのショートリードとリファレンスゲノムを使い、RNA seqのシーケンスデータをアセンブルする。複数種のシミュレーションデータセットとリアルデータセットの両方でテストされ、StringTie、Cufflinks、Bayesembler、およびTraphなどの主要アセンブラに対してリコールと精度が向上していることが述べられている。 動作は高速で、メモリ使用量も少ないとされる。

 

インストール

バイナリが配布されている。

SorceForge ダウンロード

https://sourceforge.net/projects/transcriptomeassembly/files/

 

> ./TransComb

$ ./TransComb 

    

** Error: BAM input file is not provided! **

    

    

===========================================================================

    

TransComb v.1.0 usage:

    

** Required **

    

-b <string>: BAM file produced by Tophat or Tophat2.

    

-s <string>: Strand-specific RNA-Seq reads orientation.

             1) Use <unstranded> to indicate RNA-seq reads are non-strand-specific.

             2) Use <first> to indicate fr-first-stranded RNA-seq reads.

             3) Use <second> to indicate fr-second-stranded RNA-seq reads.

    

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

    

** Options **

    

-h: Output TransComb Help Information

    

-o <string>: Output path/file name of the assembled transcripts GTF, default: ./TransComb.gtf

    

-f <string>: Minimum expression level estimated by abundance analysis for output, default: 0.

    

-l <string>: Minimum assembled transcript length, default: 500.

    

-d <string>: Minimum junction coverage fraction by maximum junction coverage, default: 0.02.

    

-D <string>: Minimum farction of unbalanced junction, default: 0.1.

    

-g <string>: Minimum  gene length, default: 200.

    

-t: Disable trimming of predicted transcripts based on coverage, default: coverage trimming is enabled.

    

-e <string>: Minimum gap length between two exons, default: 200.

    

-F <string>: Minimum seed coverage used for generate a new transcript, default: 0.

    

-a <string>: Minimum anchor length for junctions, default: 1.

    

-m <string>: Fraction of exon allowed to be covered by multi-hit reads, default: 1.

    

-v: Report the current version of TransComb and exit.

    

** Note **

    

A typical command of TransComb might be:

    

TransComb -b file.bam -s unstranded

    

===========================================================================

パスの通ったディレクトリに移動しておく。 

 

ラン

テストランを行う (unstranded)。

cd sample_test/
TransComb -b test.bam -s unstranded -l 200
  • -b <string>   BAM file produced by Tophat or Tophat2.
  • -s <string>   Strand-specific RNA-Seq reads orientation. 1) Use <unstranded> to indicate RNA-seq reads are non-strand-specific. 2) Use <first> to indicate fr-first-stranded RNA-seq reads. 3) Use <second> to indicate fr-second-stranded RNA-seq reads.

 

引用

TransComb: genome-guided transcriptome assembly via combing junctions in splicing graphs.

Genome Biol. 2016 Oct 19;17(1):213.

Liu J, Yu T, Jiang T, Li G.

 

Reference-assisted assembly CSAR

 

CSARは関連する生物のリファレンスゲノムに基づいて、対象のコンティグを効率的に順序付けする方法論。リファンレンスゲノムは、必ずしも完全でなくても動作する。

 

インストール

依存

  • MUMmer whole genome alignment package

無ければbrewで導入しておく。

本体 Github

https://github.com/ablab-nthu/CSAR

git clone https://github.com/ablab-nthu/CSAR.git
cd CSAR/
php csar.php -h #ヘルプ

$ php csar.php -h

Usage: php csar.php [option] -t <target_contigs.fna file> -r <reference.fna file> [--nuc, --pro]

Option:

-t <string>   Target genome (i.e., draft genome to be scaffolded)

 

-r <string>   Reference genome

 

--nuc         Use NUCmer to identify conserved genetic markers between target and reference genomes

 

--pro         Use PROmer to identify conserved genetic markers between target and reference genomes

 

-o <string>   Output folder to contain all the files returned by running CSAR (default: ./csar_out)

 

-h            Show help message

 

ラン

gitでcloneすればテストデータもダウンロードされる。

テストランを実行する。

php csar.php -t example/M.luteus_contigs.fna -r example/GCA_001691605.1_reference.fna --nuc -o example_output
  • -t <string>   Target genome (i.e., draft genome to be scaffolded)
  • -r <string>   Reference genome
  • --nuc    Use NUCmer to identify conserved genetic markers between target and reference genomes
  • -o <string>   Output folder to contain all the files returned by running

10秒程度でランは終わる。出力ディレクト

f:id:kazumaxneo:20180206001727j:plain

scaffolds.nuc.csarにcontigのorder結果がプリントされている。gi|240114495|ref|NZ_CABC01000027.1|などはcontigの名前、右端の0と1はコンティグの向きを表す。

f:id:kazumaxneo:20180206001837j:plain

scaffolds.nuc.csar.fnaはアセンブル結果のFASTAとなる。NNNで繋がれた領域も残っている。

 

 

引用

CSAR: a contig scaffolding tool using algebraic rearrangements.

Chen KT, Liu CL, Huang SH, Shen HT, Shieh YK, Chiu HT, Lu CL.

Bioinformatics. 2018 Jan 1;34(1):109-111.