999精品在线视频,手机成人午夜在线视频,久久不卡国产精品无码,中日无码在线观看,成人av手机在线观看,日韩精品亚洲一区中文字幕,亚洲av无码人妻,四虎国产在线观看 ?

基于流式細胞術和基因組Survey的繡球基因組大小及特征分析

2021-08-02 11:53:51陳雙雙齊香玉馮景王華娣鄧衍明
江蘇農業科學 2021年12期

陳雙雙 齊香玉 馮景 王華娣 鄧衍明

摘要:繡球花序飽滿、花型多樣、花色豐富,是重要的觀賞植物。很多繡球品種種植在酸性含鋁的土壤中時,花序可由紅色變為藍色、紫色或藍紫色,使其更受消費者喜愛。但繡球基因組信息缺乏,基因資源開發和分子機制研究相對滯后。利用流式細胞術和基于K-mer分析的基因組Survey技術估測繡球主栽品種“Bailer”(商品名:Endless SummerTM,中文名譯為無盡夏)的基因組大小。流式細胞術分析結果表明,以豌豆(Pisum sativum L.)為對照的流式細胞術估測“Bailer”基因組大小為2.17 pg(約2.12 Gb);Illumina HiSeq PE測序產生212.41 Gb的有效數據,K-mer 分析結果顯示,修正后的基因組大小為1.96 Gb;根據K-mer分析,該基因組具有明顯的雜合峰,雜合度和重復率分別為1.97%、72.80%。此外,經組裝,該基因Contigs總數為7 418 139(N50為249 bp),總讀長為1.49 Gb,Scaffolds 總數為7 071 975(N50為289 bp),總讀長為1.52 Gb,GC含量為38.75%。綜上所述,推測繡球“Bailer”基因組大小約為1.96~2.12 Gb,屬雜合度較高的復雜基因組。結果可為繡球全基因組測序及相關研究提供了理論依據。

關鍵詞:繡球;基因組大小;流式細胞術;K-mer分析;Survey技術

中圖分類號:S685.990.1 文獻標志碼: A

文章編號:1002-1302(2021)12-0039-06

收稿日期:2020-10-10

基金項目:國家自然科學基金(編號:31901359);江蘇省農業科學院探索性顛覆性創新計劃[編號:ZX(17)2004];江蘇省昆山市科技項目(編號:KN1809)。

作者簡介:陳雙雙(1986—),女,河南安陽人,博士,助理研究員,主要從事觀賞植物遺傳育種與分子生物學研究。E-mail:chenshuang7876@sina.com。

通信作者:鄧衍明,博士,研究員,主要從事花卉遺傳育種與分子生物學研究。E-mail:nksdym@163.com。

繡球(Hydrangea macrophylla)別稱八仙花,屬于虎耳草科繡球屬,因其花序色彩艷麗、開花期長、適應性廣和抗逆性強而被廣泛種植[1]。繡球是觀賞植物中的強耐鋁(Al)植物,其成熟葉片和萼片可以積累高于3 000 μg/g Al而不呈現任何毒性癥狀[2]。部分繡球品種花色易變,調節栽培土壤pH值及Al3+含量能使其呈現紅色、紫色或藍色等不同顏色[3]。大量研究表明,萼片中飛燕草素-3-O-葡萄糖苷與Al3+的絡合物調控萼片顏色,而顏色的變化取決于萼片中Al3+含量[3-4]。因此,繡球不僅是非常值得開發的觀賞植物資源,而且是研究觀賞植物藍色花形成機制的特色材料,在花卉育種領域具有較大的應用價值。

長期以來,對繡球的研究多集中在栽培和育種領域[5-6];近年來,人們開展了繡球分子生物學特性相關研究,如分子標記的開發[7]、花色變化的分子機制[8-10]、鋁轉運相關基因克隆與功能研究[11-13]等。基因組特征研究是植物基因資源開發和分子機制研究的前提,然而,關于繡球不同品種基因組大小及特征研究的報道尚少,無法為分子生物學研究提供足夠的理論基礎。早在1957年,繡球體細胞染色體數首次被鑒定為2n=2x=36[14],并在Funamoto等的研究中得到證實[15];此后,Cerbah等也發現,多數繡球種和亞種的染色體數為36,為二倍體(2n=2x=36)[16]。繡球不同品種具有不同大小的基因組,大花繡球[H. macrophylla ssp. macrophylla (Thunb.) Ser.]和粗齒繡球(H. macrophylla ssp. serrata Ser.)的1C值分別為2.15、1.92 pg,相應地GC含量分別高達41.0%、40.6%[16];進一步研究發現,大花繡球二倍體品種的1C值在1.98~2.49 pg之間,但也有研究表明,大花繡球Taube 1C值為(3.75±0.05) pg,而Veitchii為(2.35±0.05) pg[17-19]。這些結果均是基于流式細胞術估測的結果。

流式細胞術和基于K-mer分析的基因組Survey測序是估測基因組大小的常用方法,前者用于評估基因組大小和倍性水平,操作簡單,分辨率高,準確性高,后者具有高通量測序、速度快和數據量大等優點[20-21]。本研究采用流式細胞術和基于K-mer分析的基因組Survey測序對大花繡球主栽品種“Bailer”基因組大小進行估測和相互驗證,并對測序質量、雜合度、序列重復率等進行評估測定,旨在為繡球全基因組測序提供重要依據,為繡球進化分析、基因挖掘、遺傳改良以及結構基因組研究奠定基礎。

1 材料與方法

1.1 試驗材料

供試材料大花繡球“Bailer”(商品名:Endless SummerTM,中文名譯為無盡夏)種植于江蘇省農業科學院繡球種質資源圃內,該品種具有新老枝都能開花和花色極易改變等優點,是繡球的主栽品種。豌豆(Pisum sativum L.,2C=9.09 pg)作為對照[22],種子由江蘇省農業科學院江蘇省農業種質資源保護與利用平臺提供。試驗于2019年1月至3月在江蘇省農業科學院院內進行。

1.2 試驗方法

1.2.1 流式細胞分析 采集“Bailer”和豌豆幼嫩新鮮葉片,分別稱取3份,每份約1 g,置于盛有細胞裂解液LB01[15 mmol/L Tris,2 mmol/L Na2乙二胺四乙酸(EDTA),0.5 mmol/L 四鹽酸精胺,80 mmol/L KCl,20 mmol/L NaCl,0.1% Triton X-100,pH值為7.0~8.0]的培養皿中,并迅速用鋒利的刀片將其切碎。將培養皿中的解離液過濾至1.5 mL離心管中,1 000 r/min 離心 5 min 后棄上清,收集細胞沉淀,加入100 μL解離液懸浮。將待測樣品“Bailer”與對照樣品豌豆細胞懸液等比例混合,加入熒光染料PI(碘化丙啶),于4 ℃冰箱放置15 min后進行上機檢測[20]。

使用美國BD公司的流式細胞儀對“Bailer”基因組大小進行評估,變異系數控制在5%以內。每個樣品進行3次重復試驗。測試數據利用儀器自帶軟件FACSTM Software 1.0.0.650進行分析,數據統計分析利用SPSS 17.0完成。繡球基因組(核DNA含量)估算方法如下:待測樣品核DNA含量=對照樣本核DNA含量×[(待測樣品G0/G1 峰熒光強度)/(對照樣品G0/G1峰熒光強度)][23]。

1.2.2 樣品DNA制備與基因組測序 采集“Bailer”幼嫩葉片,立即用液氮速凍后保存于 -80 ℃。使用天根植物基因組DNA提取試劑盒進行DNA提取,Nanodrop 2000分光光度計(德國,eppendorf)和1%瓊脂糖凝膠電泳檢測DNA濃度及質量。檢測合格的DNA樣品通過Covaris超聲波破碎儀隨機打斷成長度為350 bp的片段,經末端修復、加A尾、加測序接頭、純化、PCR擴增等步驟完成整個文庫制備,構建好的文庫通過Illumina Hiseq PE進行測序。獲得的原始數據經過濾后利用Vsearch軟件(https://github.com/torognes/vsearch)進行評估,包括數據產量、測序錯誤率、Q20含量、Q30含量、GC含量等。同時,利用Blast工具將部分有效測序數據(clean data)與NCBI核苷酸數據庫(NT)進行比對,選取比對結果中的前3個,分析樣本是否存在外源污染。

1.2.3 K-mer分析估算基因組大小、雜合度及重復率 采用基于K-mer的分析方法來估計基因組大小和雜合率等,即從一段連續序列中迭代地選取長度為K個堿基的序列,取K=17進行分析。假設K-mer的深度頻率服從泊松分布,并且從堿基讀數中獲取的所有K-mer都可以遍歷整個基因組,則可以計算K-mer頻數分布和K-mer深度分布曲線。根據曲線獲得K-mer深度估計值,用于估計基因組大小。通過K-mer曲線的分布情況,估計基因組特征,計算雜合峰值與純合峰值的比值,從而得到雜合率;利用純合峰1.8倍后的K-mer個數所占的比例來估計重復序列比例[24]。

1.2.4 基因組組裝 利用SOAPdenovo(https://github.com/aquaskyline/SOAPdenovo2)對測序數據進行拼接組裝,選擇K-mer=41組裝Contigs和Scaffolds。將讀取長度小于100 bp的低質量序列刪除,并將過濾后的有效數據組裝為Contigs,最后利用reads之間的連接關系和插入片段大小信息,將Contigs組裝成Scaffolds[24]。

2 結果與分析

2.1 基于流式細胞術繡球基因組大小分析

當變異系數控制在5%以內時,待測樣品與對照樣品的離子清晰集中,并且可以同時檢測到混合樣品的顆粒簇(圖1),因此,解離液LB01和對照樣品豌豆對于估算繡球基因組比較適宜,保證了利用豌豆作為對照估算繡球基因組大小的準確性。利用豌豆(2C=9.09 pg)作為對照樣品,根據混合樣品PI熒光強度以及待測樣品與對照樣品峰值的倍數關系,計算繡球基因組1C值為2.17 pg(圖2、表1)。

2.2 基因組Survey與質量評估

利用Illumina Hiseq PE測序平臺對構建的3個插入片段為350 bp的DNA文庫(DES01984、DES02239和DES02238)進行雙末端測序,經過濾和修正后,獲得212.41 Gb的測序數據。測序質量評估顯示,這3個文庫Q20值分別為96.04%、96.64%、96.99%,Q30值分別為90.8%、92.12%、92.83%(表2),測序錯誤率均小于0.03%,表明測序質量較好。此外,取每個文庫的部分有效測序數據比對到NCBI核苷酸數據庫(NT)進行比對,在DES01984、DES02239和DES02238文庫中,排名前3的物種分別為粗齒繡球/珙桐(Davidia involucrata Baill.)/大花繡球、大花繡球/粗齒繡球/山茱萸(Cornus officinalis Sieb. et Zucc.)、大花繡球/粗齒繡球/葡萄(Vitis vinifera L.)。根據這些結果,認為序列比對結果為同源比對,樣本不存在污染現象(表3)。

2.3 K-mer分析

將獲得的有效數據(212.41 Gb)用于“Bailer”K-mer分析,并估測基因組大小、雜合度及序列重復率。從圖3可以看出,橫坐標表示K-mer出現的次數,縱坐標表示出現的頻率;K-mer分布曲線成峰情況較好,呈雙峰分布,在50和85附近各有1個峰值。根據survey分析,推測depth=85附近為主峰值,即K-mer的期望深度,K-mer總數為168、564、374、328(圖3、表4)。結合公式(基因組大小=K-mer-number/depth)估算得出“Bailer”基因組為1 983.11 Mb,修正后的基因組大小為 1 964.82 Mb(約1.96 Gb,表4)。從圖3還可以看出,在主峰的1/2位置處有明顯的雜合峰,而且K-mer分布曲線有嚴重的拖尾,這可能是“Bailer”基因組雜合度和重復序列含量均較高所致。根據公式計算,Bailer基因組雜合度和重復率分別高達1.97%、72.8%。

2.4 基因組組裝

利用Soap denovo軟件對所獲得的有效數據(clean data)進行denovo拼接組裝,選擇K-mer=41得到最優的拼接效果。組裝結果顯示,Contigs 總數為7 418 139,N50長度為249 bp,總讀長為 1.49 Gb;進一步組裝后,Scaffolds 總數為7 071 975,N50和N90長度分別為289 bp和113 bp,總讀長為1.52 Gb(表5)。

2.5 GC含量分布

根據GC含量與測序深度相關性分析,Bailer的GC含量為38.75%(圖4)。GC-depth分析顯示,測序無明顯偏向趨勢;GC深度分布被分為2層,這可能與基因組雜合度高有關。此外, 取低深度區域的序列與NT數據庫比對,未比對上細菌序列,因此可以排除此區域受細菌污染的可能性。

3 討論與結論

植物基因組大小估測不僅對于物種本身基因組測序具有重要意義,而且可為相關蛋白組、轉錄組以及代謝產物研究提供不可或缺的基礎資料[25-26]。本研究采用流式細胞術和基于K-mer分析的基因組Survey測序估測繡球主栽品種Bailer的基因組大小,為進一步進行基因組測序和相關基因功能研究提供了更準確的信息。

自1983年以來,流式細胞術已成為分析植物基因組大小的標準方法,應用于多種植物基因組大小(C值)的預測[16,24,27-29]。本研究中預測Bailer的1C值為2.17 pg,這與已有報道大花繡球1C值(2.15 pg)接近[16]。但與Jones等所報道的結果不同,該研究中大花繡球二倍體品種的1C值在2.23~2.49 pg之間,其中Bailmer(即Bailer)1C值為 2.37 pg[17]。流式細胞術中采用材料類型、熒光染料、裂解液的化學成分、內標的選擇、處理環境等均有可能影響同一物種基因組大小的測定結果[13]。結合已有報道,造成繡球基因組大小差異的原因可能是材料類型和熒光染料不同造成的。本研究利用PI作為熒光染料估測Bailer基因組大小,PI可通過簡單的化學計量插入到雙鏈核酸中,無堿基偏好。Cerbah等研究中對繡球1C值估測采用的是溴化乙錠(與PI類似,對細胞計數具有相同特異性)、AT堿基特異性染料Hoechst和GC堿基特異性染料mithramycin作為細胞核熒光染料,但是該研究中沒有明確所測品種[16];而Jones等研究采用的是具有AT堿基高度親和力的染料4′,6′-diamidino-2-phenylindole (DPAI)[17]。

本研究還利用了基于Illumina Hiseq測序平臺的K-mer分析技術進一步驗證并估測繡球基因組大小。根據3個不同文庫的測序數據,使用K-mer分析,修正后的Bailer的基因組大小約為 1 964.82 Mb(1.96 Gb)。根據細胞核DNA含量與基因組大小的關系(1 pg=978 Mb),流式細胞術估測繡球基因組大小為2.17 pg,換算為 2 122.26 Mb,即2.12 Gb。顯然,以上2種方法估測的繡球基因組大小不一致,流式細胞術比K-mer分析估測結果大約157.44 Mb。流式細胞術估測基因組較大可能是由于估測結果是參考已知基因組大小的對照樣品,而基于K-mer分析的基因組信息是基于數學計算,可能更加全面、準確[20,24,27]。

據報道,蘭科蝴蝶蘭(Phalaenopsis equestris)基因組雜合度為1.2%[28],康乃馨(Dianthus caryophyllus L.)基因組雜合度為0.2%[30],鐵皮石斛(Dendrobium officinale)基因組雜合度為0.48%[31],而本研究基因組雜合度評估顯示“Bailer”基因組雜合度則高達1.97%;Bailer基因組序列重復率估測為72.8%,與黃石斛(Dendrobium catenatum Lindl.)接近(78.1%)[32],然而,蝴蝶蘭、鐵皮石斛和牽牛(Ipomoea nil)基因組僅含有62.00%、63.30%和59.74%的重復序列,康乃馨基因組序列僅有33%的重復率。這些結果可能與基因組大小相關,如黃石斛基因組大小為1.1 Gb[33],蝴蝶蘭基因組大小為1.2 Gb[32],鐵皮石斛基因組大小為1.3 Gb[31],而康乃馨基因組大小只有Bailer的2/3,為622 Mb。因此,可以推測基因組大小與基因組雜合度和重復率相關。此外,與其他觀賞植物如木槿(Hibiscus syriacus,34.04%)[34]、蝴蝶蘭(Phalaenopsis cultivar ‘KHM190,30.7%)[28]、康乃馨(36%)[30]相比,Bailer基因組還有較高的GC含量(38.75%)且GC深度分布被分為2層,可能是由基因組的高雜合率引起。根據測序數據組裝得到的基因組中,Contigs、Scaffolds總長度分別為1.49、1.52 Gb,其中N50長度分別為249、289 bp,基因組雜合度高可能導致N50值較低[22,31]。顯然,繡球基因組是雜合度較高的復雜基因組。

本研究首次利用流式細胞術和K-mer分析的方法測定了繡球基因組大小,并對相應的參數進行了初步評估,主要結論如下:(1)大花繡球Bailer基因組大小估測為1.96~2.12 Gb;(2)繡球基因組具有較高的雜合度(1.97%)和重復率(72.8%),屬于復雜基因組范疇;(3)由于繡球基因組雜合率較高,導致有較高的GC含量(38.75%)。本研究結果可為繡球全基因組測序及相關研究提供理論依據。

參考文獻:

[1]鄧衍明,韓 勇,齊香玉,等. 繡球屬植物種質資源分析及其花色可調性和葉斑病抗性比較[J]. 植物資源與環境學報,2018,27(4):90-100.

[2]Ma J F,Ryan P R,Delhaize E. Aluminium tolerance in plants and the complexing role of organic acids[J]. Trends in Plant Science,2001,6(6):273-278.

[3]Schreiber H D,Jones A H,Lariviere C M,et al. Role of aluminum in red-to-blue color changes in Hydrangea macrophylla sepals[J]. Biometals,2011,24(6):1005-1015.

[4]Schreiber H D,Swink A M,Godsey T D. The chemical mechanism for Al3+ complexing with delphinidin:a model for the bluing of hydrangea sepals[J]. Journal of Inorganic Biochemistry,2010,104(7):732-739.

[5]孫曉波,蘇家樂,陳雙雙,等. 大花繡球‘無盡夏組培苗葉片再生植株的研究[J]. 中國農學通報,2020,36(16):67-72.

[6]吳可鵬,徐煒杰,謝麗麗,等. 繡球屬植物引種及適應性試驗[J]. 浙江農業科學,2020,61 (3):549-551.

[7]賈新平,孫曉波,梁麗建,等. 繡球SSR-PCR反應體系的建立與優化[J]. 華北農學報,2016,31(4):68-73.

[8]Negishi T,Oshima K,Hattori M,et al. Tonoplast-and plasma membrane-localized aquaporin-family transporters in blue Hydrangea sepals of aluminum hyperaccumulating plant[J]. PLoS One,2012,7(8):e43189.

[9]Negishi T,Oshima K,Hattori M,et al. Plasma membrane-localized Al-transporter from blue Hydrangea sepals is a member of the anion permease family[J]. Genes to Cells,2013,18(5):341-352.

[10]龔仲幸,何 勇,楊 靜,等. 外源硫酸鋁調節八仙花花青苷組成和含量變化的分子生物學機制[J]. 植物營養與肥料學報,2017,23(3):821-826.

[11]Chen H X,Lu C P,Jiang H,et al. Global transcriptome analysis reveals distinct aluminum-tolerance pathways in the Al-accumulating species Hydrangea macrophylla and marker identification[J]. PLoS One,2015,10(12):e0144927.

[12]陳海霞,范 適,胡春梅. 鋁脅迫下八仙花HmVALT基因表達的差異[J]. 西北農林科技大學學報(自然科學版),2017,45(5):147-152.

[13]李志奇,陳海霞. 八仙花NRAT1基因的克隆及表達分析[J]. 福建農業學報,2019,34(6):646-651.

[14]Ishih A,Miyase T,Suzuki T,et al. Seasonal variation in the content of a febrifugine and isofebrifugine alkaloid mixture in aerial parts of Hydrangea macrophylla var. Otaksa,with special reference to its antimalarial activity[J]. Journal of Natural Medicines,2006,61(2):213-216.

[15]Funamoto T,Tanaka R. Karyomorphological studies in some taxa of Hydrangea from Japan[J]. Kromosomo,1988(2):1583-1594.

[16]Cerbah M,Mortreau E,Brown S,et al. Genome size variation and species relationships in the genus Hydrangea[J]. Theoretical and Applied Genetics,2001,103(1):45-51.

[17]Jones K D,Reed S M. Analysis of ploidy level and its effects on guard cell length,pollen diameter,and fertility in Hydrangea macrophylla[J]. HortScience,2007,42(3):483-488.

[18]Rinehart T A,Scheffler B E,Reed S M. Ploidy variation and genetic diversity in Dichroa[J]. HortScience,2010,45(2):208-213.

[19]Peter H,Annette H,Conny T. Molecular reconstruction of an old pedigree of diploid and triploid Hydrangea macrophylla genotypes[J]. Frontiers in Plant Science,2018,9:429-442.

[20]Doleel J,Greilhuber J,Suda J. Estimation of nuclear DNA content in plants using flow cytometry[J]. Nature Protocols,2007,2(9):2233-2244.

[21]Yagi M. Recent progress in whole genome sequencing,high-density linkage maps,and genomic databases of ornamental plants[J]. Breeding Science,2018,68(1):62-70.

[22]Doleel J,Barto J. Plant DNA flow cytometry and estimation of nuclear genome size[J]. Annals of Botany,2005,95(1):99-110.

[23]弓 娜,田新民,周香艷,等. 流式細胞術在植物學研究中的應用——檢測植物核DNA含量和倍性水平[J]. 中國農學通報,2011,27(9):21-27.

[24]Wang C R,Yan H D,Li J,et al. Genome survey sequencing of purple elephant grass (Pennisetum purpureum Schum ‘Zise) and identification of its SSR markers[J]. Molecular Breeding,2018,38(7):94-104.

[25]Bennett M D,Leitch I J. Nuclear DNA amounts in angiosperms:targets,trends and tomorrow[J]. Annals of Botany,2011,107(3):467-590.

[26]Bourge M,Brown S C,Siljak-yakovlev S. Flow cytometry as tool in plant sciences,with emphasis on genome size and ploidy level assessment[J]. Genetics and Applications,2018,2(2):1-12.

[27]霍愷森,趙冬蘭,陳艷麗,等. 甘薯屬耐鹽植物馬鞍藤基因組大小及特征分析[J]. 植物遺傳資源學報,2019,20(3):728-735.

[28]Cai J,Liu X,Vanneste K,et al. The genome sequence of the orchid Phalaenopsis equestris[J]. Nature Genetics,2014,47(1):65-72.

[29]齊香玉,陳雙雙,馮 景,等. 3種瓣型茉莉基因組大小測定與比較[J]. 江蘇農業科學,2020,48(19):40-44.

[30]Yagi M,Kosugi S,Hirakawa H,et al. Sequence analysis of the genome of carnation (Dianthus caryophyllus L.)[J]. DNA Research,2014,21(3):231-241.

[31]Yan L,Wang X,Liu H,et al. The genome of Dendrobium officinale illuminates the biology of the important traditional chinese orchid herb[J]. Molecular Plant,2015,8(6):922-934.

[32]Zhang G Q,Xu Q,Bian C,et al. The Dendrobium catenatum Lindl. genome sequence provides insights into polysaccharide synthase,floral development and adaptive evolution[J]. Scientific Reports,2016,6(1):1-10.

[33]Huang J Z,Lin C P,Cheng T C,et al. The genome and transcriptome of Phalaenopsis yield insights into floral organ development and flowering regulation[J]. PeerJ,2016,4:e2017.

[34]Kim Y M,Kim S,Koo N,et al. Genome analysis of Hibiscus syriacus provides insights of polyploidization and indeterminate flowering in woody plants[J]. DNA Research,2017,24(1):71-80.

主站蜘蛛池模板: 色噜噜狠狠狠综合曰曰曰| 国产91九色在线播放| 91精品免费高清在线| 谁有在线观看日韩亚洲最新视频| 亚洲国产天堂久久综合226114| 草逼视频国产| 美女被狂躁www在线观看| 99在线小视频| 99re66精品视频在线观看| 国产毛片片精品天天看视频| julia中文字幕久久亚洲| 香蕉在线视频网站| 福利片91| 97影院午夜在线观看视频| 亚洲色图在线观看| 国产精品欧美亚洲韩国日本不卡| 中文字幕欧美日韩高清| 欧美精品v欧洲精品| 亚洲av日韩av制服丝袜| 26uuu国产精品视频| 欧美综合区自拍亚洲综合绿色 | 亚洲男人天堂2018| 亚洲福利网址| 好紧太爽了视频免费无码| 无码又爽又刺激的高潮视频| 老司机久久99久久精品播放| 国产欧美中文字幕| 亚洲永久视频| 久久99热这里只有精品免费看 | 国产成人无码播放| 97久久免费视频| 日韩欧美色综合| 日韩天堂网| 国产精品一区二区不卡的视频| 国产视频自拍一区| 免费99精品国产自在现线| 亚洲无码A视频在线| 美女国内精品自产拍在线播放 | 色哟哟色院91精品网站| 色婷婷亚洲十月十月色天| 国产无码制服丝袜| 亚洲国产一区在线观看| 免费国产一级 片内射老| 色综合久久无码网| 丝袜亚洲综合| 午夜啪啪网| 国产丝袜无码精品| 东京热一区二区三区无码视频| 亚洲欧美日韩成人在线| 亚洲综合色吧| 国产精品免费福利久久播放| 天天躁夜夜躁狠狠躁图片| 亚洲精品福利网站| 日韩人妻无码制服丝袜视频| 她的性爱视频| 韩国福利一区| 国产不卡一级毛片视频| 久久黄色影院| 青草精品视频| 欧美黄网在线| 久久久久久久久久国产精品| 亚洲日本在线免费观看| 久久精品最新免费国产成人| 亚洲国产成人精品一二区| 国产精品无码影视久久久久久久| 久久青青草原亚洲av无码| 夜夜爽免费视频| 日韩性网站| www亚洲精品| 亚洲高清资源| 国产丝袜无码一区二区视频| 日本精品视频| 国产网友愉拍精品视频| 高清乱码精品福利在线视频| 青青青国产视频手机| 成人午夜天| 女同久久精品国产99国| 永久免费av网站可以直接看的| 亚洲制服中文字幕一区二区| 成人午夜视频网站| 免费一级毛片在线观看| 精品一区国产精品|