转录组测序结题报告1．mRNA纯化：抽提得到的总RNA首先利用10U的DNaseI（Ambion，美国）在37℃消化1小时；然后利用Micropoly(A)PuristTM mRNA purification kit（Ambion，美国），进行mRNA纯化：把RNA稀释到250μl的体积，按照Kit的操作步骤（Cat.No:1919）进行；最后得到的mRNA用100μl预热的THE缓冲液洗脱，利用NanoDrop 进行定量。
2．cDNA合成：cDNA合成是在Ng等2005年发表的方法基础上改进而成（文献1，图1）。
第一链cDNA合成利用GsuI-oligo dT作为反转录引物，10μg的mRNA作为模板，用1000 单位的Superscript II reverse transcriptase (Invitrogen，美国)在42℃作用1小时完成；随后利用NaIO4（Sigma，美国）氧化mRNA的5’帽子结构，并连接生物素；通过Dynal M280磁珠（Invitrogen，美国）筛选连接了生物素的mRNA/cDNA，并通过碱裂解释放第一链cDNA；然后通过DNA ligase（TaKaRa，日本）在第一链cDNA的5’末端加上接头，然后通过Ex Taq polymerase (TaKaRa，日本)合成第二链cDNA。
最后通过GsuI酶切去除polyA和5’端接头。
图1. 全长cDNA合成示意图3．cDNA测序：合成的cDNA利用超声仪（Fisher）打断到300-500bp的范围，利用Ampure beads（Agencourt，美国）进行纯化。
随后纯化的cDNA利用TruSeq TM DNA XXmple Prep Kit – Set A (illumina，美国)制备文库，并利用TruSeq PE Cluster Kit (illumina，美国)进行扩增。
最后在illumina机器上进行测序反应。
测序得到的数据统计见表1.表1. Solexa测序统计样品对照 1 2Reads数目(对) 5,500,000 10,254,848 11,160,428Clean data 5,442,815（98.96%）10,160,130（99.08%）10,998,951（98.55%）平均长度100 100 1005．EST拼装：利用trinity进行拼装。
共得到45,308个EST cluster(contigs)。
具体拼装结果见表2和图2。
表2. 拼装统计样品XXContig数目45,308Contig平均长度698Contig长度范围201-16,169图2. Contigs长度分布(横坐标为基因长度分布，纵坐标为基因数量分布)6．数据分析：6.1 基因预测：采用EMBOSS工具包(参考文献2)中的’GetORF’对拼装得到的contigs进行基因预测，从不同contigs中找到蛋白编码序列。
6.2 基因注释：将预测得到的蛋白编码序列与GenBank的NR、GO、KEGG、KOG等数据库利用blastp进行比对，条件为E value<1e-5，选择匹配最好的一项作为注释信息。
详细结果见annotation.xls，由左至右分别为拼接软件产生的contig名称、基因功能注释、ORF起始与终止位点坐标、正反义链、氨基酸长度、KOG分类。
6.3 GO分析：GO分析利用GoPipe（参考文献3）进行，预测蛋白首先与Swiss-Prot 和TrEMBL数据库进行比对，条件为blastp，E value<1e-5，然后比对结果利用GoPipe程序，根据gene2go，得到预测蛋白的GO信息。
共有4,823个预测蛋白，匹配28,168项GO terms,如图3所示。
详细结果见annotation.xls中“GO”sheet栏。
图3. GO分布6.4 代谢通路构建：利用KEGG数据库（参考文献3），将预测蛋白与KEGG数据库进行比对，条件为双向blast，E value <1e-3；得到预测蛋白的KO number，再根据KO number，获得预测蛋白参与的代谢通路信息。
结果共有2,706个蛋白获得了KO number，它们参与的代谢通路如如图4所示。
详细结果见annotation.xls中“KEGG pathway ”sheet栏。
图4. 编码蛋白所参与的代谢通路类别6.5 表达丰度分析：首先去除低值序列得到clean reads（图5），然后mapping 到拼接的contig上（图6，图7显示mapping的结果），统计每个conig中分别来自2个样品的reads数目，接着转换成RPKM（参考文献4），最后利用DEGseq 程序包中的MARS (MA-plot-based method with Random XXmpling model)模型（参考文献5），计算每个contig代表的基因在2个样品中的表达丰度差异，FDR 值小于0.001的即被认定为具有显著性差异。
详细结果见annotation.xls中“DGE”sheet栏或”express.xlsx”。
图5. 序列质量分析（clean reads为不含N且质量大于5的碱基数至少占全长的一半）图6. 测序饱和度分析(横坐标为reads number，纵坐标为gene number)图7.基因覆盖率统计样本间差异统计详见annotation.xls中“DGE”sheet栏：由左至右分别是基因名称、基因长度、样品A 统计reads数、样品A RPKM值、样品B统计reads 数、样品B RPKM 值、样品A相对样品B表达差异倍数（取Log值）、q-value、显著性判断。
表3. 样品间显著性差差异基因统计样品上调基因数（p<0.001）下调基因数(p<0.001)1/对照2,961 1,0052/对照2,257 362/1 3,352 2,541图8.上下调基因变化(横坐标为gene，纵坐标为统计值)6.6 富集分析：对于每一个代谢通路和GO类别，我们利用超几何分布统计，计算具有显著性表达差异的基因相对全部基因的显著富集情况。
结果在2个代谢通路和7个GO terms 中差异基因具有明显的富集（FDR<0.01） 详细结果见chayi-GO.xlsx 或者chayi-KEGG.xlsx表4. GO term 富集分析结果代谢通路P value 1/对照Carbohydrate Metabolism 0.003961 2/对照Translation1.02E-10 Cell Communication9.12E-06 Signaling Molecules and Interaction 0.002699 Cardiovascular Diseases 0.002838 Immune System 0.002838 2/1Translation0.000185 Energy Metabolism0.000377表5.代谢通路富集分析结果GO TermP value1/对照cell0.028283 metabolism0.028283 2/对照structural molecule activity 2.88E-34 biosynthesis 2.66E-11 cell0.000863 motor activity0.031182 2/1structural molecule activity 1.87E-33 biosynthesis 1.30E-28 cell0.000868 electron transport 0.01974 metabolism0.0313336.7 客户定制分析：6.7.1 调控途径构建于分析Carotenoid biosynthesis 代谢途径相关基因的富集整理。
表6 Carotenoid biosynthesis 代谢途径整理6.7.2 SSR分子标记筛选详细结果见SSR.xlsx6.7.3 SNP鉴定与筛选详细结果见SNP.xlsx7. FTP文件说明所有分析结果都在FTP的对应文件夹中，具体的解释详见“RNA-Seq相关说明”。
8．参考文献：1. Rice, P., I. Longden, and A. Bleasby, EMBOSS: the European Molecular Biology OpenSoftware Suite. Trends Genet, 2000. 16(6): p. 276-7.2. Chen, Z.-Z.X., C.-H. Zhu, S., GoPipe: streamlined gene ontology annotation for batchanonymous sequences with statistics. PROGRESS IN BIOCHEMISTRY AND BIOPHYSICS, 2005. 32(2): p. 187-190.3. KanehiXX, M., et al., KEGG for representation and analysis of molecular networks involvingdiseases and drugs. Nucleic Acids Res. 38(Database issue): p. D355-60.4. Mortazavi, A., et al., Mapping and quantifying mammalian transcriptomes by RNA-Seq. NatMethods, 2008. 5(7): p. 621-8.5. Wang, L., et al., DEGseq: an R package for identifying differentially expressed genes fromRNA-seq data. Bioinformatics. 26(1): p. 136-8.。