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

秋海棠(Begonia grandis)的種內表型多樣性

2021-03-02 03:53:00李雁瓷付乃峰孫加芝曹建國田代科
植物研究 2021年5期
關鍵詞:分類分析

李雁瓷 付乃峰 孫加芝 肖 艷 曹建國 田代科*

(1. 上海師范大學生命科學學院,上海 200234;2. 上海辰山植物園/中國科學院上海辰山植物科學研究中心,上海 201602;3. 上海市資源植物功能基因組學重點實驗室,上海 201602)

秋海棠屬(BegoniaL.)是植物多樣性最豐富的類群之一,為被子植物第六大屬,目前已有1 900多個種[1~2]。秋海棠(Begonia grandisDryander)隸屬于秋海棠科(Begoniaceae)秋海棠屬(Begonia)多年生球莖種類,具有明顯的地上莖,冬季地上部分干枯死亡,地下莖休眠。廣泛分布于我國南北各地,南起云南屏邊,最北可至遼寧凌源,東起浙江寧波,最西可達西藏察隅,也是世界秋海棠科分布最北端、抗寒性最強,分布范圍最廣之一的種類[3]。該種最早的文獻記載可追溯到宋代,其花葉優美、歷史悠久、文化深厚,具有較高的觀賞、藥用和食用價值[4]。

秋海棠的種內多樣性十分豐富,不同居群在植株大小、分枝、葉色等某些性狀上的變異顯著,鑒定困難,從而造成其種內分類比較混亂[5]。以往,基于形態多樣性對該種的種內分類沒有開展系統的比較分析。種內亞種、變種的發表都是建立在簡單的形態觀察上,或者僅僅依據幾份臘葉標本,因此帶有一定的主觀性。《中國植物志》[6]及Flora of China(FOC)[7]有關該種的種內單元劃分均不太合理,在種內關系的處理上存在分歧,依據現有的性狀描述檢索無法準確定到亞種和變種,在鑒定實踐中缺乏操作性,不利于該種的合理利用和信息交流,因此有待進一步地開展定量分析。而數量分類法能對大量的生物學性狀進行全面綜合的分析,并對所有的分類性狀進行等權處理,再根據性狀間的親緣性和相似性進行歸類,從而使植物分類學研究從定性描述水平走向精確、可重復的定量分析水平,擺脫傳統分類的主觀性[8~9]。本研究以39 個秋海棠居群(包括33 個國產野生居群和6 個日本自然化居群)為研究材料,通過數量性狀變異分析、質量性狀遺傳多樣性分析、主成分分析和聚類分析對其33個表型性狀進行多樣性和數量分類研究,旨在初步探討秋海棠的主要分類特征,揭示主要性狀之間的內在聯系及各性狀所占有的信息比例,試圖為秋海棠的種內單元分類及多樣性研究提供客觀可行的依據。

1 材料與方法

1.1 試驗材料

試驗材料取自上海辰山植物園觀賞植物資源與種質創新課題組的秋海棠苗圃,以2012~2016年從全國范圍內采集和從日本引種的秋海棠居群為觀測對象。塑料盆種植(口徑18.5 cm,高15.2 cm),每盆種植一株(球)。栽培基質為腐葉土和珍珠巖和樹皮按5∶1∶1比例混合,盆底墊上一層陶粒以利排水。澆水頻度大體為夏季一天1~2 次,春秋季3~5 天1 次,冬季一周1 次,遵循基質單次澆透、保持潤而不干的原則。大棚鋪設“寒冷紗”進行遮蔭,遮光率約為40%,夏季用風機進行降溫(表1,圖1)。

表1 試驗材料來源及編碼Table1 The origin and code of B.grandis for experiment

1.2 形態觀察及數據測量

根據田代科修訂的《秋海棠野外性狀記錄表》(未發表)對苗圃收集栽培的秋海棠植株的完整生長階段進行觀察、記錄形態數據和拍照。每個居群選8~10 株健康的植株,于盛花期觀測花部等性狀,葉片測量從地上莖第一節間開始,從下往上的10 片成熟葉片中隨機選取6~8 片進行觀測,2018~2019 年連續觀察兩年。株高、葉柄、葉片、雄花被片等數量性狀用直尺測量,每株每個指標測量10次,并充分考慮最大和最小的極端類型,最終取其平均值。其中,葉片測量方法(圖2):葉長為中脈軸平行線兩端點最長距離;葉寬為同中脈平行線垂直最寬距離;葉寬側寬為葉片寬側同中脈平行線垂直最寬處;葉窄側寬為葉片窄側同中脈平行線垂直最寬處;下延葉長為從葉柄著生處開始、平行中脈,向下到寬側葉片底部邊緣頂端為止。質量性狀采用目測法,結合《中國植物志》[6]和FOC[7]中有關秋海棠的植物資料,對植株進行考察后歸納總結該種在質量性狀方面的特征,包括葉斑紋有無、葉脈和葉面顏色、花柄和內外被片顏色等。

1.3 分類性狀的選取和編碼

對原始數據進行初步統計計算后,將39 個供試居群作為分類單元(Operation taxonomic unit,OUT),以其外部形態為依據,篩選出33 個性狀進行后續編碼分析(2),其中對立的二元性狀(Binary character,簡稱為“B”)用1 和2 來表示有無;多態性狀(Multistate character,簡稱為“M”)取連續的正整數1、2、3、4、5…進行編碼;數量性狀(Numerical character,簡稱為“N”)無需編碼,直接取其測量值;導出性狀無須編碼(Derived character,簡稱為“D”),取其計算測量值[10-11]。將所有數據測定后錄入Excel,匯總為原始數據矩陣。

1.4 數據處理方法

數量性狀分析:對于數量性狀進行描述性統計,計算其平均值、標準差、極值、極差、變異系數,其中變異系數VC=標準差/均值[12~13]。

質量性狀分析:對于二元和多態性狀進行形態變異的頻數統計,計算遺傳多樣性指數(Shannon-Wiener diversity index,H′):

式中:pi表示第i種變異類型出現頻率;ln為自然對數;H′值可以描述質量性狀的總體多樣性[14],遺傳多樣性指數的計算由Microsoft Excel 2016完成[15]。

聚類及主成分分析:使用SPSS 26.0 分析數據,首先對原始數值矩陣進行數據標準化(STD)處理,以消除不同量綱對數據分析產生影響;然后Pearson 相關性為度量指標,采用組間聯接對31 個性狀進行R 型聚類分析,做出性狀指標聚類樹形圖;同時利用降維對標準化后數值矩陣進行主成分分析,將多項指標轉化為少數幾項綜合指標[16];最后根據R 型聚類分析和主成分分析結果,對居群進行Q 型聚類分析,采用平方歐式距離為度量標準,運用Ward聚類方法做出分類結果樹狀圖。

2 結果與分析

2.1 數量性狀的變異分析

對秋海棠的25 個數量性狀進行了變異分析,結果見表3。從表3 中可以看出,39 個居群數量性狀變異系數為0.077~0.441,變異系數超過0.3的性狀有10個,分別為雄蕊群高、葉柄粗、葉窄側寬、葉寬、雄花內輪被片寬、葉柄長、下延葉長、株高、雄蕊數、雄蕊柱長度,說明這些性狀在居群間的差異較大,性狀表現不穩定;變異系數低于0.12 的性狀有4個,分別為葉片長寬比、葉寬窄側比、雄花外輪被片長寬比、雄花花冠長寬比,說明秋海棠的葉片形狀、外被片及花冠的形狀在居群間差異較小,遺傳穩定性較高。

表3 秋海棠數量性狀變異分析Table 3 Variation analysis for quantitative traits of B.grandis

2.2 質量性狀多樣性分析

對秋海棠的8 個二元及多態性狀進行頻率分布和多樣性統計,從結果表4 中可以看出:39 個居群數量性狀多樣性指數變化范圍為0.287~1.483。其中,葉斑紋的多樣性指數最低,葉片無斑出現頻率最高達91.67%,說明葉片出現斑紋的情況較少;地上莖分枝多少、雄花花柄、內外輪被片顏色4 個性狀的平均多樣性指數均在0.95 以上,說明秋海棠的地上莖分枝及花部顏色在居群間的多樣性比較豐富。

表4 秋海棠質量性狀頻率分布及多樣性Table 4 Frequency distribution and diversity index on qualitative traits of B.grandis

表2 性狀及編碼Table 2 Taxonomic traits and codes

2.3 主成分分析

主成分分析是將多個性狀指標化為少數指標的一種統計方法,能從眾多性狀因素中找出主要的綜合性狀因素,使它們盡可能多地保留原始變量的信息,且彼此不相關[17]。對33 個性狀進行主成分分析,求得各主成分的特征值、方差貢獻率和累計貢獻率,結果顯示見表5。第1 主成分的特征值為11.29,方差貢獻率最大為34.22%,是最重要的主成分,其中,雄花外輪被片長和寬、雄蕊數、雄花花冠長、雄蕊群高和雄蕊群直徑的特征向量絕對值較大,說明第1主成分基本反映了這些指標的信息。第2 主成分特征值為4.26,方差貢獻率為12.91%,主要由下延葉長、葉寬側寬、葉長、葉窄側寬、葉寬等性狀決定。第3 主成分特征值為2.85,方差貢獻率為8.64%,主要由腹脈顏色、雄花外輪被片顏色、葉脈顏色、葉背顏色等性狀決定。前3個主成分累計貢獻率為55.77%,是秋海棠分類最主要的區分性狀。第4 主成分方差貢獻率為6.21%,特征向量絕對值最大的性狀為下延葉長同葉長比、雄花內輪被片長寬比和株高。第5主成分貢獻率為5.74%,雄花外輪被片長寬比、葉寬窄側比是特征向量絕對值較大的性狀。第6 主成分貢獻率為4.57%,特向量絕對值較大的性狀為葉長寬比、葉斑紋有無。一般情況下,累計貢獻率達到70%以上,事物的基本面貌可以得到反映[17]。前6個主成分的累積方差貢獻率為72.28%,包含了39個秋海棠居群的大部分性狀信息,貢獻率較大的形態學性狀是造成秋海棠居群之間表型差異的主要因素,由此可知,秋海棠居群間的差異主要表現為花大小、葉大小、葉色、花色和株高的不同。

表5 主成分的特征值,貢獻率及特征向量Table 5 The eigenvalue,contribution rate and eigenvector of selected principal components

2.4 聚類分析

2.4.1 R型聚類分析

數量分類中性狀的選取直接影響著分類的結果,因而有必要對所選性狀間的關系進行研究[18]。R 型聚類是對研究對象的性狀指標進行分類,可以反映出各性狀之間的關系及性狀選取的合理性[19]。對秋海棠的33 個性狀指標進行R 型聚類,結果從樹狀圖(圖3)可知,在D=3處做閾值等級分界線,所選性狀大多數是獨立的,可用于下一步的Q 型聚類分析,但部分性狀存在關聯性很大的現象。其中,葉長、葉寬、葉寬側寬、葉窄側寬和下延葉長這5個性狀,雄花花冠寬、雄花內被片長、雄花花冠長、雄花外輪被片長和雄花外輪被片寬這5個性狀的相關性較大。因為主成分分析結果表明葉形和花形大小的有關性狀是區分秋海棠居群的重要指標,所以都應予以保留。從R 型聚類的結果看,所選的33 個性狀可以做到彼此獨立影響居群的變異和進化,不會互相造成不必要的干擾。

2.4.2 Q型聚類分析

Q 型聚類是以OTU 為分類基本單位進行的分類運算,可以分析各OUT 間的親緣關系。對39 個秋海棠居群進行Q 型聚類分析,結果由樹狀圖(圖4)可知:在等級結合線L1(D=25)處,39 個秋海棠居群被分為2 大類(圖5),第一類包括湖北隨州(HBSZ)、山東泰山(SDTS)、湖北宜昌(HBYC)、天津薊縣(TJJX)、湖北秭歸(HBZG)、湖南張家界(HNZJJ)、河南焦作(HENJZ2、HENJZ1)、河北興隆(HBXL)、遼寧凌源(LNLY)、湖北紅安(HBHA)、湖北岳西(HBYX)、湖南桑植(HNSZ1、HNSZ2)、湖南龍山(HNLS1、HNLS2、HNLS3)河北井陘(HEBJJ)、山西靈石(SXLS)、湖北神農架(HBSNJ)、北京香山(BJXS)、湖南衡山(HNHS)、浙江寧波(ZJNB)、重慶金佛山(CQJFS)、湖南永順(HNYS1、HNYS2、HNYS3)21 個居群,主要位于長江以北和華東地區,它們的共同特征是雄花較小,其花冠長寬、內外被片長寬、雄蕊群高與直徑、雄蕊數基本上都小于平均值;第二類包括湖南常德(HNCD)、云南昆明(YNKM)、福建永泰(FJYT)、四川雅安(SCYA)、湖北利川(HBLC)、云南麗江(YNLJ)、安徽休寧(AHXN)、福 建 武 夷 山(FJWYS)、福 建 泰 寧(FJTN)、江西廬山(JXLS)、四川峨眉山(SCEMS)、浙江臨安(ZJLA1、ZJLA2)、日本富山(RBFS)、日本鐮 倉(RBLC1、RBLC2)、日 本 大 分(RBDF1、RBDF2)、日本群馬(RBQM)、日本開成町(RBKCT)、日本富山(RBFS)、18 個居群,主要位于西南、華東地區和日本,它們的共同特征是雄花較大,其花冠長寬、內外被片長寬、雄蕊群高與直徑、雄蕊數基本上都大于平均值。

在等級結合線L2(D=20)處又可將分出的大類分為4組(圖6~7)。第一組包括:湖北隨州(HBSZ)、山東泰山(SDTS)、湖北宜昌(HBYC)、天津薊縣(TJJX)、湖 北 秭 歸(HBZG)、湖 南 張 家 界(HNZJJ)、河南焦作(HENJZ2、HENJZ1)、河北興隆(HBXL)、遼寧凌源(LNLY)、湖北紅安(HBHA)、湖北岳西(HBYX)、湖南桑植(HNSZ1、HNSZ2)、湖南龍山(HNLS1、HNLS2、HNLS3)河北井陘(HEBJJ)、山西靈石(SXLS)、湖北神農架(HBSNJ)、北京香山(BJXS)、湖南衡山(HNHS)、浙江寧波(ZJNB)、重慶金佛山(CQJFS),主要位于長江以北和華中地區,它們的共同特征是植株較矮、分枝較多、葉片較小,其株高、葉長、葉寬基本上都小于平均值,葉表面沒有斑紋,大部分葉片的腹脈、葉背及背脈顏色為綠色。第二組僅湖南永順的一個特殊居群(HNYS1、HNYS2、HNYS3),位于華中地區,它們的共同特征是植株較高、分枝少、葉片較大,其株高、葉長、葉寬基本上都大于平均值,葉面有白斑,葉片的腹脈、葉背及背脈顏色皆為紅色。第三組包括:湖南常德(HNCD)、云南昆明(YNKM)、福建永泰(FJYT)、四川雅安(SCYA)、湖北利川(HBLC)、云南麗江(YNLJ)、安徽休寧(AHXN)、福建武夷山(FJWYS)、福建泰寧(FJTN)、江西廬山(JXLS)、四川峨眉山(SCEMS),主要位于西南和華東地區,它們的共同特征是植株偏矮、分枝較少、葉片偏小,其株高、葉長、葉寬大部分小于平均值,葉表面沒有斑紋,葉片的腹脈、葉背及背脈顏色皆為紅色。第四組包括:浙江臨安(ZJLA1、ZJLA2)、日本富山(RBFS)、日本鐮倉RBLC1、RBLC2)、日本大分(RBDF1、RBDF2)、日本群馬(RBQM)、日本開成町(RBKCT)、日本富山(RBFS),主要位于華東地區和日本,它們的共同特征是植株較高、分枝較少、葉片較大,其株高、葉長、葉寬基本上都大于平均值,為4組之最,葉表面沒有斑紋,大部分葉片的腹脈、葉背及背脈顏色為綠色。

3 討論

3.1 秋海棠表型性狀多樣性

變異系數在某種程度上反映出群體數量性狀的差異性和多樣性,變異系數(CV 值)大則分化程度高、多樣性顯著,相反則群體相對穩定、一致性較高[20]。選取39 個秋海棠居群的33 個形態學性狀,對其中的25個數量性狀進行了變異分析,結果顯示數量性狀變異系數為0.077~0.441,平均變異系數為0.26,表明這些性狀在秋海棠居群間變異較大,表型多樣性水平較為豐富。同時,對8 個質量性狀進行了頻率分布和多樣性統計,除葉斑紋外,其余7 個質量性狀的平均多樣性指數均較高(>0.657),表明秋海棠的地上莖分枝、花被片顏色及葉色在居群間的多樣性比較豐富。此外,前期性狀觀測中發現同一居群同時存在葉表面無斑、有白斑,葉背綠色、背脈紅色、葉背紅色,以及少數開白花個體類型,說明秋海棠的表型在居群內的多樣性也很豐富。

3.2 秋海棠分類性狀選取

通過對性狀指標的主成分分析,共提取了6個主成分。前3 個主成分累計貢獻率為55.77%,直到第6 個主成分累計貢獻率才達到72.28%,表明各性狀的貢獻率比較分散,累計貢獻率增長比較緩慢,說明在秋海棠的演化過程中性狀變異具有多樣性[21]。根據各主成分的貢獻率選出了影響較大的性狀指標,綜合而言,花大小、葉大小、葉色、花色、植株大小是秋海棠分類最主要的標準和依據。

R 型聚類分析表明,除葉片和花大小指標有較明顯的相關性外,其余性狀相對獨立,說明性狀選取較為合理。但是,由于前期記錄的數據不全,雌花(特別是被片數目、柱頭形態)和果實(特別是果實及果翅形態)的性狀因為數據采集不夠而未納入統計分析,因此,結果無法更全面地揭示各性狀間關系。在以后研究中,應進一步補充其他相關重要性狀。此外,從湖北宜昌引種的一個居群屬于軟枝、近匍匐型,花小,綜合性狀看應獨立成枝,也曾被分類學家列為不同變種,但由于栽培管理不當材料死亡,遺憾缺乏相關數據統計,今后有待補充研究。

3.3 秋海棠種內單位劃分

秋海棠的種內分類單元的分類學處理一直存在爭議。《中國植物志》[6]根據植株高矮、葉片大小及毛被情況、雄蕊柱長度和雄蕊數目、柱頭形態、果翅形態等差異將秋海棠處理為1 原亞種、2 亞種、3 變種;FOC[7]根據葉片大小及顏色、雄蕊柱長度雄蕊數目、柱頭形態的差異將秋海棠處理為1原亞種、2 亞種;李行娟[22]用三個葉綠體片段聯合分析將44 個秋海棠居群劃分為三個支系;本研究根據花大小、葉大小、植株大小、葉色、葉斑紋情況將39 個秋海棠居群分為四組,可分成4 個變種。在《中國植物志》和FOC 中,產云南(麗江)的全柱秋海棠(B.grandissubsp.holostylaIrmscher)因花柱分離、不分枝、柱頭呈頭狀或腎狀、雄蕊柱長度短于1 mm 的特征區別于其他亞種,而在本研究前期數據統計中,云南麗江(YNLJ)居群雄蕊柱長度有長達3 mm 的,又因花大小、葉片大小、植株大小、葉色等主要性狀與其他居群較為相似聚在一起,并沒有獨立成為一個亞種,同以前的分類相矛盾。研究還發現,日本的幾個居群同浙江臨安天目山的居群聚類成一枝,說明其親緣關系最近,結合形態觀察和以往的基因大小比較分析,進一步說明,日本的自然化居群的祖先可能來自浙江天目山(田代科此前如此預測過)。

造成以上秋海棠種內單位分類差異的原因有很多方面:一是《中國植物志》和FOC 依據形態學的定性描述對秋海棠進行種內劃分,觀察樣本較少會造成對不重要分類指標的錯誤側重,從而導致人為主觀設立的性狀邊界。而數量分類法則是將大量形態學數據進行定量分析,能避免過多次要性狀對分類的干擾。二是李行娟所用的葉綠體片段聚類關注的是親緣關系,而形態學性狀的聚類關注的是表型性狀,而表現型和基因型之間受到基因表達調控、個體發育調控、生長環境等因素的影響,導致某些表觀特征不能完全將遺傳變異反映出來[23~24];三是秋海棠居群的性狀多樣性比較豐富,居群間很可能發生了基因流,存在一些性狀過渡類群,因此導致亞種和變種間的界限界定困難,給種內分類帶來了挑戰。如同一居群的植株葉片顏色有腹脈綠色、腹脈紅色和葉背綠色、葉背紅色四種個體類型,葉表面有無斑和白斑兩種類型,還有一些過渡形態,這說明依據葉色和葉斑紋作為分類標準是不合適的,這同野外觀察的實際情況一致。 但是,湖南永順的一個居群所有個體葉片都具有穩定遺傳的白斑,葉片大、植株高但分枝少,即使通過種子播種繁殖的后代也同母本性狀一致,性狀十分穩定,獨立成為一個變種是否比較合理,尚需結合分子證據探討。綜上,《中國植物志》將秋海棠處理為3 個亞種及3 個變種和FOC 秋海棠處理為3 個亞種的處理方式均有待商榷。我們的研究結果和考慮到沒有納入分析的湖北宜昌居群,秋海棠可能分成5個變種更適合。然而,單純依靠少量形態學性狀不能完全反映植物的本質規律,該種居群的種內分類需要綜合考慮大量表型性狀和分子標記分析結果,今后還需要建立在大量、更全面的采樣基礎上,結合適當的分子方法開展綜合分析才能真正解決問題。

此外,在實際調查過程中,因為多人參與數據記錄,發現部分質量性狀的記錄存在一定的主觀性,使得在判斷和測量時存在誤差。例如葉毛被疏密情況、葉緣齒類型等性狀可能存在不同人不同標準的情況,故本次研究沒有納入這些性狀指標,但從實際觀察來看,這些指標的貢獻率不會太高,因此對結果影響不大。在數量分類中,分類性狀選取不同會造成不同結果,而本研究沒有將雌花和果實納入分析中,可能影響結果的準確性。因此,今后研究時需要進一步完善形態性狀鑒定標準和方法,數據記錄應盡可能由一人為主完成,同時使用一些生理指標,如植株萌發出土時間及居群的始花期、盛花期、末花期,并使用多種分類方法進行綜合分析以及從群體演化的角度進行系統分類討論。

猜你喜歡
分類分析
分類算一算
垃圾分類的困惑你有嗎
大眾健康(2021年6期)2021-06-08 19:30:06
隱蔽失效適航要求符合性驗證分析
分類討論求坐標
電力系統不平衡分析
電子制作(2018年18期)2018-11-14 01:48:24
數據分析中的分類討論
教你一招:數的分類
電力系統及其自動化發展趨勢分析
給塑料分分類吧
中西醫結合治療抑郁癥100例分析
主站蜘蛛池模板: 欧美人在线一区二区三区| 国产人成在线观看| 国产午夜小视频| 毛片久久网站小视频| 国产毛片片精品天天看视频| 亚洲欧美成人在线视频| 国产在线无码av完整版在线观看| 欧洲熟妇精品视频| 国产不卡国语在线| 午夜国产在线观看| 国产一级毛片网站| 久久无码高潮喷水| 欧美一区二区人人喊爽| 欧美翘臀一区二区三区| av在线5g无码天天| 亚洲无码37.| 亚洲欧洲日韩综合色天使| 国产一级在线播放| 日韩人妻无码制服丝袜视频| 99精品免费在线| 久久午夜夜伦鲁鲁片不卡| 久久一本日韩精品中文字幕屁孩| 国产精品极品美女自在线网站| 真实国产乱子伦视频| 国产又色又爽又黄| 亚洲丝袜第一页| 特级aaaaaaaaa毛片免费视频| 91免费国产高清观看| 中文天堂在线视频| 欧美高清国产| 在线观看国产精品日本不卡网| 亚洲区第一页| 日韩av手机在线| 青青操视频免费观看| 日韩高清欧美| 久久男人视频| 国产激情无码一区二区免费| 亚洲无码日韩一区| 囯产av无码片毛片一级| 亚洲天堂网视频| 欧美日韩一区二区在线播放| 欧美日韩免费| 国产乱子伦手机在线| 91香蕉视频下载网站| 麻豆精品在线| 国产综合在线观看视频| 在线免费亚洲无码视频| 黄色网页在线观看| 午夜成人在线视频| 国产后式a一视频| 污网站免费在线观看| 色成人亚洲| 日本一区二区不卡视频| 爽爽影院十八禁在线观看| jijzzizz老师出水喷水喷出| 激情综合激情| 久久国产拍爱| 欧美日韩国产综合视频在线观看| 成人在线天堂| 91精品国产自产在线观看| 中文字幕永久在线观看| 精品欧美日韩国产日漫一区不卡| 91麻豆精品国产91久久久久| 久久久久青草线综合超碰| 四虎国产永久在线观看| 国产欧美日韩综合在线第一| 91娇喘视频| 久久99精品久久久久纯品| 色欲色欲久久综合网| 99久久性生片| 欧美日韩另类在线| 亚洲黄网在线| 激情综合图区| 高清无码一本到东京热| 国产理论精品| 日韩欧美高清视频| 国产一区成人| 欧美狠狠干| 色偷偷av男人的天堂不卡| 精品亚洲麻豆1区2区3区| 中国国产高清免费AV片| 成人午夜网址|