第32回日本植物細胞分子生物学会(盛岡)大会・シンポジウム

シンポジウム3「バイオインフォマティクス講習会II(2014)」

8月22日(金) 14:00~17:00 E会場

 

 

講習会のパワーポイントスライド

 

Windowsを使用される方は、TeraTermをインストールしてください。

 

1)以下ダウンロードサイトの「今すぐダウンロード」をクリックし、ダウンロードします

http://www.vector.co.jp/soft/dl/win95/net/se320973.html

2)ダウンロードされたファイルをダブルクリックして、インストールします

 

1.sshを用いた遺伝研スパコンへのログイン方法

 

Windows (TeraTerm)

1)TeraTermを起動します

2)「ホスト」に「gw.ddbj.nig.ac.jp」または「gw2.ddbj.nig.ac.jp」と入力し、「OK」をクリックします

3)初めて接続する場合、以下のような「セキュリティ警告」が出ますが、気にせず「続行」をクリックします

4)「ユーザ名」と「パスフレーズ」にご自身のユーザ名とパスワードを入力し、「OK」をクリックします

5)ログインに成功すると、以下のような画面(コマンドライン)が表示されます。

 

Mac

1)「アプリケーション」>「ユーティリティ」>「ターミナル」を実行します

2)ターミナルウィンドウをできるだけ大きくしてください

3)以下のように"ssh"コマンドを実行します
      ("ssh"、" "(スペース)、"ユーザ名"の後ろに"@gw2.ddbj.nig.ac.jp")

$ ssh mkobayashi@gw2.ddbj.nig.ac.jp  # ユーザー名@gw2.ddbj.nig.ac.jp

4)以下のように問われたら、"yes"を入力し"Enter"キーを押します

Are you sure you want to continue connecting (yes/no)? yes

5)パスワードを聞かれるので、パスワードを入力し"Enter"キーを押します
(入力したパスワードは画面には表示されません)

mkobayashi@gw2.ddbj.nig.ac.jp's password:

6)以下のように表示されれば、接続成功です

---------------------------------------------------------------------
Thank you for using supercomputer system.
This node is in use for login service only. Please use 'qlogin'.
---------------------------------------------------------------------

 

2.qlogin

本ページでは、コマンドをコピー&ペーストすることで、マッピングとアセンブルの実習を行えます
TeraTermにコマンドをコピー&ペーストするには、本ページのコマンドを「ctrl + c」などでコピーし、TeraTermのコマンドライン上でマウスの右クリックをしてペーストします
Macのターミナルへのコピー&ペーストは、「cmd + c」と「cmd + v」で行うことができます

1)"qlogin"とタイプし、"Enter"キーを押します

[mkobayashi@gw2 ~]$ qlogin

2)パスワードを聞かれるので、パスワードを入力し"Enter"キーを押します

mkobayashi@nt097i's password:

 

3.データ

 

・データファイル

morioka.tar.gz

 

1)データのダウンロード

"wget"コマンドでデータをダウンロードします

[mkobayashi@nt096 ~]$ wget http://bioinf.mind.meiji.ac.jp/morioka/morioka.tar.gz

 

2)データの解凍

"tar"コマンドでデータを解凍します

[mkobayashi@nt096 ~]$ tar xzvf morioka.tar.gz

 

3)データの確認

"ls"コマンドでデータを確認します

[mkobayashi@nt096 ~]$ ls -lh morioka
合計 68M
-rw-r--r-- 1 mkobayashi ky-muagrl 993K  8月 10 22:00 2014 IRGSP-1.0_genome_chr01_1-1000020.fasta  # Nipponbareゲノム(染色体1番の1bp~1000020bp)
-rw-r--r-- 1 mkobayashi ky-muagrl  227  8月 17 21:36 2014 assemble.sh                             # アセンブル実行シェルスクリプト
drwxr-xr-x 2 mkobayashi ky-muagrl 4.0K  8月 17 21:38 2014 assemble_results                        # アセンブル結果(実行できなかった時のため)
-rw-r--r-- 1 mkobayashi ky-muagrl  475  8月 17 18:40 2014 map.sh                                  # マッピング実行シェルスクリプト
drwxr-xr-x 2 mkobayashi ky-muagrl 4.0K  8月 17 20:35 2014 map_results                             # マッピング結果(実行できなかった時のため)
-rw-r--r-- 1 mkobayashi ky-muagrl  34M  8月 10 22:00 2014 read1.fastq                             # ペアエンドリード1(DRR008446の一部)
-rw-r--r-- 1 mkobayashi ky-muagrl  34M  8月 10 22:00 2014 read2.fastq                             # ペアエンドリード2(DRR008446の一部)

 

4.マッピング

1)ディレクトリの移動

"cd"コマンドで"morioka"ディレクトリへ移動します

[mkobayashi@nt096 ~]$ cd morioka

 

2)シェルスクリプトの確認

"cat"コマンドでシェルスクリプトの内容を確認します
(bwaでマッピングし、samtoolsでマッピング結果を整理)

[mkobayashi@nt096 morioka]$ cat map.sh
#!/bin/sh
#$ -S /bin/sh
#$ -cwd
bwa index \                              # bwaのindexコマンド:リファレンス配列のインデックスを作成する
 IRGSP-1.0_genome_chr01_1-1000020.fasta  # インデックスを作成するリファレンスゲノム配列(fasta形式)

bwa aln \                                  # bwaのalnコマンド:リードをリファレンス配列にマッピングする
 IRGSP-1.0_genome_chr01_1-1000020.fasta \  # マッピングするリファレンスゲノム
 read1.fastq \                             # マッピングするリード配列(fastq形式)
 > read1.sai                               # マッピング結果を"read1.sai"へ出力

bwa aln \                                  # bwaのalnコマンド:リードをリファレンス配列にマッピングする
 IRGSP-1.0_genome_chr01_1-1000020.fasta \  # マッピングするリファレンスゲノム
 read2.fastq \                             # マッピングするリード配列(fastq形式)
 > read2.sai                               # マッピング結果を"read2.sai"へ出力
 
bwa sampe \                                # bwaのsampeコマンド:ペアエンドリードのマッピング結果をsam形式で出力する
 IRGSP-1.0_genome_chr01_1-1000020.fasta \  # マッピングしたリファレンスゲノム
 read1.sai \                               # リード1のマッピング結果ファイル
 read2.sai \                               # リード2のマッピング結果ファイル
 read1.fastq \                             # リード1のファイル
 read2.fastq \                             # リード2のファイル
 > map.sam                                 # マッピング結果ファイル(sam形式)

samtools faidx \                         # samtoolsのfaidxコマンド:リファレンス配列のインデックスを作成する
 IRGSP-1.0_genome_chr01_1-1000020.fasta  # マッピングしたリファレンスゲノム

samtools view \  # samtoolsのviewコマンド:samファイルやbamファイルを入出力する(形式を変換できる)
 -bS \           # sam形式ファイルを入力し、bam形式で出力することを指定
 map.sam \       # 入力するsam形式ファイル
 > map.bam       # 出力するbam形式ファイル

samtools sort \  # samtoolsのsortコマンド:bam形式ファイルをソートする
 map.bam \       # 入力するbam形式ファイル
 map_sort        # 出力するbam形式ファイル(ここで指定した名前に拡張子".bam"が付加されたファイルが出力される)

samtools index \ # samtoolsのindexコマンド:bam形式ファイルのインデックスを作成する
 map_sort.bam    # インデックスを作成するbam形式ファイル

bwa, samtoolsのコマンドに関する詳しい情報は、各ツールのHPからマニュアルを参照ください

bwaのHP

samtoolsのHP(version 0.1.19まで)

samtoolsのHP(version 1以降)

 

3)シェルをキューに登録

"qsub"コマンドでシェルをキューに登録します

[mkobayashi@nt096 morioka]$ qsub map.sh
Your job 498722 ("map.sh") has been submitted

 

4)キューの確認

"qstat"コマンドでキューの状況を確認します

[mkobayashi@nt096 morioka]$ qstat
job-ID  prior   name       user         state submit/start at     queue                          jclass                         slots ja-task-ID
------------------------------------------------------------------------------------------------------------------------------------------------
 498719 0.00000 QLOGIN     mkobayashi   r     08/17/2014 18:40:53 login.q@nt096i                                                    1
 498722 0.00000 map.sh     mkobayashi   qw    08/17/2014 18:41:37                                                                   1

"qstat"コマンドで"map.sh"の表示が無くなったらマッピング終了です

 

5)マッピング結果の確認

"ls"コマンドで結果ファイルを確認します

[mkobayashi@nt096 morioka]$ ls -lh
...
-rw-r--r-- 1 mkobayashi ky-muagrl 489K  8月 17 22:05 2014 IRGSP-1.0_genome_chr01_1-1000020.fasta.sa   # ゲノムのデータベースファイル
-rw-r--r-- 1 mkobayashi ky-muagrl 245K  8月 17 22:05 2014 IRGSP-1.0_genome_chr01_1-1000020.fasta.pac  # ゲノムのデータベースファイル
-rw-r--r-- 1 mkobayashi ky-muagrl 977K  8月 17 22:05 2014 IRGSP-1.0_genome_chr01_1-1000020.fasta.bwt  # ゲノムのデータベースファイル
-rw-r--r-- 1 mkobayashi ky-muagrl   40  8月 17 22:05 2014 IRGSP-1.0_genome_chr01_1-1000020.fasta.ann  # ゲノムのデータベースファイル
-rw-r--r-- 1 mkobayashi ky-muagrl   21  8月 17 22:05 2014 IRGSP-1.0_genome_chr01_1-1000020.fasta.amb  # ゲノムのデータベースファイル
-rw-r--r-- 1 mkobayashi ky-muagrl 5.0M  8月 17 22:05 2014 read1.sai                                   # リード1のインデックスファイル
-rw-r--r-- 1 mkobayashi ky-muagrl 4.8M  8月 17 22:05 2014 read2.sai                                   # リード1のインデックスファイル
-rw-r--r-- 1 mkobayashi ky-muagrl  99M  8月 17 22:06 2014 map.sam                                     # マッピング結果ファイル(SAMフォーマット)
-rw-r--r-- 1 mkobayashi ky-muagrl   22  8月 17 22:06 2014 IRGSP-1.0_genome_chr01_1-1000020.fasta.fai  # ゲノムのインデックスファイル
-rw-r--r-- 1 mkobayashi ky-muagrl  22M  8月 17 22:06 2014 map.bam                                     # マッピング結果ファイル(BAMフォーマット)
-rw-r--r-- 1 mkobayashi ky-muagrl 2.8K  8月 17 22:06 2014 map_sort.bam.bai                            # ソートしたマッピング結果のインデックスファイル
-rw-r--r-- 1 mkobayashi ky-muagrl  22M  8月 17 22:06 2014 map_sort.bam                                # ソートしたマッピング結果ファイル(BAMフォーマット)
...

 

"less"コマンドでsamファイルを確認します

[mkobayashi@nt093 morioka]$ less -S map.sam
@SQ     SN:chr01        LN:1000020
@PG     ID:bwa  PN:bwa  VN:0.6.1-r104
DRR008446.37648089      81      chr01   996     37      102M    =       999     -99     TCTCCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAACCCTAAACCCTAACCCTAAACCCTAAACCCTAAACCC        9>A<:CCBBDB@>DB?ABAAC?;;=@;(?;>FEFDB>AHA=;;=7@HHGFB8AIIGEFBIGGGDBGIHGIHGFEGIGFIHCHCGFGGJIHHFFBDDDDA@@@        XT:A:U  NM:i:4  XN:i:5  SM:i:37 AM:i:25 X0:i:1  X1:i:0  XM:i:4  XO:i:0  XG:i:0  MD:Z:0G0G1A0A97
DRR008446.37648089      161     chr01   999     25      102M    =       996     99      CCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAACCCTAAACCCTAACCCTAAACCCTAAACCATAAACCCAGA        @@@DDADDFFHHDFBFEGGAHAEHJIEEGGHGHGIDIGGIIB@B;FD89=BF2C7=@7@G@3?2=A??DE3=9A@AAB?==AC9<9<9:>C###########        XT:A:U  NM:i:5  XN:i:2  SM:i:25 AM:i:25 X0:i:1  X1:i:0  XM:i:5  XO:i:0  XG:i:0  MD:Z:0A0A89C7T0A1
DRR008446.28363971      97      chr01   1004    37      102M    =       1002    100     AACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAACCCTAAACCCTAACCCTAAACCCTAAACCCAAAACCCTAAACCAG        ?@@DDDFFHHHHHHGGHHJGEDHIIIJJJJIGHGEHHEHICHGGIGGCGBHGFHJJJDHEGHCA;.;>BBEDDB@?ACCBB<9CC?@15?ADDDA93399<B        XT:A:U  NM:i:3  SM:i:37 AM:i:37 X0:i:1  X1:i:0  XM:i:3  XO:i:0  XG:i:0  MD:Z:87T12C0T0
DRR008446.28363971      145     chr01   1002    37      102M    =       1004    -100    CTAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAACCCTAAACCCTAACCCTAAACCCTAAACCCAAAACCCTAAACC        ???9BAADBA>?C>>BDCDCCCBAADCDCBB@CEFFFFHGHFFHD=GGFDF@BIIHCHFGJJJJHHJIJIGDGJHGIHHGGHHHFEHIIDFFAHFEDFF@@@        XT:A:U  NM:i:3  SM:i:37 AM:i:37 X0:i:1  X1:i:0  XM:i:3  XO:i:0  XG:i:0  MD:Z:0T0A87T12
DRR008446.32762433      99      chr01   1001    29      102M    =       1006    102     CTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAACCCTAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAACC        CCCFFFFFHHHHHJJJJJJJJJJIIJIJJIIEGFHIIGEGHJJIFJHGFEGIGGFHHAE7@ECCCCBBFFCECCEDBC?CC?BB9<>??B?<<(93?<39?A        XT:A:R  NM:i:4  SM:i:0  AM:i:0  X0:i:2  X1:i:0  XM:i:4  XO:i:0  XG:i:0  MD:Z:66A2C0T29A1        XA:Z:chr01,+1008,102M,4;
DRR008446.32762433      147     chr01   1006    37      19M5I78M        =       1001    -102    CCCTAAACCCTAAACCCTAACCCTAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAACCCTAAACCCTAACCCTAAACCCTAAACCCTAAACCCTAAAC        1>CCBBADDCC?<>DC<5?:DCABBCCC@@BDCBEBHHEHGGGAIHHIGGEGGDJIIHHGJJJIIHDJIHGHGGIHDIIHHIHEIIHIIHDHFCFFFFFBB@        XT:A:U  NM:i:5  SM:i:37 AM:i:0  X0:i:1  X1:i:0  XM:i:0  XO:i:1  XG:i:5  MD:Z:97
DRR008446.97140311      99      chr01   1004    60      81M2I19M        =       1006    102     AACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAACCCTAAACCCTAACCCTAAACCCTAAAAACCCTAAACCCTAAACC        CCCFFFFFHHHHHIIJIJJIIIIJJJJIJGIJJJJJJHIJIIHGHIJJFHGGIJFGHIIJIG?EEHFFFFFBCCC?ACCDD<5:ACDDDB??:C<B@<9:>A        XT:A:U  NM:i:2  SM:i:37 AM:i:37 X0:i:1  X1:i:0  XM:i:0  XO:i:1  XG:i:2  MD:Z:100
DRR008446.97140311      147     chr01   1006    60      79M2I21M        =       1004    -102    CCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAACCCTAAACCCTAACCCTAAACCCTAAAAACCCTAAACCCTAAACCCT        ?CACDBA?::9?B?B@ACBAAB@CCDACA@C>CB@FCEGHACIIGEIHHHGBHHFIIIEIHFJIIGHFGJIGGIHIIGHJJIHHIIGEHFCHHHFDDDFCB@        XT:A:U  NM:i:2  SM:i:37 AM:i:37 X0:i:1  X1:i:0  XM:i:0  XO:i:1  XG:i:2  MD:Z:100

1列目:マップされたリード名
2列目:FLAG(マッピングされたか、リードの向きなどの情報)
3列目:染色体名
4列目:リードがマップされたポジション
5列目:マッピングクオリティー
6列目:CIGAR(どのようにアライメントされたかの情報)
7列目:ペアのもう一方がマップされた染色体名("="なら、同じ染色体)
8列目:ペアのもう一方がマップされたポジション
9列目:シーケンスされたゲノムフラグメントの長さ("-"はフラグメントの配列が逆向き)
10列目:リードの塩基配列
11列目:リードの塩基のクオリティー
12列目:その他の情報(リードがほかのゲノムポジションにもマップされたかどうかなど)

samフォーマットの詳細は、bwaのマニュアルの"SAM ALIGNMENT FORMAT"を参照ください

 

"samtools"の"tview"コマンドでマッピング結果を可視化します

[mkobayashi@nt096 morioka]$ samtools tview map_sort.bam IRGSP-1.0_genome_chr01_1-1000020.fasta

"→"キーで1000bp付近へ移動

または、"g"キーをタイプ後、"chr01:1000"と入力して1000bp付近へ移動します

 1001      1011      1021           1031       1041       1051       1061      1071       1081        1091      1101       1111      1121      1131
NCTAAACCCTAAACCCTAAACCCTA*****AACCCTAAACCCT*AAACCCT*AAACCCTAAACCCT*AACCCTAAACCCT*AACCCTAAACCCT**AAACCCTAAACCCTAAACCCTAAA*CCCTAAACAGCTGACAGTACGATAGATCCACGCGAG
C........................     ............. ....... .............. ............. .............  ........................ ..........        ....CG............
c,,,,,,,,,,,,,,,,,,,,,,,,*****,,,,,,,,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,,,,,,,,,,,,*,,,,,,,,,,,,,**,,,,,,,,,,,,,                              ....CG............
C........................*****.............*.......*..............*.............*.............**.....A.......AG.                           ....CG............
C........................*****.............*.......*..............*.............*.............**.................                          ....CG............
C........................*****.............*.......A..............A......*......*.............**................                           ....CG............
 ........................*****.............*.......*..............*........C..TA*.............**................C.                         ....CG............
 ........................*****.............*.......*..............*.............*.............**..................                          ...CG............
 ........................*****......*......*.......*..............*.............*.............**...................                                    ......
 ........................*****.............*.......*..............*..A..CT......*.............**..................                                       ....
  ct,,,,,,,,,,,,,,,,,,,,,*****,,,,,,,,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,,,,,,,,,,,,*,,,,,,,,,,,,,**,,,,,,a,,,,,,,,,,,,
  ,,,,,,,,,,,,,,,,,,,,,,,*****,,,,,,*,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,,,,,,,,,,,,*,,,,,,,,,,,,,**,,,,,,,,,,,,,,,,,,,,
    .....................*****.............*.......*..............*.............*.............**......A............AG
    .....................*****.............*.......*..............*.............*.............AA...................
     ....................*****.............*.......*........C.....*.............*.............**.C....................
      ,,,,,,,,,,,,,,,,,,,accct,,,,,,,,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,,,,,,,,,,,,*,,,,,,,,,,,,,**,,,,,,,,,,,,,,,,,,
      ,,,,,,,,,,,,,,,,,,,*****,,,,,,,,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,,,,,,,,,,,,*,,,,,,,,,,,,,aa,,,,,,,,,,,,,,,,,,,,,
      ...................*****.............*.......*..............*.............*.............**.......................
      t,,g,tct,,,,,,,,,,,*****,,,,,,,,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,,,,,,,,,,,,*,,,,,,,,,,,,,**,,,,,,,,,,,,,,,,,,,,,,,
      ...................*****.............*.......*..............*.............*.............**.......................
       ..................*****.............*.......*..............*.............*.............**....................C..C*
        ,,,,,,,,,,,,,,,,,*****,,,,,,,,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,,,,,,,,,,,,*,,,,,,,,,,,,,**,,,,,,,,,,,,,,,,,,,,,,,,*,
         c,,,,,,,,,,,,,,,*****,,,,,,,,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,,,,,,,,,,,,*,,,,,,,,,,,,,**,,,,,,,,,,,,,,,,,,,,,,,,*,,
           ,,,,,,,,,,,,,,*****,,,,,,,,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,,,,,,,,,,,,*,,,,,,,,,,,,,**,,,,,,,,,,,,,,,,,,,,,,,,*a,,c
           ,,,,,,,,,,,,,,*****,,,,,,,,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,,,,,,,,,,,,*,,,,,,,,,,,,,**,,,,,,,,,,,,,,,,,,,,,,,,*,,,,
           ,,,,,,,,,,,,,,*****,,,,,,*,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,,,,,,,,,,,,*,,,,,,,,,,,,,**,,,,,,,,,,,,,,,,,,,,,,,,*,,,,,
           ..............*****......*......*.......*..............*.............*.............**........................*.....
            ,t,t,,,,,,,,,*****,,,,,,,,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,,,,,*,,,,,,*,,,,,,,,,,,,,**,,,,,,,,,,,,,,,,,,,,,,,,*,,,,,,
             ............*****.............*.......*..............*.............*.............**........................*......
             a,,,,,,,,,c,*****,,,,,,,,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,a,,ct,,,,,,*,,,,,,,,,,,,,**,,,,,,,,,,,,,,,,,,,,,,,,*,,,,,,
             ,,,,,,,,,,,,*****,,,,,,,,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,,,,,,,,,,,,*,,,,,,,,,,,,,**,,,,,,,,,,,,,,,,,,,,,,,,*,,,,,,
              ...........*****.............*.......*..............*.............*.............**........................GAT.GG.
               ,,,,,,,,,,*****,,,,,,,,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,,,,,,,,,,,,*,,,,,,,,,,,,,**,,,,,,,,,,,,,,,,,,,,,,,,*,,,,,,,,
               ,,,,,,,,,,*****,,,,,,,,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,,,,,,,c,,ta*,,,,,,,,,,,,,**,,,,,,,,,,,,,,,,,,,,,,,,*,,,,,,,,
               t,,,,,,,,,*****,,,,,,,,,,,,,*,,,,,,,*,,,,,,,,,,,,,,*,,,,,,,,,,,,,*,,,,,,,,,,,,,**,,,,,,,,,,,,,,,,,,,,,,,,*,,,,,,,,

"tview"の見方

  1261      1271        1281        1291 # 1行目はゲノム上のポジション
CGCCATCGACAGTAGGTC**TGCGCGAGA*GGGG*CACCG # 2行目はリファレンスゲノムの配列
..................  CA....... .... ..... # 3行目はマッピング結果から得られるコンセンサス配列
.......   ........**CA.......G....*..... # 4行目以降がマップされたリード
.................      ......G....*..... # "."は→向きにマップされたリードの、リファレンスと同じ塩基
,,,,,,,,,,,,,,,,,,**ca,,,,  ,*,,,,a,,,,, # "C"のような大文字の塩基は、→向きにマップされたリードの、リファレンスと異なる塩基
..................**CA.......*.    ,,,,, # ","は←向きにマップされたリードの、リファレンスと同じ塩基
..................**CA.......G....*... , # "c"のような小文字の塩基は、←向きにマップされたリードの、リファレンスと異なる塩基
..................**CA.......G....*..... # "*"のポジションには、どれかのリードにinsertあり

以下のようにエラーが表示されたら、"esc"キーをタイプしてください

     +------------------------------------------------+
     | Goto: chr1:1000[bam_parse_region] fail to determine the sequence name.
     +------------------------------------------------+

"tview"を終了するには、"q"キーをタイプしてください

 

5.de novo アセンブル

1)シェルの確認

"cat"コマンドでシェルを確認します

[mkobayashi@nt097 morioka]$ cat assemble.sh
#!/bin/sh
#$ -S /bin/sh
#$ -cwd
/usr/local/pkg/velvet/velvet_1.2.10/velveth \ # velvethコマンド:配列をk-merに分割する
 out_k51 \                                    # 第一引数:結果出力ディレクトリ
 51 \                                         # 第二引数:k-merサイズ;
 -fastq \                                     # リードファイルのフォーマットがfastqであることを指定
 -shortPaired \                               # リードがショートリードのペアエンドであることを指定
 -separate \                                  # ペアエンドリードが二つのファイルからなることを指定
 read1.fastq \                                # リード1ファイル
 read2.fastq                                  # リード2ファイル

/usr/local/pkg/velvet/velvet_1.2.10/velvetg \ # velvetgコマンド:velveth結果配列の類似性からde Bruijnグラフを作成し、コンティグを作成する
 out_k51 \                                    # velveth結果のあるディレクトリ
 -scaffolding no \                            # scaffoldingを行わない事を指定
 -cov_cutoff auto                             # リードのカバレッジが低い場合に配列を信用しない。そのカバレッジの閾値。"auto"で自動計算する

velvetのコマンドに関する詳しい情報は、velvetのHPからマニュアルを参照ください

velvetのHP

 

2)シェルをキューに登録

"qsub"コマンドでシェルをキューに登録します

[mkobayashi@nt097 morioka]$ qsub assemble.sh
Your job 498785 ("assemble.sh") has been submitted

 

3)結果ファイル確認

"ls"コマンドで結果ファイルを確認します

[mkobayashi@nt097 morioka]$ ls -lh out_k51/
合計 83M
-rw-r--r-- 1 mkobayashi ky-muagrl 1.4M  8月 17 21:38 2014 Graph
-rw-r--r-- 1 mkobayashi ky-muagrl 1.3M  8月 17 21:38 2014 LastGraph
-rw-r--r-- 1 mkobayashi ky-muagrl  996  8月 17 21:38 2014 Log        # ログファイル
-rw-r--r-- 1 mkobayashi ky-muagrl 983K  8月 17 21:38 2014 PreGraph
-rw-r--r-- 1 mkobayashi ky-muagrl  43M  8月 17 21:38 2014 Roadmaps
-rw-r--r-- 1 mkobayashi ky-muagrl  40M  8月 17 21:38 2014 Sequences
-rw-r--r-- 1 mkobayashi ky-muagrl 659K  8月 17 21:38 2014 contigs.fa # コンティグファイル
-rw-r--r-- 1 mkobayashi ky-muagrl  22K  8月 17 21:38 2014 stats.txt

 

4)コンティグの統計値確認

"tail"コマンドでLogファイルを確認します

[mkobayashi@nt097 morioka]$ tail -2 out_k51/Log
Median coverage depth = 21.539821
Final graph has 347 nodes and n50 of 4731, max 39635, total 636872, using 0/304360 reads

配列数は347個、n50値は4,731bp、配列の最大長は39,635bp、総塩基長は636,872bpであることがわかります
      (ただし、これら統計値には、短くてコンティグとなっていないような配列も含まれています)

 

6.SSH接続の終了

SSH接続を終了するには、"exit"コマンドを実行しログアウトしてから、もう一度"exit"コマンドを実行してください

[mkobayashi@nt097 morioka]$ exit
logout
Connection to nt097i closed.
/home/geadmin2/UGER/utilbin/lx-amd64/qlogin_wrapper exited with exit code 0
[mkobayashi@gw2 ~]$ exit