項目數據：
- kongyu_131_PCRfree_.CCAAT_L006_R1_001.fastq.gz （100X）（19G）
kongyu_131_PCRfree_.CCAAT_L006_R2_001.fastq.gz （100X）（20G） - Y255_PCRfree_.TCCGC_L005_R1_001.fastq.gz （30X）（5.4G）
Y255_PCRfree_.TCCGC_L005_R2_001.fastq.gz （30X）（6.0G） - all.chrs.con.fasta （364M） 日本晴 ( /public/genome/rice/msu/version_7.0/all.dir/all.chrs.con.fasta)
工具：
- BWA 參照：比對工具之 BWA 使用方法
- SAMtools 參照：SAM格式 及 比對工具之 samtools 使用方法
- IGV 參照：可視化工具之 IGV 使用方法
- SOAPdenovo 參照： 基因組組裝工具之 SOAPdenovo 使用方法
策略：
- 將測序的二代reads使用BWA比對到參考基因組，分成不同的窗口，按窗口進行局部組裝，然后合並。
預備知識：
- 能用熟練使用 Perl 和 shell 寫腳本
- 會熟練使用 PBS 提交任務
- BWA使用方法
- IGV使用方法
- SOAPdenovo使用方法
具體步驟:
1.前期准備
NP=`cat $PBS_NODEFILE | wc -l` NN=`cat $PBS_NODEFILE | sort -u | wc -l` cd $PBS_O_WORKDIR export LD_LIBRARY_PATH=\$LD_LIBRARY_PATH:/public/software/htslib-1.3/lib date samtoolsdir=/public/software/samtools-1.3/bin bwadir=/public/software/bwa-0.7.12-intel dir=/public/scripts/liyan/scripts2016 sample=Y255 out=/public/home/zxli/Project/Project_Sliced_Assembly ref=/public/genome/rice/msu/version_7.0/all.dir/all.chrs.con.fasta fq1=/public/home/zxli/data/rice/Y255_PCRfree_.TCCGC_L005_R1_001.fastq.gz fq2=/public/home/zxli/data/rice/Y255_PCRfree_.TCCGC_L005_R2_001.fastq.gz
2.比對
/public/software/bwa-0.7.12-intel/bwa index $ref $bwadir/bwa mem -t $NP -f -M -R "@RG\tID:$sample\tLB:$sample\tSM:$sample\tPL:illumina\tPU:$sample" $ref $fq1 $fq2 > $out/bwamem_$sample.sam date
3.可視化的前處理
samtools view -@ 10 -bS -F 4 ./contigs_sequence_align_to_public_genome.sam > ./contigs_sequence_align_to_public_genome.bam
samtools sort -@ 10 ./contigs_sequence_align_to_public_genome.bam contigs_sequence_align_to_public_genome.sorted
samtools index contigs_sequence_align_to_public_genome.sorted.bam contigs_sequence_align_to_public_genome.sorted.bai
samtools depth contigs_sequence_align_to_public_genome.sorted.bam >depth_reads.txt
wc -l depth_reads.txt > Coverage-aln_reads.txt
$samtoolsdir/samtools view -@ $NP -Sb $out/bwamem_$sample.sam -o $out/bwamem_$sample.bam $samtoolsdir/samtools sort -@ $NP $out/bwamem_$sample.bam -o $out/bwamem_$sample.sorted.bam $samtoolsdir/samtools index $out/bwamem_$sample.sorted.bam
4.按窗口分類reads
寫perl腳本，提取SAM中的reads名稱，去fastQ里提取原始reads，按窗口分類，為下一步的局部組裝做准備。
5.局部組裝
局部組裝的問題：
已經有兩批人沒組出來了，局部組裝大多不可能組裝出完整的100K窗口，因為二代序列reads太短，重復序列太多，重復序列會導致連接中斷，一個窗口會出現很多片段，而且也沒有方法將其繼續連接起來，所以他們都半途而廢了。
后續可能會遇到的情況，必須借助后期的分析手段，將諸多片段連接成完整的序列。
杜發的文章，完全是在無參考基因組的情況下，denovo組裝，利用多種手段，才將零碎的序列組裝成完整的基因組。
其他:
PCRfree,就是DNA樣品不是通過PCR進行擴增的,防止PCR中的錯誤造成影響.
map,就是確定基因在染色體中的位置.
眾多軟件都可以利用 fastq.gz 文件, less也可以查看此類型的文件.
問題:
- fastq文件里都有些什么? 參見 fastQ格式
簡介:fastQ是二代測序下機數據的格式, 它存儲了核酸序列和相應質量評價的信息,該格式包含四行:
第一行: @HWUSI-EAS100R:6:73:941:1973#0/1 ; 以@開頭，后面是reads的ID以及其他信息，例如上例中 HWUSI-EAS100R代表Illmina設備名稱，6代表flowcell中的第六個lane，73代表第六個lane中的第73個 tile，941:1973代表該read在該tile中的x：y坐標信息；#0，若為多樣本的混合作為輸入樣本，則該標志代表樣本的編號，用來區分個樣本中的reads；/1代表paired end中的前一個read。
第二行: GATTTGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTT ; read序列
第三行: +HWUSI-EAS100R:6:73:941:1973#0/1 ; 以“+”開頭，跟隨者該read的名稱（一般於@后面的內容相同），但有時可以省略，但“+”一定不能省。
第四行: !''*((((***+))%%%++)(%%%%).1***-+*''))**55CCF>>>>>>CCCCCCC6 ; 以“+”開頭，跟隨者該read的名稱（一般於@后面的內容相同），但有時可以省略，但“+”一定不能省。
本項目中的原始reads格式如下: (每條read均為150 bp)
@ST-E00126:176:HL7H5CCXX:1:1101:23368:6319 1:N:0:GTCCGC TATCGCTTGCTCAAACGCTGGGCAACTAGTCTCTAGTCTTGGTCGAGTGTGTGGTGGACTTGGTCGTGGACATGCTTGGTTCTTAGATCGTGTTTTGTATTCAGGGTTGCTTGTACCCTTTCTTTCTTGCACTGTCCATATTCTAATGCA + AAAAFKKFKKKKKKKKKKKKKKKKKKKKKKKFFKKKKKKKKKKKKFAFKK,,FKK,A<F<KAFKKKFKKKFKKKKKKKFKF<,,,AKKKKFK,FFKKKKKKFA<KK7<,<<,,7<AFFFKKKAFFKKKKKKK77FFA,,<FKKKKAAFAF
補充：雙端測序時，fastq文件中，R1 和 R2 兩個文件中的reads是如何排列的？
R1 @ST-E00126:176:HL7H5CCXX:1:1101:7638:1221 1:N:0:NTCCGC @ST-E00126:176:HL7H5CCXX:1:1101:8572:1221 1:N:0:NTCCGC R2 @ST-E00126:176:HL7H5CCXX:1:1101:7638:1221 2:N:0:NTCCGC @ST-E00126:176:HL7H5CCXX:1:1101:8572:1221 2:N:0:NTCCGC
- 參考基因組（）是個什么玩意?
>Chr1
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNN
CTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCCTAAACCC
TAAACCCTAACCCTAAACCCTAACCCTAAACCCTAAACCCTAAACCCTAA
ACCCTAAACCCTAAACAGCTGACAGTACGATAGATCCACGCGAGAGGAAC
CGGAGAGACAACGGGATCCAGGCGCCAGCGACGGATCCGGGCGAGAGGGG
AGTGGAGATCATGGATCCGTGCGGGAGGGGAAGAAGTCGCCGAATCCGAC
......
chr1一共有865419行, 每一行50個鹼基, 最后一行23個鹼基, 一共43270923個鹼基, 約為43Mb.
該染色體的頭部 尾部 以及中間有少量的NNNN組成的gap, 是沒有組裝出來的區域.
- 水稻基因組基本常識? 參見 水稻基因組計划(百科)
分秈、粳稻兩個亞種, 一共12對染色體, 基因組大小:430Mb, 預測有32000～56000個基因,
- BWA 比對的原理, 以及之后生成的 SAM 文件該如何解讀, SAM格式是如何存儲比對信息的?
額外參考：