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

    基于高通量測序的慈竹筍轉(zhuǎn)錄組分析與基因功能注釋

    2016-12-12 05:12:36陳宇鵬曹穎胡尚連黃艷盧學琴徐剛龍治堅
    生物工程學報 2016年11期
    關(guān)鍵詞:慈竹毛竹竹筍

    陳宇鵬,曹穎,胡尚連,黃艷,盧學琴,徐剛,龍治堅

    西南科技大學 生命科學與工程學院 四川省生物質(zhì)資源利用與改性工程技術(shù)研究中心,四川 綿陽 621010

    基于高通量測序的慈竹筍轉(zhuǎn)錄組分析與基因功能注釋

    陳宇鵬,曹穎,胡尚連,黃艷,盧學琴,徐剛,龍治堅

    西南科技大學 生命科學與工程學院 四川省生物質(zhì)資源利用與改性工程技術(shù)研究中心,四川 綿陽 621010

    陳宇鵬, 曹穎, 胡尚連, 等. 基于高通量測序的慈竹筍轉(zhuǎn)錄組分析與基因功能注釋. 生物工程學報, 2016, 32(11): 1610-1623.

    Chen YP, Cao Y, Hu SL, et al. Transcriptome analysis and gene function annotation of Bambusa emeiensis shoots based on high-throughput sequencing technology. Chin J Biotech, 2016, 32(11): 1610-1623.

    慈竹是我國四川當?shù)氐膬?yōu)勢叢生竹種之一,其纖維長度和質(zhì)量較優(yōu)異,是造紙、紡織等工業(yè)的良好原料。本文利用Illumina HiSeqTM2000平臺,對10、50、100和150 cm高的慈竹筍進行轉(zhuǎn)錄組分析,共得到69.28 M條讀長 (Reads),經(jīng)從頭拼接、組裝和聚類后得到111 137條非重復(fù)序列基因Unigene,其中共有63 094條注釋到COG、GO、KEGG、Swiss-Prot和Nr數(shù)據(jù)庫中。這些Unigene不僅具有一般的功能,如轉(zhuǎn)錄和信號轉(zhuǎn)導等,還涉及到蔗糖轉(zhuǎn)運與代謝、次級代謝產(chǎn)物及細胞壁的生物合成等方面。不同高度慈竹筍的纖維素合成酶基因存在差異表達,發(fā)現(xiàn)了可能調(diào)控慈竹生長發(fā)育以及纖維素和木質(zhì)素生物合成的相關(guān)基因,為慈竹品種改良提供一定的理論基礎(chǔ)。

    慈竹,轉(zhuǎn)錄組,高通量測序,功能注釋

    竹子是重要的非木材可再生資源,以生長速度快、產(chǎn)量高、易加工、用途廣泛、經(jīng)濟價值高等特點著稱。近年來隨著分子生物學研究手段的發(fā)展,尤其是基因組和轉(zhuǎn)錄組測序技術(shù)的廣泛應(yīng)用,毛竹、麻竹等竹種的基因組和 (或)轉(zhuǎn)錄組數(shù)據(jù)已經(jīng)被公布,一些與竹筍快速生長以及竹材纖維形成的相關(guān)基因被報道。如Peng等[1-2]通過基因組和轉(zhuǎn)錄組測序發(fā)現(xiàn)了19個毛竹纖維素合酶基因,和大量與植物激素、細胞周期調(diào)控、細胞壁合成、細胞形態(tài)等相關(guān)的候選基因,為毛竹竹材特性的改良提供了依據(jù)。Liu等[3]分析了麻竹的轉(zhuǎn)錄組,篩選出105個木質(zhì)素生物合成途徑相關(guān)關(guān)鍵酶基因,并與毛竹的相關(guān)基因進行了比對,確定了潛在的與生長和發(fā)育相關(guān)的候選基因。

    慈竹Bambusa emeiensis,是四川本土大型叢生竹之一,其稈徑較小,稈壁較薄,纖維長度和質(zhì)量均優(yōu)于其他竹種,是造紙、紡織等行業(yè)的優(yōu)良原料之一,具有比較重要的應(yīng)用前景[4]。目前,對于慈竹的研究主要集中于生長狀況[5]、基因克隆[6]、遺傳多樣性[7]、生理生化指標[8]等方面。與毛竹、麻竹等竹種相比,慈竹的分子生物學研究相對滯后。有研究表明慈竹可能有70條染色體,為典型的多倍體復(fù)合體[9],其基因組測序的難度較大,而且目前慈竹轉(zhuǎn)錄組分析以及纖維素、木質(zhì)素生物合成相關(guān)的分子機理等也未見報道。

    本文以不同高度的慈竹筍為材料,比較了不同高度慈竹筍的纖維素含量,并采用Illumina HiSeqTM2000平臺對其進行了轉(zhuǎn)錄組測序。將得到的Unigene進行Cluster of Orthologous Groups of proteins (COG)、Gene Ontology (GO)、Kyoto Encyclopedia of Genes and Genomes (KEGG)、Swiss-Prot和Nonredundant protein (Nr)注釋,并與其他相似物種進行比對,旨在發(fā)現(xiàn)調(diào)控慈竹生長發(fā)育以及與纖維素、木質(zhì)素生物合成相關(guān)基因,挖掘慈竹中的珍貴基因資源。同時,本文也進一步分析了4個樣品中與纖維素合成相關(guān)的差異表達基因,為慈竹的分子遺傳育種以及品種改良等提供一定的理論基礎(chǔ)。

    1 材料與方法

    1.1 材料

    供試材料采自西南科技大學生命科學與工程學院資源圃 (年平均氣溫17.2 ℃,年平均降雨量793.5 mm)。采樣時,以露出地面10、50、100和150 cm的筍為對象,每個高度取3株獨立的筍 (測量誤差小于2 cm),并取其基部第2個節(jié)間的部分,分別切碎、混勻,其中一部分儲存于-80 ℃冰箱中備用,另一部分置于烘箱中烘干用于纖維素含量測定。

    1.2 不同高度慈竹筍纖維素含量測定

    纖維素含量的測定參照FOSS公司FibertecTMM6 1020/1021型纖維素測定儀的操作手冊進行。通過酸性洗滌纖維 (ADF)[10]法,得到纖維素含量,每個樣品生物重復(fù)3次,每組樣品運行1個空白對照。

    1.3 轉(zhuǎn)錄組測序

    1.3.1 總RNA瓊脂糖凝膠電泳檢測及測序

    RNA提取采用OMEGA BIO-TEK公司Plant RNA Kit試劑盒中的試劑及提取方法??俁NA樣品采用Nanodrop檢測和Agilent Technologies 2100 Bioanalyzer檢測。檢測合格后用NEBNext Poly (A) mRNA Magnetic Isolation Module (NEB, E7490) 富集mRNA,將mRNA片段化處理,以mRNA片段為模板,采用隨機引物法,用NEBNext mRNA Library Prep Master Mix Set for Illumina (NEB, E6110) 和NEBNext Multiplex Oligos for Illumina (NEB, E7500) 構(gòu)建上機文庫。制備好的文庫用1.8%瓊脂糖凝膠電泳進行檢測文庫插入片段大小,然后用Library Quantification Kit-Illumina GA Universal (Kapa, KK4824) 進行QPCR定量,再在cDNA片段兩端加上接頭。檢測合格的文庫在Illumina cbot上進行簇的生成,最后用Illumina HiSeqTM2000進行測序。

    1.3.2 數(shù)據(jù)分析

    轉(zhuǎn)錄組測序得到的原始數(shù)據(jù)經(jīng)過去除雜質(zhì)和冗余處理后,利用Trinity軟件[11]對經(jīng)過過濾后的高質(zhì)量數(shù)據(jù)進行從頭拼接,獲得重疊群Contig。根據(jù)其結(jié)果利用雙末端信息將來自同一轉(zhuǎn)錄本的不同Contig連在一起,做進一步的序列拼接,得到轉(zhuǎn)錄本Transcripts。在Transcripts聚類單元中選取最主要的轉(zhuǎn)錄本作為Unigene序列,并對Unigene數(shù)據(jù)進行聚類分析和進一步去冗余處理,最終得到非冗余Unigene庫。對于多樣品的組裝,由于后續(xù)的Open reading frame (ORF) 預(yù)測、Simple sequence repeats (SSR) 分析、表達豐度分析等都建立在同一套參考基因的基礎(chǔ)上,因此對各樣品得到的Unigene作進一步的聚類分析,最后得到慈竹的Unigene數(shù)據(jù)庫。

    1.3.3 功能注釋

    使用Basic Local Alignment Search Tool (Blast) 將組裝后的Unigene序列與COG、GO、KEGG、Swiss-Prot和Nr數(shù)據(jù)庫比對,得到最高序列相似性的Unigene。

    1.3.4 基因的差異表達分析

    在NR、Swiss-Rort注釋結(jié)果中進行檢索,記錄“GeneID”以及相對應(yīng)的注釋結(jié)果。根據(jù)“GeneID”記錄相應(yīng)基因的相對表達豐度值,并將值導入MeV 4.9.0分析軟件中,制作熱圖并進行聚類分析。

    1.4 差異表達基因qRT-PCR驗證

    利用primer premier 5.0對已測序正確的序列和慈竹內(nèi)參基因Tublin進行qRT-PCR引物的設(shè)計。在iQ Multicolor Real-Time PCR自動擴增儀進行反應(yīng),每個樣品重復(fù)3次,反應(yīng)體系20 μL。以內(nèi)參基因Tublin為對照,利用公式2-CT△△計算其相對表達量。

    表1 測序數(shù)據(jù)統(tǒng)計Table 1 Statistics of sequencing data

    圖1 慈竹功能基因長度分布和所占比例Fig. 1 Distribution and percentage of Unigene length from Bambusa emeiensis.

    2 結(jié)果與分析

    2.1 測序數(shù)據(jù)結(jié)果統(tǒng)計

    利用Illumina HiSeqTM2000平臺,對供試材料進行轉(zhuǎn)錄組測序,共獲得13.98 Gb的原始數(shù)據(jù),得到69.28 M條的讀長 (reads) (表1)。其中10、50和150 cm慈竹筍總堿基數(shù)均大于4 Gb,reads數(shù)大于20 M,GC含量分別為52.78%、53.12%和52.21%;而100 cm慈竹筍因測序深度較低,僅獲得7 829 716條reads,GC含量為49.49%。同時,供試材料的不確定堿基數(shù)(N%) 所占比例均≤0.03%,測序質(zhì)量值≥30的堿基(Q 30%) 所占比例均>80%,這表明測序質(zhì)量可靠。

    通過進一步分析后得到的非冗余Unigene庫,含有111 137條慈竹筍的Unigene。其中200 bp-300 bp之間的Unigene,為41 814條,占37.62%;大于1 000 bp的Unigene為21 518條,占19.35% (圖1A,1B)。

    2.2 基因功能注釋

    2.2.1 基因功能注釋

    通過數(shù)據(jù)庫比對,得到最高序列相似性的Unigene,其中有16 151條Unigene被匹配到COG數(shù)據(jù)庫中;52 872條Unigene被匹配到GO數(shù)據(jù)庫中;11 032條Unigene被匹配到KEGG數(shù)據(jù)庫中;44 552條Unigene被匹配到Swiss-Port數(shù)據(jù)庫中;而62 795條Unigene被匹配到Nr數(shù)據(jù)庫中 (表2)。

    表2 慈竹轉(zhuǎn)錄組功能注釋Table 2 Transcriptome function annotation of Bambusa emeiensis

    圖2 Unigene的COG功能分類Fig. 2 Unigene of COG function categories.

    2.2.2 功能基因的COG注釋

    COG的功能注釋與基因產(chǎn)物有直接的聯(lián)系[12]。在慈竹筍Unigene數(shù)據(jù)庫的COG功能注釋中,有16 151條Unigene具有蛋白質(zhì)功能定義。在COG功能分類中包含復(fù)制、重組、修復(fù)、信號轉(zhuǎn)導、新陳代謝、分子過程等25個分類。其中一般功能基因代表最大的一類,有4 221條Unigene;其次是復(fù)制、重組和修復(fù),有2 701條。此外,還涉及到多個植物生長發(fā)育相關(guān)以及纖維素、木質(zhì)素生物合成相關(guān)的生理生化過程,如碳水化合物轉(zhuǎn)運與代謝、氨基酸轉(zhuǎn)運與代謝、輔酶轉(zhuǎn)運與代謝、次級代謝產(chǎn)物生物合成、運輸和分解代謝、細胞壁與細胞膜的生物合成、信號轉(zhuǎn)導機制等過程 (圖2)。

    2.2.3 功能基因的GO分類

    GO數(shù)據(jù)庫可對基因和蛋白質(zhì)進行分類和注釋,其對很多物種都是適用的[13]。GO數(shù)據(jù)庫包含3大類,即細胞組分、分子功能和生物學過程。在慈竹筍Unigene數(shù)據(jù)庫中,有52 872條基因注釋到GO數(shù)據(jù)庫中。其中第一大類中的細胞部分、細胞和細胞器所占比例較高,分別達到85.71%、84.85%和80.39%;第二大類中的結(jié)合和催化活性所占比例較高,分別達到68.90%和50.82%;而第三大類中的細胞過程和代謝過程所占比例較高,分別達到75.86%和73.31% (圖3)。

    2.2.4 功能基因的KEGG注釋

    KEGG數(shù)據(jù)庫是有關(guān)基因產(chǎn)物代謝通路的主要數(shù)據(jù)庫,基于代謝通路的分析對于進一步分析解讀基因的功能非常有幫助[14]。在慈竹筍Unigene數(shù)據(jù)庫中,有11 032條Unigene注釋到KEGG數(shù)據(jù)庫中。其中注釋序列較多的3條代謝通路為RNA轉(zhuǎn)運、植物激素信號轉(zhuǎn)導與核糖體。而與植物纖維素生物合成相關(guān)的代謝通路,如淀粉和蔗糖代謝、光合生物的固碳作用和光合作用分別有257條、189條和105條注釋到KEGG數(shù)據(jù)庫中;與木質(zhì)素生物合成相關(guān)的代謝途徑,如苯丙烷類的生物合成和苯丙氨酸代謝分別有118條和110條在KEGG代謝通路中得到注釋 (表3)。

    圖3 Unigene的GO功能分類Fig. 3 Unigene of GO function categories.

    表3 Unigene數(shù)目最多的5個代謝通路及與纖維素、木質(zhì)素生物合成相關(guān)的代謝通路Table 3 Top five metabolic pathways and the metabolic pathways of cellulose and lignin biosynthesis involving Bambusa emeiensis Unigene

    2.3 慈竹Unigene與其他物種的比較分析

    首先將慈竹的Unigene與水稻Oryza sativa、毛竹Phyllostachys edulis、短柄草Brachypodium、玉米Zea mays和高粱Sorghum的Unigene進行功能比對。所得結(jié)果如表4所示,慈竹Unigene與毛竹匹配度最高,達到52.46%,可能是因為慈竹與毛竹同屬竹類植物的原因。其次是與水稻相匹配達到50.41%,可能是因為水稻為模式植物,其基因序列信息較其他幾種禾本科更全面。

    表4 慈竹Unigene與其他物種的比較分析Table 4 Analysis on Unigene between Bambusa emeiensis and other species

    根據(jù)表4的結(jié)果,再將慈竹的Unigene與水稻轉(zhuǎn)錄本序列和毛竹編碼序列數(shù)據(jù)庫進行Blast比對,選擇E-value<1e-10。63 094個慈竹Unigene中34 127個 (54.09%) Unigene與毛竹轉(zhuǎn)錄組相匹配,33 685個 (53.39%) Unigene與水稻相匹配。在這些配對的序列中29 591個(46.90%) Unigene與毛竹和水稻都相匹配,4 536個 (7.19%) Unigene只與水稻相匹配,4 094個(6.49%) Unigene只與毛竹相匹配。還有許多預(yù)測的Unigene與水稻和毛竹都不匹配,它們編碼的蛋白質(zhì)主要與生理過程、結(jié)合活性和催化功能有關(guān),其中有大量預(yù)測的Unigene可能是慈竹特有的 (圖4)。

    圖4 慈竹Unigene與水稻和毛竹的比較分析Fig. 4 Analysis on Unigene between Bambusa emeiensis, Oryza sativa and Phyllostachys edulis.

    表5 慈竹纖維素合成的功能基因與水稻、短柄草、玉米、高粱和毛竹的比較Table 5 Bambusa emeiensis functional Unigene of cellulose biosynthesis comparison with Oryza sativa, Brachypodium, Zea mays, Sorghum and Phyllostachys edulis

    表6 慈竹木質(zhì)素合成的功能基因與水稻、短柄草、玉米、高粱和毛竹的比較Table 6 Bambusa emeiensis functional Unigene of lignin biosynthesis comparison with Oryza sativa, Brachypodium, Zea mays and Sorghum and Phyllostachys edulis

    2.4 纖維素、木質(zhì)素生物合成相關(guān)基因與其他物種比對分析

    2.4.1 纖維素生物合成相關(guān)酶基因與其他物種比對分析

    將鑒定出的慈竹纖維素生物合成途徑中的關(guān)鍵酶基因與水稻、短柄草、玉米、高粱和毛竹進行比對分析,毛竹的相關(guān)數(shù)據(jù)引用自Peng等[1]。由表5可知,在Swiss-Prot數(shù)據(jù)庫中,慈竹纖維素合酶Cellulose synthase (CesA) 的數(shù)量與其他幾個物種相比差異較大,CesA僅有9條。纖維素合酶相似蛋白Cellulose synthase-like protein (Csl) 則是毛竹的最多有76條,慈竹僅有11條,與毛竹相比差異比較明顯。而短柄草中最少,只在Swiss-Prot數(shù)據(jù)庫中鑒定到7條Csl基因。蔗糖合酶Sucrose synthase (SuSy) 數(shù)量最多的是水稻 (21),其次是短柄草 (11),慈竹只在Swiss-Prot數(shù)據(jù)庫中鑒定出7條。

    2.4.2 木質(zhì)素生物合成相關(guān)酶基因與其他物種比對分析

    經(jīng)Swiss-Prot數(shù)據(jù)庫功能注釋,篩選出188個木質(zhì)素生物合成途徑相關(guān)關(guān)鍵酶基因。與其他禾本科植物相比,慈竹苯丙氨酸解氨酶Phenylalanine ammonia-lyase (PAL)、4-香豆酸輔酶A連接酶4-coumarate-CoA ligase (4CL)、肉桂酰輔酶A還原酶Cinnamoyl-CoA reductase (CCR)、肉桂醇脫氫酶Cinnamoyl alcohol dehydrogenase (CAD) 等7個關(guān)鍵酶基因均略少于水稻 (表6)。在Swiss-Prot數(shù)據(jù)庫中,慈竹、短柄草、玉米和高粱的漆酶Laccase (LAC) 基因是數(shù)量最豐富的,分別有17、32、20和20條,毛竹則未檢測到。水稻的CCR,咖啡酰輔酶A-O-甲基轉(zhuǎn)移酶Caffeoyl-CoA O-methyltransferase (CCoAOMT) 和咖啡酸/5-羥基阿魏酸-O-甲基轉(zhuǎn)移酶Caffeic acid-O-methyltransferase (COMT) 基因數(shù)量較其他物種豐富,分別為45條、10條和10條。在慈竹中數(shù)量較少的是COMT,為2條。在毛竹中PAL最多,為8條,其余則較少,與慈竹相比差異較大。毛竹的相關(guān)數(shù)據(jù)引用自Peng等[1]。

    2.4.3 轉(zhuǎn)錄因子相關(guān)基因與其他物種比對分析

    轉(zhuǎn)錄因子是由核基因編碼的一類蛋白質(zhì),能選擇性地激活或抑制某些基因的轉(zhuǎn)錄與表達[15]。從Swiss-Prot數(shù)據(jù)庫中得到的轉(zhuǎn)錄因子如表7所示。慈竹MYB和bHLH轉(zhuǎn)錄因子數(shù)量較水稻、短柄草、玉米和高粱多,且轉(zhuǎn)錄因子整體數(shù)量上與毛竹相比差異較明顯。

    2.5 差異表達基因分析

    2.5.1 差異表達基因的篩選與統(tǒng)計

    為了解4個不同高度慈竹筍的差異表達基因,采用DESeq軟件,設(shè)定FDR<0.01、Fold Change (FC) ≥2作為篩選差異表達基因的關(guān)鍵指標。將10、50、100和150 cm四個樣品進行兩兩比對后,得到以下結(jié)果,10 vs 50 cm組有1 828條差異表達基因;10 vs 100 cm組有1 542條差異表達基因;10 vs 150 cm組有2 601條差異表達基因;50 vs 100 cm組有1 224條差異表達基因;50 vs 150 cm組有1 745條差異表達基因;100 vs 150 cm組有1 231條差異表達基因。將同一組合在不同的數(shù)據(jù)庫間進行比較,提取并統(tǒng)計篩選得到的差異表達基因的注釋信息(表8)。

    表7 慈竹轉(zhuǎn)錄因子相關(guān)基因與水稻、短柄草、玉米、高粱和毛竹的比較Table 7 Transcription factor of Bambusa emeiensis comparison with Oryza sativa, Brachypodium, Zea mays and Sorghum and Phyllostachys edulis

    表8 注釋的差異表達基因數(shù)目統(tǒng)計Table 8 Summary of differential expression gene annotated

    2.5.2 纖維素生物合成相關(guān)差異表達基因分析

    慈竹纖維素含量直接影響著造紙效率。通過酸性洗滌纖維 (ADF) 法分析了4個不同高度的慈竹筍基部的纖維素含量。結(jié)果表明,在筍的早期生長階段,隨著筍的伸長,慈竹筍基部的纖維素含量逐漸增加 (圖5)。通過進一步的轉(zhuǎn)錄組數(shù)據(jù)分析,可以從分子機理上了解慈竹筍在發(fā)育過程中纖維素含量的變化。

    圖5 不同高度慈竹筍纖維素含量Fig. 5 The content of cellulose in different height of Bambusa emeiensis.

    CesA在纖維素生物合成過程中起著關(guān)鍵的調(diào)節(jié)作用。轉(zhuǎn)錄組數(shù)據(jù)分析表明,CesA基因的RPKM值在10和50 cm的樣品中聚在一簇;而100和150 cm的樣品聚為另一簇,表明CesA基因隨慈竹筍高度的增長呈現(xiàn)不同的表達模式(圖6A,6B)。其中,CesA1、CesA5、CesA6和CesA8同源的Unigene有著較高的相對表達量,但卻隨著筍的伸長,呈現(xiàn)出降低的趨勢(圖6A,6B)。在4個不同生長高度的慈竹筍中,與CesA4、CesA7、CesA9同源的10條Unigene,雖然呈現(xiàn)著較低的表達,卻隨著筍的伸長,呈現(xiàn)出增高的趨勢 (圖6A,6B)。這些Unigene在GO數(shù)據(jù)庫的功能注釋中都涉及多個GO功能條目,包括細胞壁、初生細胞壁生物合成、次生細胞壁生物合成、纖維素合酶活性、纖維素生物合成過程、木質(zhì)部發(fā)育、細胞壁增厚等,而在KEGG注釋中,10條Unigene都參與的KEGG通路是K10999 (纖維素合酶) (表9)。

    圖6 不同生長高度的慈竹筍中CesA的相對表達水平 (RPKM值) 和差異表達水平 (Log2值)Fig. 6 The relative expression level (RPKM value) and the differential expression level (Log2 value) of CesA in the shoots of different height from Bambusa emeiensis.

    表9 纖維素合酶基因的功能注釋Table 9 Annotation of cellulose synthase genes

    2.6 差異表達基因qRT-PCR驗證

    根據(jù)纖維素生物合成相關(guān)差異表達基因熱圖分析,篩選出CesA9 (T1_Unigene_BMK.22002)進行qRT-PCR驗證。如圖7所示,在100 cm的慈竹筍中,CesA9的相對表達量大約是10 cm筍的9.5倍,表明其表達水平可能直接影響著筍的纖維素含量。

    圖7 CesA9 (T1_Unigene_BMK.22002) 基因qRT-PCR驗證Fig. 7 Quantitative RT-PCR verification of CesA9 (T1_Unigene_BMK.22002) gene.

    3 討論與結(jié)論

    測序技術(shù)日趨成熟,水稻[16]和小麥[17]等禾本科比較有代表性植物的基因組測序已經(jīng)完成,并且已經(jīng)轉(zhuǎn)向?qū)嶋H應(yīng)用方面的研究。在竹子中毛竹[18]基因組測序也已經(jīng)完成,為竹類植物的分子生物學研究奠定了良好基礎(chǔ)。轉(zhuǎn)錄組測序已經(jīng)不再僅限于擬南芥等模式植物,越來越多非模式植物的轉(zhuǎn)錄組也開始得到研究,其對應(yīng)的生長發(fā)育狀態(tài)與機理也逐漸被人們所了解。本文利用轉(zhuǎn)錄組測序技術(shù),對處于早期生長的4個不同高度的慈竹筍進行轉(zhuǎn)錄組從頭測序,得到了非冗余的111 137條慈竹筍的Unigene。COG和KEGG功能注釋表明,除了基本的功能外,這些Unigene涉及到蔗糖轉(zhuǎn)運與代謝、次級代謝產(chǎn)物 (如與木質(zhì)素相關(guān)的苯丙烷類的生物合成) 及細胞壁的生物合成等方面(圖2,表3)。

    值得一提的是竹類植物MYB類轉(zhuǎn)錄因子數(shù)量明顯多于另外4個物種。有研究指出MYB轉(zhuǎn)錄因子家族能通過調(diào)控參與苯丙烷代謝途徑的相關(guān)基因表達,進而影響木質(zhì)部細胞的形成[19]。此外,MYB類轉(zhuǎn)錄因子也會調(diào)控部分纖維素生物合成相關(guān)的基因表達,如MYB46可直接調(diào)控與次生壁相關(guān)的部分纖維素合酶基因,從而影響纖維素的合成[20]。因此,進一步探討MYB類轉(zhuǎn)錄因子在慈竹中的功能,對慈竹分子遺傳改良具有重要的意義。

    慈竹纖維素含量豐富,是較好的造紙原料,其與硬頭黃[21]等竹種在四川已形成了以其為龍頭的竹漿造紙工業(yè)產(chǎn)業(yè)鏈。竹材中的纖維素含量直接影響竹漿造紙的效率。10、50、100和150 cm高度的筍中,纖維素含量具有明顯差異。因此本文進一步深入分析了4個樣品中與纖維素合成相關(guān)的差異表達基因。纖維素的生物合成是一個高度復(fù)雜的生物過程,有多個酶參與其生物合成,其中包括CesA、Csl、SuSy等[22]。CesA是目前研究最多的纖維素生物合成基因,它編碼的纖維素合酶不僅是纖維素生物合成的場所,而且在纖維素生物合成過程中起著關(guān)鍵的調(diào)節(jié)作用。纖維素合酶本身具有多基因現(xiàn)象,在纖維素生物合成過程中,多個CesA基因彼此相關(guān),協(xié)同作用[23]。在慈竹Unigene數(shù)據(jù)庫中,分別與Swiss-port和Nr數(shù)據(jù)庫比對時,與CesA基因同源的Unigene分別有49條和33條;與Csl基因同源的Unigene分別有36條和27條;與SuSy基因同源的Unigene分別有14條和7條。這些基因序列與毛竹、水稻、玉米等物種的CesA、Csl和SuSy基因具有很高的同源性。進一步分析表明,不同高度慈竹筍的纖維素合酶CesA4、CesA7和CesA9基因存在差異表達且隨著筍的伸長,總體呈上升的趨勢 (圖6)。有研究表明,竹類植物筍期的高速生長與筍的居間分生組織中細胞的快速生長與伸長有關(guān),初生壁的發(fā)育確保了細胞的正常生長,而次生壁的發(fā)育所起到的支撐作用也確保了筍的快速生長[24-25]。另有研究表明,不同的CesA基因分別在細胞初生壁和次生壁合成時期起不同的作用[23]。水稻的OsCesA4、OsCesA7和OsCesA9是細胞次生壁合成所必需的,在木質(zhì)部的細胞中高度表達。OsCesA4、OsCesA7或OsCesA9突變后會導致水稻莖桿變脆且纖維素含量降低,進而影響水稻的正常生長[26]。因此推測慈竹的纖維素合酶基因CesA4、CesA7和CesA9與細胞次生壁的合成有關(guān)且可能影響慈竹筍的生長發(fā)育,這有待進一步研究。

    REFERENCES

    [1] Peng ZH, Lu Y, Li LB, et al. The draft genome of the fast-growing non-timber forest species moso bamboo (Phyllostachys heterocycla). Nat Genet, 2013, 45(4): 456-463.

    [2] Peng ZH, Zhang CL, Zhang Y, et al. Transcriptome sequencing and analysis of the fast growing shoots of moso bamboo (Phyllostachys edulis). PLoS ONE, 2013, 8(11): e78944.

    [3] Liu MY, Qiao GR, Jiang J, et al. Transcriptome sequencing and De novo analysis for ma bamboo (Dendrocalamus latiflorus Munro) using the illumina platform. PLoS ONE, 2012, 7(10): e46766.

    [4] Qi JQ, Chi B, Xie JL, et al. Study on variations of fiber morphology and tissue proportion of Neosinocalamus affinis culm. Trans China Pulp Paper, 2013, 28(3): 1-4 (in Chinese).齊錦秋, 池冰, 謝九龍, 等. 慈竹纖維形態(tài)及組織比量的研究. 中國造紙學報, 2013, 28(3): 1-4.

    [5] Liu Q, He H, Shen ZP. Study on the growth status and correlative environment factors of bamboo, Neosinocalamus affinis from Chengdu region. Sichuan Environ, 2001, 20(4): 43-46 (in Chinese).劉慶, 何海, 沈昭萍. 成都地區(qū)慈竹生長狀況及其與環(huán)境因子關(guān)系的初步分析. 四川環(huán)境, 2001, 20(4): 43-46.

    [6] Hu SL, Cao Y, Huang SX, et al. Cloning and bioinformation analysis of 4CL gene in Neosinocalamus affinis. J Northwest A&F Univ: Nat Sci Ed, 2009, 37(8): 204-210 (in Chinese).胡尚連, 曹穎, 黃勝雄, 等. 慈竹4CL基因的克隆及其生物信息學分析. 西北農(nóng)林科技大學學報: 自然科學版, 2009, 37(8): 204-210.

    [7] Chen QB, Jiang Y, Lu XQ, et al. Assessment of genetic diversity in Neosinocalamus affinis from different regions of Sichuan province. J Northwest A&F Univ: Nat Sci Ed, 2009, 37(6): 187-193 (in Chinese).陳其兵, 蔣瑤, 盧學琴, 等. 四川不同地區(qū)慈竹的遺傳多樣性研究. 西北農(nóng)林科技大學學報: 自然科學版, 2009, 37(6): 187-193.

    [8] Hu SL, Jiang Y, Chen QB, et al. Physical and chemical properties of 2 species in bamboos from the different regions in Sichuan province. Bull Botan Res, 2010, 30(6): 708-712 (in Chinese).胡尚連, 蔣瑤, 陳其兵, 等. 四川2種叢生竹理化特性及纖維形態(tài)研究. 植物研究, 2010, 30(6): 708-712.

    [9] Li XL, Lin RS, Fung HL, et al. Chromosome numbers of some caespitose bamboos native in or introduced to China. Acta Phytotaxonom Sin, 2001, 39(5): 433-442 (in Chinese).李秀蘭, 林汝順, 馮學琳, 等. 中國部分叢生竹類染色體數(shù)目報道. 植物分類學報, 2001, 39(5): 433-442.

    [10] Tang WW, Wu XL, Zhou J. Comparison of neutral detergent fiber, acid detergent fiber and acid detergent lignin contents of Dongxiang wild rice and cultivated rice. Agric Sci Technol, 2010, 11(3): 11-14.

    [11] Grabherr MG, Haas BJ, Yassour M, et al. Full-length transcriptome assembly from RNA-Seq data without a reference genome. Nat Biotechnol, 2011, 29(7): 644-652.

    [12] Tatusov RL, Galperin MY, Natale DA, et al. The COG database: a tool for genome-scale analysis of protein functions and evolution. Nucleic Acids Res, 2000, 28(1): 33-36.

    [13] Ashburner M, Ball CA, Blake JA, et al. Gene Ontology: tool for the unification of biology. Nat Genet, 2000, 25(1): 25-29.

    [14] Kanehisa M, Goto S, Kawashima S, et al. The KEGG resource for deciphering the genome. Nucleic Acids Res, 2004, 32(Suppl 1): D277-D280.

    [15] Zhang CY, Long Y, Feng J, et al. Transcriptional regulation of plant genes and it’s significance in biology. Hereditas (Beijing), 2007, 29(7): 793-799 (in Chinese).張椿雨, 龍艷, 馮吉, 等. 植物基因在轉(zhuǎn)錄水平上的調(diào)控及其生物學意義. 遺傳, 2007, 29(7): 793-799.

    [16] Lin CY, Trinh NN, Lin CW, et al. Transcriptome analysis of phytohormone, transporters and signaling pathways in response to vanadium stress in rice roots. Plant Physiol Biochem, 2013, 66: 98-104.

    [17] Zhang J, Liu W, Han H, et al. De novo transcriptome sequencing of Agropyron cristatum to identify available gene resources for the enhancement of wheat. Genomics, 2015, 106(2): 129-136.

    [18] Zhao GZ, Sun HY, Zhao HS, et al. A study of genome sequencing of Phyllostachys edulis and its data applications. World Bamb Ratt, 2015, 13(3): 8-13 (in Chinese).趙廣枝, 孫化雨, 趙韓生, 等. 毛竹基因組測序及數(shù)據(jù)應(yīng)用研究現(xiàn)狀. 世界竹藤通訊, 2015, 13(3): 8-13.

    [19] Zhao CS, Craig JC, Petzold HE, et al. The xylem and phloem transcriptomes from secondary tissues of the Arabidopsis root-hypocotyl. Plant Physiol, 2005, 138(2): 803-818.

    [20] Kim WC, Ko JH, Kim JY, et al. MYB46 directly regulates the gene expression of secondary wall associated cellulose synthases in Arabidopsis. Plant J, 2013, 73(1): 26-36.

    [21] Chen QB, Gao SP, Liu L. The selection of paper-pulp bamboo species and the development of bamboo paper sector in Sichuan. J Bamb Res, 2002, 21(4): 47-51 (inChinese).陳其兵, 高素萍, 劉麗. 四川省優(yōu)良紙漿竹種選擇與竹紙產(chǎn)業(yè)化發(fā)展. 竹子研究匯刊, 2002, 21(4): 47-51.

    [22] Cheng X, Hao HQ, Peng L. Recent progresses on cellulose synthesis in cell wall of plants. J Trop Subtr Bot, 2011, 19(3): 283-290 (in Chinese).程曦, 郝懷慶, 彭勵. 植物細胞壁中纖維素合成的研究進展. 熱帶亞熱帶植物學報, 2011, 19(3): 283-290.

    [23] McFarlane HE, D?ring A, Persson S. The cell biology of cellulose synthesis. Annu Rev Plant Biol, 2014, 65: 69-94.

    [24] Gan XH, Ding YL. Investigation on the variation of fiber wall in Phyllostachys edulis culms. For Res, 2006, 19(4): 457-462 (in Chinese).甘小洪, 丁雨龍. 毛竹莖稈纖維發(fā)育過程中細胞壁的變化規(guī)律研究. 林業(yè)科學研究, 2006, 19(4): 457-462.

    [25] Cui K, He CY, Zhang JG, et al. Characteristics of temporal and spatial tissue development during the rapidly growing stage of moso bamboo culms. For Res, 2012, 25(4): 425-431 (in Chinese).崔凱, 何彩云, 張建國, 等. 毛竹莖稈組織速生的時空發(fā)育特征. 林業(yè)科學研究, 2012, 25(4): 425-431.

    [26] Tanaka K, Murata K, Yamazaki M, et al. Three distinct rice cellulose synthase catalytic subunit genes required for cellulose synthesis in the secondary wall. Plant Physiol, 2003, 133(1): 73-83.

    (本文責編 陳宏宇)

    March 9, 2016; Accepted: July 20, 2016

    Shanglian Hu. Tel/Fax: +86-138-81194095; E-mail: hushanglian@126.com

    Transcriptome analysis and gene function annotation of Bambusa emeiensis shoots based on high-throughput
    sequencing technology

    Yupeng Chen, Ying Cao, Shanglian Hu, Yan Huang, Xueqin Lu, Gang Xu, and Zhijian Long
    Engineering Research Center for Biomass Resource Utilization and Modification of Sichuan Province, School of Life Science and Engineering, Southwest University of Science and Technology, Mianyang 621010, Sichuan, China

    Bambusa emeiensis is one of the preponderant species of sympodial bamboos in Sichuan province of China, and has excellent fiber length and quality as raw materials for papermaking, textile and other industries. In this study, with the application of Illumina HiSeqTM2000 platform, we analyzed transcriptome in B. emeiensis with different heights of 10, 50, 100 and 150 cm. A total of 69.28 M reads were obtained, and a sum up of 111 137 bands of Unigenes were acquired following de novo stitching, assembly and clustering, among which there were 63 094 bands that had been integrated in the COG, GO, KEGG, Swiss-Prot and Nr databases using annotated methods. These Unigenes not only had general functions, such as transcription and signal transduction, but were also involved in sucrose transport and metabolism, secondary metabolites and cell wall biosynthesis. There was significant difference regarding the expression of cellulose synthase gene in B. emeiensis at different heights, relevant genes were found that might be responsible for the regulation of the growth and development of B. emeiensis as well as the biosynthesis of cellulose and lignin. Our findings could provide some elementary theories for breed improvement of B. emeiensis.

    Bambusa emeiensis, transcriptome, high-throughput sequencing, function annotation

    Supported by: National Natural Science Foundation of China (Nos. 31400257, 31400333), Breeding Program Fund Project by the 13th Five-Year Plan of Sichuan Province, Fund of Engineering Research Center for Biomass Resource Utilizaiton and Modification of Sichuan Province (Nos. 12zxsk07, 13zxsk01, 14tdgc05), Major Project of Education Department in Sichuan Province (No. 16ZA0145), Postgraduate Innovation Fund Project by Southwest University of Science and Technology (No. 15ycx089).

    國家自然科學基金 (Nos. 31400257, 31400333),四川省“十三五”育種公關(guān)項目,四川省生物質(zhì)資源利用與改性工程技術(shù)研究中心基金 (Nos. 12zxsk07, 13zxsk01, 14tdgc05),四川省教育廳科研項目 (No. 16ZA0145),西南科技大學研究生創(chuàng)新基金 (No. 15ycx089) 資助。

    猜你喜歡
    慈竹毛竹竹筍
    金城寺
    文史雜志(2023年2期)2023-05-30 14:21:32
    長足大竹象取食脅迫下慈竹bZIP轉(zhuǎn)錄因子的鑒定及生物信息學分析
    節(jié)骨草和毛竹苗
    竹筍
    竹筍
    乙酰化慈竹的制備及其尺寸穩(wěn)定性研究
    竹子學報(2017年2期)2017-09-15 16:04:31
    小竹筍的夢
    小布老虎(2017年4期)2017-08-10 08:22:40
    挖竹筍
    毛竹護筍養(yǎng)竹技術(shù)探討
    毛竹幼林四季管理技術(shù)
    少妇人妻精品综合一区二区 | 国内少妇人妻偷人精品xxx网站| 日韩欧美在线乱码| 亚洲国产日韩欧美精品在线观看| 99热精品在线国产| 国产又黄又爽又无遮挡在线| 欧美极品一区二区三区四区| 亚洲av电影不卡..在线观看| 精品午夜福利视频在线观看一区| 亚洲精品乱码久久久v下载方式| 日本在线视频免费播放| 国内揄拍国产精品人妻在线| 色综合色国产| 日本成人三级电影网站| 国语自产精品视频在线第100页| 国产高清三级在线| 精品久久久久久久久久免费视频| 大又大粗又爽又黄少妇毛片口| 国产在线男女| 色5月婷婷丁香| 国产探花极品一区二区| 在线a可以看的网站| 中文亚洲av片在线观看爽| 亚洲电影在线观看av| 日本黄色视频三级网站网址| 深夜a级毛片| 搡女人真爽免费视频火全软件 | 午夜免费激情av| 别揉我奶头~嗯~啊~动态视频| 亚洲欧美中文字幕日韩二区| 日本撒尿小便嘘嘘汇集6| 欧美极品一区二区三区四区| 日韩欧美国产在线观看| 成年av动漫网址| 精品人妻一区二区三区麻豆 | 欧美区成人在线视频| 国产老妇女一区| 国产乱人视频| 97人妻精品一区二区三区麻豆| 国内揄拍国产精品人妻在线| 五月玫瑰六月丁香| 久久久精品欧美日韩精品| 在线观看美女被高潮喷水网站| 舔av片在线| 亚洲五月天丁香| 一边摸一边抽搐一进一小说| 综合色丁香网| 中国美女看黄片| 国产精品一区www在线观看| 最新在线观看一区二区三区| 免费看日本二区| 又黄又爽又刺激的免费视频.| 久久久欧美国产精品| 小说图片视频综合网站| 亚洲最大成人中文| 国产探花极品一区二区| 国产精品不卡视频一区二区| 蜜桃久久精品国产亚洲av| 夜夜夜夜夜久久久久| 亚洲av熟女| 国产真实乱freesex| 日本黄色片子视频| 久久久久久久久久久丰满| 国产乱人偷精品视频| a级一级毛片免费在线观看| 精品久久久久久成人av| 国产高清三级在线| 日韩一区二区视频免费看| 国产精品亚洲美女久久久| 人妻制服诱惑在线中文字幕| 亚洲性久久影院| 国产精品国产高清国产av| 永久网站在线| 此物有八面人人有两片| 日韩精品有码人妻一区| 97人妻精品一区二区三区麻豆| 国产精品永久免费网站| 亚洲无线观看免费| 高清午夜精品一区二区三区 | 人人妻,人人澡人人爽秒播| 草草在线视频免费看| 精品久久国产蜜桃| 免费大片18禁| 国产成人影院久久av| 国产中年淑女户外野战色| 成熟少妇高潮喷水视频| 精品久久久久久久久久免费视频| 国产激情偷乱视频一区二区| 日韩欧美在线乱码| 日本成人三级电影网站| 男女啪啪激烈高潮av片| 一区福利在线观看| 亚洲av二区三区四区| 麻豆久久精品国产亚洲av| 欧美三级亚洲精品| www日本黄色视频网| 久久久久久久久大av| 欧美一区二区亚洲| 久久精品国产亚洲av天美| 免费av不卡在线播放| 22中文网久久字幕| 国产综合懂色| 亚洲精品乱码久久久v下载方式| 久久久久久久久久成人| av中文乱码字幕在线| 欧美性猛交黑人性爽| 最近的中文字幕免费完整| 精品久久久久久久久久免费视频| 色综合色国产| 神马国产精品三级电影在线观看| 嫩草影视91久久| 国产一区二区在线av高清观看| 久久久精品欧美日韩精品| 亚洲成人久久爱视频| 免费观看精品视频网站| 欧美不卡视频在线免费观看| 精品一区二区三区人妻视频| 精品久久久久久久久亚洲| 免费大片18禁| 非洲黑人性xxxx精品又粗又长| 国产大屁股一区二区在线视频| av在线亚洲专区| 国产精品电影一区二区三区| 久久精品久久久久久噜噜老黄 | 国产乱人偷精品视频| 国产欧美日韩精品一区二区| 国产精品1区2区在线观看.| 九九在线视频观看精品| 午夜福利18| 国产成人福利小说| 男女视频在线观看网站免费| 久久草成人影院| 天天一区二区日本电影三级| 男人舔奶头视频| 我要搜黄色片| 亚洲天堂国产精品一区在线| 国产av在哪里看| 两个人视频免费观看高清| 尤物成人国产欧美一区二区三区| 久久精品人妻少妇| 免费搜索国产男女视频| 国产大屁股一区二区在线视频| 国产高清视频在线播放一区| 欧美在线一区亚洲| 99精品在免费线老司机午夜| 日韩欧美精品v在线| 美女内射精品一级片tv| 亚洲三级黄色毛片| 亚洲精品成人久久久久久| 午夜福利在线观看吧| 中文字幕久久专区| 欧美bdsm另类| av女优亚洲男人天堂| 中文字幕av在线有码专区| 日本色播在线视频| 国产亚洲精品久久久com| 免费一级毛片在线播放高清视频| 日日撸夜夜添| 亚洲国产欧美人成| 免费看日本二区| 国产精品亚洲一级av第二区| 乱人视频在线观看| 国产精品一区二区免费欧美| 2021天堂中文幕一二区在线观| 久久韩国三级中文字幕| 国内精品久久久久精免费| 午夜福利在线在线| 久久国内精品自在自线图片| 国产免费男女视频| 亚洲自偷自拍三级| 最近2019中文字幕mv第一页| 看免费成人av毛片| 色综合亚洲欧美另类图片| 成人av在线播放网站| 伦精品一区二区三区| 菩萨蛮人人尽说江南好唐韦庄 | 色哟哟·www| 国产69精品久久久久777片| 国产毛片a区久久久久| 丰满人妻一区二区三区视频av| 最近2019中文字幕mv第一页| 白带黄色成豆腐渣| 黑人高潮一二区| aaaaa片日本免费| 人人妻人人澡人人爽人人夜夜 | 99久久成人亚洲精品观看| 一级黄色大片毛片| 国产精品女同一区二区软件| 国产精品永久免费网站| 男人舔奶头视频| av卡一久久| 一级毛片aaaaaa免费看小| 精品人妻视频免费看| 男插女下体视频免费在线播放| 亚洲人成网站在线观看播放| 日韩精品有码人妻一区| 秋霞在线观看毛片| 成年女人看的毛片在线观看| 亚洲欧美日韩卡通动漫| 日韩一区二区视频免费看| 欧美日韩精品成人综合77777| 在现免费观看毛片| 国产成人freesex在线 | 亚洲国产精品合色在线| 少妇的逼好多水| 91精品国产九色| 色哟哟哟哟哟哟| 免费看美女性在线毛片视频| 又爽又黄a免费视频| 午夜福利在线在线| 欧美成人一区二区免费高清观看| 一个人看视频在线观看www免费| 亚州av有码| 一级av片app| 成年女人看的毛片在线观看| 一个人观看的视频www高清免费观看| 国产欧美日韩精品一区二区| 亚洲一级一片aⅴ在线观看| 国产日本99.免费观看| 亚洲人成网站在线播放欧美日韩| or卡值多少钱| 啦啦啦观看免费观看视频高清| 亚洲国产精品sss在线观看| 老司机影院成人| 国产精品乱码一区二三区的特点| 99在线人妻在线中文字幕| 免费看日本二区| 亚洲成人久久性| 欧美激情在线99| 啦啦啦韩国在线观看视频| 五月伊人婷婷丁香| a级毛色黄片| 在线天堂最新版资源| 欧美日韩乱码在线| 在线国产一区二区在线| 免费高清视频大片| 尤物成人国产欧美一区二区三区| 日韩中字成人| 午夜福利视频1000在线观看| 免费人成在线观看视频色| 最新中文字幕久久久久| 少妇熟女欧美另类| 久久久久久九九精品二区国产| 久久韩国三级中文字幕| 伊人久久精品亚洲午夜| 99九九线精品视频在线观看视频| 精品久久久久久久人妻蜜臀av| 国产精品久久久久久精品电影| 国产亚洲欧美98| 不卡一级毛片| 丝袜美腿在线中文| 搞女人的毛片| 俄罗斯特黄特色一大片| 日韩大尺度精品在线看网址| 午夜亚洲福利在线播放| 寂寞人妻少妇视频99o| 两个人视频免费观看高清| 午夜视频国产福利| 床上黄色一级片| 尤物成人国产欧美一区二区三区| 午夜a级毛片| 国产日本99.免费观看| 国产极品精品免费视频能看的| av在线播放精品| 看黄色毛片网站| av免费在线看不卡| 精品无人区乱码1区二区| 久久鲁丝午夜福利片| www日本黄色视频网| 99久国产av精品国产电影| 麻豆一二三区av精品| 自拍偷自拍亚洲精品老妇| 两个人视频免费观看高清| 干丝袜人妻中文字幕| 久久久久性生活片| 男女那种视频在线观看| 亚洲乱码一区二区免费版| 亚洲精品国产av成人精品 | 蜜臀久久99精品久久宅男| 国产高清有码在线观看视频| 久久人妻av系列| 亚洲av中文av极速乱| 岛国在线免费视频观看| av黄色大香蕉| 天天躁夜夜躁狠狠久久av| av在线天堂中文字幕| 观看免费一级毛片| 日韩一本色道免费dvd| 亚洲av美国av| 99热全是精品| h日本视频在线播放| 麻豆国产97在线/欧美| 免费看日本二区| 在线天堂最新版资源| 18禁在线播放成人免费| 日韩欧美国产在线观看| 此物有八面人人有两片| 99久久无色码亚洲精品果冻| 国产综合懂色| 亚洲精品粉嫩美女一区| 国产中年淑女户外野战色| 大香蕉久久网| 高清日韩中文字幕在线| 久久久久九九精品影院| 国产欧美日韩精品一区二区| 18禁黄网站禁片免费观看直播| 国模一区二区三区四区视频| 99精品在免费线老司机午夜| 18禁黄网站禁片免费观看直播| 亚洲精品国产av成人精品 | 国产精品免费一区二区三区在线| 亚洲欧美日韩高清专用| 91在线观看av| 婷婷六月久久综合丁香| 久久九九热精品免费| 亚洲欧美成人精品一区二区| 国产一区二区在线观看日韩| 久久99热这里只有精品18| 人人妻人人澡人人爽人人夜夜 | 校园春色视频在线观看| 欧美日本亚洲视频在线播放| 亚洲成a人片在线一区二区| 国产高清视频在线播放一区| 麻豆成人午夜福利视频| 1000部很黄的大片| 人妻少妇偷人精品九色| 欧美一区二区精品小视频在线| 免费观看的影片在线观看| 欧美区成人在线视频| 六月丁香七月| 变态另类成人亚洲欧美熟女| 国产精品久久视频播放| 在线观看66精品国产| 91久久精品国产一区二区成人| 熟妇人妻久久中文字幕3abv| 男女之事视频高清在线观看| 国产免费男女视频| 岛国在线免费视频观看| 看十八女毛片水多多多| videossex国产| 欧美不卡视频在线免费观看| 性插视频无遮挡在线免费观看| 欧美激情久久久久久爽电影| 日日摸夜夜添夜夜爱| 看免费成人av毛片| 久久6这里有精品| 一个人看视频在线观看www免费| 欧美性感艳星| 淫妇啪啪啪对白视频| 2021天堂中文幕一二区在线观| 99热6这里只有精品| 亚洲精华国产精华液的使用体验 | 国产午夜精品久久久久久一区二区三区 | av.在线天堂| 国产精品三级大全| 一级av片app| 99久久久亚洲精品蜜臀av| 我要搜黄色片| 淫妇啪啪啪对白视频| 色噜噜av男人的天堂激情| 可以在线观看的亚洲视频| 欧美激情国产日韩精品一区| 99热网站在线观看| 在线看三级毛片| 伦理电影大哥的女人| 成年女人永久免费观看视频| 精品人妻一区二区三区麻豆 | 免费观看人在逋| 99热这里只有是精品50| 亚洲av中文av极速乱| 国产精品久久久久久av不卡| 黄色欧美视频在线观看| 国产不卡一卡二| 麻豆成人午夜福利视频| 婷婷亚洲欧美| 97在线视频观看| 中文字幕人妻熟人妻熟丝袜美| 色在线成人网| 99在线视频只有这里精品首页| 亚洲最大成人中文| 亚洲精品国产成人久久av| 久久久精品大字幕| 亚洲av中文字字幕乱码综合| 美女大奶头视频| 日本欧美国产在线视频| 久久久久国产网址| 九九热线精品视视频播放| 亚洲第一电影网av| 国产午夜精品久久久久久一区二区三区 | 日韩欧美免费精品| 国产探花在线观看一区二区| 婷婷六月久久综合丁香| 美女高潮的动态| ponron亚洲| 一本久久中文字幕| 国产精品乱码一区二三区的特点| 久久精品综合一区二区三区| 国产大屁股一区二区在线视频| 精品人妻视频免费看| 国产麻豆成人av免费视频| 国产老妇女一区| 男插女下体视频免费在线播放| 亚洲欧美精品自产自拍| 日韩 亚洲 欧美在线| 精品国内亚洲2022精品成人| 欧美一区二区国产精品久久精品| 久久亚洲国产成人精品v| 日本免费a在线| 国产精品乱码一区二三区的特点| 美女内射精品一级片tv| 国内精品美女久久久久久| 久久久久九九精品影院| 亚洲色图av天堂| 日韩成人av中文字幕在线观看 | 亚洲最大成人av| 99热这里只有是精品50| 18禁黄网站禁片免费观看直播| 最近2019中文字幕mv第一页| 女同久久另类99精品国产91| 乱人视频在线观看| 精品久久久噜噜| 国内揄拍国产精品人妻在线| 国产亚洲欧美98| 天天躁夜夜躁狠狠久久av| 亚洲va在线va天堂va国产| 白带黄色成豆腐渣| 成人漫画全彩无遮挡| 精品国内亚洲2022精品成人| 久久久久久伊人网av| 欧美成人免费av一区二区三区| 成人毛片a级毛片在线播放| 日本一二三区视频观看| 色综合色国产| 日韩成人av中文字幕在线观看 | 91久久精品国产一区二区三区| 国产综合懂色| 草草在线视频免费看| 久久草成人影院| 高清毛片免费看| 少妇丰满av| 直男gayav资源| 国产精品,欧美在线| 欧美一级a爱片免费观看看| 久久久精品94久久精品| 两性午夜刺激爽爽歪歪视频在线观看| 亚洲精品影视一区二区三区av| 亚洲电影在线观看av| 日韩 亚洲 欧美在线| 亚洲乱码一区二区免费版| 欧美另类亚洲清纯唯美| 亚洲18禁久久av| 小蜜桃在线观看免费完整版高清| 亚洲七黄色美女视频| 精品不卡国产一区二区三区| 国产精品美女特级片免费视频播放器| 亚洲国产精品sss在线观看| 搡老熟女国产l中国老女人| 男女边吃奶边做爰视频| 国产精品99久久久久久久久| 一边摸一边抽搐一进一小说| 久久久久久九九精品二区国产| 成年免费大片在线观看| 久久久成人免费电影| av黄色大香蕉| 色综合站精品国产| 国产精品国产高清国产av| 亚洲va在线va天堂va国产| 一进一出抽搐gif免费好疼| 久久精品国产亚洲网站| 草草在线视频免费看| 黄色欧美视频在线观看| 一级毛片久久久久久久久女| 在线a可以看的网站| 别揉我奶头~嗯~啊~动态视频| 午夜影院日韩av| 黄色欧美视频在线观看| 免费看美女性在线毛片视频| 国产老妇女一区| 精品人妻熟女av久视频| 美女免费视频网站| 级片在线观看| 欧美一区二区国产精品久久精品| 亚洲精品亚洲一区二区| 极品教师在线视频| 国产精品人妻久久久影院| 日韩在线高清观看一区二区三区| 欧美成人a在线观看| 亚洲成人久久性| 女人被狂操c到高潮| 日韩制服骚丝袜av| 国产激情偷乱视频一区二区| 嫩草影院新地址| 成人特级黄色片久久久久久久| 国产麻豆成人av免费视频| www日本黄色视频网| 十八禁国产超污无遮挡网站| 成人漫画全彩无遮挡| av福利片在线观看| 久99久视频精品免费| 国产真实伦视频高清在线观看| 午夜久久久久精精品| 国产黄a三级三级三级人| 内射极品少妇av片p| 亚洲精品国产av成人精品 | 免费无遮挡裸体视频| 好男人在线观看高清免费视频| 2021天堂中文幕一二区在线观| 国内久久婷婷六月综合欲色啪| 色哟哟·www| 日韩,欧美,国产一区二区三区 | 嫩草影院入口| 日本黄大片高清| 麻豆精品久久久久久蜜桃| 久久久成人免费电影| 国产在线男女| 美女内射精品一级片tv| 97碰自拍视频| 国产乱人偷精品视频| 亚洲无线在线观看| 免费av不卡在线播放| 51国产日韩欧美| 一a级毛片在线观看| 啦啦啦观看免费观看视频高清| 日本免费a在线| 免费观看在线日韩| 久久精品国产99精品国产亚洲性色| 亚洲一级一片aⅴ在线观看| 中国美白少妇内射xxxbb| 性色avwww在线观看| 蜜臀久久99精品久久宅男| 尤物成人国产欧美一区二区三区| 国产一区二区三区在线臀色熟女| 亚洲成人av在线免费| 最近的中文字幕免费完整| 久久这里只有精品中国| 午夜免费男女啪啪视频观看 | 嫩草影视91久久| 国产精品永久免费网站| 桃色一区二区三区在线观看| 99久久中文字幕三级久久日本| 日韩成人伦理影院| 最后的刺客免费高清国语| 久久久午夜欧美精品| 又爽又黄a免费视频| 亚洲性夜色夜夜综合| 欧美日韩综合久久久久久| 欧美色欧美亚洲另类二区| av免费在线看不卡| 欧美人与善性xxx| 久久久色成人| 色噜噜av男人的天堂激情| 99久久精品一区二区三区| 99久久无色码亚洲精品果冻| 男人舔奶头视频| 精品人妻一区二区三区麻豆 | av专区在线播放| 最近中文字幕高清免费大全6| 三级国产精品欧美在线观看| 一卡2卡三卡四卡精品乱码亚洲| 久久精品影院6| 99久久九九国产精品国产免费| 在线免费十八禁| 搡女人真爽免费视频火全软件 | 国产欧美日韩精品亚洲av| 久久久久国产网址| 日韩 亚洲 欧美在线| www.色视频.com| 国产午夜精品论理片| 久久久色成人| 久久精品久久久久久噜噜老黄 | 校园春色视频在线观看| 最近手机中文字幕大全| 国产真实乱freesex| 精品午夜福利视频在线观看一区| 99久久九九国产精品国产免费| 有码 亚洲区| 日韩欧美国产在线观看| 亚洲成人中文字幕在线播放| 99热这里只有是精品在线观看| 你懂的网址亚洲精品在线观看 | 亚洲欧美日韩高清在线视频| 一级a爱片免费观看的视频| 成人性生交大片免费视频hd| 又黄又爽又免费观看的视频| 美女xxoo啪啪120秒动态图| 看片在线看免费视频| 97超视频在线观看视频| 精品99又大又爽又粗少妇毛片| 欧美日本亚洲视频在线播放| 久99久视频精品免费| 亚洲av二区三区四区| 激情 狠狠 欧美| 成人一区二区视频在线观看| 亚洲av成人av| 一级毛片aaaaaa免费看小| 日本精品一区二区三区蜜桃| 天堂网av新在线| 中国国产av一级| 免费看光身美女| 深爱激情五月婷婷| 看片在线看免费视频| 有码 亚洲区| 午夜视频国产福利| 成年女人永久免费观看视频| 亚洲欧美成人精品一区二区| 色哟哟·www| 欧美潮喷喷水| 精品午夜福利在线看| 真实男女啪啪啪动态图| 免费电影在线观看免费观看| 亚洲av电影不卡..在线观看| 婷婷色综合大香蕉| 美女高潮的动态| 日韩欧美精品v在线| 亚洲成人久久爱视频|