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

    基于VTCI和分位數(shù)回歸模型的冬小麥單產(chǎn)估測方法

    2017-07-31 20:54:09王鵬新張樹譽(yù)
    關(guān)鍵詞:估產(chǎn)關(guān)中平原位數(shù)

    王 蕾 王鵬新 李 俐 張樹譽(yù)

    (1.中國農(nóng)業(yè)大學(xué)信息與電氣工程學(xué)院,北京100083;2.農(nóng)業(yè)部農(nóng)業(yè)災(zāi)害遙感重點(diǎn)實(shí)驗(yàn)室,北京100083; 3.陜西省氣象局,西安710014)

    基于VTCI和分位數(shù)回歸模型的冬小麥單產(chǎn)估測方法

    王 蕾1,2王鵬新1,2李 俐1,2張樹譽(yù)3

    (1.中國農(nóng)業(yè)大學(xué)信息與電氣工程學(xué)院,北京100083;2.農(nóng)業(yè)部農(nóng)業(yè)災(zāi)害遙感重點(diǎn)實(shí)驗(yàn)室,北京100083; 3.陜西省氣象局,西安710014)

    條件植被溫度指數(shù)(VTCI)是一種綜合了歸一化植被指數(shù)(NDVI)與地表溫度(LST)的遙感干旱監(jiān)測方法,在關(guān)中平原的近實(shí)時(shí)干旱監(jiān)測中具有其適用性。分位數(shù)回歸能全面反映因變量的條件分布在不同分位數(shù)處的特征,回歸結(jié)果穩(wěn)健可靠。為了進(jìn)一步研究VTCI干旱監(jiān)測結(jié)果與小麥單產(chǎn)之間的關(guān)系及提高冬小麥單產(chǎn)估測精度,構(gòu)建了不同分位數(shù)τ(0.1,0.3,0.5,0.7,0.9)下關(guān)中平原各市2008—2014年的冬小麥主要生育期VTCI與單產(chǎn)之間的線性回歸模型,并基于中位數(shù)(τ=0.5)回歸模型對研究區(qū)域的冬小麥單產(chǎn)進(jìn)行了估測。結(jié)果表明,分位數(shù)回歸模型比較全面地反映了不同分位數(shù)下冬小麥單產(chǎn)分布與VTCI之間的相關(guān)程度,彌補(bǔ)了最小二乘估產(chǎn)模型回歸結(jié)果單一、易受異常值影響等的不足。中位數(shù)回歸模型的單產(chǎn)估測結(jié)果與實(shí)際單產(chǎn)之間的相對誤差和均方根誤差的最小值及平均值均低于最小二乘回歸模型,估測精度較高。此外,中位數(shù)單產(chǎn)估測模型獲取的冬小麥估產(chǎn)結(jié)果在年際變化規(guī)律與空間分布特征上與實(shí)際產(chǎn)量均較相符,說明分位數(shù)回歸在研究VTCI與產(chǎn)量之間的關(guān)系及冬小麥單產(chǎn)估測中具有其適用性與可靠性。

    冬小麥;分位數(shù)回歸;條件植被溫度指數(shù);遙感;估產(chǎn)

    引言

    作物估產(chǎn)信息是國民經(jīng)濟(jì)宏觀調(diào)控的重要信息,及時(shí)、準(zhǔn)確、大范圍地獲取作物估產(chǎn)信息對農(nóng)業(yè)經(jīng)濟(jì)發(fā)展與糧食政策的制定具有重要的現(xiàn)實(shí)意義[1-3]。遙感技術(shù)可充分利用地物表面的光譜、時(shí)間、空間和方向信息,具有準(zhǔn)確客觀、宏觀快速的特點(diǎn),成為大面積農(nóng)業(yè)資源調(diào)查、作物估產(chǎn)、精準(zhǔn)農(nóng)業(yè)中應(yīng)用最為廣泛的技術(shù)手段之一[4]。常用的基于遙感信息的作物產(chǎn)量估測方法主要有機(jī)理模型法、經(jīng)驗(yàn)?zāi)P头鞍霗C(jī)理半經(jīng)驗(yàn)?zāi)P头ǖ龋?]。其中,經(jīng)驗(yàn)?zāi)P头ㄍㄟ^分析遙感數(shù)據(jù)與作物長勢指標(biāo)之間的線性關(guān)系及作物長勢指標(biāo)與產(chǎn)量之間的相關(guān)性實(shí)現(xiàn)作物產(chǎn)量的間接估測,是一種簡單易行的估產(chǎn)方法[5];機(jī)理模型法參數(shù)與過程復(fù)雜,目前多利用遙感數(shù)據(jù)與作物模型的耦合進(jìn)行作物產(chǎn)量的估測與預(yù)測[6];半經(jīng)驗(yàn)半機(jī)理模型中以光能利用效率模型應(yīng)用最廣泛,但存在一些關(guān)鍵參數(shù)難以定量模擬的問題[7]。

    在全球變暖的背景下,農(nóng)業(yè)干旱已經(jīng)成為影響作物生長發(fā)育與產(chǎn)量形成的主要因素[8-9]。目前,通過構(gòu)建遙感干旱指數(shù)與作物單產(chǎn)之間的關(guān)系可初步實(shí)現(xiàn)農(nóng)業(yè)干旱對糧食安全影響的定量分析[10-11]。常見的遙感干旱監(jiān)測指數(shù)主要有基于能量平衡、微波、植被狀態(tài)和熱紅外的遙感監(jiān)測方法等[12]。其中,條件植被溫度指數(shù)(Vegetation temperature condition index,VTCI)是基于歸一化植被指數(shù)(NDVI)和地表溫度(LST)的散點(diǎn)圖呈三角形區(qū)域分布的空間特征提出的一種干旱監(jiān)測方法[13-14],適用于監(jiān)測一特定年內(nèi)某一時(shí)期區(qū)域級的干旱程度,可通過構(gòu)建VTCI與作物長勢指標(biāo)之間的定量關(guān)系進(jìn)行干旱影響評估、作物產(chǎn)量估測及預(yù)測等研究[15-17]。

    在以往的作物估產(chǎn)研究中,大多應(yīng)用最小二乘法(Ordinary least squares,OLS)的經(jīng)典線性回歸進(jìn)行模型構(gòu)建與數(shù)據(jù)分析,以及遙感干旱監(jiān)測指數(shù)與作物產(chǎn)量之間的關(guān)系研究[18-20]。然而,基于經(jīng)典最小二乘法的線性回歸模型只能基于因變量與自變量在均值水平上的相關(guān)關(guān)系得到一條回歸直線,挖掘的信息量有限。此外,利用最小二乘法進(jìn)行回歸分析時(shí)對隨機(jī)誤差的分布特征要求嚴(yán)格,對于一些實(shí)際問題,最小二乘法的線性假設(shè)可能過于苛刻,難以得到無偏、有效的參數(shù)估計(jì)值[21]。為了彌補(bǔ)經(jīng)典最小二乘法在回歸分析中的缺陷,KOENKER和BASSETT于1978年提出了分位數(shù)回歸(Quantile regression,QR)的理論[22]。分位數(shù)回歸利用因變量的條件分位數(shù)建模,較之經(jīng)典的最小二乘法應(yīng)用條件更加寬松,挖掘信息更加豐富,可以更精確地描述自變量對于因變量的變化范圍以及條件分布形狀的影響,在一定程度上比最小二乘的回歸結(jié)果更加穩(wěn)?。?1,23]。這一方法最初主要應(yīng)用于經(jīng)濟(jì)領(lǐng)域,隨著計(jì)算機(jī)的普及被迅速應(yīng)用于醫(yī)學(xué)、教育、社會(huì)、環(huán)境科學(xué)等多個(gè)領(lǐng)域[24-26],但在遙感干旱影響評估與作物單產(chǎn)估測方面的應(yīng)用還未見報(bào)道。本文基于分位數(shù)回歸理論構(gòu)建關(guān)中平原冬小麥主要生育期的條件植被溫度指數(shù)與冬小麥單產(chǎn)之間的分位數(shù)回歸模型,研究不同分位數(shù)下的冬小麥單產(chǎn)受干旱的影響程度,并通過對比最小二乘法估產(chǎn)模型與分位數(shù)估產(chǎn)模型的精度,驗(yàn)證分位數(shù)回歸模型在干旱影響評估與冬小麥單產(chǎn)估測中的適用性。

    1 材料與方法

    1.1 研究區(qū)域

    研究區(qū)域?yàn)槲挥陉兾魇≈胁康年P(guān)中平原,總面積5.55×104km2,東西長約300 km,西窄東寬,地勢西高東低,土壤肥沃,墾植指數(shù)達(dá)70%以上,主要包括渭南市、西安市、咸陽市、寶雞市、銅川市5個(gè)地級市及楊凌農(nóng)業(yè)高新技術(shù)產(chǎn)業(yè)示范區(qū)。該區(qū)域地處季風(fēng)區(qū)的邊緣,氣候溫暖,年平均氣溫6~13℃,年均降水量550~700mm,降水量偏少,且降水分布西部優(yōu)于東部,在空間分布上不均勻,年際變化與年內(nèi)變化均較大[27]。尤其自20世紀(jì)90年代以來,關(guān)中平原整體呈現(xiàn)暖干化趨勢,春旱、伏旱顯著,關(guān)中平原的農(nóng)作物在全生育期呈現(xiàn)不同程度的水分不足現(xiàn)象,進(jìn)而造成該區(qū)域的糧食減產(chǎn),干旱已成為該地區(qū)農(nóng)業(yè)增產(chǎn)增收的主要制約因素[17]。

    1.2 數(shù)據(jù)處理

    VTCI已被證實(shí)是一種近實(shí)時(shí)的干旱監(jiān)測方法,其定義為[13-14]

    式中 Ni——某一像素的NDVI

    LNi——某一像素的NDVI值為Ni時(shí)的地表溫度

    LNimax——NDVI為Ni時(shí)研究區(qū)域內(nèi)所有像素地表溫度的最大值

    LNimin——NDVI為Ni時(shí)研究區(qū)域內(nèi)所有像素地表溫度的最小值

    a、b、a'、b'——待定系數(shù),由研究區(qū)域的NDVI和LST散點(diǎn)圖近似獲得

    采用關(guān)中平原 2008—2014年 3—5月份的Aqua-MODIS的日地表溫度產(chǎn)品(MYD11A1)和日地表反射率產(chǎn)品(MYD09GA)反演得到日LST和日NDVI產(chǎn)品,應(yīng)用最大值合成技術(shù),分別生成旬NDVI和LST最大值合成產(chǎn)品,通過對多年某一旬的LST最大值合成產(chǎn)品逐像素取最小值,生成多年的旬LST最大-最小值合成產(chǎn)品。經(jīng)上述方法確定冷、熱邊界后,根據(jù)VTCI的計(jì)算公式(式(1)),生成以旬為單位的VTCI時(shí)間序列數(shù)據(jù)。依據(jù)關(guān)中平原冬小麥的生長情況,將冬小麥越冬后的主要生育期劃分為返青期(3月上旬—3月中旬)、拔節(jié)期(3月下旬—4月中旬)、抽穗-灌漿期(4月下旬—5月上旬)和乳熟期(5月中旬—5月下旬)4個(gè)生育時(shí)期[19],取每個(gè)生育時(shí)期內(nèi)多旬VTCI的均值作為該生育時(shí)期的VTCI值,計(jì)算出關(guān)中平原各市2008—2014年各生育時(shí)期VTCI。然后運(yùn)用基于熵值法與因子權(quán)重排序法的歸一組合賦權(quán)法[28]確定關(guān)中平原五市(除銅川市)的冬小麥各個(gè)生育時(shí)期的權(quán)重,進(jìn)而得到關(guān)中平原冬小麥每年的加權(quán)VTCI,即冬小麥的主要生育期VTCI。

    關(guān)中平原各市2008—2014年的冬小麥單產(chǎn)數(shù)據(jù)主要來源于研究年份的陜西省統(tǒng)計(jì)年鑒。

    1.3 分位數(shù)回歸

    1.3.1 分位數(shù)回歸理論

    分位數(shù)回歸方法是將回歸方法與條件分位數(shù)進(jìn)行結(jié)合,通過最小化離差絕對值的加權(quán)和,依據(jù)因變量的條件分位數(shù)對自變量進(jìn)行回歸,進(jìn)而得到所有分位數(shù)下的回歸模型[23]。其定義為:設(shè)隨機(jī)變量Y的一組隨機(jī)變量為{y1,y2,…,yn},Y的分布函數(shù)為F(y)=P(Y≤y),則對于任意的τ(0<τ<1),Y的第τ分位數(shù)函數(shù)為

    對于任意的τ(0<τ<1),定義損失函數(shù)ρτ(z)為

    式中 I(·)——示性函數(shù)

    由式(5)可知,損失函數(shù)為分段函數(shù)(圖1),且ρτ(z)≥0。

    圖1 分位數(shù)回歸的損失函數(shù)ρτ(z)Fig.1 Loss function of quantile regression

    為方便積分,將ρτ(z)表示為

    1.3.2 線性分位數(shù)回歸

    線性分位數(shù)回歸模型的方程為

    式中 Y——因變量 β——未知參數(shù)

    ε——誤差項(xiàng)

    對于一般的線性條件分位數(shù)函數(shù),為了分析自變量X對于因變量Y在其各分位數(shù)τ上的影響,需要求解滿足

    的β(τ)。目前對式(10)的算法主要包括單純形算法、內(nèi)點(diǎn)算法與平滑算法等[23]。進(jìn)而得到 Y的τ(0<τ<1)分位數(shù)函數(shù)為

    式中 β(τ)——式(10)中極小化問題的解

    由此,在不同的τ(0<τ<1)下就能得到不同的分位數(shù)函數(shù),即所有y在x上的條件分布的一簇直線。

    2 結(jié)果與分析

    2.1 分位數(shù)回歸模型的構(gòu)建

    為了研究干旱對不同地區(qū)的小麥產(chǎn)量分布產(chǎn)生的影響,本文對關(guān)中平原四市(除銅川市)2008—2014年的冬小麥主要生育期VTCI與冬小麥單產(chǎn)的28組數(shù)據(jù)進(jìn)行分析,以小麥單產(chǎn)作為因變量(Y),以主要生育期VTCI的干旱監(jiān)測結(jié)果作為自變量(X),建立VTCI與小麥單產(chǎn)在不同分位數(shù)τ(0.1,0.3,0.5,0.7,0.9)下的線性分位數(shù)回歸模型,通過參數(shù)估計(jì)算法得到常數(shù)項(xiàng)與X系數(shù)(β(τ))(表1)。由表1可知,不同分位數(shù)下的常數(shù)項(xiàng)與X系數(shù)的估計(jì)值均落在了置信區(qū)間內(nèi),且以中位數(shù)(τ= 0.5)回歸與高分位數(shù)(τ=0.7與τ=0.9)回歸的置信區(qū)間較小,說明中位數(shù)與高分位數(shù)參數(shù)估計(jì)的置信水平高于低分位數(shù)處,回歸模型的擬合優(yōu)度較好。

    表1 分位數(shù)回歸的系數(shù)Tab.1 Coefficients of quantile regression

    以0.05為步長,估計(jì)各個(gè)分位數(shù)(0.05,0.10,0.15,…,0.95)下的常數(shù)項(xiàng)與X系數(shù),及其對應(yīng)的置信區(qū)間,得到了分位數(shù)回歸在不同分位數(shù)下的系數(shù)(圖2)。由圖2可知,常數(shù)項(xiàng)的估計(jì)值(圖2a)主要介于(400,2 000)之間,在最小二乘回歸常數(shù)項(xiàng)估計(jì)值的附近上下波動(dòng),并隨著分位數(shù)的增加呈現(xiàn)出先降低后增加的變化趨勢,置信區(qū)間在低分位數(shù)較大,在高分位處逐漸縮小。自變量X系數(shù)(圖2b)在低分位數(shù)處(τ<0.4)隨著分位數(shù)的增加呈現(xiàn)先降低后顯著增加的趨勢,并在大約τ=0.35處,其取值與最小二乘線性回歸相同;在τ∈(0.4,0.6)時(shí),X系數(shù)的估計(jì)值較大且發(fā)展趨勢平穩(wěn);當(dāng)τ>0.6時(shí),隨著分位數(shù)的增加,X系數(shù)的估計(jì)值緩慢下降,且τ=0.75處的X系數(shù)估計(jì)值與最小二乘法相同;其置信區(qū)間的變化規(guī)律與常數(shù)項(xiàng)類似,低分位數(shù)處較大,隨著分位數(shù)的增加逐漸縮小。上述分析表明,最小二乘法的線性回歸使用簡單的線性回歸參數(shù)估計(jì)值對所有的數(shù)據(jù)進(jìn)行估算,只能得到VTCI干旱監(jiān)測結(jié)果與冬小麥單產(chǎn)的條件分布在均值位置上的一條擬合直線,這樣勢必會(huì)造成產(chǎn)量估測結(jié)果的偏高或偏低。而分位數(shù)回歸在每個(gè)分位數(shù)處的系數(shù)估計(jì)值及其置信區(qū)間不盡相同,說明在冬小麥減產(chǎn)、正?;蜇S產(chǎn)的年份或地區(qū),VTCI值與冬小麥產(chǎn)量之間的線性關(guān)系不同,由此可依據(jù)研究年份間關(guān)中平原各市冬小麥的生長與干旱狀況,構(gòu)建不同分位數(shù)下的小麥估產(chǎn)模型,得到更加切合實(shí)際情況的冬小麥單產(chǎn)估測結(jié)果。

    圖2 不同分位數(shù)下的常數(shù)項(xiàng)及X系數(shù)Fig.2 Coefficients of intercept and X under different quantiles

    圖3 不同分位數(shù)下的擬合線Fig.3 Fitted lines under different quantiles

    不同分位數(shù)下的回歸直線(圖3)中,最小二乘回歸、中位數(shù)回歸(τ=0.5)以及τ為0.9、0.7、0.3、0.1下的分位數(shù)回歸直線的變化趨勢相同,均揭示了冬小麥單產(chǎn)在旱情嚴(yán)重時(shí)偏低,而旱情較輕或無旱時(shí)冬小麥單產(chǎn)較高的規(guī)律。但不同分位數(shù)下回歸直線的斜率卻不盡相同,中位數(shù)回歸時(shí)的斜率最大,高于中高分位(τ=0.7、0.9)與中低分位(τ=0.3、0.1)時(shí)的直線斜率,說明在小麥產(chǎn)量正常的年份或地區(qū),小麥產(chǎn)量與干旱之間的線性關(guān)系最為密切;中位數(shù)回歸與最小二乘法的回歸直線的位置不同,說明條件密度的不對稱性,且最小二乘回歸受到異常值(VTCI值較大時(shí)產(chǎn)量卻偏低的點(diǎn))的影響,回歸的穩(wěn)健度受到削弱,基于最小二乘法的估產(chǎn)結(jié)果會(huì)在無旱情發(fā)生(VTCI大于0.57)時(shí)偏低,而在有旱情發(fā)生(VTCI小于0.57)時(shí)略微偏高,這樣會(huì)導(dǎo)致有旱災(zāi)發(fā)生的年份或區(qū)域的冬小麥產(chǎn)量被高估,而無旱年份或區(qū)域的冬小麥產(chǎn)量被低估;分位數(shù)回歸可彌補(bǔ)最小二乘回歸的局限性,基于不同的小麥單產(chǎn)分布條件,構(gòu)建了不同的單產(chǎn)估測模型,比經(jīng)典的最小二乘回歸反映更多的局部信息,從而實(shí)現(xiàn)對干旱與冬小麥產(chǎn)量之間關(guān)系的深入分析及小麥估產(chǎn)精度的提高。

    2.2 分位數(shù)回歸模型的精度評價(jià)

    分位數(shù)回歸模型中,0.5分位點(diǎn)處的中位數(shù)回歸模型可以較好地解決最小二乘法中某些“離群值”影響回歸顯著性的問題,是一種較為穩(wěn)健的線性回歸模型[23]。同時(shí),由于0.5分位點(diǎn)處于因變量的中間位置,在對所有的數(shù)據(jù)進(jìn)行擬合時(shí)較為適宜。故本文運(yùn)用分位數(shù)回歸模型中的0.5分位點(diǎn)進(jìn)行研究區(qū)域2008—2014年的冬小麥單產(chǎn)估測,并結(jié)合實(shí)際單產(chǎn)對線性分位數(shù)估產(chǎn)模型的估測精度進(jìn)行評價(jià)。在τ=0.5處構(gòu)建的分位數(shù)回歸模型為

    中位數(shù)回歸結(jié)果t檢驗(yàn)的P值小于0.001,達(dá)到極顯著水平。類似地,令變量Ylm表示小麥單產(chǎn),變量Xlm表示VTCI的干旱監(jiān)測結(jié)果,用最小二乘法求得回歸方程,建立的最小二乘回歸模型為

    結(jié)果表明,最小二乘回歸結(jié)果的F顯著性檢驗(yàn)的P值小于0.001,達(dá)到極顯著水平。

    為了檢驗(yàn)上述兩種回歸模型對小麥單產(chǎn)估測精度的相對優(yōu)劣,分別應(yīng)用最小二乘法及分位數(shù)回歸模型計(jì)算得到關(guān)中平原四市2008—2014年的冬小麥單產(chǎn)估測值。統(tǒng)計(jì)得到了兩模型的冬小麥估測單產(chǎn)與實(shí)際單產(chǎn)之間的相對誤差和均方根誤差的最大值、最小值及平均值(表2)。由表2可知,應(yīng)用分位數(shù)回歸模型獲取的2008—2014年的冬小麥單產(chǎn)估測結(jié)果與實(shí)際單產(chǎn)之間相對誤差的最小值(0)與平均值(5.6%)均低于最小二乘回歸,分位數(shù)回歸的相對誤差最大值(26.7%)受異常值的影響略高于最小二乘的最大值(25.5%);類似地,分位數(shù)回歸的均方根誤差最大值為623.0 kg/hm2,高于最小二乘的最大值(594.9 kg/hm2),最小值僅為0,低于最小二乘回歸的最小值(25.0 kg/hm2),分位數(shù)回歸的平均值為160.3 kg/hm2,低于最小二乘的平均值(172.1 kg/hm2)。上述結(jié)果表明,分位數(shù)回歸的單產(chǎn)估測結(jié)果的相對誤差與均方根誤差整體上低于最小二乘回歸,基于分位數(shù)回歸的估產(chǎn)精度優(yōu)于最小二乘法,故分位數(shù)回歸模型可用于VTCI與冬小麥單產(chǎn)之間的相關(guān)關(guān)系及冬小麥產(chǎn)量估測等的研究。

    表2 最小二乘與分位數(shù)回歸的估測單產(chǎn)與實(shí)際單產(chǎn)之間的相對誤差與均方根誤差Tab.2 Relative errors and rootmean square errors between estimated and actual yields using the least squares and quantile regressions

    2.3 基于分位數(shù)回歸的冬小麥單產(chǎn)估測

    基于τ=0.5處的分位數(shù)回歸模型估測關(guān)中平原2008—2014年的冬小麥單產(chǎn)(圖4)。從冬小麥單產(chǎn)的空間分布特征來看,冬小麥單產(chǎn)以中部最高,西部次之,東部最低;從年際變化規(guī)律來看,關(guān)中平原2008—2014年關(guān)中平原冬小麥的單產(chǎn)呈現(xiàn)個(gè)別年份波動(dòng),總體增長的趨勢,2013年關(guān)中平原冬小麥較之前幾年減產(chǎn)嚴(yán)重,這是由于2013年關(guān)中平原降水較少,旱情偏重,造成冬小麥生育期內(nèi)的水分不足,進(jìn)而形成一定程度的減產(chǎn)。由此,基于分位數(shù)回歸模型的冬小麥單產(chǎn)估測結(jié)果與實(shí)際情況相符,模型的估測精度較高。

    3 討論

    分位數(shù)回歸較之經(jīng)典的最小二乘回歸具有其獨(dú)特的優(yōu)勢,原因在于:最小二乘回歸系數(shù)的估計(jì)為最佳線性無偏估計(jì),對誤差的分布要求嚴(yán)格,且對于存在重尾或離群值的數(shù)據(jù),其穩(wěn)健性會(huì)被削弱;此外,最小二乘回歸模型的回歸結(jié)果為基于自變量與因變量在均值水平上的相關(guān)關(guān)系得到的一條擬合直線,所反映的信息不夠豐富。分位數(shù)回歸采用加權(quán)的最小一乘回歸法,在某些情況下比經(jīng)典的最小二乘回歸法更加穩(wěn)健,且分位數(shù)回歸可以得到不同分位數(shù)下因變量在自變量上的條件分布的軌跡,有效突出局部之間的相關(guān)關(guān)系,故可以作為經(jīng)典最小二乘法線性回歸的有益補(bǔ)充。本文中,分位數(shù)回歸的方法能夠針對不同的產(chǎn)量分布條件建立各個(gè)分位點(diǎn)的回歸模型,分位數(shù)回歸結(jié)果不僅包含VTCI干旱監(jiān)測數(shù)據(jù)與冬小麥單產(chǎn)中心位置分布的相關(guān)關(guān)系,還可以度量干旱與冬小麥單產(chǎn)在各個(gè)不同位置(產(chǎn)量偏高或偏低)下的線性相關(guān)程度,因而,分位數(shù)回歸模型對不同分位點(diǎn)上干旱對冬小麥單產(chǎn)作用差異的刻畫更加全面,為干旱影響評估方面的研究提供了大量寶貴信息,這比最小二乘法表現(xiàn)出的平均水平所包含的信息量更有價(jià)值。

    圖4 關(guān)中平原2008—2014年的冬小麥單產(chǎn)估測結(jié)果Fig.4 Estimated winter wheat yields in Guanzhong Plain from 2008 to 2014

    不同分位數(shù)下的小麥單產(chǎn)與主要生育期VTCI之間回歸直線的整體發(fā)展趨勢一致,即冬小麥單產(chǎn)隨著旱情的加劇而降低,隨著旱情的減緩而增加,這是由于研究區(qū)域冬小麥的生長和最終產(chǎn)量主要受到干旱因素的制約,一般而言,某區(qū)域某一年份的干旱程度與該年份的冬小麥產(chǎn)量的高低密切相關(guān),故利用VTCI干旱監(jiān)測結(jié)果與冬小麥單產(chǎn)之間的線性關(guān)系可實(shí)現(xiàn)小麥的單產(chǎn)估測。分位數(shù)回歸模型中,τ=0.5分位數(shù)處于所有因變量的中間位置,利用該分位數(shù)對關(guān)中四市多年的冬小麥單產(chǎn)進(jìn)行估測最為可靠,其估測精度相對最小二乘法有所提高。然而,冬小麥產(chǎn)量的影響因素繁多,除干旱因素外,冬小麥產(chǎn)量的形成還受到多種非干旱因素(病蟲害、凍害、田間管理措施等)的綜合作用,只有對冬小麥產(chǎn)量的各種影響要素進(jìn)行全面考慮,才能得到更加科學(xué)、合理的產(chǎn)量估測結(jié)果,這些問題是今后通過分位數(shù)回歸模型進(jìn)行小麥單產(chǎn)估測研究的重點(diǎn),可嘗試通過多元的分位數(shù)回歸模型等來解決。

    4 結(jié)論

    以關(guān)中平原四市2008—2014年的冬小麥主要生育期的VTCI干旱監(jiān)測結(jié)果作為自變量,冬小麥單產(chǎn)作為因變量,分別建立最小二乘回歸模型與不同分位點(diǎn)(0.1,0.3,0.5,0.7,0.9)下的分位數(shù)回歸模型,對比分析中位數(shù)回歸、最小二乘回歸與其他分位數(shù)下的回歸結(jié)果及模型精度,并基于中位數(shù)回歸與最小二乘回歸構(gòu)建的單產(chǎn)估測模型對研究區(qū)域的冬小麥單產(chǎn)進(jìn)行估測。主要結(jié)論如下:

    (1)分位數(shù)回歸得到了冬小麥單產(chǎn)在VTCI干旱監(jiān)測結(jié)果下的條件分布的一簇直線,能夠比較全面的度量不同分位數(shù)下干旱與冬小麥單產(chǎn)之間的相關(guān)關(guān)系。而且,分位數(shù)回歸采用加權(quán)的最小一乘回歸法,對異常值不敏感,建立的估產(chǎn)模型更具可靠性。由此,分位數(shù)回歸克服了普通線性回歸在回歸結(jié)果單一、易受異常值的影響等方面的不足,適用于關(guān)中平原冬小麥干旱影響評估的研究中。

    (2)基于中位數(shù)回歸構(gòu)建的關(guān)中平原冬小麥的單產(chǎn)估測模型的精度較高,估測單產(chǎn)與實(shí)際單產(chǎn)之間的相對誤差均值為5.6%,最小值為0;均方根誤差的均值為160.3 kg/hm2,最小值僅為0,估產(chǎn)精度整體上優(yōu)于最小二乘回歸的結(jié)果。此外,應(yīng)用中位數(shù)回歸模型獲取的冬小麥估產(chǎn)結(jié)果在年際變化規(guī)律與空間分布特征上與實(shí)際產(chǎn)量均較相符,說明中位數(shù)回歸及分位數(shù)回歸在研究干旱與產(chǎn)量之間的關(guān)系及冬小麥單產(chǎn)估測中具有其適用性。

    1 KOWALIKW,DABROWSKA-ZIELINSKA K,MERONIM,et al.Yield estimation using SPOT-VEGETATION products:a case study of wheat in European countries[J].International Journal of Applied Earth Observation and Geoinformation,2014,32:228-239.

    2 趙春江.農(nóng)業(yè)遙感研究與應(yīng)用進(jìn)展[J/OL].農(nóng)業(yè)機(jī)械學(xué)報(bào),2014,45(12):277-293.http:∥www.j-csam.org/jcsam/ch/ reader/view_abstract.a(chǎn)spx?flag=1&file_no=20141241&journal_id=jcsam.DOI:10.6041/j.issn.1000-1298.2014.12.041.ZHAO Chunjiang.Advances of research and application in remote sensing for agriculture[J/OL].Transactions of the Chinese Society for Agricultural Machinery,2014,45(12):277-293.(in Chinese)

    3 DORAISWAMY PC,SINCLAIR T R,HOLLINGER S,et al.Application of MODIS derived parameters for regional crop yield assessment[J].Remote Sensing of Environment,2005,97(2):192-202.

    4 楊鶴松,王鵬新,孫威.條件植被溫度指數(shù)在華北平原干旱監(jiān)測中的應(yīng)用[J].北京師范大學(xué)學(xué)報(bào):自然科學(xué)版,2007,43(3):314-318.YANG Hesong,WANG Pengxin,SUN Wei.Application of the vegetation temperature condition index to drought monitoring in North China Plain[J].Journal of Beijing Normal University:Natural Science,2007,43(3):314-318.(in Chinese)

    5 李軍玲,郭其樂,彭記永.基于MODIS數(shù)據(jù)的河南省冬小麥產(chǎn)量遙感估算模型[J].生態(tài)環(huán)境學(xué)報(bào),2012,21(10):1665-1669.LI Junling,GUOQile,PENG Jiyong.Remote sensing estimationmodel of Henan Provincewinterwheat yield based on MODISdata[J].Ecology and Environmental Sciences,2012,21(10):1665-1669.(in Chinese)

    6 DEWIT A,DUVEILLER G,DEFOURNY P.Estimating regional winter wheat yield with WOFOST through the assimilation of green area index retrieved from MODIS observations[J].Agricultural and Forest Meteorology,2012,164:39-52.

    7 BASTIAANSSENW G M,ALIS.A new crop yield forecasting model based on satellite measurements applied across the Indus Basin,Pakistan[J].Agriculture,Ecosystems&Environment,2003,94(3):321-340.

    8 WU J J,ZHOU L,MO X Y,etal.Droughtmonitoring and analysis in China based on the integrated surface drought index(ISDI)[J].International Journal of Applied Earth Observation and Geoinformation,2015,41:23-33.

    9 李芬,于文金,張建新,等.干旱災(zāi)害評估研究進(jìn)展[J].地理科學(xué)進(jìn)展,2011,28(7):891-898.LIFen,YUWenjin,ZHANG Jianxin,et al.Review of drought disaster evaluation[J].Progress in Geography,2011,28(7): 891-898.(in Chinese)

    10 史舟,梁宗正,楊媛媛,等.農(nóng)業(yè)遙感研究現(xiàn)狀與展望[J/OL].農(nóng)業(yè)機(jī)械學(xué)報(bào),2015,46(2):247-260.http:∥www.jcsam.org/jcsam/ch/reader/view_abstract.a(chǎn)spx?flag=1&file_no=20150237&journal_id=jcsam.DOI:10.6041/j.issn.1000-1298.2015.02.037.SHIZhou,LIANG Zongzheng,YANG Yuanyuan,et al.Status and prospect of agriculture remote sensing[J/OL].Transactions of the Chinese Society for Agricultural Machinery,2015,46(2):247-260.(in Chinese)

    11 李文娟,覃志豪,林綠.農(nóng)業(yè)旱災(zāi)對國家糧食安全影響程度的定量分析[J].自然災(zāi)害學(xué)報(bào),2010,19(3):111-118.LIWenjuan,QIN Zhihao,LIN Lü.Quantitative analysis of agro-drought impact on food security in China[J].Journal of Nature Disasters,2010,19(3):111-118.(in Chinese)

    12 周磊,武建軍,張潔.以遙感為基礎(chǔ)的干旱監(jiān)測方法研究進(jìn)展[J].地理科學(xué),2014,28(1):85-91.ZHOU Lei,WU Jianjun,ZHANG Jie.Remote sensing based droughtmonitoring approach and research progress[J].Scientia Geographica Sinica,2014,28(1):85-91.(in Chinese)

    13 王鵬新,龔健雅,李小文.條件植被溫度指數(shù)及其在干旱監(jiān)測中的應(yīng)用[J].武漢大學(xué)學(xué)報(bào):信息科學(xué)版,2001,26(5): 412-418.WANG Pengxin,GONG Jianya,LIXiaowen.Vegetation temperature condition index and its application for droughtmonitoring[J].Geomatics and Information Science ofWuhan University,2001,26(5):412-418.(in Chinese)

    14 SUNW,WANG P X,ZHANG SY,et al.Using the vegetation temperature condition index for time series drought occurrence monitoring in the Guanzhong Plain,PR China[J].International Journal of Remote Sensing,2008,29(17):5133-5144.

    15 陳陽,范建容,郭芬芬,等.條件植被溫度指數(shù)在云南干旱監(jiān)測中的應(yīng)用[J].農(nóng)業(yè)工程學(xué)報(bào),2011,27(5):231-236.CHEN Yang,F(xiàn)AN Jianrong,GUO Fenfen,et al.Application of the vegetation temperature condition index to droughtmonitoring in Yunnan Province[J].Transactions of the CSAE,2011,27(5):231-236.(in Chinese)

    16 李艷,王鵬新,劉峻明,等.基于條件植被溫度指數(shù)的冬小麥主要生育時(shí)期干旱監(jiān)測效果評價(jià)I:因子權(quán)重排序法和熵值法組合賦權(quán)[J].干旱地區(qū)農(nóng)業(yè)研究,2013,31(6):159-163.LIYan,WANG Pengxin,LIU Junming,etal.Evaluation of droughtmonitoring effects in themain growth and developmentstages ofwinterwheatusing vegetation temperature condition index.I:factorweight sortingmethod and entropymethod[J].Agricultural Research in the Arid Areas,2013,31(6):159-163.(in Chinese)

    17 田苗,王鵬新,張樹譽(yù),等.基于條件植被溫度指數(shù)的冬小麥產(chǎn)量預(yù)測[J/OL].農(nóng)業(yè)機(jī)械學(xué)報(bào),2014,45(2):239-245.http:∥www.j-csam.org/jcsam/ch/reader/view_abstract.a(chǎn)spx?flag=1&file_no=20140240&journal_id=jcsam.DOI:10.6041/j.issn.1000-1298.2014.02.040.TIAN Miao,WANG Pengxin,ZHANG Shuyu,et al.Winter wheat yield forecasting based on vegetation temperature condition index[J/OL].Transactions of the Chinese Society for Agricultural Machinery,2014,45(2):239-245.(in Chinese)

    18 張順謙,卿清濤,侯美亭,等.基于溫度植被干旱指數(shù)的四川伏旱遙感監(jiān)測與影響評估[J].農(nóng)業(yè)工程學(xué)報(bào),2007,23(9):141-146.ZHANG Shunqian,QING Qingtao,HOU Meiting,et al.Remote sensing and impact estimation for Sichuan hot-drought based on temperature vegetation dryness index[J].Transactions of the CSAE,2007,23(9):141-146.(in Chinese)19 黃弘,王鵬新,李俐.關(guān)中平原小麥生育期VTCI加權(quán)估算及其與產(chǎn)量的相關(guān)性研究[J].干旱地區(qū)農(nóng)業(yè)研究,2011,29(6):173-178.HUANG Hong,WANG Pengxin,LILi.Correlations between weighted VTCIin key growth and developmentstages ofwinterwheat and wheat yields in the Guanzhong Plain[J].Agricultural Research in the Arid Areas,2011,29(6):173-178.(in Chinese)

    20 黃健熙,羅倩,劉曉暄,等.基于時(shí)間序列MODISNDVI的冬小麥產(chǎn)量預(yù)測方法[J/OL].農(nóng)業(yè)機(jī)械學(xué)報(bào),2016,47(2): 295-301.http:∥www.j-csam.org/jcsam/ch/reader/view_abstract.a(chǎn)spx?flag=1&file_no=20160239&journal_id=jcsam.DOI:10.6041/j.issn.1000-1298.2016.02.039.HUANG Jianxi,LUO Qian,LIU Xiaoxuan,et al.Winter wheat yield forecasting based on time series of MODISNDVI[J/OL].Transactions of the Chinese Society for Agricultural Machinery,2016,47(2):295-301.(in Chinese)

    21 羅玉波.分位數(shù)回歸模型及其應(yīng)用[M].北京:知識產(chǎn)權(quán)出版社,2009.

    22 KOENKER R,BASSETTG J.Regression quantiles[J].Journal of the Econometric Society,1987,46(1):33-50.

    23 陳建寶,丁軍軍.分位數(shù)回歸技術(shù)綜述[J].統(tǒng)計(jì)與信息論壇,2008,23(3):89-96.CHEN Jianbao,DING Junjun.A review of technologies on quantile regression[J].Statistics&Information Forum,2008,23 (3):89-96.(in Chinese)

    24 FOURNIER JM,KOSKE I.Public employment and earnings inequality:an analysis based on conditional and unconditional quantile regressions[J].Economics Letters,2013,121:263-266.

    25 ARUNRAJN S,AHRENSD.A hybrid seasonal autoregressive integrated moving average and quantile regression for daily food sales forecasting[J].International Journal of Production Economics,2015,170:321-335.

    26 ROCCHINID,CADE B S.Quantile regression applied to spectral distance decay[J].IEEE Geoscience and Remote Sensing Letters,2008,5(4):640-643.

    27 賀音,張聰娥,張黎.基于SPEI的陜西近40年干旱時(shí)空特性分析[J].陜西氣象,2014(5):26-32.HE Yin,ZHANG Conge,ZHANG Li.Analyzing the temporal and spatial characteristics of drought in Shaanxi Province in recent 40 years based on standardized precipitation evapotranspiration index[J].Journalof ShaanxiMeteorology,2014(5):26-32.(in Chinese)

    28 王蕾.基于條件植被溫度指數(shù)的冬小麥主要生育期干旱影響評估研究[D].北京:中國農(nóng)業(yè)大學(xué),2015.WANG Lei.Estimation of drought impacton themain growth stage ofwinterwheatbased on vegetation temperature condition index[D].Beijing:China Agricultural University,2015.(in Chinese)

    Winter Wheat Yield Estimation Method Based on Quantile Regression Model and Remotely Sensed Vegetation Tem perature Condition Index

    WANG Lei1,2WANG Pengxin1,2LILi1,2ZHANG Shuyu3

    (1.College of Information and Electrical Engineering,China Agricultural University,Beijing 100083,China 2.Key Laboratory of Remote Sensing for Agri-Hazards,Ministry of Agriculture,Beijing 100083,China 3.Shaanxi Provincial Meteorological Bureau,Xi’an 710014,China)

    Vegetation temperature condition index(VTCI)combines normalized difference vegetation index(NDVI)and land surface temperature(LST),and is applicable to amore accuratemonitoring of droughts in Guanzhong Plain,Shaanxi Province,China.Quantile regression is a tool for comprehensively reflecting the conditional distribution characters under different quantiles,and its regression results are steady and reliable.In order to achieve a better correlation between winter wheat yield and the weighted VTCIaswell as a higher yield estimation accuracy,linear regressionmodels between the weighted VTCI and yields in the cities of Guanzhong Plain in the years from 2008 to 2014 were analyzed by using the quantile regression whose quantiles were set to be 0.1,0.3,0.5,0.7 and 0.9,respectively.These quantile regression results roundly reflected the distribution of the yields under different drought conditions and were beneficial supplement of the linear regression from which the single fitted line and impressionable results from outlierswere obtained.Thewheat yield estimationmodel based on themedian regression(quantile equalled to 0.5)was used to monitor the wheat yields in the cities of Guanzhong Plain from 2008 to 2014,the average andminimum values of the relative errors and the rootmean square errors(RMSE)between the estimated yields and the actual yieldswere all lower than those derived from the ordinary least square method.Additionally,the characteristics of inter-annual evolution and spatial distribution of the estimated yields using the median regression model were in good agreement with theactual situation,which indicated that the quantile regression was feasible and reliable in the research of winter wheat yield estimation and the relationship between yield and drought.

    winter wheat;quantile regression;vegetation temperature condition index;remote sensing; yield estimation

    S127;TP79

    A

    1000-1298(2017)07-0167-07

    2016-11-21

    2016-12-14

    國家自然科學(xué)基金項(xiàng)目(41371390)

    王蕾(1988—),女,博士生,主要從事定量遙感及其在干旱預(yù)測中的應(yīng)用研究,E-mail:409118258@qq.com

    王鵬新(1965—),男,教授,博士生導(dǎo)師,主要從事定量遙感及其在農(nóng)業(yè)中的應(yīng)用研究,E-mail:wangpx@cau.edu.cn

    10.6041/j.issn.1000-1298.2017.07.021

    猜你喜歡
    估產(chǎn)關(guān)中平原位數(shù)
    關(guān)中平原人為土形成的歷史探析
    基于無人機(jī)多光譜遙感數(shù)據(jù)的煙草植被指數(shù)估產(chǎn)模型研究
    五次完全冪的少位數(shù)三進(jìn)制展開
    遙感技術(shù)在大豆種植情況監(jiān)測中的應(yīng)用
    基于三角模型的關(guān)中5市土地生態(tài)安全區(qū)域差異分析
    《關(guān)中平原城市群發(fā)展規(guī)劃》獲批發(fā)布
    新西部(2018年3期)2018-03-21 10:09:08
    在廣仁寺
    基于地級市的區(qū)域水稻遙感估產(chǎn)與空間化研究
    基于SAR技術(shù)的高原山區(qū)煙草估產(chǎn)模型
    遙感衛(wèi)星CCD相機(jī)量化位數(shù)的選擇
    午夜免费观看网址| 窝窝影院91人妻| 桃红色精品国产亚洲av| 少妇的逼好多水| 久久精品综合一区二区三区| 99久久精品热视频| 午夜免费成人在线视频| 日韩精品中文字幕看吧| 国产色婷婷99| 久久久久久大精品| 黄色视频,在线免费观看| 亚洲av日韩精品久久久久久密| 久久久国产精品麻豆| 一个人免费在线观看电影| 一边摸一边抽搐一进一小说| 成人av在线播放网站| 亚洲成a人片在线一区二区| 亚洲avbb在线观看| 看免费av毛片| 国产伦一二天堂av在线观看| 国产探花极品一区二区| 久久精品91蜜桃| 日韩人妻高清精品专区| 日韩欧美国产一区二区入口| www日本黄色视频网| 久久久成人免费电影| 国产欧美日韩精品亚洲av| 国产欧美日韩一区二区三| 国产精品香港三级国产av潘金莲| 国产不卡一卡二| 欧美bdsm另类| 国产精品久久久人人做人人爽| 草草在线视频免费看| 国产主播在线观看一区二区| 久久亚洲精品不卡| 丰满人妻一区二区三区视频av | a级毛片a级免费在线| 观看美女的网站| 欧美在线一区亚洲| 真人一进一出gif抽搐免费| 国产精品亚洲av一区麻豆| 一区二区三区国产精品乱码| 午夜免费成人在线视频| 午夜亚洲福利在线播放| 国产精品 欧美亚洲| 欧美不卡视频在线免费观看| 久久精品国产亚洲av香蕉五月| 欧美在线一区亚洲| 色尼玛亚洲综合影院| 国产精品一区二区三区四区免费观看 | 一a级毛片在线观看| 日韩欧美免费精品| 午夜福利成人在线免费观看| 变态另类成人亚洲欧美熟女| 欧美黑人欧美精品刺激| 啪啪无遮挡十八禁网站| 最近视频中文字幕2019在线8| 久久精品国产综合久久久| 亚洲欧美日韩高清在线视频| 国产中年淑女户外野战色| 最近最新中文字幕大全免费视频| 精品免费久久久久久久清纯| 亚洲aⅴ乱码一区二区在线播放| 18禁美女被吸乳视频| 成年人黄色毛片网站| 久9热在线精品视频| 久久性视频一级片| 国产美女午夜福利| 麻豆国产av国片精品| 国产视频内射| 亚洲不卡免费看| 制服人妻中文乱码| 内射极品少妇av片p| 亚洲成人免费电影在线观看| 日本与韩国留学比较| 丁香欧美五月| 久久久久亚洲av毛片大全| 国产aⅴ精品一区二区三区波| 亚洲国产精品久久男人天堂| 特级一级黄色大片| 两个人的视频大全免费| 亚洲精品日韩av片在线观看 | 日本精品一区二区三区蜜桃| 久久精品影院6| 免费看美女性在线毛片视频| 最近最新免费中文字幕在线| tocl精华| 午夜激情福利司机影院| 亚洲第一欧美日韩一区二区三区| av黄色大香蕉| 免费一级毛片在线播放高清视频| 国产精品乱码一区二三区的特点| 国产午夜福利久久久久久| 精品国产美女av久久久久小说| 黄片大片在线免费观看| 国产精品野战在线观看| 午夜免费成人在线视频| 天天躁日日操中文字幕| 亚洲人成网站高清观看| 国产野战对白在线观看| 欧美zozozo另类| www.熟女人妻精品国产| 12—13女人毛片做爰片一| 国产亚洲欧美在线一区二区| 欧美中文日本在线观看视频| 搡老熟女国产l中国老女人| 一区二区三区高清视频在线| 国产一区在线观看成人免费| 午夜免费观看网址| 午夜亚洲福利在线播放| 国产成人av教育| 精品人妻一区二区三区麻豆 | 亚洲色图av天堂| 午夜影院日韩av| 88av欧美| 男人舔女人下体高潮全视频| 在线观看av片永久免费下载| 国产三级中文精品| 天堂动漫精品| 99riav亚洲国产免费| 国产精品亚洲av一区麻豆| 欧美最黄视频在线播放免费| h日本视频在线播放| 嫩草影视91久久| 国产私拍福利视频在线观看| av视频在线观看入口| 又爽又黄无遮挡网站| 9191精品国产免费久久| 精品国产亚洲在线| 国产91精品成人一区二区三区| 亚洲七黄色美女视频| aaaaa片日本免费| 脱女人内裤的视频| 国产视频一区二区在线看| 国产不卡一卡二| 宅男免费午夜| 国产精品综合久久久久久久免费| 日本 av在线| 亚洲av一区综合| 99在线视频只有这里精品首页| 久久人妻av系列| 亚洲精品一卡2卡三卡4卡5卡| 3wmmmm亚洲av在线观看| 国产一区二区在线av高清观看| 中国美女看黄片| 国产精品一及| 亚洲精品国产精品久久久不卡| 日本 欧美在线| 熟女少妇亚洲综合色aaa.| 91九色精品人成在线观看| 国产熟女xx| 99riav亚洲国产免费| 99久久无色码亚洲精品果冻| 99在线视频只有这里精品首页| 国产69精品久久久久777片| 啦啦啦韩国在线观看视频| 亚洲精品一卡2卡三卡4卡5卡| 精品久久久久久,| 88av欧美| 久久九九热精品免费| 国产亚洲欧美98| 亚洲天堂国产精品一区在线| 国产午夜精品论理片| 免费av不卡在线播放| 国产爱豆传媒在线观看| 久久99热这里只有精品18| 叶爱在线成人免费视频播放| 免费看十八禁软件| 夜夜夜夜夜久久久久| 亚洲人成网站在线播| 久久久久免费精品人妻一区二区| 欧美区成人在线视频| 久久精品国产自在天天线| av视频在线观看入口| 最新在线观看一区二区三区| 美女被艹到高潮喷水动态| 国产亚洲精品久久久久久毛片| 人人妻人人澡欧美一区二区| 亚洲av不卡在线观看| 欧美3d第一页| 日本撒尿小便嘘嘘汇集6| 中文在线观看免费www的网站| www日本在线高清视频| www日本在线高清视频| 日本成人三级电影网站| 不卡一级毛片| 国产免费一级a男人的天堂| 日韩欧美在线乱码| 欧美乱色亚洲激情| 最好的美女福利视频网| 国产精品自产拍在线观看55亚洲| 久久久久国内视频| 欧美极品一区二区三区四区| 99国产精品一区二区三区| 神马国产精品三级电影在线观看| 国产av不卡久久| 性色avwww在线观看| 特大巨黑吊av在线直播| 真人做人爱边吃奶动态| 亚洲色图av天堂| 91字幕亚洲| 香蕉av资源在线| 久久精品91无色码中文字幕| 中出人妻视频一区二区| 日日摸夜夜添夜夜添小说| 中文亚洲av片在线观看爽| av欧美777| 男女之事视频高清在线观看| 三级毛片av免费| a级一级毛片免费在线观看| 99久久久亚洲精品蜜臀av| 九色成人免费人妻av| 午夜影院日韩av| 国产真人三级小视频在线观看| 在线观看舔阴道视频| 国内精品久久久久精免费| 老司机深夜福利视频在线观看| 老熟妇仑乱视频hdxx| 久久精品国产亚洲av香蕉五月| 亚洲av第一区精品v没综合| 亚洲性夜色夜夜综合| 日本免费a在线| 欧美激情在线99| 一进一出抽搐动态| 中文字幕高清在线视频| 国产野战对白在线观看| 一夜夜www| 欧美色欧美亚洲另类二区| 亚洲在线观看片| 老汉色∧v一级毛片| 男女床上黄色一级片免费看| 啪啪无遮挡十八禁网站| 少妇人妻精品综合一区二区 | 宅男免费午夜| av在线天堂中文字幕| 99riav亚洲国产免费| 老司机福利观看| 色综合婷婷激情| 国产精品99久久99久久久不卡| 国产毛片a区久久久久| 少妇高潮的动态图| netflix在线观看网站| 夜夜躁狠狠躁天天躁| 高潮久久久久久久久久久不卡| 丁香六月欧美| 淫妇啪啪啪对白视频| 99久国产av精品| 真人一进一出gif抽搐免费| 美女 人体艺术 gogo| 国产精品 欧美亚洲| 真实男女啪啪啪动态图| 九色国产91popny在线| 国产淫片久久久久久久久 | 欧美一级a爱片免费观看看| 99久久成人亚洲精品观看| 国产亚洲精品一区二区www| 亚洲av成人不卡在线观看播放网| 变态另类丝袜制服| 国产成人aa在线观看| 99热这里只有是精品50| 欧美三级亚洲精品| 亚洲精品色激情综合| 日韩亚洲欧美综合| 天堂网av新在线| 我的老师免费观看完整版| 国产国拍精品亚洲av在线观看 | 国产精品野战在线观看| ponron亚洲| 嫩草影院精品99| 欧美黄色淫秽网站| 少妇熟女aⅴ在线视频| 搡老岳熟女国产| 久久久久久久精品吃奶| 中出人妻视频一区二区| 亚洲最大成人手机在线| 成人18禁在线播放| 亚洲人成网站在线播放欧美日韩| 国产精品免费一区二区三区在线| a在线观看视频网站| 国产高清videossex| 色哟哟哟哟哟哟| 亚洲精品国产精品久久久不卡| 少妇的丰满在线观看| 黄色日韩在线| 亚洲,欧美精品.| 精品熟女少妇八av免费久了| 欧美绝顶高潮抽搐喷水| 人人妻,人人澡人人爽秒播| 99国产综合亚洲精品| 99久国产av精品| 麻豆成人av在线观看| 美女黄网站色视频| 精品熟女少妇八av免费久了| 中文资源天堂在线| 一a级毛片在线观看| 高清毛片免费观看视频网站| 怎么达到女性高潮| av天堂在线播放| 日日干狠狠操夜夜爽| 国产精品 欧美亚洲| 午夜福利在线观看免费完整高清在 | 国产精品,欧美在线| 成年人黄色毛片网站| 狂野欧美激情性xxxx| 一本一本综合久久| 亚洲第一电影网av| 波野结衣二区三区在线 | 免费搜索国产男女视频| 国产成人av激情在线播放| 婷婷精品国产亚洲av在线| 丝袜美腿在线中文| 午夜久久久久精精品| 久久香蕉精品热| 岛国视频午夜一区免费看| 午夜激情福利司机影院| 人人妻人人看人人澡| 午夜免费成人在线视频| 搡老妇女老女人老熟妇| 黄色丝袜av网址大全| 悠悠久久av| 免费大片18禁| 国产淫片久久久久久久久 | 国产乱人伦免费视频| 免费无遮挡裸体视频| 波野结衣二区三区在线 | 婷婷亚洲欧美| 国产爱豆传媒在线观看| 一a级毛片在线观看| 噜噜噜噜噜久久久久久91| 午夜精品在线福利| 黄色丝袜av网址大全| 欧美日本视频| 亚洲成人免费电影在线观看| 午夜福利欧美成人| 最近视频中文字幕2019在线8| 久久精品91无色码中文字幕| 老司机午夜十八禁免费视频| 伊人久久精品亚洲午夜| 亚洲美女视频黄频| 69av精品久久久久久| 男女那种视频在线观看| 好男人电影高清在线观看| 国产亚洲欧美在线一区二区| 亚洲avbb在线观看| 国产午夜福利久久久久久| 成年女人永久免费观看视频| 国产欧美日韩一区二区精品| 免费观看精品视频网站| 久久久久久人人人人人| 一个人免费在线观看电影| 国产精品亚洲av一区麻豆| 精品福利观看| 免费电影在线观看免费观看| 欧美日韩福利视频一区二区| 国产熟女xx| 男女做爰动态图高潮gif福利片| 男插女下体视频免费在线播放| 国产欧美日韩精品一区二区| 欧美又色又爽又黄视频| 日本a在线网址| 免费观看精品视频网站| 中文字幕人成人乱码亚洲影| 欧美zozozo另类| 露出奶头的视频| 亚洲人与动物交配视频| 亚洲欧美日韩东京热| 日韩欧美一区二区三区在线观看| 男人和女人高潮做爰伦理| 99久久九九国产精品国产免费| 久久精品综合一区二区三区| 99久久成人亚洲精品观看| 每晚都被弄得嗷嗷叫到高潮| 国产一区二区三区在线臀色熟女| 亚洲国产中文字幕在线视频| 三级国产精品欧美在线观看| 一个人免费在线观看电影| 日本一本二区三区精品| 国产一区二区在线观看日韩 | 成年免费大片在线观看| 精品一区二区三区人妻视频| 我要搜黄色片| 国产色爽女视频免费观看| 人人妻人人澡欧美一区二区| 午夜精品在线福利| 手机成人av网站| 亚洲性夜色夜夜综合| 怎么达到女性高潮| 99久久精品热视频| 久久久国产成人精品二区| 午夜福利高清视频| 桃红色精品国产亚洲av| 国产精品一区二区三区四区免费观看 | 日本黄大片高清| 亚洲欧美一区二区三区黑人| 99久久精品一区二区三区| 黄色女人牲交| 欧美一级毛片孕妇| 日本免费一区二区三区高清不卡| 午夜免费激情av| 老熟妇乱子伦视频在线观看| 亚洲狠狠婷婷综合久久图片| 国产色婷婷99| 亚洲激情在线av| 亚洲第一欧美日韩一区二区三区| 亚洲av成人不卡在线观看播放网| 久久九九热精品免费| 国产精品98久久久久久宅男小说| 欧美性猛交╳xxx乱大交人| 99热精品在线国产| 精品国产美女av久久久久小说| 色综合欧美亚洲国产小说| 婷婷精品国产亚洲av在线| 禁无遮挡网站| 亚洲色图av天堂| 叶爱在线成人免费视频播放| 男女那种视频在线观看| av中文乱码字幕在线| 精品熟女少妇八av免费久了| 一区二区三区激情视频| 69av精品久久久久久| 亚洲国产日韩欧美精品在线观看 | 免费在线观看日本一区| 91麻豆精品激情在线观看国产| 精品99又大又爽又粗少妇毛片 | 韩国av一区二区三区四区| 变态另类丝袜制服| 国产高清三级在线| 国产精品久久久人人做人人爽| 99久久综合精品五月天人人| 男女午夜视频在线观看| 欧美日韩一级在线毛片| 精品一区二区三区视频在线观看免费| 午夜福利在线观看吧| 国产精品一区二区三区四区免费观看 | 成人午夜高清在线视频| 欧美日韩国产亚洲二区| a级毛片a级免费在线| 国产激情欧美一区二区| 天堂动漫精品| 国产色爽女视频免费观看| www.色视频.com| 99热这里只有精品一区| 亚洲成人精品中文字幕电影| 很黄的视频免费| 久久精品人妻少妇| 精品欧美国产一区二区三| 国产亚洲精品久久久久久毛片| 嫩草影院入口| 亚洲色图av天堂| 69av精品久久久久久| 亚洲va日本ⅴa欧美va伊人久久| 国产精品女同一区二区软件 | 欧美性猛交黑人性爽| 久久精品国产99精品国产亚洲性色| 午夜免费成人在线视频| 非洲黑人性xxxx精品又粗又长| 精品电影一区二区在线| 1024手机看黄色片| 夜夜看夜夜爽夜夜摸| 最后的刺客免费高清国语| av福利片在线观看| 欧美成人a在线观看| 白带黄色成豆腐渣| а√天堂www在线а√下载| 国产精品99久久久久久久久| 国产精品爽爽va在线观看网站| 真人做人爱边吃奶动态| 成年女人永久免费观看视频| 成人欧美大片| 国产精品三级大全| 免费一级毛片在线播放高清视频| 国产成人福利小说| 人妻丰满熟妇av一区二区三区| 欧美黑人欧美精品刺激| 国产在线精品亚洲第一网站| 99热只有精品国产| 久久久久久久亚洲中文字幕 | 乱人视频在线观看| 在线播放无遮挡| 最近视频中文字幕2019在线8| 亚洲午夜理论影院| 极品教师在线免费播放| 99久久99久久久精品蜜桃| 国产一区二区三区在线臀色熟女| 欧美高清成人免费视频www| 国产一区二区三区在线臀色熟女| 亚洲成av人片在线播放无| 天天添夜夜摸| 亚洲国产日韩欧美精品在线观看 | 内射极品少妇av片p| 性色av乱码一区二区三区2| 国产探花极品一区二区| 日日干狠狠操夜夜爽| 99久久九九国产精品国产免费| 国产精品久久久久久人妻精品电影| 国产黄片美女视频| 一区二区三区免费毛片| 九九在线视频观看精品| 免费观看的影片在线观看| 亚洲中文字幕一区二区三区有码在线看| 亚洲最大成人手机在线| 老汉色∧v一级毛片| 久久性视频一级片| 精品一区二区三区视频在线观看免费| tocl精华| ponron亚洲| 99热6这里只有精品| 2021天堂中文幕一二区在线观| av片东京热男人的天堂| 国产三级在线视频| 美女高潮的动态| 丁香欧美五月| 精品久久久久久成人av| 一区福利在线观看| 国产在视频线在精品| 国产成年人精品一区二区| 精品国内亚洲2022精品成人| 久久久久久久午夜电影| 日韩欧美在线乱码| 免费观看人在逋| 又黄又爽又免费观看的视频| 91久久精品国产一区二区成人 | 99在线人妻在线中文字幕| 国产激情偷乱视频一区二区| 久久久精品欧美日韩精品| 国产成+人综合+亚洲专区| av视频在线观看入口| av天堂中文字幕网| 嫁个100分男人电影在线观看| 老汉色∧v一级毛片| 黄色片一级片一级黄色片| 亚洲精品一区av在线观看| 老司机午夜福利在线观看视频| 久久精品夜夜夜夜夜久久蜜豆| 国产av一区在线观看免费| 亚洲,欧美精品.| 99精品在免费线老司机午夜| 欧美一区二区国产精品久久精品| 欧美乱码精品一区二区三区| 亚洲欧美一区二区三区黑人| 波野结衣二区三区在线 | 色综合婷婷激情| 一个人免费在线观看电影| 国产乱人伦免费视频| 日本一二三区视频观看| 日韩免费av在线播放| 国产亚洲欧美98| 午夜a级毛片| 国产高清三级在线| 国产精品日韩av在线免费观看| 亚洲成av人片在线播放无| 91在线精品国自产拍蜜月 | 午夜免费激情av| 欧美一区二区精品小视频在线| av片东京热男人的天堂| 色噜噜av男人的天堂激情| 国产综合懂色| 99在线人妻在线中文字幕| av专区在线播放| 久久精品综合一区二区三区| 香蕉久久夜色| 狠狠狠狠99中文字幕| 波野结衣二区三区在线 | 婷婷丁香在线五月| 国产真实伦视频高清在线观看 | 两个人视频免费观看高清| 日本成人三级电影网站| 午夜精品久久久久久毛片777| 人人妻人人澡欧美一区二区| 精品99又大又爽又粗少妇毛片 | 欧美区成人在线视频| 麻豆成人午夜福利视频| 午夜影院日韩av| 欧美国产日韩亚洲一区| 一二三四社区在线视频社区8| 久久精品国产综合久久久| 丝袜美腿在线中文| 国产极品精品免费视频能看的| 一级黄片播放器| 国产精品香港三级国产av潘金莲| xxx96com| 国产成人欧美在线观看| 他把我摸到了高潮在线观看| 丁香六月欧美| а√天堂www在线а√下载| 99视频精品全部免费 在线| 亚洲人成电影免费在线| 亚洲国产精品999在线| 日韩欧美国产在线观看| 成年女人永久免费观看视频| 免费看光身美女| 午夜免费激情av| 人人妻,人人澡人人爽秒播| 精品久久久久久成人av| 亚洲av成人av| 国产成+人综合+亚洲专区| 免费看十八禁软件| 大型黄色视频在线免费观看| 亚洲一区高清亚洲精品| 国产精品一及| 欧美高清成人免费视频www| 欧美黑人巨大hd| 久久久久久九九精品二区国产| 欧美黄色淫秽网站| 久久久久久久久大av| 国产精品野战在线观看| 午夜亚洲福利在线播放| 亚洲精品一卡2卡三卡4卡5卡| 天美传媒精品一区二区| 黄色日韩在线| 亚洲欧美日韩高清专用| 又粗又爽又猛毛片免费看| 国产成人欧美在线观看| 国产精品久久久人人做人人爽| 精品久久久久久久久久免费视频|