• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    利用WGCNA進(jìn)行玉米花期基因共表達(dá)模塊鑒定

    2019-01-17 01:25:02楊宇昕桑志勤代文雙
    作物學(xué)報(bào) 2019年2期
    關(guān)鍵詞:共表達(dá)開(kāi)花生物學(xué)

    楊宇昕 桑志勤,2 許 誠(chéng) 代文雙 鄒 棖,*

    ?

    利用WGCNA進(jìn)行玉米花期基因共表達(dá)模塊鑒定

    楊宇昕1桑志勤1,2許 誠(chéng)1代文雙1鄒 棖1,*

    1中國(guó)農(nóng)業(yè)科學(xué)院作物科學(xué)研究所, 北京 100081;2新疆農(nóng)墾科學(xué)院, 新疆石河子 832000

    權(quán)重基因共表達(dá)網(wǎng)絡(luò)分析(weighted gene co-expression network analysis, WGCNA)是系統(tǒng)生物學(xué)的一種研究方法, 在挖掘生物學(xué)數(shù)據(jù)與特定性狀之間的生物學(xué)關(guān)系方面具有十分重要的作用。本研究利用玉米(L.)自交系B73的14份不同發(fā)育階段的轉(zhuǎn)錄組數(shù)據(jù), 篩選掉低表達(dá)豐度的基因, 最終得到了22,426個(gè)高表達(dá)的基因用于創(chuàng)建基因表達(dá)矩陣; 利用不同組織作為性狀, 創(chuàng)建表型矩陣。然后利用R軟件中的WGCNA包建立了共表達(dá)網(wǎng)絡(luò), 共得到20個(gè)模塊。本研究將與組織相關(guān)性高于0.65的模塊定義為組織特異性模塊, 最終鑒定到14個(gè)組織特異性模塊。利用在線網(wǎng)站Agrigo對(duì)組織特異性模塊中的基因進(jìn)行GO (gene ontology)富集分析, 發(fā)現(xiàn)14個(gè)模塊中均可以得到富集種類(lèi)。開(kāi)花作為玉米生育周期中的一個(gè)重要生理過(guò)程, 不僅代表著植物從營(yíng)養(yǎng)生長(zhǎng)到生殖生長(zhǎng)的轉(zhuǎn)變, 也關(guān)系到產(chǎn)量、株高和抗逆性等農(nóng)藝性狀。本研究發(fā)現(xiàn)8個(gè)組織特異性模塊中的基因可以富集到與開(kāi)花調(diào)控的代謝通路。此外, 有17個(gè)已經(jīng)報(bào)道過(guò)的開(kāi)花時(shí)間調(diào)控基因存在于共表達(dá)模塊中, 并且主要分布在Blue模塊和Darkmagenta模塊, 因此本研究重點(diǎn)關(guān)注了這2個(gè)模塊內(nèi)部的基因調(diào)控網(wǎng)絡(luò)。本研究通過(guò)計(jì)算不同組織中的基因表達(dá)豐度, 并聯(lián)合權(quán)重基因共表達(dá)網(wǎng)絡(luò)分析的方法, 鑒定到了具有生物學(xué)意義的共表達(dá)基因模塊, 挖掘到了數(shù)個(gè)開(kāi)花相關(guān)的模塊, 有助于揭示玉米開(kāi)花調(diào)控的遺傳機(jī)制。

    玉米; 權(quán)重基因共表達(dá)網(wǎng)絡(luò); 發(fā)育調(diào)控; 開(kāi)花; 轉(zhuǎn)錄組

    近年來(lái), 隨著高通量測(cè)序技術(shù)的飛速發(fā)展與成本不斷的降低, 越來(lái)越多的研究人員采用多樣本的轉(zhuǎn)錄組測(cè)序進(jìn)行生物學(xué)問(wèn)題的研究。用多樣本的轉(zhuǎn)錄組數(shù)據(jù)可以系統(tǒng)地研究復(fù)雜條件下(例如多個(gè)組織、多個(gè)環(huán)境、多試驗(yàn)變量)的生物學(xué)問(wèn)題, 然而傳統(tǒng)的兩個(gè)樣本的相關(guān)性比較方法已經(jīng)不能有效的處理海量的生物學(xué)數(shù)據(jù)。面對(duì)這個(gè)問(wèn)題, 網(wǎng)絡(luò)研究的方法脫穎而出, 當(dāng)前生物學(xué)研究中廣泛應(yīng)用的網(wǎng)絡(luò)研究方法有基因網(wǎng)絡(luò)[1]、蛋白質(zhì)調(diào)控網(wǎng)絡(luò)[2]和代謝網(wǎng)絡(luò)[3]等。借助于基因組學(xué)和轉(zhuǎn)錄組學(xué)的快速發(fā)展, 基因調(diào)控網(wǎng)絡(luò)在理論和實(shí)踐層面都取得了較快的發(fā)展。在基因調(diào)控網(wǎng)絡(luò)中, 應(yīng)用較廣泛的是權(quán)重基因共表達(dá)網(wǎng)絡(luò)。前人研究指出, 權(quán)重基因共表達(dá)網(wǎng)絡(luò)內(nèi)部的基因聯(lián)通符合無(wú)尺度的網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu), 即網(wǎng)絡(luò)中少數(shù)基因作為網(wǎng)絡(luò)中的核心節(jié)點(diǎn)與大部分基因具有聯(lián)系, 該過(guò)程的實(shí)現(xiàn)依賴(lài)于對(duì)基因間的連接系數(shù)進(jìn)行冪函數(shù)處理, 使網(wǎng)絡(luò)逐漸趨向無(wú)尺度分布[4]。因此借助網(wǎng)絡(luò)分析, 可以有效地利用海量的轉(zhuǎn)錄組數(shù)據(jù), 將眾多基因劃分為相似的基因表達(dá)模塊進(jìn)行研究。權(quán)重基因共表達(dá)網(wǎng)絡(luò)分析(weighted gene co-expression network analysis, WGCNA)根據(jù)基因表達(dá)水平的不同將基因聚類(lèi), 具有相似表達(dá)模式的基因被劃分到同一個(gè)模塊中。其目的在于研究共表達(dá)模塊和目標(biāo)性狀之間的生物學(xué)相關(guān)性、鑒定組織特異性模塊并研究模塊內(nèi)的核心基因。WGCNA作為一種系統(tǒng)生物學(xué)方法, 已經(jīng)廣泛應(yīng)用于多樣本的生物學(xué)問(wèn)題研究當(dāng)中。例如在動(dòng)物中, 研究人員利用四個(gè)年齡段的恒河猴[5](L.), 分別選取8個(gè)組織進(jìn)行轉(zhuǎn)錄組測(cè)序, 通過(guò)將不同的lncRNA和mRNA分為不同的模塊, 闡釋了lncRNA-lncRNA, lncRNA-mRNA基因間的共表達(dá)關(guān)系。WGCNA也廣泛應(yīng)用在植物研究中, 為了研究油菜(L.)在干旱處理?xiàng)l件下的mRNA的表達(dá)模式[6], 研究人員利用干旱和對(duì)照處理下的48個(gè)轉(zhuǎn)錄組數(shù)據(jù)進(jìn)行了共表達(dá)網(wǎng)絡(luò)構(gòu)建, 進(jìn)一步鑒定到與干旱處理相關(guān)的生物學(xué)模塊, 并且挖掘了目標(biāo)模塊內(nèi)的核心基因。在對(duì)草莓(L.)花的研究中[7], 研究者選取了6個(gè)組織, 12個(gè)發(fā)育時(shí)期的樣品, 進(jìn)行了轉(zhuǎn)錄組測(cè)序進(jìn)而得到了共表達(dá)網(wǎng)絡(luò), 在特異性模塊中鑒定到了大量的轉(zhuǎn)錄因子作為模塊中的核心基因, 對(duì)促進(jìn)花芽分化具有重要的作用, 此外該研究還表明WGCNA比K-MEAN聚類(lèi)更具有特異性, 這是因?yàn)閃GCNA構(gòu)建的網(wǎng)絡(luò)符合無(wú)尺度拓?fù)浞植? 因此可以得到更多的具有生物學(xué)意義的模塊。

    玉米開(kāi)花是一個(gè)多基因控制的數(shù)量性狀, 代表了玉米從生殖生長(zhǎng)到營(yíng)養(yǎng)生長(zhǎng)的轉(zhuǎn)變。開(kāi)花是玉米生育周期中的重要農(nóng)藝性狀, 關(guān)系到葉片數(shù)量、產(chǎn)量、營(yíng)養(yǎng)價(jià)值以及抗逆性等性狀。近年來(lái), 研究人員利用分子生物學(xué)的方法, 已經(jīng)在玉米中鑒定到了數(shù)十個(gè)開(kāi)花調(diào)控基因。例如, 通過(guò)數(shù)量性狀位點(diǎn)定位, 已經(jīng)克隆了數(shù)個(gè)開(kāi)花調(diào)控基因, 包括[8]、[9]和[10]。研究人員結(jié)合突變體分析和轉(zhuǎn)基因技術(shù)鑒定到了在玉米中的同源基因[11], 研究結(jié)果表明是一個(gè)開(kāi)花促進(jìn)因子, 該基因在葉片中表達(dá)并通過(guò)韌皮部運(yùn)送到玉米的頂端分生組織, 在頂端分生組織中與轉(zhuǎn)錄因子結(jié)合, 進(jìn)而激活開(kāi)花促進(jìn)因子[12], 而的表達(dá)則標(biāo)志著玉米開(kāi)始從營(yíng)養(yǎng)生長(zhǎng)轉(zhuǎn)變到生殖生長(zhǎng)。此外, 有研究使用5000份重組玉米自交系進(jìn)行巢式關(guān)聯(lián)分析, 鑒定到數(shù)十個(gè)參與光周期調(diào)控的候選基因[13]。隨著玉米參考基因組以及單倍型圖譜的相繼構(gòu)建完成, 利用重測(cè)序、轉(zhuǎn)錄組和代謝組等生物信息學(xué)手段研究玉米的農(nóng)藝性狀已經(jīng)成為一個(gè)有力的手段, 然而前人的研究著眼于單個(gè)調(diào)控因子(如順式表達(dá)元件, 反式作用因子), 鮮有利用網(wǎng)絡(luò)分析的方法解析玉米開(kāi)花遺傳機(jī)制的報(bào)道。

    本研究利用玉米14個(gè)不同組織的RNA-seq數(shù)據(jù), 計(jì)算得到全基因組水平上的基因表達(dá)量, 進(jìn)而構(gòu)建了具有高度拓?fù)湫再|(zhì)的共表達(dá)網(wǎng)絡(luò), 得到了基因共表達(dá)模塊, 并進(jìn)一步篩選到了與組織高度關(guān)聯(lián)的共表達(dá)模塊。通過(guò)對(duì)組織特異性模塊的基因富集分析發(fā)現(xiàn)了具有生物學(xué)功能的模塊。此外, 結(jié)合富集分析的結(jié)果, 鑒定到了數(shù)個(gè)與開(kāi)花調(diào)控相關(guān)的模塊。本研究利用轉(zhuǎn)錄組數(shù)據(jù)得到了開(kāi)花相關(guān)的共表達(dá)模塊, 挖掘了模塊內(nèi)的開(kāi)花相關(guān)的調(diào)控基因, 為玉米開(kāi)花的遺傳機(jī)制解析提供了新的研究思路和實(shí)驗(yàn)結(jié)果。

    1 材料與方法

    1.1 玉米轉(zhuǎn)錄組數(shù)據(jù)獲取與分析

    玉米不同發(fā)育階段的轉(zhuǎn)錄組數(shù)據(jù)來(lái)源于NCBI(national center for biotechnology information)的SRA (Sequence Read Archive)數(shù)據(jù)庫(kù), 獲取碼為PRJNA171684[14]。轉(zhuǎn)錄組數(shù)據(jù)來(lái)源于葉片、頂端分生組織、花絲、花藥、果穗等組織。首先利用SRAtoolkit軟件的fast-dump參數(shù)將高通量測(cè)序數(shù)據(jù)轉(zhuǎn)變?yōu)镕astq文件, 利用FastQC[15]軟件對(duì)原始測(cè)序數(shù)據(jù)進(jìn)行質(zhì)量評(píng)估, 接著利用Trimmomatic軟件去除接頭、低質(zhì)量和未知堿基數(shù)目過(guò)多的reads, 得到clean reads。參考Pertea流程[16], 利用Hisat2軟件將clean reads比對(duì)到玉米B73參考基因組。玉米參考基因組序列和注釋信息下載自ensemble數(shù)據(jù)庫(kù)發(fā)布的玉米基因組信息(ftp://ftp.ensemblgenomes.org/ pub/plants/release-27/genbank/zea_mays/), 利用Stringtie軟件對(duì)轉(zhuǎn)錄本組裝, 之后采用R軟件的Ballgown包計(jì)算基因在每一個(gè)組織中的轉(zhuǎn)錄本表達(dá)量。使用FPKM (fragments per kilobase of exon per million fragments mapped, 每千堿基外顯子百萬(wàn)片段數(shù))值來(lái)衡量基因的表達(dá)水平。

    1.2 共表達(dá)網(wǎng)絡(luò)構(gòu)建

    基因的表達(dá)譜矩陣來(lái)自所有樣本的基因表達(dá)量??紤]到在各個(gè)組織中均低表達(dá)的基因不具有生物學(xué)意義, 在后續(xù)的計(jì)算中還會(huì)增加運(yùn)算量, 因此為了保證共表達(dá)網(wǎng)絡(luò)的特異性和快捷性, 基因在各個(gè)組織中的最大FPKM值小于5的將被過(guò)濾[17-18]。首先, 構(gòu)建共表達(dá)網(wǎng)絡(luò)需要計(jì)算每2個(gè)基因間的相關(guān)系數(shù), 進(jìn)而得到相似矩陣。使用公式表示如下: Smn=cor(xm, xn), S=[Smn], 式中Smn表示基因m和基因n之間的皮爾森相關(guān)系數(shù), S表示相似矩陣。

    R軟件中的WGCNA包提供了一系列函數(shù)用來(lái)構(gòu)建權(quán)重基因共表達(dá)網(wǎng)絡(luò)[19-20]。為了使網(wǎng)絡(luò)符合無(wú)尺度網(wǎng)絡(luò)分布, 需要選擇合適的權(quán)重值, 利用WGCNA包中的函數(shù)pickSoftThreshold計(jì)算權(quán)重值。根據(jù)圖1的結(jié)果, 選擇軟閾值β=14用來(lái)構(gòu)建共表達(dá)網(wǎng)絡(luò), 在該參數(shù)下每個(gè)基因的平均鄰接系數(shù)是280, 中位數(shù)是165, 最大連接系數(shù)是1140。鄰接矩陣的計(jì)算公式為Amn=[(1+Smn)/2]β, 為了清除背景噪音和偽關(guān)聯(lián)帶來(lái)的誤差, 將鄰接矩陣轉(zhuǎn)為拓?fù)渲丿B矩陣(Topological Overlap Matrix, TOM), 同時(shí)利用函數(shù)dissTom=1–TOM對(duì)拓?fù)渲丿B矩陣取逆得到相異性矩陣。利用函數(shù)hclust對(duì)相異矩陣進(jìn)行層次聚類(lèi)。針對(duì)產(chǎn)生的聚類(lèi)樹(shù), 采用動(dòng)態(tài)切割 (Dynamic Tree Cut) 法切割基因聚類(lèi)樹(shù), 相關(guān)參數(shù)為minModuleSize=30, 該過(guò)程可以把表達(dá)模式相近的基因合并在同一個(gè)分支, 每一個(gè)分支代表一個(gè)共表達(dá)模塊(圖2)。

    圖1 軟閾值確定

    圖中的橫軸均代表軟閾值(β)。A: 縱坐標(biāo)對(duì)應(yīng)的是無(wú)尺度網(wǎng)絡(luò)模型指數(shù); B: 每一個(gè)軟閾值對(duì)應(yīng)的網(wǎng)絡(luò)平均連接程度。

    The abscissa represents the soft threshold (β). A: ordinate corresponds to the index of scale free network model; B: the average link degree of each soft threshold.

    1.3 組織特異性模塊鑒定

    具有相似表達(dá)模式的基因被認(rèn)為是具有類(lèi)似的生物學(xué)功能。模塊是一組表達(dá)模式高度相似的基因集合。為了深入研究與組織高度相關(guān)聯(lián)的模塊, 計(jì)算了模塊特征向量(module eigengene, ME)和不同組織之間的相關(guān)系數(shù), 相關(guān)系數(shù)越大, 該模塊與組織間關(guān)聯(lián)度越高, 而相關(guān)系數(shù)值越小, 模塊與組織間關(guān)聯(lián)度低。為了研究具有高度生物學(xué)特性的模塊, 本研究定義相關(guān)系數(shù)的閾值為0.65, 即任何模塊與組織的相關(guān)系數(shù)高于0.65將被定義為組織特異性模塊[21]。

    1.4 模塊基因的富集分析

    首先提取出組織特異性模塊中的基因, 再將模塊基因提交到Agrigo[22](GO Analysis Toolkit and Database for Agricultural Community)(http://bioinfo. cau.edu.cn/agriGO/index.php)在線分析工具, 采用奇異富集分析(singular enrichment analysis, SEA)法進(jìn)行富集分析。使用玉米參考基因組B73作為參考數(shù)據(jù)庫(kù), 利用Fisher’s校驗(yàn)和Bonferronni多重校驗(yàn)來(lái)進(jìn)行顯著GO類(lèi)型的篩選。

    1.5 玉米開(kāi)花調(diào)控通路的相關(guān)模塊

    開(kāi)花作為一個(gè)復(fù)雜的農(nóng)藝性狀, 不僅與外界環(huán)境信號(hào)有關(guān), 同時(shí)也涉及到內(nèi)源性的信號(hào)調(diào)控[23]。根據(jù)已有的文獻(xiàn)報(bào)道, 玉米開(kāi)花受到日照、赤霉素以及DNA的甲基化修飾等因素的調(diào)控。利用開(kāi)花調(diào)控相關(guān)的一些關(guān)鍵詞例如開(kāi)花(flower)、光(light)、赤霉素(gibberellin)、甲基化(methylation)、分生組織(meristem)等, 在組織特異性模塊中查找參與調(diào)控開(kāi)花的模塊與基因。

    1.6 基因互作網(wǎng)絡(luò)的構(gòu)建

    根據(jù)拓?fù)渲丿B矩陣計(jì)算出模塊中不同基因間的權(quán)重值, 權(quán)重值越高, 代表基因間的關(guān)聯(lián)程度越高。篩選出模塊間相關(guān)系數(shù)排名較高的基因, 同時(shí)挑選出模塊中已經(jīng)報(bào)道過(guò)的開(kāi)花基因, 利用這些基因進(jìn)行互作網(wǎng)絡(luò)的構(gòu)建, 網(wǎng)絡(luò)中每一個(gè)節(jié)點(diǎn)代表一個(gè)基因, 連接線兩端的基因具有相似的生物學(xué)功能。

    2 結(jié)果與分析

    2.1 基因共表達(dá)網(wǎng)絡(luò)構(gòu)建

    通過(guò)數(shù)據(jù)分析共獲得39,625個(gè)基因的表達(dá)譜數(shù)據(jù), 過(guò)濾低表達(dá)量的基因, 獲得22,426個(gè)高表達(dá)基因。計(jì)算得到了基因表達(dá)水平的樣本聚類(lèi)和性狀關(guān)聯(lián), 圖3的結(jié)果表明每一個(gè)組織的基因聚類(lèi)樹(shù)可以和組織很好地對(duì)應(yīng)。通過(guò)對(duì)權(quán)重值篩選, 最終選擇β=14來(lái)構(gòu)建網(wǎng)絡(luò), 采用動(dòng)態(tài)剪切樹(shù)法合并表達(dá)相似的模塊, 共獲得20個(gè)共表達(dá)模塊(圖2), 不同的顏色代表不同的模塊。模塊中的基因數(shù)目分配根據(jù)其表達(dá)量進(jìn)行一個(gè)相關(guān)度的聚類(lèi), 聚類(lèi)度較高的基因被分配到一個(gè)模塊中, 其中Black模塊包含的基因數(shù)目最多, 是8492個(gè)基因, Floralwhite模塊內(nèi)的基因數(shù)量最少, 含87個(gè)基因, 平均每個(gè)模塊包含的基因數(shù)目是1123, Grey模塊是一組未分配到其他模塊的基因集合(圖2)。

    圖2 基因聚類(lèi)樹(shù)和模塊切割

    A: 基于拓?fù)渲丿B構(gòu)建的基因聚類(lèi)樹(shù)。B: 動(dòng)態(tài)混合切割法得到的基因模塊, 顏色代表模塊。C: 合并相似表達(dá)模式的基因模塊。

    A: clustering of genes based on the topological overlap. B: the gene modules obtained from the dynamic tree cut. C: the merged modules with similar expression pattern.

    圖3 樣本層次聚類(lèi)樹(shù)及對(duì)應(yīng)的組織信息

    A: 基于歐氏距離得到的基因聚類(lèi)樹(shù)。橫軸表示不同的玉米組織, 縱軸代表基因間的聚類(lèi)高度。B: 性狀與基因聚類(lèi)樹(shù)關(guān)聯(lián)熱圖。白色的代表低關(guān)聯(lián)度, 紅色的代表高關(guān)聯(lián)度。

    A: the gene cluster tree based on the Euclidean distance. The abscissa represents the different tissues. The ordinate represents the cluster height of genes. B: the association trait heatmap, the white means a low correlation, the red means a high correlation.

    圖4 共表達(dá)模塊中基因數(shù)目分布

    橫軸代表模塊, 縱軸代表模塊中的基因數(shù)目。

    The abscissa represents the modules and the vertical ordinate represents the gene number of each module.

    2.2 組織特異性模塊鑒定

    在20個(gè)模塊中有14個(gè)與組織存在高度特異性(圖5)。大多數(shù)組織都有與其高度相關(guān)聯(lián)的模塊, 例如Blue模塊與13葉期葉片存在高度相關(guān)性(=0.96,=6E–08), Darkmagenta與v18雄穗連接度較高(=0.98,=2E–09)。14個(gè)組織特異性模塊中沒(méi)有任何一個(gè)模塊和v18未成熟果穗與授粉前果穗兩個(gè)組織存在高度相關(guān)性。我們推測(cè)這兩個(gè)時(shí)期組織中的基因主要進(jìn)行常規(guī)的生命活動(dòng), 因此沒(méi)有組織特異性的模塊與之對(duì)應(yīng)。

    圖5 模塊與性狀關(guān)聯(lián)熱圖

    橫軸表示不同的性狀, 縱軸表示每一個(gè)模塊的特征向量。紅色的格子代表性狀與模塊具有正相關(guān)性, 綠色的格子代表性狀與模塊具有負(fù)相關(guān)性。

    Each column represents the co-expression module. Red color of each box represents the positive correlation between module and trait. Green color of each box represents the negative relationships between module and trait.

    2.3 組織特異性模塊的功能富集分析

    具有相似表達(dá)模式的基因可能具有相同的生物學(xué)功能, 通過(guò)共表達(dá)網(wǎng)絡(luò)構(gòu)建的模塊在很大程度上與特定的生物組織相聯(lián)系, 而特定的生物組織又與特定生命活動(dòng)存在高度協(xié)同性[4]。為了進(jìn)一步解析組織特異性模塊的生物學(xué)功能, 本研究利用Agrigo工具對(duì)模塊基因進(jìn)行富集分析, 表明這14個(gè)模塊都可以富集得到顯著的GO通路, 這些GO通路涉及到生物學(xué)過(guò)程(biological process, BP), 分子功能 (molecular function, MF)以及細(xì)胞組分(cellular component, CC)。與根發(fā)育高度關(guān)聯(lián)Darkgrey模塊可以富集得到水轉(zhuǎn)運(yùn)(GO:0006833), 流體運(yùn)送(GO:0042044)等調(diào)控通路, 表明WGCNA可以構(gòu)建到具有生物學(xué)意義的共表達(dá)模塊。通過(guò)對(duì)開(kāi)花發(fā)育與調(diào)控的關(guān)鍵字進(jìn)行檢索發(fā)現(xiàn), Blue、Darkmagenta、Darkred、Darkslateblue、Orange、Plum2、Yellow和Bisuqe4模塊都富集到了與開(kāi)花相關(guān)的調(diào)控通路(圖6), 這些通路主要是響應(yīng)光照刺激的調(diào)控通路(GO:0009416)、頂端分生組織形態(tài)保持(GO:0010073)和花器官發(fā)育(GO:0048437)等。Blue模塊與13葉期葉片存在最高的相關(guān)性(=0.96,=6E–08), Darkmagenta模塊和v18雄穗具有最高的相關(guān)性(=0.98,=2E–09), 這兩個(gè)模塊對(duì)應(yīng)的組織是葉片和雄穗, 是和開(kāi)花密切相關(guān)的組織, 此外, Blue模塊包含了7個(gè)已經(jīng)報(bào)道的開(kāi)花基因, Darkmagenta模塊也包含6個(gè)已經(jīng)報(bào)道的開(kāi)花基因, 雖然富集得到了開(kāi)花調(diào)控通路, 但是開(kāi)花基因數(shù)目較少或者沒(méi)有已知開(kāi)花基因存在于這些模塊, 因此本研究重點(diǎn)關(guān)注了Blue模塊(圖7)和Darkmagenta (圖8)模塊, 篩選出模塊中高連通性的基因, 利用Cytoscape[24]軟件對(duì)網(wǎng)絡(luò)進(jìn)行了可視化處理, 處于網(wǎng)絡(luò)中的基因是該模塊中的高連接度的基因。

    2.4 開(kāi)花基因互作網(wǎng)絡(luò)構(gòu)建

    考慮到開(kāi)花時(shí)間調(diào)控基因主要分布在Blue和Darkmagenta模塊, 因此利用這2個(gè)模塊中已經(jīng)報(bào)道的開(kāi)花基因構(gòu)建基因互作網(wǎng)絡(luò)。Blue模塊中共發(fā)現(xiàn)7個(gè)報(bào)道過(guò)的開(kāi)花基因, 其中、[25]和[26]3個(gè)基因的連通性在模塊排名為前1%, 可以作為模塊中的核心基因,和作為開(kāi)花促進(jìn)因子, 是開(kāi)花調(diào)控通路中的重要基因。此外[27]、[28-29]、[30]和[31]基因在模塊中聯(lián)通性排名靠后, 但是這些基因被證明是玉米開(kāi)花調(diào)控通路中關(guān)鍵基因, 推測(cè)WGCNA構(gòu)建的模塊中, 基因的連通性基于基因間表達(dá)量, 而基因間表達(dá)量受多種因素影響, 導(dǎo)致關(guān)鍵調(diào)控基因在模塊中的相關(guān)性較低。因此, 本研究以、、、、、、等基因作為樞紐基因, 篩選與開(kāi)花基因高連通性的基因, 將其定義為候選開(kāi)花基因, 利用開(kāi)花基因和候選花期基因構(gòu)建局部調(diào)控網(wǎng)絡(luò)(圖9), 在該網(wǎng)絡(luò)中, 每一個(gè)節(jié)點(diǎn)代表一個(gè)基因, 節(jié)點(diǎn)間通過(guò)連接線相聯(lián)系, 處于連接線兩端的基因被認(rèn)為具有相同生物功能。與v18雄穗高度關(guān)聯(lián)的Darkmagenta模塊包含6個(gè)開(kāi)花調(diào)控基因, 即[32]、[33]、[32]、[34]、[34]和[35], 在該模塊中有2個(gè)基因的連通性排名前10%, 分別是和, 參照Blue模塊的研究思路, 利用Darkmagenta已經(jīng)報(bào)道過(guò)的開(kāi)花基因構(gòu)建Darkmagenta模塊內(nèi)的局部調(diào)控網(wǎng)絡(luò)(圖10)。此外, 檢索了候選開(kāi)花基因在擬南芥中的同源基因, 并借助TAIR (https://www.arabidopsis.org/)網(wǎng)站注釋了擬南芥同源基因的功能(表1)。

    圖6 開(kāi)花富集通路

    每一行代表開(kāi)花富集通路, 每一列代表組織特異性模塊。點(diǎn)的大小代表多重校驗(yàn)的值大小, 點(diǎn)的顏色表示輸入基因與背景基因的比值。

    Each row corresponds to a GO enrichment pathway, column to a tissue specific module. The point size is calculated by the-value of multiple check. The Richfactor represent the ratio of input genes to the background genes.

    圖7 Blue模塊內(nèi)的基因共表達(dá)網(wǎng)絡(luò)

    圖8 Darkmagenta模塊內(nèi)基因共表達(dá)網(wǎng)絡(luò)

    3 討論

    傳統(tǒng)的生物學(xué)研究側(cè)重于從分子水平闡釋單個(gè)功能元件(如DNA、mRNA和蛋白質(zhì))對(duì)生命活動(dòng)的影響, 雖然該方法對(duì)于揭示具體性狀的遺傳機(jī)理具有非常重要的意義, 但只能局部地解釋某一生命活動(dòng)發(fā)生的原因。隨著測(cè)序技術(shù)的快速發(fā)展, 傳統(tǒng)的生物學(xué)研究不能充分有效地挖掘海量數(shù)據(jù)中蘊(yùn)含的生物學(xué)意義。網(wǎng)絡(luò)作為系統(tǒng)生物學(xué)的一種研究手段, 借助基因組、轉(zhuǎn)錄組和代謝組的數(shù)據(jù), 廣泛地應(yīng)用于生命科學(xué)的探索中。相比于其他調(diào)控網(wǎng)絡(luò), WGCNA可以特異地篩選出與性狀相關(guān)的基因, 并進(jìn)行模塊化分類(lèi), 得到具有高度生物學(xué)意義的共表達(dá)模塊, 已經(jīng)被證明是一種高效的數(shù)據(jù)挖掘手段[36]。

    本研究通過(guò)對(duì)14個(gè)組織特異性模塊進(jìn)行富集分析, 發(fā)現(xiàn)這14個(gè)模塊均可以得到具有生物學(xué)意義的調(diào)控通路(附表1)。例如在與根高度關(guān)聯(lián)的Darkgrey模塊, 富集分析的GO term包含水轉(zhuǎn)運(yùn)(GO:0006833)、流體運(yùn)輸(GO:0042044)等通路, 這些通路和根部所具有的生物功能高度吻合, 該結(jié)果證明了共表達(dá)網(wǎng)絡(luò)構(gòu)建模塊的特異性。本研究直接利用轉(zhuǎn)錄組數(shù)據(jù)對(duì)應(yīng)的組織, 例如根、葉片和花藥等作為研究的性狀, 該方法在前人的研究中已有過(guò)相關(guān)報(bào)道, 如在黃瓜(L.) 10個(gè)組織的共表達(dá)網(wǎng)絡(luò)構(gòu)建中, 研究人員成功鑒定到數(shù)個(gè)與黃瓜苦味合成代謝相關(guān)的模塊[18]。

    表1 Darkmagenta模塊和Blue模塊中候選開(kāi)花基因的功能注釋。

    本研究重點(diǎn)關(guān)注了開(kāi)花這個(gè)重要的農(nóng)藝性狀。開(kāi)花作為玉米生命活動(dòng)的中心環(huán)節(jié), 對(duì)其遺傳機(jī)制的解析雖然取得了一定的進(jìn)展, 例如通過(guò)數(shù)量性狀位點(diǎn)定位和關(guān)聯(lián)分析等方法鑒定了許多重要的開(kāi)花基因, 但是這些研究手段都有一定的局限性, 例如數(shù)量性狀定位不能反映參考群體較為廣泛的遺傳異質(zhì)性[37], 此外, 許多復(fù)雜的農(nóng)藝性狀是一個(gè)與群體結(jié)構(gòu)存在高度相關(guān)的性狀, 例如開(kāi)花, 而且利用關(guān)聯(lián)分析鑒定控制目標(biāo)性狀相關(guān)的調(diào)控基因容易受到群體結(jié)構(gòu)的影響[38]。權(quán)重基因共表達(dá)網(wǎng)絡(luò)則通過(guò)共表達(dá)模塊化處理, 將成千上萬(wàn)的基因分配到模塊中, 通過(guò)研究模塊蘊(yùn)含的生物學(xué)意義進(jìn)一步挖掘相應(yīng)基因的功能, 因此較好地解決了對(duì)復(fù)雜性狀解析能力不足的問(wèn)題, 被廣泛地應(yīng)用于基因組學(xué)的研究之中。

    圖9 Blue模塊開(kāi)花相關(guān)的基因網(wǎng)絡(luò)

    紅色節(jié)點(diǎn)是參與開(kāi)花調(diào)控通路的基因。

    Each node corresponding to a gene, the red node is the hub gene of the network and is the flowering time gene which had reported.

    在本研究中的富集分析結(jié)果中, 8個(gè)組織特異性模塊包含開(kāi)花的調(diào)控通路(圖6), 例如Blue模塊中的響應(yīng)光刺激調(diào)控通路(GO:0009416), Darkolivegreen模塊中的頂端分生組織形態(tài)維持(GO:0010022)和Plum2模塊的花器官發(fā)育(GO:0048437)。此外, 結(jié)合已經(jīng)報(bào)道過(guò)的開(kāi)花基因, 進(jìn)一步發(fā)現(xiàn)Blue模塊和Darkmagenta模塊含有較多開(kāi)花相關(guān)的調(diào)控基因, 而且這些開(kāi)花基因的連通性在模塊中很高, 處于樞紐地位, 例如Blue模塊中的和基因是開(kāi)花激活因子, Darkmagenta模塊中的和基因, 是時(shí)鐘通路的重要參與基因。通過(guò)對(duì)模塊中基因的局部網(wǎng)絡(luò)化, 挖掘到了與開(kāi)花基因連接程度較高的基因, 這些基因可以作為開(kāi)花候選基因被進(jìn)一步研究。例如在Blue模塊中, 與8高連通性的基因, 其在擬南芥中的同源基因是,基因編碼一個(gè)調(diào)節(jié)開(kāi)花的小蛋白質(zhì), 并且參與赤霉素調(diào)控通路, 在光照誘導(dǎo)后該基因在頂端分生組織表達(dá),的高連接度基因在擬南芥中的同源基因是, 通過(guò)功能注釋發(fā)現(xiàn)這個(gè)基因編碼葉綠素結(jié)合蛋白D1, 屬于光系統(tǒng)II反應(yīng)中心, 參與光周期調(diào)控通路。在Darkmagenta模塊中,的高連接度基因是, 其在擬南芥中的同源基因是, 該基因編碼一個(gè)F-box蛋白, 響應(yīng)光照刺激, 參與光周期途徑;是的高連接度基因, 該基因在擬南芥中的同源基因是, 主要參與花粉管生長(zhǎng)和受精。以上結(jié)果表明, 在WGCNA中, 與目標(biāo)基因存在高連通性的基因具有與其相似的生物學(xué)意義, 這可以為挖掘目標(biāo)基因提供新的研究思路。在本研究中, 重點(diǎn)關(guān)注了Blue和Darkmagenta兩個(gè)開(kāi)花相關(guān)的調(diào)控模塊, 其余6個(gè)開(kāi)花相關(guān)的基因模塊雖然沒(méi)有被詳細(xì)的論述, 但是其包含開(kāi)花富集通路, 可進(jìn)一步挖掘其蘊(yùn)含的生物學(xué)意義。

    圖10 Darkmagenta模塊的開(kāi)花相關(guān)的基因網(wǎng)絡(luò)

    紅色節(jié)點(diǎn)基因參與開(kāi)花調(diào)控通路。

    Each node corresponding to a gene, the red node is the hub gene of the network and is the flowering time gene which had reported.

    本研究鑒定到的組織特異性模塊對(duì)應(yīng)的組織都是玉米開(kāi)花調(diào)控和發(fā)育的重要參與器官, 因此, 利用共表達(dá)網(wǎng)絡(luò)解析具體的農(nóng)藝性狀, 分析模塊的生物功能, 挖掘模塊中的目標(biāo)基因, 將為解析復(fù)雜的農(nóng)藝性狀提供重要的參考依據(jù)。

    4 結(jié)論

    構(gòu)建了與性狀相關(guān)聯(lián)的權(quán)重基因共表達(dá)網(wǎng)絡(luò), 得到了14個(gè)組織特異性模塊。揭示了組織特異性模塊蘊(yùn)含的生物學(xué)意義, 鑒定了開(kāi)花相關(guān)的調(diào)控基因, 并用其構(gòu)建了局部的具有生物學(xué)意義的網(wǎng)絡(luò)。本研究結(jié)果可為后續(xù)解析玉米發(fā)育的遺傳機(jī)理提供重要的實(shí)驗(yàn)數(shù)據(jù)和理論基礎(chǔ)。

    [1] Stuart J M, Segal E, Koller D, Kim S K. A gene-coexpression network for global discovery of conserved genetic modules., 2003, 302: 249–255.

    [2] Jeong H, Mason S P, Barabási A L, Oltvai Z N. Lethality and centrality in protein networks., 2001, 411: 41.

    [3] Gille C, Hoffmann S, Holzhütter H G. METANNOGEN: compiling features of biochemical reactions needed for the reconstruction of metabolic networks., 2007, 1: 5.

    [4] Barabási A L, Oltvai Z N. Network biology: understanding the cell's functional organization., 2004, 5: 101–113.

    [5] Liu S, Wang Z, Chen D, Zhang B, Tian R R, Wu J, Zhang Y, Xu K Y, Yang L M, Cheng C, Ma J, Lv L B, Zheng Y T, Hu X T, Yi Z, Wang X T, Li J L. Annotation and cluster analysis of spatiotemporal- and sex-related lncRNA expression in rhesus macaque brain., 2017, 27: 1608–1620.

    [6] Greenham K, Guadagno C R, Gehan M A, Mockler T C, Weinig C, Ewers B E, McClung C R. Temporal network analysis identifies early physiological and transcriptomic indicators of mild drought in., 2017, 6: e29655.

    [7] Hollender C A, Kang C, Darwish O, Geretz A, Matthews B F, Slovin J, Alkharouf N, Liu Z. Floral transcriptomes in woodland strawberry uncover developing receptacle and anther gene networks., 2014,165: 1062–1075.

    [8] Vl?du?u C, McLaughlin J, Phillips R L. Fine mapping and characterization of linked quantitative trait loci involved in the transition of the maize apical meristem from vegetative to generative structures., 1999, 153: 993–1007.

    [9] Wong A Y, Colasanti J. Maize floral regulator protein INDETERMINATE1 is localized to developing leaves and is not altered by light or the sink/source transition., 2007, 58: 403–414.

    [10] Muszynski M G, Dam T, Li B, Shirbroun D M, Hou Z, Bruggemann E, Archibald R, Ananiev E V, Danilevskaya O N. Delayed flowering1 encodes a basic leucine zipper protein that mediates floral inductive signals at the shoot apex in maize., 2006, 142: 1523–1536.

    [11] Meng X, Muszynski M G, Danilevskaya O N. The-likegene functions as a floral activator and is involved in photoperiod sensitivity in maize., 2011, 23: 942–960.

    [12] Danilevskaya O N.encodes a basic leucine zipper protein that mediates floral inductive signals at the shoot apex in maize., 2006, 142: 1523–1536.

    [13] Coles N D, McMullen M D, Balint-Kurti P J, Pratt R C, Holland J B. Genetic control of photoperiod sensitivity in maize revealed by joint multiple population analysis., 2010, 184: 799–812.

    [14] Sekhon R S, Lin H, Childs K L, Hansey C N, Buell C R, de Leon N, Kaeppler S M. Genome-wide atlas of transcription during maize development., 2011, 66: 553–563.

    [15] Kroll K W, Mokaram N E, Pelletier A R, Frankhouser D E, Westphal M S, Stump P A, Stump C L, Bundschuh R, Blachly J S, Yan P. Quality Control for RNA-Seq (QuaCRS): an integrated quality control pipeline., 2014, 13: 7–14.

    [16] Pertea M, Kim D, Pertea G M, Leek J T, Salzberg S L. Transcript-level expression analysis of RNA-seq experiments with HISAT, StringTie and Ballgown., 2016, 11: 1650–1667.

    [17] 魏凱, 張婷婷, 馬磊. 豬基因共表達(dá)網(wǎng)絡(luò)模塊的構(gòu)建及功能分析. 畜牧獸醫(yī)學(xué)報(bào), 2017, 48: 2205–2215.Wei K, Zhang T T, Ma L. Construction and functional analysis of gene co-expression network modules., 2017, 48: 2205–2215 (in Chinese with English abstract).

    [18] 林行眾, 張忠華, 楊清, 黃三文. 黃瓜共表達(dá)基因模塊的識(shí)別及其特點(diǎn)分析. 農(nóng)業(yè)生物技術(shù)學(xué)報(bào), 2017, 23: 1121–1130. Lin X Z, Zhang Z H, Yang Q, Huang S W. Identification and characterization analysis of co-expression gene modules in cucumber (L.)., 2017, 23: 1121–1130 (in Chinese with English abstract).

    [19] Zhang B, Kirov S, Snoddy J. WebGestalt: an integrated system for exploring gene sets in various biological contexts., 2005, 33: W741–W748.

    [20] Langfelder P, Horvath S. WGCNA: an R package for weighted correlation network analysis., 2008, 9: 559.

    [21] Downs G S, Bi Y M, Colasanti J, Wu W, Chen X, Zhu T, Rothstein S J, Lukens L N. A developmental transcriptional network for maize defines coexpression modules., 2013, 161: 1830–1843.

    [22] Du Z, Zhou X, Ling Y, Zhang Z H, Su Z. agriGO: a GO analysis toolkit for the agricultural community., 2010, 38: 64–70

    [23] Dong Z S, Danilevskaya O, Abadie T, Messina C, Coles N, Cooper M. A gene regulatory network model for floral transition of the shoot apex in maize and its dynamic modeling., 2012, 7: e43450.

    [24] Shannon P, Markiel A, Ozier O, Baliga N S, Wang J T, Ramage D, Amin N, Schwikowski B, Ideker T. Cytoscape: a software environment for integrated models of biomolecular interaction networks., 2003, 13: 2498–2504.

    [25] Mascheretti I, Turner K, Brivio R S, Hand A, Colasanti J, Rossi V. Florigen-Encoding genes of day-neutral and photoperiod- sensitive maize are regulated by different chromatin modifications at the floral transition., 2015, 168: 1351–1363.

    [26] Khan S, Rowe S C, Harmon F G. Coordination of the maize transcriptome by a conserved circadian clock., 2010, 10: 126.

    [27] Sheehan M J, Kennedy L M, Costich D E, Brutnell T P. Subfunctionalization ofandin the control of seedling and mature plant traits in maize., 2007, 49: 338–353.

    [28] Thornsberry J M, Goodman M M, Doebley J, Kresovich S, Nielsen D, Buckler E S. Dwarf8 polymorphisms associate with variation in flowering time., 2001, 28: 286.

    [29] Larsson S J, Lipka A E, Buckler E S. Lessons fromon the strengths and weaknesses of structured association mapping., 2013, 9: e1003246.

    [30] Lawit S J, Wych H M, Xu D, Kundu S, Tomes D T. Maize DELLA proteins dwarf plant8 and dwarf plant9 as modulators of plant development., 2010, 51: 1854–1868.

    [31] Wang X, Wu L, Zhang S, Wu L, Ku L, Wei X, Xie L, Chen Y. Robust expression and association ofwith circadian rhythms in maize., 2011, 30: 1261–1272.

    [32] Miller T A, Muslin E H, Dorweiler J E. A maize CONSTANS-like gene,, exhibits distinct diurnal expression patterns in varied photoperiods., 2008, 227: 1377–1388.

    [33] Sheehan M J, Farmer P R, Brutnell T P. Structure and expression of maize phytochrome family homeologs., 2004, 167: 1395–1405

    [34] Hayes K R, Beatty M, Meng X, Simmons C R, Habben J E, Danilevskaya O N. Maize global transcriptomics reveals pervasive leaf diurnal rhythms but rhythms in developing ears are largely limited to the core oscillator., 2010, 5: e12887.

    [35] Chardon F, Virlon B, Moreau L, Falque M, Joets J, Decousset L, Murigneux A, Charcosset A. Genetic architecture of flowering time in maize as inferred from quantitative trait loci meta- analysis and synteny conservation with the rice genome., 2004, 168: 2169–2185.

    [36] Zhao W, Langfelder P, Fuller T, Dong J, Li A, Hovarth S. Weighted gene coexpression network analysis: state of the art., 2010, 20: 281–300.

    [37] Holland J B. Genetic architecture of complex traits in plants., 2007, 10: 156–161.

    [38] Camus-Kulandaivelu L, Veyrieras J B, Madur D, Combes V, Fourmann M, Barraud S, Dubreuil P, Gouesnard B, Manicacci D, Charcosset A. Maize adaptation to temperate climate: relationship between population structure and polymorphism in thegene., 2006, 172: 2449–2463.

    附表1 組織特異性模塊的GO富集分析結(jié)果(部分)

    Supplementary table 1 GO enrichment analysis results of tissue specific module (part)

    模塊ModuleGO條目GO term本體Ontology描述DescriptionP值P-value Darkorange2GO: 0009889BP生物合成過(guò)程調(diào)節(jié)Regulation of biosynthetic process6.90E–07 Darkorange2GO: 0006355BPDNA依賴(lài)的轉(zhuǎn)錄調(diào)節(jié)Regulation of transcription, DNA-dependent4.40E–07 Darkorange2GO: 0015267MF通道活性Channel activity0.00068 Darkorange2GO: 0022838MF底物特異性通道活性Substrate-specific channel activity0.00068 BlueGO: 0009628BP非生物刺激的反應(yīng)Response to abiotic stimulus3.60E–10 BlueGO: 0009416BP光刺激響應(yīng)Response to light stimulus1.60E–08 BlueGO: 0003700MF轉(zhuǎn)錄因子活性Transcription factor activity0.00017 BlueGO: 0009535CC葉綠體類(lèi)囊體膜Chloroplast thylakoid membrane0.00022 DarkredGO: 0051186BP輔因子代謝過(guò)程Cofactor metabolic process1.20E–19 DarkredGO: 0004252MF絲氨酸型肽鏈內(nèi)切酶活性Serine-type endopeptidase activity0.00015 DarkredGO: 0009543CC葉綠體類(lèi)囊體腔 Chloroplast thylakoid lumen4.40E–13 DarkslateblueGO: 0016051BP碳水化合物生物合成過(guò)程Carbohydrate biosynthetic process5.60E–05 DarkslateblueGO: 0016830MF碳-碳裂解酶活性Carbon-carbon lyase activity3.40E–06 DarkslateblueGO: 0042651CC類(lèi)囊體膜Thylakoid membrane0.00013 TurquoiseGO: 0060560BP形態(tài)發(fā)生發(fā)育Developmental growth involved in morphogenesis1.70E–11 TurquoiseGO: 0015299MF溶質(zhì): 氫反向轉(zhuǎn)運(yùn)蛋白活性Solute: hydrogen antiporter activity4.60E–05 TurquoiseGO: 0031224CC內(nèi)膜Intrinsic to membrane8.50E–06 DarkmagentaGO: 0010927BP涉及形態(tài)發(fā)生的內(nèi)膜組裝Cellular component assembly involved in morphogenesis4.10E–09 DarkmagentaGO: 0004553MF水解酶活性, 水解O-糖基化合物Hydrolase activity, hydrolyzing O-glycosyl compounds7.60E–06 DarkmagentaGO: 0030312CC高爾基體Golgi apparatus0.0004 Bisque4GO: 0042180BP細(xì)胞酮代謝過(guò)程Cellular ketone metabolic process3.00E–08 Bisque4GO: 0044281BP小分子代謝過(guò)程Small molecule metabolic process9.90E–06 Bisque4GO: 0042221BP相應(yīng)化學(xué)刺激Response to chemical stimulus0.00032 DarkgreyGO: 0010033BP對(duì)有機(jī)物質(zhì)的反應(yīng)Response to organic substance1.60E–07 DarkgreyGO: 0004553MF水解酶活性, 水解O-糖基化合物Hydrolase activity, hydrolyzing O-glycosyl compounds4.40E–10 DarkgreyGO: 0005740CC線粒體包膜Mitochondrial envelope7.40E–05 FloralwhiteGO: 0006412BP蛋白質(zhì)翻譯Translation9.10E–07 FloralwhiteGO: 0005198MF結(jié)構(gòu)分子活性Structural molecule activity5.80E–08 FloralwhiteGO: 0043227CC膜有界細(xì)胞器Membrane-bounded organelle0.0085 DarkolivegreenGO: 0010022BP分生組織決定Meristem determinacy6.10E–07

    (附續(xù)表1)

    模塊ModuleGO條目GO term本體Ontology描述DescriptionP值P-value DarkolivegreenGO: 0003700MF轉(zhuǎn)錄因子活性Transcription factor activity5.50E–05 OrangeGO: 0042542BP過(guò)氧化氫反應(yīng)Response to hydrogen peroxide2.90E–05 OrangeGO: 0009526CC質(zhì)體Plastid9.60E–09 YellowGO: 0009605BP響應(yīng)外界刺激Response to external stimulus8.40E–10 YellowGO: 0022892MF底物特異性轉(zhuǎn)運(yùn)體活性Substrate-specific transporter activity0.00013 Plum2GO: 0048437BP花器官發(fā)育Floral organ development9.20E–05 Plum2GO: 0005576CC胞外區(qū)Extracellular region9.70E–05 GreenyellowGO: 0034220BP離子跨膜轉(zhuǎn)運(yùn)Ion transmembrane transport5.00E–13 GreenyellowGO: 0032561MF鳥(niǎo)苷酸核糖核酸結(jié)合Guanyl ribonucleotide binding2.40E–10 GreenyellowGO: 0043231CC內(nèi)細(xì)胞器膜Organelle inner membrane1.20E–17

    Identification of maize flowering gene co-expression modules by WGCNA

    YANG Yu-Xin1, SANG Zhi-Qin1,2, XU Cheng1, DAI Wen-Shuang1, and ZOU Cheng1,*

    1Institute of Crop Sciences, Chinese Academy of Agricultural Sciences, Beijing 100081, China;2Xinjiang Academy of Agricultural and Reclamation Science, Shihezi 832000, Xinjiang, China

    Weighted gene co-expression network analysis (WGCNA) is one of the research methods in systematic biology. It can effectively analyze the complex samples, and has been extensively used in the analysis of complicated traits for many samples. Weighted gene co-expression network has the characteristics of scale-free distribution and could construct the scale free network. The genes with similar expression level can be clustered and assigned to a module, then the relationships between co-expression modules and specific tissues can be furtherly analyzed. Our research utilized the transcriptome data of 14 different tissues of maize (L.) inbred line B73, and calculated the gene expression level of the whole genome. Through filtering out the genes with low expression level we finally got 22,426 genes with high expression level to construct the gene expression matrix. We utilized the different tissues as the trait to construct the trait matrix. The weighted gene co-expression network analysis packages of R software was used to perform the co-expression network analysis, and 20 co-expression modules were identified. We finally obtained 14 tissue specific modules which were highly correlated with traits (> 0.65). The enrichment analysis tool Agrigo was taken to perform the GO enrichment of the tissue specific module genes, all the 14 tissues could be enriched in GO terms. Flowering is one of the important agronomic traits in the life cycle of maize controlled by external environment signals and genetic factors. Maize flowering not only represents the transition from the vegetative growth to reproductive growth, also relates to grain yield, plant height and resistance. In our research, we detected eight tissue specific modules, which could be obtained within flowering time related pathways. In addition, 17 flowering genes which have been reported in the literatures were assigned to the co-expression modules, and mainly assigned to the Blue and Darkmagenta modules. Therefore, we focused on the network of Blue and Darkmagenta modules. Our research calculated the gene expression abundances, and detected several flowering time related modules, which will contribute to revealing the genetic mechanism of maize flowering time regulation.

    maize; weighted gene co-expression network; development regulation; flowering; transcriptome

    2018-07-17;

    2018-10-08;

    2018-11-08.

    10.3724/SP.J.1006.2019.83053

    鄒棖,E-mail: zoucheng@caas.cn

    E-mail: yyx0719@126.com

    本研究由國(guó)家重點(diǎn)研發(fā)計(jì)劃項(xiàng)目(2016YFD0100303)和國(guó)家自然科學(xué)基金項(xiàng)目(31371638)資助。

    This study was supported by the National Key Research and Development Program of China (2016YFD0100303), the National Natural Science Foundation of China (31371638).

    URL:http://kns.cnki.net/kcms/detail/11.1809.S.20181105.0933.004.htm

    猜你喜歡
    共表達(dá)開(kāi)花生物學(xué)
    谷稗的生物學(xué)特性和栽培技術(shù)
    侵襲性垂體腺瘤中l(wèi)ncRNA-mRNA的共表達(dá)網(wǎng)絡(luò)
    初中生物學(xué)糾錯(cuò)本的建立與使用
    《一棵開(kāi)花的樹(shù)》
    金橋(2019年12期)2019-08-13 07:16:40
    雨開(kāi)花
    文苑(2019年14期)2019-08-09 02:14:06
    一雙開(kāi)花的舊靴
    大灰狼(2019年5期)2019-05-29 17:45:26
    初中生物學(xué)糾錯(cuò)本的建立與使用
    膀胱癌相關(guān)lncRNA及其共表達(dá)mRNA的初步篩選與功能預(yù)測(cè)
    中國(guó)流行株HIV-1gag-gp120與IL-2/IL-6共表達(dá)核酸疫苗質(zhì)粒的構(gòu)建和實(shí)驗(yàn)免疫研究
    国产精品98久久久久久宅男小说| 中国美女看黄片| 极品教师在线免费播放| 午夜精品在线福利| 亚洲avbb在线观看| 两性夫妻黄色片| 黄色毛片三级朝国网站| www国产在线视频色| 麻豆国产av国片精品| 亚洲精品久久国产高清桃花| 成人高潮视频无遮挡免费网站| 欧美国产日韩亚洲一区| 国产真实乱freesex| 中文字幕熟女人妻在线| 欧美中文综合在线视频| 嫁个100分男人电影在线观看| 亚洲av日韩精品久久久久久密| 亚洲18禁久久av| av中文乱码字幕在线| 熟女电影av网| 九色国产91popny在线| 欧美中文日本在线观看视频| 中文字幕精品亚洲无线码一区| 十八禁人妻一区二区| 久久中文字幕人妻熟女| 国产精品久久久人人做人人爽| 亚洲 欧美 日韩 在线 免费| 国产一级毛片七仙女欲春2| 国产探花在线观看一区二区| 国产精华一区二区三区| 亚洲av电影在线进入| 日本精品一区二区三区蜜桃| 国产片内射在线| 免费在线观看成人毛片| 亚洲人与动物交配视频| 国产片内射在线| 可以在线观看毛片的网站| 亚洲午夜精品一区,二区,三区| 老司机靠b影院| www日本黄色视频网| АⅤ资源中文在线天堂| 亚洲一码二码三码区别大吗| 99riav亚洲国产免费| 制服诱惑二区| 国产久久久一区二区三区| 免费一级毛片在线播放高清视频| 国产精品98久久久久久宅男小说| 亚洲成人精品中文字幕电影| 国产欧美日韩精品亚洲av| 欧美zozozo另类| 欧美日韩福利视频一区二区| 亚洲av成人精品一区久久| 中文资源天堂在线| 日韩免费av在线播放| 午夜a级毛片| 99在线人妻在线中文字幕| 亚洲av熟女| 亚洲精华国产精华精| 中文字幕精品亚洲无线码一区| 无人区码免费观看不卡| 国产成人av激情在线播放| 亚洲欧美日韩东京热| 精品一区二区三区视频在线观看免费| 97人妻精品一区二区三区麻豆| 亚洲精品国产精品久久久不卡| 老鸭窝网址在线观看| 国产真实乱freesex| 日本三级黄在线观看| 久久国产精品人妻蜜桃| 免费无遮挡裸体视频| ponron亚洲| 男女做爰动态图高潮gif福利片| 久久精品aⅴ一区二区三区四区| 久久久国产精品麻豆| 在线观看免费日韩欧美大片| 狂野欧美激情性xxxx| 看免费av毛片| 亚洲熟女毛片儿| 亚洲国产精品sss在线观看| 亚洲熟妇熟女久久| 一边摸一边抽搐一进一小说| 夜夜夜夜夜久久久久| 中文字幕久久专区| 一边摸一边做爽爽视频免费| 日韩精品免费视频一区二区三区| 黄频高清免费视频| 91国产中文字幕| 人妻久久中文字幕网| 天堂影院成人在线观看| 国产伦一二天堂av在线观看| 成人av一区二区三区在线看| 亚洲精品美女久久久久99蜜臀| 一夜夜www| 亚洲一卡2卡3卡4卡5卡精品中文| 淫妇啪啪啪对白视频| 国产精品国产高清国产av| 婷婷精品国产亚洲av在线| 欧美黑人欧美精品刺激| 嫁个100分男人电影在线观看| 精品欧美一区二区三区在线| 两性午夜刺激爽爽歪歪视频在线观看 | a级毛片a级免费在线| 人成视频在线观看免费观看| 欧美中文日本在线观看视频| 国产精品98久久久久久宅男小说| 婷婷丁香在线五月| 国产探花在线观看一区二区| 看免费av毛片| 国产伦在线观看视频一区| 国产成人精品久久二区二区免费| 9191精品国产免费久久| 黑人操中国人逼视频| 一进一出抽搐gif免费好疼| 欧美黑人欧美精品刺激| 97超级碰碰碰精品色视频在线观看| 久久精品国产综合久久久| 午夜影院日韩av| 国产亚洲精品久久久久5区| 国产v大片淫在线免费观看| 欧美激情久久久久久爽电影| 免费在线观看完整版高清| 香蕉国产在线看| 啦啦啦韩国在线观看视频| 午夜福利高清视频| 黄频高清免费视频| 最新美女视频免费是黄的| 欧美三级亚洲精品| 亚洲国产欧美人成| 校园春色视频在线观看| 久久久久久免费高清国产稀缺| 天堂动漫精品| 成年免费大片在线观看| 黄色成人免费大全| 床上黄色一级片| bbb黄色大片| 9191精品国产免费久久| 在线永久观看黄色视频| 操出白浆在线播放| 操出白浆在线播放| 很黄的视频免费| 日本在线视频免费播放| 一级毛片精品| 国内精品一区二区在线观看| 午夜福利高清视频| 久久中文看片网| 免费看美女性在线毛片视频| 亚洲av电影不卡..在线观看| 国产精品1区2区在线观看.| 国产蜜桃级精品一区二区三区| 久久性视频一级片| 国产激情欧美一区二区| 成人三级黄色视频| 午夜精品在线福利| 两性午夜刺激爽爽歪歪视频在线观看 | 精品久久久久久久毛片微露脸| e午夜精品久久久久久久| 天天一区二区日本电影三级| 亚洲av成人av| 久久草成人影院| av免费在线观看网站| 亚洲av美国av| av片东京热男人的天堂| 国产精品久久久久久久电影 | 亚洲美女黄片视频| 男人的好看免费观看在线视频 | 视频区欧美日本亚洲| 久久久国产成人精品二区| 国产免费男女视频| 亚洲美女视频黄频| 三级毛片av免费| 两性午夜刺激爽爽歪歪视频在线观看 | 午夜精品久久久久久毛片777| 一级毛片高清免费大全| 国产精品av视频在线免费观看| 俺也久久电影网| 欧美最黄视频在线播放免费| 久久精品人妻少妇| 岛国视频午夜一区免费看| 男插女下体视频免费在线播放| 久久 成人 亚洲| 国产亚洲精品第一综合不卡| a级毛片a级免费在线| 黄色片一级片一级黄色片| 男女那种视频在线观看| 久久香蕉精品热| 欧美日韩中文字幕国产精品一区二区三区| 国产激情久久老熟女| 一边摸一边做爽爽视频免费| 日韩欧美免费精品| 久久久精品大字幕| 午夜福利18| 母亲3免费完整高清在线观看| 九色国产91popny在线| 男人舔奶头视频| 男女下面进入的视频免费午夜| 丁香欧美五月| 女生性感内裤真人,穿戴方法视频| 国产高清有码在线观看视频 | 免费在线观看日本一区| 91在线观看av| 国产亚洲精品综合一区在线观看 | 精品国产美女av久久久久小说| 久久精品影院6| 久久精品国产亚洲av高清一级| 亚洲欧美激情综合另类| 天堂影院成人在线观看| 免费无遮挡裸体视频| avwww免费| 国产一区二区在线av高清观看| 久久亚洲精品不卡| 午夜a级毛片| 国产高清激情床上av| 琪琪午夜伦伦电影理论片6080| 亚洲国产欧美人成| 国产伦人伦偷精品视频| 久久久久国内视频| 免费在线观看黄色视频的| 亚洲五月婷婷丁香| 欧美日韩一级在线毛片| 亚洲熟妇熟女久久| 香蕉丝袜av| 成年人黄色毛片网站| bbb黄色大片| 色在线成人网| 免费电影在线观看免费观看| 国产成人精品无人区| 最近在线观看免费完整版| 亚洲一区高清亚洲精品| av片东京热男人的天堂| 精品国内亚洲2022精品成人| 中文亚洲av片在线观看爽| 日本 av在线| a级毛片a级免费在线| 91国产中文字幕| 久久精品综合一区二区三区| 亚洲av美国av| 丝袜美腿诱惑在线| 久久久国产成人免费| 免费在线观看亚洲国产| 法律面前人人平等表现在哪些方面| 亚洲 欧美一区二区三区| 国内精品一区二区在线观看| 999精品在线视频| 国产一区在线观看成人免费| 色噜噜av男人的天堂激情| 亚洲av日韩精品久久久久久密| 最近在线观看免费完整版| 中文字幕最新亚洲高清| av免费在线观看网站| 无遮挡黄片免费观看| 欧美性猛交╳xxx乱大交人| 亚洲人成电影免费在线| 欧美一级a爱片免费观看看 | 岛国视频午夜一区免费看| 黄片大片在线免费观看| 日本 欧美在线| 免费搜索国产男女视频| 久久性视频一级片| 中文字幕久久专区| 成人永久免费在线观看视频| e午夜精品久久久久久久| 色尼玛亚洲综合影院| 日本熟妇午夜| 国产三级在线视频| 亚洲性夜色夜夜综合| 三级毛片av免费| 夜夜躁狠狠躁天天躁| 国产精品一区二区精品视频观看| 日韩精品青青久久久久久| 欧美一级a爱片免费观看看 | 99久久精品热视频| 最近视频中文字幕2019在线8| 1024视频免费在线观看| 欧美日本视频| 88av欧美| 一本精品99久久精品77| 久久国产乱子伦精品免费另类| 成年人黄色毛片网站| 欧美日韩瑟瑟在线播放| 欧美日韩亚洲综合一区二区三区_| 级片在线观看| 亚洲欧美日韩无卡精品| 中文字幕精品亚洲无线码一区| 亚洲中文字幕日韩| 1024视频免费在线观看| 手机成人av网站| 国产精品日韩av在线免费观看| 成人18禁在线播放| 精华霜和精华液先用哪个| 在线观看免费午夜福利视频| 亚洲一区二区三区色噜噜| 狠狠狠狠99中文字幕| 在线国产一区二区在线| 免费在线观看完整版高清| 成熟少妇高潮喷水视频| 天天躁狠狠躁夜夜躁狠狠躁| 久9热在线精品视频| 黄色成人免费大全| 母亲3免费完整高清在线观看| 好看av亚洲va欧美ⅴa在| 亚洲精品美女久久久久99蜜臀| 亚洲国产精品久久男人天堂| 村上凉子中文字幕在线| 久9热在线精品视频| 小说图片视频综合网站| 久久久精品欧美日韩精品| 亚洲av成人精品一区久久| 亚洲美女视频黄频| 亚洲国产中文字幕在线视频| 美女免费视频网站| 狂野欧美激情性xxxx| 久久久久久九九精品二区国产 | 亚洲一区二区三区色噜噜| 国产精品 国内视频| 国产午夜精品论理片| 亚洲男人天堂网一区| 国产成人一区二区三区免费视频网站| 欧美一级a爱片免费观看看 | 桃色一区二区三区在线观看| av福利片在线观看| 国产精品日韩av在线免费观看| 性欧美人与动物交配| 欧美绝顶高潮抽搐喷水| 国产一区二区在线观看日韩 | 亚洲一区二区三区不卡视频| 啦啦啦免费观看视频1| 国产成人精品无人区| 亚洲一区高清亚洲精品| 精品国产亚洲在线| 亚洲片人在线观看| 香蕉丝袜av| 国产成人av激情在线播放| 老熟妇仑乱视频hdxx| 欧美在线一区亚洲| 怎么达到女性高潮| 午夜免费成人在线视频| 少妇被粗大的猛进出69影院| 成人手机av| 亚洲国产欧美网| 久久久久久人人人人人| 老汉色av国产亚洲站长工具| 国产成人精品久久二区二区91| 国产99久久九九免费精品| avwww免费| 给我免费播放毛片高清在线观看| 天堂影院成人在线观看| 色在线成人网| 亚洲av熟女| 女人被狂操c到高潮| 小说图片视频综合网站| 午夜亚洲福利在线播放| 天堂av国产一区二区熟女人妻 | 色噜噜av男人的天堂激情| 99riav亚洲国产免费| 中文字幕av在线有码专区| 中文字幕人成人乱码亚洲影| 色精品久久人妻99蜜桃| 国产精品一区二区免费欧美| 国产一区二区三区视频了| 久久中文字幕一级| 国产又黄又爽又无遮挡在线| 又粗又爽又猛毛片免费看| 人成视频在线观看免费观看| 中文字幕久久专区| 人妻夜夜爽99麻豆av| www.999成人在线观看| 精品欧美国产一区二区三| 激情在线观看视频在线高清| 19禁男女啪啪无遮挡网站| 亚洲天堂国产精品一区在线| 亚洲国产精品sss在线观看| 日本 欧美在线| 国产探花在线观看一区二区| 婷婷精品国产亚洲av在线| 精品久久久久久久末码| 亚洲 欧美一区二区三区| 波多野结衣高清作品| 在线播放国产精品三级| 亚洲成人国产一区在线观看| 两个人看的免费小视频| 免费在线观看完整版高清| 成人18禁高潮啪啪吃奶动态图| 久久欧美精品欧美久久欧美| 最近视频中文字幕2019在线8| 中出人妻视频一区二区| 手机成人av网站| 亚洲精品美女久久久久99蜜臀| 国产v大片淫在线免费观看| 国产视频一区二区在线看| 中亚洲国语对白在线视频| 亚洲专区国产一区二区| www.www免费av| 日本成人三级电影网站| 成人手机av| 精品一区二区三区视频在线观看免费| 国产熟女午夜一区二区三区| 老司机靠b影院| 别揉我奶头~嗯~啊~动态视频| 三级国产精品欧美在线观看 | 伦理电影免费视频| 宅男免费午夜| 高潮久久久久久久久久久不卡| 欧美高清成人免费视频www| 麻豆国产97在线/欧美 | 国产成年人精品一区二区| 91老司机精品| 亚洲中文av在线| 成人精品一区二区免费| 黄色视频,在线免费观看| 黑人操中国人逼视频| 国产精品av久久久久免费| 一本大道久久a久久精品| 亚洲精品中文字幕一二三四区| 欧美日韩瑟瑟在线播放| 国模一区二区三区四区视频 | 国产一区在线观看成人免费| √禁漫天堂资源中文www| 国产精品久久久人人做人人爽| 精品久久久久久成人av| 国产私拍福利视频在线观看| 久久久精品大字幕| 欧美性猛交黑人性爽| 日韩欧美在线乱码| 舔av片在线| 国产亚洲精品av在线| 少妇裸体淫交视频免费看高清 | 亚洲中文字幕一区二区三区有码在线看 | 精品国产超薄肉色丝袜足j| 亚洲精品久久国产高清桃花| 国语自产精品视频在线第100页| 搡老妇女老女人老熟妇| 亚洲精品一卡2卡三卡4卡5卡| 丝袜人妻中文字幕| 老汉色∧v一级毛片| 1024手机看黄色片| 日韩欧美国产一区二区入口| 中出人妻视频一区二区| 他把我摸到了高潮在线观看| 中文字幕久久专区| 舔av片在线| 亚洲欧美日韩东京热| 午夜福利在线在线| 日韩大码丰满熟妇| 久久久久久人人人人人| 无人区码免费观看不卡| 成人三级黄色视频| 国产高清videossex| 草草在线视频免费看| 国产精品 欧美亚洲| 国产主播在线观看一区二区| 精品一区二区三区av网在线观看| 琪琪午夜伦伦电影理论片6080| 中文资源天堂在线| 夜夜爽天天搞| 亚洲精品中文字幕在线视频| 欧美人与性动交α欧美精品济南到| 国产精品久久久久久久电影 | 香蕉丝袜av| 曰老女人黄片| 成年人黄色毛片网站| 91字幕亚洲| 国产区一区二久久| 一区福利在线观看| 精品第一国产精品| 色综合站精品国产| 色综合婷婷激情| 亚洲精品中文字幕一二三四区| av福利片在线观看| 欧美日韩一级在线毛片| 久久国产乱子伦精品免费另类| 国产精品久久电影中文字幕| 国产高清视频在线观看网站| 日韩欧美在线乱码| 亚洲黑人精品在线| 亚洲成av人片在线播放无| 狠狠狠狠99中文字幕| 麻豆成人av在线观看| 亚洲一卡2卡3卡4卡5卡精品中文| 久久国产乱子伦精品免费另类| 日本三级黄在线观看| 亚洲成人精品中文字幕电影| 亚洲国产中文字幕在线视频| 免费高清视频大片| 啦啦啦韩国在线观看视频| 最近视频中文字幕2019在线8| 又黄又爽又免费观看的视频| 在线观看免费日韩欧美大片| 国产精品影院久久| 午夜成年电影在线免费观看| 亚洲精品美女久久久久99蜜臀| 黑人操中国人逼视频| 日韩av在线大香蕉| 在线观看免费视频日本深夜| 狂野欧美白嫩少妇大欣赏| 一级a爱片免费观看的视频| 高清在线国产一区| 国产精品日韩av在线免费观看| 99国产精品一区二区三区| 亚洲精华国产精华精| 超碰成人久久| 香蕉丝袜av| 99热这里只有精品一区 | 美女午夜性视频免费| 非洲黑人性xxxx精品又粗又长| 国产爱豆传媒在线观看 | www.自偷自拍.com| 午夜福利视频1000在线观看| 身体一侧抽搐| 后天国语完整版免费观看| 精品一区二区三区视频在线观看免费| 麻豆av在线久日| 天天一区二区日本电影三级| 99国产精品一区二区蜜桃av| 老司机在亚洲福利影院| 这个男人来自地球电影免费观看| 后天国语完整版免费观看| 老司机靠b影院| 精品一区二区三区四区五区乱码| 一区二区三区高清视频在线| 1024手机看黄色片| 国产精品久久久久久亚洲av鲁大| 99久久精品国产亚洲精品| 亚洲熟女毛片儿| 好男人在线观看高清免费视频| 欧美一区二区精品小视频在线| 国产真实乱freesex| 三级男女做爰猛烈吃奶摸视频| 久久天躁狠狠躁夜夜2o2o| 小说图片视频综合网站| 国产精华一区二区三区| 老汉色∧v一级毛片| 麻豆国产av国片精品| 男人舔奶头视频| 一级作爱视频免费观看| 午夜福利在线在线| 久久久国产成人精品二区| 人人妻人人澡欧美一区二区| 国产97色在线日韩免费| 婷婷精品国产亚洲av在线| 日本一二三区视频观看| 亚洲avbb在线观看| 一a级毛片在线观看| 制服人妻中文乱码| 欧美黄色片欧美黄色片| 国产99久久九九免费精品| 一个人免费在线观看的高清视频| 亚洲中文日韩欧美视频| 黄色视频,在线免费观看| 在线十欧美十亚洲十日本专区| 色播亚洲综合网| 欧美激情久久久久久爽电影| 免费人成视频x8x8入口观看| 欧美日韩瑟瑟在线播放| 亚洲五月天丁香| 在线观看66精品国产| 一级片免费观看大全| 成人18禁在线播放| 成人特级黄色片久久久久久久| 丝袜人妻中文字幕| 精品国产超薄肉色丝袜足j| 久久性视频一级片| 麻豆久久精品国产亚洲av| 欧美日韩一级在线毛片| 嫩草影视91久久| 久久久久久亚洲精品国产蜜桃av| 成人永久免费在线观看视频| 欧美一区二区国产精品久久精品 | 欧美极品一区二区三区四区| 欧美黑人精品巨大| 久久午夜综合久久蜜桃| 亚洲精品美女久久久久99蜜臀| 免费观看精品视频网站| 90打野战视频偷拍视频| 特级一级黄色大片| 国产精品精品国产色婷婷| 91在线观看av| 久久中文字幕人妻熟女| 精品久久蜜臀av无| 成人国产综合亚洲| 一个人观看的视频www高清免费观看 | 中文字幕人妻丝袜一区二区| 精品久久久久久久人妻蜜臀av| 久久中文字幕人妻熟女| 国产精品99久久99久久久不卡| 国产av又大| 男人舔女人的私密视频| 91麻豆av在线| 国产亚洲精品第一综合不卡| 日本在线视频免费播放| www.999成人在线观看| 亚洲国产日韩欧美精品在线观看 | 亚洲av成人不卡在线观看播放网| 精品无人区乱码1区二区| 一级a爱片免费观看的视频| 丰满人妻熟妇乱又伦精品不卡| 亚洲精品美女久久av网站| 9191精品国产免费久久| 大型av网站在线播放| 精品国产超薄肉色丝袜足j| 久久这里只有精品19| 男女午夜视频在线观看| 国产精品永久免费网站| 国产精品九九99| 99国产精品一区二区三区| 一本久久中文字幕| 97碰自拍视频| 成年免费大片在线观看| 丰满的人妻完整版| 波多野结衣巨乳人妻| 免费人成视频x8x8入口观看| 999精品在线视频| 热99re8久久精品国产| 午夜福利欧美成人| 91老司机精品| av在线天堂中文字幕|