陆地棉基因组的纤维发育基因共表达网络构建与分析(3)
首先，我们在这里做一些数据描述：这包括报道的FRGs的表达模式，其基因组区域和与lncRNA的进化关系，接下来我们描述了我们注释的lncRNA的一些表达特征。我们发现FGRs在基因组上呈随机分布，并不全部偏好在纤文中表达。但是FRGs的周围是否存在lncRNA是亚基因组特异的，并且在进化过程中lncRNA都出现在早期的基因分支中，在演化的后期lncRNA逐渐消失在基因附近(上下游10kb)。尽管FGRs在基因组上随机分布，我们发现它们在At和Dt的第12染色体上有偏向性的分布。
其次，我们做一些数据解释：根据构建的基因共表达网络，可以获得基因模块内高度连接的转录本，模块内的基因都是高度正相关或者高度负相关的。通过整合已报道的基因功能、四个阶段中模块表达量随时间变化的表达趋势、GO富集分析和表达模式等数据，发现了可能与纤文发育相关的11个模块。根据这一分析结果，我们可以用“模块（一簇基因）驱动的纤文发育”的假说来解释纤文发育”。
第三，我们在这里进行了一些数据探索：通过预测miRNA靶点来预测lncRNA和mRNA之间的可能的竞争关系；lncRNA作为 miRNA前体的可能性检测；为共表达网络中的三个重要的基因簇构建了三个贝叶斯网络，以推断一个模块内基因之间的因果效应。
我们希望能够在转录组的水平上，系统性的观察纤文发育这个复杂的生物形状，并预测可能的组成成分，相互作用机制和调节模式，以期节省生物实验成本。我们希望我们探究的这些基本的知识可以为最终实现植物改良提供一定的借鉴。
1. 材料与方法
1.1. 数据集和分析工具
报道的基因，mRNA-seq文件和显著的纤文发育相关的GWAS基因座来自本实验室。小RNA测序的清洁数据(clean reads)[7] 也来自本实验室。另外一部分棉花纤文和胚珠中小RNA的测序数据来自National Center for Biotechnology Information (NCBI) ，登录号为GSE61774。
mRNA序列样品为：
种子发芽-0h /种子萌发-5h /种子发芽-10h /子叶种子萌发-24h /子叶种子萌发-48h /子叶种子萌发-72h /子叶种子萌发-96h /子叶种子萌发-120h /种子萌发根 - 24h /种子萌发根 - 48h /种子萌发根 - 72h /种子发芽根 - 96h /种子萌发根/ 120h /根/茎/叶/花/环/花瓣/雄蕊/雌蕊/ calycle / -3dpa胚珠/ -1dpa胚珠/ 0dpa胚珠/ 1dpa胚珠/ 3dpa胚珠/ 5dpa胚珠/ 10dpa胚珠/ 20dpa胚珠/ 25dpa胚珠/ 35dpa胚珠/ 5dpa纤文/ 10dpa纤文/ 20dpa纤文/ 25dpa纤文和几次重复。
我们运用R (R Development Core Team (2016))，R studio( RStudio Team (2016)) 和Perl 来进行数据的处理与统计。
1.2. 得到纤文相关基因(FRGs)
在获得报道的纤文相关基因后，我们发现83个基因可以是棉花纤文相关基因，通过搜索原始发表的文章，我们获得了NCBI基因和其他特定基因的登录号和其它特征信息。然后我们通过Perl脚本根据登录号下载了具体的基因序列。比我们采用本地比对工具blast +的blastn子程序，将83个基因序列作为查询条件 (query)，实验室cds库文件作为参考序列 (reference)，进行局部序列比对。通过分析结果，通过身份、比对长度、不匹配、评估和比特分数存在于TM-1中来确定这些基因为纤文相关基因。对于一些难以区分的序列，我们使用emboss软件中的swith-water子程序根据比对质量分数来确定。因为swith-water是一种经典的全局比对软件，它可以为mRNA-mRNA等长度相差无几的序列比对提供重要的参考信息。最后，我们得到了52个可以对应于TM-1的棉纤文相关基因。中间结果和最终结果文件可以在补充资料中找到[表S1]。
1.3. 得到差异表达基因(DEGs)