• <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ù)的選擇
    中文字幕色久视频| 亚洲第一青青草原| 免费看美女性在线毛片视频| 很黄的视频免费| 脱女人内裤的视频| 国产精品免费视频内射| 男人舔女人下体高潮全视频| 一个人观看的视频www高清免费观看 | 久久久久久国产a免费观看| 69av精品久久久久久| 极品人妻少妇av视频| 国产aⅴ精品一区二区三区波| 超碰成人久久| 久久国产精品男人的天堂亚洲| av在线天堂中文字幕| 国产成人av教育| 啦啦啦观看免费观看视频高清 | 国产成人av教育| 久久青草综合色| 亚洲人成网站在线播放欧美日韩| 精品国产美女av久久久久小说| 夜夜爽天天搞| 国产精品一区二区在线不卡| 精品一区二区三区视频在线观看免费| 欧美老熟妇乱子伦牲交| 制服诱惑二区| 可以在线观看毛片的网站| 欧美日韩黄片免| 亚洲成人国产一区在线观看| 后天国语完整版免费观看| 久久午夜综合久久蜜桃| 亚洲精品国产区一区二| 久久久久久国产a免费观看| 中文字幕人妻熟女乱码| 亚洲av成人一区二区三| 99国产极品粉嫩在线观看| 成人免费观看视频高清| 亚洲免费av在线视频| 亚洲午夜理论影院| 亚洲欧美日韩无卡精品| 黄频高清免费视频| 欧美国产日韩亚洲一区| 久久久久久免费高清国产稀缺| 国产私拍福利视频在线观看| 在线观看免费日韩欧美大片| 国产激情欧美一区二区| av中文乱码字幕在线| 亚洲av电影在线进入| 午夜老司机福利片| 老司机深夜福利视频在线观看| 亚洲少妇的诱惑av| 久久久久久久精品吃奶| 亚洲国产中文字幕在线视频| 久久香蕉精品热| 丰满的人妻完整版| 一边摸一边抽搐一进一小说| 天堂影院成人在线观看| 国产精品综合久久久久久久免费 | 黄色丝袜av网址大全| 国产精品99久久99久久久不卡| 一级a爱片免费观看的视频| 国产伦人伦偷精品视频| 亚洲精品国产一区二区精华液| 一边摸一边抽搐一进一出视频| 国产真人三级小视频在线观看| 两个人视频免费观看高清| 自拍欧美九色日韩亚洲蝌蚪91| 国产精品亚洲美女久久久| 久久久国产欧美日韩av| 亚洲黑人精品在线| 久久人妻av系列| 亚洲精品在线美女| 免费高清在线观看日韩| 极品人妻少妇av视频| 女人爽到高潮嗷嗷叫在线视频| 波多野结衣高清无吗| 嫩草影院精品99| av片东京热男人的天堂| 777久久人妻少妇嫩草av网站| 国产主播在线观看一区二区| 精品欧美国产一区二区三| 国内精品久久久久精免费| 国产高清视频在线播放一区| 国产av精品麻豆| 久久久国产成人精品二区| 欧美日韩瑟瑟在线播放| 久久天堂一区二区三区四区| 制服诱惑二区| 一本大道久久a久久精品| 丁香欧美五月| 禁无遮挡网站| 亚洲午夜理论影院| 久久精品国产亚洲av香蕉五月| 最近最新免费中文字幕在线| 亚洲色图综合在线观看| 国产成人系列免费观看| 桃红色精品国产亚洲av| 黄色毛片三级朝国网站| 巨乳人妻的诱惑在线观看| 色播亚洲综合网| 多毛熟女@视频| av免费在线观看网站| 午夜免费观看网址| 久99久视频精品免费| 色综合婷婷激情| 亚洲国产精品999在线| 欧美成人性av电影在线观看| 久久香蕉激情| 欧美乱色亚洲激情| 黄色女人牲交| 国产成人av教育| 又黄又粗又硬又大视频| 亚洲第一青青草原| 国产成人一区二区三区免费视频网站| 亚洲一区中文字幕在线| 天堂动漫精品| 女人被狂操c到高潮| 国产精品乱码一区二三区的特点 | 99精品在免费线老司机午夜| 亚洲第一av免费看| 久久精品成人免费网站| 色播在线永久视频| 亚洲精品久久成人aⅴ小说| 高清毛片免费观看视频网站| 9191精品国产免费久久| 男女做爰动态图高潮gif福利片 | 亚洲熟女毛片儿| 午夜福利,免费看| 亚洲黑人精品在线| 精品国产一区二区久久| 一卡2卡三卡四卡精品乱码亚洲| av天堂久久9| 亚洲一区二区三区不卡视频| 欧美av亚洲av综合av国产av| 一区二区三区精品91| 国产精品野战在线观看| 婷婷六月久久综合丁香| 制服诱惑二区| 国产精品国产高清国产av| 精品国产美女av久久久久小说| 欧美色视频一区免费| 亚洲av片天天在线观看| 久久精品国产清高在天天线| 精品一区二区三区四区五区乱码| 一级a爱视频在线免费观看| 日本精品一区二区三区蜜桃| 女人被狂操c到高潮| 日韩精品中文字幕看吧| 最好的美女福利视频网| 男女下面进入的视频免费午夜 | 亚洲 国产 在线| 无遮挡黄片免费观看| 欧美日韩乱码在线| 国产成人欧美在线观看| 免费在线观看日本一区| 国产乱人伦免费视频| 伊人久久大香线蕉亚洲五| 欧美午夜高清在线| 黄片小视频在线播放| 亚洲专区国产一区二区| 国产成年人精品一区二区| 两性夫妻黄色片| 天堂影院成人在线观看| 一级a爱视频在线免费观看| 午夜福利欧美成人| av网站免费在线观看视频| 在线播放国产精品三级| 亚洲av五月六月丁香网| 欧美黑人精品巨大| 久久精品91无色码中文字幕| 国产精品永久免费网站| 一级a爱视频在线免费观看| 亚洲第一青青草原| 久久香蕉激情| 亚洲欧美日韩高清在线视频| 国产一区二区三区视频了| 免费在线观看完整版高清| 国产熟女午夜一区二区三区| 国产成人免费无遮挡视频| 色在线成人网| 丝袜美腿诱惑在线| 性少妇av在线| 欧美性长视频在线观看| 亚洲专区中文字幕在线| 黄色女人牲交| 欧美乱妇无乱码| 日韩一卡2卡3卡4卡2021年| e午夜精品久久久久久久| 中出人妻视频一区二区| 亚洲国产日韩欧美精品在线观看 | 最新美女视频免费是黄的| 国产视频一区二区在线看| 99国产精品免费福利视频| 两性夫妻黄色片| 50天的宝宝边吃奶边哭怎么回事| 久久国产精品影院| 黄色女人牲交| a在线观看视频网站| 国产高清videossex| 一二三四在线观看免费中文在| 欧美黑人欧美精品刺激| 久久人妻福利社区极品人妻图片| 久久草成人影院| 免费一级毛片在线播放高清视频 | 一级作爱视频免费观看| 亚洲第一电影网av| 91九色精品人成在线观看| 神马国产精品三级电影在线观看 | 美女国产高潮福利片在线看| 国产视频一区二区在线看| 人人妻人人澡欧美一区二区 | 国产精品野战在线观看| 黄色 视频免费看| 91精品三级在线观看| 99riav亚洲国产免费| 欧洲精品卡2卡3卡4卡5卡区| 国产免费男女视频| 日韩精品青青久久久久久| 亚洲中文av在线| 午夜精品久久久久久毛片777| av网站免费在线观看视频| 久久久久久国产a免费观看| 亚洲欧美激情综合另类| 757午夜福利合集在线观看| 欧美av亚洲av综合av国产av| 亚洲av熟女| 国产高清videossex| 亚洲一区二区三区不卡视频| 亚洲欧美日韩无卡精品| 欧美成人性av电影在线观看| 亚洲国产中文字幕在线视频| 欧美日韩精品网址| 18禁国产床啪视频网站| 久久久久九九精品影院| 亚洲国产精品久久男人天堂| 成人国产综合亚洲| 精品一品国产午夜福利视频| 午夜亚洲福利在线播放| 国产精华一区二区三区| 亚洲欧美激情综合另类| 欧美大码av| 高潮久久久久久久久久久不卡| 亚洲精品在线美女| 免费高清视频大片| 久久香蕉激情| 亚洲色图综合在线观看| 18禁观看日本| 亚洲中文av在线| 午夜精品久久久久久毛片777| 久久国产亚洲av麻豆专区| 亚洲伊人色综图| 久久人妻福利社区极品人妻图片| 成人精品一区二区免费| 久久中文看片网| 久久草成人影院| 午夜精品国产一区二区电影| 禁无遮挡网站| 成人国产综合亚洲| 狠狠狠狠99中文字幕| 一级毛片高清免费大全| 亚洲男人天堂网一区| 成人亚洲精品av一区二区| www.自偷自拍.com| 可以在线观看的亚洲视频| 国产精品免费视频内射| 亚洲国产精品成人综合色| 亚洲国产日韩欧美精品在线观看 | 高清黄色对白视频在线免费看| 一区二区三区国产精品乱码| 19禁男女啪啪无遮挡网站| 丰满人妻熟妇乱又伦精品不卡| 亚洲自拍偷在线| 在线观看免费午夜福利视频| 亚洲一区二区三区色噜噜| 国产精品 国内视频| 亚洲国产欧美网| 岛国视频午夜一区免费看| 国产亚洲欧美精品永久| 国产精品影院久久| 99久久国产精品久久久| 18禁观看日本| 日韩有码中文字幕| 女人被躁到高潮嗷嗷叫费观| 给我免费播放毛片高清在线观看| 一级a爱视频在线免费观看| 男女之事视频高清在线观看| 国产精品亚洲美女久久久| 国产亚洲欧美在线一区二区| 欧美最黄视频在线播放免费| 国产欧美日韩精品亚洲av| 国产av一区在线观看免费| 国产精品久久久av美女十八| 在线av久久热| 亚洲自拍偷在线| 国产成人av教育| 久久婷婷成人综合色麻豆| 午夜福利在线观看吧| 中文字幕精品免费在线观看视频| 俄罗斯特黄特色一大片| 国产亚洲精品一区二区www| 国产精品日韩av在线免费观看 | 亚洲九九香蕉| 国产亚洲精品一区二区www| 国产三级黄色录像| 日日爽夜夜爽网站| 国产精品九九99| 国产精品一区二区在线不卡| 久久天堂一区二区三区四区| 国产人伦9x9x在线观看| 国产一卡二卡三卡精品| 91大片在线观看| 在线观看一区二区三区| 亚洲第一欧美日韩一区二区三区| 国产av在哪里看| 国产极品粉嫩免费观看在线| 午夜福利高清视频| 中文字幕人妻熟女乱码| 午夜福利高清视频| 国产欧美日韩综合在线一区二区| 天天躁夜夜躁狠狠躁躁| bbb黄色大片| 免费在线观看影片大全网站| 波多野结衣一区麻豆| 久久香蕉精品热| 亚洲精华国产精华精| 欧美日韩瑟瑟在线播放| 大型av网站在线播放| 超碰成人久久| 大型av网站在线播放| 搞女人的毛片| 最近最新中文字幕大全免费视频| 757午夜福利合集在线观看| 少妇裸体淫交视频免费看高清 | 国产区一区二久久| 国产av在哪里看| 非洲黑人性xxxx精品又粗又长| 老汉色∧v一级毛片| 国产色视频综合| 精品少妇一区二区三区视频日本电影| 国产精品 国内视频| 欧美不卡视频在线免费观看 | 黄色视频,在线免费观看| 久久久久久久久中文| 久久久久久国产a免费观看| 亚洲色图综合在线观看| 老汉色av国产亚洲站长工具| 丁香六月欧美| 久久久久久大精品| 国产亚洲精品av在线| 欧美日韩中文字幕国产精品一区二区三区 | 18禁美女被吸乳视频| 韩国av一区二区三区四区| 亚洲人成77777在线视频| 中文字幕高清在线视频| or卡值多少钱| 搞女人的毛片| 免费人成视频x8x8入口观看| 一级a爱视频在线免费观看| 一级片免费观看大全| 久久久久国产精品人妻aⅴ院| 男女之事视频高清在线观看| 老熟妇乱子伦视频在线观看| 久久精品成人免费网站| 国产亚洲精品久久久久久毛片| 久久久久久久久久久久大奶| 一边摸一边做爽爽视频免费| 久久午夜亚洲精品久久| 国产欧美日韩一区二区三区在线| 很黄的视频免费| 神马国产精品三级电影在线观看 | 欧美人与性动交α欧美精品济南到| 999精品在线视频| 麻豆国产av国片精品| 国产99白浆流出| 午夜精品在线福利| 少妇粗大呻吟视频| 很黄的视频免费| 国产男靠女视频免费网站| 国产亚洲欧美98| 国产三级在线视频| 99久久综合精品五月天人人| 亚洲精品久久成人aⅴ小说| www.精华液| 在线永久观看黄色视频| 国产aⅴ精品一区二区三区波| 国产精品免费一区二区三区在线| 精品国产乱子伦一区二区三区| 国产精品一区二区三区四区久久 | 不卡av一区二区三区| 亚洲七黄色美女视频| www日本在线高清视频| 亚洲色图综合在线观看| 午夜福利免费观看在线| 99国产精品99久久久久| 国产精品久久久人人做人人爽| 很黄的视频免费| 亚洲最大成人中文| 女性被躁到高潮视频| 琪琪午夜伦伦电影理论片6080| 国产99久久九九免费精品| 一进一出抽搐gif免费好疼| 18美女黄网站色大片免费观看| 在线国产一区二区在线| 97人妻精品一区二区三区麻豆 | 欧美精品亚洲一区二区| 长腿黑丝高跟| 一卡2卡三卡四卡精品乱码亚洲| 少妇粗大呻吟视频| 久久久精品欧美日韩精品| 欧美大码av| 久久久久亚洲av毛片大全| 一边摸一边抽搐一进一小说| 午夜日韩欧美国产| 成人亚洲精品一区在线观看| 9191精品国产免费久久| 91av网站免费观看| 在线十欧美十亚洲十日本专区| 香蕉丝袜av| 美国免费a级毛片| 欧美日韩亚洲综合一区二区三区_| 韩国av一区二区三区四区| 日韩视频一区二区在线观看| 国产亚洲精品一区二区www| 国内毛片毛片毛片毛片毛片| 欧美乱色亚洲激情| 亚洲一码二码三码区别大吗| 手机成人av网站| 制服丝袜大香蕉在线| 欧美大码av| 黑人巨大精品欧美一区二区蜜桃| 国产精品久久视频播放| 精品卡一卡二卡四卡免费| 国产精品久久久人人做人人爽| 久久久精品国产亚洲av高清涩受| 国产亚洲精品av在线| 精品久久久久久,| 禁无遮挡网站| 麻豆久久精品国产亚洲av| 非洲黑人性xxxx精品又粗又长| 天天躁狠狠躁夜夜躁狠狠躁| 国产片内射在线| 十八禁网站免费在线| 国产成人免费无遮挡视频| 免费看美女性在线毛片视频| 国产精品1区2区在线观看.| 久久性视频一级片| 韩国精品一区二区三区| 亚洲欧美精品综合久久99| 国产免费av片在线观看野外av| 91老司机精品| 国产精品自产拍在线观看55亚洲| 女性被躁到高潮视频| 亚洲精品av麻豆狂野| 国产欧美日韩精品亚洲av| 亚洲av片天天在线观看| 啪啪无遮挡十八禁网站| 亚洲一区二区三区色噜噜| x7x7x7水蜜桃| 一区福利在线观看| 看免费av毛片| 少妇裸体淫交视频免费看高清 | 色尼玛亚洲综合影院| 天堂动漫精品| 国产伦一二天堂av在线观看| 亚洲精品久久国产高清桃花| 欧美一区二区精品小视频在线| 亚洲一区高清亚洲精品| 国产欧美日韩一区二区三| xxx96com| 精品乱码久久久久久99久播| 首页视频小说图片口味搜索| 久久人人爽av亚洲精品天堂| 亚洲精品一区av在线观看| 亚洲国产精品成人综合色| 欧美绝顶高潮抽搐喷水| 国产精品九九99| 亚洲欧美日韩高清在线视频| 亚洲成人国产一区在线观看| 丰满人妻熟妇乱又伦精品不卡| 久久亚洲精品不卡| 自拍欧美九色日韩亚洲蝌蚪91| av中文乱码字幕在线| 亚洲人成伊人成综合网2020| 亚洲 国产 在线| 久久久精品国产亚洲av高清涩受| 午夜视频精品福利| 成年人黄色毛片网站| 精品国产亚洲在线| 91国产中文字幕| 国产午夜福利久久久久久| 国产亚洲av嫩草精品影院| 亚洲成人免费电影在线观看| 十八禁人妻一区二区| 黑人巨大精品欧美一区二区mp4| 1024视频免费在线观看| 国产精品免费视频内射| 一级片免费观看大全| 国产精品99久久99久久久不卡| 老鸭窝网址在线观看| 黄色 视频免费看| 亚洲成国产人片在线观看| 18禁国产床啪视频网站| 咕卡用的链子| 丝袜人妻中文字幕| av中文乱码字幕在线| 成人18禁在线播放| 久久精品人人爽人人爽视色| 精品福利观看| 无人区码免费观看不卡| 欧美黑人精品巨大| 男男h啪啪无遮挡| 国产成人精品无人区| 久久婷婷成人综合色麻豆| 狂野欧美激情性xxxx| 人人妻人人爽人人添夜夜欢视频| 国产成+人综合+亚洲专区| 男人的好看免费观看在线视频 | 亚洲国产精品合色在线| 欧美色视频一区免费| 亚洲国产毛片av蜜桃av| xxx96com| 精品免费久久久久久久清纯| 国产成人欧美| 国产亚洲欧美精品永久| 美女高潮到喷水免费观看| 手机成人av网站| 在线观看免费午夜福利视频| 国产主播在线观看一区二区| 一区二区三区国产精品乱码| 国产精品久久视频播放| 国产av在哪里看| 最新美女视频免费是黄的| 久久中文字幕一级| 国产99白浆流出| 精品久久蜜臀av无| 亚洲熟妇中文字幕五十中出| 一本综合久久免费| 亚洲色图 男人天堂 中文字幕| 久久国产乱子伦精品免费另类| 欧美久久黑人一区二区| 国产亚洲精品久久久久5区| 国产午夜福利久久久久久| 中文字幕色久视频| 99在线人妻在线中文字幕| 美女国产高潮福利片在线看| 国产麻豆成人av免费视频| 婷婷六月久久综合丁香| 黑人操中国人逼视频| 桃红色精品国产亚洲av| 亚洲国产精品999在线| 黄片播放在线免费| а√天堂www在线а√下载| 色哟哟哟哟哟哟| 国产亚洲欧美精品永久| 欧美日韩中文字幕国产精品一区二区三区 | 村上凉子中文字幕在线| 精品第一国产精品| а√天堂www在线а√下载| 少妇粗大呻吟视频| 日本在线视频免费播放| 免费高清在线观看日韩| 夜夜看夜夜爽夜夜摸| 最新在线观看一区二区三区| 日本精品一区二区三区蜜桃| 亚洲一区中文字幕在线| 国产成人一区二区三区免费视频网站| 久久精品国产清高在天天线| 熟妇人妻久久中文字幕3abv| 久久久国产成人免费| 久久精品aⅴ一区二区三区四区| 亚洲av成人一区二区三| 给我免费播放毛片高清在线观看| 国产真人三级小视频在线观看| 国产1区2区3区精品| 一边摸一边抽搐一进一小说| 精品久久久久久,| 国产精品电影一区二区三区| 久久天躁狠狠躁夜夜2o2o| 国产精品久久久久久人妻精品电影| 女人被躁到高潮嗷嗷叫费观| 成人18禁在线播放| 日韩欧美免费精品| 免费在线观看影片大全网站| 欧美国产日韩亚洲一区| 国产av精品麻豆| 亚洲一区二区三区色噜噜| 精品久久久久久成人av| 精品久久久久久久毛片微露脸| 亚洲精品在线美女| 国产激情久久老熟女| 成在线人永久免费视频| 国产成人欧美| 亚洲第一电影网av| 性色av乱码一区二区三区2| 18禁观看日本| 久久国产精品人妻蜜桃| 亚洲av熟女| 人人妻,人人澡人人爽秒播| 国产成人av激情在线播放| av网站免费在线观看视频| 涩涩av久久男人的天堂| 99精品在免费线老司机午夜| 日韩av在线大香蕉| 亚洲 国产 在线| 久久人妻av系列| 日韩av在线大香蕉| 午夜影院日韩av| 黄色毛片三级朝国网站| 亚洲五月婷婷丁香| av在线播放免费不卡| 91九色精品人成在线观看| 久久香蕉激情| 宅男免费午夜| 亚洲一区中文字幕在线| 桃色一区二区三区在线观看|