項(xiàng)目數(shù)據(jù):
工具:
策略:
預(yù)備知識(shí):
具體步驟: 1.前期準(zhǔn)備 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.比對(duì) /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,按窗口分類,為下一步的局部組裝做準(zhǔn)備。
5.局部組裝
局部組裝的問(wèn)題: 已經(jīng)有兩批人沒(méi)組出來(lái)了,局部組裝大多不可能組裝出完整的100K窗口,因?yàn)槎蛄衦eads太短,重復(fù)序列太多,重復(fù)序列會(huì)導(dǎo)致連接中斷,一個(gè)窗口會(huì)出現(xiàn)很多片段,而且也沒(méi)有方法將其繼續(xù)連接起來(lái),所以他們都半途而廢了。 后續(xù)可能會(huì)遇到的情況,必須借助后期的分析手段,將諸多片段連接成完整的序列。 杜發(fā)的文章,完全是在無(wú)參考基因組的情況下,denovo組裝,利用多種手段,才將零碎的序列組裝成完整的基因組。
其他: PCRfree,就是DNA樣品不是通過(guò)PCR進(jìn)行擴(kuò)增的,防止PCR中的錯(cuò)誤造成影響. map,就是確定基因在染色體中的位置. 眾多軟件都可以利用 fastq.gz 文件, less也可以查看此類型的文件.
問(wèn)題:
簡(jiǎn)介:fastQ是二代測(cè)序下機(jī)數(shù)據(jù)的格式, 它存儲(chǔ)了核酸序列和相應(yīng)質(zhì)量評(píng)價(jià)的信息,該格式包含四行: 第一行: @HWUSI-EAS100R:6:73:941:1973#0/1 ; 以@開頭,后面是reads的ID以及其他信息,例如上例中 HWUSI-EAS100R代表Illmina設(shè)備名稱,6代表flowcell中的第六個(gè)lane,73代表第六個(gè)lane中的第73個(gè) tile,941:1973代表該read在該tile中的x:y坐標(biāo)信息;#0,若為多樣本的混合作為輸入樣本,則該標(biāo)志代表樣本的編號(hào),用來(lái)區(qū)分個(gè)樣本中的reads;/1代表paired end中的前一個(gè)read。 第二行: GATTTGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTT ; read序列 第三行: +HWUSI-EAS100R:6:73:941:1973#0/1 ; 以“+”開頭,跟隨者該read的名稱(一般于@后面的內(nèi)容相同),但有時(shí)可以省略,但“+”一定不能省。 第四行: !''*((((***+))%%%++)(%%%%).1***-+*''))**55CCF>>>>>>CCCCCCC6 ; 以“+”開頭,跟隨者該read的名稱(一般于@后面的內(nèi)容相同),但有時(shí)可以省略,但“+”一定不能省。 本項(xiàng)目中的原始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 補(bǔ)充:雙端測(cè)序時(shí),fastq文件中,R1 和 R2 兩個(gè)文件中的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個(gè)堿基, 最后一行23個(gè)堿基, 一共43270923個(gè)堿基, 約為43Mb. 該染色體的頭部 尾部 以及中間有少量的NNNN組成的gap, 是沒(méi)有組裝出來(lái)的區(qū)域.
分秈、粳稻兩個(gè)亞種, 一共12對(duì)染色體, 基因組大小:430Mb, 預(yù)測(cè)有32000~56000個(gè)基因,
額外參考: |
|
來(lái)自: 生物_醫(yī)藥_科研 > 《基因組組裝》