2014年11月3日月曜日

BEDファイルのソート

BEDファイルがでたらめにソートしていると、
chr1
chr10
chr11
chr12
…
chr19
chr2
chr20
chr21
chr22
chrX
chrY
という感じになっている。つまり、chr~ごと文字として認識しているから。
sort -k 1,1 -k2,2n hoge.bed 

#とすれば、以下のように理想的な形にソートされる。
chr1
chr2
chr3
…
chr21
chr22
chrX
chrY

2014年9月2日火曜日

MacでホストOSとゲストOSでファイルを共有する方法 on VirtualBox

だいたいググって出てくる結果と変わらない。

  1. 正常にMacにVirtualBoxがインストールされていれば、 ~/VirtualBox VMsができてる。
  2. ここに適当にVM-Shareという適当なディレクトリを作成する。
  3. VirtualBoxを起動して、設定の共有フォルダで右側の+をクリック
  4. 先ほど作ったディレクトリのパスを入力して、適当にshareと名前を付ける。
  5. 一度VirtualBoxを落として再起動、適宜イメージを起動。
  6. ターミナルで、sudo mkdir /mac
  7. パス聞かれるから入力。入力しても文字でないけど。気にしない。
  8.  sudo mount -t vboxsf share macとすれば、ゲスト側のmacというディレクトリとホスト側のVM-Shareが共有される。

2014年8月11日月曜日

liftover

既に解析済みのbedやgif, positionデータはあるものの、現在自分が解析しているゲノムのバージョンとは異なる場合があるが、UCSCにあるliftoverでそろえる事が可能。

http://hgdownload.cse.ucsc.edu/admin/exe/で環境に合わせて
私の場合はlinux.x86_64/以下にあるliftoverをダウンロードchmod +x liftoverでlocal/binなどパスが通っている場所に置いておく。さらに、例えば、mm8からmm10に変換したい場合は、http://hgdownload.soe.ucsc.edu/goldenPath/mm8/liftOver/mm8ToMm10.over.chain.gzをダウンロードして解凍。適当に共有ディレクトリ等に置いておく。

#変換したいファイルがbedなら(chr1 1000099 1000200)
liftOver Oct4_mm8.bed mm8ToMm10.over.chain Oct4_mm10.bed unmapped.txt
Oct4_mm10.bed>>変換されたbed 
unmapped.txt>>変換できなかったもの
#変換したいファイルがpositionなら(chr1:1000100-1000200)
liftOver -positions Oct4_mm8.txt mm8ToMm10.over.chain Oct4_mm10.txt unmapped.txt
Oct4_mm10.txt>>変換されたposition
unmapped.txt>>変換できなかったのもの
liftOver_XXXXXX.bed>>変換前のbedファイル
liftOver_XXXXXX.bedmapped>>変換後のbedファイル
positionでもliftover出来て、そのbedを得られるのは便利

GEOのRAWデータをターミナルからダウンロード

あるアクセッション番号のRAWデータをダウンロードしたい時に、webのリンクからだと直接ダウンロードできないので、以下の様に自分のダウンロードしたい番号を入れる。親ディレクトリは下三桁をnnnとして、それより上の桁をダウンロードしたい番号と同じにする。
以下は例
wget ftp://ftp.ncbi.nlm.nih.gov/geo/series/GSE11nnn/GSE11431/suppl/GSE11431_RAW.tar

2014年7月25日金曜日

サーバーの特定のディレクトリをマウントする

サーバー側で作成したpdfやpngなどのファイルをいちいちscpするのはめんどくさいのでマウントして手元のPCから直接覗けるようにする。 基本的にhttp://blog.offline-net.com/2014/03/02/mavericks-osxfuse-sshfs-local2remoteをみればわかる 私も上記解説と同じ環境である。 デスクトップにhogeという名前のフォルダを先に作成する。 以下の様にターミナルで叩けばマウントできる。
sshfs ユーザー名@サーバー名(hoge.co.jp):/home/ユーザー名/PATH/TO/指定ディレクトリ /Users/ユーザー名/Desktop/hoge
指定ディレクトリ以下にディレクトリがあってもアクセスできる。

fastaファイルの特定の位置をマスクする方法

bedtoolsの中のmaskFastaFromBedというツールを使う。 使い方は
maskFastaFromBed -fi [マスクかけたいfasta ファイル] -bed  [位置を指定した bedファイル] -fo [ 出力fa staファイル名] -mc N(Nでマスクする)
マスクかけたいfasta ファイルはmulti fastaでOK bedも複数行でOK 例として
test.fa
>chr1
TCTCATGCTAATCAGGATCTGGGGCACAGTACAAACCTCAAATGCTGTAGTTGAG
>chr2
GCTATTTCAGCCCTAGTGGAAGGGGGTGGGGCGGCACTAGGGCAGCAGCCCGTGCCCGCTGCACTTCCCTTATCGGATTACTGCTTAACTCCCACCCCCTCCTAGGTCCAAAAAGTATACTAATTACATGGA
>chr3
GCAGTACATTCCACAGCATCCTTAGGTGTTGCTCAGAGTCTACCTAACACCTAAGTACAAAAATGAAGCCCTAGGGCAGCCCTAAATGGTTTGTTTCCGCAAAATAGTGATCAACTCAAGTAAACAAGGGTT
>chr4
CATAAGACACATTGAAGACACCGCCACAGAACTTTACTGCCCTCAAGAGTCTTCTCTCATCTCTGGGATGCCTGTGAGGC

test.bed
chr1 5 10
chr2 10 20
chr2 21 23
chr5 10 20

maskFastaFromBed -fi test.fa -bed  test.bed -fo out.fa -mc N

out.fa
>chr1
TCTCANNNNNATCAGGATCTGGGGCACAGTACAAACCTCAAATGCTGTAGTTGAG
>chr2
GCTATTTCAGNNNNNNNNNNANNGGGTGGGGCGGCACTAGGGCAGCAGCCCGTGCCCGCTGCACTTCCCTTATCGGATTACTGCTTAACTCCCACCCCCTCCTAGGTCCAAAAAGTATACTAATTACATGGA
>chr3
GCAGTACATTCCACAGCATCCTTAGGTGTTGCTCAGAGTCTACCTAACACCTAAGTACAAAAATGAAGCCCTAGGGCAGCCCTAAATGGTTTGTTTCCGCAAAATAGTGATCAACTCAAGTAAACAAGGGTT
>chr4
CATAAGACACATTGAAGACACCGCCACAGAACTTTACTGCCCTCAAGAGTCTTCTCTCATCTCTGGGATGCCTGTGAGGC
となるはず。 bedは定義がstartの数字を含まないので5と書いてあれば6文字目を示しendは10と書いてあれば10文字目を示す事に注意する。

2014年7月24日木曜日

マップしなかったリード・マップしたリード・ユニークなリードの数え方

bowtieが終わったあとに表示される以下の様なリード数が表示されるがこれをあとでコマンドから確認する方法
135582741 reads; of these:
135582741 (100.00%) were unpaired; of these:
22383842 (16.51%) aligned 0 times
61518065 (45.37%) aligned exactly 1 time
51680834 (38.12%) aligned >1 times
83.49% overall alignment rate 
基本的にフラグが立っているのでそれを利用する。
#reads; of these:これはfastQCした結果のhtmlを見るかbowtieにかける前のfastqをwc -lして4で割れば良い。
#aligned 0 times:マップされなかったリード数の表示は
samtools view -cf 4 filename.bam
#マップされたリード数の表示は
samtools view -cF 4 filename.bam
#aligned >1 times:マルチマップされたリード数の表示は
grep 'XS:' ./filename.sam | wc -l
#aligned exactly 1 time:ユニークにマップされたリード数の表示は
マップされたリード数からマルチマップされたリード数を引くか
sed '/XS:/d' Enh1.sam | wc -lとして、全リード数からマルチマップのリードをsamから取り除いた行カウントを行う。その値からsamのheaderのライン数とマップされなかったリード数を引く。
ちなみに、pairをきちんと組んでいるものを抽出するには この組み合わせ
samtools view -S -f 99 -S filename.sam > ./filename.99.sam
samtools view -S -f 147 -S filename.sam > ./filename.147.sam
この組み合わせ
samtools view -S -f 83 -S filename.sam > ./filename.83.sam
samtools view -S -f 163 -S filename.sam > ./filename.163.sam
http://ppotato.wordpress.com/2010/08/25/samtool-bitwise-flag-paired-reads/がわかりやすい。
<<<<<追記>>>>>>
 XS:による方法は非常に限られた状況でのみ成立するっぽいです.
 bowtie2でSEでやった時は少なくともbowtieを走らせたあとの結果と一致しました。
ちなみにPEでやったものを同様に行うと、想定よりも大きい値が出たため、PEにはこの方法は適応できません。
XS:i:<N>
Alignment score for the best-scoring alignment found other than the alignment reported. Can be negative. Can be greater than 0 in --local mode (but not in --end-to-end mode). Only present if the SAM record is for an aligned read and more than one alignment was found for the read. Note that, when the read is part of a concordantly-aligned pair, this score could be greater than AS:i. ってbowtie2のmanualに記載があるけれどよくわかりません( TДT)

SE/PEどちらの場合も、よくいわれるNH:i:のフラグはbowtie2では立ちませんでした。MAPQの値で行う事も出来ませんでした。

確実にmulti alignを取り除く方法知っている方いたら教えてくださいm(_ _ )m

Trimmomatic

だいぶ慣れてきたので、インストールの仕方よりもツールの使い方をメインに記載するかな。 アダプタートリマーツールであり、クリーニングも出来るTrimmomatic (http://www.usadellab.org/cms/?page=trimmomatic) 使い方としては、デフォルトにadoptersにいくつか入っているが、適当に自作してよい。例えば、Truseq4-SE.faとして
>TruSeq_v1v2_Index12Adapter
GATCGGAAGAGCACACGTCTGAACTCCAGTCACCTTGTAATCTCGTATGCCGTCTTCTGCTTG
TruSeq_v1v2_UniversalAdapter
AATGATACGGCGACCACCGAGATCTACACTCTTTCCCTACACGACGCTCTTCCGATCT
というファイルを作って、以下のコマンドを叩く。
java -jar Trimmomatic-0.32/trimmomatic-0.32.jar SE \
-threads 8 -phred33 \
-trimlog ./DNseqE14.log.txt \
 ./DNseqE14.fastq \
 ./trimmomatic/DNseqE14.trim.fastq \
 ILLUMINACLIP:../../../../../local/Trimmomatic-0.32/adapters/TruSeq4-SE.fa:2:40:15 \
 LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:20
ILLUMINACLIP自体がアダプタートリム。その他はクオリティ値でトリム。 ペアエンドの時は少し書き方が異なるがhttp://www.usadellab.org/cms/?page=trimmomaticに参考例の記載がある。
このツールはシーケンスするフラグメントが非常に短い時にインサートを突き抜けてアダプターもシーケンスしてしまった時の除去に非常に有効なツール

2014年5月18日日曜日

ChIPpeakAnnoのTSSデータを追加する

ChIPpeakAnnoはデフォルトでいくつかのTSSデータが入っている (R-version/library/ChIPpeakAnno/data)
今回はデフォルトにはない、酵母のTSSの情報を得るのが目的

既にbiomaRtがインストールされている前提で、R上で以下のコマンドを入力
mart1 = useMart(biomart="ensembl")
listDatasets(mart1)
とするとリストがでてくる。
13番目にscerevisiae_gene_ensemblがある。

sacCer3= useMart(biomart="ensembl", datset="scerevisiae_gene_ensembl")
TSS.scerevisiae.sacCer3 = getAnnotation(mart=sacCer3, featureType="TSS")
として、
ls()とすると
mart1
TSS.mouse.sacCer3
sacCer3
となるので
rm('mart1')
rm('sacCer3')
として、TSS.mouse.sacCer3のみの状態にしてから
save(list=ls(), file = "./TSS.scerevisiae.sacCer3.rda")
Rを終了させて、カレントディレクトリでlsするとTSS.scerevisiae.sacCer3.rdaがあるのでこれを
/PATH/TO/R-version/library/data/
以下にコピーすれば、デフォルトに入っているマウスやヒトのデータを扱うのと同じように扱えるようになる。

2014年5月16日金曜日

シーケンスデータの一括ダウンロード

論文を見るとGEOのaccsession番号が書かれているが複数のシーケンスデータが置かれている事がほとんどで、前回は一つだけをダウンロードする方法を書いた。深く解析しようと思って全てのシーケンスデータをダウンロードしたいが、一個ずつやってるとめんどくさいので。シェルスクリプトを記述する事でそれを解決する。
シーケンスデータは多くの場合連番でSRAにupされているが、SRRだとかSRXだとか色々番号があって対応関係を把握するのが先決。
例えば以下GSMの番号は連番。
一番上のリンクへ行くと
最下部にとあるのでftpへ行く(chromeで見た方がいい)
となっている。SRR824704にSRR824704.sraという目的のファイルが入っている。
同様にGMS1111729もみてみると
となっていて、SRR824727にSRR824727.sraという目的のファイルが入っている。
なので今回はSRX263794-263817とSRR824704-824727がそれぞれ対応しているという事がわかる。シーケンスデータによっては全てのデータが連番ではないこともあるので注意が必要。
シェルスクリプトの作法はよくわからないのでググるヽ(゚∀゚)メ(゚∀゚)メ(゚∀゚)ノ
vimを使っているので、
vim SRA_download.sh
としてファイルを作る。
シェルスクリプト作る際は思った挙動するかわからないので、echoでおもっているようなスクリプトがかかれているか確認する。
#!/bin/sh
# 全てのSRX Number
srx_list=({263794..263817})
# 全てのSRA Accussion Number
sra_list=({824704..824727})
# ダウンロードする FTP サイトの共通部分までのパス
base_url=ftp://ftp-trace.ncbi.nlm.nih.gov/sra/sra-instant/reads/ByExp/sra/SRX%2FSRX263/
for i in {0..23}
do
        echo "wget ${base_url}/SRX${srx_list[i]}/SRR${sra_list[i]}/SRR${sra_list[i]}.sra"
        #sleep 10m
done
と記述して上書きする(esc :wq)
chmod u+x SRA_download.sh
としてから
./SRA_download.sh
とすると確認が出来るので出てきたURLからftp://以下から.sraの一つ上流までをコピーしてブラウザーで確認する。chromeでftp://ftp-trace.ncbi.nlm.nih.gov/sra/sra-instant/reads/ByExp/sra/SRX%2FSRX263/SRX263809/SRR824719とすれば、先ほどのリンクにたどれる。
これで一応問題ないので
#!/bin/sh
# 全てのSRX Number
srx_list=({263794..263817})
# 全てのSRA Accussion Number
sra_list=({824704..824727})
# ダウンロードする FTP サイトの共通部分までのパス
base_url=ftp://ftp-trace.ncbi.nlm.nih.gov/sra/sra-instant/reads/ByExp/sra/SRX%2FSRX263/
for i in {0..23}
do
       wget ${base_url}/SRX${srx_list[i]}/SRR${sra_list[i]}/SRR${sra_list[i]}.sra
       sleep 10m
done
とecho" "と#を取り払って、上書きする(esc :wq)
./SRA_download.sh
でダウンロードが始まる。

説明しとくと
#!/bin/sh---シェルスクリプトのおまじない
# 全てのSRX Number
srx_list=({263794..263817})先に調べておいたSRXの番号のリスト
# 全てのSRA Accussion Number
sra_list=({824704..824727})先に調べておいたSRRの番号のリスト
# ダウンロードする FTP サイトの共通部分までのパス
base_url=ftp://ftp-trace.ncbi.nlm.nih.gov/sra/sra-instant/reads/ByExp/sra/SRX%2FSRX263/
#全データが24個あるので24回ループさせる。(1番目は0らしい。)
for i in {0..23}
#なんかやらせる時-do
do
#URL共通部分/SRX+番号(リストのi番目)/SRA+番号(リストのi番目)/SRA+番号(リストのi番目).sraってこと
${base_url}/SRX${srx_list[i]}/SRR${sra_list[i]}/SRR${sra_list[i]}.sra
#ftpって色々制限あったりするらしいのでダウンロード終わってから10分待ってから次のダウンロード始める。
sleep 10m
#やらせるの終わったら-done
done


2014年4月19日土曜日

Rを使ったChIP-seq解析(近傍遺伝子や距離等のアノテーション)


高次解析に必要な関数が色々作られていて、その関数集団をパッケージと表現するらしい。ライフサイエンスのパッケージ集にbioconductorというものがありHPは http://www.bioconductor.org


ターミナルを開いて、local/work/以下にRディレクトリをさくせいして、Rにて作業する。ちなみにMacsまでのデータはlocal/workにある。
cd local/work
mkdir R
cd R
R
とすると以下のようになる。
R version 3.1.0 (2014-04-10) -- "Spring Dance"
Copyright (C) 2014 The R Foundation for Statistical Computing
Platform: x86_64-unknown-linux-gnu (64-bit)

R は、自由なソフトウェアであり、「完全に無保証」です。
 一定の条件に従えば、自由にこれを再配布することができます。
 配布条件の詳細に関しては、'license()' あるいは 'licence()' と入力してください。

R は多くの貢献者による共同プロジェクトです。
 詳しくは 'contributors()' と入力してください。
 また、R や R のパッケージを出版物で引用する際の形式については
 'citation()' と入力してください。

 'demo()' と入力すればデモをみることができます。
 'help()' とすればオンラインヘルプが出ます。
 'help.start()' で HTML ブラウザによるヘルプがみられます。
 'q()' と入力すれば R を終了します。
>
となってRが使えるようになる。Rが起動している時は>が左についている状態にある。簡単なコマンドについてはググれば色々出てくる。計算や変数への代入などにある程度慣れてから、bioconductorをインストールしてみた。とりあえず書いてある通りやってみる。

>source("http://bioconductor.org/biocLite.R")

URL 'http://www.bioconductor.org/packages/2.14/bioc/src/contrib/BiocInstaller_1.13.3.tar.gz' を試しています
Content type 'application/x-gzip' length 14185 bytes (13 Kb)

 開かれた URL

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

downloaded 13 Kb

* installing *source* package ‘BiocInstaller’ ...

** R

** inst

** preparing package for lazy loading

** help

*** installing help indices

** building package indices

** testing if installed package can be loaded

Bioconductor version 2.14 (BiocInstaller 1.13.3), ?biocLite for help

* DONE (BiocInstaller)

The downloaded source packages are in

     ‘/tmp/RtmpIpHodt/downloaded_packages’

Updating HTML index of packages in '.Library'

Making 'packages.html' ... done

Bioconductor version 2.14 (BiocInstaller 1.13.3), ?biocLite for help
>
となる。数分かかる。ChIPpeakAnnoというパッケージをダウンロードする。
>biocLite("ChIPpeakAnno")

BioC_mirror: http://bioconductor.org

Using Bioconductor version 2.14 (BiocInstaller 1.13.3), R version 3.1.0.

Installing package(s) 'ChIPpeakAnno’

以下色々続く
結構かかる5分くらい。つぎにChIPpeakAnnoを使うにはRを起動するたびに以下のコマンドを叩く必要があるらしい
>library("ChIPpeakAnno")

要求されたパッケージ biomaRt をロード中です

要求されたパッケージ IRanges をロード中です

要求されたパッケージ Biostrings をロード中です

要求されたパッケージ XVector をロード中です
10秒くらいでおわる。
次に解析手順。
MACSから出てきたピークのbedファイルをread.tableで読み込む。
それをRのパッケージで使いやすい形式にBED2RangedDataで変換する。
>SoxDay0.df <- read.table("../SoxDay0_peaks.bed", header= FALSE) 
> head(SoxDay0.df)
    V1      V2      V3          V4     V5
1 chr1 3062864 3063136 MACS_peak_1  81.45
2 chr1 3482902 3483152 MACS_peak_2  65.14
3 chr1 4150841 4151273 MACS_peak_3 222.85
4 chr1 4660129 4660269 MACS_peak_4  53.63
5 chr1 4802472 4802926 MACS_peak_5 193.69
6 chr1 4855364 4855775 MACS_peak_6 270.83 
> SoxDay0.gr <- BED2RangedData(SoxDay0.df, header= FALSE)
どんな中身かを少しだけ見るにはhead()とする。
>head(SoxDay0.gr)
RangedData with 6 rows and 2 value columns across 23 spaces
               space             ranges |      strand     score
            <factor>          <IRanges> | <character> <numeric>
MACS_peak_1        1 [3062864, 3063136] |           +     81.45
MACS_peak_2        1 [3482902, 3483152] |           +     65.14
MACS_peak_3        1 [4150841, 4151273] |           +    222.85
MACS_peak_4        1 [4660129, 4660269] |           +     53.63
MACS_peak_5        1 [4802472, 4802926] |           +    193.69
MACS_peak_6        1 [4855364, 4855775] |           +    270.83 
このフォーマットにしておくと今後色んな解析を行うのに便利らしい。ここまでのデータをセーブしておく。(save())
>save(list=ls(), file = "./bed2rangedData.rdat”)
次にTSSのデータがパッケージ内に入っているのでdata()としてそれを読み込む。annotatePeakInBatchを使ってそのTSSデータと先ほど作成したRangedDataを結びつけてEnsEMBL Gene ID、遺伝子までの距離、上流・下流の情報をアノテーションする。
>data(TSS.mouse.GRCm38)
>SoxDay0.anno <- annotatePeakInBatch(RangedData(SoxDay0.gr), AnnotationData = TSS.mouse.GRCm38,output="both") 
>head(SoxDay0.anno)
RangedData with 6 rows and 9 value columns across 23 spaces
                                    space               ranges |          peak
                                 <factor>            <IRanges> |   <character>
MACS_peak_1 ENSMUSG00000090025          1 [ 3062864,  3063136] |   MACS_peak_1
MACS_peak_10 ENSMUSG00000033740         1 [ 6458167,  6458317] |  MACS_peak_10
MACS_peak_100 ENSMUSG00000089964        1 [74824905, 74825354] | MACS_peak_100
MACS_peak_101 ENSMUSG00000082012        1 [76036406, 76036593] | MACS_peak_101
MACS_peak_106 ENSMUSG00000094946        1 [83844615, 83844958] | MACS_peak_106
MACS_peak_107 ENSMUSG00000062590        1 [86187886, 86188233] | MACS_peak_107
                                      strand            feature start_position
                                 <character>        <character>      <numeric>
MACS_peak_1 ENSMUSG00000090025             + ENSMUSG00000090025        3054233
MACS_peak_10 ENSMUSG00000033740            + ENSMUSG00000033740        6487231
MACS_peak_100 ENSMUSG00000089964           + ENSMUSG00000089964       74850895
MACS_peak_101 ENSMUSG00000082012           + ENSMUSG00000082012       75636390
MACS_peak_106 ENSMUSG00000094946           + ENSMUSG00000094946       83795912
MACS_peak_107 ENSMUSG00000062590           + ENSMUSG00000062590       86154780
                                 end_position insideFeature distancetoFeature
                                    <numeric>   <character>         <numeric>
MACS_peak_1 ENSMUSG00000090025        3054733    downstream              8631
MACS_peak_10 ENSMUSG00000033740       6860940      upstream            -29064
MACS_peak_100 ENSMUSG00000089964     74879953      upstream            -25990
MACS_peak_101 ENSMUSG00000082012     75637179    downstream            400016
MACS_peak_106 ENSMUSG00000094946     83796032    downstream             48703
MACS_peak_107 ENSMUSG00000062590     86278284        inside             33106
                                 shortestDistance fromOverlappingOrNearest
                                        <numeric>              <character>
MACS_peak_1 ENSMUSG00000090025               8131             NearestStart
MACS_peak_10 ENSMUSG00000033740             28914             NearestStart
MACS_peak_100 ENSMUSG00000089964            25541             NearestStart
MACS_peak_101 ENSMUSG00000082012           399227             NearestStart
MACS_peak_106 ENSMUSG00000094946            48583             NearestStart
MACS_peak_107 ENSMUSG00000062590            33106             NearestStart 
上流や下流の情報でまとめたテーブルを作ってみる。
>SoxDay0.anno.table <- table(as.data.frame(SoxDay0.anno)$insideFeature) 
>head(SoxDay0.anno.table)

    downstream includeFeature         inside     overlapEnd   overlapStart
          1197              1           1228             14             14
      upstream
          1209 
これを画像として出力させるには以下のようにする。
>png("SoxDay0.anno.barplot.png")
>barplot(SoxDay0.anno.table)
>dev.off()
ちょっと待つ
null device
          1 
とでる。local/work/R内にSoxDay0.anno.barplot.pngが出来てる これをscpでPCに転送
scp kosugi@サーバー名:/home/kosugi/local/work/R/SoxDay0.anno.barplot.png /Users/TK/Analysis/ChIPpeakAnno/
などとすれば見れる。こんな感じ。
次は org.Mm.eg.dbを使ってGene Ontologyを実践してみる。

2014年4月15日火曜日

ChIP-seq解析の流れ

大体の解析の流れはこんな感じ。
  1. 解析に必要なツールのインストール、コントロールを含む.fastqまたは.SRAファイル (SRA Toolkitでfastqへ変換する必要有)をダウンロード
  2. fastqファイルの中身をFastQCで確認して、FASTX-Toolkitのfastq_quality_trimmer, fastq_quality_filterでクリーニング
  3. Bowtie2でマッピング(リファレンスゲノムをダウンロードしておく必要有)
  4. Samtoolsでバイナリー化
  5. 必要が有ればbedtoolsのintersectBedを使ってリピート配列の除去
  6. MACSでピークコール
  7. R+Bioconductorを用いて高次解析
    • ChIPpeakAnno
      • アノテーション付加(近傍TSSまでの距離・EnsEMBL Gene ID・上流下流)
    • org.Mm.eg.db
      • Gene ontology
    • BSgenome
      • ピーク領域の配列を得る
    • rGADEM
      • de novo motif search
    • MotIV
      • 既存データとの類似性チェック
Rを使った高次解析については二階堂さんのHPが参考になる。とりあえず、1-6に関して実際にやった流れを記す。長いのと、Rの使い勝手がすこしちがうので7に関しては別に記す。

1については既に終わっている。

2から始める。(workディレクトリにfastqファイルがある前提でworkで作業を行う。)

とりあえずfastqc
fastqc -t 8 Sox2_ChIP-seq_Oct3_4_KODay0.fastq -o ./
しばらく待つと、Sox2_ChIP-seq_Oct3_4_KODay0.fastqcというディレクトリとそのzipファイルが出来る。それをディレクトリごとローカルに転送する。ローカル側から
scp -r kosugi@サーバーアカウント:/home/kosugi/local/work/Sox2_ChIP-seq_Oct3_4_KODay0.fastqc.zip /Users/TK/fastqc/  
すると指定した場所にディレクトリごと転送される。その中にfastqc_report.htmlがあるのでブラウザーにドラッグアンドドロップする。
こんな感じになる。3'末のクオリティーが低いものが多い事がわかる。FASTQ形式についてwikiで調べればこの数値についてある程度解説が見られる。とりあえず今回は3'末端からクオリティー値を20未満のものをトリミングして、30bp未満のreadを破棄。さらに、80%以上がクオリティー値が20以上ものを抽出する事にした。オプションの詳細はFASTX-Toolkitのusageを確認する。パイプ(|)でつなげてそのまま実行する。
fastq_quality_trimmer -t 20 -l 30 -Q 33 -i Sox2_ChIP-seq_Oct3_4_KODay0.fastq | fastq_quality_filter -q 20 -p 80 -Q 33 -o ./Sox2_ChIP-seq_Oct3_4_KODay0.clean.fastq
このfastqファイルを先ほどと同様にfastqcしてfastqc_report.htmlをみる。
クリーニングされた事が見た目でわかる。
それぞれのfastqファイルをwc -lすると行数がわかる
wc -l Sox2_ChIP-seq_Oct3_4_KODay0.fastq
62931980 Sox2_ChIP-seq_Oct3_4_KODay0.fastq
wc -l Sox2_ChIP-seq_Oct3_4_KODay0.clean.fastq
59910120 Sox2_ChIP-seq_Oct3_4_KODay0.clean.fastq
fastqは4行で1readなので大体0.8M read程除去されたことがわかる。

3のBowite2によるマッピング

とりあえずコマンド。10分程かかる。read数やマシーン性能に依存する。
bowtie2 -t -p 8 -x ../bowtie2-2.2.1/indexes/mm10 -U Sox2_ChIP-seq_Oct3_4_KODay0.clean.fastq -S SoxDay0.sam


オプションの意味は 時間:-t
コア数:-p
インデックスのパス:-x
fastqのパス:-U
出力名:-S

4のSamtoolsによるバイナリ化

とりあえずコマンド。
samtools view -Sb SoxDay0.sam > SoxDay0.bam
[samopen] SAM header is present: 66 sequences.等と出て、数分待つ。 オプションの意味は
-S:インプットがsamファイル 
-b:アウトプットがbamファイル
bowtieで吐き出されるのは人間が理解できるsamファイルであり、これ以降のステップではバイナリー化しないと認識してくれないらしくこのステップが必要になる。つぎにソートを行う。igvというツールを使う時にindexを作る必要が有り、そのためには順序よく並んでないと怒られるらしい。
samtools sort SoxDay0.bam SoxDay0.sorted
すると、SoxDay0.sorted.bamというのがつくられる。次にindexを作る。
samtools index SoxDay0.sorted.bam
数分待つとSoxDay0.sorted.bam.baiファイルが出来る。今回はigvでマップされたものを確認しない。次にマップ率を出す。flagstatを使う。
samtools flagstat SoxDay0.sorted.bam > SoxDay0.sorted.bam.summary.txt 
less SoxDay0.sorted.bam.summary.txt
14977530 + 0 in total (QC-passed reads + QC-failed reads)
0 + 0 duplicates
14529891 + 0 mapped (97.01%:-nan%)
0 + 0 paired in sequencing
0 + 0 read1
0 + 0 read2
0 + 0 properly paired (-nan%:-nan%)
0 + 0 with itself and mate mapped
0 + 0 singletons (-nan%:-nan%)
0 + 0 with mate mapped to a different chr
0 + 0 with mate mapped to a different chr (mapQ>=5)

97%がマップされている。

5のintersectBedを使ってリピート配列の除去

必要が有ればrepeat region にマップされたreadsを除去する。これには準備が必要でUCSC Table Browserからデータをダウンロードしておく。以下のように設定して、get outputしたらダウンロードできる。

local/shareに転送した。これで準備完了
scp /Users/TK/Downloads/repeats_mm10.bed kosugi@サーバー名:/home/kosugi/local/share
以下がrepeat regionを除去するコマンド。
intersectBed -abam SoxDay0.sorted.bam -b ../share/repeats_mm10.bed -v > SoxDay0.rmRepeat.bam
これで除去できた。 6のMACSによるピークコール コントロールのサンプルについても同様の作業を行って、bamファイルが2つある状況にする。準備ができたら以下のコマンドを叩く。数分待つ。
macs14 -t SoxDay0.sorted.bam -c SoxKODay2.sorted.bam -f bam -g mm -n SoxDay0
オプションの意味は -t:sampleファイル -c:controlファイル -f:ファイル形式 -g:ゲノムサイズ -n:出力名 終わると、色んなファイルが出力される。
SoxDay0_negative_peaks.xls
SoxDay0_summits.bed
SoxDay0_peaks.bed
SoxDay0_peaks.xls
SoxDay0_model.r
例えば、SoxDay0_peaks.bedの中身はこんな感じ。
head SoxDay0_peaks.bed
chr1 3062864 3063136 MACS_peak_1 81.45
chr1 3482902 3483152 MACS_peak_2 65.14
chr1 4150841 4151273 MACS_peak_3 222.85
chr1 4660129 4660269 MACS_peak_4 53.63
chr1 4802472 4802926 MACS_peak_5 193.69
chr1 4855364 4855775 MACS_peak_6 270.83
chr1 4972272 4972486 MACS_peak_7 71.27
chr1 4972497 4972701 MACS_peak_8 84.01
chr1 6448821 6449289 MACS_peak_9 126.84
chr1 6458167 6458317 MACS_peak_10 62.50
以下のように、ローカル側からxlsファイルを転送すればエクセルで見られる。3000ヶ所くらいピークが見つかった。
scp kosugi@サーバー名:/home/kosugi/local/work/SoxDay0_peaks.xls /Users/TK/MACS/ 
wiggleファイルが欲しい場合はMACSを実行する時に -wとオプションをつける。(結構時間がかかる数十分)

2014年4月13日日曜日

Rのインストール

解析をするにあたって、Rという統計解析ツールがあると便利そうなのでダウンロード、インストールを行う。root権限無しの画面無しの条件にインストールする情報がすくない状況だったので、やや試行錯誤した。
HPはここ。Source Code for all Platformsの最新版である、 R-3.1.0.tar.gzをダウンロードする。
ターミナルを開いてsrcに移動し、wget "と記述し、R-3.1.0.tar.gzをドラッグアンドドロップしたら"と記述する。
cd src/
wget "http://cran.ism.ac.jp/src/base/R-3/R-3.1.0.tar.gz"
次に解凍。.tar.gzなのでtar zxvfで解凍。-Cで解凍先にlocalを指定する。
tar zxvf R-3.1.0.tar.gz -C ../local/
localのR-3.1.0に移動してREADMEやINSTALL、doc/manual/R-admin.texiをみてもなんかいまいちわからない。なのでここを参照する。

2.4 Instrationをみると、root権限が無い状況で/usr/local/bin/にインストールできない時には、prefixが使える。./configure optionの項目を見ると、計算機を使うような画面無しの条件だと--with-x=noオプションを追加して./configureをすればいいらしい。
./configure --prefix=/home/kosugi/local --with-x=no
makeする。
make
makeに結構時間かかる。5分くらい。そのご以下のコマンドを叩く。
make check
make install
これが終わると、R-3.1.0/bin以下に様々なプログラムがある事がわかる。 なのでパスを通す。
export PATH=/home/kosugi/local/R-3.1.0/bin/:$PATH
R -hと叩くと
R -h

Usage: R [options] [< infile] [> outfile]
   or: R CMD command [arguments]

Start R, a system for statistical computation and graphics, with the
specified options, or invoke an R tool via the 'R CMD' interface.

Options:
以下オプションの説明
となったので問題なく使えそう。vimで.bashrcにパスを記述しておく。ファイルが開いたらiを押すと書き込める。以下を一番下に追記する。追記したらescキーを押して:wq
vim ~/.bashrc
##R
export PATH=/home/kosugi/local/R-3.1.0/bin/:$PATH
再ログインして
echo $PATH
/home/kosugi/local/R-3.1.0/bin/:~その他のパス
となる。R -hとコマンドを叩くと先ほどと同様にusageがでるのでこれでRが使えるようになったぽい。他のツールと同様でprefixでどうにかなるが、画面無しのオプション--with-x=noこれをつける事がエラーを出さないポイントか・・・。 ただ全然使い方わからない。統合TVで検索したらR+Bioconductorを使ったNGS解析1R+Bioconductorを使ったNGS解析2というものが動画でUPされていてここで紹介されている二階堂さんのLearning Rのページが非常に参考になる。まずはこれを参考にしてRというツールを使う事に慣れて、ChIP-seqのピーク検出後の高次解析を行えるようにする。

bedtoolsのインストール

リピート配列の部分を除去したほうがde novo motif searchをする時に繰り返し配列に含まれる配列があまり出てこなくなるらしい。
その為にbedtoolsのintersectBedというプログラムを使うと簡単に除去できるらしい
とりあえずダウンロードする。HPはここ
release can found hereと書かれてる右側のリンクへ移動する。


ターミナルを開いてsrcに移動しwget"と記述してbedtools-2.19.1.tar.gzのリンクをドラッグアンドドロップして"とする
cd src/
wget "https://github.com/arq5x/bedtools2/releases/download/v2.19.1/bedtools-2.19.1.tar.gz"
localに解凍する。.tar.gzなのでtar zxvfで解凍。-Cで解凍先を指定する。
tar zxvf bedtools-2.19.1.tar.gz -C ../local/

localのbedtools-2.19.1に移動してREADME.mdをみてもインストールに関する事が書いていない。のでここにあるマニュアルをダウンロードしてみるとinstallという項目がある。
make clean
make all
とすればよいらしい
make clean
* Cleaning-up BamTools API
* Cleaning up.
となって
make all
Building BEDTools:
=========================================================
DETECTED_VERSION = v2.19.1
CURRENT_VERSION  = v2.19.0-15-g733c84b
Updating version file.
 * Creating BamTools API
...以下100行ぐらいコンパイルが続く
これが終わると、bedtools-2.19.1/bin以下に様々なプログラムがある事がわかる。 なのでパスを通す。
export PATH=/home/kosugi/local/bedtools2-2.19.1/bin/:$PATH
intersectBedと叩くと
intersectBed

Tool:    bedtools intersect (aka intersectBed)
Version: v2.19.1
Summary: Report overlaps between two feature files.

Usage:   bedtools intersect [OPTIONS] -a  -b 
以下オプションの説明
となったので問題なく使えそう。vimで.bashrcにパスを記述しておく。ファイルが開いたらiを押すと書き込める。以下を一番下に追記する。追記したらescキーを押して:wq
vim ~/.bashrc
##bedtools
export PATH=/home/kosugi/local/bedtools2-2.19.1/bin/:$PATH
再ログインして
echo $PATH
/home/kosugi/local/bedtools2-2.19.1/bin/:~その他のパス
となっているので問題なくつかえそうだ。

2014年4月11日金曜日

FASTX-Toolkitのインストール

ちょっと試していくうちにシーケンスリードをある一定基準でトリミングしたり、破棄したりするツールが必要だという事に気づいた。とりあえずFASTX-Toolkitというのが便利らしいのでダウンロードする。HPはここ
Download & Installationをクリックする。
Fastx-toolkit version 0.0.13 requires libgtextutils-0.6 (available here for download)とかかれているのでlibgtextutilsも必要。とりあえず最新版であるlibgtextutils-0.7.tar.gzをまず先にインストールする。 ターミナルを開いて
cd src/
wget "https://github.com/agordon/fastx_toolkit/releases/download/0.0.14/fastx_toolkit-0.0.14.tar.bz2"
localに解凍する。
tar zxvf libgtextutils-0.7.tar.gz -C ../local/
ディレクトリを移動して、INSTALLファイルを見るとIt is recommended to use the following options: ./configure --prefix=/boot/common といった記述が有るので、./configureとコンパイル、インストールをする。
./configure --prefix=/home/kosugi/local/
make
make install
とすると/local/lib以下にlibgtextutils-0.7関連のファイルの存在が確認できる。これでいいのだろうか・・・ とりあえず、本題の最新のfastx_toolkit-0.0.14.tar.bz2をダウンロードして、localに解凍する。
wget "https://github.com/agordon/fastx_toolkit/releases/download/0.0.14/fastx_toolkit-0.0.14.tar.bz2"
tar jxvf fastx_toolkit-0.0.14.tar.bz2 -C ../local/
localにできたディレクトリに移動して./configureとコンパイル、インストールをする。
./configure --prefix=/home/kosugi/local/
とすると
checking for GTEXTUTILS... configure: error: Package requirements (gtextutils) were not met:

No package 'gtextutils' found

Consider adjusting the PKG_CONFIG_PATH environment variable if you
installed software in a non-standard prefix.

Alternatively, you may set the environment variables GTEXTUTILS_CFLAGS
and GTEXTUTILS_LIBS to avoid the need to call pkg-config.
See the pkg-config man page for more details.
と怒られる。どうやらさきほど出来た/home/kosugi/local/lib/にpkgconfigディレクトリがあるのでそれにPKG_CONFIG_PATHのパスを通せという事のようだ。 vimで.bashrcに記述。ファイルが開いたらiを押すと書き込める。以下を一番下に追記する。追記したらescキーを押して:wq
vim ~/.bashrc
export PKG_CONFIG_PATH=/home/kosugi/local/lib/pkgconfig/:$PKG_CONFIG_PATH
一度ターミナルを切って再接続して確認する。
echo $PKG_CONFIG_PATH
/home/kosugi/local/lib/pkgconfig/:
多分これでよいはず。fastx_toolkit-0.0.14ディレクトリにもどって
./configure --prefix=/home/kosugi/local/
make
make install
とすると、local/bin以下にたくさんのツールが有る事が確認できる。すでに、/home/kosugi/local/binにはパスが通っているので、たとえばホームディレクトリで
fastq_quality_trimmer -h
usage: fastq_quality_trimmer [-h] [-v] [-t N] [-l N] [-z] [-i INFILE] [-o OUTFILE]
Part of FASTX Toolkit 0.0.14 by A. Gordon (assafgordon@gmail.com)
以下オプションの説明
無事インストールがすんだ。./configure時に特にprefixする必要が無い状況であれば大して困らないと思うがroot権限が無いような状況でもなんとかインストールできる事がわかった。

2014年4月10日木曜日

ChIP-seqデータのダウンロード

Adachi K, Nikaido I, Ohta H, Ohtsuka S et al. Context-dependent wiring of Sox2 regulatory networks for self-renewal of embryonic and trophoblast stem cells. Mol Cell 2013 Nov 7;52(3):380-92. PMID: 24120664のデータを用いて、解析を行ってみる。論文を見るとデータがある場所はGSE28455ということで、ページの下部にChIPのデータはGSE28453にまとまっている模様

リンク先は以下のようになっていてSamplesのmoreをクリックすると全体がでてくる。

今回は練習がてら、GSM703186 Sox2 ChIP-seq, Oct3/4KO Day0とコントロールにGSM703191 Sox2 ChIP-seq, Sox2KO Day2を使う。GSM703186 Sox2 ChIP-seq, Oct3/4KO Day0のデータをダウンロードするのでGSM703186をクリックし、下部へ移動すると

すでに著者らが解析済みのファイル(BW : BigWig形式)とシーケンサーから吐き出されただけのデータがあり、今回はSRX057755をダウンロードする。ただ、NCBIからこのままダウロードする事も可能だが、前回も言った通り(SRA Toolkitのインストール).SRA形式からいちいち変換しなくてはならない。しかし、DDBJのDRAならfastq形式であると言う事を聞いたので、DRAからダウンロードする。SRX057755とググればトップにでてきた。
リンクをクリックすると以下のようなページに飛ぶ。
FASTQと書いてあるところからダウンロードできるのだが、safariを使っている場合は、うまく出来ないので、ここからはターミナルで行う。ターミナルを開いてsrcに移動する。lftpと打ってから、FASTQのリンクをドラッグアンドドロップする。
lftp ftp://ftp.ddbj.nig.ac.jp/ddbj_database/dra/fastq/SRA035/SRA035367/SRX057755
とすると
cd 成功、cwd=/ddbj_database/dra/fastq/SRA035/SRA035367/SRX057755

lftp ftp.ddbj.nig.ac.jp:/ddbj_database/dra/fastq/SRA035/SRA035367/SRX057755>
となるのでlsで中身を確認する。
lftp ftp.ddbj.nig.ac.jp:/ddbj_database/dra/fastq/SRA035/SRA035367/SRX057755>ls
-rw-r--r--   1 51005    51005    723342854 Nov  1 05:09 SRR185892.fastq.bz2
fastqファイルが見つかったのでこれを転送する。
lftp ftp.ddbj.nig.ac.jp:/ddbj_database/dra/fastq/SRA035/SRA035367/SRX057755>get SRR185892.fastq.bz2
とすると転送が始まる。wgetでは無い事に注意。転送がおわったらexitすればftpから戻って来れる。 あとは解凍するだけ。~/local/shareというディレクトリを作成してここに解凍する事にする。だけど、bunzip2で解凍先を指定する方法がわからなかったので、とりあえずsrcで解凍ほんで~/local/shareにコピーしてSox2_ChIP-seq_Oct3_4_KODay0.fastqにリネーム
cd src/
bunzip2 -d SRR185892.fastq.bz2
cp ../local/share
cd ../local/share/
mv SRR185892.fastq Sox2_ChIP-seq_Oct3_4_KODay0.fastq
ls
Sox2_ChIP-seq_Oct3_4_KODay0.fastq
同様のステップをSox2 ChIP-seq, Sox2KO Day2のデータも行った。 これでfastqcにかけられるところまできた。

2014年4月6日日曜日

SRA Toolkitのインストール

NCBIにupされているシーケンスデータは.sraフォーマットであり.fastqフォーマットに変換する必要があるらしい。なのでそれに必要なツールをダウンロード、インストールする。
HPはここhttp://trace.ncbi.nlm.nih.gov/Traces/sra/sra.cgi?view=software
環境に合わせてCentOS Linux 64 bit architectureをダウンロードする。
ターミナルを開いてsrcに移動してドラッグアンドドロップ
cd src/
wget "http://ftp-trace.ncbi.nlm.nih.gov/sra/sdk/2.3.5/sratoolkit.2.3.5-centos_linux64.tar.gz"
--2014-04-06 17:13:46--  http://ftp-trace.ncbi.nlm.nih.gov/sra/sdk/2.3.5/sratoolkit.2.3.5-centos_linux64.tar.gz
ftp-trace.ncbi.nlm.nih.gov をDNSに問いあわせています... 130.14.250.13, 2607:f220:41e:250::12
ftp-trace.ncbi.nlm.nih.gov|130.14.250.13|:80 に接続しています... 接続しました。
HTTP による接続要求を送信しました、応答を待っています... 200 OK
長さ: 64933999 (62M) [application/x-gzip]
`sratoolkit.2.3.5-centos_linux64.tar.gz' に保存中

100%[=================================================>] 64,933,999  3.74M/s 時間 17s

2014-04-06 17:14:06 (3.56 MB/s) - `sratoolkit.2.3.5-centos_linux64.tar.gz' へ保存完了 [64933999/64933999]
次に解凍 .tar.gzなのでtar zxvfで解凍。-Cで解凍先を指定する。
tar zxvf sratoolkit.2.3.5-centos_linux64.tar.gz -C ../local/
解凍されたsratoolkit.2.3.5-centos_linux64ディレクトリのbin以下に実行ファイルが入っている。その中のfastq-dumpが目的の実行ファイル。あとは/home/kosugi/local/sratoolkit.2.3.5-centos_linux64/binにパスを通す。
export PATH=/home/kosugi/local/sratoolkit.2.3.5-centos_linux64/bin/:$PATH
ほんで、home ディレクトリに戻ってfastq-dumpでうまくいくか確認
fastq-dump
Usage:
  fastq-dump [options]  [...]
  fastq-dump [options] 

Use option --help for more information

fastq-dump : 2.3.5
あとは同じ事を./bashrcに記述して再ログインしてもホームディレクトリで実行できるか確認する。
vim .bashrc
ファイルが開いたらiを押すと書き込める。以下を一番下に追記する。追記したらescキーを押して:wq
## sratoolkit
export PATH=/home/kosugi/local/sratoolkit.2.3.5-centos_linux64/bin/:$PATH  :$PATHを後ろに持ってくると今まで書かれたパスの前に記述という事らしい。
念のため再ログインしてwhichでどのパスでプログラムが実行されている確認する。
which fastq-dump
~/local/sratoolkit.2.3.5-centos_linux64/bin/fastq-dump
これでSRAからダウンロードしたファイルをfastqファイルに変換できる。

MACSのインストール

HPはここhttp://liulab.dfci.harvard.edu/MACS/
右上のdownloadをクリック
ターミナルを開いてsrcに移動、Download version 1.4.2をドラッグアンドドロップ
cd src/
wget "https://github.com/downloads/taoliu/MACS/MACS-1.4.2-1.tar.gz"
--2014-04-06 13:55:08--  https://github.com/downloads/taoliu/MACS/MACS-1.4.2-1.tar.gz
github.com をDNSに問いあわせています... 192.30.252.128
github.com|192.30.252.128|:443 に接続しています... 接続しました。
HTTP による接続要求を送信しました、応答を待っています... 302 Found
場所: http://cloud.github.com/downloads/taoliu/MACS/MACS-1.4.2-1.tar.gz [続く]
--2014-04-06 13:55:09--  http://cloud.github.com/downloads/taoliu/MACS/MACS-1.4.2-1.tar.gz
cloud.github.com をDNSに問いあわせています... 54.230.27.179, 205.251.209.90, 205.251.209.218, ...
cloud.github.com|54.230.27.179|:80 に接続しています... 接続しました。
HTTP による接続要求を送信しました、応答を待っています... 200 OK
長さ: 66569 (65K) [application/gzip]
`MACS-1.4.2-1.tar.gz' に保存中

100%[======================================>] 66,569      66.4K/s 時間 1.0s

2014-04-06 13:55:11 (66.4 KB/s) - `MACS-1.4.2-1.tar.gz' へ保存完了 [66569/66569]
次に解凍 .tar.gzなのでtar jxvfで解凍。-Cで解凍先を指定する。
tar zxvf MACS-1.4.2-1.tar.gz -C ../local/
解凍されたMACS-1.4.2-1ディレクトリに移動して、INSTALLファイルをみる
cd ../local/MACS-1.4.2-1
less INSTALL

* Prerequisite
Python version must be equal to 2.6 or 2.7 to run MACS. We recommend
using the version 2.6.5.

* Install under Debian or Ubuntu Linux system
The most convenient way to install MACS is through Debian APT system,
so that it can be perfectly integrated in the Python environment of
your operation system. You can easily manage the package, and the
uninstall is much easier. Download the deb package from MACS download
page, and type this in the commend line:

$ dpkg -i macs_1.4.0.deb

To uninstall, type:
$ dpkg -r macs_1.4.0

MACS uses Python's distutils tools for source installations. To
install a source distribution of MACS, unpack the distribution tarball
and open up a command terminal. Go to the directory where you unpacked
MACS, and simply run the install script :

$ python setup.py install

By default, the script will install python library and executable
codes globally, which means you need to be root or administrator of
the machine so as to complete the installation. Please contact the
administrator of that machine if you want their help. If you need to
provide a nonstandard install prefix, or any other nonstandard
options, you can provide many command line options to the install
script. Use the –help option to see a brief list of available options:

$ python setup.py --help

For example, if I want to install everything under my own HOME
directory, use this command:

$ python setup.py install --prefix /home/taoliu/

ようするにpython2.6か2.7が必要。Debian or Ubuntu Linux systemではdpkgというコマンドが使えるようだが、今の環境はCentOSなのでpython setup.py installを使う。prefixを使えば指定先にインストールできる。という事らしい。
python --version
Python 2.7.3
となったので問題ないかな。
python setup.py install --prefix /home/kosugi/local
と叩くと以下のようになり、Python2.7を使っているという事がわかる。
/usr/local/Python-2.7.3/lib/python2.7/distutils/dist.py:267: UserWarning: Unknown distribution option: 'console'
  warnings.warn(msg)
/usr/local/Python-2.7.3/lib/python2.7/distutils/dist.py:267: UserWarning: Unknown distribution option: 'app'
  warnings.warn(msg)
running install
running build
running build_py
creating build
creating build/lib
creating build/lib/MACS14
copying lib/__init__.py -> build/lib/MACS14
copying lib/Constants.py -> build/lib/MACS14
copying lib/OptValidator.py -> build/lib/MACS14
copying lib/OutputWriter.py -> build/lib/MACS14
copying lib/PeakDetect.py -> build/lib/MACS14
copying lib/PeakModel.py -> build/lib/MACS14
copying lib/Prob.py -> build/lib/MACS14
creating build/lib/MACS14/IO
copying lib/IO/__init__.py -> build/lib/MACS14/IO
copying lib/IO/bedGraphIO.py -> build/lib/MACS14/IO
copying lib/IO/BinKeeper.py -> build/lib/MACS14/IO
copying lib/IO/FeatIO.py -> build/lib/MACS14/IO
copying lib/IO/Parser.py -> build/lib/MACS14/IO
copying lib/IO/WiggleIO.py -> build/lib/MACS14/IO
running build_scripts
creating build/scripts-2.7
copying and adjusting bin/macs14 -> build/scripts-2.7
copying and adjusting bin/elandmulti2bed -> build/scripts-2.7
copying and adjusting bin/elandresult2bed -> build/scripts-2.7
copying and adjusting bin/elandexport2bed -> build/scripts-2.7
copying and adjusting bin/sam2bed -> build/scripts-2.7
copying and adjusting bin/wignorm -> build/scripts-2.7
changing mode of build/scripts-2.7/macs14 from 644 to 755
changing mode of build/scripts-2.7/elandmulti2bed from 644 to 755
changing mode of build/scripts-2.7/elandresult2bed from 644 to 755
changing mode of build/scripts-2.7/elandexport2bed from 644 to 755
changing mode of build/scripts-2.7/sam2bed from 644 to 755
changing mode of build/scripts-2.7/wignorm from 644 to 755
running install_lib
creating /home/kosugi/local/lib
creating /home/kosugi/local/lib/python2.7
creating /home/kosugi/local/lib/python2.7/site-packages
creating /home/kosugi/local/lib/python2.7/site-packages/MACS14
copying build/lib/MACS14/__init__.py -> /home/kosugi/local/lib/python2.7/site-packages/MACS14
copying build/lib/MACS14/Constants.py -> /home/kosugi/local/lib/python2.7/site-packages/MACS14
copying build/lib/MACS14/OptValidator.py -> /home/kosugi/local/lib/python2.7/site-packages/MACS14
copying build/lib/MACS14/OutputWriter.py -> /home/kosugi/local/lib/python2.7/site-packages/MACS14
copying build/lib/MACS14/PeakDetect.py -> /home/kosugi/local/lib/python2.7/site-packages/MACS14
copying build/lib/MACS14/PeakModel.py -> /home/kosugi/local/lib/python2.7/site-packages/MACS14
copying build/lib/MACS14/Prob.py -> /home/kosugi/local/lib/python2.7/site-packages/MACS14
creating /home/kosugi/local/lib/python2.7/site-packages/MACS14/IO
copying build/lib/MACS14/IO/__init__.py -> /home/kosugi/local/lib/python2.7/site-packages/MACS14/IO
copying build/lib/MACS14/IO/bedGraphIO.py -> /home/kosugi/local/lib/python2.7/site-packages/MACS14/IO
copying build/lib/MACS14/IO/BinKeeper.py -> /home/kosugi/local/lib/python2.7/site-packages/MACS14/IO
copying build/lib/MACS14/IO/FeatIO.py -> /home/kosugi/local/lib/python2.7/site-packages/MACS14/IO
copying build/lib/MACS14/IO/Parser.py -> /home/kosugi/local/lib/python2.7/site-packages/MACS14/IO
copying build/lib/MACS14/IO/WiggleIO.py -> /home/kosugi/local/lib/python2.7/site-packages/MACS14/IO
byte-compiling /home/kosugi/local/lib/python2.7/site-packages/MACS14/__init__.py to __init__.pyc
byte-compiling /home/kosugi/local/lib/python2.7/site-packages/MACS14/Constants.py to Constants.pyc
byte-compiling /home/kosugi/local/lib/python2.7/site-packages/MACS14/OptValidator.py to OptValidator.pyc
byte-compiling /home/kosugi/local/lib/python2.7/site-packages/MACS14/OutputWriter.py to OutputWriter.pyc
byte-compiling /home/kosugi/local/lib/python2.7/site-packages/MACS14/PeakDetect.py to PeakDetect.pyc
byte-compiling /home/kosugi/local/lib/python2.7/site-packages/MACS14/PeakModel.py to PeakModel.pyc
byte-compiling /home/kosugi/local/lib/python2.7/site-packages/MACS14/Prob.py to Prob.pyc
byte-compiling /home/kosugi/local/lib/python2.7/site-packages/MACS14/IO/__init__.py to __init__.pyc
byte-compiling /home/kosugi/local/lib/python2.7/site-packages/MACS14/IO/bedGraphIO.py to bedGraphIO.pyc
byte-compiling /home/kosugi/local/lib/python2.7/site-packages/MACS14/IO/BinKeeper.py to BinKeeper.pyc
byte-compiling /home/kosugi/local/lib/python2.7/site-packages/MACS14/IO/FeatIO.py to FeatIO.pyc
byte-compiling /home/kosugi/local/lib/python2.7/site-packages/MACS14/IO/Parser.py to Parser.pyc
byte-compiling /home/kosugi/local/lib/python2.7/site-packages/MACS14/IO/WiggleIO.py to WiggleIO.pyc
running install_scripts
copying build/scripts-2.7/macs14 -> /home/kosugi/local/bin
copying build/scripts-2.7/elandmulti2bed -> /home/kosugi/local/bin
copying build/scripts-2.7/elandresult2bed -> /home/kosugi/local/bin
copying build/scripts-2.7/elandexport2bed -> /home/kosugi/local/bin
copying build/scripts-2.7/sam2bed -> /home/kosugi/local/bin
copying build/scripts-2.7/wignorm -> /home/kosugi/local/bin
changing mode of /home/kosugi/local/bin/macs14 to 755
changing mode of /home/kosugi/local/bin/elandmulti2bed to 755
changing mode of /home/kosugi/local/bin/elandresult2bed to 755
changing mode of /home/kosugi/local/bin/elandexport2bed to 755
changing mode of /home/kosugi/local/bin/sam2bed to 755
changing mode of /home/kosugi/local/bin/wignorm to 755
running install_egg_info
Writing /home/kosugi/local/lib/python2.7/site-packages/MACS-1.4.2-py2.7.egg-info
今回は、/home/kosugi/local/binにemacs14という実行ファイルとその他の実行ファイルが作られている。
elandexport2bed
elandmulti2bed
elandresult2bed
macs14
sam2bed
wignorm
INSTALLファイルにはOn Linux, using bash, I include the new value in my PYTHONPATH by adding this line to my ~/.bashrc となっているので
export PYTHONPATH=/home/kosugi/local/lib/python2.7/site-packages:$PYTHONPATH
と記述し再ログインして、macs14とコマンドを叩くと
Usage: macs14 <-t data-blogger-escaped-tfile=""> [-n name] [-g genomesize] [options]
Example: macs14 -t ChIP.bam -c Control.bam -f BAM -g h -n test -w --call-subpeaks
macs14 -- Model-based Analysis for ChIP-Sequencing
Options:
以下オプションの説明
となったのでインストール出来た。ヽ(゚∀゚)メ(゚∀゚)メ(゚∀゚)ノ 今回は、インストール先を自分の好きな場所にする事が出来た。--prefix 指定先とすると自動的にbinとlibというディレクトリが作られるようだ。 前もってbinを作ってPATHも通しておいた(ディレクトリの作成)ので、順調に進んだかな。ちゃんとreadmeやinstallといったファイルを見ながら確認してやれば結構敷居低いと感じた。これでChIP-seqの解析できるかな

2014年4月5日土曜日

Samtoolsのインストール

HPはここhttp://samtools.sourceforge.net

真ん中らへんにあるdownload pageをクリック
 tabixの上にあるsamtoolsをクリック
0.1.19をクリック
ターミナルを起動して、srcに移動する。 samtools-0.1.19.tar.bz2をターミナルにドラッグアンドドロップ。末尾に/downloadとついていたら消す。

cd src/
wget "http://sourceforge.net/projects/samtools/files/samtools/0.1.19/samtools-0.1.19.tar.bz2"
--2014-04-05 17:03:23--  http://sourceforge.net/projects/samtools/files/samtools/0.1.19/samtools-0.1.19.tar.bz2
sourceforge.net をDNSに問いあわせています... 216.34.181.60
sourceforge.net|216.34.181.60|:80 に接続しています... 接続しました。
HTTP による接続要求を送信しました、応答を待っています... 302 Found
場所: http://sourceforge.net/projects/samtools/files/samtools/0.1.19/samtools-0.1.19.tar.bz2/download [続く]
--2014-04-05 17:03:24--  http://sourceforge.net/projects/samtools/files/samtools/0.1.19/samtools-0.1.19.tar.bz2/download
sourceforge.net|216.34.181.60|:80 に接続しています... 接続しました。
HTTP による接続要求を送信しました、応答を待っています... 302 Found
場所: http://downloads.sourceforge.net/project/samtools/samtools/0.1.19/samtools-0.1.19.tar.bz2?r=&ts=1396685026&use_mirror=cznic [続く]
--2014-04-05 17:03:24--  http://downloads.sourceforge.net/project/samtools/samtools/0.1.19/samtools-0.1.19.tar.bz2?r=&ts=1396685026&use_mirror=cznic
downloads.sourceforge.net をDNSに問いあわせています... 216.34.181.59
downloads.sourceforge.net|216.34.181.59|:80 に接続しています... 接続しました。
HTTP による接続要求を送信しました、応答を待っています... 302 Found
場所: http://cznic.dl.sourceforge.net/project/samtools/samtools/0.1.19/samtools-0.1.19.tar.bz2 [続く]
--2014-04-05 17:03:25--  http://cznic.dl.sourceforge.net/project/samtools/samtools/0.1.19/samtools-0.1.19.tar.bz2
cznic.dl.sourceforge.net をDNSに問いあわせています... 217.31.202.30, 2001:1488:ffff::30
cznic.dl.sourceforge.net|217.31.202.30|:80 に接続しています... 接続しました。
HTTP による接続要求を送信しました、応答を待っています... 200 OK
長さ: 514507 (502K) [application/octet-stream]
`samtools-0.1.19.tar.bz2' に保存中

100%[======================================>] 514,507      108K/s 時間 4.6s

2014-04-05 17:03:30 (108 KB/s) - `samtools-0.1.19.tar.bz2' へ保存完了 [514507/514507]
次に解凍 .tar.bz2なのでtar jxvfで解凍。-Cで解凍先を指定する。
tar jxvf samtools-0.1.19.tar.bz2 -C ../local
解凍されたsamtools-0.1.19ディレクトリに移動。
cd ../local/samtools-0.1.19
lsで確認するといろんな種類のファイルが入っている。INSTALLファイルを見るとSAMtools depends on the zlib library . Version 1.2.3+ is preferred and with 1.2.3+ you can compile razip and use it to compress a FASTA file.とかなんとか書かれているのでよくわからんzlib libraryがサーバー内にインストールされてるか確認する。
rpm -qa | grep lib
RPMのパッケージについてインストールされているものを確認できるらしい・・・| grep libはその中でlibを含んでいるものを抽出という意味。
jzlib-1.0.7-7.5.el6.x86_64
zlib-devel-1.2.3-29.el6.x86_64
jzlib-javadoc-1.0.7-7.5.el6.x86_64
zlib-1.2.3-29.el6.x86_64
jzlib-demo-1.0.7-7.5.el6.x86_64
zlib-static-1.2.3-29.el6.x86_64
zlib-1.2.3-29.el6.i686
となるので要求のものはインストールされている模様。とりあえずINSTALLをみるとmakeってコマンド打てばいいみたい。
make
lsでみると、さらにファイルが増えていて、samtoolsという実行権限がついたファイルを確認できる。samtoolsは-hとしても意味わからんと怒られるので 単純にsamtoolsと叩く
samtools

Program: samtools (Tools for alignments in the SAM format)
Version: 0.1.19-44428cd

Usage:   samtools  [options]
以下オプションの説明
後はパスを指定する。
export PATH=/home/kosugi/local/samtools-0.1.19/:$PATH
としてディレクトリ移動してもコマンドが動く事を確認できた。 あとは.bashrcに記述。
vim .bashrc
ファイルが開いたらiを押すと書き込める。以下を一番下に追記する。追記したらescキーを押して:wq
##Samtools
export PATH=/home/kosugi/local/samtools-0.1.19/:$PATH :$PATHを後ろに持ってくると今まで書かれたパスの前に記述という事らしい。
念のため再ログインしてwhichでどのパスでプログラムが実行されている確認する。
which samtools
~/local/samtools-0.1.19/samtools
コンパイルが必要なパターンで初めてのmake・・・エラー起きなくてよかった。


Bowtie2に使うゲノムダウンロード

HPはここhttp://bowtie-bio.sourceforge.net/bowtie2/index.shtml
右下にindexesっていうとこがある。
ここからhg19とmm10のダウンロードを行う。
インデックスがついたゲノムファイルを置いておくディレクトリを~/local/bowtie2-2.2.1につくる
cd ../local/bowtie2-2.2.1
mkdir indexes
早速ダウンロード H. sapiens, UCSC hg19をクリックするとPCにダウンロードされるので、ターミナルにリンク先をコピペする。ターミナルへのリンクのコピペはリンクをドラッグアンドドロップすればできる。
ターミナルを起動して、srcに移動してダウンロード。
cd src/
wget "ftp://ftp.ccb.jhu.edu/pub/data/bowtie2_indexes/hg19.zip"
次に解凍。zipなのでunzipを使う。 localに解凍する -dで解凍先に../local/bowtie2-2.2.1を指定
unzip hg19.zip -d ../local/bowtie2-2.2.1/indexes/
mm10も同様に行う。
wget "ftp://ftp.ccb.jhu.edu/pub/data/bowtie2_indexes/mm10.zip"
次に解凍。zipなのでunzipを使う。 localに解凍する -dで解凍先に../local/bowtie2-2.2.1を指定
unzip mm10.zip -d ../local/bowtie2-2.2.1/indexes/
一応indexesディレクトリがどうなったか確認すると以下のようになっている。
ls -1 ../local/bowtie2-2.2.1/indexes/
hg19.1.bt2
hg19.2.bt2
hg19.3.bt2
hg19.4.bt2
hg19.rev.1.bt2
hg19.rev.2.bt2
make_hg19.sh
make_mm10.sh
mm10.1.bt2
mm10.2.bt2
mm10.3.bt2
mm10.4.bt2
mm10.rev.1.bt2
mm10.rev.2.bt2
一つのゲノムファイルに付き6つの.bt2ファイルがある。.shファイルは何に使うか不明。ひとまずこれでOK。fasta形式のゲノムファイルのみが手元に有る場合は、bowtie2のbowtie-buildというコマンドでindexを作る事ができるらしい。そういう機会がでてきたらusageで使い方を確認して行ってみる

2014年4月4日金曜日

Bowtie 2のインストール

HPはここhttp://bowtie-bio.sourceforge.net/bowtie2/index.shtml
Lastest Releaseにあるリンクをクリック
bowtie2-2.2.1-linux-x86_64.zipをクリックするとPCにダウンロードされるので、ターミナルにリンク先をコピペする。ターミナルへのリンクのコピペはリンクをドラッグアンドドロップすればできる。ただ今回はhttp/~/~/~/bowtie2-2.2.1-linux-x86_64.zip/downloadとなっているとうまくいかないので最後の/downloadを消す必要が有る。
ターミナルを起動して、srcに移動する。
cd src/
wget "http://sourceforge.net/projects/bowtie-bio/files/bowtie2/2.2.1/bowtie2-2.2.1-linux-x86_64.zip"
--2014-04-04 11:36:17--  http://sourceforge.net/projects/bowtie-bio/files/bowtie2/2.2.1/bowtie2-2.2.1-linux-x86_64.zip
sourceforge.net をDNSに問いあわせています... 216.34.181.60
sourceforge.net|216.34.181.60|:80 に接続しています... 接続しました。
HTTP による接続要求を送信しました、応答を待っています... 302 Found
場所: http://sourceforge.net/projects/bowtie-bio/files/bowtie2/2.2.1/bowtie2-2.2.1-linux-x86_64.zip/download [続く]
--2014-04-04 11:36:17--  http://sourceforge.net/projects/bowtie-bio/files/bowtie2/2.2.1/bowtie2-2.2.1-linux-x86_64.zip/download
sourceforge.net|216.34.181.60|:80 に接続しています... 接続しました。
HTTP による接続要求を送信しました、応答を待っています... 302 Found
場所: http://downloads.sourceforge.net/project/bowtie-bio/bowtie2/2.2.1/bowtie2-2.2.1-linux-x86_64.zip?r=&ts=1396579000&use_mirror=jaist [続く]
--2014-04-04 11:36:18--  http://downloads.sourceforge.net/project/bowtie-bio/bowtie2/2.2.1/bowtie2-2.2.1-linux-x86_64.zip?r=&ts=1396579000&use_mirror=jaist
downloads.sourceforge.net をDNSに問いあわせています... 216.34.181.59
downloads.sourceforge.net|216.34.181.59|:80 に接続しています... 接続しました。
HTTP による接続要求を送信しました、応答を待っています... 302 Found
場所: http://jaist.dl.sourceforge.net/project/bowtie-bio/bowtie2/2.2.1/bowtie2-2.2.1-linux-x86_64.zip [続く]
--2014-04-04 11:36:18--  http://jaist.dl.sourceforge.net/project/bowtie-bio/bowtie2/2.2.1/bowtie2-2.2.1-linux-x86_64.zip
jaist.dl.sourceforge.net をDNSに問いあわせています... 150.65.7.130, 2001:df0:2ed:feed::feed
jaist.dl.sourceforge.net|150.65.7.130|:80 に接続しています... 接続しました。
HTTP による接続要求を送信しました、応答を待っています... 200 OK
長さ: 26374304 (25M) [application/octet-stream]
`bowtie2-2.2.1-linux-x86_64.zip' に保存中

100%[======================================>] 26,374,304  2.43M/s 時間 10s

2014-04-04 11:36:29 (2.41 MB/s) - `bowtie2-2.2.1-linux-x86_64.zip' へ保存完了 [26374304/26374304]
次に解凍。zipなのでunzipを使う。 localに解凍する -dで解凍先を指定 ちなみにコマンドやファイル名はtabキーで補完されるので使うと便利
unzip bowtie2-2.2.1-linux-x86_64.zip -d ../local/
解凍すると、bowtie2-2.2.1ディレクトリができる。先ほどのHPのPATHのところにBy adding your new Bowtie 2 directory to your PATH environment variable, you ensure that whenever you run bowtie2, bowtie2-build or bowtie2-inspect from the command line, you will get the version you just installed without having to specify the entire path. This is recommended for most users. To do this, follow your operating system's instructions for adding the directory to your PATH.と記述されているので適当にパスを通す。
export PATH=/home/kosugi/local/bowtie2-2.2.1:$PATH
ほんで、home ディレクトリに戻ってbowtie2 -hでusageがでるか確認
bowtie2 -h
Bowtie 2 version 2.2.1 by Ben Langmead (langmea@cs.jhu.edu, www.cs.jhu.edu/~langmea)
以下Usage
動いた!これで場所を選ばず実行できる。
再ログインしても良いようにパスを.bashrcに記述
vim .bashrc
ファイルが開いたらiを押すと書き込める。以下を一番下に追記する。追記したらescキーを押して:wq
## bowtie2
export PATH=/home/kosugi/local/bowtie2-2.2.1:$PATH  :$PATHを後ろに持ってくると今まで書かれたパスの前に記述という事らしい。
念のため再ログインしてwhichでどのパスでプログラムが実行されている確認する。
which bowtie2
~/local/bowtie2-2.2.1/bowtie2
順調順調^^

FastQCのインストール

HPはここhttp://www.bioinformatics.babraham.ac.uk/projects/fastqc/
Download Nowをクリック
FastQC v0.10.1 (Win/Linux zip file)をクリックしてしまうとそのままPCにダウンロードされるのでターミナルにリンク先をコピペする。ターミナルへのリンクのコピペはリンクをドラッグアンドドロップすればできる
ターミナルを起動して、srcに移動する。
cd src/
wget "http://www.bioinformatics.babraham.ac.uk/projects/fastqc/fastqc_v0.10.1.zip"
--2014-04-03 23:36:07--  http://www.bioinformatics.babraham.ac.uk/projects/fastqc/fastqc_v0.10.1.zip
www.bioinformatics.babraham.ac.uk をDNSに問いあわせています... 149.155.132.143
www.bioinformatics.babraham.ac.uk|149.155.132.143|:80 に接続しています... 接続しました。
HTTP による接続要求を送信しました、応答を待っています... 200 OK
長さ: 596632 (583K) [application/zip]
`fastqc_v0.10.1.zip' に保存中

100%[======================================>] 596,632      103K/s 時間 5.9s

2014-04-03 23:36:14 (99.2 KB/s) - `fastqc_v0.10.1.zip' へ保存完了 [596632/596632]
次に解凍。zipなのでunzipを使う。 localに解凍する -dで解凍先を指定 ちなみにコマンドやファイル名はtabキーで補完されるので使うと便利
unzip unzip fastqc_v0.10.1.zip ../local
解凍すると、FastQCディレクトリができるので移動してInstall.txtを見るとActually installing FastQC is as simple as unzipping the zip file it comes in into a suitable location. That's it. Once unzipped it's ready to go.fastqcとなっているのでそのまま使えるらしい。実際にやってみると実行権限が与えられていないのでfastqcファイルに実行権限を与える。与えたら-hでusageがでるか確認
cd FastQC/
chmod +x fastqc
fastqc -h
versionやoptionの説明がでればOK /home/kosugi/local/FastQCにパスを通す。
export PATH=/home/kosugi/local/FastQC:$PATH
ほんで、home ディレクトリに戻ってfastqc -hでうまくいくか確認
[kosugi@~]$ fastqc -h
FastQC - A high throughput sequence QC analysis tool~
〜以下usage〜
動いた!これで場所を選ばず実行できる。
再ログインしても良いようにパスを.bashrcに記述
vim .bashrc
ファイルが開いたらiを押すと書き込める。以下を一番下に追記する。追記したらescキーを押して:wq
## FastQC
export PATH=/home/kosugi/local/FastQC:$PATH  :$PATHを後ろに持ってくると今まで書かれたパスの前に記述という事らしい。
念のため再ログインしてwhichでどのパスでプログラムが実行されている確認する。
which fastqc
~/local/FastQC/fastqc
とりあえず、これでいいのかな・・・

2014年4月3日木曜日

ChIP-seq解析に必要なダウンロードするもの

ChIP-seqの解析に使う代表的なツールとしてなにが必要かな・・・
論文を参考にすると
シーケンスの配列クオリティチェックツール
FastQC : http://www.bioinformatics.babraham.ac.uk/projects/fastqc/
マッピングツール
Bowtie 2 : http://bowtie-bio.sourceforge.net/bowtie2/index.shtml
SAM/BAMファイルの操作ツール
Samtools : http://samtools.sourceforge.net
ピークコールツール
MACS : http://liulab.dfci.harvard.edu/MACS/index.html
もちろんゲノムも必要。Bowtie 2のサイトの右下にリンクからおとす。
(indexが含まれているようなので色々都合が良いと思われる。)
それぞれ、ダウンロード-解凍-インストール-パスを通す等という順番でやっていく。
ツール関係で超絶参考になるサイトは以下のNGS Surfer's Wiki
NGS Surfer's Wikiのツール情報


syntaxhighlighterを使ってコマンドを表示する事にした。

解析の話とは無関係だけど、見た目をすこしみやすくしたいのでsyntaxhighlighterをBloggerに対応させてみた。コマンドラインはsyntaxhighlighterを使って表示したほうが文章と区別できてよいかな。結構、はまった。
Syntax Highlighter Scripts Generatorで全てにチェックを入れてgenerateしてコピー
それをBloggerのテンプレート→カスタマイズの隣のHTMLの編集での直前にペースト-保存
以下にある一行目を削除
<link href='http://alexgorbatchev.com/pub/sh/current/styles/shCore.css' rel='stylesheet' type='text/css'/>
一番下から二行目のSyntaxHighlighter.all();の上に以下を追加すると右上にでる、はてなマーク消える。
SyntaxHighlighter.defaults['toolbar'] = false;

次にcss
shCore.css
にあるcssをメモ帳にでもコピペして以下の場所を適宜変更。
  margin: 0 !important;
  outline: 0 !important;
  overflow: visible !important;
  padding: 2px !important;//0から2pxに変更
  .syntaxhighlighter.ie {
  font-size: .9em !important;
  padding: 1px 0 20px 0 !important; // 1px->20pxに変更
}

これをBloggerのテンプレート>カスタマイズ>上級者向け>cssにコピペ-適用する
使い方はHTML編集でコマンドのところをpreでくくる感じ
<pre class="brush:言語名(コマンドならbash HTMLならhtml)">
//コマンド
</pre>

コマンドだとこんな感じになる
mkdir hoge hogeディレクトリを作る
cd ../     一つ上のディレクトリへ移動
less hoge  hogeの中身を見る
ls -la     ディレクトリの中身確認
wget "URL" ダウンロード
wgetの色変わらんのね。なんでやろ。
NGSの解析とは全く関係ないところにやや時間を費やしてしまった・・・。

2014年4月1日火曜日

ディレクトリの作成(2014/04/04追記パスの通し方)

無事にサーバーに接続できたので、マッピングツールやらゲノムやらダウンロードしたい!
けどその前に、ssh接続してすぐのディレクトリにはなにがあるやろ?
確認。という事でls
[kosugi@~]$ ls

何も無いw(実は隠しファイルはある。-a で.付きのものが見える。)

ダウンロード用、作業フォルダのようなものを作っておこう。
mkdir src     ダウンロード
mkdir work     作業
とした。
[kosugi@~]$ ls 
src work
ディレクトリできた。

srcにツールやゲノムをとりあえずダウンロードしていこう。
<2014/04/03追記>
一般的にはlocalとするようなのでworkはlocal以下に作る。
binもlocal以下に作っておく。
[kosugi@~]$ ls
src  work
[kosugi@~]$ rmdir work
[kosugi@~]$ mkdir local
[kosugi@~]$ ls
local  src
[kosugi@]$ cd local
[kosugi@local]$ mkdir work
[kosugi@local]$ mkdir bin
[kosugi@local]$ ls
bin work
とりあえずこれでOK。
binにはパスも通しておく。パスについてはよくわからないけれどネットに色んな情報が書かれているのでそれを参考にする。パスを通した場所にプログラムを置いていないと自由にコマンドが使えない事だけは理解できた。さらに、サーバーで既通っているパスの上流にパスを通したい。普通はサーバーの/usr/local/binにパスが通っていて、そこにインストールされるかファイルをコピーする。今回は権限が与えられていないので自分のhomeディレクトリに作ったbinにパスを通しておく。
export PATH=/home/kosugi/local/bin:$PATH
とすれば通るはず。
echo $PATH
/home/kosugi/local/bin:/他の色んなパス/usr/local/bin:/bin:/usr/bin:/usr/local/sbin:/usr/sbin:/sbin
と言う事で通った。 ただこのままだと、ログオフしたら消える。
echo $PATH
/他の色んなパス/usr/local/bin:/bin:/usr/bin:/usr/local/sbin:/usr/sbin:/sbin
なので毎回打ち直す必要が有るが、そんなめんどくさい事していられないのでひとまず.bashrcに書き込むといいらしい。するとこのログインした時の初期化問題が解決する。書き込むにはvimというコマンドを使うようだ。(emacsというのもあるらしい。)書き込むにはi 保存終了はescキーを押して:wqと打ち込む
vim .bashrc
ファイルが開いたら書き込める。以下を適当な位置に追加する。
## 自分のlocal/binにパスを通す。
export PATH=/home/kosugi/local/bin:$PATH  :$PATHを後ろに持ってくると今まで書かれたパスの前に記述という事らしい。
と書いて終了。ログインし直しても
echo $PATH
/home/kosugi/local/bin:/他の色んなパス/usr/local/bin:/bin:/usr/bin:/usr/local/sbin:/usr/sbin:/sbin
先頭に/home/kosugi/local/binのパスが通っている。よかったよかった。

サーバーに接続するまでの流れ


サーバーで解析するにあたって
・サーバー構築(サーバー構築のLinux本を一通り読む必要あり
  1. 初期設定
  2. ユーザアカウント作成
・自分自身のPCからサーバーに接続する必要が有る。
  1. ssh接続のための秘密鍵を作成する必要がある。

そこらへんの設定が終われば、コマンド一行でサーバにssh接続できるようになるとのこと。
今回は既にサーバーが準備されている状態であり、秘密鍵の作成までは、管理者にやってもらった。ssh-keygen -t dsaっていうのを使ってやるみたい。一般的らしいのでググると結構でてくる。

PCはmac book pro を使っていて、ターミナルを起動すると
kosugi:~ TK$
となる。そこでls -l(自分がいるディレクトリのファイルの情報をみるコマンド)を叩く。
kosugi:~ TK$ ls -l .ssh
total 16
-rw-------   1 TK  staff   672  3 31 15:22 id_dsa
-rw-r--r--   1 TK  staff   801  3 31 22:38 known_hosts
とid_dsaがちゃんとある事がわかる。ここにssh接続に必要な秘密鍵がある。

kosugi:~ TK$ ssh -i ~/.ssh/id_dsa kosugi@サーバー名とコマンドを叩くと
kosugi:~ TK$ ssh -i ~/.ssh/id_dsa kosugi@サーバー名
Last login: Tue Apr  1 02:11:06 2014 from 〜
[kosugi@サーバー名~]$
というわけでサーバーにssh接続されるようになったヽ(゚∀゚)メ(゚∀゚)メ(゚∀゚)ノ

ただ、Linuxをあつかった事がない人がまず行き詰まるのは、おそらくこのステップ。
よーわからん。って言うのが私の感想。やってもらわなかったらなかなか難しいと思う。時間がある時に、virtual boxなどの仮想環境で秘密鍵の作成はやってみようと思う。けどMac単体で解析する場合や、仮想環境を立ち上げている場合なら基本的にssh接続する必要がないと思われるので、ここの詳細は今はとりあえず省略。

Linuxのお勉強


ほとんどコマンドも知らないので、覚える必要あり。
以下の資料を教科書にしてみる。
LPIC:Linux標準教科書(Ver2.0.0)
http://www.lpi.or.jp/linuxtext/text.shtml
PDF版はアンケートに答える必要あり。
Android版もあるみたい。
https://play.google.com/store/apps/details?id=com.lpijapan.linux

コマンドをひたすら頭に叩き込む。
めっちゃ使うやつは少したってからまとめる。

解析するプラットフォーム


やり方はいろいろあるようです。



  1. GUIソフトで解析(例えば、CLCbio Genome Work bench )
    • 比較的簡単。見た目にわかりやすい。パソコン普通に扱えればできそう。
      • お金かかる。1ライセンス、60万程度。
      • 導入されているプログラムしか使えない。
  2. Macで解析 (Unix)
    • OS Xを生かして、仮想環境無しで始められる。参考資料転がってそう。
    • 緒方さんが作成した、お家でできるMac Bookでやる次世代シーケンスデータ解析がわかりやすい(http://www.ipad-zine.com/b/1520
      • PC性能に解析スピードが依存するので、比較的解析時間かかる。
  3. Linuxで解析
    • 仮想マシーン(Linuxに慣れるにはちょうどいい。virtual box・VM等)
      • PC性能に解析スピードが依存するので、比較的解析時間かかる。
    • サーバー → いわゆる本番環境(PC性能高いので解析スピードが早い。
      • サーバーが必要。→環境がないならお金かかる。
      • セキュリティやサーバーの構築に知識が必要。
  4. Galaxyで解析
    • webベースのツールで必要なツールが大体そろっている。
      • webでやるとuploadにやたら時間かかるときがある。
      • 使いたいツールが使えない。
    • サーバーにシステムごと導入できる。
      • サーバーに一から導入するのは結構大変そう。
2, 3 (4も?)については以下の知識が必要と思われる。
  • Unix/Linux
  • ネットワーク(インターネット)
  • ハードウェア
これらは、必要に応じて、少しずつ勉強していく。まずはLinuxのコマンドをざっと覚える。

幸いにもラボにサーバーが存在する事、基本的な構築が済んでおり管理者にUserアカウントを作成してもらえたので、すぐサーバーを使える!

サーバーを使うデメリットがだいぶ減ったので3-2の本番環境で始める事にした。

0から始めるっていう割には環境整ってるやん!
環境というよりも、解析した経験がゼロってことですね。