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

    基于重測(cè)序的陸地棉InDel標(biāo)記開(kāi)發(fā)與評(píng)價(jià)

    2019-01-17 01:28:14溫天旺林忠旭
    作物學(xué)報(bào) 2019年2期
    關(guān)鍵詞:棉區(qū)種質(zhì)多態(tài)性

    吳 迷 汪 念 沈 超 黃 聰 溫天旺 林忠旭

    ?

    基于重測(cè)序的陸地棉InDel標(biāo)記開(kāi)發(fā)與評(píng)價(jià)

    吳 迷 汪 念 沈 超 黃 聰 溫天旺 林忠旭*

    華中農(nóng)業(yè)大學(xué)植物科學(xué)技術(shù)學(xué)院 / 作物遺傳改良國(guó)家重點(diǎn)實(shí)驗(yàn)室, 湖北武漢 430070

    堿基插入/缺失(InDel)是基因組中豐富的遺傳變異形式。InDel以其密度高、易于基因型分型等優(yōu)點(diǎn)成為分子標(biāo)記開(kāi)發(fā)的理想來(lái)源。本研究利用262份陸地棉品系重測(cè)序數(shù)據(jù)鑒定的InDel位點(diǎn), 在全基因組范圍內(nèi)設(shè)計(jì)了3206個(gè)InDel標(biāo)記并挑選均勻分布的320個(gè)標(biāo)記進(jìn)行驗(yàn)證。320個(gè)標(biāo)記篩選出87個(gè)多態(tài)性標(biāo)記, 多態(tài)性率為26.88%。利用多態(tài)性標(biāo)記對(duì)不同地理來(lái)源的262份陸地棉種質(zhì)資源進(jìn)行基因分型, 共檢測(cè)到160個(gè)等位位點(diǎn); 多態(tài)性信息含量(PIC)為0.0836~0.3750, 平均值為0.3073; 基因多樣性指數(shù)變異范圍為0.0874~0.5000, 平均值為0.3876, 表明我國(guó)陸地棉遺傳基礎(chǔ)相對(duì)狹窄。群體結(jié)構(gòu)分析將262份陸地棉品系大致劃分為2個(gè)亞群, 聚類分析和主成分分析的結(jié)果與之基本一致。采用混合線性模型(Mixed linear model)對(duì)6個(gè)纖維品質(zhì)性狀的關(guān)聯(lián)分析檢測(cè)到65個(gè)關(guān)聯(lián)位點(diǎn)(< 0.01), 各位點(diǎn)對(duì)表型變異貢獻(xiàn)率為2.57%~8.12%。本研究旨在利用重測(cè)序數(shù)據(jù)開(kāi)發(fā)全基因組范圍的可用于凝膠檢測(cè)的InDel標(biāo)記, 為棉花種質(zhì)資源研究和分子標(biāo)記輔助選擇育種提供便捷工具。

    陸地棉; InDel標(biāo)記; 遺傳多樣性; 群體結(jié)構(gòu); 關(guān)聯(lián)分析

    棉花作為世界上重要的經(jīng)濟(jì)作物, 是紡織制造業(yè)中天然纖維的主要來(lái)源。陸地棉()是棉花的主要栽培種。我國(guó)棉花育種史較短, 前期改良主要依賴于國(guó)外引種[1], 當(dāng)前選育和種植的大部分陸地棉品種多由美國(guó)引進(jìn)的岱字棉、斯字棉、金字棉等少數(shù)種質(zhì)資源衍生而來(lái), 這也導(dǎo)致我國(guó)陸地棉種質(zhì)資源相對(duì)匱乏、遺傳基礎(chǔ)較為狹窄、多樣性水平低[2]。為有效拓寬陸地棉的遺傳基礎(chǔ)、加快不同育種目標(biāo)下的棉花新品種選育及遺傳改良進(jìn)程, 將傳統(tǒng)育種方法和分子標(biāo)記輔助選擇育種方法結(jié)合是有效可行的。目前, 多種分子標(biāo)記如簡(jiǎn)單序列重復(fù)(SSR)、單核苷酸多態(tài)性(SNP)等廣泛地應(yīng)用于棉花種質(zhì)資源遺傳多樣性評(píng)價(jià)、親緣關(guān)系鑒定、遺傳圖譜的構(gòu)建、重要農(nóng)藝性狀QTL定位及圖位克隆和種質(zhì)鑒定等方面[3-4]。

    InDel標(biāo)記是根據(jù)不同個(gè)體基因組同源序列發(fā)生的核苷酸片段插入或缺失而開(kāi)發(fā)的, 具有基因組內(nèi)分布廣、密度高、變異穩(wěn)定、可重復(fù)性好、價(jià)格低廉等優(yōu)點(diǎn), 在水稻[5]、玉米[6]、油菜[7]等作物中被廣泛應(yīng)用。2015年Zhang等[8]完成了異源四倍體陸地棉品種‘TM-1’全基因組測(cè)序及組裝工作。2017年Wang等[9]完成了352份棉花材料重測(cè)序, 同時(shí)發(fā)掘大量InDel變異位點(diǎn), 為全基因組范圍的InDel標(biāo)記開(kāi)發(fā)提供了條件。

    棉花許多重要的農(nóng)藝性狀如產(chǎn)量、纖維品質(zhì)、生育期等均為數(shù)量性狀, 解析復(fù)雜數(shù)量性狀的遺傳基礎(chǔ)主要有連鎖分析和關(guān)聯(lián)分析兩種方法。關(guān)聯(lián)分析是一種以連鎖不平衡為基礎(chǔ), 鑒定自然群體目標(biāo)性狀與遺傳標(biāo)記關(guān)系的分析方法。它具有無(wú)需構(gòu)建作圖群體、省時(shí)省力、解析精度高等優(yōu)點(diǎn), 也是傳統(tǒng)連鎖分析QTL定位的有益補(bǔ)充。在棉花中, 國(guó)內(nèi)外廣泛開(kāi)展了關(guān)聯(lián)分析相關(guān)的研究工作[9-11], 然而這些研究大部分是基于SNP的, 利用InDel進(jìn)行關(guān)聯(lián)分析的研究相對(duì)較少。

    本研究基于Wang等[9]的重測(cè)序結(jié)果, 篩選InDel變異位點(diǎn), 開(kāi)發(fā)覆蓋全基因組范圍的、可用于凝膠檢測(cè)的InDel標(biāo)記, 并對(duì)重測(cè)序的262份陸地棉材料進(jìn)行遺傳多樣性、親緣關(guān)系和群體結(jié)構(gòu)的分析, 以及通過(guò)關(guān)聯(lián)分析挖掘與纖維品質(zhì)性狀相關(guān)聯(lián)的InDel標(biāo)記, 從而評(píng)價(jià)這些標(biāo)記在棉花遺傳和育種中的應(yīng)用潛力。

    1 材料與方法

    1.1 供試材料及性狀表型數(shù)據(jù)

    采用262份陸地棉種質(zhì)資源材料, 根據(jù)地理來(lái)源分類, 包括中國(guó)長(zhǎng)江流域棉區(qū)的75份、中國(guó)黃河流域棉區(qū)的118份、中國(guó)西北內(nèi)陸棉區(qū)的33份、中國(guó)北方特早熟棉區(qū)的17份、中國(guó)南方棉區(qū)的2份, 以及美國(guó)的11份和前蘇聯(lián)的6份[9]。依據(jù)Nie等[11]的SSR標(biāo)記在503份種質(zhì)材料中的擴(kuò)增結(jié)果, 挑選了24份多態(tài)性高的材料用于后續(xù)多態(tài)性標(biāo)記篩選, 這些材料比較全面地覆蓋了種質(zhì)材料的多態(tài)性信息, 利用24份材料作為參考, 可以快速篩選出具有代表性的多態(tài)性標(biāo)記。24份種質(zhì)材料編號(hào)分別為ZY40、ZY43、ZY54、ZY91、ZY121、ZY122、ZY211、ZY238、ZY240、ZY262、ZY320、ZY358、ZY361、ZY362、ZY376、ZY410、ZY413、ZY423、ZY434、ZY452、ZY468、ZY481、ZY495和ZY507。

    種質(zhì)材料的表型數(shù)據(jù)在Nie等[11]的文章中已報(bào)道, 從中獲取了2年4點(diǎn)共8個(gè)環(huán)境下的表型數(shù)據(jù)及最佳線性無(wú)偏預(yù)測(cè)(BLUP)后的表型值, 8個(gè)環(huán)境分別為2012年湖北黃岡(E1)、2012年河南原陽(yáng)(E2)、2012年新疆庫(kù)爾勒(E3)、2012年新疆石河子(E4)、2013年湖北黃岡(E5)、2013年河南原陽(yáng)(E6)、2013年新疆庫(kù)爾勒(E7)和2013年新疆石河子(E8)。選取與纖維品質(zhì)相關(guān)的6個(gè)表型數(shù)據(jù), 分別為纖維上半部平均長(zhǎng)度(FUHML)、纖維比強(qiáng)度(FS)、纖維伸長(zhǎng)率(FE)、纖維整齊度(FU)、馬克隆值(MV)和短纖維率(SF)。

    1.2 InDel標(biāo)記設(shè)計(jì)

    基于Wang等[9]檢測(cè)的InDel位點(diǎn)集合, 篩選出測(cè)序深度大于7、插入/缺失堿基數(shù)大于30 bp且最小等位基因頻率大于0.05的InDel位點(diǎn)。提取InDel位點(diǎn)上下游150 bp序列, 將其中含有“N”的位點(diǎn)過(guò)濾。以Zhang等[8]發(fā)表的陸地棉標(biāo)準(zhǔn)品系‘TM-1’基因組序列為參考進(jìn)行比對(duì), 閾值設(shè)置為1×10–10, 選取完全匹配的InDel位點(diǎn)作為候選標(biāo)記開(kāi)發(fā)位點(diǎn)。使用Primer 3.0[12]設(shè)計(jì)引物, 引物長(zhǎng)度范圍為18~24 bp, GC含量范圍在40%~60%之間, 退火溫度范圍為54~60℃。最后, 去除位于scaffold上的引物及SSR引物。引物由金斯瑞生物科技有限公司合成。

    1.3 基因組DNA的提取、PCR擴(kuò)增及電泳檢測(cè)

    從田間取幼嫩葉片, 參照Li等[13]改良的CTAB法提取基因組DNA。以微量分光光度計(jì)Nanodrop 2000檢測(cè)其質(zhì)量和濃度, 并稀釋至終濃度20~40 ngμL–1, 于–20℃保存?zhèn)溆谩?/p>

    采用10 μL PCR擴(kuò)增體系, 其中包括1×reaction buffer、0.8 mmol L–1MgCl2、0.25 mmol L–1dNTPs、2.5 nmol L–1primers、25 ng DNA、0.5 UDNA polymerase (上海博彩生物科技有限公司)。擴(kuò)增程序?yàn)?5℃預(yù)變性5 min; 然后95℃變性50 s, 退火45 s (根據(jù)不同引物設(shè)置退火溫度), 72℃延伸60 s, 34個(gè)循環(huán); 最后72℃延伸5 min, 15℃保存。用6%非變性聚丙烯酰胺凝膠電泳分離擴(kuò)增產(chǎn)物, 銀染檢測(cè)分離條帶。

    以人工讀帶, 同一位置上清晰且重復(fù)性好的條帶記為“1”, 無(wú)帶記為“0”, 缺失數(shù)據(jù)記為“–”。

    1.4 數(shù)據(jù)統(tǒng)計(jì)分析

    使用POWERMARKER V3.25軟件[14]計(jì)算多態(tài)性信息含量(polymorphic information content, PIC)和基因多樣性指數(shù); 計(jì)算各種質(zhì)材料間Nei’s遺傳距離,依據(jù)距離矩陣構(gòu)建系統(tǒng)發(fā)生樹(shù), 進(jìn)行Neighbor- Joining聚類分析, iTOL[15]用于繪制和編輯聚類圖。

    利用STRUCTURE 2.3.4軟件[16]計(jì)算各種質(zhì)材料相應(yīng)的Q值, 用于分析自然群體遺傳結(jié)構(gòu)。設(shè)置組群數(shù)值為1~10, 每個(gè)值獨(dú)立重復(fù)運(yùn)算7次, 將MCMC (Markov Chain Monte Carlo)開(kāi)始時(shí)的不作數(shù)迭代(Length of burn-in period)設(shè)為100,000次, 再將不作數(shù)迭代后的MCMC設(shè)為100,000次; 利用STRUCTURE HARVESTER[17]分析輸出數(shù)據(jù), 依據(jù)ln()和Δ值來(lái)確定最佳K值[18]。

    利用TASSEL 5.0軟件[19]進(jìn)行主成分分析(Principal component analysis, PCA), 估算PCA矩陣對(duì)表型變異的解釋率。同時(shí), 以群體結(jié)構(gòu)(Q)和親緣關(guān)系()為協(xié)變量的混合線性模型(MLM)進(jìn)行性狀-基因型的關(guān)聯(lián)分析, 計(jì)算各標(biāo)記位點(diǎn)在<0.01時(shí)對(duì)表型變異的貢獻(xiàn)率(2)。

    2 結(jié)果與分析

    2.1 InDel標(biāo)記開(kāi)發(fā)及多態(tài)性分析

    基于Wang等[9]重測(cè)序檢測(cè)的InDel位點(diǎn), 篩選并開(kāi)發(fā)了3206個(gè)InDel標(biāo)記, 可在網(wǎng)頁(yè)(http://cotton. hzau.edu.cn/EN/download.php)中下載標(biāo)記的相關(guān)信息, 文件名稱為“ZY_InDels data”。在棉花基因組中平均每隔5.88 Mb挑選一個(gè)標(biāo)記, 最終挑選出均勻分布于26條染色體的320個(gè)標(biāo)記進(jìn)行合成。

    利用24份多態(tài)性高的種質(zhì)材料進(jìn)行多態(tài)性標(biāo)記的初步篩選, 檢測(cè)到87個(gè)多態(tài)性標(biāo)記, 占總標(biāo)記數(shù)的26.88%, 詳細(xì)信息見(jiàn)附表1。多態(tài)性標(biāo)記涵蓋了整個(gè)棉花基因組, 每條染色體至少一個(gè)標(biāo)記, 平均每條染色體上3.35個(gè)。圖1為標(biāo)記HAU_ID_ D11-01在部分材料中的擴(kuò)增結(jié)果。其他230個(gè)標(biāo)記無(wú)多態(tài)性, 3個(gè)標(biāo)記沒(méi)有得到擴(kuò)增產(chǎn)物。

    圖1 標(biāo)記HAU_ID_D11-01在部分品系中擴(kuò)增結(jié)果

    2.2 遺傳多樣性分析

    利用87個(gè)多態(tài)性標(biāo)記對(duì)262份棉花種質(zhì)資源材料進(jìn)行基因型分析, 共檢測(cè)到160個(gè)等位變異位點(diǎn), 平均每對(duì)標(biāo)記1.84個(gè)。基因多樣性指數(shù)變幅為0.0874 (HAU_ID_D07-04)~0.5000 (HAU_ID_D05- 08), 平均值為0.3876; 多態(tài)性信息含量(PIC)變幅為0.0836 (HAU_ID_D07-04)~0.3750 (HAU_ID_A07- 02), 平均值為0.3073 (附表1)。

    2.3 群體結(jié)構(gòu)分析

    ln()值隨假定亞群數(shù)值(1~10)的增大而持續(xù)增大。在=2時(shí)Δ達(dá)到最大(圖2), 由此, 在遺傳結(jié)構(gòu)上將262份棉花種質(zhì)材料劃分為POP-1和POP-2兩個(gè)亞群, 繪制供試材料的群體結(jié)構(gòu)圖(圖3)。POP-1包含175份材料, 黃河流域棉區(qū)(84份)和長(zhǎng)江流域棉區(qū)(71份)的種質(zhì)材料占該組的88.57%; POP-2包含87份材料, 黃河流域棉區(qū)(34份)、西北內(nèi)陸棉區(qū)(23份)和北方特早熟棉區(qū)(13份)的種質(zhì)材料占該組的80.46%。每個(gè)亞群都包含不同種植區(qū)域的種質(zhì)材料, 但同時(shí)又以同一種植區(qū)域的種質(zhì)材料為主。

    主成分分析(PCA)中PC1、PC2分別解釋了個(gè)體變異的7.58%和4.55%。以PC1為x軸、PC2為y軸構(gòu)建一個(gè)二維坐標(biāo)圖(圖4)。種質(zhì)材料大致上可分為2個(gè)集群, 一個(gè)以長(zhǎng)江流域棉區(qū)和黃河流域棉區(qū)的種質(zhì)材料為主, 另一個(gè)以黃河流域棉區(qū)、西北內(nèi)陸棉區(qū)和北方特早熟棉區(qū)的種質(zhì)材料為主。不同來(lái)源的材料在一定程度上相對(duì)獨(dú)立, 而伴隨育種過(guò)程中的基因交換和基因滲入, 材料互相混合。PCA觀察結(jié)果與STRUCTURE軟件的分析結(jié)果較為一致。

    圖2 K值與ln P(D)、ΔK值折線圖

    A: ln()與值變化折線圖; B: Δ值隨值變化折線圖。

    A: Line chart of ln() value with change of; B: Line chart of Δvalue with change of.

    2.4 基于遺傳距離的聚類分析

    262份材料的Nei’s遺傳距離(Genetic distance, GD)平均為0.499, 最高為1.075, 最低為0.012。遺傳距離最高的材料均來(lái)自黃河流域棉區(qū), 分別是ZY31 (晉棉13)和ZY297 (晉棉9)。基于距離矩陣構(gòu)建系統(tǒng)發(fā)生樹(shù), Neighbor-Joining聚類圖(圖5)顯示262份陸地棉材料可以劃分為G1 (紅色線條)和G2 (綠色線條)兩大類群, G1包含182份種質(zhì)材料(中國(guó)黃河流域棉區(qū)93份、中國(guó)長(zhǎng)江流域棉區(qū)70份、中國(guó)西北內(nèi)陸棉區(qū)7份、中國(guó)北方特早熟棉區(qū)5份、美國(guó)5份、前蘇聯(lián)2份); G2包含80份種質(zhì)材料(中國(guó)黃河流域棉區(qū)25份、中國(guó)長(zhǎng)江流域棉區(qū)5份、中國(guó)西北內(nèi)陸棉區(qū)26份、中國(guó)北方特早熟棉區(qū)12份、中國(guó)南方棉區(qū)2份、美國(guó)6份、前蘇聯(lián)4份), 該組材料來(lái)源更為豐富, 其平均遺傳距離也更高。總體來(lái)看, 262份種質(zhì)材料的聚類分析結(jié)果與系譜信息相吻合。

    圖3 基于InDel標(biāo)記的262份陸地棉品系群體結(jié)構(gòu)圖

    圖4 基于InDel標(biāo)記的262份陸地棉品系的PCA圖

    NIR: 中國(guó)西北內(nèi)陸棉區(qū); NSEMR: 中國(guó)北方特早熟棉區(qū); SCR: 中國(guó)南方棉區(qū); SU: 前蘇聯(lián); USA: 美國(guó); YRR: 中國(guó)黃河流域棉區(qū); YtRR: 中國(guó)長(zhǎng)江流域棉區(qū)。

    NIR: Northwestern Inland Region of China; NSEMR: Northern Specific Early Maturation Region of China; SCR: South China Region; SU: Former Soviet Union; USA: American; YRR: Yellow River Region of China; YtRR: Yangtze River Region of China.

    圖5 基于遺傳距離的系統(tǒng)發(fā)生樹(shù)

    紅色線條: G1, 包含182份種質(zhì)材料; 綠色線條: G2, 包含80份種質(zhì)材料。

    Red lines: group 1, including 182 accessions; Green lines: group 2, including 80 accessions.

    2.5 性狀-標(biāo)記的關(guān)聯(lián)分析

    采用TASSEL 5.0軟件中MLM(Q+K)模型進(jìn)行性狀-標(biāo)記的關(guān)聯(lián)分析。表型性狀包括E1~E8和最佳線性無(wú)偏預(yù)測(cè)(BLUP)的表型值。在<0.01顯著條件下, 共檢測(cè)到65個(gè)與纖維品質(zhì)性狀相關(guān)聯(lián)的位點(diǎn)(附表2)。與上半部纖維長(zhǎng)度相關(guān)聯(lián)的位點(diǎn)有12個(gè), 在不同環(huán)境中對(duì)表型變異解釋率范圍為2.97% (HAU_ID_D12-10)~5.92% (HAU_ID_D09-10), 平均值為3.43%。對(duì)其余5個(gè)纖維品質(zhì)性狀, FS、FE、FU、MV和SF分別檢測(cè)到7、13、10、9和14個(gè)關(guān)聯(lián)位點(diǎn),表型變異的解釋率范圍分別為2.57%~7.36%、2.61%~5.66%、2.74%~8.12%、2.72%~5.58%和2.66%~ 6.99%。65個(gè)關(guān)聯(lián)位點(diǎn)中, 有23個(gè)至少能在2個(gè)環(huán)境中被檢測(cè)到, 有11個(gè)至少能在4個(gè)環(huán)境中被檢測(cè)到。其中, 標(biāo)記HAU_ID_D02-06與FE在8個(gè)環(huán)境中均被檢測(cè)到關(guān)聯(lián), 標(biāo)記HAU_ID_D12-10與FU在6個(gè)環(huán)境中被檢測(cè)到關(guān)聯(lián), 標(biāo)記HAU_ID_A01-15與FL在5個(gè)環(huán)境中被檢測(cè)到關(guān)聯(lián)。

    同一標(biāo)記位點(diǎn)與多種性狀相關(guān)聯(lián)可能是多效性引起的。統(tǒng)計(jì)發(fā)現(xiàn)19個(gè)標(biāo)記位點(diǎn)為多效應(yīng)位點(diǎn), 這些位點(diǎn)至少與2個(gè)性狀相關(guān)聯(lián)(表1)。其中, 標(biāo)記HAU_ID_D07-09同時(shí)與5個(gè)性狀(FUHML、FS、FE、FU、SF)相關(guān)聯(lián); 標(biāo)記HAU_ID_D12-10同時(shí)與4個(gè)性狀(FUHML、FS、FU、SF)相關(guān)聯(lián); 5個(gè)標(biāo)記同時(shí)與3個(gè)性狀相關(guān)聯(lián), 分別為HAU_ID_A01-15 (FUHML、FU和SF)、HAU_ID_A03-09 (FE、FU和SF)、HAU_ID_A10-09 (FS、FU和MV)、HAU_ID_ D06-06 (FUHML、FU和SF)和HAU_ID_D06-07 (FUHML、FU和SF); 其余12個(gè)位點(diǎn)與2個(gè)性狀相關(guān)聯(lián)。

    表1 多效應(yīng)標(biāo)記位點(diǎn)

    3 討論

    3.1 InDel標(biāo)記開(kāi)發(fā)及應(yīng)用

    InDel變異廣泛存在于基因組中, 其多態(tài)性頻率僅次于SNP。InDel標(biāo)記以其密度大、準(zhǔn)確性高、具共顯性、檢測(cè)的經(jīng)濟(jì)快速和高效等特點(diǎn)成為分子標(biāo)記開(kāi)發(fā)的寶貴資源。2017年Wang等[9]對(duì)352份種質(zhì)材料進(jìn)行重測(cè)序, 我們挑選了其中262份陸地棉材料進(jìn)行InDel標(biāo)記開(kāi)發(fā), 由于中等大小的插入/缺失位點(diǎn)易于利用凝膠電泳進(jìn)行基因分型, 參考Wu等[20]的InDel篩選標(biāo)準(zhǔn), 篩選出插入/缺失片段長(zhǎng)度大于30 bp的InDel位點(diǎn), 并將之轉(zhuǎn)化為基于PCR的可用于凝膠檢測(cè)的InDel標(biāo)記, 最終設(shè)計(jì)了3206個(gè)InDel標(biāo)記并從中挑選了均勻分布于基因組的320個(gè)InDel標(biāo)記用于多態(tài)性驗(yàn)證。

    為了快速篩選多態(tài)性標(biāo)記, 我們根據(jù)Nie等[11]的SSR標(biāo)記在503份種質(zhì)材料中的擴(kuò)增結(jié)果, 挑選出具有代表性的、材料間多態(tài)性豐富的24份種質(zhì)用于標(biāo)記篩選。以320個(gè)InDel篩選到87個(gè)多態(tài)性標(biāo)記, 多態(tài)性比例為26.88%。在陸地棉中, Wang等[21]基于RAD-seq測(cè)序數(shù)據(jù)開(kāi)發(fā)了121對(duì)InDel標(biāo)記, 兩親本間多態(tài)性率為25.62%, 與本研究結(jié)果相似, 均呈現(xiàn)出一個(gè)較低的多態(tài)性水平。從實(shí)驗(yàn)方法上來(lái)看, 本研究?jī)H采用24份種質(zhì)材料篩選多態(tài)性標(biāo)記, 會(huì)導(dǎo)致部分標(biāo)記的多態(tài)性不能被檢測(cè)。同時(shí), 采用聚丙烯酰胺凝膠電泳檢測(cè)標(biāo)記的多態(tài)性, 受凝膠的分辨率及分離效果等因素影響, 導(dǎo)致多態(tài)性標(biāo)記的識(shí)別能力不夠。此外, 陸地棉是異源四倍體物種, 基因組較為復(fù)雜, A亞組和D亞組間同源性程度較高, 在目前的基因組序列質(zhì)量還不夠完善的情況下, 重測(cè)序的序列同基因組序列比對(duì)時(shí)可能將來(lái)自不同亞組的序列比對(duì)到同一位置, 造成假陽(yáng)性, 因而無(wú)法檢測(cè)到標(biāo)記的多態(tài)性。但與以往在陸地棉種內(nèi)開(kāi)發(fā)的SSR標(biāo)記多態(tài)性率(2.5%~6.3%)相比[21-23], 基于重測(cè)序開(kāi)發(fā)的InDel標(biāo)記多態(tài)性率有所提高, 該研究結(jié)果為將來(lái)陸地棉引物的開(kāi)發(fā)提供了參考。

    目前, SSR標(biāo)記及SNP標(biāo)記在棉花的遺傳資源評(píng)價(jià)、種質(zhì)鑒定、遺傳圖譜構(gòu)建和基因定位與圖位克隆[22,24-25]等方面廣泛應(yīng)用, 而InDel標(biāo)記則報(bào)道較少。新一代測(cè)序技術(shù)的發(fā)展和應(yīng)用, 為棉花基因組測(cè)序及高密度InDel標(biāo)記的開(kāi)發(fā)提供了良好的條件[8-9],可以預(yù)見(jiàn)InDel標(biāo)記將有廣闊的應(yīng)用前景。

    3.2 遺傳多樣性分析

    種質(zhì)資源是作物遺傳改良的基礎(chǔ), 由于陸地棉品種同質(zhì)性較強(qiáng), 遺傳基礎(chǔ)相對(duì)狹窄, 分析陸地棉種質(zhì)資源材料的遺傳多樣性、群體結(jié)構(gòu)及親緣關(guān)系對(duì)拓寬棉花種質(zhì)基礎(chǔ)及棉花育種起到重要作用。

    本研究篩選的87個(gè)InDel標(biāo)記均勻分布于棉花基因組上。基因多樣性指數(shù)分布在0.087,4到0.500,0之間; 多態(tài)性信息含量分布在0.083,6到0.387,6之間, 均值為0.307,3, 能夠較為全面地分析陸地棉的遺傳多樣性。此前, Huang等[10]以63,058個(gè)SNPs對(duì)503份種質(zhì)資源材料進(jìn)行分析, PIC值變幅為0.279~0.369, 均值為0.332。Fang等[24]挑選448對(duì)SSR標(biāo)記分析了193份地理來(lái)源不同的陸地棉材料, PIC均值為0.29, 根據(jù)資料考證, 他們指出自20世紀(jì)以來(lái)美國(guó)陸地棉品種的遺傳多樣性逐漸下降。Ai等[26]和Tyagi等[27]的結(jié)果也表明陸地棉遺傳多樣性水平普遍偏低。本研究與前人研究結(jié)果基本一致。

    3.3 聚類分析和群體結(jié)構(gòu)分析

    基于InDel標(biāo)記的聚類分析將262份種質(zhì)材料大致劃分為兩類(圖5)。聚類結(jié)果往往受到群體大小及群體內(nèi)成員構(gòu)成差異等因素的影響, 2017年Wang等[9]基于全基因組SNP標(biāo)記將352份種質(zhì)材料劃分為兩類, 能有效地將中國(guó)與國(guó)外的種質(zhì)及野生種質(zhì)區(qū)分開(kāi)。本研究所采用的262份種質(zhì)材料是從352份材料中挑選得來(lái)的陸地棉, 剔除了大量野生種質(zhì)及國(guó)外種質(zhì), 僅與Wang等[9]的國(guó)內(nèi)種質(zhì)材料聚類模式相比, 聚類結(jié)果較為一致, 均能將多數(shù)西北內(nèi)陸棉區(qū)、北方特早熟棉區(qū)的材料與其他材料區(qū)分開(kāi)。同時(shí), 我們發(fā)現(xiàn)G1組89.6%的材料來(lái)源于長(zhǎng)江流域和黃河流域, 由于地理位置相近、基礎(chǔ)種質(zhì)資源相同或相似及育種過(guò)程中的人工定向選擇, 使得這兩大棉區(qū)的材料遺傳距離相近, 遺傳差異較小, 更加傾向于聚在一簇。G2組中不同地理來(lái)源的材料互相混合, 多數(shù)材料遺傳距離較遠(yuǎn)。后續(xù)可考慮通過(guò)加強(qiáng)不同棉區(qū)的種質(zhì)交流、加強(qiáng)國(guó)外引種力度來(lái)有效豐富我國(guó)棉花品種的遺傳多樣性。

    對(duì)群體結(jié)構(gòu)的正確評(píng)估是進(jìn)行復(fù)雜性狀關(guān)聯(lián)分析的前提[28]。群體結(jié)構(gòu)分析結(jié)果表明262份材料可大致分為2個(gè)亞群, 與聚類分析結(jié)果相對(duì)比, POP-1與G1、POP-2與G2共有的種質(zhì)材料分別為161份和66份, 重疊率分別為88.46%和82.50%。本研究結(jié)合了STRUCTURE、聚類分析及主成分分析(PCA)這3種方法對(duì)群體結(jié)構(gòu)進(jìn)行分析, 結(jié)果一致性較好, 可有效避免由群體結(jié)構(gòu)分層、等位基因分布不均等造成的標(biāo)記與性狀間的虛假關(guān)聯(lián)。

    3.4 關(guān)聯(lián)分析

    產(chǎn)量和纖維品質(zhì)是陸地棉主要的育種目標(biāo)性狀, 多是復(fù)雜的數(shù)量性狀, 受環(huán)境影響較大, 遺傳基礎(chǔ)復(fù)雜。近年來(lái), 關(guān)聯(lián)分析廣泛應(yīng)用到棉花各性狀尤其是纖維品質(zhì)性狀QTL定位中, 并挖掘到很多與陸地棉纖維品質(zhì)性狀相關(guān)的QTL[4, 9-11]。不同的關(guān)聯(lián)分析方法得到的結(jié)果不盡相同。以群體結(jié)構(gòu)(Q)和親緣關(guān)系(K)作為協(xié)變量的MLM (Q+K)優(yōu)于GLM (Q)和MLM(K), 能有效降低關(guān)聯(lián)的假陽(yáng)性, 提高檢測(cè)結(jié)果的準(zhǔn)確度和真實(shí)性[29]。

    本研究利用多環(huán)境的表型數(shù)據(jù)進(jìn)行標(biāo)記-性狀的關(guān)聯(lián)分析。36個(gè)標(biāo)記共檢測(cè)到65個(gè)關(guān)聯(lián)位點(diǎn)(<0.01), 各位點(diǎn)對(duì)表型變異的解釋率為2.57%~ 8.12%。標(biāo)記的開(kāi)發(fā)是基于Wang等[9]的重測(cè)序數(shù)據(jù), 故參照其結(jié)果設(shè)置LD, 將置信區(qū)間設(shè)置為±296 kb。與前人的研究結(jié)果比較后發(fā)現(xiàn)8個(gè)位點(diǎn)已有報(bào)道, 其中4個(gè)位點(diǎn)與本研究關(guān)聯(lián)到相同性狀。其余28個(gè)位點(diǎn)可能是一些新的標(biāo)記位點(diǎn), 有待進(jìn)一步的實(shí)驗(yàn)驗(yàn)證。

    4個(gè)具有相同關(guān)聯(lián)性狀的標(biāo)記中, 一個(gè)標(biāo)記(HAU_ID_D05-02)與馬克隆值相關(guān), 且在5個(gè)環(huán)境中檢測(cè)到關(guān)聯(lián), 該標(biāo)記位點(diǎn)(Pos: D05_13953192)和Ma等[30]檢測(cè)到的位點(diǎn)(Pos: D05_13941003)十分相近, 表型變異解釋率達(dá)5.58%。其他3個(gè)標(biāo)記(HAU_ ID_A07-02、HAU_ID_D02-01、HAU_ID_D04-07)都與纖維伸長(zhǎng)率相關(guān), 其中, 標(biāo)記HAU_ID_A07-02在4個(gè)環(huán)境中檢測(cè)到關(guān)聯(lián), 表型變異解釋率為3.95%,該位點(diǎn)被Zhang等[31]和Fang等[32]共同檢測(cè); 標(biāo)記HAU_ID_D02-01表型變異解釋率為3.10%, 與Ma等[30]共同檢測(cè); 標(biāo)記HAU_ID_D04-07在5個(gè)環(huán)境中檢測(cè)到關(guān)聯(lián), 表型變異解釋率達(dá)5.22%, 該關(guān)聯(lián)位點(diǎn)在Wang等[9]、Huang等[10]以及其他兩篇報(bào)道[33-34]中共同檢測(cè)。2013年Liu等[35]鑒定到一個(gè)參與赤霉素信號(hào)響應(yīng)的基因GhGASL, 組織表達(dá)模式分析表明該基因在棉花纖維中特異性表達(dá), 可能參與纖維發(fā)育早期的纖維細(xì)胞的伸長(zhǎng), 將該基因的特異性引物比對(duì)到‘TM-1’基因組[8]上, 引物位點(diǎn)(Pos: D04_47893998)與本研究的標(biāo)記位點(diǎn)HAU_ID_ D04-07 (Pos: D04_47451029)緊密連鎖。

    在不同報(bào)道中能夠重復(fù)檢測(cè)到的位點(diǎn)支撐了本研究的關(guān)聯(lián)結(jié)果, 這些關(guān)聯(lián)位點(diǎn)具有可重復(fù)性和穩(wěn)定性, 可用于QTL精細(xì)定位中候選區(qū)間的確定, 同時(shí)有助于后期分子設(shè)計(jì)育種中材料的選擇及鑒定。此外, 本研究檢測(cè)到19個(gè)多效應(yīng)位點(diǎn)同樣值得關(guān)注, 例如位于D12染色體的標(biāo)記HAU_ID_D12-10, 與之關(guān)聯(lián)的4個(gè)性狀(FUHML、FS、FU、SF)能在多個(gè)環(huán)境中穩(wěn)定檢測(cè), 且標(biāo)記對(duì)性狀表型變異解釋率處于相對(duì)較高的水平。這些多效應(yīng)位點(diǎn)為解析不同性狀之間的基因網(wǎng)絡(luò)及遺傳聯(lián)系提供參考依據(jù), 對(duì)指導(dǎo)育種中多性狀協(xié)同改良具有重要價(jià)值。

    4 結(jié)論

    基于全基因組重測(cè)序數(shù)據(jù)設(shè)計(jì)了3206個(gè)InDel標(biāo)記, 挑選320個(gè)標(biāo)記進(jìn)行合成并檢測(cè)到87個(gè)多態(tài)性標(biāo)記, 在陸地棉種內(nèi)標(biāo)記開(kāi)發(fā)中呈現(xiàn)較高的多態(tài)性水平。262份陸地棉種質(zhì)資源的遺傳多樣性分析表明種質(zhì)材料的遺傳基礎(chǔ)較為狹窄。群體結(jié)構(gòu)分析將種質(zhì)材料大致劃分為2個(gè)亞群, 與基于遺傳距離的聚類分析及主成分分析的結(jié)果相同。關(guān)聯(lián)分析檢測(cè)到65個(gè)與纖維品質(zhì)性狀顯著相關(guān)的標(biāo)記位點(diǎn), 其中23個(gè)標(biāo)記位點(diǎn)能在多個(gè)環(huán)境中被重復(fù)穩(wěn)定檢測(cè)到。結(jié)果表明了基于重測(cè)序開(kāi)發(fā)InDel標(biāo)記是有效可行的, 本研究開(kāi)發(fā)的標(biāo)記為進(jìn)一步解析陸地棉纖維品質(zhì)的遺傳機(jī)理提供了理論依據(jù), 為棉花分子標(biāo)記輔助育種提供了高效便捷的基因組工具。

    附表 請(qǐng)見(jiàn)網(wǎng)絡(luò)版: 1) 本刊網(wǎng)站http://zwxb. chinacrops.org/; 2)中國(guó)知網(wǎng)http://www.cnki.net/; 3) 萬(wàn)方數(shù)據(jù)http://c.wanfangdata.com.cn/Periodical-zuowxb. aspx。

    [1] 喻樹(shù)迅. 中國(guó)棉花產(chǎn)業(yè)百年發(fā)展歷程. 農(nóng)學(xué)學(xué)報(bào), 2018, (1): 85–91. Yu S X. The development of cotton production in the recent hundred years of China., 2018, (1): 85–91 (in Chinese with English abstract).

    [2] 董承光, 王娟, 周小鳳, 馬曉梅, 李生秀, 余渝, 李保成. 基于表型性狀的陸地棉種質(zhì)資源遺傳多樣性分析. 植物遺傳資源學(xué)報(bào), 2016, 17: 438–446. Dong C G, Wang J, Zhou X F, Ma X M, Li S X, Yu Y, Li B C. Evaluation on genetic diversity of cotton germplasm resources (L.) on morphological characters., 2016, 17: 438–446 (in Chinese with English abstract).

    [3] Dong C G, Wang J, Yu Y, Li B C, Chen Q J. Association mapping and favourable QTL alleles for fibre quality traits in Upland cotton (L.).,2018, 97: e1–e12.

    [4] Huang C, Shen C, Wen T W, Gao B, Zhu D, Li X, Ahmed M M, Li D, Lin Z X. SSR-based association mapping of fiber quality in upland cotton using an eight-way MAGIC population., 2018, 293: 793–805.

    [5] Sahu P K, Mondal S, Sharma D, Vishwakarma G, Kumar V, Das B K. InDel marker based genetic differentiation and genetic diversity in traditional rice (L.) landraces of Chhattisgarh, India., 2017, 12: e0188864.

    [6] Liu J, Qu J T, Yang C, Tang D G, Li J W, Lan H, Rong T Z. Development of genome-wide insertion and deletion markers for maize, based on next-generation sequencing data., 2015, 16: 601.

    [7] Liu B, Wang Y, Zhai W, Deng J, Wang H, Cui Y, Cheng F, Wang X W, Wu J. Development of InDel markers forbased on whole-genome re-sequencing., 2013, 126: 231–239.

    [8] Zhang T, Hu Y, Jiang W, Fang L, Guan X, Chen J, Zhang J, Saski C A, Scheffler B E, Stelly D M.Sequencing of allotetraploid cotton (L. acc. TM-1) provides a resource for fiber improvement.,2015, 33: 531–537.

    [9] Wang M J, Tu L L, Min L, Lin Z X, Wang P C, Yang Q Y, Ye Z X, Shen C, Li J Y, Zhang X L. Asymmetric subgenome selection and cis-regulatory divergence during cotton domestication., 2017, 49: 579–587.

    [10] Huang C, Nie X H, Shen C, You C Y, Li W, Zhao W X, Zhang X L, Lin Z X. Population structure and genetic basis of the agronomic traits of upland cotton in China revealed by a genome-wide association study using high-density SNPs., 2017, 15: 1374–1386.

    [11] Nie X H, Huang C, You C Y, Li W, Zhao W X, Shen C, Zhang B B, Wang H T, Yan Z H, Dai B S, Wang M J, Zhang X L, Lin Z X. Genome-wide SSR-based association mapping for fiber quality in nation-wide upland cotton inbreed cultivars in China., 2016, 17: 352.

    [12] Rozen S, Skaletsky H. Primer3 on the WWW for general users and for biologist programmers., 2000, 132: 365–386.

    [13] Li X M, Yuan D J, Wang H T, Chen X M, Wang B, Lin Z X, Zhang X L. Increasing cotton genome coverage with polymorphic SSRs as revealed by SSCP., 2012, 55: 459–470.

    [14] Botstein D, White R L, Skolnick M, Davis R W. Construction of a genetic linkage map in man using restriction fragment length polymorphisms., 1980, 32: 314–331.

    [15] Letunic I, Bork P. Interactive tree of life (iTOL) v3: an online tool for the display and annotation of phylogenetic and other trees., 2016, 44: W242–W245.

    [16] Pritchard J K, Stephens M, Donnelly P. Inferences of population structure using multilocus genotype data., 2000, 155: 945–959.

    [17] Earl D A, Vonholdt B M. STRUCTURE HARVESTER: a website and program for visualizing STRUCTURE output and implementing the Evanno method., 2012, 4: 359–361.

    [18] Evanno G, Regnaut S, Goudet J. Detecting the number of clusters of individuals using the software STRUCTURE: a simulation study., 2010, 14: 2611–2620.

    [19] Bradbury P J, Zhang Z, Kroon D E, Casstevens T M, Ramdoss Y, Buckler E S. TASSEL: software for association mapping of complex traits in diverse samples., 2007, 23: 2633–2635.

    [20] Wu D H, Wu H P, Wang C S, Tseng H Y, Hwu K K. Genome-wide InDel marker system for application in rice breeding and mapping studies., 2013, 192: 131–143.

    [21] Wang H, Jin X, Zhang B, Shen C, Lin Z. Enrichment of an intraspecific genetic map of upland cotton by developing markers using parental RAD sequencing., 2015, 22: 147–160.

    [22] Lin Z X, Zhang Y X, Zhang X L, Guo X P. A high-density integrative linkage map for., 2009, 166: 35–45.

    [23] Liu R Z, Wang B H, Guo W Z, Qin Y S, Wang L G, Zhang Y M, Zhang T Z. Quantitative trait loci mapping for yield and its components by using two immortalized populations of a heterotic hybrid inL., 2012, 29: 297–311.

    [24] Fang D D, Li P, Thyssen G, Hinze L L, Percy R G. A microsatellite-based genome-wide analysis of genetic diversity and linkage disequilibrium in Upland cotton (L.)., 2013, 191: 391–401.

    [25] Wen T W, Wu M, Shen C, Gao B, Zhu D, Zhang X L, You C Y, Lin Z X. Linkage and association mapping reveals the genetic basis of brown fibre ()., 2018, 16: 1654–1666.

    [26] Ai X, Liang Y, Wang J, Zheng J, Gong Z, Guo J, Li X, Qu Y. Genetic diversity and structure of elite cotton germplasm (L.) using genome-wide SNP data., 2017, 145: 409–416.

    [27] Tyagi P, Gore M A, Bowman D T, Campbell B T, Udall J A, Kuraparthy V. Genetic diversity and population structure in the US Upland cotton (L.)., 2014, 127: 283–295.

    [28] Flint-Garcia S A, Thornsberry J M, Buckler E S. Structure of linkage disequilibrium in plants., 2003, 54: 357–374.

    [29] Zhang Z W, Ersoz E, Lai C Q, Todhunter R J, Tiwari H K, Gore M A, Bradbury P J, Yu J M, Arnett D K, Ordovas J M, Buckler E S. Mixed linear model approach adapted for genome-wide association studies., 2010, 4: 355–360.

    [30] Ma Z Y, He S P, Wang X F, Sun J L, Zhang Y, Zhang G Y, Wu L Q, Li Z K, Liu Z H, Sun G F, Du X M. Resequencing a core collection of upland cotton identifies genomic variation and loci influencing fiber quality and yield., 2018, 50: 803–813.

    [31] Zhang S W, Feng L C, Xing L T, Yang B, Gao X, Zhu X F, Zhang T Z, Zhou B L. New QTLs for lint percentage and boll weight mined in introgression lines from two feral landraces intoacc TM-1., 2016, 135: 90–101.

    [32] Fang L, Wang Q, Hu Y, Jia Y H, Chen J D, Liu B L, Zhang Z Y, Guan X Y, Chen S Q, Zhou B L, Du X M. Genomic analyses in cotton identify signatures of selection and loci associated with fiber quality and yield traits., 2017, 49: 1089–1098.

    [33] Islam M S, Thyssen G N, Jenkins J N, Zeng L, Delhom C D, McCarty J C, Deng D D, Hinchliffe D J, Jones D C, Fang D D. A MAGIC population-based genome-wide association study reveals functional association ofgene with superior fiber quality in cotton., 2016, 17: 903.

    [34] Sun Z W, Wang X F, Liu Z W, Gu Q S, Zhang Y, Li Z K, Ke H F, Yang J, Wu J H, Wu L Q, Zhang G Y, Zhang C Y, Ma Z Y. Genome-wide association study discovered genetic variation and candidate genes of fibre quality traits inL., 2017, 15: 982–996.

    [35] Liu Z H, Zhu L, Shi H Y, Chen Y, Zhang J M, Zheng Y, Li X B. Cotton GASL genes encoding putative gibberellin- regulated proteins are involved in response to GA signaling in fiber development., 2013, 40: 4561–4570.

    Development and evaluation of InDel markers in cotton based on whole-genome re-sequencing data

    WU Mi, WANG Nian, SHEN Chao, HUANG Cong, WEN Tian-Wang, and LIN Zhong-Xu*

    National Key Laboratory of Crop Genetic Improvement / College of Plant Science and Technology, Huazhong Agricultural University, Wuhan 430070, Hubei, China

    Insertion and deletion (InDel) are abundant forms of genetic variation in the genome. InDel has been recognized as an ideal source for marker development due to its high-density distribution and genotyping efficiency. In this study, the whole genome re-sequencing data of 262 upland cotton accessions were applied to identify 3206 InDel markers, and 320 markers with uniform distribution across the genome were selected to be evaluated. Eighty-seven polymorphic markers were identified, accounting for 26.88% of screened markers. A total of 160 allelic loci were detected using the 87 polymorphic markers in the 262 upland cotton accessions with an average polymorphic information content (PIC) of 0.3073 (ranging from 0.0836 to 0.3750) and an average genetic diversity of 0.3876 (ranging from 0.0874 to 0.5000), indicating a relatively low genetic diversity. Population structure analysis revealed extensive admixture and identified two subgroups, clustering analysis and principal component analysis supported the subgroups identified by STRUCTURE. Association analysis were performed by MLM (Mixed linear model), and 65 marker loci were associated with fiber quality traits (< 0.01), explaining 2.57%–8.12% of the phenotypic variation. Genome-wide and gel based InDel markers developed based on re-sequencing data in this study provide a facile tool for cotton germplasm resources research and molecular marker assisted selection breeding.

    Upland cotton; InDel marker; genetic diversity; population structure; association analysis

    2018-07-19;

    2018-10-08;

    2018-11-16.

    10.3724/SP.J.1006.2019.84100

    林忠旭, E-mail: linzhongxu@mail.hzau.edu.cn

    E-mail: 656542680@qq.com

    本研究由湖北省技術(shù)創(chuàng)新專項(xiàng)(2018ABA082)資助。

    This study was supported by the Technology Innovation Program of Hubei Province (2018ABA082).

    URL: http://kns.cnki.net/kcms/detail/11.1809.S.20181115.1104.002.html

    猜你喜歡
    棉區(qū)種質(zhì)多態(tài)性
    華南地區(qū)最大農(nóng)作物種質(zhì)資源保護(hù)庫(kù)建成
    單核苷酸多態(tài)性與中醫(yī)證候相關(guān)性研究進(jìn)展
    美國(guó)農(nóng)業(yè)部公布各產(chǎn)區(qū)2018年度棉花產(chǎn)量情況
    亞麻抗白粉病種質(zhì)資源的鑒定與篩選
    貴州玉米種質(zhì)資源遺傳多樣性及核心種質(zhì)庫(kù)構(gòu)建
    馬鈴薯cpDNA/mtDNA多態(tài)性的多重PCR檢測(cè)
    紅錐種質(zhì)早期生長(zhǎng)表現(xiàn)
    GlobalFiler~? PCR擴(kuò)增試劑盒驗(yàn)證及其STR遺傳多態(tài)性
    洞庭湖棉區(qū)油茬棉花直播一播全苗技術(shù)
    作物研究(2014年6期)2014-03-01 03:39:12
    淺議兵團(tuán)北疆棉區(qū)種子生產(chǎn)存在的問(wèn)題及解決對(duì)策
    女性被躁到高潮视频| av又黄又爽大尺度在线免费看| 国产在线一区二区三区精| 啦啦啦视频在线资源免费观看| 99香蕉大伊视频| 90打野战视频偷拍视频| 飞空精品影院首页| 一级毛片女人18水好多| 免费高清在线观看视频在线观看| h视频一区二区三区| 在线观看免费高清a一片| 亚洲欧美成人综合另类久久久| 中文字幕另类日韩欧美亚洲嫩草| 欧美精品高潮呻吟av久久| 亚洲午夜精品一区,二区,三区| 伦理电影免费视频| 精品福利永久在线观看| 韩国高清视频一区二区三区| 欧美国产精品一级二级三级| 久久99一区二区三区| 黄片小视频在线播放| 天天躁夜夜躁狠狠躁躁| 久久久久网色| 亚洲中文日韩欧美视频| 久久久久精品国产欧美久久久 | 搡老岳熟女国产| 久热这里只有精品99| av视频免费观看在线观看| 人妻一区二区av| 99久久人妻综合| 一本久久精品| 脱女人内裤的视频| 国产片内射在线| 99热网站在线观看| 精品免费久久久久久久清纯 | 久久久久国产一级毛片高清牌| 99国产精品免费福利视频| 久久精品国产a三级三级三级| av不卡在线播放| 青青草视频在线视频观看| 亚洲黑人精品在线| 大香蕉久久网| 国产精品麻豆人妻色哟哟久久| 欧美 亚洲 国产 日韩一| 自拍欧美九色日韩亚洲蝌蚪91| 欧美日韩亚洲综合一区二区三区_| 久久国产精品影院| av一本久久久久| 国产熟女午夜一区二区三区| 丁香六月天网| 伦理电影免费视频| 国产精品 国内视频| 久久精品熟女亚洲av麻豆精品| 国产日韩欧美视频二区| 国产成人免费无遮挡视频| 免费黄频网站在线观看国产| 999精品在线视频| 日日爽夜夜爽网站| 午夜福利在线观看吧| 国产亚洲av片在线观看秒播厂| 久久国产精品人妻蜜桃| 亚洲一区二区三区欧美精品| 久久久久久久国产电影| 亚洲欧洲精品一区二区精品久久久| 99精品欧美一区二区三区四区| 大片免费播放器 马上看| 国产精品秋霞免费鲁丝片| 人人妻人人澡人人爽人人夜夜| 国产精品麻豆人妻色哟哟久久| 欧美一级毛片孕妇| 交换朋友夫妻互换小说| 日韩熟女老妇一区二区性免费视频| 精品少妇一区二区三区视频日本电影| 777米奇影视久久| 亚洲国产精品一区二区三区在线| 一二三四社区在线视频社区8| 精品一区二区三卡| 蜜桃在线观看..| 手机成人av网站| 欧美大码av| 久久久久久久久免费视频了| 热99re8久久精品国产| 久久九九热精品免费| 免费看十八禁软件| 久久人妻福利社区极品人妻图片| 十八禁高潮呻吟视频| 亚洲成国产人片在线观看| 精品国产一区二区三区久久久樱花| 宅男免费午夜| 午夜福利视频在线观看免费| 亚洲中文字幕日韩| 又黄又粗又硬又大视频| 窝窝影院91人妻| 国产97色在线日韩免费| 国产在线视频一区二区| 久久精品国产亚洲av香蕉五月 | 久久国产精品大桥未久av| 国产高清videossex| 亚洲天堂av无毛| 亚洲精品自拍成人| 国产一级毛片在线| 午夜老司机福利片| av网站在线播放免费| 一进一出抽搐动态| 国产成人啪精品午夜网站| 首页视频小说图片口味搜索| 午夜91福利影院| 三级毛片av免费| 成人手机av| 久久久久国内视频| 久久天堂一区二区三区四区| 久久久欧美国产精品| 亚洲欧美色中文字幕在线| 精品人妻1区二区| 免费在线观看视频国产中文字幕亚洲 | 欧美日韩黄片免| 欧美精品亚洲一区二区| 亚洲精品一卡2卡三卡4卡5卡 | 丁香六月天网| 国产欧美日韩一区二区三区在线| 三上悠亚av全集在线观看| 国产精品一区二区在线不卡| 少妇的丰满在线观看| 国产三级黄色录像| 国产精品久久久av美女十八| 国产野战对白在线观看| 国产一区二区三区综合在线观看| 国产欧美日韩精品亚洲av| 成年人免费黄色播放视频| 人人妻人人澡人人爽人人夜夜| 亚洲人成77777在线视频| 亚洲欧美精品自产自拍| 亚洲欧美一区二区三区黑人| 99国产精品免费福利视频| 欧美老熟妇乱子伦牲交| 两个人免费观看高清视频| 亚洲欧美精品自产自拍| 色婷婷av一区二区三区视频| 国产男女内射视频| 欧美性长视频在线观看| 极品少妇高潮喷水抽搐| 精品一区二区三卡| 日韩三级视频一区二区三区| 狠狠婷婷综合久久久久久88av| 日韩人妻精品一区2区三区| 又大又爽又粗| 亚洲熟女精品中文字幕| 首页视频小说图片口味搜索| 色综合欧美亚洲国产小说| 老司机午夜十八禁免费视频| 一本一本久久a久久精品综合妖精| 男女无遮挡免费网站观看| 午夜福利乱码中文字幕| 黄片小视频在线播放| 美女高潮喷水抽搐中文字幕| 国产不卡av网站在线观看| 欧美午夜高清在线| av在线老鸭窝| 欧美av亚洲av综合av国产av| 男人舔女人的私密视频| 国产欧美日韩一区二区三区在线| 国产欧美日韩一区二区三区在线| 国产高清视频在线播放一区 | 99精国产麻豆久久婷婷| 手机成人av网站| av一本久久久久| 高清在线国产一区| 人成视频在线观看免费观看| av网站免费在线观看视频| 亚洲激情五月婷婷啪啪| 91国产中文字幕| 国产国语露脸激情在线看| 爱豆传媒免费全集在线观看| 青春草亚洲视频在线观看| 亚洲精品久久午夜乱码| 欧美日韩精品网址| 少妇粗大呻吟视频| 久久狼人影院| 一二三四社区在线视频社区8| 啪啪无遮挡十八禁网站| 波多野结衣一区麻豆| 老司机福利观看| 久久久久精品国产欧美久久久 | 黄色毛片三级朝国网站| 老司机亚洲免费影院| 97在线人人人人妻| 欧美激情 高清一区二区三区| 人人妻人人澡人人爽人人夜夜| cao死你这个sao货| 亚洲熟女精品中文字幕| 午夜激情久久久久久久| 日日夜夜操网爽| 亚洲欧美精品综合一区二区三区| 精品久久久精品久久久| 亚洲va日本ⅴa欧美va伊人久久 | 日韩精品免费视频一区二区三区| 精品人妻熟女毛片av久久网站| 美女高潮到喷水免费观看| 麻豆国产av国片精品| 亚洲人成电影免费在线| 老司机午夜十八禁免费视频| 大片电影免费在线观看免费| 久久久久精品国产欧美久久久 | 岛国在线观看网站| 少妇 在线观看| 18禁观看日本| 午夜影院在线不卡| 亚洲视频免费观看视频| 欧美午夜高清在线| 肉色欧美久久久久久久蜜桃| 午夜激情久久久久久久| 看免费av毛片| a级毛片黄视频| 另类精品久久| 99国产极品粉嫩在线观看| 欧美老熟妇乱子伦牲交| 中文字幕人妻丝袜制服| 成年美女黄网站色视频大全免费| 91成年电影在线观看| 精品人妻熟女毛片av久久网站| 啦啦啦 在线观看视频| videosex国产| 欧美老熟妇乱子伦牲交| 成在线人永久免费视频| 久久性视频一级片| 中文精品一卡2卡3卡4更新| 亚洲全国av大片| 男女边摸边吃奶| 亚洲黑人精品在线| 香蕉国产在线看| 精品一区在线观看国产| 老司机亚洲免费影院| av线在线观看网站| 久久精品熟女亚洲av麻豆精品| 久久av网站| 男女免费视频国产| 亚洲伊人久久精品综合| 脱女人内裤的视频| 老司机亚洲免费影院| 成人国产一区最新在线观看| 国产一区二区三区综合在线观看| 超色免费av| av网站在线播放免费| 午夜91福利影院| 亚洲欧美日韩另类电影网站| 久久99热这里只频精品6学生| 美女福利国产在线| 一二三四在线观看免费中文在| 亚洲精品国产区一区二| 在线亚洲精品国产二区图片欧美| 丝袜喷水一区| 亚洲精品一二三| 亚洲一码二码三码区别大吗| 免费高清在线观看视频在线观看| 午夜成年电影在线免费观看| 午夜福利影视在线免费观看| 女人被躁到高潮嗷嗷叫费观| 18禁国产床啪视频网站| 成年美女黄网站色视频大全免费| 国产区一区二久久| 真人做人爱边吃奶动态| 男人爽女人下面视频在线观看| 久久人人97超碰香蕉20202| 日韩精品免费视频一区二区三区| 日韩三级视频一区二区三区| 国产av又大| 久久午夜综合久久蜜桃| 国产精品99久久99久久久不卡| 一个人免费看片子| 久久中文字幕一级| 交换朋友夫妻互换小说| 熟女少妇亚洲综合色aaa.| 91老司机精品| 久久影院123| 精品高清国产在线一区| xxxhd国产人妻xxx| 窝窝影院91人妻| 黄色毛片三级朝国网站| av免费在线观看网站| 少妇精品久久久久久久| 日韩视频一区二区在线观看| 欧美日韩一级在线毛片| 制服人妻中文乱码| 精品福利观看| 男女床上黄色一级片免费看| 亚洲精品自拍成人| 美女扒开内裤让男人捅视频| 老司机深夜福利视频在线观看 | 叶爱在线成人免费视频播放| 亚洲人成电影观看| 男女国产视频网站| 国产在视频线精品| 一个人免费看片子| 国内毛片毛片毛片毛片毛片| 亚洲人成电影免费在线| 50天的宝宝边吃奶边哭怎么回事| 日韩 亚洲 欧美在线| 日韩欧美一区视频在线观看| 精品免费久久久久久久清纯 | 女人被躁到高潮嗷嗷叫费观| 女性生殖器流出的白浆| 男女边摸边吃奶| 国产男女内射视频| 久久精品亚洲熟妇少妇任你| 丝瓜视频免费看黄片| 日本黄色日本黄色录像| 免费人妻精品一区二区三区视频| 男女国产视频网站| 国产精品国产三级国产专区5o| a级片在线免费高清观看视频| 亚洲av国产av综合av卡| 亚洲国产中文字幕在线视频| 在线 av 中文字幕| 免费少妇av软件| 亚洲性夜色夜夜综合| 老司机午夜十八禁免费视频| 午夜91福利影院| 精品亚洲乱码少妇综合久久| 欧美激情极品国产一区二区三区| 国产av又大| 脱女人内裤的视频| 黄色 视频免费看| 久久精品aⅴ一区二区三区四区| 成年美女黄网站色视频大全免费| 国产精品亚洲av一区麻豆| 亚洲精品av麻豆狂野| 欧美日韩黄片免| 国产欧美亚洲国产| 一级,二级,三级黄色视频| 久久午夜综合久久蜜桃| 亚洲中文字幕日韩| 在线精品无人区一区二区三| 最新的欧美精品一区二区| 久久精品国产综合久久久| 精品人妻在线不人妻| 国产国语露脸激情在线看| 老司机深夜福利视频在线观看 | 日本wwww免费看| e午夜精品久久久久久久| 日日夜夜操网爽| 精品人妻一区二区三区麻豆| 欧美日韩中文字幕国产精品一区二区三区 | 久久女婷五月综合色啪小说| 在线观看免费日韩欧美大片| 男女高潮啪啪啪动态图| 国产成人欧美| 啦啦啦 在线观看视频| 性色av乱码一区二区三区2| 一级片'在线观看视频| 亚洲一区中文字幕在线| 午夜免费观看性视频| 一级片免费观看大全| 国产精品免费大片| 又紧又爽又黄一区二区| 欧美激情极品国产一区二区三区| 三级毛片av免费| 50天的宝宝边吃奶边哭怎么回事| 在线观看www视频免费| 国产一卡二卡三卡精品| 久久久久国内视频| 欧美乱码精品一区二区三区| 国产视频一区二区在线看| 在线精品无人区一区二区三| 久久午夜综合久久蜜桃| 黑人巨大精品欧美一区二区mp4| 99re6热这里在线精品视频| 国产成人av教育| 在线观看舔阴道视频| videos熟女内射| 久久人妻熟女aⅴ| 国产亚洲午夜精品一区二区久久| 亚洲精品成人av观看孕妇| 中文字幕人妻丝袜一区二区| av福利片在线| 9191精品国产免费久久| 搡老岳熟女国产| 久久久久国内视频| av国产精品久久久久影院| 黑人巨大精品欧美一区二区蜜桃| 久久女婷五月综合色啪小说| 日韩一卡2卡3卡4卡2021年| 老司机亚洲免费影院| 国产欧美日韩一区二区三 | 久久天躁狠狠躁夜夜2o2o| 精品久久蜜臀av无| 欧美 亚洲 国产 日韩一| h视频一区二区三区| 一区二区三区精品91| 真人做人爱边吃奶动态| 欧美av亚洲av综合av国产av| 一级,二级,三级黄色视频| 黄色怎么调成土黄色| videos熟女内射| 十八禁人妻一区二区| 青青草视频在线视频观看| 三上悠亚av全集在线观看| 90打野战视频偷拍视频| 午夜成年电影在线免费观看| 欧美 日韩 精品 国产| 99精国产麻豆久久婷婷| 777米奇影视久久| 一边摸一边做爽爽视频免费| kizo精华| 国产精品秋霞免费鲁丝片| 欧美精品啪啪一区二区三区 | 中文字幕另类日韩欧美亚洲嫩草| 亚洲国产精品一区三区| 超碰97精品在线观看| 国产精品.久久久| 少妇的丰满在线观看| 国产精品国产三级国产专区5o| 性色av乱码一区二区三区2| www.熟女人妻精品国产| 侵犯人妻中文字幕一二三四区| 日本黄色日本黄色录像| 亚洲免费av在线视频| 男人爽女人下面视频在线观看| 18禁国产床啪视频网站| 一进一出抽搐动态| 美国免费a级毛片| 精品人妻在线不人妻| 亚洲熟女精品中文字幕| 国产在视频线精品| 欧美午夜高清在线| 国产欧美日韩综合在线一区二区| 高清欧美精品videossex| 精品人妻在线不人妻| 亚洲天堂av无毛| 一本久久精品| 久久毛片免费看一区二区三区| 麻豆av在线久日| 国产av国产精品国产| 久久国产精品男人的天堂亚洲| 中国美女看黄片| 少妇猛男粗大的猛烈进出视频| h视频一区二区三区| 最新在线观看一区二区三区| 免费在线观看影片大全网站| 日本a在线网址| 国产99久久九九免费精品| 天天影视国产精品| 一进一出抽搐动态| 黑人巨大精品欧美一区二区mp4| 欧美日韩亚洲综合一区二区三区_| kizo精华| 久久中文看片网| 中文字幕色久视频| 一级a爱视频在线免费观看| 黄色视频不卡| 欧美另类亚洲清纯唯美| 亚洲国产欧美日韩在线播放| 亚洲七黄色美女视频| 淫妇啪啪啪对白视频 | 91麻豆精品激情在线观看国产 | 亚洲国产欧美在线一区| 久久久精品区二区三区| 久久国产精品大桥未久av| www.精华液| 亚洲国产成人一精品久久久| 狠狠狠狠99中文字幕| 黄片播放在线免费| 无遮挡黄片免费观看| 久久性视频一级片| 亚洲欧美清纯卡通| 欧美老熟妇乱子伦牲交| 日本vs欧美在线观看视频| 亚洲五月色婷婷综合| 女人久久www免费人成看片| 成人三级做爰电影| 老司机靠b影院| 纵有疾风起免费观看全集完整版| 久久国产亚洲av麻豆专区| 国产主播在线观看一区二区| 操出白浆在线播放| a级毛片黄视频| 五月开心婷婷网| 国产片内射在线| 99久久人妻综合| 亚洲天堂av无毛| 秋霞在线观看毛片| 日本a在线网址| 亚洲欧美一区二区三区久久| 深夜精品福利| 久久久国产一区二区| 婷婷成人精品国产| 欧美精品一区二区大全| 在线观看免费高清a一片| 亚洲精品美女久久av网站| 一级片'在线观看视频| 午夜免费成人在线视频| 成年人黄色毛片网站| 美女国产高潮福利片在线看| 亚洲欧美日韩另类电影网站| 亚洲av日韩在线播放| 国产在视频线精品| 亚洲av电影在线进入| 久久综合国产亚洲精品| 国产精品香港三级国产av潘金莲| 久久精品久久久久久噜噜老黄| 动漫黄色视频在线观看| 啪啪无遮挡十八禁网站| 精品久久久精品久久久| 亚洲精品粉嫩美女一区| 成人18禁高潮啪啪吃奶动态图| 国产真人三级小视频在线观看| 日韩制服骚丝袜av| 99国产极品粉嫩在线观看| 久久精品久久久久久噜噜老黄| 国产精品99久久99久久久不卡| 久久人妻福利社区极品人妻图片| 日日摸夜夜添夜夜添小说| 精品乱码久久久久久99久播| 国产精品欧美亚洲77777| 国产亚洲av片在线观看秒播厂| 最黄视频免费看| 黄网站色视频无遮挡免费观看| 高清黄色对白视频在线免费看| 在线观看舔阴道视频| 一二三四社区在线视频社区8| 精品一区二区三区av网在线观看 | 久久精品亚洲熟妇少妇任你| 大陆偷拍与自拍| 午夜福利视频精品| 黄色视频不卡| 天堂8中文在线网| 91麻豆精品激情在线观看国产 | 精品第一国产精品| 少妇被粗大的猛进出69影院| 老司机午夜福利在线观看视频 | 免费在线观看黄色视频的| 侵犯人妻中文字幕一二三四区| 黄色视频不卡| 首页视频小说图片口味搜索| 伊人久久大香线蕉亚洲五| 在线观看免费午夜福利视频| 另类精品久久| 在线观看免费日韩欧美大片| 99热全是精品| 啪啪无遮挡十八禁网站| 国产精品熟女久久久久浪| 亚洲九九香蕉| 99国产精品一区二区三区| 69av精品久久久久久 | 欧美日韩国产mv在线观看视频| 十八禁网站网址无遮挡| 一级毛片女人18水好多| 高清视频免费观看一区二区| 欧美另类一区| 两个人免费观看高清视频| 亚洲欧洲日产国产| avwww免费| 老司机福利观看| 亚洲第一av免费看| 在线av久久热| 一个人免费看片子| 首页视频小说图片口味搜索| 新久久久久国产一级毛片| 99久久精品国产亚洲精品| 午夜免费成人在线视频| 国产麻豆69| 成年av动漫网址| 久久久精品国产亚洲av高清涩受| 亚洲 国产 在线| 女警被强在线播放| 韩国精品一区二区三区| 亚洲国产av影院在线观看| 人人妻,人人澡人人爽秒播| 亚洲伊人久久精品综合| 精品人妻一区二区三区麻豆| 两人在一起打扑克的视频| 亚洲伊人色综图| 国产成人av教育| 久久精品国产a三级三级三级| 亚洲av男天堂| 在线观看人妻少妇| 国产精品亚洲av一区麻豆| tube8黄色片| 国产亚洲欧美精品永久| 在线观看一区二区三区激情| 精品欧美一区二区三区在线| 天天影视国产精品| 色综合欧美亚洲国产小说| 精品国产一区二区三区四区第35| 国产激情久久老熟女| 人妻 亚洲 视频| 国产精品一区二区在线不卡| 久久天堂一区二区三区四区| 成人影院久久| 99精品欧美一区二区三区四区| 久久99热这里只频精品6学生| 国产精品.久久久| 在线观看舔阴道视频| 性色av乱码一区二区三区2| 国产精品熟女久久久久浪| 男女免费视频国产| 9热在线视频观看99| 色精品久久人妻99蜜桃| 天堂俺去俺来也www色官网| 国产av又大| 午夜精品国产一区二区电影| bbb黄色大片| 久久精品国产a三级三级三级| 国产97色在线日韩免费| 久久国产精品男人的天堂亚洲| 国产亚洲欧美精品永久| 考比视频在线观看| 亚洲国产精品一区三区| 欧美日本中文国产一区发布| 日日夜夜操网爽| 国产一区二区三区在线臀色熟女 | 不卡一级毛片| 亚洲一卡2卡3卡4卡5卡精品中文| 亚洲全国av大片| 18禁裸乳无遮挡动漫免费视频| 精品一区二区三卡|