亚洲免费av电影一区二区三区,日韩爱爱视频,51精品视频一区二区三区,91视频爱爱,日韩欧美在线播放视频,中文字幕少妇AV,亚洲电影中文字幕,久久久久亚洲av成人网址,久久综合视频网站,国产在线不卡免费播放

        ?

        基于SPCA和遙感指數(shù)的干旱內(nèi)陸河流域生態(tài)脆弱性時(shí)空演變及動(dòng)因分析
        ——以石羊河流域?yàn)槔?/h1>
        2019-05-13 07:16:28郭澤呈龐素菲李振亞周俊菊頡斌斌
        生態(tài)學(xué)報(bào) 2019年7期
        關(guān)鍵詞:石羊河脆弱性流域

        郭澤呈,魏 偉,*,龐素菲,李振亞,周俊菊,頡斌斌

        1 西北師范大學(xué)地理與環(huán)境科學(xué)學(xué)院, 蘭州 730070 2 蘭州城市學(xué)院城市經(jīng)濟(jì)與旅游文化學(xué)院, 蘭州 730070

        生態(tài)脆弱性是指在特定的時(shí)間和空間尺度上,相對(duì)于外界擾動(dòng)影響具有感知性和恢復(fù)力,是生態(tài)系統(tǒng)自身固有的屬性,受外界干擾因素的影響下表現(xiàn)出來(lái)[1]。近年來(lái),隨著社會(huì)經(jīng)濟(jì)的高速發(fā)展,人類社會(huì)活動(dòng)同生態(tài)環(huán)境之間的矛盾日益突出,人類對(duì)自然資源的不合理利用,工農(nóng)業(yè)發(fā)展帶來(lái)的污染等給周邊生態(tài)平衡造成破壞[2],有關(guān)生態(tài)脆弱性以及脆弱性評(píng)估的研究逐漸成為全球性研究熱點(diǎn),受到國(guó)內(nèi)外眾多學(xué)者的關(guān)注[3- 4]。尤其在國(guó)內(nèi),眾多學(xué)者從不同區(qū)域、不同尺度和不同方法的角度展開(kāi)了大量的研究工作,如王讓會(huì)等[1]從景觀生態(tài)學(xué)角度,選取10個(gè)綜合性敏感因子,采用生態(tài)脆弱性指數(shù)評(píng)價(jià)法對(duì)塔里木河流域生態(tài)敏感性及恢復(fù)力進(jìn)行了分析;陳群利等[5]從自然、社會(huì)和水土流失景觀格局三方面著手,選取13個(gè)評(píng)價(jià)指標(biāo)構(gòu)建生態(tài)脆弱性評(píng)價(jià)指標(biāo)體系,采用集對(duì)分析(SPA)方法對(duì)貴州省畢節(jié)地區(qū)生態(tài)脆弱性進(jìn)行了分析;徐慶勇等[6]利用遙感技術(shù),結(jié)合空間主成分分析和層次分析法,對(duì)晉北地區(qū)進(jìn)行生態(tài)脆弱性評(píng)價(jià),為晉北地區(qū)合理調(diào)控人類活動(dòng),保護(hù)和治理生態(tài)環(huán)境提供一定的依據(jù);王志杰等[7]則基于“壓力-狀態(tài)-響應(yīng)”評(píng)價(jià)模型框架,利用空間主成分分析對(duì)南水北調(diào)中線漢中市水源地生態(tài)脆弱性進(jìn)行了定量評(píng)價(jià)。從上述研究成果來(lái)看,區(qū)域生態(tài)脆弱性評(píng)價(jià)不僅從宏觀層面上了解區(qū)域的生態(tài)脆弱現(xiàn)狀以及分布特征,而且有助于生態(tài)治理措施的合理實(shí)施和脆弱生態(tài)環(huán)境的長(zhǎng)期保護(hù)。但從評(píng)價(jià)的指標(biāo)體系數(shù)據(jù)源來(lái)看,多以氣象、遙感,數(shù)字高程模型等數(shù)據(jù)為主,部分研究也涉及了遙感數(shù)據(jù)與地理空間數(shù)據(jù)相結(jié)合的方法計(jì)算生態(tài)脆弱性,但充分利用遙感模型和遙感指數(shù)從流域?qū)用嫔?尤其在干旱地區(qū)利用遙感指標(biāo)對(duì)其進(jìn)行深入分析,這部分成果相對(duì)較少。此外,利用傳統(tǒng)方法,由于主觀因素影響,在構(gòu)建指標(biāo)體系時(shí)將選取的指標(biāo)全部納入評(píng)價(jià)模型中,從而忽視評(píng)價(jià)指標(biāo)間是否具有高相關(guān)性,造成評(píng)價(jià)結(jié)果的不確定性,其可靠度并不高[8]。因此,從目前關(guān)于生態(tài)脆弱性的綜合評(píng)價(jià)來(lái)看,借助遙感數(shù)據(jù)構(gòu)建多指標(biāo)綜合評(píng)價(jià)模型,不但能夠客觀、快速的分析其現(xiàn)狀空間分布特征,而且亦可分析其時(shí)空變化規(guī)律,其結(jié)果不受主觀影響,可信度較高。

        石羊河流域地處黃土、青藏和蒙新三大高原的交匯過(guò)渡帶,是生態(tài)脆弱和環(huán)境變化敏感的流域之一[9],今年來(lái)生態(tài)問(wèn)題日益突出,成為我國(guó)干旱區(qū)內(nèi)陸河流域生態(tài)退化的典型區(qū)域[10],受到了社會(huì)各界的廣泛關(guān)注?;诖?本文重點(diǎn)在于以典型的西北干旱內(nèi)陸河流域——石羊河流域?yàn)檠芯繉?duì)象,以徐涵秋提出的新型遙感生態(tài)指數(shù)(RSEI)[11]為基礎(chǔ),選取其基于遙感的綠度、濕度、熱度和干度指標(biāo)作為研究流域的評(píng)價(jià)指標(biāo),采用空間主成分分析便捷、快速地對(duì)石羊河流域2000—2016年的生態(tài)脆弱性進(jìn)行定量的評(píng)價(jià),以可視化形式展現(xiàn)該流域生態(tài)脆弱性的時(shí)空分布特征;同時(shí)利用這些指標(biāo)來(lái)研究影響石羊河流域生態(tài)脆弱性演變的重要因素,用以揭示該流域生態(tài)脆弱性時(shí)空演變規(guī)律,旨在為該流域的生態(tài)脆弱性保護(hù)及治理提供理論基礎(chǔ)和決策依據(jù)。

        1 研究區(qū)與數(shù)據(jù)源

        1.1 研究區(qū)概況

        石羊河流域位于甘肅省河西走廊東部,烏鞘嶺以西,祁連山北麓,東南與甘肅省白銀、蘭州兩市相接,西北與甘肅省張掖市毗鄰,西南緊靠青海省,東北與內(nèi)蒙古自治區(qū)接壤,屬河西走廊三大內(nèi)陸河流域之一,介于東經(jīng)101°41′—104°16′,北緯36°29′—39°27′之間。流域的行政區(qū)劃包括張掖市肅南裕固族自治縣和山丹縣的部分地區(qū),武威市的古浪縣、涼州區(qū)、民勤縣全部及天??h部分,金昌市的永昌縣及金川區(qū)全部以及白銀市景泰縣的少部分地區(qū),共涉及4市9縣,流域總面積約40578平方公里。石羊河流域深居大陸腹地,海拔1000—5000m,屬大陸性溫帶干旱氣候,具有太陽(yáng)輻射較強(qiáng)、日照充足、溫差大、降水少、蒸發(fā)強(qiáng)烈以及空氣干燥的氣候特點(diǎn)。流域內(nèi)土壤、植被類型分布因受氣候、水文和地形等自然條件的影響,形成了明顯的土壤-植被垂直帶譜[12]。

        1.2 數(shù)據(jù)源與預(yù)處理

        論文采用的數(shù)據(jù)包括:數(shù)字高程模型(DEM)數(shù)據(jù)和遙感影像數(shù)據(jù)。數(shù)字高程模型(DEM)數(shù)據(jù)的空間分辨率為30m×30m,用于提取石羊河流域的海拔等信息;遙感影像數(shù)據(jù)為地理空間數(shù)據(jù)云網(wǎng)站(http://www.gscloud.cn/)提供的2000年Landsat TM影像和2016年Landsat OLI影像,空間分辨率均為30m。在數(shù)據(jù)處理時(shí)利用ENVI5.1軟件分別對(duì)各時(shí)期每幅影像進(jìn)行輻射定標(biāo),將各幅影像的像元值轉(zhuǎn)換成傳感器處的反射率;使用FLAASH絕對(duì)大氣校正對(duì)各幅影像進(jìn)行校正,以減少不同時(shí)期影像在光照、大氣等條件產(chǎn)生的差異,并分別對(duì)不同時(shí)期影像進(jìn)行無(wú)縫鑲嵌,囊括研究區(qū)域的所有范圍;最后通過(guò)石羊河流域邊界對(duì)所有影像進(jìn)行裁剪,并采用二次多項(xiàng)式和最鄰近像元法對(duì)不同時(shí)相影像間配準(zhǔn),配準(zhǔn)的均方根誤差應(yīng)小于0.5個(gè)像元內(nèi),滿足本次研究對(duì)數(shù)據(jù)精度的要求。

        2 研究方法

        2.1 評(píng)價(jià)指標(biāo)計(jì)算

        2.1.1濕度指標(biāo)(Wet)

        纓帽變換的亮度、綠度、濕度分量與生態(tài)質(zhì)量密切相關(guān)[13- 14]。其中,濕度分量可反映研究區(qū)域土壤和植被的濕度狀況。因此,本文采用濕度分量來(lái)代表濕度指標(biāo)?;贚andsat TM和Landsat OLI反射率影像數(shù)據(jù)的公式分別為[15-16]:

        WetTM=0.0315ρB+0.2021ρG+0.3102ρR+0.1594ρNIR-0.6806ρSWIR1-0.6109ρSWIR2

        (1)

        WetOLI=0.1511ρB+0.1972ρG+0.3283ρR+0.3407ρNIR-0.7117ρSWIR1-0.4559ρSWIR2

        (2)

        式中:ρB、ρG、ρR、ρNIR、ρSWIR1、ρSWIR2分別為TM和OLI的藍(lán)、綠、紅、近紅外、短波紅外1、短波紅外2波段的反射率數(shù)據(jù)。

        2.1.2綠度指標(biāo)(NDVI)

        植被是反映區(qū)域生態(tài)質(zhì)量好壞極其重要的因素。歸一化植被指數(shù)(Normalized difference vegetation index, NDVI)是應(yīng)用最廣泛的植被指數(shù),能良好地反映植物生物量、葉面積指數(shù)與植被覆蓋度之間的關(guān)系[17]。因此,可選用NDVI作為本文的綠度指標(biāo),其公式為:

        NDVI=(ρNIR-ρR)/(ρNIR+ρR)

        (3)

        式中:ρNIR、ρNIR分別為TM5與OLI8的紅和近紅外波段的反射率數(shù)據(jù)。

        2.1.3干度指標(biāo)(NDSI)

        研究流域的城市用地會(huì)對(duì)周邊生態(tài)產(chǎn)生一定的影響,造成地表的“干化”,對(duì)生態(tài)產(chǎn)生負(fù)面影響,這里選用建筑指數(shù)(IBI)來(lái)表示[18]; 而研究流域內(nèi)有大面積的荒漠區(qū),其大面積的裸露土地是整個(gè)研究流域的重要的生態(tài)弊病,這里選用裸土指數(shù)(SI)來(lái)表達(dá)研究流域的裸露狀態(tài)[19]。最終干度指標(biāo)是由建筑指數(shù)(IBI)和裸土指數(shù)(SI)合成,記為NDSI[11],其計(jì)算公式為:

        IBI={2ρSWIR1/(ρSWIR1+ρNIR)-[ρNIR/(ρNIR+ρR)+ρG/(ρG+ρSWIR1)]}/{2ρSWIR1/(ρSWIR1+ρNIR)+

        [ρNIR/(ρNIR+ρR)+ρG/(ρG+ρSWIR1)]}

        (4)

        SI=[(ρSWIR1+ρR)-(ρB+ρNIR)/(ρSWIR1+ρR)+(ρB+ρNIR)]

        (5)

        NDSI=(IBI+SI)/2

        (6)

        式中:ρG、ρB、ρR、ρNIR、ρSWIR1分別為TM和OLI綠、藍(lán)、紅、近紅外和短波紅外1波段的反射率數(shù)據(jù)。

        2.1.4熱度指標(biāo)(LST)

        熱度指標(biāo)采用地表溫度[11]表達(dá),地表溫度與植被和水資源關(guān)系密切。本文選取Landsat用戶手冊(cè)的模型[20-21]計(jì)算出亮溫,再通過(guò)地表比輻射率[22]校正后,獲得地表溫度,其公式為:

        L=gain×DN+bias

        (7)

        Tb=K2/ln(K1/L+1)

        (8)

        LST=Tb/[1+(λTb/ρ)lnε]

        (9)

        式中:L為TM熱紅外6波段和OLI熱紅外10波段輻射亮度值;DN為TM熱紅外6波段和OLI熱紅外10波段像元值;gain為熱紅外波段增益值,bias為熱紅外波段偏置值,而gain和bias值可從影像的頭文件中獲取[23];Tb為傳感器處溫度值,即亮溫;K1和K2為傳感器處的定標(biāo)參數(shù)。其中,TM處,K1=607.76W/(m2·sr·μm),K2=1260.56K,OLI處,K1=774.89W/(m2·sr·μm),K2=1321.08K;λ為熱紅外波段波長(zhǎng);ρ=1.4380×104μm K;ε為地表比輻射率。

        2.2 評(píng)價(jià)指標(biāo)的標(biāo)準(zhǔn)化

        各評(píng)價(jià)指標(biāo)具有不同的性質(zhì),其量綱不同,無(wú)法直接進(jìn)行生態(tài)脆弱性的評(píng)價(jià),因此,必須對(duì)各評(píng)價(jià)指標(biāo)進(jìn)行標(biāo)準(zhǔn)化處理,以解決指標(biāo)間無(wú)法直接比較的矛盾。本文選取的四個(gè)指標(biāo)根據(jù)對(duì)生態(tài)脆弱性的貢獻(xiàn)可分為正向指標(biāo)和逆向指標(biāo)。正向指標(biāo)代表指標(biāo)值越大,其生態(tài)脆弱性越高; 逆向指標(biāo)則反之。其中,正向指標(biāo)包括干度和熱度,逆向指標(biāo)包括綠度和濕度。指標(biāo)標(biāo)準(zhǔn)化則采用極差標(biāo)準(zhǔn)化的方法,其公式為[7-8]:

        正向指標(biāo):

        (10)

        逆向指標(biāo):

        (11)

        式中:SIi表示第i指標(biāo)的標(biāo)準(zhǔn)化值,其值域范圍為0—10;Ii為第i指標(biāo)的實(shí)際值;Imax為第i指標(biāo)的最大值;Imin為第i指標(biāo)的最小值。

        2.3 生態(tài)脆弱性綜合評(píng)價(jià)模型

        在眾多的生態(tài)脆弱性評(píng)價(jià)方法中,空間主成分分析法(SPCA)是通過(guò)對(duì)特征光譜空間坐標(biāo)軸的旋轉(zhuǎn),將相關(guān)的多變量空間數(shù)據(jù)轉(zhuǎn)化為少數(shù)幾個(gè)不相關(guān)的綜合指標(biāo),實(shí)現(xiàn)用較少的綜合指標(biāo)最大限度的保留原來(lái)較多變量所反映的信息。當(dāng)累計(jì)方差貢獻(xiàn)率大于或等于85%,就能代表絕大多數(shù)的有關(guān)信息[24]。采用空間主成分分析法并不需要人為的確定各個(gè)指標(biāo)的權(quán)重,可以避免因人為因素而造成最終結(jié)果的偏差。本文在ENVI 5.1軟件平臺(tái)上,將評(píng)價(jià)指標(biāo)體系中標(biāo)準(zhǔn)化后的綠度、濕度、干度和熱度四個(gè)評(píng)價(jià)指標(biāo)進(jìn)行空間主成分分析,計(jì)算出生態(tài)脆弱性評(píng)價(jià)指數(shù)(EVI)。其計(jì)算公式為[7-8]:

        EVI=r1Y1+r2Y2+r3Y3+···+rnYn

        (12)

        式中:EVI為生態(tài)脆弱性指數(shù);ri為第i個(gè)主成分;Yi為第i個(gè)主成分的貢獻(xiàn)率。

        為進(jìn)行生態(tài)脆弱性指數(shù)(EVI)的對(duì)比和度量,因此對(duì)EVI進(jìn)行標(biāo)準(zhǔn)化處理,其處理公式為:

        (13)

        式中:SIEVI為生態(tài)脆弱性指數(shù)的標(biāo)準(zhǔn)化值,值域范圍為0—10;EVIi為生態(tài)脆弱性指數(shù)實(shí)際的值;EVImax和EVImin為生態(tài)脆弱性指數(shù)的最大值和最小值。

        根據(jù)空間主成分分析原理,由公式(12)和表2得出石羊河流域生態(tài)脆弱性評(píng)價(jià)反演模型:

        EVI2000=0.8854×PC1+0.0859×PC2+0.0220×PC3

        (14)

        EVI2016=0.8592×PC1+0.1162×PC2+0.0209×PC3

        (15)

        式中:EVI2000和EVI2016分別為2000年和2016年的生態(tài)脆弱性指數(shù);PC1—PC3為原始空間變量進(jìn)行主成分變換后的前3個(gè)主成分因子。2000年和2016年前3個(gè)主成分因子的累計(jì)貢獻(xiàn)率均達(dá)到99%,而兩年的第4主成分因子經(jīng)對(duì)比后發(fā)現(xiàn),其表述的信息與實(shí)際生態(tài)脆弱性現(xiàn)狀并不符合,大部分信息為噪聲,因此可以忽略不計(jì)。

        表1 空間主成分分析結(jié)果

        2.4 生態(tài)脆弱性分級(jí)及生態(tài)脆弱性整體指數(shù)

        通過(guò)ArcGIS 10.2 中的Raster Calculator工具獲取研究流域2000和2016年生態(tài)脆弱性指數(shù)的空間分布。根據(jù)現(xiàn)有的生態(tài)脆弱性評(píng)價(jià)的分級(jí)標(biāo)準(zhǔn)[7-8],結(jié)合研究流域的流域特征,將研究流域生態(tài)脆弱性(EVI)分成5個(gè)等級(jí): I:微度脆弱(0—2), II:輕度脆弱(2—4),II:中度脆弱(4—6),IV:強(qiáng)度脆弱(6—8),V:重度脆弱(8—10)。同時(shí)為了研究不同年份在不同空間單元下的生態(tài)脆弱性狀況的整體差異狀況,本文采用生態(tài)脆弱性整體指數(shù)(Ecological vulnerability body index, EVBI)進(jìn)行分析,其計(jì)算公式[8]為:

        (16)

        式中:EVBI為生態(tài)脆弱性差異指數(shù),Pi為第i類生態(tài)脆弱等級(jí)值,這里為1—5;Ai為第i類脆弱性等級(jí)的面積;S為研究區(qū)域總面積。

        2.5 生態(tài)脆弱性時(shí)空演變格局提取

        利用ArcGIS 10.2 中的Raster Calculator工具進(jìn)行運(yùn)算,將2000年和2016年的生態(tài)脆弱性分級(jí)圖進(jìn)行空間疊加,提取生態(tài)脆弱性變化動(dòng)態(tài)圖斑[25],運(yùn)算公式為:

        CodeClassification_2000—2016=10×CodeClassification_2000+CodeClassification_2016

        (17)

        式中:CodeClassification_2000—2016為生態(tài)脆弱性分級(jí)變化類型代碼;CodeClassification_2000和CodeClassification_2016分別為5種生態(tài)脆弱性等級(jí)類型代碼(1—5),其中1—5分別代表微度、輕度、中度、強(qiáng)度和重度脆弱。因此,在運(yùn)算后的生態(tài)脆弱性分級(jí)變化類型代碼中,十位數(shù)為2000年生態(tài)脆弱性分級(jí)類型,而個(gè)位數(shù)為2016年生態(tài)脆弱性分級(jí)類型,CodeClassification_2000—2016所代表的含義為2000年生態(tài)脆弱性分級(jí)類型轉(zhuǎn)變成2016年生態(tài)脆弱性分級(jí)類型(如12指2000年的微度脆弱類型轉(zhuǎn)變?yōu)?016年的輕度脆弱類型)。

        2.6 生態(tài)脆弱性不同海拔梯度下的空間分布及時(shí)空演變格局提取

        石羊河流域地勢(shì)南高北低,自西南向東北傾斜,地貌類型多樣,從南到北,依次經(jīng)過(guò)祁連山地、平原走廊、低山丘陵和荒漠區(qū),相對(duì)高差達(dá)3900m以上[12]。不同海拔梯度上,其溫度特征、植被覆蓋、水源涵養(yǎng)以及地表裸露等都具有明顯的差異,生態(tài)脆弱程度也存在著差異。因此,有必要研究石羊河流域在不同海拔梯度下的生態(tài)脆弱性垂直地帶性差異。

        基于石羊河流域海拔分布特征,依據(jù)陳志明等[26]根據(jù)國(guó)家DTM數(shù)據(jù)進(jìn)行高程頻數(shù)統(tǒng)計(jì)得到的地貌類型劃分標(biāo)準(zhǔn),即可將石羊河流域地貌類型劃分成3種:中山(海拔1000—2000m)、高中山(海拔2000—3000m)和高山(海拔>3000m)。將石羊河流域生態(tài)脆弱性分級(jí)評(píng)價(jià)結(jié)果與海拔分級(jí)圖疊加[25],即可獲得不同海拔梯度下各脆弱性等級(jí)的空間分布狀況。其運(yùn)算公式為:

        EEVI=10×ElevationClassification+EVIClassification

        (18)

        式中:EEVI為不同海拔梯度下各脆弱性等級(jí)分布; ElevationClassification為海拔分級(jí)結(jié)果;EVIClassification為生態(tài)脆弱性分級(jí)結(jié)果。

        為研究2000年和2016年兩期在不同海拔梯度下生態(tài)脆弱性分級(jí)變化狀況,利用ArcGIS 10.2 中的Raster Calculator工具進(jìn)行運(yùn)算,將不同海拔梯度下2000年和2016年生態(tài)脆弱性等級(jí)分布圖進(jìn)行空間疊加[25],其運(yùn)算公式為:

        CodeEEVI_2000—2016=100×CodeEEVI_2000+CodeEEVI_2016

        (19)

        式中: CodeEEVI_2000—2016為不同海拔梯度下生態(tài)脆弱性分級(jí)類型變化代碼,CodeEEVI_2000和CodeEEVI_2016分別為不同海拔梯度下15種生態(tài)脆弱性分級(jí)原代碼(11—35)。其中,千位數(shù)和十位數(shù)代表3種海拔分級(jí)代碼(1—3分別代表中山、高中山和高山;百位數(shù)和個(gè)位數(shù)為生態(tài)脆弱性等級(jí)類型代碼(1—5分別代表微度、輕度、中度、強(qiáng)度和重度脆弱)。CodeEEVI_2000—2016代表不同海拔梯度下2000年生態(tài)脆弱性類型轉(zhuǎn)變成2016年生態(tài)脆弱性類型 (如1112指地貌類型為中山,2000年的微度脆弱類型轉(zhuǎn)變成2016年的輕度脆弱類型)。

        3 結(jié)果與分析

        3.1 生態(tài)脆弱性指標(biāo)變化特征

        表2是各年份4個(gè)指標(biāo)的統(tǒng)計(jì)值。統(tǒng)計(jì)結(jié)果表明,17年間,石羊河流域的濕度和綠度指標(biāo)呈上升的趨勢(shì),均值分別從2000年的-0.29、0.22上升到2016年的-0.15、0.27,增幅分別為48.3%、22.7%,證明該流域水源涵養(yǎng)能力變好;植被覆蓋也呈逐年增加的趨勢(shì);干度指標(biāo)有所下降,均值由2000年的0.09下降到2016年的0.06,減幅為33.3%,表明該流域地表裸露程度有所降低;而與植被和水資源關(guān)系密切的地表溫度呈逐年上升趨勢(shì),其均值從2000年的37.40上升到2016年的39.33,增幅為0.05%,說(shuō)明該流域水熱平衡差異進(jìn)一步增加,對(duì)未來(lái)生態(tài)脆弱性影響顯著。從生態(tài)脆弱性指標(biāo)空間分布圖可大致看出其空間分布(圖1),17年間,人工綠洲及周邊地區(qū)的植被覆蓋和土壤濕度狀況有所好轉(zhuǎn),而祁連山區(qū)有一定的萎縮現(xiàn)象;地表裸露狀況在城鎮(zhèn)和祁連山區(qū)有所擴(kuò)張,而荒漠區(qū)地表裸露狀況有所減緩;流域的地表溫度在整體空間分布上呈逐年上升的趨勢(shì)。

        3.2 石羊河流域生態(tài)脆弱性整體特征

        通過(guò)對(duì)上述指標(biāo)進(jìn)行生態(tài)脆弱性評(píng)價(jià)模型的計(jì)算,反演出不同年份生態(tài)脆弱性等級(jí)空間分布,以此來(lái)進(jìn)行研究流域生態(tài)脆弱性的時(shí)空演變。從整體上看,該流域的生態(tài)脆弱性整體指數(shù)(EVBI)從2000年的3.47下降到2016年的3.39,增幅不大。整個(gè)流域生態(tài)脆弱性以中度脆弱和強(qiáng)度脆弱為主。圖2是2000年和2016年各生態(tài)脆弱等級(jí)面積分布狀況。2016年,石羊河流域微度、輕度、中度、強(qiáng)度和重度脆弱區(qū)的面積分別為2493.90 km2、5845.84 km2、6496.70 km2、24661.61km2和1080.79 km2,微度、輕度脆弱區(qū)面積比2000年分別增加815.13 km2和61.35 km2;強(qiáng)度脆弱區(qū)面積比2000年增加2962.50 km2;中度、重度脆弱區(qū)面積比2000年分別減少1615.67 km2、2223.30 km2。在空間分布上(圖3), 2000年和2016年的石羊河流域大部分區(qū)域生態(tài)脆弱性處于強(qiáng)度脆弱水平,人工綠洲以及祁連山地區(qū)多為中度脆弱及以下水平,而重度脆弱區(qū)主要位于中東部的荒漠區(qū)。這主要因?yàn)槭蚝恿饔蛏罹哟箨懜沟?氣候干燥,水資源匱乏,植被覆蓋較少,大部分區(qū)域處于沙漠化狀態(tài),其生態(tài)脆弱程度較高。為進(jìn)一步探究石羊河流域生態(tài)脆弱等級(jí)17年間隨時(shí)間推移的空間變化狀況,由公式(15)計(jì)算獲得生態(tài)脆弱性時(shí)空演變格局(圖4),統(tǒng)計(jì)出17年間不同生態(tài)脆弱等級(jí)面積轉(zhuǎn)移矩陣(表3)。從表3中可知,17年間各類生態(tài)脆弱性等級(jí)面積轉(zhuǎn)移的總和為10117.30 km2。2000年,微度脆弱主要轉(zhuǎn)向輕度脆弱,轉(zhuǎn)移面積為582.05 km2;輕度脆弱主要轉(zhuǎn)向微度和中度脆弱,轉(zhuǎn)移面積分別為1023.61 km2和1275.20 km2;中度脆弱主要轉(zhuǎn)向輕度和強(qiáng)度脆弱,轉(zhuǎn)移面積分別為1505.46 km2和2070.34 km2;強(qiáng)度脆弱主要轉(zhuǎn)向中度脆弱,轉(zhuǎn)移面積為932.26 km2;重度脆弱主要轉(zhuǎn)向強(qiáng)度脆弱,轉(zhuǎn)移面積為2218.63 km2;2016年新增的微度脆弱面積主要由輕度脆弱轉(zhuǎn)化而來(lái),轉(zhuǎn)化面積占微度脆弱新增面積的70.2%;輕度脆弱面積主要由中度脆弱轉(zhuǎn)化而來(lái),轉(zhuǎn)化面積占比為60.8%;中度脆弱面積主要由輕度脆弱面積轉(zhuǎn)化而來(lái),轉(zhuǎn)化面積占比為55.6%;強(qiáng)度脆弱面積主要由重度脆弱面積轉(zhuǎn)化而來(lái),轉(zhuǎn)化面積占比為50.3%;重度脆弱面積主要由強(qiáng)度脆弱轉(zhuǎn)化而來(lái),轉(zhuǎn)化面積占比為99.5%。17年間,微度、輕度脆弱區(qū)的新增面積分別為1458.31km2和2474.88 km2。在空間變化上,人工綠洲以及周邊地區(qū)生態(tài)脆弱性向低脆弱性轉(zhuǎn)移,生態(tài)有所改善,沙漠化治理取得一定成效[27];但是,因城鎮(zhèn)擴(kuò)張、祁連山區(qū)人類過(guò)度開(kāi)采資源等導(dǎo)致城鎮(zhèn)及祁連山地區(qū)生態(tài)有所萎縮,生態(tài)脆弱性向高脆弱性轉(zhuǎn)移。研究區(qū)內(nèi)生態(tài)脆弱性總體上降低面積大于增高面積,從側(cè)面反映出該流域生態(tài)往良好方向發(fā)展。

        表2 各年份4個(gè)指標(biāo)的統(tǒng)計(jì)值

        圖2 石羊河流域生態(tài)脆弱性等級(jí)面積分布 Fig.2 Area distribution of ecological vulnerability levels in Shiyang River Basin

        圖3 石羊河流域生態(tài)脆弱性等級(jí)空間分布Fig.3 Spatial distribution of ecological vulnerability levels in Shiyang River Basin

        圖4 石羊河流域生態(tài)脆弱性時(shí)空演變格局Fig.4 Spatiotemporal evolution pattern of EVI in Shiyang River Basin

        脆弱性等級(jí)Vulnerability level2016微度脆弱Slight vulnerability 輕度脆弱Lightvulnerability中度脆弱Moderatevulnerability強(qiáng)度脆弱Strongvulnerability重度脆弱Heavyvulnerability2000年總計(jì)Total of 20002000年轉(zhuǎn)移量Reduction of 20002000微度脆弱Slight vulnerability1035.58 582.05 49.81 11.33 0.00 1678.77 643.19 輕度脆弱Light vulnerability1023.61 3370.96 1275.20 114.66 0.06 5784.49 2413.53 中度脆弱 Moderate vulnerability332.91 1505.46 4203.33 2070.34 0.34 8112.37 3909.05 強(qiáng)度脆弱Strong vulnerability90.51 354.33 932.26 20246.64 75.37 21699.11 1452.47 重度脆弱Heavy vulnerability11.29 33.04 36.11 2218.63 1005.02 3304.09 2299.07 2016年總計(jì)Total of 20062493.90 5845.84 6496.70 24661.61 1080.79 40578.83 -2016年新增量Increment of 20161458.31 2474.88 2293.37 4414.97 75.77 -10117.30

        3.3 石羊河流域不同海拔梯度上生態(tài)脆弱性時(shí)空演變特征

        圖5 不同海拔梯度下生態(tài)脆弱性等級(jí)面積分布Fig.5 Area distribution of ecological vulnerability levels in different elevations

        圖5和圖6為17年間石羊河流域在不同海拔梯度上生態(tài)脆弱性等級(jí)的面積及空間分布情況。該流域生態(tài)脆弱性最高的海拔梯度為中山區(qū)(1000—2000m),2000和2016年生態(tài)脆弱性整體指數(shù)(EVBI)分別為3.88和3.72,高于各年份EVBI值的整體狀況(2000年為3.47,2016年為3.39),生態(tài)脆弱性等級(jí)空間分布以強(qiáng)度脆弱為主(面積占比分別為71.13%和77.69%),主要分布在大面積的荒漠區(qū);高中山區(qū)(2000—3000m)次之,2000年和2016年的EVBI值分別為2.98和3.04,生態(tài)脆弱性等級(jí)空間分布以中度脆弱和強(qiáng)度脆弱為主(面積占比分別為53.80%、24.52%和37.76%、37.15%),這是由于在該海拔梯度下的祁連山山區(qū)受人類活動(dòng)的干擾,地表裸露程度有所增大,抗外界干擾能力和自我恢復(fù)能力有所下降,生態(tài)脆弱性有所增高;而高山區(qū)(>3000m)生態(tài)脆弱性等級(jí)較低,2000年和2016年的EVBI值分別為1.88和2.03,生態(tài)脆弱性等級(jí)空間分布以輕度脆弱為主(面積占比分別為60.54%和52.95%),明顯小于其他海拔梯度和整體的生態(tài)脆弱水平。

        圖6 不同海拔梯度下生態(tài)脆弱性等級(jí)空間分布Fig.6 Spatial distribution of EVI in different elevations

        為獲取在不同海拔梯度下生態(tài)脆弱性等級(jí)的變化情況,由公式(19)得到17年間不同海拔梯度下生態(tài)脆弱性等級(jí)變化的空間分布(圖7),統(tǒng)計(jì)其面積轉(zhuǎn)移矩陣(表4)。統(tǒng)計(jì)結(jié)果表明,中山區(qū)呈現(xiàn)生態(tài)脆弱性從較高等級(jí)向較低等級(jí)轉(zhuǎn)換的趨勢(shì),以微度、輕度脆弱增長(zhǎng),重度脆弱降低為主,2016年微度、輕度脆弱新增面積(分別為747.20 km2和1293.07 km2)大于2000年轉(zhuǎn)移面積(分別為87.96 km2和913.11 km2),而2016年重度脆弱新增面積(75.22 km2)小于2000年轉(zhuǎn)移面積(2299.92 km2);高中山區(qū)和高山區(qū)則呈現(xiàn)生態(tài)脆弱性等級(jí)向高脆弱等級(jí)方向發(fā)展的趨勢(shì),高中山區(qū)以強(qiáng)度脆弱增長(zhǎng)為主,2016年強(qiáng)度脆弱新增面積(1390.81 km2)大于2000年轉(zhuǎn)移面積(275.78 km2);高山區(qū)以中度脆弱增長(zhǎng),微度、輕度脆弱降低為主。2016年中度脆弱新增面積(613.50 km2)大于2000年轉(zhuǎn)移面積(199.71 km2);2016年微度、輕度脆弱新增面積(311.67 km2和558.12 km2)小于2000年轉(zhuǎn)移面積(447.57 km2和886.35 km2)。從空間分布上也可得到相同的規(guī)律,中山區(qū)生態(tài)脆弱性降低的主要原因在于綠洲擴(kuò)張和沙漠化的有效治理,而其中部分地區(qū)由于城鎮(zhèn)化進(jìn)程加快導(dǎo)致輕度、中度脆弱往中度、強(qiáng)度脆弱轉(zhuǎn)移;高中山區(qū)和高山區(qū)生態(tài)脆弱性升高與人類過(guò)度的開(kāi)發(fā)礦產(chǎn)資源、部分水電設(shè)施違法建設(shè)和運(yùn)行等因素有關(guān)。

        圖7中變化類型為無(wú)變化,其為三種地貌類型在2000—2016年各生態(tài)脆弱性等級(jí)未發(fā)生變化的總和,面積為29865.02km2

        3.4 石羊河流域各縣區(qū)生態(tài)脆弱性差異

        圖8 石羊河流域各縣區(qū)生態(tài)脆弱性整體指數(shù)Fig.8 EVBI of each county in Shiyang River Basin

        根據(jù)石羊河流域各縣區(qū)的生態(tài)脆弱性整體指數(shù)(EVBI)進(jìn)一步比較2000年和2016年生態(tài)脆弱性的空間差異性(圖8)。結(jié)果表明,2016年較2000年相比,金川區(qū)、涼州區(qū)、永昌縣、民勤縣和古浪縣的EVBI值呈下降趨勢(shì),其中民勤縣的降幅最大,從2000年的3.95下降到2016年的3.83,下降了0.12,說(shuō)明從2001年起的一系列生態(tài)治理的政策和措施使得民勤縣生態(tài)脆弱性有所降低,植被覆蓋增多,水源涵養(yǎng)變好,沙漠化程度得到一定的遏止。但兩個(gè)年份中金川區(qū)、涼州區(qū)、永昌縣、民勤縣和古浪縣的EVBI值均大于3,生態(tài)脆弱性處于中度偏高水平,除了各縣區(qū)內(nèi)有一定的荒漠區(qū)影響生態(tài)脆弱水平外,各縣區(qū)建設(shè)用地增多、人口流動(dòng)大、人類活動(dòng)頻繁等也是導(dǎo)致一定程度上區(qū)域內(nèi)生態(tài)脆弱性增高的原因。而流域內(nèi)的天??h和肅南縣生態(tài)脆弱性處于輕度偏低水平,EVBI值呈上升趨勢(shì),分別由2000年的1.96、2.25上升到2016年的2.05、2.43。天??h和肅南縣位于祁連山區(qū),與人類過(guò)多的干擾活動(dòng)有關(guān),導(dǎo)致植被覆蓋,水源涵養(yǎng)能力和礦產(chǎn)資源儲(chǔ)存等受到一定的影響,生態(tài)問(wèn)題日漸突出,因此今后應(yīng)需加強(qiáng)生態(tài)脆弱性的治理和恢復(fù)。

        4 討論

        4.1 生態(tài)脆弱性評(píng)價(jià)指標(biāo)體系與評(píng)價(jià)方法選取的合理性

        由于多光譜影像具有多波段和相關(guān)性較強(qiáng)的特征,導(dǎo)致多光譜影像內(nèi)部信息的冗余度較高[28]。因此基于遙感影像計(jì)算的指數(shù)間可能存在信息冗余,不僅計(jì)算量增多,而且會(huì)直接影響評(píng)價(jià)結(jié)果的精確性。本文在構(gòu)建指標(biāo)體系時(shí),想在不受人為因素和主觀條件約束情況下,通過(guò)遙感模型和指數(shù),客觀和快速的評(píng)價(jià)研究流域在17年間生態(tài)脆弱性變化及空間演變規(guī)律,經(jīng)過(guò)慎重篩選最終選擇了徐涵秋提出的新型遙感生態(tài)指數(shù)(RSEI)的指標(biāo)體系[11],其可以良好的反映研究流域植被覆蓋、土壤濕度、地表裸露和地表溫度狀況,而且該指標(biāo)體系間不存在較明顯的相關(guān)性[29];而本文所選用的空間主成分分析法本身就具有去除各指標(biāo)間一定的相關(guān)性,降低數(shù)據(jù)冗余的作用[11,30]。為驗(yàn)證各指標(biāo)的相關(guān)性,本文以2016年為例,對(duì)各指標(biāo)進(jìn)行共線性診斷[8]。常用的共線性診斷指標(biāo)主要有兩個(gè):方差膨脹因子(VIF)和容忍度(TOL)。這兩個(gè)指標(biāo)互為倒數(shù),當(dāng)VIF>10 (即TOL<0.1) 時(shí),表明所選指標(biāo)的多元共線性較為嚴(yán)重。在ArcGIS 10.2中,采用3km×3km的格網(wǎng)貫穿全影像的方法,共均勻生成研究區(qū)4523個(gè)點(diǎn),然后分別利用這些點(diǎn)讀取各指標(biāo)與EVI的值,利用SPSS 21.0計(jì)算出各指標(biāo)的VIF與TOL(表5)。 從計(jì)算結(jié)果可以看出,各指標(biāo)的VIF均小于10,TOL均大于0.1,表明各指標(biāo)間不存在明顯的相關(guān)性。綜上所述,本文選取綠度、濕度、干度和熱度作為評(píng)價(jià)指標(biāo)體系,空間主成分分析法作為評(píng)價(jià)方法是可取的。

        表5 多元共線性診斷結(jié)果

        4.2 石羊河流域生態(tài)脆弱性的演變動(dòng)因

        為進(jìn)一步探究石羊河流域生態(tài)脆弱性的時(shí)空演變規(guī)律,本文引入地理探測(cè)器[31]作為分析工具,用以診斷出生態(tài)脆弱性的主導(dǎo)影響因素。地理探測(cè)器包括4個(gè)探測(cè)器:風(fēng)險(xiǎn)探測(cè)器、因子探測(cè)器、生態(tài)探測(cè)器與交互作用探測(cè)器。而利用其中的因子探測(cè)器可探測(cè)某因子是否是形成生態(tài)脆弱性時(shí)空分布格局的原因以及在多大程度上解釋了生態(tài)脆弱性的空間分異機(jī)理[31]。具體做法是以生態(tài)脆弱性指數(shù)(EVI)作為因變量,將選取的4個(gè)指標(biāo)作為自變量因子,將自變量進(jìn)行分層,由數(shù)值量轉(zhuǎn)為類型量,2000和2016年各指標(biāo)均采用自然斷點(diǎn)法分成5類,代表不同程度的生態(tài)脆弱性類型;然后在ArcGIS 10.2中,采用3km×3km的格網(wǎng)貫穿全影像的方法,共均勻生成研究區(qū)4523個(gè)點(diǎn),將因變量值和自變量值過(guò)格網(wǎng)點(diǎn)匹配起來(lái),進(jìn)行因子探測(cè)分析[32],得出各因子對(duì)生態(tài)脆弱性指數(shù)(EVI)的影響力值(q值,q值越大表示該因子對(duì)生態(tài)脆弱性指數(shù)(EVI)的影響越大)和因子解釋力值 (p值,p值越大表示該因子對(duì)生態(tài)脆弱性指數(shù)(EVI)的解釋力越小)[33]。分析結(jié)果顯示(表6): 2000和2016年的p值均為0,4個(gè)因子對(duì)石羊河流域生態(tài)脆弱性的解釋力都很充足;q值均大于0.5,4個(gè)因子對(duì)石羊河流域生態(tài)脆弱性的影響均為顯著。從生態(tài)脆弱性隨時(shí)間的演變過(guò)程來(lái)看,2016年較2000年相比,4個(gè)因子的q值均有不同程度的變化,說(shuō)明17年間,4個(gè)因子對(duì)石羊河流域生態(tài)脆弱性演變的影響程度存在波動(dòng)變化的趨勢(shì);從生態(tài)脆弱性主導(dǎo)影響因子角度來(lái)看,2016年較2000年相比,在q排序中,干度較其他因子對(duì)生態(tài)脆弱性變化的影響程度有所提升,而濕度有所下降。其主要原因是由于石羊河流域濕度含量的增長(zhǎng),促使植被覆蓋增多,植被對(duì)流域生態(tài)脆弱性的影響變大;植被覆蓋的多少在一定程度上抑制地表溫度對(duì)生態(tài)脆弱性的影響,因此地表溫度的q值有所減??;干度q值變大主要受城市擴(kuò)張和祁連山生態(tài)破壞導(dǎo)致在部分空間分布上地表裸露程度變大的影響。針對(duì)影響因子存在的復(fù)雜耦合關(guān)系對(duì)石羊河流域生態(tài)脆弱性在不同的地理單元上產(chǎn)生的明顯空間分異特征,當(dāng)?shù)卣畱?yīng)從政策制度上完善健全,落實(shí)各方責(zé)任,協(xié)調(diào)當(dāng)?shù)鼐用窠?jīng)濟(jì)發(fā)展與生態(tài)脆弱性治理的平衡,采取正確、適當(dāng)?shù)拇胧┡c方法抑制生態(tài)脆弱性向更為脆弱的方向發(fā)展,加快石羊河流域生態(tài)治理的步伐,鞏固生態(tài)治理建設(shè)的成果,才能從根源上去治理石羊河流域脆弱的生態(tài)環(huán)境。

        4.3 研究結(jié)果的不足

        考慮到生態(tài)脆弱性評(píng)價(jià)具有的綜合性、復(fù)雜性和不明確性等特點(diǎn),而目前的評(píng)價(jià)方法尚不能做到全面、科學(xué)和客觀的評(píng)價(jià)。因此,本研究也只是在總結(jié)前人研究成果的基礎(chǔ)上,以干旱內(nèi)陸河流域獨(dú)特的自然條件為依據(jù),“綠”、“濕”作為干旱內(nèi)陸河流域生態(tài)脆弱性的主要決定因素,而“干”、“熱”為干旱內(nèi)陸河流域主要生態(tài)脆弱特征,利用這四個(gè)方面結(jié)合遙感和評(píng)價(jià)模型對(duì)石羊河流域生態(tài)脆弱性進(jìn)行宏觀、快速和客觀評(píng)價(jià),揭示其生態(tài)脆弱性的演變動(dòng)因,不僅為石羊河流域生態(tài)脆弱性治理提供一定的參考價(jià)值,而且對(duì)干旱內(nèi)陸河流域生態(tài)脆弱性評(píng)價(jià)提供一定的思路和借鑒。本文的主要目的是通過(guò)遙感指數(shù)思想,試圖探索出一種利用客觀評(píng)價(jià)模型對(duì)干旱內(nèi)陸地區(qū)生態(tài)脆弱性評(píng)價(jià)的方法。但在指標(biāo)體系的確定上,本文仍存在一些不足:干旱內(nèi)陸河流域生態(tài)脆弱性變化是一個(gè)十分復(fù)雜的問(wèn)題,涉及自然、生態(tài)、社會(huì)經(jīng)濟(jì)、人類活動(dòng)等各個(gè)方面,不可能用1個(gè)或幾個(gè)指標(biāo)來(lái)完全表征,而本研究只考慮了自然因素對(duì)石羊河流域生態(tài)脆弱性的影響。要想更加科學(xué)、全面的反映石羊河流域生態(tài)脆弱性時(shí)空變化規(guī)律仍需進(jìn)一步的探索與研究。

        表6 石羊河流域4個(gè)影響因子的地理探測(cè)結(jié)果

        5 結(jié)論

        本文基于流域的生態(tài)特征和遙感的快速、客觀、宏觀和大面積觀測(cè)特點(diǎn),選取濕度、綠度、干度和熱度等4個(gè)指標(biāo)構(gòu)建生態(tài)脆弱性的評(píng)價(jià)指標(biāo)體系,運(yùn)用空間主成分分析法(SPCA),對(duì)石羊河流域生態(tài)脆弱性狀況和時(shí)空特征進(jìn)行分析。得到以下結(jié)論:

        (1) 從各遙感指數(shù)空間分布來(lái)看,濕度和綠度指標(biāo)在17年間呈增加趨勢(shì),證明該流域水源涵養(yǎng)能力變好,植被覆蓋率變大;干度指標(biāo)值有所下降,表明該流域地表裸露程度有所降低;而與植被和水資源關(guān)系密切的地表溫度呈逐年上升趨勢(shì),說(shuō)明該流域水熱平衡差異進(jìn)一步增加,對(duì)未來(lái)生態(tài)脆弱性影響顯著;

        (2) 從全流域生態(tài)脆弱性時(shí)空演變特征來(lái)看,該區(qū)域主要以強(qiáng)度和中度脆弱為主,高脆弱區(qū)主要分布在荒漠和城鎮(zhèn)區(qū),低脆弱區(qū)主要位于人工綠洲和祁連山地區(qū)。17年間生態(tài)脆弱性整體上呈緩慢降低趨勢(shì);

        (3) 從不同的海拔生態(tài)脆弱性分布來(lái)看,中山區(qū)(1000—2000m)最高,以強(qiáng)度脆弱為主;高中山區(qū)(2000—3000m)次之,以中度和強(qiáng)度脆弱為主;高山區(qū)(>3000m)最低,以輕度脆弱為主。17年間中山區(qū)生態(tài)脆弱性有所下降,而高中山與高山區(qū)卻呈上升的趨勢(shì);

        (4) 從不同的行政區(qū)劃生態(tài)脆弱性分布來(lái)看,金川區(qū)、涼州區(qū)、永昌縣、民勤縣和古浪縣處于中度和強(qiáng)度脆弱水平,17年間生態(tài)脆弱性呈逐漸下降的趨勢(shì);而天祝縣和肅南縣整體處于輕度和微度脆弱水平,17年間生態(tài)脆弱性呈逐漸上升趨勢(shì)。

        (5) 從生態(tài)脆弱性的演變動(dòng)因來(lái)看,4個(gè)指標(biāo)對(duì)石羊河流域生態(tài)脆弱性影響均為顯著。2000年生態(tài)脆弱性的主導(dǎo)影響因子依次為熱度>濕度>綠度>干度;2016年生態(tài)脆弱性的主導(dǎo)影響因子依次為熱度>干度>綠度>濕度。

        猜你喜歡
        石羊河脆弱性流域
        人民黃河(2023年7期)2023-08-27 15:41:53
        壓油溝小流域
        基于不同旱情指數(shù)的石羊河流域春旱監(jiān)測(cè)研究
        堡子溝流域綜合治理
        羅堰小流域
        石羊河流域永昌縣地下水及水資源供需平衡分析
        水利規(guī)劃與設(shè)計(jì)(2018年1期)2018-01-31 01:53:37
        煤礦電網(wǎng)脆弱性評(píng)估
        電子制作(2017年10期)2017-04-18 07:23:09
        殺毒軟件中指令虛擬機(jī)的脆弱性分析
        基于攻擊圖的工控系統(tǒng)脆弱性量化方法

        国产熟女精品一区二区三区| 可以直接在线看国产在线片网址| 久久国内精品自在自线| 亚洲中字慕日产2020| 中文字幕精品一区二区2021年| 91在线在线啪永久地址| 丝袜美女美腿一区二区| 国产亚洲成人精品久久| 国产午夜精品一区二区| 国产亚洲蜜芽精品久久| 日本一区二区在线播放观看| 人妻熟女翘屁股中文字幕| 丰满熟妇乱又伦精品| 少妇做爰免费视频网站| 亚洲成人免费无码| 久久精品国产亚洲av专区| 日本亚洲系列中文字幕| 波多野结衣久久精品99e| 人人妻人人澡人人爽人人精品| 欧美1区二区三区公司| 日本岛国一区二区三区四区| 日本做受120秒免费视频| 狠狠色狠狠色综合久久第一次| 精品视频在线观看一区二区有| 日韩美女亚洲性一区二区| 国产乱码精品一区二区三区四川人| 国产美女在线精品亚洲二区| 蜜桃一区二区三区在线视频| 国产特级毛片aaaaaa高潮流水| 色婷婷久久一区二区三区麻豆| 国产性一交一乱一伦一色一情| 亚洲第一页在线观看视频网站| 日韩人妻无码精品一专区二区三区| 中文字幕亚洲欧美日韩2019| 色综合色综合久久综合频道| av在线播放中文专区| 波多野结衣的av一区二区三区| 国产一区二区精品在线观看| 精品国产污黄网站在线观看| 亚洲熟女综合色一区二区三区| 国产成人精品一区二区不卡|