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

基于無人機多光譜影像的不同施氮量水稻LAI反演方法研究

2024-12-31 00:00:00曾廣泉馬韜張孟希戴妍陳凱丁繼輝俞雙恩王中文
江蘇農業科學 2024年20期
關鍵詞:無人機水稻

doi:10.15889/j.issn.1002-1302.2024.20.006

摘要:為探明水稻生育期內葉面積指數(LAI)的變化情況,建立可快速準確估測不同生育期水稻LAI的模型。在蒸滲測坑內進行不同施氮量下的水稻栽培試驗,基于無人機采集不同時期水稻測坑多光譜數據,對計算得出的植被指數進行相關性分析,篩選出與5個生育期相關性最高的前5種植被指數,利用多元線性回歸、偏最小二乘回歸、隨機森林、貝葉斯嶺、梯度提升回歸這5種方法構建預測模型;模型構建后,引入施氮量作為變量對模型進行優化,比較各生育期最優的預測模型。結果表明,機器學習算法下,隨機森林回歸模型對揚花期的水稻LAI預測精度和穩定性最好(r2=0.85,MSE=0.33);施氮量對LAI有顯著影響,引入施氮量作為變量對各時期、各模型的預測精度都有所提高(r2平均值提高0.15)。從整體看,機器學習算法下,隨機森林模型可對各時期的水稻LAI指數進行較好的預測;在分蘗后期,各模型的預測精度和穩定性最好。期待本研究建立的預測模型可為水稻LAI遙感監測和田間精細管理提供參考,并為精準農業定量化研究提供技術支持。

關鍵詞:水稻;無人機;施氮量;葉面積指數;多光譜影像;模型計算

中圖分類號:S127" 文獻標志碼:A

文章編號:1002-1302(2024)20-0041-08

收稿日期:2023-10-23

基金項目:國家自然科學基金(編號:51879074、52109051);江蘇省水利科技項目(編號:2020048)。

作者簡介:曾廣泉(1998—),男,甘肅蘭州人,碩士,主要從事節水灌溉與作物信息方面的研究。E-mail:zengqhalo@qq.com。

通信作者:張孟希,副教授,高級工程師,主要從事水土資源高效利用研究。E-mail:221310030004@hhu.edu.cn。

水稻是我國重要的糧食作物,2022年水稻產量達到2.08億t,占我國糧食總產量的32.92%[1]。水稻生產關系到全國糧食安全[2],了解水稻生長狀況對于實現稻米的高產穩產至關重要。葉面積指數(LAI)是單位面積土地上植株單側葉片平鋪面積的總和,與作物的色素含量、碳循環、生物量、物候等生化參數密切相關,是表現作物生長狀況的的重要指標[3-4]。了解水稻葉面積指數對于掌握水稻生長狀況與產量估算有重要意義。

傳統的LAI測量方法需要進行田間破壞性取樣和人工測量分析,存在工作效率低、工作量大等弊端,無法滿足大片區域的自動化植物表型分析需求[5]。為進行水稻LAI的反演,于強等結合干物質重量和Logistic模型對LAI變化進行研究[6];葉宏寶等通過建立葉齡模型對LAI進行模擬[7]。隨著遙感影像和地理信息系統技術的發展,無人機遙感已成為農業監測中常用的技術手段之一;該技術具備機動靈活、作業周期短、受地形條件影響小等特點[8]。近年來,無人機遙感估測已被運用于各類作物的指標反演中[9],不僅減少了工作量,還提高了指標測量的效率和精度。通過衛星、無人機遙感等技術,王福民等進行了植被指數對LAI估算的研究[10];杭艷紅等結合光譜與紋理特征對水稻葉面積進行估算[11]。在引入了機器學習、神經網絡算法等工具后,孫玉婷等研究了支持向量機對水稻葉面積的預測[12];蘇中濱等改進QGA-ELM算法應用于水稻LAI反演[13]。但以上研究在對不同施氮量下的LAI進行反演時得出的精度不夠高。

在水稻生育期內,各時期LAI值的變化會呈現“先升高后回落”的現象[14-15],這是因為水稻在進入孕穗期后,干物質分配由葉、莖等營養器官轉向穗等生殖器官,且在進入黃熟期后,葉片逐漸變黃、變干,導致LAI在整個生育期達到峰值后便開始回落。在水稻的生長過程中,氮元素有至關重要的作用,合理增施氮肥,可使葉片面積、形態生長達到最適合水稻生長和營養吸收的狀態[16]。在研究中加入施氮量,對LAI預測的影響因素具有重要意義。

結合以上分析,本研究使用無人機拍攝得到多光譜影像,經過后續計算處理得到光譜指數,運用多元回歸、機器學習算法等方法構建模型,并通過交叉驗證和加入施氮量作為變量優化、修正模型,得到適用性更高、精度更為準確的模型,旨在為無人機遙感監測在智慧農業中的應用提供一定的基礎與參考。

1" 材料與方法

1.1" 試驗區概況

試驗于2022年6—10月在江蘇省南京市河海大學江寧校區節水園區(118°47′E,31°55′N)進行。該地屬亞熱帶季風性氣候,降水充沛,年平均降水量約1 090.4 mm,年平均氣溫約15.4℃。

試驗區內設有32個固定式蒸滲測坑,每個測坑長2.5 m、寬2 m、深2 m。測坑內土壤質地為黏壤土,0~30 cm深度體積質量為1.46 g/cm3,飽和土壤含水量為30.4%,田間持水量為26.5%,總孔隙度為45.8%。

1.2" 試驗設計

試驗水稻品種為南粳9108,種植密度為 34萬穴/hm2,每穴3 株籽苗。于6月2日育秧,6月26日移栽,10月21號收割。試驗設置5種氮肥梯度處理N1、N2、N3、N4、N5,對應施氮量分別為0、150、225、300、375 kg/hm2,按照基肥 ∶蘗肥 ∶穗肥=4 ∶3 ∶3的配比施用氮肥(以純氮計算)。基肥、蘗肥、穗肥的施用時間分別為6月25日、7月6日、8月18日。磷肥(以P2O5計算)按照75 kg/hm2在施基肥時施用,鉀肥(以K2O計算)按照120 kg/hm2分基肥和穗肥各1/2施用。每種施氮處理設重復3次。

1.3" 數據采集

1.3.1" 葉面積指數數據采集

LAI數據采集于2022年水稻分蘗前期(7月11日)、分蘗中期(7月23日)、分蘗后期(8月7日)、拔節期(8月17日)、揚花期(9月2日)、乳熟期(9月16日)、黃熟期(10月3日)。數據采集通過破壞性取樣方法,在各試驗小區內隨機選取4穴,每穴選取5張代表性葉片,使用EPSON 4990 Photo掃描儀進行掃描,使用FIJI軟件分析測定葉片面積,結合比葉重法[17]測得LAI。

1.3.2" 圖像數據采集

圖像數據由大疆精靈4多光譜版無人機拍攝獲取。于2022年7—9月中午11:00—13:00,選擇無云、光照條件好的天氣進行拍攝。無人機飛行高度設定為8 m,飛行旁向重疊率為80%,航向重疊率為70%。

拍攝用無人機配備6個1/2.9英寸互補金屬氧化物半導體傳感器:其中1個彩色傳感器用于可見光成像;5個單色傳感器用于多光譜成像,分別獲取藍、綠、紅、紅邊、近紅外波段。單個傳感器的有效像素為208萬px,多光譜傳感器的具體參數見表1。

1.4" 植被指數計算

植被指數是用來衡量地表植被狀況的一種指標[18],它通常采用遙感技術,通過測量所接收到的反射輻射能力,計算出植被與非植被區域的光譜反射率差異,從而得到相應的植被指數值;各植被指數能夠在一定條件下用來定量表明植被的生長狀況。本研究選取20種植被指數(表2)進行分析,使用皮爾遜相關系數選取相關性較好的植被指數參與模型構建。

1.5" 模型構建與評價

選取水稻分蘗后期、拔節期、揚花期、乳熟期、黃熟期的試驗數據構建模型,選取2種多元線性回歸方法及3種機器學習方法建立模型,分別是多元線性回歸模型、偏最小二乘回歸模型、隨機森林模型、貝葉斯嶺模型、梯度提升回歸模型。訓練機器學習相關模型時,將訓練中70%的試驗數據作為訓練集、30%的數據作為驗證集,并進行5折交叉驗證[37]。

模型構建完成后,采用決定系數(r2)和均方根誤差(RMSE)評價模型。一般來說,r2越高,表明模型越穩定;RMSE越低,表示模型精度越高。

1.6" 數據處理與分析

數據統計分析采用Microsoft Excel 2021、IBM SPSS 25,繪圖采用Origin 2022。顯著性分析采用Duncans新復極差法,相關性分析采用Person相關系數法。

2" 結果與分析

2.1" 不同施氮量下水稻LAI的變化規律

由圖1表明,不同施氮量下的水稻LAI在整個生育期內總體呈先升后降的趨勢。在結束返青期進入分蘗期時,施氮量對水稻LAI的影響并不顯著;在進入分蘗后期后,施氮量對水稻LAI的影響逐漸增大。整個生育期內,LAI值從分蘗初期開始不斷增加,在揚花期達到峰值,進入乳熟期后下降。這是由于進入乳熟期后,葉片光合作用減弱,養分向穗部轉移,植株下部的葉片開始衰老死亡,從而造成LAI逐漸減小[38]。

其中,N5處理的LAI值在分蘗前期的均值為0.75,到揚花期逐漸增加到6.47,到黃熟期后回落到4.32。在各生育期內,LAI值隨施氮量的增加而增大;其中N2處理相比N1處理的增幅最為顯著,N5相比N4的增幅較小。揚花期,N2處理的LAI均值比N1處理增加了1.64,而N5處理的LAI均值比N4處理僅增加0.37。

2.2" 不同生育期植被指數與LAI的相關性

選取水稻分蘗后期、拔節期、揚花期、乳熟期、黃熟期相應的LAI,使用選取的20種植被指數與其

進行Person相關性分析,得到圖2所示結果。圖2顯示,在分蘗后期,WDRVI、NDVI、OSAVI這3項指數與LAI的相關性最好,相關系數值(r值)分別為0.932、0.929、0.927;在拔節期,RVI、CIg、CIre這3項指數與LAI的相關性最好,r值分別為0.874、0.873、0.868;在揚花期,CVI、CIg、CIre這3項指數與LAI的相關性最好,r值分別為0.836、0.835、0.833;在乳熟期,TVI、DVI、EVI這3項指數與LAI的相關性最好,r值分別為0.718、0.715、0.711;在黃熟期,CIg、CVI、GNDVI這3項指數與LAI的相關性最好,r值分別為0.720、0.709、0.707。由此說明,不同生育期下植被指數與LAI的相關性不同;在分析不同生育期的LAI時,可分別選取相關性較好的植被指數進行分析, 用優選的植被指數構建水稻的

LAI反演模型。

根據圖2得出的結果,選取各生育期與LAI相關性最好的5個植被指數,具體情況見表3。

2.3" 基于優選植被指數的水稻LAI模型構建

2.3.1" 基于植被指數的模型構建

根據各時期所選植被指數情況,對各生育期進行植被指數與LAI的模型構建。分別將各時期的相關植被指數放入多元回歸模型(MLR)、偏最小二乘回歸模型(PLS)、隨機森林模型(RF)、貝葉斯嶺模型(BR)、梯度提升回歸模型(GBR)中進行訓練。由表4可知,在各模型的訓練結果中,對水稻LAI值的模擬預測以分蘗后期到揚花期之間較好(r2值均在0.6~0.9之間);其中最好的是RF模型對揚花期的模擬結果,r2=0.85,RMSE=0.33。

2.3.2" 模型的改進與驗證

基于前人的研究[39]以及“2.1”節中對本試驗區水稻做出的研究結論,可得出施氮量對水稻的LAI變化有顯著影響,在建立模型時可加入施氮量這一自變量對模型進行改進。由圖3結果可知,在加入施氮量這一自變量后,各生育期模型的精度均有不同程度的提升。梯度提升回歸模型中,5個生育期的r2值由0.78、0.61、0.86、0.63、0.56分別提高到0.88、0.83、0.95、0.79、0.69,各時期的r2值平均提高了0.15。

從5種預測模型在不同生育期的預測結果來看(圖3),大部分模型的預測值和實測值之間存在著較好的線性關系。其中,黃熟期在不同模型下5個時期中的預測效果最差,其5種預測模型的r2值分別為0.52、0.50、0.52、0.56、0.56,均低于各模型r2的平均值(0.68);改進后黃熟期模型的r2值分別為0.51、0.51、0.58、0.64、0.69,均低于改進后各模型r2的平均值(0.74)。

3" 討論與結論

本研究中,計算整理了不同時期與水稻LAI相

關性最好的植被指數,將以上植被指數帶入多種預

測模型中對LAI進行預測,得到各時期最適宜的模擬結果;加入施氮量這一控制量對模型進行改進,使各模型的精度均有不同程度的提高。

在前人的研究中,多運用光譜指數以及由光譜指數進行線性運算得出的植被指數進行建模;但無人機獲取的光譜信息較容易受到日光條件的干擾,而杭艷紅等的研究中也發現僅用光譜指數與植被指數對LAI進行估計存在一定的局限性[11]。同時,氮元素在水稻的生長過程中有至關重要的作用,其對水稻的LAI變化也有顯著影響[40-41]。研究發現,在分蘗后期到揚花期這段葉片生長最繁盛的時期,施氮量對水稻LAI的大小有顯著影響。因此,本研究結合了植被指數和施氮量來提高水稻LAI的估算準確度。

在本研究中,綜合不同時期所選擇的植被指數來看,CVI、CIg、CIre等與近紅外波段關聯較高的植被指數與LAI的相關性最好,相關性r值分別達到

了0.87、0.86、0.83。結合植株的生理特性,植株葉片構造中的柵狀組織和冠層葉片會對紅外太陽輻射進行散射,在近紅外波段形成高反射平臺,從而使得水稻LAI指數與近紅外波段的光譜反射率表現出良好的相關性[42]。這一特性也可很好地解釋同一模型在水稻不同生長時期得出不同預測準確度情況。研究發現,水稻在分蘗期到揚花期葉片不斷生長的這一階段內,各模型的預測精度(r2均值為0.76)要好于進入乳熟期后,植株養分向穗部轉移、植株下部的葉片開始衰老枯黃的這一階段(r2均值為0.56)。

在本研究中,通過比較各時期建立的預測模型,可看出機器學習算法下的隨機森林模型、貝葉斯嶺回歸模型、梯度提升模型的穩定性、預測精度均好于多元線性模型、偏最小二乘回歸模型。這與水稻的生長特性相關。在分蘗期后,水稻葉片的生長速度逐漸變緩,從而導致LAI指數與植被指數的變化為非線性關系[28]。因此,在構建水稻LAI預測模型時,可以考慮采用更多非線性模型進行預測。但機器學習算法下的模型存在黑箱模型的特性[43],使得模型不易理解和解釋。故計劃在后續研究中把不同品種、不同播種方式之間的差異性作為考慮因素,以更充分地了解植被指數與水稻LAI之間的關系。

綜上所述,得出以下結論:(1)各生育期內植被指數與水稻LAI的相關性普遍較強,但不同生育期內的與LAI相關性最好的植被指數類型不同,其中CVI、CIg、CIre這3項植被指數與LAI的相關性最好。

(2)各模型的橫向對比中,機器學習算法下的隨機森林模型對揚花期的預測精度和穩定性最好;各時期的橫向對比中,各模型對分蘗后期的預測精度和穩定性最好。

(3)施氮量對水稻LAI有顯著影響,加入施氮量這一變量對模型進行優化,可使各時期的模型預測精度有所提高。

參考文獻:

[1]國家統計局. 國家統計局關于2022年糧食產量數據的公告[J]. 現代面粉工業,2023,37(1):18.

[2]章秀福,王丹英,方福平,等. 中國糧食安全和水稻生產[J]. 農業現代化研究,2005,26(2):85-88.

[3]常好雪,蔡曉斌,陳曉玲,等. 基于實測光譜的植被指數對水稻葉面積指數的響應特征分析[J]. 光譜學與光譜分析,2018,38(1):205-211.

[4]Gitelson A A,Via A,Arkebauer T J,et al. Remote estimation of leaf area index and green leaf biomass in maize canopies[J]. Geophysical Research Letters,2003,30(5):1248.

[5]韓煥豪,崔遠來,時元智,等. SunScan冠層分析儀在水稻葉面積指數測量中的應用[J]. 灌溉排水學報,2015,34(8):44-48.

[6]于" 強,傅抱璞,姚克敏.水稻葉面積指數的普適增長模型[J]. 中國農業氣象,1995(2):6-8.

[7]葉宏寶,孟亞利,湯" 亮,等. 水稻葉齡與葉面積指數動態的模擬研究[J]. 中國水稻科學,2008,22(6):625-630.

[8]劉" 楊,馮海寬,黃" 玨,等. 基于無人機高光譜影像的馬鈴薯株高和地上生物量估算[J]. 農業機械學報,2021,52(2):188-198.

[9]陳仲新,任建強,唐華俊,等. 農業遙感研究應用進展與展望[J]. 遙感學報,2016,20(5):748-767.

[10]王福民,黃敬峰,唐延林,等. 新型植被指數及其在水稻葉面積指數估算上的應用[J]. 中國水稻科學,2007,21(2):159-166.

[11]杭艷紅,蘇" 歡,于滋洋,等. 結合無人機光譜與紋理特征和覆蓋度的水稻葉面積指數估算[J]. 農業工程學報,2021,37(9):64-71.

[12]孫玉婷,楊紅云,王映龍,等. 基于支持向量機的水稻葉面積測定[J]. 江蘇農業學報,2018,34(5):1027-1035.

[13]蘇中濱,陸藝偉,谷俊濤,等. 改進的QGA-ELM算法水稻葉面積指數反演模型[J]. 光譜學與光譜分析,2021,41(4):1227-1233.

[14]申關望,祁玉良,魯偉林,等. 施氮量對雜交水稻D優3138干物質轉運與稻米品質的影響[J]. 貴州農業科學,2017,45(2):23-25.

[15]劉振波,葛海嘯,葛云健,等. 冠層集聚指數對水稻LAI測量精度的影響[J]. 江蘇農業科學,2018,46(24):352-354.

[16]陳惠哲,朱德峰,林賢青,等. 穗肥施氮量對水稻劍葉生長及披垂的影響[J]. 西南農業學報,2007,20(6):1246-1249.

[17]劉镕源,王紀華,楊貴軍,等. 冬小麥葉面積指數地面測量方法的比較[J]. 農業工程學報,2011,27(3):220-224.

[18]田慶久,閔祥軍. 植被指數研究進展[J]. 地球科學進展,1998,13(4):327-333.

[19]Ahamed T,Tian L,Zhang Y,et al. A review of remote sensing methods for biomass feedstock production[J]. Biomass and Bioenergy,2011,35(7):2455-2469.

[20]Hunt E R Jr,Daughtry C S T,Eitel J U H,et al. Remote sensing leaf chlorophyll content using a visible band index[J]. Agronomy Journal,2011,103(4):1090-1099.

[21]Vincini M,Frazzi E,DAlessio P. A broad-band leaf chlorophyll vegetation index at the canopy scale[J]. Precision Agriculture,2008,9(5):303-319.

[22]Xu L Y,Xie X D,Li S. Correlation analysis of the urban heat island effect and the spatial and temporal distribution of atmospheric particulates using TM images in Beijing[J]. Environmental Pollution,2013,178:102-114.

[23]Dorigo W A,Zurita-Milla R,de Wit A J W,et al. A review on reflective remote sensing and data assimilation techniques for enhanced agroecosystem modeling[J]. International Journal of Applied Earth Observation and Geoinformation,2007,9(2):165-193.

[24]Jiang Z,Huete A,Didan K,et al. Development of a two-band enhanced vegetation index without a blue band[J]. Remote Sensing of Environment,2008,112(10):3833-3845.

[25]Gobron N,Pinty B,Verstraete M M,et al. Advanced vegetation indices optimized for up-coming sensors:design,performance,and applications[J]. IEEE Transactions on Geoscience and Remote Sensing,2000,38(6):2489-2505.

[26]le Maire G,Franois C,Dufrêne E. Towards universal broad leaf chlorophyll indices using PROSPECT simulated database and hyperspectral reflectance measurements[J]. Remote Sensing of Environment,2004,89(1):1-28.

[27]Pu R L,Gong P,Yu Q. Comparative analysis of EO-1 ALI and Hyperion,and landsat ETM+ data for mapping forest crown closure and leaf area index[J]. Sensors,2008,8(6):3744-3766.

[28]Haboudane D. Hyperspectral vegetation indices and novel algorithms for predicting green LAI of crop canopies:modeling and validation in the context of precision agriculture[J]. Remote Sensing of Environment,2004,90(3):337-352.

[29]Becker-Reshef I,Vermote E,Lindeman M,et al. A generalized regression-based model for forecasting winter wheat yields in Kansas and Ukraine using MODIS data[J]. Remote Sensing of Environment,2010,114(6):1312-1323.

[30]Wu C Y,Niu Z,Tang Q,et al. Estimating chlorophyll content from hyperspectral vegetation indices:modeling and validation[J]. Agricultural and Forest Meteorology,2008,148(8/9):1230-1241.

[31]Roujean J L,Breon F M. Estimating PAR absorbed by vegetation from bidirectional reflectance measurements[J]. Remote Sensing of Environment,1995,51(3):375-384.

[32]Birth G S,McVey G R. Measuring the color of growing turf with a reflectance Spectrophotometer1[J]. Agronomy Journal,1968,60(6):640-643.

[33]Gitelson A A. Wide Dynamic Range Vegetation Index for remote quantification of biophysical characteristics of vegetation[J]. Journal of Plant Physiology,2004,161(2):165-173.

[34]Broge N H,Leblanc E. Comparing prediction power and stability of broadband and hyperspectral vegetation indices for estimation of green leaf area index and canopy chlorophyll density[J]. Remote Sensing of Environment,2001,76(2):156-172.

[35]Gitelson A A,Kaufman Y J,Stark R,et al. Novel algorithms for remote estimation of vegetation fraction[J]. Remote Sensing of Environment,2002,80(1):76-87.

[36]Herrmann I,Pimstein A,Karnieli A,et al. LAI assessment of wheat and potato crops by VENμS and Sentinel-2 bands[J]. Remote Sensing of Environment,2011,115(8):2141-2151.

[37]Karal . Performance comparison of different kernel functions in SVM for different k value in k-fold cross-validation[C]//2020 Innovations in Intelligent Systems and Applications Conference (ASYU).Istanbul:IEEE,2020:1-5.

[38]秦占飛,申" 健,謝寶妮,等. 引黃灌區水稻葉面積指數的高光譜估測模型[J]. 武漢大學學報(信息科學版),2017,42(8):1159-1166.

[39]李艷大,孫濱峰,曹中盛,等. 基于作物生長監測診斷儀的雙季稻葉面積指數監測模型[J]. 農業工程學報,2020,36(10):141-149.

[40]田廣麗. 氮水平及栽培密度影響水稻生產能力的機制研究[D]. 南京:南京農業大學,2016.

[41]楊" 洪,李旭毅,卿發紅,等. 不同產量水平水稻群體光合特性和產量構成差異[J]. 江蘇農業學報,2023,39(5):1089-1096.

[42]Liu W D,Xiang Y Q,Zheng L F,et al. Relationship between rice LAI,CH.D,and hyperspectral data[C]//Multispectral and Hyperspectral Remote Sensing Instruments and Applications.Hangzhou,China:SPIE,2003:352-359.

[43]張笑宇,沈" 超,藺琛皓,等. 面向機器學習模型安全的測試與修復[J]. 電子學報,2022,50(12):2884-2918.

猜你喜歡
無人機水稻
什么是海水稻
有了這種合成酶 水稻可以耐鹽了
今日農業(2021年21期)2021-11-26 05:07:00
水稻種植60天就能收獲啦
軍事文摘(2021年22期)2021-11-26 00:43:51
油菜可以像水稻一樣實現機插
今日農業(2021年14期)2021-10-14 08:35:40
一季水稻
文苑(2020年6期)2020-06-22 08:41:52
水稻花
文苑(2019年22期)2019-12-07 05:29:00
高職院校新開設無人機專業的探討
人間(2016年26期)2016-11-03 17:52:40
利用無人機進行航測工作的方式方法
一種適用于輸電線路跨線牽引無人機的飛行方案設計
科技視界(2016年22期)2016-10-18 14:30:27
淺析無人機技術在我國的發展前景
企業導報(2016年9期)2016-05-26 20:58:26
主站蜘蛛池模板: 国产欧美一区二区三区视频在线观看| 亚洲天堂精品在线观看| 欧美成人国产| 国产91导航| 色亚洲成人| 亚洲天堂自拍| 亚洲欧洲免费视频| 久久99精品国产麻豆宅宅| 青青热久麻豆精品视频在线观看| 亚洲高清中文字幕| 在线观看无码av五月花| A级毛片高清免费视频就| 美女啪啪无遮挡| 亚洲一区二区三区香蕉| 无码电影在线观看| 欧美一区二区三区不卡免费| 亚洲欧美精品在线| 国产成人在线无码免费视频| 国产真实乱人视频| 国产熟睡乱子伦视频网站| 91美女视频在线| 88av在线| 久久永久精品免费视频| 亚洲精品视频免费| 亚洲国产精品久久久久秋霞影院| 亚洲婷婷丁香| 青草娱乐极品免费视频| 999福利激情视频| 91色在线观看| 一本综合久久| 欧美色图久久| 欧美日韩国产成人在线观看| 丰满人妻久久中文字幕| 91亚洲视频下载| 97人人做人人爽香蕉精品| 国产精品不卡片视频免费观看| 国产一级在线观看www色| 午夜国产理论| 日本午夜影院| 亚洲人成网7777777国产| av免费在线观看美女叉开腿| 婷婷色丁香综合激情| 国产激情无码一区二区三区免费| 亚洲二区视频| 青青热久免费精品视频6| 亚洲视频黄| 露脸真实国语乱在线观看| 国产在线视频欧美亚综合| 国产精品99一区不卡| 一级毛片免费观看不卡视频| 国产97区一区二区三区无码| 精品第一国产综合精品Aⅴ| 九九热这里只有国产精品| 天天干天天色综合网| 国产真实乱人视频| 国内99精品激情视频精品| 97青青青国产在线播放| 台湾AV国片精品女同性| 97视频在线精品国自产拍| V一区无码内射国产| 丝袜美女被出水视频一区| 十八禁美女裸体网站| 国产一区二区三区在线观看免费| 亚洲精品大秀视频| 浮力影院国产第一页| 中文精品久久久久国产网址| 国产精品lululu在线观看| 久久久国产精品无码专区| 中国成人在线视频| 青青青国产视频| 亚洲最大福利视频网| 久久亚洲天堂| 中文纯内无码H| 国产福利免费视频| 97国产在线观看| a欧美在线| 伊人精品视频免费在线| 狼友视频国产精品首页| 国产精品污视频| 国产一在线观看| 国禁国产you女视频网站| 日韩精品无码免费专网站|