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

基于地統(tǒng)計學(xué)的高原縣域耕地土壤養(yǎng)分最優(yōu)插值方法研究

2015-03-23 05:46:57唐希賢張英明
中國水土保持 2015年11期
關(guān)鍵詞:耕地

唐希賢,張英明

(1.臨夏縣農(nóng)村合作經(jīng)濟經(jīng)營管理站,甘肅 臨夏 731800;2.臨夏縣韓集鎮(zhèn)農(nóng)業(yè)技術(shù)推廣站,甘肅 臨夏 731800)

基于地統(tǒng)計學(xué)的高原縣域耕地土壤養(yǎng)分最優(yōu)插值方法研究

唐希賢1,張英明2

(1.臨夏縣農(nóng)村合作經(jīng)濟經(jīng)營管理站,甘肅 臨夏 731800;2.臨夏縣韓集鎮(zhèn)農(nóng)業(yè)技術(shù)推廣站,甘肅 臨夏 731800)

地統(tǒng)計學(xué);耕地;土壤養(yǎng)分;空間插值方法;甘肅

選擇適宜區(qū)域土壤屬性特征的空間插值方法是高效揭示土壤養(yǎng)分空間分異的基礎(chǔ)。以地處青藏高原與黃土高原過渡區(qū)的甘肅省臨夏縣為案例區(qū),利用該縣2007年農(nóng)業(yè)部測土施肥采集的5 475個樣點數(shù)據(jù),以及地統(tǒng)計學(xué)的交叉驗證和樣點驗證兩種評價方法,篩選出適合研究區(qū)耕地土壤養(yǎng)分(有機質(zhì)、堿解氮、有效磷和速效鉀)分布預(yù)測的高效插值方法。從交叉驗證結(jié)果來看,有機質(zhì)和堿解氮采用析取克里格插值精度較高,而有效磷和速效鉀分別采用普通克里格和簡單克里格插值精度較高。從樣點驗證結(jié)果來看,有機質(zhì)和速效鉀分別以析取克里格和簡單克里格插值精度較高,這也與交叉驗證相吻合;而堿解氮和有效磷分別采用簡單克里格和泛克里格插值精度較高,與交叉驗證結(jié)果并不一致。由于交叉驗證是具有一定“欺騙性”的自我驗證,本研究認為基于樣點驗證的精度評價方法更可靠,因此臨夏縣耕地堿解氮和有效磷適宜的插值方法分別為簡單克里格和泛克里格。

土壤養(yǎng)分是植物生長發(fā)育所必需的,準確獲取它的空間分布規(guī)律是實現(xiàn)精準農(nóng)業(yè)變量施肥的前提條件,而采樣與分析是目前掌握土壤養(yǎng)分狀況的最主要途徑[1-3]。然而,無論采樣密度多大,均不可能覆蓋所有區(qū)域,因此如何利用有限的樣點數(shù)據(jù)來獲得更為詳盡的土壤屬性空間分布信息一直是國內(nèi)外研究的熱點??臻g插值模型可以將離散的數(shù)據(jù)點轉(zhuǎn)化為連貫的數(shù)據(jù)曲面,從而對未知樣點的養(yǎng)分狀況進行有效預(yù)測,是表征土壤養(yǎng)分空間分布特征的重要手段[4]。但不同模型的插值原理和計算方法有所差異,其預(yù)測結(jié)果反映的土壤養(yǎng)分空間變異性可能不一樣[5-6]。因此,選擇適宜不同土壤屬性特征的插值方法對揭示其空間分布規(guī)律具有重要意義。

近年來,隨著“3S”技術(shù)的迅猛發(fā)展,基于地統(tǒng)計學(xué)的插值技術(shù)被廣泛應(yīng)用到土壤屬性空間變異性研究中。地統(tǒng)計學(xué)不僅考慮調(diào)查樣點本身的數(shù)值,還考慮與鄰近樣點的空間位置關(guān)系,經(jīng)常被用來研究既有一定隨機性又有一定結(jié)構(gòu)性的各種變量空間分布。但很多研究表明,不同土壤養(yǎng)分適宜的插值方法可能有所差異。劉吉平等[7]利用吉林省榆樹市弓棚鎮(zhèn)的386個采樣數(shù)據(jù)分析了不同插值方法對農(nóng)田土壤堿解氮空間變異性預(yù)測結(jié)果的影響,結(jié)果表明當采樣密度較大時,克里格模型插值精度高于BP神經(jīng)網(wǎng)絡(luò)模型。孫義祥等[8]基于1 975個調(diào)查樣點的研究結(jié)果表明,不同插值方法對縣域土壤有效磷空間變異性的預(yù)測效果差異較大,且表現(xiàn)為析取克里格、普通克里格、簡單克里格和局部多項式插值精度高于全局多項式、反距離加權(quán)法、泛克里格和徑向基函數(shù)法,以析取克里格法插值效果最好。李曉燕等[9]通過364個采樣數(shù)據(jù)分析不同插值方法的預(yù)測精度,結(jié)果表明,吉林省德惠市土壤速效鉀的空間變異性預(yù)測適宜采用普通克里格指數(shù)模型,其插值效果優(yōu)于反距離權(quán)重法和球面多項式法。

從以上研究也可以看出,目前土壤養(yǎng)分空間插值方法預(yù)測精度研究主要以小規(guī)模樣本數(shù)據(jù)集為基礎(chǔ),且有關(guān)同一研究區(qū)不同土壤養(yǎng)分對插值模型的適宜性的研究相對較少,這導(dǎo)致土壤養(yǎng)分變異性研究存在較大不確定性。青藏高原位于我國西南部,約占我國陸地總面積的25%,是我國面積最大的高原;黃土高原位于我國中部偏北,面積約為40萬km2,是世界最大的黃土堆積區(qū),暴雨徑流沖刷疏松黃土導(dǎo)致水土流失嚴重[10]。因此,確定適用于高原地貌類型區(qū)的耕地土壤養(yǎng)分點面拓展模型,對合理設(shè)計養(yǎng)分流失防控措施和肥力提升計劃皆具有十分重要的意義[11]。為此,本研究以地處青藏高原與黃土高原過渡地區(qū)的甘肅省臨夏縣為案例區(qū),利用該縣2007年農(nóng)業(yè)部測土施肥大規(guī)模采集的5 475個樣點數(shù)據(jù),研究地統(tǒng)計學(xué)中常用的普通克里格、簡單克里格、泛克里格和析取克里格4種插值方法對該地區(qū)耕地土壤有機質(zhì)、堿解氮、有效磷和速效鉀空間變異性預(yù)測精度的影響,以期找到不同土壤養(yǎng)分最適宜的插值模型,為我國青藏高原與黃土高原過渡區(qū)縣域耕地養(yǎng)分預(yù)測選擇高效插值方法提供科學(xué)依據(jù)。

1 研究區(qū)概況

臨夏縣位于甘肅省中部,地理坐標為北緯34°57′~36°12′、東經(jīng)102°41′~103°40′(圖1)[12],地處青藏高原與黃土高原的過渡地帶,屬溫帶半濕潤氣候區(qū)。臨夏縣地形較為復(fù)雜,海拔差距較大,總體呈西南高東北低趨勢??h域南部和西部均是山地,海拔一般在3 000 m 以上;中東部以黃土高原為主,海拔介于1 750~2 300 m之間,地貌青藏、黃土高原兼有,屬甘肅省典型地貌類型。根據(jù)第二次土壤普查統(tǒng)計結(jié)果,臨夏縣耕地土壤類型包括壚土、紅土、黃綿土、黑土和山地棕壤,且以黑土、壚土和紅土面積分布最廣,占全縣耕地總面積的95%以上[13]。

圖1 臨夏縣地理位置圖

2 數(shù)據(jù)來源及研究方法

2.1 數(shù)據(jù)來源

本研究使用的耕地土壤調(diào)查樣點數(shù)據(jù)來源于農(nóng)業(yè)部測土配方施肥項目。2007年8—9月臨夏縣共采集土樣5 475份,取樣采用多點混合土樣采集法,取土深度為20 cm,并利用GPS定位獲取樣點經(jīng)緯度。樣點屬性包括采樣地塊基本情況(統(tǒng)一編號、經(jīng)緯度、土壤類型等)和土壤測試結(jié)果(有機質(zhì)、堿解氮、有效磷、速效鉀等),其中有機質(zhì)采用K2Cr2O7-H2SO4溶液油浴法測試,堿解氮采用堿解擴散法,有效磷采用NaHCO3浸提-鉬銻抗比色法,速效鉀采用NH4OAc浸提-火焰光度法。利用ArcGIS 9.3軟件將樣點經(jīng)緯度坐標投影轉(zhuǎn)換成北京1954直角坐標系坐標,并與縣域耕地利用現(xiàn)狀圖相匹配,然后根據(jù)編碼規(guī)則將樣點屬性值導(dǎo)入圖層并保存入庫,最終形成臨夏縣耕地土壤調(diào)查樣點數(shù)據(jù)庫。

2.2 數(shù)據(jù)處理軟件與地統(tǒng)計學(xué)插值方法

本研究利用SPSS 19.0分析土壤有機質(zhì)、堿解氮、有效磷和速效鉀的基本統(tǒng)計特征,借助ArcGIS 9.3的地統(tǒng)計分析工具對縣域耕地土壤養(yǎng)分分別進行普通克里格(Ordinary Kriging,OK)、簡單克里格(Simple Kriging,SK)、泛克里格(Universal Kriging,UK)和析取克里格(Disjunctive Kriging,DK)空間插值。這4種插值模型均以變異函數(shù)理論和結(jié)構(gòu)分析為基礎(chǔ),能最大限度地利用樣本信息對區(qū)域化變量進行無偏最優(yōu)估計,在考慮樣點距離的同時,還考慮已知樣點的空間分布及與未知樣點的空間方位關(guān)系,是目前應(yīng)用較為廣泛的地統(tǒng)計學(xué)插值方法[14]。

本研究從兩個方面對插值精度進行評估:一是地理信息系統(tǒng)軟件自帶的交叉驗證功能,具體評定標準為平均預(yù)測誤差(ME)和標準預(yù)測誤差(MSE)越接近于0越好,均方根誤差(RMSE)越小且越接近于平均標準差(ASE)越好,均方根標準誤差(RMSSE)越接近于1越好,且優(yōu)先考慮RMSE[15]。二是從調(diào)查樣點中隨機抽取20%(1 095個采樣點)作為驗證樣點(圖2),通過預(yù)測值與實測值的相關(guān)系數(shù)r及RMSE兩個指標進行預(yù)測精度評定,當RMSE越小、r值越大時,預(yù)測精度越高。

(1)

式中:N為樣點數(shù)量;Xoi為樣點實測值,Xpi為樣點預(yù)測值。

圖2 基于耕地土壤的樣點分布

3 結(jié)果與討論

3.1 基于交叉驗證的空間插值方法精度評價

本研究利用ArcGIS 9.3對臨夏縣2007年采集的5 475個耕地土壤調(diào)查樣點數(shù)據(jù)進行了普通克里格、簡單克里格、泛克里格和析取克里格4種不同地統(tǒng)計學(xué)方法的空間插值,并通過該軟件的交叉驗證功能進行精度評價(見表1)。從表1可以看出,縣域耕地土壤有機質(zhì)和堿解氮的RMSE分別為5.14~5.21 g/kg和20.76~21.32 mg/kg,均以析取克里格的RMSE最小,表明析取克里格插值效果優(yōu)于普通克里格、簡單克里格和泛克里格,能更好地模擬臨夏縣耕地土壤有機質(zhì)和堿解氮的空間分布特征。有效磷和速效鉀的RMSE分別為17.06~17.12 mg/kg和32.18~32.38 mg/kg,且有效磷以普通克里格的RMSE最小,速效鉀以簡單克里格的RMSE最小,說明臨夏縣耕地土壤有效磷和速效鉀分別采用普通克里格和簡單克里格的插值精度較高。從各種土壤養(yǎng)分不同插值方法交叉驗證的誤差大小來看,泛克里格的RMSE均高于普通克里格、簡單克里格和析取克里格,表明臨夏縣耕地土壤養(yǎng)分插值以泛克里格預(yù)測精度最低,因此不宜采用該方法進行縣域土壤養(yǎng)分空間變異特征分析。從不同土壤養(yǎng)分的歸一化均方根誤差(RMSE與平均值的比值)來看,有效磷歸一化均方根誤差高于59.00%,而有機質(zhì)、堿解氮和速效鉀歸一化均方根誤差均低于30.00%,表明基于地統(tǒng)計學(xué)的臨夏縣耕地土壤有效磷插值效果不夠理想,預(yù)測精度明顯低于有機質(zhì)、堿解氮和速效鉀。

表1 臨夏縣耕地土壤養(yǎng)分地統(tǒng)計學(xué)插值方法交叉驗證精度評價

3.2 基于驗證樣點的空間插值方法精度評價

基于驗證樣點的空間插值結(jié)果表明,基于析取克里格插值的臨夏縣耕地土壤有機質(zhì)含量預(yù)測值與實測值相關(guān)系數(shù)r為0.588 7,高于普通克里格(0.579 0)、簡單克里格(0.582 0)和泛克里格(0.579 1),RMSE為5.177 1 g/kg,小于普通克里格(5.288 2)、簡單克里格(5.208 6)和泛克里格(5.282 2),表明縣域耕地土壤有機質(zhì)含量預(yù)測值析取克里格插值精度較高,這與交叉驗證得出的結(jié)果相一致;另外,基于簡單克里格插值的縣域耕地土壤速效鉀含量預(yù)測值與實測值相關(guān)系數(shù)r為0.747 3,高于其他3種插值方法,且RMSE最小,說明臨夏縣耕地土壤速效鉀含量預(yù)測值簡單克里格插值效果較好,這也與交叉驗證結(jié)果相一致。但從空間插值結(jié)果也可知,臨夏縣耕地土壤堿解氮和有效磷含量預(yù)測值分別以簡單克里格和泛克里格插值精度最高,這與上述的交叉驗證結(jié)果不一致。交叉驗證是ArcGIS軟件利用樣點數(shù)據(jù)進行的自我驗證,得到的樣點預(yù)測值已經(jīng)受到該樣點實測值的影響,從理論上說具有一定的“欺騙性”,而基于實測樣點所進行的精度評價不同于交叉驗證的自我驗證,其已剔除了驗證樣點實測數(shù)據(jù)對預(yù)測值的影響,是完全通過已知點來推算未知點的養(yǎng)分含量,并利用驗證點進行精度分析的一種檢驗方法[16],因此基于實測樣點的精度分析結(jié)論比交叉驗證有更大的可靠性,故本研究認為臨夏縣耕地土壤堿解氮和有效磷含量預(yù)測分別采用簡單克里格和泛克里格插值模型相對更優(yōu)。

從臨夏縣不同地統(tǒng)計學(xué)插值方法對各個養(yǎng)分指標預(yù)測的適宜性來看,有效磷含量1 095個驗證樣點的預(yù)測值與實測值相關(guān)系數(shù)r均低于0.40,明顯小于有機質(zhì)、堿解氮和速效鉀,表明臨夏縣耕地土壤有效磷含量的地統(tǒng)計學(xué)空間插值精度低于有機質(zhì)、堿解氮和速效鉀。該結(jié)論與交叉驗證結(jié)果是一致的,這可能與有效磷本身的空間分布特征有關(guān)。臨夏縣耕地土壤有效磷含量極大值與極小值相差200多倍,變異系數(shù)為62.86%,是有機質(zhì)和堿解氮變異系數(shù)的2倍、速效鉀變異系數(shù)的1.5倍。有研究表明,強烈的空間變異性對插值效果具有一定的影響,所以有效磷含量插值精度整體較低[17]。有機質(zhì)、堿解氮和速效鉀1 095個驗證樣點的預(yù)測值與實測值相關(guān)系數(shù)均高于0.57,其中速效鉀的相關(guān)系數(shù)r均高于0.74;另一方面,有機質(zhì)、堿解氮和速效鉀1 095個驗證樣點的歸一化均方根誤差均低于29.00%,其中有機質(zhì)歸一化均方根誤差均低于25.00%,而有效磷歸一化均方根誤差均高于57.00%,約是其他3種養(yǎng)分的2倍,表明基于地統(tǒng)計學(xué)的臨夏縣耕地土壤有機質(zhì)、堿解氮和速效鉀插值效果整體較好,預(yù)測精度明顯高于有效磷,這與交叉驗證結(jié)果相一致。造成這一現(xiàn)象的主要原因可能是有機質(zhì)、堿解氮和速效鉀的空間變異程度均相對較低,有利于土壤屬性值的空間預(yù)測,該結(jié)論也符合很多學(xué)者認為的“較低的空間變異性在一定程度上能保證較高的空間預(yù)測精度”的一般規(guī)律[17-19]。

3.3 縣域耕地土壤養(yǎng)分空間分布特征

4種不同地統(tǒng)計學(xué)空間插值方法預(yù)測的臨夏縣耕地養(yǎng)分空間分布特征基本一致,即:有機質(zhì)和堿解氮含量以中部及西北部較高,東北部較低;有效磷含量以東北部較高,中部較低;速效鉀含量以東南部較高,其他區(qū)域分布較為均勻,且以東北部含量最低。上述分布規(guī)律可能與該地區(qū)復(fù)雜的地形特征有關(guān),臨夏縣地勢西南高而東北低,西部和南部多為山地,是青藏高原東北的邊緣隆起部分,海拔一般在3 000 m以上,氣溫相對較低,微生物活性較差,土壤有機質(zhì)和全氮分解較慢;而中東部地區(qū)屬黃土高原,海拔介于1 750~2 300 m之間,明顯低于西南部,氣溫相對較高,土壤微生物生理活動旺盛,土壤有機質(zhì)和全氮分解速度較快[20]。此外,臨夏縣降水量分布極為不均,呈東北部少而西南部多的規(guī)律,而很多研究表明降水量大有利于土壤有機質(zhì)和全氮積累,這是東北部有機質(zhì)和全氮含量低的另一原因[21]。有效磷含量分布規(guī)律正好與有機質(zhì)和全氮相反,可能是因為高海拔地區(qū)的土壤磷素受雨水沖刷而通過地表徑流的方式聚集在了地勢較低的地區(qū)。臨夏縣東南部速效鉀含量明顯高于其他區(qū)域,這可能與該區(qū)域耕地土壤質(zhì)地以黏土和黏壤土為主有很大關(guān)系,有研究表明,土壤黏粒含量與速效鉀含量呈極顯著的正相關(guān)[22]。

4 結(jié) 論

本研究以大規(guī)模采集的5 475個樣點數(shù)據(jù)為基礎(chǔ),利用兩種驗證方法綜合評價了普通克里格、簡單克里格、泛克里格和析取克里格4種不同地統(tǒng)計學(xué)插值方法對臨夏縣耕地土壤養(yǎng)分含量空間預(yù)測精度的影響,結(jié)果表明:基于驗證樣點的空間插值方法精度評價比交叉驗證法具有更高的可靠性;臨夏縣耕地土壤有機質(zhì)含量預(yù)測適宜采用析取克里格插值方法,堿解氮和速效鉀含量預(yù)測適宜采用簡單克里格插值方法,有效磷含量預(yù)測適宜采用泛克里格插值方法。

[1] 趙明松,張甘霖,王德彩,等.徐淮黃泛平原土壤有機質(zhì)空間變異特征及主控因素分析[J].土壤學(xué)報,2013,50(1):1-11.

[2] 李翔,潘瑜春,馬景宇,等.基于多種土壤養(yǎng)分的精準管理分區(qū)方法研究[J].土壤學(xué)報,2007,44(1):14-20.

[3] 李啟全,岳天祥,范澤孟,等.中國表層土壤全氮的空間模擬分析[J].地理研究,2010,29(11):1981-1992.

[4] 楊子清,陳平留,劉健,等.基于Kriging法的森林土壤養(yǎng)分空間插值[J].福建農(nóng)林大學(xué)學(xué)報:自然科學(xué)版,2012,41(3):296-300.

[5] 王秀,苗孝可,孟志軍,等.插值方法對GIS土壤養(yǎng)分插值結(jié)果的影響[J].土壤通報,2005,36(6):826-830.

[6] 趙永存,史學(xué)正,于東升,等.不同方法預(yù)測河北省土壤有機碳密度空間分布特征的研究[J].土壤學(xué)報,2005,42(3):379-385.

[7] 劉吉平,于佳,于洋,等.采樣密度和插值對農(nóng)田土壤堿解氮空間變異性的影響[J].土壤通報,2014,45(4):789-794.

[8] 孫義祥,吳傳洲,朱克保,等.插值方法與樣點數(shù)對縣域土壤有效磷空間變異特征評價的影響[J].應(yīng)用生態(tài)學(xué)報,2009,20(3):673-678.

[9] 李曉燕,張樹文.吉林省德惠市土壤速效鉀的空間分異及不同插值方法的比較[J].水土保持學(xué)報,2004,18(4):97-100.

[10] 信忠保,許炯心,余新曉.近50年黃土高原水土流失的時空變化[J].生態(tài)學(xué)報,2009,29(3):1129-1139.

[11] 陳濤,常慶瑞,劉京,等.黃土高原南麓縣域耕地土壤速效養(yǎng)分時空變異[J].生態(tài)學(xué)報,2013,33(2):554-564.

[12] 劉興華.耕地資源社會保障問題研究——以甘肅省臨夏縣為例[D].蘭州:甘肅農(nóng)業(yè)大學(xué),2012.

[13] 甘肅省土壤普查辦公室.甘肅土壤[M].北京:農(nóng)業(yè)出版社,1990.

[14] 張峰,杜群,葛宏立,等.基于地統(tǒng)計學(xué)和CFI樣地的浙江省森林碳空間分布研究[J].生態(tài)學(xué)報,2012,32(16):5275-5286.

[15] ESRI.Using Geostatistical Analyst[M].United States of America:ESRI Inc,2004.

[16] 范勝龍,黃炎和,林金石.表征土壤有機碳區(qū)域分布的優(yōu)化空間插值模型研究——以福建省龍海市為例[J].水土保持研究,2011,18(6):1-5.

[17] 趙永存.土壤屬性表征的空間尺度效應(yīng)和不確定性研究——以河北省土壤有機碳為例[D].南京:中國科學(xué)院南京土壤研究所,2005.

[18] 趙春生,張佳寶.土壤特性空間變異研究的定量方法[J].土壤學(xué)報,1995,23(5):45-53.

[19] 張忠啟.樣點布置模式及密度對揭示土壤有機碳空間變異的影響[D].南京:中國科學(xué)院南京土壤研究所,2010.

[20] 張金屯.全球氣候變化對自然土壤碳、氮循環(huán)的影響[J].地理科學(xué),1998,18(5):463-471.

[21] 陳建國,田大倫,閆文德,等.鄂西南區(qū)域土壤有機碳儲量、密度及其影響因子[J].中南林業(yè)科技大學(xué)學(xué)報,2011,31(5):57-62.

[22] 龐夙,陶曉秋,張英,等.會理縣新植煙區(qū)土壤速效鉀含量空間變異特征及其影響因子[J].中國煙草科學(xué),2012,33(1):32-42.

(責(zé)任編輯 徐素霞)

S158.2

A

1000-0941(2015)11-0049-04

唐希賢(1987—),男,甘肅臨夏縣人,農(nóng)業(yè)經(jīng)濟師,主要從事農(nóng)業(yè)科技推廣、農(nóng)業(yè)產(chǎn)業(yè)化及農(nóng)村經(jīng)濟管理工作;通信作者張英明(1985—),男,甘肅臨夏縣人,助理農(nóng)藝師,主要從事農(nóng)業(yè)科技推廣工作。

2015-03-25

猜你喜歡
耕地
自然資源部:加強黑土耕地保護
我國將加快制定耕地保護法
堅決落實耕地保護“軍令狀” 牢牢掌握糧食安全主動權(quán)
浙江人大(2022年4期)2022-04-28 21:37:09
保護耕地
北京測繪(2021年12期)2022-01-22 03:33:36
新增200億元列入耕地地力保護補貼支出
耕地保護需要強化系統(tǒng)觀念
耕地種田也能成為風(fēng)景
當代陜西(2019年10期)2019-06-03 10:12:12
耕地時節(jié)
扎緊四個“口袋” 打造耕地保護新常態(tài)
耕地質(zhì)量
主站蜘蛛池模板: 四虎国产精品永久一区| 露脸国产精品自产在线播| 亚洲欧美另类视频| 亚洲国产精品一区二区高清无码久久| 毛片大全免费观看| 乱人伦99久久| 91无码人妻精品一区| 91精品国产一区自在线拍| 午夜精品国产自在| 亚洲AV电影不卡在线观看| 国产一区三区二区中文在线| 五月激情婷婷综合| 日韩欧美中文| 国产香蕉97碰碰视频VA碰碰看| 国产福利大秀91| 真人高潮娇喘嗯啊在线观看| 91青青草视频| 99久久这里只精品麻豆| 91小视频在线观看| 456亚洲人成高清在线| 久久久久亚洲AV成人网站软件| 色网在线视频| 国产在线精品香蕉麻豆| 国产激情影院| a亚洲视频| 国产性生交xxxxx免费| 国产精品自在自线免费观看| 亚洲精品午夜无码电影网| 亚洲永久色| 亚洲综合天堂网| 亚洲综合色婷婷中文字幕| 男人的天堂久久精品激情| 天堂中文在线资源| 中文字幕亚洲另类天堂| 午夜日b视频| 欧美日韩第三页| 青青草原国产| 日日摸夜夜爽无码| 国产精品成人观看视频国产| 天天躁夜夜躁狠狠躁图片| 色国产视频| 亚洲最新地址| 东京热av无码电影一区二区| 亚洲综合国产一区二区三区| 91小视频版在线观看www| 午夜无码一区二区三区在线app| 亚洲一区二区成人| 日本一区中文字幕最新在线| 亚洲国产在一区二区三区| 欧美无专区| 亚洲男人的天堂在线观看| 高清无码一本到东京热| 无码综合天天久久综合网| 成年A级毛片| 婷婷六月综合网| 国产一区二区三区精品欧美日韩| 国产欧美精品一区二区| 9久久伊人精品综合| 亚洲精品无码在线播放网站| 在线精品亚洲一区二区古装| 香蕉在线视频网站| 欧美综合区自拍亚洲综合绿色| 91精品日韩人妻无码久久| 国产性生大片免费观看性欧美| 狠狠色综合网| 在线观看的黄网| 成人在线天堂| 日韩无码真实干出血视频| 国产视频a| 激情六月丁香婷婷四房播| 婷婷久久综合九色综合88| 精品视频一区二区三区在线播| 大学生久久香蕉国产线观看| 国产成人精品男人的天堂下载 | 精品国产免费观看| 国产91视频观看| 在线免费亚洲无码视频| 综合亚洲网| 色哟哟精品无码网站在线播放视频| 狠狠色香婷婷久久亚洲精品| 伊人福利视频| 全部免费特黄特色大片视频|