macでインフォマティクス

macでインフォマティクス

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

bamからbigWigとWiggle Formatに変換するツール

 

bamからwiggleファイルに変換して、viewerに取り込み、カバレッジtrackとして表示できる。ただしそれにはsamtoolsのpileupを使いbamからwiggleファイルを作る必要があり、作り方がやや面倒だった。現在では、コマンド一発でwiggleファイルを作るツールが公開されている。このツールを使ってbigwig/wiggleに変換し、IGVにカバレッジtrackを追加する流れを書いておく。

 

インストール

wiggleとbigwigに変換するツール

https://github.com/MikeAxtell/bam2wig

バイナリなのでpathを通して実行権を与えるだけでランできる。ただしbigwigに変換するには、下のリンクからwigtobigwigツールをダウンロードしてパスを通しておく必要がある。wigtobigwigがないと、wiggleファイルだけ作られる。

 

wigtobigwigのダウンロード - unix

Index of /admin/exe/linux.x86_64

 

wigtobigwigのダウンロード - mac OSX

bam2wig (linux)

 > ./wigToBigWig

$ wigToBigWig

wigToBigWig v 4 - Convert ascii format wig file (in fixedStep, variableStep

or bedGraph format) to binary big wig format.

usage:

   wigToBigWig in.wig chrom.sizes out.bw

Where in.wig is in one of the ascii wiggle formats, but not including track lines

and chrom.sizes is a two-column file/URL: <chromosome name> <size in bases>

and out.bw is the output indexed big wig file.

If the assembly <db> is hosted by UCSC, chrom.sizes can be a URL like

  http://hgdownload.cse.ucsc.edu/goldenPath//bigZips/.chrom.sizes

or you may use the script fetchChromSizes to download the chrom.sizes file.

If not hosted by UCSC, a chrom.sizes file can be generated by running

twoBitInfo on the assembly .2bit file.

options:

   -blockSize=N - Number of items to bundle in r-tree.  Default 256

   -itemsPerSlot=N - Number of data points bundled at lowest level. Default 1024

   -clip - If set just issue warning messages rather than dying if wig

                  file contains items off end of chromosome.

   -unc - If set, do not use compression.

   -fixedSummaries - If set, use a predefined sequence of summary levels.

   -keepAllChromosomes - If set, store all chromosomes in b-tree.

> ./bam2wig

r$ bam2wig 

bam2wig

Version 1.4

 

USAGE:

bam2wig [options] [.bam]

 

DEPENDENCIES:

samtools [in PATH]

wigToBigWig [in PATH -- optional]

 

OPTIONS:

-t [float] : Threshold for reporting. Defaults to 0 (all depths reported)

-s [top/bottom/both] : Strand to report. Defaults to both. 'bottom' will quantify in negative numbers. Used in conjunction with option -g will subset the loci to those on the indicated strand.

-r [string] : Limit analysis to specified read group

-l [integer] : Minimum size of read to report. Default: no minimum.

-L [integer] : Maximum size of read to report. Default: no maximum.

-c [chr:start-stop] : Limit analysis to the indicated locus

-g [string] : Path to a GFF3 file. Loci in the GFF3 will be reported. Negates any setting of -c.

-F [integer] : Argument to be passed to samtools view's -F option ... which alignments will be ignored. Defaults to 3844

-m : Normalize reads to reads per million before calculating depths

-d : Treat as degradome/PARE data .. tally only the depths of 5' ends. Option -s must be either 'bottom' or 'top'

-D : Name of output directory. Defaults to [/YourBamsPath/YourBamsBasename_bam2wig/].

-i : Instead of abundance, instead report fraction of reads meeting the size requirements given by -l and/or -L

-h : Get help (print this message)

-v : Print version number and quit

 

DOCUMENTATION:

type 'perldoc bam2wig'

  

ラン

bam2wig input.bam 

使いそうなオプション

-r  : Limit analysis to specified read group

-l : Minimum size of read to report. Default: no minimum.

-L : Maximum size of read to report. Default: no maximum.

-c  : Limit analysis to the indicated locus.

-m : Count each as (1/(n mapped reads)) * 1E6. In other words, normalize to reads per million.

ランが終わると、bamのあるディレクトリにbamファイル名のディレクトリが追加される。その中にbigwigとwig形式のファイルが作られる。

 

IGVで表示する。まずbamを取り込む。

igv -g input.fasta input.bam

f:id:kazumaxneo:20170720151839j:plain

IGVが起動してbamが表示された。さらにFIle -> Load From filesからbigwigファイルを追加する。

 

カバレッジが深い部位があると変化が見えにくい。その場合、パネルの左側で右クリックし、scaleをautoに変更する。

f:id:kazumaxneo:20170720151202j:plain

 

ついでにtrack のheightも変更する。ここでは250にした。

f:id:kazumaxneo:20170720151301j:plain

カバレッジ変化が見やすくなった。

f:id:kazumaxneo:20170720151927j:plain

 

 bigwigはbainaryだが、wiggleファイルはテキストファイルである。下はwiggle変換ツールで作ったwiggleファイルの中身を見たものである。

variableStep chrom=AP017303.1 span=90

390 2

variableStep chrom=AP017303.1 span=2

522 2

variableStep chrom=AP017303.1 span=87

524 4

variableStep chrom=AP017303.1 span=1

611 5

variableStep chrom=AP017303.1 span=2

612 3

variableStep chrom=AP017303.1 span=4

614 1

variableStep chrom=AP017303.1 span=83

618 3

variableStep chrom=AP017303.1 span=7

IGVでこの部位のリードを見てみる。

f:id:kazumaxneo:20170720144929j:plain矢印部分がwiuggleファイルで表記されている部位である。

 

左側を拡大する。

f:id:kazumaxneo:20170720145200j:plain

wiglgleでは

variableStep chrom=AP017303.1 span=90

390 2

となっている。つまり390からspanの90ntだけカバレッジ1の領域が続くことを表している。

  

wiggle formatの詳細はUCSCのリンクを確認してください。

Genome Browser Wiggle Track Format

 

他にdeeptoolsのbamCoverageを使う方法もある。

引用

https://github.com/MikeAxtell/bam2wig

Biostars

https://webcache.googleusercontent.com/search?q=cache:VrSkb-8q6SMJ:https://www.biostars.org/p/2699/+&cd=1&hl=ja&ct=clnk&gl=jp