• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    基于TRMM 3B43數(shù)據(jù)的川西高原月降水量空間降尺度模擬*

    2016-05-27 02:58:50馮文蘭成都信息工程大學(xué)資源環(huán)境學(xué)院成都610225
    中國(guó)農(nóng)業(yè)氣象 2016年2期
    關(guān)鍵詞:滯后性多元線性回歸

    鄭 杰,閭 利,馮文蘭,涂 坤(成都信息工程大學(xué)資源環(huán)境學(xué)院,成都 610225)

    ?

    基于TRMM 3B43數(shù)據(jù)的川西高原月降水量空間降尺度模擬*

    鄭杰,閭利,馮文蘭**,涂坤
    (成都信息工程大學(xué)資源環(huán)境學(xué)院,成都 610225)

    摘要:利用2001-2013年TRMM 3B43、MODIS-NDVI、DEM、氣象觀測(cè)等數(shù)據(jù),在分析植被對(duì)降水響應(yīng)滯后性的基礎(chǔ)上,構(gòu)建了TRMM 3B43數(shù)據(jù)中月降水量與經(jīng)緯度、海拔、坡向和NDVI因子間的多元線性回歸方程式,作為川西高原月降水量資料的降尺度計(jì)算模型,采用“回歸方程+殘差”的插值方法獲取研究區(qū)2001-2013年1km空間分辨率的月降水量空間數(shù)據(jù),并利用區(qū)內(nèi)16個(gè)氣象站點(diǎn)的觀測(cè)數(shù)據(jù)與模擬結(jié)果進(jìn)行了相關(guān)分析和誤差檢驗(yàn)。結(jié)果表明:(1)各氣象觀測(cè)站點(diǎn)基于TRMM 3B43資料的降尺度模擬降水量的數(shù)據(jù)均具有很高的精度,其中,精度最高的稻城站模擬結(jié)果與站點(diǎn)觀測(cè)值的相關(guān)系數(shù)高達(dá)0.9839,精度最低的小金站相關(guān)系數(shù)亦高達(dá)0.8781;(2)在月、年尺度上,降尺度模擬降水量的數(shù)據(jù)亦具有很高的精度,其中,5-10月的精度明顯高于其它月份,濕潤(rùn)年份精度總體高于干旱年份;(3)降尺度模擬降水量與站點(diǎn)實(shí)測(cè)降水量整體上相關(guān)系數(shù)為0.9499,偏差為0.0866,兩者吻合度較高,但降尺度模擬降水量值略偏高;(4)降尺度在月尺度上能基本保證TRMM 3B43原始數(shù)據(jù)的精度,而在年尺度上能有效提高原始數(shù)據(jù)的精度,加之對(duì)空間分辨率的提高,可為獲得更加全面、精細(xì)的降水分布數(shù)據(jù)提供有效方法。

    關(guān)鍵詞:TRMM3B43數(shù)據(jù);降水?dāng)?shù)據(jù);空間降尺度;滯后性;多元線性回歸

    鄭杰,閭利,馮文蘭,等.基于TRMM 3B43數(shù)據(jù)的川西高原月降水量空間降尺度模擬[J].中國(guó)農(nóng)業(yè)氣象,2016,37(2):245-254

    降水是全球水循環(huán)的基本組成物質(zhì),具有顯著的時(shí)空變異性,影響一個(gè)區(qū)域的水分和熱量分布,因此,獲取精確度較高的降水?dāng)?shù)據(jù)對(duì)了解降水量的時(shí)空變化特征和掌握區(qū)域能量平衡起著重要作用[1-3]。

    傳統(tǒng)的降水資料主要由地面氣象站直接觀測(cè)獲取,并通過(guò)實(shí)測(cè)點(diǎn)插值計(jì)算出研究區(qū)的空間降水?dāng)?shù)據(jù),然而,站點(diǎn)觀測(cè)數(shù)據(jù)僅能反映站點(diǎn)周邊一定范圍內(nèi)地區(qū)的降水量,很難有效反映整個(gè)區(qū)域的空間降水信息[4]。近年來(lái),隨著RS和GIS技術(shù)的快速發(fā)展,通過(guò)衛(wèi)星遙感數(shù)據(jù)對(duì)降水進(jìn)行探測(cè)、反演,覆蓋范圍廣、空間分辨率較高且受氣候和地形限制較小的衛(wèi)星降水觀測(cè)數(shù)據(jù)得以實(shí)現(xiàn)[5-6]。熱帶降水測(cè)量衛(wèi)星(Tropical Rainfall Measuring Mission,TRMM)于1997年成功發(fā)射,是由美國(guó)宇航局(NASA)和日本國(guó)家空間發(fā)展局(JAXA)共同研制的第一顆專門用于定量測(cè)量熱帶和亞熱帶地區(qū)降水的氣象衛(wèi)星,由于融合了多個(gè)衛(wèi)星的數(shù)據(jù)以及地面觀測(cè)數(shù)據(jù),該數(shù)據(jù)空間分辨率和數(shù)據(jù)質(zhì)量均較高,且能保證降水?dāng)?shù)據(jù)的均一性[7-8]。然而,盡管TRMM降水?dāng)?shù)據(jù)的空間分辨率在同類衛(wèi)星降水?dāng)?shù)據(jù)中最高,但仍難以滿足流域的水文和氣候模型等研究對(duì)較高空間分辨率降水?dāng)?shù)據(jù)的需求,因此,進(jìn)一步提高TRMM數(shù)據(jù)的空間分辨率及其精度非常必要。

    目前,國(guó)內(nèi)外已有學(xué)者開展了氣象觀測(cè)站點(diǎn)數(shù)據(jù)的多元回歸統(tǒng)計(jì)插值研究以及TRMM 3B43數(shù)據(jù)的空間降尺度研究。周鎖銓等[9]表明,可用多元線性回歸方程表示降水量和高程、坡度、坡向等地形信息之間的關(guān)系,并在此基礎(chǔ)上建立了降水量和地形因子關(guān)系的模型,獲取了長(zhǎng)江中上游地區(qū)4km分辨率的空間降水信息;馬金輝等[10]通過(guò)建立高程、坡度、坡向等地形因子和TRMM降水之間的多元線性回歸方程,將TRMM 3B43數(shù)據(jù)的空間分辨率由原來(lái)的0.25°提高至1km。然而,降水是一個(gè)復(fù)雜的天氣過(guò)程,且近年來(lái),人類活動(dòng)對(duì)地表過(guò)程的影響顯著,影響氣候的不確定因素增多[11]。植被是連結(jié)氣候和人文因素的重要紐帶[12],因而,Jia等[13-15]將降水量與歸一化植被指數(shù)(NDVI)的相關(guān)性引入降水?dāng)?shù)據(jù)的空間插值以及降尺度研究中,既不完全依賴于地形因子與降水?dāng)?shù)據(jù)的相關(guān)性,又對(duì)降水?dāng)?shù)據(jù)的空間變異信息進(jìn)行了補(bǔ)充,可減小空間插值和降尺度的誤差。嵇濤等[16]綜合考慮了地形和植被因子,分別從不同空間尺度(0.25°、0.5°、0.75°和1°)對(duì)TRMM 3B43數(shù)據(jù)進(jìn)行回歸模擬,結(jié)果指出,川渝地區(qū)0.25°分辨率下降水?dāng)?shù)據(jù)的空間降尺度效果最好,且能進(jìn)一步提高降水估算的精度。然而,目前已有的成果多針對(duì)TRMM 3B43的年降水?dāng)?shù)據(jù)進(jìn)行空間降尺度研究,對(duì)月尺度降水?dāng)?shù)據(jù)的空間降尺度研究相對(duì)較少。

    川西高原地區(qū)是長(zhǎng)江上游重要的生態(tài)屏障,是地理、氣象學(xué)者研究的重點(diǎn)區(qū)域。鑒于該區(qū)域氣象站點(diǎn)稀少且分布不均勻,本研究充分借鑒已有研究成果的思路,進(jìn)一步考慮植被生長(zhǎng)對(duì)降水響應(yīng)的滯后性特征,結(jié)合地形和植被因子構(gòu)建TRMM 3B43降水?dāng)?shù)據(jù)與經(jīng)緯度、高程、坡向和NDVI的多元線性回歸模型,開展TRMM 3B43月降水?dāng)?shù)據(jù)在地形復(fù)雜地區(qū)的空間降尺度研究,以期為地形復(fù)雜的山區(qū)獲取高時(shí)空分辨率的降水?dāng)?shù)據(jù)提供有效方法。

    1 資料與方法

    1.1研究區(qū)概況

    川西高原地區(qū)(97 °26′-104°27′E,27°57′-34°21′N)處于橫斷山脈東段,青藏高原的東部邊緣,包括甘孜藏族自治州和阿壩藏族羌族自治州,面積約20萬(wàn)km2。平均海拔3900m,地勢(shì)高差懸殊,是四川省海拔最高的地區(qū),總體上呈現(xiàn)西部高東部低的態(tài)勢(shì)(圖1a)。

    川西高原屬高山高原高寒氣候區(qū),以寒溫帶氣候?yàn)橹?,全區(qū)氣候垂直變化顯著,從河谷到山脊依次可分為亞熱帶、暖溫帶、中溫帶、寒溫帶、亞寒帶、寒帶和永凍帶。區(qū)內(nèi)氣候特征的水平差異也很顯著,河谷干暖,山地冷濕。全區(qū)冬寒夏涼,太陽(yáng)輻射強(qiáng)烈,光能資源豐富,年均溫4~12℃,極端最低氣溫-20℃以下,年降水量500~900mm,多集中在5-10月[17]。研究區(qū)地形地貌復(fù)雜,山脈縱橫,地形起伏大(圖1b),坡向等地形因子對(duì)植被覆蓋和降水等的空間分布有較大影響。

    圖1 研究區(qū)高程和氣象站點(diǎn)分布(a)及坡向(b)Fig. 1 Elevation and distribution of meteorological stations(a),aspect(b) in study area

    1.2數(shù)據(jù)來(lái)源及預(yù)處理

    研究使用的TRMM 3B43數(shù)據(jù)為TRMM第6版3級(jí)產(chǎn)品(V6-3B43)的月降水資料(http://precip. gsfc.nasa.gov/),空間分辨率為0.25°×0.25°,時(shí)間分辨率為1個(gè)月。利用IDL程序讀取HDF格式的2001-2013年的TRMM 3B43數(shù)據(jù),并將其轉(zhuǎn)化為矢量格式的點(diǎn)數(shù)據(jù)。

    MODIS NDVI數(shù)據(jù)來(lái)源于NASA Terra衛(wèi)星提供的MODIS 13A2 級(jí)產(chǎn)品中16d最大值合成NDVI數(shù)據(jù)集(http://www.gscloud.cn/),數(shù)據(jù)空間分辨率為1km×1km,時(shí)間覆蓋范圍為2001年1月-2014年1月,該數(shù)據(jù)在制備過(guò)程中經(jīng)過(guò)輻射校正、大氣校正、幾何校正等預(yù)處理。為反映植被變化對(duì)降水的響應(yīng),對(duì)該數(shù)據(jù)進(jìn)行投影轉(zhuǎn)換和拼接等處理,并采用最大值合成法(MVC)獲取月NDVI值。

    SRTM-DEM數(shù)據(jù)來(lái)源于地理空間數(shù)據(jù)云平臺(tái)(http://www.gscloud.cn/),該數(shù)據(jù)由美國(guó)宇航局(NASA)和國(guó)防部國(guó)家測(cè)繪局(NIMA)聯(lián)合測(cè)量,空間分辨率為90m×90m。經(jīng)檢查,研究區(qū)內(nèi)沒(méi)有數(shù)據(jù)空洞。為建立地形因子與降水的相關(guān)關(guān)系,首先對(duì)DEM數(shù)據(jù)進(jìn)行投影轉(zhuǎn)換和裁剪,再采用最臨近插值法將其重采樣至空間分辨率為1km×1km。

    氣象資料來(lái)源于中國(guó)地面氣候資料月值數(shù)據(jù)集(http://cdc.nmic.cn/home.do),包括甘孜藏族自治州和阿壩藏族羌族自治州境內(nèi)的16個(gè)氣象站點(diǎn)(圖1a)2001-2013年的月降水量數(shù)據(jù)。為檢驗(yàn)研究結(jié)果的精度,利用ArcGIS軟件生成氣象站點(diǎn)降水量的矢量點(diǎn)圖層。

    1.3研究方法

    1.3.1植被對(duì)降水的滯后性分析

    為了說(shuō)明植被NDVI因子在降尺度模型構(gòu)建中的可行性,本文引入植被NDVI的局部Moran’s I指數(shù),該指數(shù)用來(lái)反映NDVI分布的空間自相關(guān)性,即揭示當(dāng)前像元與其鄰近的空間像元NDVI值的相似性,其計(jì)算式為

    式中,Ii為局部Moran’s I指數(shù),xi和xj分別為空間像元i和j的NDVI值,wij為權(quán)重矩陣,n為像元總數(shù);若Moran’s I指數(shù)負(fù)值太大,說(shuō)明植被NDVI會(huì)影響降尺度模型的精度,需排除這些點(diǎn)的NDVI值。本文運(yùn)用OpenGeoDa軟件計(jì)算局部Moran’s I指數(shù)。

    植被生長(zhǎng)受降水的影響較大,本文在降尺度模型中考慮了植被因子,但降水從到地表被植被吸收并反映到植被的生長(zhǎng)狀況需要一定的時(shí)間[18],即當(dāng)前月份的植被NDVI值與當(dāng)月降水的關(guān)系并不很密切,在當(dāng)月降水量降尺度模型中選擇當(dāng)月NDVI值作為植被因子其模擬效果并不理想,因此,本研究通過(guò)對(duì)比當(dāng)月、前一月和前兩月降水變化與植被NDVI的相關(guān)關(guān)系,揭示植被NDVI對(duì)降水響應(yīng)的滯后性。

    1.3.2降水?dāng)?shù)據(jù)的空間統(tǒng)計(jì)降尺度方法

    1.3.2.1降尺度原理

    統(tǒng)計(jì)降尺度方法是指采用統(tǒng)計(jì)經(jīng)驗(yàn)方法確立不同空間尺度影像之間某一特征量的線性和非線性函數(shù)關(guān)系,從而針對(duì)柵格影像進(jìn)行空間尺度的轉(zhuǎn)換,是低分辨率氣候資料轉(zhuǎn)換為高分辨率氣候資料的重要手段之一[19]。本文在統(tǒng)計(jì)降尺度方法的基礎(chǔ)上,結(jié)合經(jīng)緯度、高程、坡向和植被因子,構(gòu)建適于該地區(qū)的TRMM月降水量多元線性回歸模型,然后加上回歸預(yù)測(cè)的殘差值,得到降尺度后的降水空間分布數(shù)據(jù),其模型為

    P= F(X,Y,H,V,A)+D TRMMHR(2)

    式中,F(xiàn)(X,Y,H,V,A)表示由地形和地理因子線性回歸得到的降水預(yù)測(cè)值,X為經(jīng)度(o),Y為緯度(o),H為高程(m),V為植被因子,A為坡向因子,ΔTRMMHR為由地形、地理和植被因子引起的降水分布的殘差(mm)。馬金輝等[10]在建立降水量與地形因子之間的回歸方程時(shí),采用坡向的余弦值來(lái)計(jì)算降水量,取得了較理想的結(jié)果,因此,本文在月降水量降尺度模型的建立時(shí)也采用坡向的余弦值來(lái)量化坡向因子。

    1.3.2.2降尺度步驟

    (1)對(duì)研究區(qū)邊界向外做25km的緩沖區(qū)(約向外增加一個(gè)TRMM像元),以緩沖區(qū)邊界裁取研究區(qū)內(nèi)2001-2013年TRMM 3B43月降水量格點(diǎn)數(shù)據(jù)。

    (2)在0.25°分辨率下提取每一個(gè)對(duì)應(yīng)位置上TRMM像元的逐月降水量值,以及該點(diǎn)對(duì)應(yīng)的經(jīng)緯度信息(X、Y)、海拔值(H)、坡向的余弦值(A)、植被因子的值(V),構(gòu)建多元線性回歸模型F(X,Y,H,V,A)0.25°。

    (3)在0.25°分辨率下,計(jì)算模擬預(yù)測(cè)值與TRMM原始降水量數(shù)據(jù)的殘差ΔTRMMLR,并運(yùn)用樣條函數(shù)法插值到1km空間分辨率殘差影像,即

    ΔTRMMLR= TRMM3B43-F(X,Y,H,V,A)0.25o(3)

    (4)將1km分辨率的經(jīng)緯度、DEM、坡向因子和植被因子等影像代入多元線性回歸模型F(X,Y,H,V,A)0.25°中計(jì)算1km分辨率的預(yù)測(cè)降水量F(X,Y,H,V,A)。再利用式(2)通過(guò)預(yù)測(cè)降水量數(shù)據(jù)和殘差數(shù)據(jù)疊加得到1km分辨率的降尺度模擬降水量。

    1.3.3模擬結(jié)果檢驗(yàn)方法

    分別計(jì)算降尺度模型模擬的降水量數(shù)據(jù)、TRMM 3B43原始降水量數(shù)據(jù)與地面站點(diǎn)實(shí)測(cè)降水量數(shù)據(jù)的相關(guān)系數(shù)(R)、偏差(Bias)和均方根誤差(RMSE),對(duì)兩組數(shù)據(jù)的一致性和誤差進(jìn)行檢驗(yàn)分析,計(jì)算式為

    其中

    式中,iP和Mi分別表示站點(diǎn)實(shí)測(cè)降水量和降尺度模擬降水量或TRMM 3B43降水量數(shù)據(jù),n為樣本容量;相關(guān)系數(shù)R用來(lái)反映兩組數(shù)據(jù)的密切程度,R的取值范圍為[0,1],數(shù)值越接近1,表示數(shù)據(jù)的一致性越好;偏差反映站點(diǎn)實(shí)測(cè)數(shù)據(jù)與降尺度模擬數(shù)據(jù)或TRMM 3B43原始數(shù)據(jù)的偏離程度,其中Bias>0,說(shuō)明估算值大于實(shí)測(cè)值,Bias越接近0,兩組數(shù)據(jù)的吻合度越高,本研究以|Bias|表示相對(duì)誤差;均方根誤差也用來(lái)表示降尺度模擬降水量與站點(diǎn)實(shí)測(cè)降水量之間的偏差,數(shù)值越小,表示兩組數(shù)據(jù)誤差越小。

    2 結(jié)果與分析

    2.1月降水量降尺度模型的建立

    2.1.1植被指數(shù)對(duì)月降水量響應(yīng)的時(shí)間滯后性

    像元NDVI的Moran’s I指數(shù)為負(fù),主要是由河流、湖泊和人類活動(dòng)干擾等因素造成的,NDVI值的空間變異性較大,會(huì)影響降尺度模型的構(gòu)建??紤]到不同地物NDVI的大小有所差異,本研究將Moran’s I指數(shù)<-1的區(qū)域認(rèn)為是水體的干擾。計(jì)算研究區(qū)多年年最大NDVI的Moran’s I指數(shù)發(fā)現(xiàn),研究區(qū)平均Moran’s I指數(shù)為0.17,其中,65%像元的NDVI值的Moran’s I指數(shù)>0,Moran’s I指數(shù)<-1(空間變異較大)的僅占0.04%,綜合考慮Moran’s I指數(shù),故引入植被NDVI作為因子之一進(jìn)行TRMM 3B43數(shù)據(jù)的降尺度研究是可行的。同時(shí),為了排除水體NDVI值對(duì)降尺度模型精度的影響,本研究在構(gòu)建降尺度模型時(shí),未將水體的NDVI值參與方程的構(gòu)建。

    川西高原逐月NDVI與當(dāng)月、前一月和前兩月降水量的相關(guān)系數(shù)見(jiàn)表1。由表可見(jiàn),(1)從年內(nèi)的對(duì)比來(lái)看,多數(shù)月份的月NDVI值與前一月降水量有更高的相關(guān)系數(shù),尤其在4、5、6和10月,兩者相關(guān)性均通過(guò)顯著性檢驗(yàn);(2)無(wú)論是在植被生長(zhǎng)期(5-10月)還是全年平均,植被月NDVI與前一月降水量的相關(guān)系數(shù)(R前一月)>當(dāng)月(R當(dāng)月)>前兩月(R前兩月)。綜上可知,研究區(qū)植被對(duì)降水的響應(yīng)具有一定的滯后性,整體來(lái)看滯后期約為1個(gè)月,這與相關(guān)研究成果一致[20-21]。因此,本研究選取下一個(gè)月的植被NDVI值作為當(dāng)月降水的降尺度模型的植被因子。

    表 1 川西高原逐月NDVI與當(dāng)月、前一月及前兩月降水量的相關(guān)系數(shù)Table 1 Correlation coefficients between NDVI of current month and precipitation of current month, previous month, month before previous one in the Western Sichuan Plateau

    2.1.2月降水量數(shù)據(jù)的降尺度模型

    考慮到川西高原地區(qū)植被對(duì)降水響應(yīng)的滯后性,同時(shí)考慮山區(qū)存在最大或最小降水高程帶[22],將式(2)展開為一次多項(xiàng)式,并進(jìn)一步對(duì)高程和NDVI進(jìn)行三次多項(xiàng)式的展開,得到最終的降尺度模型為

    式中,Pi為第i月的降尺度模擬降水量(mm);X為經(jīng)度(°);Y為緯度(°);A為坡向的余弦值;NDVIi+1為第i+1月的NDVI值,當(dāng)i=12時(shí),NDVIi+1為下一年1月的NDVI值;H為海拔(m);ΔTRMMHR為1km空間分辨率的殘差影像;a,b,c,d,…,j,k為多元線性回歸模型的回歸系數(shù)。

    2.2月降水量降尺度模型的有效性檢驗(yàn)

    2.2.1各站點(diǎn)精度檢驗(yàn)

    將2001-2013年16個(gè)氣象觀測(cè)站的月降水量,與各站點(diǎn)對(duì)應(yīng)地理位置的降尺度模型模擬的月降水量進(jìn)行對(duì)比,其散點(diǎn)分布和相關(guān)系數(shù)見(jiàn)圖2。由圖可見(jiàn),(1)色達(dá)、甘孜、理塘、稻城和九龍站,其模擬值與實(shí)測(cè)值的相關(guān)系數(shù)R>0.97(P<0.001),且散點(diǎn)多集中在線性趨勢(shì)線周邊,說(shuō)明這些氣象站點(diǎn)的降尺度模擬結(jié)果與實(shí)測(cè)降水量的相關(guān)性顯著且偏離程度較小,數(shù)據(jù)精度高;(2)若爾蓋、馬爾康和石渠等(圖2 f-圖2 l)7個(gè)站點(diǎn),其模擬值與實(shí)測(cè)值的相關(guān)系數(shù)R在0.95~0.97(P<0.001),散點(diǎn)分布相對(duì)較集中,僅少部分點(diǎn)偏離線性趨勢(shì)線較遠(yuǎn),數(shù)據(jù)精度相對(duì)較高;(3)紅原、小金、巴塘和康定站點(diǎn)模擬值與實(shí)測(cè)值的相關(guān)系數(shù)偏低(R<0.95,P<0.001),其中小金站的相關(guān)系數(shù)最低,僅0.8781,明顯低于研究區(qū)其它站點(diǎn),且其散點(diǎn)分布較分散,有較多的點(diǎn)偏離線性趨勢(shì)線的程度較高,說(shuō)明該站點(diǎn)降尺度結(jié)果的精度相對(duì)偏低??傮w上(對(duì)比圖1a),西部海拔相對(duì)較高的甘孜藏族自治州各站的精度相對(duì)較高,海拔相對(duì)較低的阿壩藏族羌族自治州的精度相對(duì)偏低。

    以研究區(qū)氣象站點(diǎn)所在1km范圍內(nèi)的平均高程為自變量,分別以各站點(diǎn)降尺度模型模擬的降水量與實(shí)測(cè)值之間的相關(guān)系數(shù)R和相對(duì)誤差|Bias|為因變量作三次多項(xiàng)式回歸分析,高程與R和|Bias|的相關(guān)系數(shù)分別為0.7409和0.6998,均通過(guò)0.05水平的顯著性檢驗(yàn),呈現(xiàn)較強(qiáng)的三次函數(shù)關(guān)系特征。隨著海拔的升高,相關(guān)系數(shù)也隨之增加,且增加的速度不斷減緩,|Bias|值則隨著海拔的升高呈現(xiàn)急劇減小—緩慢增大的波動(dòng)變化趨勢(shì)。由圖2和圖3可見(jiàn),R值最低的站點(diǎn),其|Bias|最大,即降尺度模擬數(shù)據(jù)精度最差的站點(diǎn)是海拔高度最低的小金站,其原因是:(1)由于數(shù)據(jù)記錄的不連續(xù)性以及融合地面雨量統(tǒng)計(jì)資料時(shí)算法本身的不足,導(dǎo)致TRMM 3B43原始降水?dāng)?shù)據(jù)的精度存在一定問(wèn)題;(2)小金站點(diǎn)地處河谷區(qū),海拔高度2273m,但周圍高山脊的海拔高達(dá)4500m(圖1a),高程落差2000m左右,地形起伏較大,氣候垂直變化顯著,且?guī)X谷地貌對(duì)大氣流通的阻隔作用較大,導(dǎo)致該地區(qū)降水誤差明顯[23]。

    圖2 2001-2013年利用降尺度模型模擬各站點(diǎn)月降水量與實(shí)測(cè)值的散點(diǎn)圖Fig. 2 Scatter diagram of monthly precipitation between simulated by downscaling model and observed during 2001-2013

    圖3 川西高原地區(qū)高程與相關(guān)系數(shù)R(a)、相對(duì)偏差|Bias|(b)的散點(diǎn)圖Fig.3 Scatter diagram of elevation, R(a) and |Bias|(b) in the Western Sichuan Plateau

    2.2.2各月精度檢驗(yàn)

    分別統(tǒng)計(jì)歷年降尺度模型模擬的各月降水量和TRMM 3B43數(shù)據(jù)中各月降水量資料,與16個(gè)氣象觀測(cè)站的月降水量進(jìn)行對(duì)比,其相關(guān)系數(shù)、偏差和均方根誤差統(tǒng)計(jì)結(jié)果見(jiàn)表2。由表可見(jiàn),(1)各月降尺度模型模擬的降水量與實(shí)測(cè)站點(diǎn)降水量的相關(guān)系數(shù)均通過(guò)0.001水平的顯著性檢驗(yàn),1月和12月的相關(guān)系數(shù)相對(duì)偏低(R<0.6),其余月份兩組數(shù)據(jù)均具有較強(qiáng)的相關(guān)性(R>0.7)。其中降水量相對(duì)集中的5-10月,降尺度模型模擬的降水量與站點(diǎn)實(shí)測(cè)數(shù)據(jù)的偏差較低(Bias<0.1),同時(shí)偏差值相對(duì)較高的1-4月和11-12月降尺度模擬結(jié)果的均方根誤差較低,可見(jiàn),對(duì)川西高原地區(qū)2001-2013年的TRMM 3B43月降水?dāng)?shù)據(jù)進(jìn)行降尺度操作取得很好的效果,降尺度模擬結(jié)果比站點(diǎn)實(shí)測(cè)數(shù)據(jù)略偏高;(2)降尺度模型模擬的各月降水量、TRMM 3B43原始數(shù)據(jù)各月降水量數(shù)據(jù)2-10月的R、Bias和RMSE的值基本相近,其中4、5、6和9月降尺度模型模擬的降水量的精度比TRMM 3B43原始數(shù)據(jù)的精度略偏高;相比而言,1、11和12月降尺度模型模擬的降水量的精度較之TRMM 3B43原始數(shù)據(jù)明顯偏低,說(shuō)明多數(shù)情況下降尺度結(jié)果能夠保證原始數(shù)據(jù)的精度,但降水量偏低的情況下降尺度效果不佳,其原因一是降水量偏低的月份,降水量與植被指數(shù)的相關(guān)性偏低(表1),多元線性回歸模型的精度相對(duì)偏低,導(dǎo)致降尺度效果不佳;二是降水量偏低的月份,TRMM 3B43原始數(shù)據(jù)相對(duì)于站點(diǎn)實(shí)測(cè)值的誤差也較大,精度較低,導(dǎo)致降尺度模擬數(shù)據(jù)的精度較低。

    2.2.3典型干濕年份精度檢驗(yàn)

    選取2001-2013年間降水量最多的年份2012年(776.64mm)作為濕潤(rùn)年份,降水量最少的年份2006年(598.27mm)作為干旱年份,將典型干、濕年份的降尺度模型模擬數(shù)據(jù)和TRMM 3B43數(shù)據(jù)的各月降水量資料,與16個(gè)氣象觀測(cè)站的月降水量進(jìn)行對(duì)比。由表3可見(jiàn),(1)干濕年份兩組數(shù)據(jù)的相關(guān)系數(shù)數(shù)值接近,且均通過(guò)0.001水平的顯著性檢驗(yàn),濕潤(rùn)年份降尺度模型模擬的數(shù)據(jù)的偏差和均方根誤差均低于TRMM 3B43原始數(shù)據(jù),干旱年份降尺度模型模擬的數(shù)據(jù)的偏差低于TRMM 3B43原始數(shù)據(jù),均方根誤差略高于TRMM 3B43原始數(shù)據(jù),可見(jiàn),降尺度模型模擬的數(shù)據(jù)精度較高。(2)2012年兩組數(shù)據(jù)的均方根誤差均略低于2006年,但其相關(guān)系數(shù)大于2006年,且偏差更接近0,說(shuō)明濕潤(rùn)年份的降尺度模型模擬的數(shù)據(jù)精度略高于干旱年份。

    表2 2001-2013年兩組數(shù)據(jù)各月降水量與實(shí)測(cè)值間的誤差統(tǒng)計(jì)Table 2 Statistics of error between monthly precipitation of two sets data and observed during 2001-2013

    表3 典型年兩組數(shù)據(jù)的檢驗(yàn)結(jié)果Table 3 Validation results of two sets data in two typical years

    2.2.4整體精度檢驗(yàn)

    將2001-2013年16個(gè)氣象觀測(cè)站對(duì)應(yīng)位置、年份和月份的降尺度模型模擬數(shù)據(jù)和TRMM 3B43數(shù)據(jù)的各月降水量值分別與站點(diǎn)實(shí)測(cè)值作一元線性回歸分析,其散點(diǎn)分布和誤差統(tǒng)計(jì)見(jiàn)圖4。由圖4a可見(jiàn),(1)降尺度模型模擬的月降水量和站點(diǎn)實(shí)測(cè)月降水量的相關(guān)系數(shù)達(dá)到0.9499(P<0.001),兩者存在顯著的線性相關(guān)關(guān)系;兩組數(shù)據(jù)的偏差為0.0866,表明降尺度模型模擬的月降水量值整體略偏高。(2)對(duì)比圖4b可見(jiàn),降尺度模型模擬的月降水量與站點(diǎn)實(shí)測(cè)月降水量的R值、Bias值、RMSE值雖然比TRMM 3B43原始數(shù)據(jù)與之相對(duì)應(yīng)的值略偏低,但基本接近,表明降尺度模擬的結(jié)果與站點(diǎn)實(shí)測(cè)數(shù)據(jù)誤差形成的主要原因是TRMM 3B43數(shù)據(jù)本身存在一定的誤差[10];同時(shí),川西高原地勢(shì)復(fù)雜,地形因素對(duì)大氣流動(dòng)的影響很大,而本文在利用TRMM 3B43數(shù)據(jù)結(jié)合經(jīng)緯度、高程、坡向和NDVI進(jìn)行多元線性回歸分析時(shí),高程和坡向并不能完全代表整個(gè)地形因素對(duì)降水的影響,導(dǎo)致個(gè)別氣象觀測(cè)站點(diǎn)的模擬值誤差偏大(圖2)。

    整體而言,降尺度模型模擬的月降水量的精度較高,雖略低于TRMM 3B43原始數(shù)據(jù)的精度,但相差不大,說(shuō)明本文TRMM 3B43數(shù)據(jù)降尺度的模型是可行的,同時(shí)又提高了降水?dāng)?shù)據(jù)的空間分辨率,能夠更加全面地反映川西高原地區(qū)降水量分布的空間格局。

    圖 4 2001-2013年降尺度模擬月降水量(a)和TRMM3B43月降水量(b)與實(shí)測(cè)值間的散點(diǎn)圖Fig. 4 Scatter diagram of monthly precipitation between observed and simulated by downscaling model(a) and TRMM 3B43 data(b) during 2001-2013

    3 結(jié)論與討論

    (1)在空間上,各氣象觀測(cè)站點(diǎn)降尺度模型模擬的降水量數(shù)據(jù)精度較高,其中,精度最高的稻城站模擬結(jié)果與站點(diǎn)實(shí)測(cè)值的相關(guān)系數(shù)高達(dá)0.9839,精度最低的小金站相關(guān)系數(shù)亦高達(dá)0.8781,且該數(shù)據(jù)精度的空間差異明顯,西部的甘孜藏族自治州的數(shù)據(jù)精度比東部的阿壩藏族羌族自治州偏高,這主要與TRMM 3B43原始數(shù)據(jù)自身的精度、研究區(qū)復(fù)雜的地形地貌以及多元線性回歸模型的精度有關(guān)。

    (2)在月尺度上,降尺度模型模擬的月降水量數(shù)據(jù)與站點(diǎn)實(shí)測(cè)數(shù)據(jù)的吻合度較好,能較好地反映區(qū)域降水量的年內(nèi)變化,但各月份的精度有所差異,降水相對(duì)集中的5-10月比其它月份數(shù)據(jù)精度較高。與TRMM 3B43原始降水?dāng)?shù)據(jù)相比,降尺度模型模擬的月降水量在不降低其數(shù)據(jù)精度的基礎(chǔ)上,大幅提高了數(shù)據(jù)的空間分辨率,使其能更加全面、精細(xì)地反映研究區(qū)降水的空間分布特征。

    (3)在年尺度上,典型干、濕年份降尺度模型模擬的降水量數(shù)據(jù)的精度較高,且略高于TRMM 3B43原始數(shù)據(jù)的精度。同時(shí),濕潤(rùn)年份的精度高于干旱年份,這與在降尺度過(guò)程中同樣考慮到植被NDVI因子的Immerzeel等[13,16]的研究成果不一致,其原因一是植被的影響,濕潤(rùn)年份與干旱年份相比,降水量更充足,植被對(duì)降水響應(yīng)的滯后性更明顯,降水集中的5-10月多元回歸模型的效果更好,精度更高;二是TRMM 3B43原始數(shù)據(jù)的影響,濕潤(rùn)年份TRMM原始數(shù)據(jù)的精度比干旱年份偏高,與本文降尺度模型模擬結(jié)果干濕年份的精度對(duì)比結(jié)果一致,這也進(jìn)一步說(shuō)明了在降尺度過(guò)程中考慮植被對(duì)降水響應(yīng)的滯后性,降尺度模型模擬的數(shù)據(jù)與TRMM原始數(shù)據(jù)能呈現(xiàn)更好的一致性。

    (4)研究區(qū)2001-2013年156個(gè)月的降尺度模型模擬的月降水量與16個(gè)地面觀測(cè)實(shí)測(cè)值整體上的R=0.9499,Bias=0.0866,RMSE=19.92,說(shuō)明降尺度結(jié)果與站點(diǎn)實(shí)測(cè)降水?dāng)?shù)據(jù)之間存在顯著的線性相關(guān)關(guān)系,數(shù)據(jù)誤差較小,數(shù)據(jù)精度較高,降尺度模型可行。且與TRMM 3B43原始數(shù)據(jù)的整體精度相比,降尺度模型模擬的數(shù)據(jù)精度與之非常接近,只是略偏低,說(shuō)明本文降尺度數(shù)據(jù)適于川西高原地區(qū)的降水研究。但馬金輝等[10,24]在對(duì)TRMM 3B43年降水量降尺度研究中發(fā)現(xiàn)降尺度結(jié)果的精度較之原始數(shù)據(jù)偏高,可能是因?yàn)轳R金輝等是針對(duì)TRMM年降水?dāng)?shù)據(jù)進(jìn)行降尺度研究,而本文是對(duì)TRMM月降水?dāng)?shù)據(jù)進(jìn)行降尺度研究,研究的時(shí)間尺度不同,且年內(nèi)不同月份降水量的降尺度精度有所差異,說(shuō)明降水量的年內(nèi)變化是影響年尺度和月尺度降水?dāng)?shù)據(jù)降尺度精度的主要原因。

    參考文獻(xiàn)References

    [1] Goovaerts P.Geostatistical approaches for incorporating elevation into the spatial interpolation of rainfall[J].Journal of Hydrology, 2000,228(1):113-129.

    [2] Bohnenstengel S,Schlünzen K,Beyrich F.Representativity of in situ precipitation measurements :a case study for the LITFASS area in North-Eastern Germany[J].Journal of Hydrology,2011,400(3-4):387-395.

    [3]曾麗紅,宋開山,張柏,等.1960-2008年吉林省降水量的時(shí)空演變特征[J].中國(guó)農(nóng)業(yè)氣象,2010,31(3):344-352.

    Zeng L H,Song K S,Zhang B,et al.Spatial-temporal precipitation variation over Jilin Province during 1960-2008[J].Chinese Journal of Agrometeorology,2010,31(3):344-352.(in Chinese)

    [4]朱會(huì)義,劉述林,賈紹鳳.自然地理要素空間插值的幾個(gè)問(wèn)題[J].地理研究,2004,23(4):425-432.

    Zhu H Y,Liu S L,Jia S F.Problems of the spatial interpolation of physical geographical elements[J].Geographica Research, 2004,23(4):425-432.(in Chinese)

    [5]張濤,李寶林,何元慶,等.基于TRMM訂正數(shù)據(jù)的橫斷山區(qū)降水時(shí)空分布特征[J].自然資源學(xué)報(bào),2015,30(2):261-270.

    Zhang T,Li B L,He Y Q,et al.Spatial and temporal distribution of precipitation based on corrected TRMM data in Hengduan Mountains[J].Journal of Natural Resources,2015,30(2):261-270. (in Chinese)

    [6]呂洋,楊勝天,蔡明勇,等.TRMM衛(wèi)星降水?dāng)?shù)據(jù)在雅魯藏布江流域的適用性分析[J].自然資源學(xué)報(bào),2013,28(8):1414-1425.

    Lv Y,Yang S T,Cai M Y,et al.The applicability analysis of TRMM precipitation data in the Yarlung Zangbo River Basin[J]. Journal of Natural Resources,2013,28(8):1414-1425.(in Chinese)

    [7]曾紅偉,李麗娟.瀾滄江及周邊流域TRMM3B43數(shù)據(jù)精度檢驗(yàn)[J].地理學(xué)報(bào),2011,66(7):994-1004.

    Zeng H W,Li L J.Accuracy validation of TRMM3B43data in Lancang River Basin[J].Acta Geographica Sinica,2011,66(7): 994-1004.(in Chinese)

    [8]朱國(guó)鋒,蒲燾,張濤,等.TRMM降水?dāng)?shù)據(jù)在橫斷山區(qū)的精度[J].地理科學(xué),2013,33(9):1125-1131.

    Zhu G F,Pu T,Zhang T,et al.The accuracy of TRMM precipitation data in Hengduan Mountainous region,China[J]. Scientia Geographica Sinica,2013,33(9):1125-1131.(in Chinese)

    [9]周鎖銓,薛根元,周麗峰,等.基于GIS降水空間分析的逐步插值方法[J].氣象學(xué)報(bào),2006,64(1):100-111.

    Zhou S Q,Xue G Y,Zhou L F,et al.The stepwise interpolation approach of precipitation for spatial analysis based on GIS[J]. Acta Meteopologica Sinica,2006,64(1):100-111.(in Chinese)

    [10]馬金輝,屈創(chuàng),張海筱,等.2001-2010年石羊河流域上游TRMM降水資料的降尺度研究[J].地理科學(xué)進(jìn)展,2013, 32(9):1423-1432.

    Ma J H,Qu C,Zhang H X,et al.Spatial downscaling of TRMM precipitation data based on DEM in the upstream of Shiyang River Basin during 2001-2010[J].Progress in Geography,2013,32(9):1423-1432.(in Chinese)

    [11]徐成東,孔云峰,仝文偉.線性加權(quán)回歸模型的高原山地區(qū)域降水空間插值研究[J].地球信息科學(xué),2008,10(1):14-19.

    Xu C D,Kong Y F,Tong W W.A weighted Linear Regression Model for precipitation spatial interpolation in Altiplano and Mountain in area[J].Geo-information Science,2008,10(1):14-19. (in Chinese)

    [12]李曉光,劉華民,王立新,等.鄂爾多斯高原植被覆蓋變化及其與氣候和人類活動(dòng)的關(guān)系[J].中國(guó)農(nóng)業(yè)氣象,2014,35(4): 470-476.

    Li X G,Liu H M,Wang L X,et al.Vegetation cover change and its relationship between climate and human activities in Ordos Plateau[J].Chinese Journal of Agrometeorology,2014, 35(4):470-476.(in Chinese)

    [13]Immerzeel W W,Rutten M M,Droogers P.Spatial downscaling of TRMM precipitation using vegetative response on the Iberian Peninsula[J].Remote Sensing of Environment,2009, 113(2):362-370.

    [14]王智,吳友均,梁鳳超,等.新疆地區(qū)年降水量的空間插值方法研究[J].中國(guó)農(nóng)業(yè)氣象,2011,32(3):331-337.

    Wang Z,Wu Y J,Liang F C,et al.Study on spatial interpolation method of annual precipitation in Xinjiang[J].Chinese Journal of Agrometeorology,2011,32(3):331-337.(in Chinese)

    [15]Jia S,Zhu W,Lv A,et al.A statistical spatial downscaling algorithm of TRMM precipitation based on NDVI and DEM in the Qaidam Basin of China[J].Remote Sensing of Environment,2011,115(12):3069-3079.

    [16]嵇濤,劉睿,楊華,等.多源遙感數(shù)據(jù)的降水空間降尺度研究:以川渝地區(qū)為例[J].地球信息科學(xué)學(xué)報(bào),2015,17(1):108-117.

    Ji T,Liu R,Yang H,et al.Spatial downscaling of precipitation using multi-source remote sensing data:a case study of Sichuan-Chongqing Region[J].Journal of Geo-information Science,2015,17(1):108-117.(in Chinese)

    [17]張虹嬌.川西北高原氣候變化特征研究[J].西南大學(xué)學(xué)報(bào)(自然科學(xué)版),2014,36(12):148-156.

    Zhang H J.A study on the characteristics of climate change on Northwestern Sichuan Plateau[J].Journal of Southwest University(Natural Science Edition),2014,36(12):148-156.(inChinese)

    [18]姜琳,馮文蘭,郭兵.雅魯藏布江流域近13年植被覆蓋動(dòng)態(tài)監(jiān)測(cè)及與降水因子的相關(guān)性分析[J].長(zhǎng)江流域資源與環(huán)境,2014,23(11):1610-1619.

    Jiang L,Feng W L,Guo B.Analysis of dynamic monitoring of vegetation change and the correlation with precipitation factor in Yalu Tsangpo River Basin during the past 13 years[J].Resources and Environment in the Yangtze Basin, 2014,23(11):1610-1619.(in Chinese)

    [19] Chen J,Brissette F P,Leconte R.Coupling statistical and dynamical methods for spatial downscaling of precipitation[J]. Climatic Change,2012,114(3-4):509-526.

    [20]白建軍,白江濤,王磊.2000-2010年陜北地區(qū)植被NDVI時(shí)空變化及其與區(qū)域氣候的關(guān)系[J].地理科學(xué),2014,34(7): 882-888.

    Bai J J,Bai J T,Wang L.Spatio-temporal change of vegetation NDVI and its relation with regional climate in Northern Shaanxi Province in 2000-2010[J].Scientia Geographica Sinica,2014,34(7):882-888.(in Chinese)

    [21]白淑英,王莉,史建橋.長(zhǎng)江流域NDVI對(duì)氣候變化響應(yīng)的時(shí)滯效應(yīng)[J].中國(guó)農(nóng)業(yè)氣象,2012,33(4):579-586.

    Bai S Y,Wang L,Shi J Q.Time lag effect of NDVI response to climatic change in Yangtze River Basin[J].Chinese Journal of Agrometeorology,2012,33(4):579-586.(in Chinese)

    [22]王超,趙傳燕.TRMM多衛(wèi)星資料在黑河上游降水時(shí)空特征研究中的應(yīng)用[J].自然資源學(xué)報(bào),2013,28(5):862-872.

    Wang C,Zhao C Y.A study of the spatio-temporal distribution of precipitation in upper reaches of Heihe River of China using TRMM data[J].Journal of Natural Resources, 2013,28(5):862-872.(in Chinese)

    [23]嵇濤,楊華,劉睿,等.TRMM衛(wèi)星降水?dāng)?shù)據(jù)在川渝地區(qū)的適用性分析[J].地理科學(xué)進(jìn)展,2014,33(10):1375-1386.

    Ji T,Yang H,Liu R,et al.Applicability analysis of the TRMM precipitation data in the Sichuan-Chongqing region[J].Progress in Geography,2014,33(10):1375-1386.(in Chinese)

    [24]李凈,張曉.TRMM降水?dāng)?shù)據(jù)的空間降尺度方法研究[J].地理科學(xué),2015,35(9):1164-1169.

    Li J,Zhang X.Downscaling method of TRMM satellite precipitation data[J].Scientia Geographica Sinica,2015,35(9): 1164-1169.(in Chinese)

    Spatial Downscaling Simulation of Monthly Precipitation Based on TRMM 3B43 Data in the Western Sichuan Plateau

    ZHENG Jie, LV Li, FENG Wen-lan, TU Kun
    (College of Resources and Environment, Chengdu University of Information Technology, Chengdu 610225, China)

    Abstract:Using TRMM 3B43 data, MODIS-NDVI, DEM, meteorological data of observation stations during 2001-2013, a multiple linear regression model was built among TRMM 3B43 monthly precipitation and longitude, latitude, altitude, aspect, NDVI factor as downscaling model of monthly precipitation data in the Western Sichuan Plateau based on the analysis of the lag of vegetation response to precipitation. Then, combined with regression equation and residuals, interpolation method was adopted to obtain monthly precipitation data with 1km spatial resolution. Finally, the accuracy of simulated data obtained by downscaling model was tested by the correlation analysis and error detection between the simulated results and the observation data of 16 meteorological stations in the study area. The results showed as follows: (1) Precipitation simulated by downscaling model based on TRMM 3B43 data had a high precision in all meteorological observation stations. Daocheng site displayed the highest accuracy, the correlation coefficient between simulation results and observed values attained 0.9839, while Xiaojin site displayed the lowest accuracy, the correlation coefficient was 0.8781. (2) The precipitation simulated by downscaling model displayed high accuracy in the whole study area at both monthly and yearly time scale. The accuracy of simulated results from May to October was significantly higher than other months, as well as typical wet year (2012) was higher than dry year (2006). (3) The precipitation simulated by downscaling model displayed high accuracy to the observation data on the whole(R=0.9499, Bias=0.0866), while thevalue of simulated precipitation was slightly higher. (4) Compared to the original data TRMM 3B43, the simulated data by downscaling model guaranteed the accuracy and improved the spatial resolution. So, this method could provide an effective way to produce a more sophisticated precipitation data with higher spatial resolution.

    Key words:TRMM 3B43 data; Precipitation data; Spatial downscaling; Lag; Multiple linear regression

    doi:10.3969/j.issn.1000-6362.2016.02.015

    * 收稿日期:2015-08-11**通訊作者。E-mail:fwl@cuit.edu.cn

    基金項(xiàng)目:國(guó)家自然科學(xué)基金項(xiàng)目(41301653)

    作者簡(jiǎn)介:鄭杰(1991-),碩士生,主要從事3S集成與氣象應(yīng)用研究。E-mail:zhengjie0601@sina.com

    猜你喜歡
    滯后性多元線性回歸
    基于教材 深入課堂 勇于實(shí)踐
    東方教育(2016年22期)2017-04-07 16:58:42
    電信立法若干問(wèn)題研究
    卷宗(2016年10期)2017-01-21 18:31:14
    電信立法若干問(wèn)題研究
    卷宗(2016年10期)2017-01-21 18:30:43
    淺析醫(yī)院成本核算工作中出現(xiàn)的問(wèn)題
    基于組合模型的卷煙市場(chǎng)需求預(yù)測(cè)研究
    基于多元線性回歸分析的冬季鳥類生境選擇研究
    我國(guó)上市商業(yè)銀行信貸資產(chǎn)證券化效應(yīng)實(shí)證研究
    云學(xué)習(xí)平臺(tái)大學(xué)生學(xué)業(yè)成績(jī)預(yù)測(cè)與干預(yù)研究
    全國(guó)主要市轄區(qū)的房?jī)r(jià)收入比影響因素研究
    商(2016年20期)2016-07-04 01:23:26
    利用計(jì)量工具比較東西部的經(jīng)濟(jì)狀況
    商(2016年5期)2016-03-28 12:14:30
    国产三级黄色录像| 又大又爽又粗| 50天的宝宝边吃奶边哭怎么回事| 熟女少妇亚洲综合色aaa.| 成人18禁高潮啪啪吃奶动态图| 人人澡人人妻人| 日韩免费高清中文字幕av| 亚洲成人免费电影在线观看| 在线亚洲精品国产二区图片欧美| 一级黄色大片毛片| 亚洲精品国产av蜜桃| 99久久人妻综合| 国产又爽黄色视频| 青春草视频在线免费观看| 欧美精品高潮呻吟av久久| h视频一区二区三区| 久久精品国产亚洲av香蕉五月 | 久久久精品94久久精品| 国产亚洲精品一区二区www | 别揉我奶头~嗯~啊~动态视频 | 自线自在国产av| 在线观看人妻少妇| 宅男免费午夜| 国产免费现黄频在线看| 国产福利在线免费观看视频| 在线观看免费午夜福利视频| 淫妇啪啪啪对白视频 | 成人亚洲精品一区在线观看| 热99re8久久精品国产| 黑人欧美特级aaaaaa片| 亚洲精品国产区一区二| 成年人午夜在线观看视频| 亚洲精品第二区| 在线观看免费视频网站a站| 成年美女黄网站色视频大全免费| av网站在线播放免费| 夜夜骑夜夜射夜夜干| 伊人久久大香线蕉亚洲五| 国产欧美日韩精品亚洲av| av超薄肉色丝袜交足视频| 777米奇影视久久| 亚洲专区字幕在线| 91成人精品电影| 国产亚洲av片在线观看秒播厂| 国产在线视频一区二区| 后天国语完整版免费观看| 免费在线观看黄色视频的| 女人精品久久久久毛片| 国产成人av激情在线播放| 黄色毛片三级朝国网站| 午夜精品国产一区二区电影| 性高湖久久久久久久久免费观看| 国产成人av教育| 黄色a级毛片大全视频| 亚洲专区国产一区二区| 欧美午夜高清在线| 国产黄色免费在线视频| 国产精品一区二区在线不卡| 国产精品国产三级国产专区5o| 两个人免费观看高清视频| 99国产精品一区二区蜜桃av | 欧美成狂野欧美在线观看| 亚洲精品国产色婷婷电影| 一二三四在线观看免费中文在| 国产真人三级小视频在线观看| 欧美激情高清一区二区三区| av网站在线播放免费| 色婷婷久久久亚洲欧美| 午夜福利乱码中文字幕| a级毛片在线看网站| 老司机福利观看| 免费日韩欧美在线观看| 久久狼人影院| 欧美精品一区二区大全| 精品一区在线观看国产| 美女视频免费永久观看网站| av有码第一页| e午夜精品久久久久久久| 国产成人精品久久二区二区免费| 国产成人系列免费观看| 交换朋友夫妻互换小说| 伦理电影免费视频| 少妇粗大呻吟视频| 亚洲av电影在线进入| 精品久久久久久电影网| 一本大道久久a久久精品| 精品人妻一区二区三区麻豆| 肉色欧美久久久久久久蜜桃| 免费久久久久久久精品成人欧美视频| 亚洲中文日韩欧美视频| av天堂在线播放| 天天添夜夜摸| av网站在线播放免费| 在线十欧美十亚洲十日本专区| tube8黄色片| 国产伦人伦偷精品视频| 国产精品 欧美亚洲| 欧美精品亚洲一区二区| 欧美老熟妇乱子伦牲交| 69精品国产乱码久久久| 9191精品国产免费久久| 日韩欧美国产一区二区入口| 性色av乱码一区二区三区2| 成年女人毛片免费观看观看9 | 中文字幕色久视频| 老司机影院成人| 日本欧美视频一区| 中国美女看黄片| 国产高清国产精品国产三级| 黄色 视频免费看| 亚洲精品国产精品久久久不卡| 人人妻人人澡人人看| 男男h啪啪无遮挡| 日韩制服丝袜自拍偷拍| 99国产精品一区二区蜜桃av | 国产成人一区二区三区免费视频网站| netflix在线观看网站| 无遮挡黄片免费观看| 国产亚洲精品久久久久5区| 777米奇影视久久| 亚洲专区国产一区二区| 日韩,欧美,国产一区二区三区| 波多野结衣一区麻豆| 天堂8中文在线网| 99热全是精品| 亚洲精品国产精品久久久不卡| 日韩一卡2卡3卡4卡2021年| 国产亚洲精品一区二区www | 日本91视频免费播放| 日本a在线网址| 亚洲国产欧美在线一区| 亚洲专区国产一区二区| 国产免费现黄频在线看| 2018国产大陆天天弄谢| 女人被躁到高潮嗷嗷叫费观| 性高湖久久久久久久久免费观看| 老司机福利观看| 麻豆乱淫一区二区| 欧美日本中文国产一区发布| 99久久精品国产亚洲精品| 亚洲自偷自拍图片 自拍| 国产精品久久久av美女十八| 一区二区三区精品91| 亚洲 欧美一区二区三区| 久久人人97超碰香蕉20202| 人人妻人人澡人人爽人人夜夜| 爱豆传媒免费全集在线观看| 亚洲熟女毛片儿| 亚洲成人国产一区在线观看| 69精品国产乱码久久久| 免费日韩欧美在线观看| 欧美日本中文国产一区发布| av又黄又爽大尺度在线免费看| 亚洲午夜精品一区,二区,三区| 国产欧美日韩一区二区三 | 成年人午夜在线观看视频| 在线永久观看黄色视频| 中文欧美无线码| 久久久久精品国产欧美久久久 | 亚洲视频免费观看视频| 亚洲精品美女久久av网站| 久久精品国产亚洲av香蕉五月 | 日韩一区二区三区影片| 黄色怎么调成土黄色| 国产精品免费视频内射| 极品少妇高潮喷水抽搐| 97在线人人人人妻| 国产欧美亚洲国产| 日本一区二区免费在线视频| 在线观看免费高清a一片| 一区二区三区激情视频| 久久久欧美国产精品| 久久人妻福利社区极品人妻图片| 爱豆传媒免费全集在线观看| 午夜日韩欧美国产| 999精品在线视频| 91麻豆精品激情在线观看国产 | 亚洲综合色网址| 欧美乱码精品一区二区三区| 亚洲精品国产av蜜桃| 99精品欧美一区二区三区四区| 啪啪无遮挡十八禁网站| 一本—道久久a久久精品蜜桃钙片| 丁香六月欧美| 999久久久精品免费观看国产| 欧美黑人精品巨大| 成年美女黄网站色视频大全免费| 亚洲伊人久久精品综合| 黄色毛片三级朝国网站| 国产欧美日韩一区二区三区在线| 天天躁狠狠躁夜夜躁狠狠躁| 亚洲男人天堂网一区| 中国美女看黄片| 国产成人啪精品午夜网站| 中文字幕精品免费在线观看视频| 丝瓜视频免费看黄片| 亚洲中文av在线| 成人黄色视频免费在线看| 精品一区在线观看国产| 欧美精品亚洲一区二区| avwww免费| 高清av免费在线| 夜夜骑夜夜射夜夜干| 天天躁狠狠躁夜夜躁狠狠躁| 国产精品av久久久久免费| 欧美另类一区| 高清欧美精品videossex| tube8黄色片| 国产有黄有色有爽视频| 在线av久久热| 又黄又粗又硬又大视频| 精品亚洲乱码少妇综合久久| 成人亚洲精品一区在线观看| 国产成人av教育| 亚洲国产欧美网| 亚洲少妇的诱惑av| 99精国产麻豆久久婷婷| 色婷婷av一区二区三区视频| 精品久久久精品久久久| 久久人人97超碰香蕉20202| 性高湖久久久久久久久免费观看| 日本撒尿小便嘘嘘汇集6| 女性被躁到高潮视频| 国产欧美日韩一区二区三 | 日本猛色少妇xxxxx猛交久久| 亚洲精品国产av成人精品| 久久精品熟女亚洲av麻豆精品| 免费观看人在逋| 99热国产这里只有精品6| 久久久久久久精品精品| 性色av一级| 亚洲久久久国产精品| tocl精华| 亚洲中文日韩欧美视频| 97在线人人人人妻| 91av网站免费观看| 亚洲国产欧美在线一区| 曰老女人黄片| 热99久久久久精品小说推荐| 超色免费av| 精品少妇久久久久久888优播| 欧美亚洲日本最大视频资源| 啦啦啦 在线观看视频| 久久香蕉激情| 亚洲精品久久午夜乱码| 999久久久国产精品视频| 天天添夜夜摸| 亚洲国产中文字幕在线视频| 亚洲专区中文字幕在线| 欧美日韩视频精品一区| 亚洲欧美精品综合一区二区三区| 国产有黄有色有爽视频| 亚洲av片天天在线观看| 美女高潮喷水抽搐中文字幕| 色综合欧美亚洲国产小说| 久久久久久久大尺度免费视频| 天天操日日干夜夜撸| 99国产精品99久久久久| 精品少妇黑人巨大在线播放| 国产成人啪精品午夜网站| 丝瓜视频免费看黄片| 又大又爽又粗| 午夜激情av网站| 成人亚洲精品一区在线观看| 久久国产精品人妻蜜桃| 午夜免费观看性视频| 热99re8久久精品国产| 亚洲欧美激情在线| 肉色欧美久久久久久久蜜桃| xxxhd国产人妻xxx| 久久亚洲精品不卡| 91麻豆av在线| 亚洲国产毛片av蜜桃av| 三级毛片av免费| 麻豆乱淫一区二区| 日韩人妻精品一区2区三区| 精品少妇内射三级| 91麻豆精品激情在线观看国产 | 午夜福利,免费看| 婷婷丁香在线五月| 99国产极品粉嫩在线观看| 超色免费av| 老司机影院成人| 女人久久www免费人成看片| 欧美精品av麻豆av| 久久狼人影院| 亚洲第一欧美日韩一区二区三区 | 精品一区二区三卡| 菩萨蛮人人尽说江南好唐韦庄| 亚洲精品第二区| 成年美女黄网站色视频大全免费| 国产1区2区3区精品| 久久久欧美国产精品| 首页视频小说图片口味搜索| 每晚都被弄得嗷嗷叫到高潮| 一区在线观看完整版| 亚洲专区国产一区二区| 午夜福利在线观看吧| 91大片在线观看| 欧美亚洲日本最大视频资源| 国产成人影院久久av| 色婷婷久久久亚洲欧美| 美女脱内裤让男人舔精品视频| 自拍欧美九色日韩亚洲蝌蚪91| 国产精品久久久久久精品电影小说| 色94色欧美一区二区| 亚洲精品国产一区二区精华液| 97精品久久久久久久久久精品| 久久人妻熟女aⅴ| 精品卡一卡二卡四卡免费| 91精品三级在线观看| 欧美少妇被猛烈插入视频| 国产成人系列免费观看| 十分钟在线观看高清视频www| 亚洲精品国产av蜜桃| 嫁个100分男人电影在线观看| 国产精品熟女久久久久浪| 久久女婷五月综合色啪小说| 一区二区三区四区激情视频| 18禁国产床啪视频网站| 三级毛片av免费| 国产av国产精品国产| 国产精品秋霞免费鲁丝片| 免费一级毛片在线播放高清视频 | 欧美精品一区二区免费开放| svipshipincom国产片| 亚洲熟女毛片儿| 女人久久www免费人成看片| 少妇的丰满在线观看| 欧美精品亚洲一区二区| 99re6热这里在线精品视频| 亚洲成国产人片在线观看| 手机成人av网站| 日本精品一区二区三区蜜桃| 欧美日韩一级在线毛片| 在线观看免费日韩欧美大片| 精品国产乱码久久久久久小说| 超碰成人久久| 午夜福利一区二区在线看| 在线观看免费日韩欧美大片| 亚洲欧美一区二区三区久久| 日韩欧美一区视频在线观看| av片东京热男人的天堂| a 毛片基地| 亚洲国产精品999| 亚洲天堂av无毛| 亚洲自偷自拍图片 自拍| 亚洲男人天堂网一区| 777米奇影视久久| 飞空精品影院首页| 精品福利永久在线观看| 搡老乐熟女国产| 丁香六月天网| 中文欧美无线码| 女性被躁到高潮视频| 精品人妻熟女毛片av久久网站| 国产97色在线日韩免费| 亚洲精品在线美女| 成年美女黄网站色视频大全免费| 天堂8中文在线网| 日韩视频在线欧美| 免费观看a级毛片全部| 欧美日韩精品网址| 亚洲av日韩精品久久久久久密| svipshipincom国产片| 国产深夜福利视频在线观看| 中国美女看黄片| 人人澡人人妻人| 老熟妇仑乱视频hdxx| 国产亚洲精品久久久久5区| 欧美一级毛片孕妇| 亚洲欧美日韩另类电影网站| 日本精品一区二区三区蜜桃| 亚洲第一欧美日韩一区二区三区 | 亚洲第一欧美日韩一区二区三区 | 韩国高清视频一区二区三区| 久久国产精品大桥未久av| 国产精品一区二区在线不卡| 好男人电影高清在线观看| 老熟女久久久| 亚洲国产欧美一区二区综合| 欧美另类亚洲清纯唯美| 老司机靠b影院| 国产精品av久久久久免费| 91大片在线观看| 无遮挡黄片免费观看| 香蕉国产在线看| 午夜激情av网站| 免费高清在线观看日韩| 国产精品欧美亚洲77777| 天天躁夜夜躁狠狠躁躁| 美女大奶头黄色视频| 日韩免费高清中文字幕av| 高清视频免费观看一区二区| 大码成人一级视频| 老司机靠b影院| 中文字幕人妻丝袜制服| 精品国产国语对白av| 欧美激情 高清一区二区三区| 国产精品亚洲av一区麻豆| 午夜两性在线视频| 色94色欧美一区二区| www.熟女人妻精品国产| 美女中出高潮动态图| 免费看十八禁软件| 欧美日韩国产mv在线观看视频| 国产av一区二区精品久久| 亚洲国产av影院在线观看| 97人妻天天添夜夜摸| 免费高清在线观看日韩| 亚洲精品av麻豆狂野| 另类亚洲欧美激情| 最新在线观看一区二区三区| 又紧又爽又黄一区二区| 国产在线免费精品| www.熟女人妻精品国产| 国产在线一区二区三区精| 男女边摸边吃奶| 亚洲av片天天在线观看| 精品乱码久久久久久99久播| 一本—道久久a久久精品蜜桃钙片| 精品国产超薄肉色丝袜足j| 男女之事视频高清在线观看| 中文字幕av电影在线播放| 日韩欧美一区视频在线观看| 人人妻人人澡人人看| 国产精品一二三区在线看| 免费观看av网站的网址| av在线播放精品| 久久中文看片网| av网站在线播放免费| 久久人人97超碰香蕉20202| 久久久久国内视频| 不卡av一区二区三区| videosex国产| 一级毛片女人18水好多| 久久久国产精品麻豆| 亚洲精品国产av成人精品| 肉色欧美久久久久久久蜜桃| 欧美日韩精品网址| 中文欧美无线码| 国产高清视频在线播放一区 | 日本五十路高清| 成人国语在线视频| 久久狼人影院| 久久久国产一区二区| 亚洲伊人色综图| 亚洲av男天堂| 午夜福利一区二区在线看| 丰满迷人的少妇在线观看| 中文字幕人妻丝袜一区二区| 黄色毛片三级朝国网站| 免费观看av网站的网址| 国产视频一区二区在线看| 欧美日韩视频精品一区| 国产亚洲精品第一综合不卡| 久久国产精品人妻蜜桃| 亚洲第一青青草原| 真人做人爱边吃奶动态| 狠狠婷婷综合久久久久久88av| 一边摸一边抽搐一进一出视频| 视频在线观看一区二区三区| 在线亚洲精品国产二区图片欧美| 十分钟在线观看高清视频www| 亚洲国产看品久久| 青青草视频在线视频观看| 美女高潮喷水抽搐中文字幕| 久久精品国产综合久久久| 在线看a的网站| 日韩 亚洲 欧美在线| 成人黄色视频免费在线看| 嫩草影视91久久| 老司机福利观看| 91麻豆精品激情在线观看国产 | 久久精品人人爽人人爽视色| 欧美精品av麻豆av| 免费少妇av软件| 久久久国产成人免费| 999精品在线视频| 久久九九热精品免费| 日韩一区二区三区影片| 男女边摸边吃奶| 亚洲欧美精品综合一区二区三区| 黄频高清免费视频| 真人做人爱边吃奶动态| 亚洲国产日韩一区二区| 中国国产av一级| 亚洲精品久久成人aⅴ小说| 成年女人毛片免费观看观看9 | 麻豆av在线久日| 亚洲国产精品成人久久小说| 制服诱惑二区| 黑人巨大精品欧美一区二区mp4| 性色av乱码一区二区三区2| 欧美乱码精品一区二区三区| 十八禁网站免费在线| 成人影院久久| 男女无遮挡免费网站观看| 18禁裸乳无遮挡动漫免费视频| 99久久人妻综合| 一区二区av电影网| 午夜免费成人在线视频| 午夜成年电影在线免费观看| 国产三级黄色录像| 岛国在线观看网站| 国产av又大| 国产在线一区二区三区精| 午夜成年电影在线免费观看| 黄色视频不卡| 女性被躁到高潮视频| 久久精品亚洲av国产电影网| 亚洲成人免费av在线播放| 男女国产视频网站| 免费高清在线观看日韩| 国产精品1区2区在线观看. | 国产亚洲欧美在线一区二区| 国产老妇伦熟女老妇高清| 黄片大片在线免费观看| 在线观看人妻少妇| 成人影院久久| 精品国产一区二区三区久久久樱花| 亚洲国产精品成人久久小说| 免费在线观看黄色视频的| 国产xxxxx性猛交| 黄色怎么调成土黄色| 精品国产乱码久久久久久男人| 午夜两性在线视频| 日本撒尿小便嘘嘘汇集6| 欧美变态另类bdsm刘玥| 免费一级毛片在线播放高清视频 | 亚洲综合色网址| 日本精品一区二区三区蜜桃| 无限看片的www在线观看| 亚洲精品日韩在线中文字幕| 十八禁高潮呻吟视频| 免费av中文字幕在线| 成人亚洲精品一区在线观看| 九色亚洲精品在线播放| 亚洲国产av影院在线观看| 午夜免费观看性视频| 中文字幕高清在线视频| a级毛片在线看网站| 亚洲专区中文字幕在线| 欧美一级毛片孕妇| 国产欧美亚洲国产| 下体分泌物呈黄色| 人人妻人人澡人人看| 一区福利在线观看| 免费一级毛片在线播放高清视频 | 欧美中文综合在线视频| 亚洲精品粉嫩美女一区| 97人妻天天添夜夜摸| 国产欧美日韩综合在线一区二区| 80岁老熟妇乱子伦牲交| 日日摸夜夜添夜夜添小说| 久久久久久久久久久久大奶| 久久久久久久国产电影| 欧美另类一区| 免费高清在线观看日韩| 丰满饥渴人妻一区二区三| 国产精品亚洲av一区麻豆| www.999成人在线观看| 一本大道久久a久久精品| 波多野结衣一区麻豆| 国产xxxxx性猛交| 国产1区2区3区精品| 久久久久久久久久久久大奶| 精品亚洲乱码少妇综合久久| 中国国产av一级| 天堂中文最新版在线下载| 蜜桃国产av成人99| www.精华液| 人人妻人人澡人人爽人人夜夜| 午夜两性在线视频| 亚洲av男天堂| 一二三四社区在线视频社区8| 欧美国产精品va在线观看不卡| 精品国产乱码久久久久久小说| 日韩欧美一区二区三区在线观看 | 国产国语露脸激情在线看| 精品视频人人做人人爽| 91精品国产国语对白视频| 丝袜人妻中文字幕| 91成人精品电影| 国产99久久九九免费精品| 亚洲av欧美aⅴ国产| 97在线人人人人妻| 国产精品一区二区免费欧美 | 性色av一级| 亚洲成国产人片在线观看| 一二三四社区在线视频社区8| 欧美日韩av久久| 亚洲国产av新网站| 精品乱码久久久久久99久播| 久久精品亚洲熟妇少妇任你| 视频区图区小说| 亚洲视频免费观看视频| 不卡av一区二区三区| 亚洲人成电影免费在线| 蜜桃在线观看..| 我要看黄色一级片免费的| 免费在线观看黄色视频的| 热99久久久久精品小说推荐| 高清欧美精品videossex| 99香蕉大伊视频| 亚洲视频免费观看视频| 欧美人与性动交α欧美软件| 女警被强在线播放| 宅男免费午夜| 欧美性长视频在线观看| 一区二区av电影网| 精品乱码久久久久久99久播| 亚洲国产精品成人久久小说| 日韩制服骚丝袜av| 日本精品一区二区三区蜜桃|