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

砒砂巖區植被覆蓋度環境驅動因子量化分析
——基于地理探測器

2022-05-20 09:47:22馬曉妮任宗萍謝夢瑤李占斌
生態學報 2022年8期
關鍵詞:區域

馬曉妮,任宗萍,*,謝夢瑤,李占斌,李 鵬,2,張 星

1 西安理工大學省部共建西北旱區生態水利國家重點實驗室,西安 710048 2 西安理工大學旱區生態水文與災害防治國家林業局重點實驗室,西安 710048

植被作為陸地生態系統的主要組成部分,同氣候變化關系緊密,在全球物質和能量的循環過程中起著至關重要的作用[1]。近年來,在全球氣候變化以及人類活動的影響下,陸地植被發生著不同程度和方式的變化[2]。特別是進入20世紀80年代以來,全球范圍內的植被覆蓋率一直在持續增加,成為生物圈變化最顯著的特征之一[3]。其中,近30年來中國歸一化植被指數(Normalized Difference Vegetation Index,NDVI)呈增加趨勢的面積占比達53.8%,極顯著增加的面積達29.3%[4]。黃土高原作為退耕還林還草工程實施的重點區域,是我國植被覆蓋變化最顯著的地區,約43.4%的區域植被覆蓋度呈顯著改善[5]。因此,掌握該區植被動態變化以及探索退耕后植被覆蓋度對環境因子的響應特征至關重要。

旱區的植被覆蓋率增加主要受降水、溫度、土壤類型和植被類型的影響,位于干旱半干旱地區的黃土高原,其不同區域的植被恢復對水熱條件的響應也不盡相同[6-7]。砒砂巖區位于黃土高原北部,地形破碎,降水量少且以暴雨為主,土壤類型以砒砂巖土、黃綿土、風沙土為主,土層較薄,抗沖抗蝕性較差,難以滿足大多數植物生長的需求,導致該區植被覆蓋特征以及其驅動因子有別于典型的黃土區[8]。砒砂巖區受風力、水力、凍融等作用,區域內侵蝕劇烈,為黃河粗泥沙的集中來源區,植被恢復是該區域主要的治理措施之一[9]。近年來先后實施了一系列項目對區域生態環境狀況進行治理,但由于該區生態環境波動性明顯,原生環境極易受到損壞,面對全球氣候變化,脆弱而敏感的砒砂巖地區近十幾年來植被發生的變化以及與環境因子的關系有待考究[10-12],而且現有關于砒砂巖區植被變化特征及其影響因子研究中,綜合考慮多因子對植被覆蓋度的響應研究較少。

地理探測器是一種新的空間分析模型,主要是用于探測地理要素的空間分異性,定量化多環境驅動因子對植被覆蓋度的直接和間接作用[13]。因此本文基于1999—2018年的NDVI數據,結合氣象數據、土壤數據以及DEM數據,通過Mann-Kendall趨勢檢驗、Hurst指數分析砒砂巖區1999—2018年植被覆蓋度時空變化特征,利用地理探測器方法分析氣候、土壤、地形等因子對植被覆蓋變化的解釋力,確定促進植被生長的環境因子最適宜的范圍和類型,對砒砂巖區生態效益的合理評估以及進一步的植被建設等具有重要的理論和實踐意義。

1 研究區概況

砒砂巖區位于38°10′—40°10′N,108°45′—111°31′E,是黃土高原侵蝕最劇烈的區域,主要分布在內蒙古自治區鄂爾多斯市的東勝區、準格爾旗、伊金霍洛旗、達拉特旗、杭錦旗,在陜西省的神木、府谷兩縣,山西省的河曲、保德兩縣,內蒙古的清水河縣有零星分布[14]。研究區總面積1.67萬km2,主要包括了鄂爾多斯盆地在中生代形成的地層中的砂巖、粉砂巖和泥巖等多種沉積巖,該區西北部為低丘陵地貌,東南部為高丘陵地貌,最東部黃河西岸為低山地貌[15]。屬于暖溫帶北緣半干旱大陸氣候,年平均氣溫6.5—9.6℃,年平均降雨量315—442mm,年平均蒸發量2200—2600mm。根據地被物類型及覆蓋程度將其分為覆土砒砂巖區、覆沙砒砂巖區、裸露砒砂巖區三種類型(圖1)。

圖1 砒砂巖區位置圖Fig.1 Location map of the Pisha sandstone area

2 數據與方法

2.1 數據來源與處理

遙感數據:砒砂巖區1999—2018年NDVI數據(分辨率1000m),來源于資源環境科學數據中心的中國年植被指數(NDVI)空間分布數據集(http://www.resdc.cn)。數字高程模型(分辨率30m)、土地利用數據(分辨率30m),來源于地理空間數據云(http://www.gscloud.cn)。土壤類型數據(分辨率1000m),來源于聯合國糧農組織網(http://www.fao.org)。土壤體積含水量數據(分辨率0.0625°×0.0625°),來源于中國氣象數據網CLDAS-V2.0(http://data.cma.cn)。利用Arcgis10.2軟件對以上遙感數據進行投影轉換和掩膜提取,并重采樣為1000m分辨率的空間數據。

氣象數據:1999—2018年砒砂巖區及其周邊共12個氣象站點的日降雨量以及日平均氣溫數據(圖1),來源于資源環境科學與數據中心的中國氣象要素站點觀測逐日數據集;基于反距離權重內插法生成1000m分辨率柵格數據。

地理探測器輸入數據處理:基于Arcgis10.2軟件創建漁網工具,生成1km×1km格網,共15732個格點作為采樣點,采用自然斷點法將年均降水量、年均氣溫、土壤水分、海拔、坡度分為9類,坡向分為9類(表1),土壤類型分為17類,土地利用分為6類,提取植被覆蓋度及各個環境因子數據至采樣點輸入地理探測器軟件,用于分析植被覆蓋度與環境因子之間的關系。

表1 環境因子自然斷點法分區Table 1 Classification of environmental factors by natural breaks

2.2 研究方法

2.2.1植被覆蓋度計算

根據像元二分模型原理推求砒砂巖區植被覆蓋度,計算公式如下[16]:

(1)

其中:

(2)

(3)

式中:NDVIsoil為完全是裸土或無植被覆蓋區域的NDVI值,NDVIveg則代表完全被植被所覆蓋的像元的NDVI值。NDVImax、NDVImin分別為歸一化植被指數的最大值和最小值。

由于采用的是最大值合成法生成的植被指數數據集,參考相關文獻[16],近似取FVCmax為100%,FVCmin為0%,即公式可變為

(4)

為避免噪聲影響,NDVImax、NDVImin一般取累積概率5%和95%的NDVI作為區域最小和最大的NDVI值。

2.2.2植被覆蓋度變化趨勢分析方法

本研究采用變異系數(Cv)表示砒砂巖區植被覆蓋度變化的穩定性,分別采用Mann-Kendall方法和Hurst指數逐像元分析近20年區域植被覆蓋度變化趨勢以及未來變化趨勢。

(1)變異系數(Cv)

變異系數Cv反映的是相對變異,即隨機變量的離散程度,計算公式為:

Cv=σ/μ

(5)

式中:σ為標準差;μ為均值。

根據Nielsen的劃分標準,當Cv≤10%、10%

(2)Mann-Kendall趨勢檢驗及Sen傾斜度

Mann-Kendall是水文氣象時間序列數據趨勢檢驗中使用最廣泛的非參數檢驗方法,其優點是不需要樣本遵循特定的分布,且不受異常值的干擾,用于分析近20年砒砂巖區植被覆蓋度變化趨勢,具體公式見文獻[18]。

趨勢性大小以及變化方向使用Sen傾斜度表示[19],其中β的計算公式如下:

(6)

式中:β表示植被覆蓋度上升或下降的趨勢,其中β>0代表植被覆蓋度呈上升趨勢,β<0代表植被覆蓋度呈下降趨勢。

(3)Hurst指數

Hurst指數可判斷時間序列未來變化的趨勢,多用于定量描述長時間序列信息相關性,其取值范圍為H(0

2.2.3植被覆蓋度驅動因子定量分析方法

(1)地理探測器

地理探測器是探測空間分異性以及揭示其背后驅動力的一組統計學方法,本論文應用了地理探測器中的因子探測器,交互探測器,風險探測器以及生態探測器來研究植被驅動力以及多因子交互作用[21]。因子探測器主要探測不同環境因子對FVC空間分異的影響大小,其中q值越高,影響力越大,表達式為:

(7)

(8)

式中:h為FVC值影響因子的分類或分區;Nh和N分別為層h和全區的單元數;σh2和σ2分別為層h和全區的Y值的方差。SSW和SST分別為層內方差之和及全區總方差,q的值域為[0,1]。

交互探測器主要是識別不同因子之間的交互作用,即評估兩個因子共同作用時是否會增加或減弱對FVC的解釋力或這些因子對FVC的影響是相互獨立的。評估方法是首先分別計算兩種影響因素X1和X2對Y的q值,并且計算它們交互時的q值,對三者之間的q值大小進行比較,主要分為以下幾種結果(表2)。

表2 交互探測器判斷依據Table 2 Judgment basis of interaction detector

生態探測器用于比較兩因子對FVC空間分布的影響是否有顯著差異,用F統計量來檢驗:

(9)

式中:NX1和NX2表示兩因子的樣本量;SSWX1和SSWX2表示兩因子形成分層的層內方差之和。

風險探測器用于探測不同環境因子對FVC影響的適宜范圍或類型,用t統計量來檢驗:

(10)

(2)偏相關分析

偏相關分析用于進一步探究砒砂巖區植被覆蓋度對單因子降水和氣溫的響應,計算公式:

(11)

式中:Rxy,z為自變量z固定后因變量x與自變量y的偏相關系數;Rxy,Rxz,Ryz為兩因子間的相關系數;x為植被覆蓋度,y為年降水量,z為年均氣溫。

3 結果與分析

3.1 植被覆蓋度時空特征

基于砒砂巖區NDVI數據集(1999—2018年)提取了該區植被覆蓋度,用于分析砒砂巖區植被時空變化特征。砒砂巖區平均植被覆蓋度年際變化如圖2所示,區域內植被覆蓋度總體呈波動上升趨勢,平均上升幅度為0.086/10a。根據砒砂巖區平均植被覆蓋度年際變化圖,2013年區域平均植被覆蓋度在研究時段內達到最大值(56.5%),將2013年植被覆蓋度與1999年及2018年區域植被覆蓋度進行成對樣本T檢驗,結果顯示,1999年和2018年的區域植被覆蓋度與2013年相比均存在顯著性差異(P<0.01),2013年植被覆蓋相比1999年顯著增加,而與2018年相比出現減少趨勢,表明該區植被覆蓋度敏感性較強。

圖2 砒砂巖區平均植被覆蓋度年際變化圖 Fig.2 Annual variation of mean vegetation coverage in the Pisha sandstone area陰影代表95%置信區間

統計1999—2018年平均植被覆蓋度,并根據SL190—2007《土壤侵蝕分類分級標準》[22],即10%、30%、45%、60%四個節點將砒砂巖區分為裸地、低植被覆蓋度、中低植被覆蓋度、中植被覆蓋度、高植被覆蓋度5個等級,得到砒砂巖區多年平均植被覆蓋度空間分級圖(圖3)。砒砂巖區多年平均植被覆蓋度為42.3%,空間上總體呈現出東南高、西北低的空間格局,即覆土區>覆沙區>裸露區。分別統計各砒砂巖分區不同植被覆蓋度等級的面積比例可知覆土區處于中植被覆蓋度水平,該等級面積占比57.1%;覆沙區與裸露區以中低植被覆蓋度為主,占比分別是56.1%、48.4%。

圖3 1999—2018年砒砂巖區平均植被覆蓋度等級及變異系數空間分布Fig.3 Spatial distribution of average vegetation coverage grade and coefficient of variation in the Pisha sandstone area from 1999 to 2018

砒砂巖區近20年植被覆蓋度年均空間變化范圍為4.8%—89.5%,變異系數空間變化范圍為0.10—0.90(圖3),整個砒砂巖區植被覆蓋度變化均屬于中等程度變異。由表3可知,四個分區中,覆土區植被覆蓋度最大,均值為49.0%,平均變化率為0.096/10 a,變異系數最小,最穩定,而裸露區植被覆蓋度最小,均值為30.7%,生態系統比其他類型區脆弱。

表3 不同砒砂巖區植被覆蓋度變化的統計特征Table 3 Statistical characteristics of vegetation coverage changes in different Pisha sandstone areas

3.2 植被覆蓋度變化趨勢分析

1999—2018年砒砂巖區植被覆蓋變化趨勢圖(圖4)表明,呈增加趨勢的面積為15383.58km2,約占砒砂巖區總面積的92.2%;呈減少趨勢的面積為1322.64km2,約占砒砂巖區總面積的7.8%。正值區域大于負值區域面積,表明近20年砒砂巖區植被有改善的趨勢。由表4可知,植被覆蓋度極顯著增加的面積占45.6%,主要分布在覆土區和覆沙區,兩者分別占25.6%和10.2%;顯著增加的面積占砒砂巖區總面積的16.8%,其中裸露區顯著增加的面積為4.5%。顯著與極顯著增加的區域基本相同,主要分布在砒砂巖區東部區域。顯著和極顯著減少的區域面積占比分別為0.6%和0.5%,零星分布在裸露區中部??偟膩碚f,砒砂巖區植被覆蓋度整體得到改善的面積比例遠大于退化。

圖4 砒砂巖區植被覆蓋度變化趨勢及可持續性變化等級空間分布Fig.4 Spatial distribution of vegetation coverage change trend and the sustainable change grades in the Pisha sandstone area

表4 不同砒砂巖區植被覆蓋度空間變化趨勢面積占比統計/%Table 4 Statistics on the spatial variation trend of vegetation coverage in different Pisha sandstone areas

砒砂巖區1999—2018年植被覆蓋度的Hurst指數介于0.23-0.99之間,植被覆蓋度呈正向持續性變化(Hurst>0.5)的區域占總面積的56.2%,反向持續變化(Hurst<0.5)面積比重為43.8%。說明砒砂巖區植被覆蓋度變化將在短期內保持現有的發展趨勢,即未來植被覆蓋度的變化呈現正向趨勢。Hurst指數按照砒砂巖分區進行統計可知(圖4),砒砂巖區植被持續性變化趨勢主要以持續性改善為主,面積占比50.5%,主要分布在覆土區南部,植被由改善向退化方向變化的面積約占41.6%,主要分布在覆土區東北部以及裸露區西北部,植被由退化向改善方向變化的面積約2.4%,零星分布在裸露區中部,其余5.6%的植被呈持續性退化趨勢。

3.3 植被覆蓋度變化驅動因子量化

為了探究不同環境因子對砒砂巖區植被覆蓋度空間分布的影響,選取8種因子并對其進行地理探測器分析,得到影響因子的解釋力q值(見表5)。解釋力從大到小依次為降水、土壤水分、氣溫、海拔、土壤類型、土地利用類型、坡度和坡向(P<0.01)。降水、土壤水分和氣溫的q值均大于0.5,是區域內植被覆蓋空間分布的主導環境因子。土壤類型的q值為0.333,與海拔q值(0.319)接近,即對植被覆蓋空間分布的解釋力中等,坡度和坡向對植被覆蓋的解釋力最弱,q值分別為0.054和0.003。

表5 1999—2018年砒砂巖區環境因子對植被覆蓋的解釋力Table 5 The explanatory power of environmental factors in the Pisha sandstone area on vegetation coverage from 1999 to 2018

雙因子交互作用下會加強對植被覆蓋度的解釋力(圖5),坡向同其他因子交互時均為非線性增強趨勢,其余各因子交互均為雙因子增強趨勢。其中,降水∩土壤水分(q值為0.651)、氣溫∩土壤水分(q值為0.647)、土壤類型∩土壤水分(q值為0.641)和降水∩土壤類型(q值為0.628)對砒砂巖區植被覆蓋空間分布的解釋力較大。由此可以看出,降水作為砒砂巖區植被覆蓋空間分布的主導氣候類因子,在同其他環境因子的交互作用下,對植被覆蓋度影響最大。

圖5 1999—2018年砒砂巖區各環境因子交互探測解釋力(q)以及生態檢測Fig.5 The interactive detection explanatory power (q)and ecological detection of various environmental factors in the Pisha sandstone from 1999 to 2018*雙因子間無顯著差異

土壤類型單因素對該區域植被覆蓋度的解釋力較弱,但同土壤水分的交互作用解釋力較大(q值為0.641),表明土壤類型在滿足一定的土壤水分時對植被覆蓋度產生顯著的影響。此外,為檢測雙因子間對植被覆蓋度空間分布是否存在顯著差異,對各因子做生態探測,結果表明:土壤水分和氣溫、土壤類型和海拔對植被覆蓋度空間分布無顯著差異(檢測值為N),其余各因子間均有顯著性差異(P<0.05),說明土壤水分和氣溫對植被空間分布具有相似的機理,二者的解釋力值也較相近。

土地利用、坡度以及坡向對植被覆蓋度的解釋力不超過0.1,即使在交互作用下,對區域植被覆蓋度的影響也較低,因此本研究選取年均降水量、年均氣溫、土壤類型、土壤水分、高程輸入風險探測器以確定各環境因子促進植被生長的適宜范圍或類型,其中植被覆蓋度均值越大,相對應的環境因子范圍或類型更適合植被生長。對于氣候因子來說,年均降水量劃分為9個分區,其中各分區的年均降水量范圍分別是315—343、343—364、364—379、379—390、390—400、400—409、409—419、419—429、429—442 mm。在年均降水量>429 mm時,植被覆蓋度均值最大,其中>429 mm的區域占砒砂巖區域面積9.9%,且有統計檢驗表明,該分區(429—442 mm)同其他分區具有顯著性差異(P<0.05)。同樣,年均氣溫也劃分為9個分區,數值越大表明溫度越高,在9分區(9.1—9.6℃),植被覆蓋度均值最大為58.0%。

不同土壤類型的植被覆蓋度呈波動變化(表6),其中潛育雛形土、飽和薄層土以及石灰性雛形土更適宜植被生長,植被覆蓋度均值最大分別為57.1%、56.0%、52.0%。土壤水分與植被覆蓋度呈正相關,最適宜的分區為9分區,即滿足土壤水分在0.23-0.26 m3/m3范圍下,植被覆蓋度可達55.1%。對于不同的高程范圍,植被覆蓋度均值隨著高程的增加而降低,在1分區(774—956 m)達到最大值51.0%。

表6 不同環境因子的適宜范圍或類型Table 6 Suitable ranges or types of different environmental factors

4 討論

植被對環境變化相當敏感,尤其是在半干旱地區,氣候是對植被生長具有決定性影響的自然驅動力[23]。全球半干旱地區的植被綠化受降水控制的區域約50%,受氣溫控制的約7%[24]。砒砂巖區降水量稀少,且區域雨熱同期,土壤水分主要靠降水量補給[25]。肖強[26]等研究發現黃土高原西北部的植被覆蓋度與降水呈顯著正相關,而與溫度呈顯著負相關。He[27]等發現黃土高原90.12%的區域降水和NDVI呈正相關,且通過顯著性檢驗的區域主要分布在黃土高原北部地區。黃土高原地處干旱半干旱區,降水是限制該區域植被生長的主要因子,具體表現為植被覆蓋度隨著降水的增加而增大[28]。本研究也進一步證實了降水是主要的影響因素,在砒砂巖區不同的區域,植被覆蓋度同降水以及氣溫的偏相關性方向以及程度均有不同(圖6),其中降水同植被覆蓋度呈正相關的區域面積占96.8%,且69.7%的面積顯著正相關(P<0.05),年均氣溫同植被覆蓋度呈負相關的區域面積占比68.5%,僅3.9%的區域面積通過顯著性檢驗。溫度升高會導致地面水分蒸發增加,以至于出現土壤干化現象,對植被的生長不利[29]。在半干旱地區,氣溫主要是調節了植被的年內生長規律,同降水相比來說,氣溫與植被的年際變化相關性一般較小,降水對植被覆蓋的影響更加敏感,是西北干旱區植被變化的關鍵因子[30]。

圖6 植被覆蓋度同降水與氣溫的偏相關系數空間分布Fig.6 Spatial distribution of partial correlation coefficient between vegetation coverage and precipitation and temperature

砒砂巖區1999—2018年植被生長狀況整體呈增加趨勢,該結論同之前的研究結果是一致的[31]。而本論文中進一步指出,該區植被總體雖呈現向好趨勢,但是仍舊面臨較大的挑戰。植被這種向好態勢一方面與地區氣象環境關聯較密,一方面也與區域治理措施相關。為恢復黃土高原生態環境狀況,自1999年以來實施退耕還林工程,其中約2×106hm2的耕地變為草地和森林,該區近80%的區域植被覆蓋率增加[32]。砒砂巖區地形破碎,生態環境惡劣,針對于該區域的治理,國家先后實施了沙棘資源建設,國家生態工程建設等水土保持項目[33]。雖然取得了一定成效,對當地生態環境的改善具有重要意義,但受該地區自身雨熱條件限制,這些措施的作用仍然是有限的,砒砂巖地區的生態惡化和水土流失問題依然沒有得到完全有效的解決,區域植被建設進入了一定瓶頸期。同時也有研究表明陜甘寧地區15年來植被雖有一定程度的恢復,但由于能源設施以及礦產資源的開發,北部地區植被退化的風險越來越大,植被覆蓋度未來退化趨勢面積占比達58.7%[34-35]。因此,研究該區植被覆蓋特征并量化環境驅動因子的影響,能夠針對覆土、覆沙、裸露砒砂巖區不同環境背景條件下合理優化植被格局,對于維系該區植被建設成果,實現區域生態保護同資源利用協調發展具有重要意義。

5 結論

本研究分析了近20年砒砂巖區的植被覆蓋度變化及未來變化趨勢,并利用地理探測器定量分析了不同環境因子對植被覆蓋度的影響,具體結論如下:

(1)1999—2018年砒砂巖區植被覆蓋度呈增加趨勢,平均上升幅度為0.086/10a,多年平均植被覆蓋度為42.3%。空間上,砒砂巖區植被覆蓋度呈現從東南向西北遞減的空間分布特征,其三個區的植被覆蓋度從大到小變化依次是覆土區、覆沙區、裸露區。

(2)近20年砒砂巖區植被覆蓋度得到改善的比重較大,其中45.5%的區域呈極顯著增加,主要分布在覆土區和覆沙區,顯著和極顯著減少的區域零星分布在裸露區中部。該區植被覆蓋度未來變化趨勢以持續性改善為主,面積占比50.5%,主要分布在覆土區南部,雖然區域植被總體呈向好趨勢,但是仍舊面臨較大挑戰,約41.6%的植被存在由改善向退化轉變的風險。

(3)降水、土壤水分和氣溫是影響區域內植被覆蓋空間分布的主導環境因子,坡度和坡向對植被覆蓋的解釋力最弱;環境因子間的交互作用以雙因子增強趨勢為主,其中降水∩土壤水分解釋力最高為0.651,各因子對砒砂巖區植被生長的影響均有其適宜的范圍/類型,研究結果為該區域下一步植被恢復工作提供相應科學依據。

猜你喜歡
區域
分割區域
探尋區域創新的密碼
科學(2020年5期)2020-11-26 08:19:22
基于BM3D的復雜紋理區域圖像去噪
軟件(2020年3期)2020-04-20 01:45:18
小區域、大發展
商周刊(2018年15期)2018-07-27 01:41:20
論“戎”的活動區域
敦煌學輯刊(2018年1期)2018-07-09 05:46:42
區域發展篇
區域經濟
關于四色猜想
分區域
公司治理與技術創新:分區域比較
主站蜘蛛池模板: 亚洲精品视频免费观看| 中文字幕va| 免费无码一区二区| 九九精品在线观看| 69av免费视频| 亚洲美女操| 九九热在线视频| 亚洲欧美国产高清va在线播放| 亚洲中文制服丝袜欧美精品| 国产女人18水真多毛片18精品 | 亚洲Av激情网五月天| 国产精品一区二区无码免费看片| 沈阳少妇高潮在线| 久久久久免费精品国产| 国产人前露出系列视频| 亚洲欧州色色免费AV| 99热线精品大全在线观看| 久久人人97超碰人人澡爱香蕉 | 久久一本精品久久久ー99| 四虎精品免费久久| 重口调教一区二区视频| 国产大片喷水在线在线视频| 黄色国产在线| 国产自在自线午夜精品视频| 另类综合视频| 欧美a级完整在线观看| 久久成人18免费| 亚洲黄色视频在线观看一区| 亚洲资源站av无码网址| 欧美亚洲第一页| 国产精品视频导航| 久久夜色撩人精品国产| 国产视频自拍一区| 欧美亚洲国产精品第一页| 麻豆精品在线视频| 国产精品浪潮Av| 日本尹人综合香蕉在线观看| 手机精品视频在线观看免费| 久久美女精品| 青青网在线国产| 欧美亚洲日韩不卡在线在线观看| 国产网友愉拍精品| 在线日韩一区二区| 成年人久久黄色网站| 97精品国产高清久久久久蜜芽 | 人妻精品久久无码区| 午夜福利在线观看入口| 亚洲天堂免费| 国产久操视频| 2022国产91精品久久久久久| 国产精品人莉莉成在线播放| 啦啦啦网站在线观看a毛片| 波多野结衣视频一区二区| 欧美日韩在线亚洲国产人| 国产精品污视频| 国产日本欧美在线观看| 成人在线观看不卡| 国产JIZzJIzz视频全部免费| 国产成人精品一区二区秒拍1o| 日本高清免费一本在线观看| 免费a级毛片视频| 中国精品自拍| 全部免费特黄特色大片视频| 四虎在线高清无码| 污网站免费在线观看| 日本人又色又爽的视频| 亚洲一级毛片| 国产精品自在在线午夜| 伊人欧美在线| 久久免费精品琪琪| 国产中文在线亚洲精品官网| 日韩欧美中文字幕一本| 日韩精品一区二区深田咏美| 成人欧美在线观看| 青青国产成人免费精品视频| V一区无码内射国产| 亚洲日韩国产精品无码专区| 国产一区二区福利| 激情综合婷婷丁香五月尤物| 999国内精品视频免费| 日韩国产 在线| 亚洲日韩欧美在线观看|