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

    基于貝葉斯模型平均法和逐步回歸法構(gòu)建杉木單木胸徑生長(zhǎng)模型*

    2021-12-01 01:50:50魯樂樂張雄清張建國(guó)
    林業(yè)科學(xué) 2021年9期
    關(guān)鍵詞:單木后驗(yàn)杉木

    魯樂樂 王 震 張雄清 張建國(guó)

    (1.中國(guó)林業(yè)科學(xué)研究院林業(yè)研究所 國(guó)家林業(yè)和草原局林木培育重點(diǎn)實(shí)驗(yàn)室 北京100091;2.南京林業(yè)大學(xué)南方現(xiàn)代林業(yè)協(xié)同創(chuàng)新中心 南京210037)

    單木胸徑是估算立木材積、評(píng)價(jià)森林生產(chǎn)力的重要基礎(chǔ)因子(Lhotkaetal., 2011)。影響單木胸徑生長(zhǎng)變化的因子很多,如氣候、立地、林齡、林分密度、林分結(jié)構(gòu)等(Clarketal., 2014),總體上可以分為內(nèi)部因子(主要指林分因子)和外部因子(主要指氣候因子)兩大類(Zhangetal., 2017)。

    氣候因子是單木胸徑生長(zhǎng)變化的重要驅(qū)動(dòng)因子,已引起學(xué)者廣泛關(guān)注(Allenetal., 2015; Hanewinkeletal., 2012)。Henderson 等(2009)建立長(zhǎng)葉松(Pinuspalustris)胸徑生長(zhǎng)與氣候因子之間的關(guān)系模型,結(jié)果表明在所有氣候變量中,溫度是最弱的驅(qū)動(dòng)因子。黎敬業(yè)等(2019)分析馬尾松(Pinusmassoniana)徑向生長(zhǎng)與氣候因子的關(guān)系,結(jié)果發(fā)現(xiàn)高海拔地區(qū)馬尾松徑向生長(zhǎng)對(duì)氣候因子年際波動(dòng)敏感性較強(qiáng),與生長(zhǎng)季前冬季光溫條件和生長(zhǎng)季內(nèi)7月降雨呈正相關(guān),生長(zhǎng)-氣候關(guān)系在不同樣點(diǎn)間表現(xiàn)出較強(qiáng)一致性。王延芳等(2020)基于 Vaganov-Shashkin 模型模擬青海云杉(Piceacrassifolia)徑向生長(zhǎng),結(jié)果表明降雨對(duì)祁連山中部低海拔地區(qū)青海云杉徑向生長(zhǎng)具有決定性作用。

    對(duì)于林分因子,包括林分密度(Zhangetal., 2015)、立地(Aragaoetal., 2009)、年齡(Fosteretal., 2014)等對(duì)單木胸徑生長(zhǎng)變化的驅(qū)動(dòng)性也不小。通常來說密度越大,單木胸徑生長(zhǎng)量越小,立地條件越好,單木胸徑生長(zhǎng)量越大(DeRoseetal., 2009; Zhangetal., 2015)。王冬至等(2017)分析華北落葉松(Larixprincipis-rupprechtii)-白樺(Betulaplatyphylla)針闊混交林胸徑生長(zhǎng),結(jié)果表明年齡和競(jìng)爭(zhēng)指數(shù)與胸徑生長(zhǎng)量呈負(fù)相關(guān)關(guān)系,胸徑和林分優(yōu)勢(shì)高與胸徑生長(zhǎng)量呈正相關(guān)關(guān)系。張海平(2017)構(gòu)建基于氣象因子的天然白樺林單木胸徑生長(zhǎng)模型,結(jié)果發(fā)現(xiàn)溫度和降雨有利于白樺生長(zhǎng)。

    綜上可知,氣候因子和林分因子對(duì)單木胸徑生長(zhǎng)均有一定驅(qū)動(dòng)性。在進(jìn)行模型分析選擇變量時(shí),一般大部分研究都是主觀指定包含某些變量的某個(gè)特定模型為最優(yōu)模型,由此導(dǎo)致所建立模型的預(yù)測(cè)能力減弱,這是頻率統(tǒng)計(jì)模型的一個(gè)弱點(diǎn)。如逐步回歸法(stepwise regression, SR)根據(jù)變量顯著性(P),選擇顯著的變量保留在模型中; 而實(shí)際上事先并不知道包含哪些變量的模型為最優(yōu)模型,即模型本身存在不確定性(Rafteryetal., 1991; Draper, 1995)。當(dāng)基于研究數(shù)據(jù)建立回歸模型時(shí),若忽略模型本身的不確性僅以某一特定模型進(jìn)行分析,一方面會(huì)低估模型的不確定性,高估預(yù)測(cè)結(jié)果,限制模型的適用范圍; 另一方面會(huì)比正常情況下更趨向于拒絕無效假設(shè),從而產(chǎn)生誤導(dǎo)性結(jié)果(Viallefontetal., 2001; 張志杰等, 2007)。

    貝葉斯模型平均法(Bayesian model averaging,BMA)是近些年發(fā)展較快的一種多因子和模型不確定性分析方法(Picardetal., 2012)。與傳統(tǒng)分析方法不同,BMA并非僅指定一個(gè)最優(yōu)模型,而是考慮模型空間內(nèi)所有可能的模型以及各模型的后驗(yàn)概率,將模型的后驗(yàn)概率加權(quán)平均得到貝葉斯平均模型,并根據(jù)模型中變量估計(jì)值不等于 0 的后驗(yàn)概率進(jìn)行變量重要性排序,判斷各模型出現(xiàn)的可能性(羅劍鋒, 2003)。近年來,BMA 在林業(yè)領(lǐng)域也有些研究,如Bullock等(2007)成功地將 BMA 應(yīng)用于林木胸徑分布中、Picard等(2012)利用 BMA分析了熱帶森林地面生物量、Zhang 等(2014)利用 BMA 構(gòu)建了杉木(Cunninghamialanceolata)林分?jǐn)嗝娣e生長(zhǎng)模型、Lu等(2019)利用BMA構(gòu)建了枯損率模型等; 但是,在錯(cuò)綜復(fù)雜的影響因子中如何更好地選擇因子,并構(gòu)建不確定性胸徑生長(zhǎng)模型還未見報(bào)道。

    杉木是我國(guó)速生高產(chǎn)用材樹種之一,也是我國(guó)南方亞熱帶地區(qū)最重要的鄉(xiāng)土栽培樹種。據(jù)第九次全國(guó)森林資源清查結(jié)果,杉木人工林總造林面積達(dá) 990.2萬hm2,占全國(guó)十大主要造林樹種面積的27.23%,蓄積占32.57%,均居各大造林樹種之首。本研究擬探索杉木人工林單木胸徑生長(zhǎng)量變化的驅(qū)動(dòng)因子,比較不同驅(qū)動(dòng)因子的重要性,構(gòu)建不確定性單木胸徑生長(zhǎng)模型,以期為杉木經(jīng)營(yíng)管理者科學(xué)經(jīng)營(yíng)管理杉木人工林提供參考。

    1 研究區(qū)概況與數(shù)據(jù)調(diào)查

    研究區(qū)位于福建省邵武市衛(wèi)閩林場(chǎng)(117°43′E,27°05′N),為杉木中心產(chǎn)區(qū),地貌主要為低山、高丘,海拔250~700 m,坡度25°~35°。屬亞熱帶季風(fēng)氣候,年均氣溫17.7 ℃,1月平均氣溫6.8 ℃,7月平均氣溫28 ℃,極端低溫-7.9 ℃,年均降雨量1 768 mm,年均相對(duì)濕度82%。

    研究區(qū)杉木密度試驗(yàn)林為1982年使用1年生苗木造林,采用隨機(jī)區(qū)組試驗(yàn)設(shè)計(jì),設(shè)置5種密度: A,2 m×3 m(1 667 株·hm-2); B, 2 m×1.5 m(3 333 株·hm-2); C,2 m×1 m(5 000 株·hm-2); D,1 m×1.5 m(6 667 株·hm-2);E, 1 m×1 m(10 000 株·hm-2)。每種密度重復(fù)3次,共設(shè)立樣地15塊,每塊樣地大小為20 m×30 m。造林后1984—1990年逐年調(diào)查,之后到2010年隔年調(diào)查,每次調(diào)查均在秋冬季進(jìn)行。測(cè)量樹高時(shí),每塊樣地隨機(jī)抽取50株樹,計(jì)算最高5株樹的平均高作為林分優(yōu)勢(shì)木平均高。

    2 研究方法

    2.1 變量選擇

    內(nèi)部因子主要指單木大小因子、競(jìng)爭(zhēng)因子和林分因子。單木大小因子包括胸徑(diameter at breast height, DBH,cm)、年齡(age,A,a); 競(jìng)爭(zhēng)因子包括胸高斷面積(basal area, BA,m2·hm-2)、大于對(duì)象木的斷面積和(sum of basal areas of trees larger than the subject tree, BAL,m2·hm-2); 林分因子包括林分密度(number,N,trees·hm-2)、優(yōu)勢(shì)木平均高(dominant height,Hd,m)、林分平方平均胸徑(quadratic mean diameter,Dg,cm)。內(nèi)部因子統(tǒng)計(jì)值見表1。

    外部因子主要指氣候因子,通過ClimateAP軟件獲得(Wangetal., 2012),包括年均氣溫(mean annual temperature,MAT,℃)、年均降雨量(mean annual precipitation,AP,mm)、最冷月平均溫度(mean coldest month temperature,MCMT,℃)、年均干旱指數(shù)(annual heat-moisture index,AHM)、低于0 ℃天數(shù)(degree-days below 0 ℃, DD0, d)、夏季平均最高溫度(summer mean maximum temperature,SMMT,℃)、冬季平均最低溫度(winter mean minimum temperature,WMMT,℃)、最熱月平均溫度(mean warmest month temperature,MWMT,℃)和春季平均氣溫(spring mean temperature,SMT,℃)。氣候因子統(tǒng)計(jì)值見表2。

    表2 氣候因子統(tǒng)計(jì)值①Tab.2 Summary statistics of climate factors

    2.2 模型構(gòu)建

    根據(jù)樹木生長(zhǎng)規(guī)律,胸徑生長(zhǎng)量(annual increment of DBH, AID)大于或等于0; 但是當(dāng)AID等于0時(shí),其對(duì)數(shù)變換lnAID無法計(jì)算。為避免這些誤差,Calama等(2005)和李春明(2012)在對(duì)數(shù)變換之前給胸徑生長(zhǎng)增加常數(shù)1,以胸徑生長(zhǎng)量加1的對(duì)數(shù)值ln(AID +1)為模型因變量。本研究構(gòu)建的胸徑生長(zhǎng)量模型形式如下:

    ln(AID+1)=x0+x1lnDBH+x2DBH+

    x3BA+x4Dg+x5Hd+x6lnA+x7N+

    x8BAL+x9MCMT+x10AP+x11AHM+

    x12SMMT+x13WMWT+x14SMT+x15MAT+

    x16DD0+x17MWMT。

    (1)

    其中: AID為 1984—2010 年每2年間隔的年均胸徑生長(zhǎng)量(為排除不同間隔期帶來的影響,僅研究每2年間隔的年均胸徑生長(zhǎng)量);x0~x17為參數(shù)。

    數(shù)據(jù)的對(duì)數(shù)變換會(huì)導(dǎo)致胸徑生長(zhǎng)量估計(jì)偏差,Baskerville(1972)給出了對(duì)這種偏差進(jìn)行修正估計(jì)的近似值:

    (2)

    2.3 BMA原理

    假設(shè)y為研究感興趣的量(胸徑生長(zhǎng)量),Z為調(diào)查所得數(shù)據(jù)(氣候因子和內(nèi)部因子),f={f1,…,fk}(k=1,…,n)代表所有可能模型組成的模型空間,而包含哪些變量的模型為最優(yōu)模型事先并不知道,即模型本身存在不確定性。通常情況下,可以考慮的模型數(shù)量很大,如有p個(gè)自變量,那么在僅考慮主效應(yīng)時(shí)其模型數(shù)量為2p。BMA選取可能模型的子集,并用這些模型的后驗(yàn)概率進(jìn)行所有推理和預(yù)測(cè)(Kassetal., 1995)。根據(jù)貝葉斯模型平均法,y的后驗(yàn)分布為:

    (3)

    式中:P(y|fk,Z)為給定調(diào)查數(shù)據(jù)Z和模型fk條件下y的后驗(yàn)分布;P(fk|Z)為給定調(diào)查數(shù)據(jù)Z條件下fk為最優(yōu)模型的概率。

    由方程(3)可知,y的后驗(yàn)分布實(shí)際上是以后驗(yàn)概率P(fk|Z)為權(quán)重對(duì)所有模型后驗(yàn)分布進(jìn)行加權(quán)的一個(gè)平均值。在方程(3)中,假定fk為最優(yōu)模型,則概率分布為:

    (4)

    式中:β=(β0,β1,…,βn)為模型fk回歸系數(shù)的向量。

    對(duì)于變量的先驗(yàn)分布,本研究利用無信息先驗(yàn)分布(Boxetal., 1973)。

    模型空間中所有模型數(shù)量很大,假設(shè)自變量為5個(gè),那么即使不考慮變量間的交互作用,可構(gòu)建的模型也可達(dá)32個(gè),然而在實(shí)際研究中只有滿足要求的模型才能被選進(jìn)來。Madigan等(1994)建議采用Occam窗法排除不當(dāng)?shù)哪P鸵詼p少模型數(shù)量: 當(dāng)一個(gè)模型后驗(yàn)概率小于最優(yōu)模型后驗(yàn)概率的5%時(shí),則從模型空間中移除該模型。本研究采用Occam窗法排除不當(dāng)?shù)哪P?,排序得到后?yàn)概率大的模型。

    根據(jù)貝葉斯理論,備擇假設(shè)H1為:Xk為杉木單木胸徑生長(zhǎng)量變化驅(qū)動(dòng)因子的后驗(yàn)概率有多大,即βk不等于0的后驗(yàn)概率有多大。BMA以所有包含變量Xk的杉木單木胸徑生長(zhǎng)模型后驗(yàn)概率的和作為βk不等于0的后驗(yàn)概率的估計(jì):

    (5)

    式中:Ii為0/1指示變量,當(dāng)βk在模型fi中,Ii=1,反之Ii=0;Ω為模型空間。

    Kass等(1995)認(rèn)為,BMA變量后驗(yàn)概率推斷的規(guī)則如下:P(βk≠0|Z)<0.5,表示沒有證據(jù)表明Xk是杉木單木胸徑生長(zhǎng)量變化的驅(qū)動(dòng)因子; 0.5 ≤P(βk≠0|Z)<0.75,表示有弱的證據(jù)表明Xk是杉木單木胸徑生長(zhǎng)量變化的驅(qū)動(dòng)因子; 0.75 ≤P(βk≠0|Z)<0.95,表示有強(qiáng)的證據(jù)表明Xk是杉木單木胸徑生長(zhǎng)量變化的驅(qū)動(dòng)因子;P(βk≠0|Z)≥0.95,表示有很強(qiáng)的證據(jù)表明Xk是杉木單木胸徑生長(zhǎng)量變化的驅(qū)動(dòng)因子。

    若Xk是杉木單木胸徑生長(zhǎng)量變化的驅(qū)動(dòng)因子,則根據(jù)BMA模型似然函數(shù)推導(dǎo)得出杉木單木胸徑生長(zhǎng)模型中變量因子Xk的點(diǎn)估計(jì)的后驗(yàn)概率均值為:

    (6)

    之后對(duì)該變量后驗(yàn)概率均值進(jìn)行排序,選擇杉木單木胸徑生長(zhǎng)量變化的主要驅(qū)動(dòng)因子,并構(gòu)建相應(yīng)的不確定性胸徑生長(zhǎng)模型。

    2.4 模型評(píng)價(jià)

    采用決定系數(shù)(R2)和平均絕對(duì)差(mean absolute difference, MAD)作為模型預(yù)測(cè)的評(píng)價(jià)指標(biāo)。R2越大、MAD越小,模型精度越高:

    (7)

    (8)

    為探究不同初植密度單木胸徑生長(zhǎng)響應(yīng)差異,本研究分不同密度以及所有樣本進(jìn)行分析處理。在不同處理中,隨機(jī)抽取總數(shù)據(jù)的60%作為建模數(shù)據(jù),剩余40%作為檢驗(yàn)數(shù)據(jù)。同時(shí),分別利用BMA和SR(雙向選擇)2種方法分析杉木單木胸徑生長(zhǎng)量與內(nèi)部和外部因子的關(guān)系。BMA分析通過R軟件BMA包實(shí)現(xiàn)(Rafteryetal., 2005),SR分析通過step函數(shù)實(shí)現(xiàn)。

    3 結(jié)果與分析

    3.1 BMA和SR確定的最佳模型比較

    在E密度中,SR確定的模型與BMA確定的最優(yōu)模型相似。在A和D密度中,SR確定的模型與BMA確定的第二優(yōu)模型相似。在B、C密度和全樣本數(shù)據(jù)中,SR確定的模型不在BMA確定的后驗(yàn)概率較大的前幾個(gè)模型中(表3)。B、C密度和全樣本數(shù)據(jù)3種處理在模型選擇上表現(xiàn)出明顯差異。

    表3 BMA和SR確定的單木胸徑生長(zhǎng)模型和BMA模型的后驗(yàn)概率①Tab.3 Individual tree diameter growth model by BMA and SR methods and posterior probability of BMA model

    除此之外,2種方法確定的模型R2和MAD在5種密度水平和所有密度水平的集合數(shù)據(jù)集中幾乎相同,模型預(yù)測(cè)效果沒有明顯差異(表4)。

    表4 基于BMA和SR的單木胸徑生長(zhǎng)模型評(píng)價(jià)Tab.4 Evaluation of individual tree diameter growth model based on BMA and SR methods

    3.2 BMA和SR的參數(shù)估計(jì)值和P比較

    A密度:BMA模型中,單木胸徑生長(zhǎng)量隨DBH、MWMT、MCMT、AP和AHM增加而增加,隨N、Dg、lnDBH、lnA、SMMT和WMMT增加而減小;SR模型中,單木胸徑生長(zhǎng)量隨BA、DBH、MWMT、MCMT、AP、AHM、DD0和SMT增加而增加,隨N、lnDBH、lnA、MAT和WMMT增加而減小(表5、圖1)。

    B密度:BMA模型中,單木胸徑生長(zhǎng)量與Hd、DBH、MWMT、MCMT和AP呈正相關(guān)關(guān)系,與N、Dg、lnDBH、lnA、BAL、SMMT和WMMT呈負(fù)相關(guān)關(guān)系;SR模型中,單木胸徑生長(zhǎng)量與DD0和SMT呈正相關(guān)關(guān)系(表5、圖1)。

    C密度:BMA模型中,單木胸徑生長(zhǎng)量與BA、DBH、lnDBH、MWMT、MCMT、AP和DD0呈正相關(guān)關(guān)系,與N、Dg、lnA、BAL、SMMT和WMMT呈負(fù)相關(guān)關(guān)系;SR模型中,單木胸徑生長(zhǎng)量與Hd、MAT呈負(fù)相關(guān)關(guān)系,與AHM和SMMT呈正相關(guān)關(guān)系(表5、圖1)。

    D密度:BMA模型中,單木胸徑生長(zhǎng)量與BA、Hd、DBH和MAT顯著正相關(guān),與N、Dg、lnA、BAL、AHM、DD0、WMMT和SMT顯著負(fù)相關(guān);SR模型中,單木胸徑生長(zhǎng)量與AP和SMMT顯著負(fù)相關(guān),與MWMT顯著正相關(guān)(表5、圖1)。

    E密度:BMA模型中,單木胸徑生長(zhǎng)量與BA、Dg、DBH、MWMT、MCMT、DD0和SMT呈正相關(guān)關(guān)系,與lnDBH、lnA、BAL、MAT、SMMT和WMMT呈負(fù)相關(guān)關(guān)系;SR模型中,單木胸徑生長(zhǎng)量也與AHM呈負(fù)相關(guān)關(guān)系(表5、圖1)。

    全樣本數(shù)據(jù)處理中,BMA選擇的變量與D密度相似,較D密度多選擇了MWMT、MCMT和AP,且顯著正相關(guān),少選擇了AHM;逐步回歸模型中,單木胸徑生長(zhǎng)量與lnDBH顯著負(fù)相關(guān),與AHM顯著正相關(guān)(表5、圖1)。

    在多數(shù)處理中,BMA和SR的年均胸徑生長(zhǎng)量隨DBH、BA、MWMT、MCMT、AP增加而增加,隨N、lnA、BAL、WMMT增加而減小,且隨期初胸徑增大,年均胸徑生長(zhǎng)量逐漸變小(圖2)。Dg在A、B、C、D密度和全樣本數(shù)據(jù)中與胸徑生長(zhǎng)量表現(xiàn)為負(fù)相關(guān)關(guān)系,但是在E密度中表現(xiàn)為正相關(guān)關(guān)系。與DBH不同,lnDBH只在C密度中與胸徑生長(zhǎng)量表現(xiàn)為正效應(yīng)關(guān)系,在A、B、E密度中表現(xiàn)為負(fù)效應(yīng)關(guān)系,在另外2種處理中與胸徑生長(zhǎng)量的關(guān)系不顯著。MAT在中低密度(A、B、C密度)中對(duì)胸徑生長(zhǎng)量的影響不大,在D密度中表現(xiàn)為正效應(yīng)關(guān)系,在E密度和全樣本數(shù)據(jù)中表現(xiàn)為負(fù)效應(yīng)關(guān)系。對(duì)于AHM、SMMT和SMT,BMA只在較少處理中認(rèn)為與胸徑生長(zhǎng)量顯著相關(guān),但是SR在大多數(shù)處理中都選擇了這些變量(表5)。

    BMA認(rèn)為,Hd在B、D密度和全樣本數(shù)據(jù)中與胸徑生長(zhǎng)量呈正相關(guān)關(guān)系,與SR略有不同,SR認(rèn)為在C密度中胸徑生長(zhǎng)量隨Hd增加而減小。在A、B、C密度和全樣本數(shù)據(jù)中,BMA認(rèn)為胸徑生長(zhǎng)量隨AP增加而增加;而在D密度中,SR認(rèn)為AP與胸徑生長(zhǎng)量呈負(fù)效應(yīng)關(guān)系。

    此外,SR選擇的變量比BMA更多。2種方法選擇的內(nèi)部因子幾乎相同,BMA只在A密度中多選擇了Dg,SR在A密度中多選擇了BA,在全樣本數(shù)據(jù)中多選擇了lnDBH,但是在每種處理中,SR都比BMA多選擇了氣候因子。由此可見,SR更容易選擇冗余變量到模型中(表5)。

    4 討論

    4.1 BMA與SR的比較

    本研究基于BMA與SR構(gòu)建單木胸徑生長(zhǎng)模型,SR根據(jù)自變量對(duì)因變量的意義選擇變量,但是容易忽略模型本身的不確定性; BMA的優(yōu)點(diǎn)已在幾種不同類型模型中進(jìn)行了評(píng)估,如線性回歸、生存模型等(Raftery, 1996; Volinskyetal., 1997; Murphyetal., 2001),均表明BMA可提高模型性能。本研究結(jié)果同樣表明,BMA在變量選擇上優(yōu)于SR。

    采用SR擬合的模型在3種處理下位于BMA選擇的后驗(yàn)概率高的前幾個(gè)模型中,僅在E密度中具有與BMA選擇最佳模型相同的預(yù)測(cè)變量(表3),這表明模型本身存在不確定性,SR忽略了這一點(diǎn)。

    表5 SR和BMA構(gòu)建單木胸徑生長(zhǎng)模型的參數(shù)估計(jì)值①Tab.5 The parameter estimates of individual tree diameter growth model by BMA and SR methods

    圖1 不同密度處理下BMA模型各自變量后驗(yàn)概率Fig.1 The posterior probability of variables in individual tree diameter growth model determined by BMA method

    圖2 年均胸徑生長(zhǎng)量與期初胸徑的關(guān)系Fig.2 The relationship between annual increment of DBH (AID)and diameter at the beginning of growth period

    單木胸徑生長(zhǎng)模型中,BMA認(rèn)為,Hd在B、D密度和全樣本數(shù)據(jù)中與單木胸徑生長(zhǎng)量呈正相關(guān)關(guān)系,而SR認(rèn)為,在C密度中單木胸徑生長(zhǎng)量隨Hd增加而減小。一般情況下,胸徑生長(zhǎng)量應(yīng)隨優(yōu)勢(shì)高增加而增加。BMA認(rèn)為,在A、B、C密度和全樣本數(shù)據(jù)中單木胸徑生長(zhǎng)量隨年均降雨量增加而增加,而在D密度中,SR認(rèn)為年均降雨量與胸徑生長(zhǎng)量呈負(fù)效應(yīng)關(guān)系。有研究發(fā)現(xiàn),降雨量增加可以促進(jìn)樹木生長(zhǎng)(余黎等, 2014; 張海平, 2017)。可見,BMA的結(jié)論比SR具有更大的可解釋性。

    此外,SR選擇的變量比BMA更多,即選擇冗余變量較多。單木胸徑生長(zhǎng)模型中SR選擇了更多的氣候因子,余黎等(2014)、歐強(qiáng)新等(2019)發(fā)現(xiàn),單木胸徑年均生長(zhǎng)量受氣候因子影響較小,而競(jìng)爭(zhēng)因子和單木大小因子是影響單木胸徑年均生長(zhǎng)量的主要因子,這與BMA的結(jié)果是一致的。Genell等(2010)利用模擬數(shù)據(jù)表明,BMA不選擇冗余變量的概率高于SR,且選擇真實(shí)預(yù)測(cè)因子的概率與SR相似。

    需要說明的是,2種方法在模型預(yù)測(cè)性能上的差異較小。Zapata-Cuartas等(2012)研究表明,對(duì)于小樣本,貝葉斯方法的性能優(yōu)于經(jīng)典方法。Wang等(2019)將分層貝葉斯方法與非線性混合效應(yīng)模型進(jìn)行比較,并估計(jì)地上生物量,結(jié)果發(fā)現(xiàn)2種方法在大樣本時(shí)表現(xiàn)相似,而對(duì)于小樣本貝葉斯方法要好得多。

    4.2 單木胸徑生長(zhǎng)量的影響因子

    內(nèi)部因子中,年均胸徑生長(zhǎng)量隨單木期初胸徑增加而增加,隨lnA增加而減小,這表明單木期初胸徑越大,單木生長(zhǎng)優(yōu)勢(shì)越大,且在幼齡時(shí)增長(zhǎng)迅速,隨著年齡增加增長(zhǎng)速度變緩,即隨著期初胸徑增大,年均胸徑生長(zhǎng)量變小(圖2)。一些研究表明,期初胸徑達(dá)到一定大小后,即達(dá)到大徑階后,其增長(zhǎng)逐漸緩慢(張海平, 2017)。但lnDBH只在C密度中對(duì)年均胸徑生長(zhǎng)量有正效應(yīng)影響,在3種處理(A、B、E密度)中與年均胸徑生長(zhǎng)量呈負(fù)相關(guān)關(guān)系。

    年均胸徑生長(zhǎng)量隨林分密度、林分平方平均胸徑和大于對(duì)象木的斷面積和增加而減小,與以往研究結(jié)果相同(Brooksetal., 1998; Huangetal., 2010)。大于對(duì)象木的斷面積之和越小,林分內(nèi)單木平均擁擠程度越??; 林分密度和林分平方平均胸徑越小,林分內(nèi)單木對(duì)生長(zhǎng)空間的平均占有和利用越大,故單木年均胸徑生長(zhǎng)量越大(余黎等, 2014)。這都表明競(jìng)爭(zhēng)越小,胸徑生長(zhǎng)量越大。但是在中高密度(C、D、E密度)和全樣本數(shù)據(jù)中,胸高斷面積對(duì)年均胸徑生長(zhǎng)量有正效應(yīng)影響。在B、E密度和全樣本數(shù)據(jù)中,年均胸徑生長(zhǎng)量隨優(yōu)勢(shì)木平均高增加而增加,在其他處理中優(yōu)勢(shì)木平均高對(duì)年均胸徑生長(zhǎng)量沒有影響。一般來說,樹木生長(zhǎng)與優(yōu)勢(shì)木平均高呈正相關(guān)關(guān)系,王冬至等(2017)研究華北落葉松-白樺混交林發(fā)現(xiàn),立地質(zhì)量對(duì)胸徑生長(zhǎng)量也具有一定影響,并與胸徑生長(zhǎng)量呈正相關(guān)關(guān)系,但影響較小。

    氣候因子中,最熱月平均溫度、最冷月平均溫度和年均降雨量均對(duì)年均胸徑生長(zhǎng)量有促進(jìn)作用。有研究發(fā)現(xiàn),年均降雨量對(duì)冷杉(Abiesnephrolepis)胸徑生長(zhǎng)量表現(xiàn)為正效應(yīng)關(guān)系,但達(dá)到一定閾值后表現(xiàn)為負(fù)效應(yīng)關(guān)系(余黎等, 2014),即當(dāng)水分成為限制因子時(shí),降雨往往對(duì)生長(zhǎng)表現(xiàn)出正效應(yīng),而當(dāng)水分充足或過多時(shí),逐漸表現(xiàn)出負(fù)效應(yīng)。也有研究發(fā)現(xiàn),生長(zhǎng)季水分多少可直接或間接影響細(xì)胞的分裂速率(程瑞梅等, 2015),而且還在一定程度上影響生長(zhǎng)季中樹木體內(nèi)營(yíng)養(yǎng)物質(zhì)的積累(姜倩倩等, 2012)。張海平(2017)研究發(fā)現(xiàn),生長(zhǎng)季最低溫對(duì)年均胸徑生長(zhǎng)量表現(xiàn)為正效應(yīng)關(guān)系。李文馨等(2015)發(fā)現(xiàn),年均積溫和年均降雨對(duì)胸徑生長(zhǎng)量有促進(jìn)作用,該結(jié)果與本研究結(jié)果相反。

    另外,在低密度(A密度)時(shí),干旱指數(shù)越大年均胸徑生長(zhǎng)量越大,在高密度(D密度)時(shí),干旱指數(shù)越大年均生長(zhǎng)量越小。余黎等(2014)發(fā)現(xiàn),年均胸徑生長(zhǎng)量隨年均氣溫與年均降雨量之比增大而下降,與本研究在D密度的結(jié)果是一致的,在A密度呈正相關(guān)的原因可能是低密度時(shí)降雨對(duì)每株樹都是充足的,所以干旱指數(shù)對(duì)年均胸徑生長(zhǎng)量沒有抑制作用。

    5 結(jié)論

    本研究采用逐步回歸法和貝葉斯模型平均法構(gòu)建杉木人工林單木胸徑生長(zhǎng)模型,逐步回歸法獲得模型的后驗(yàn)概率小于BMA獲得最佳模型(最高后驗(yàn)概率)或逐步回歸模型不在BMA模型空間前幾個(gè)后驗(yàn)概率高的模型中,由此表現(xiàn)出杉木單木胸徑生長(zhǎng)模型本身的不確定性。此外,逐步回歸法更容易選擇冗余變量。

    內(nèi)部因子和氣候因子均會(huì)影響杉木單木胸徑生長(zhǎng),但單木胸徑年均生長(zhǎng)量受氣候因子的影響較小,而競(jìng)爭(zhēng)因子和單木大小因子是影響單木胸徑年均生長(zhǎng)量的主要因子。單木胸徑生長(zhǎng)量隨林分密度、林分平方平均胸徑、大于對(duì)象木的斷面積和、年齡、冬季平均最低溫度增加而減小,隨期初胸徑、胸高斷面積、優(yōu)勢(shì)木平均高、最冷月平均溫度、最熱月平均溫度、年均降雨量增加而增加??偟膩碚f,杉木單木胸徑生長(zhǎng)量隨競(jìng)爭(zhēng)增加而減小,隨溫度和降雨增加而增加。

    猜你喜歡
    單木后驗(yàn)杉木
    地基與無人機(jī)激光雷達(dá)結(jié)合提取單木參數(shù)
    融合LiDAR點(diǎn)云與高分影像的單木檢測(cè)方法研究
    基于對(duì)偶理論的橢圓變分不等式的后驗(yàn)誤差分析(英)
    杉木黃化病的防治技術(shù)措施研究
    貝葉斯統(tǒng)計(jì)中單參數(shù)后驗(yàn)分布的精確計(jì)算方法
    無人機(jī)影像匹配點(diǎn)云單木識(shí)別算法
    遙感信息(2019年1期)2019-03-22 01:38:16
    基于雙尺度體元覆蓋密度的TLS點(diǎn)云數(shù)據(jù)單木識(shí)別算法
    森林工程(2018年5期)2018-05-14 13:54:30
    杉木萌芽更新關(guān)鍵技術(shù)
    杉木育苗化學(xué)防除雜草技術(shù)
    一種基于最大后驗(yàn)框架的聚類分析多基線干涉SAR高度重建算法
    欧美黑人巨大hd| 欧美极品一区二区三区四区| 国产色婷婷99| 少妇的逼水好多| 亚洲在线观看片| 好看av亚洲va欧美ⅴa在| 日韩国内少妇激情av| 女同久久另类99精品国产91| av天堂在线播放| 日韩中字成人| 免费av毛片视频| 综合色av麻豆| 一本一本综合久久| 欧美日韩亚洲国产一区二区在线观看| 国产伦精品一区二区三区四那| 精品一区二区三区视频在线| 直男gayav资源| 欧美高清成人免费视频www| 女同久久另类99精品国产91| 色噜噜av男人的天堂激情| 男女床上黄色一级片免费看| 国产一区二区三区视频了| 欧美三级亚洲精品| 久久精品国产自在天天线| 深夜精品福利| 麻豆国产97在线/欧美| 一个人免费在线观看的高清视频| 中亚洲国语对白在线视频| 午夜精品久久久久久毛片777| 国产v大片淫在线免费观看| 怎么达到女性高潮| 在线天堂最新版资源| 午夜福利在线在线| 搡老岳熟女国产| 国产欧美日韩精品亚洲av| 最好的美女福利视频网| 少妇人妻一区二区三区视频| 天堂√8在线中文| 日韩av在线大香蕉| 中出人妻视频一区二区| 中文字幕人妻熟人妻熟丝袜美| 国产精品人妻久久久久久| 全区人妻精品视频| 男人舔女人下体高潮全视频| 深夜a级毛片| 久久这里只有精品中国| 免费观看人在逋| 国产视频内射| 在线观看舔阴道视频| 欧美日韩黄片免| 午夜福利18| 婷婷色综合大香蕉| 综合色av麻豆| 在线免费观看不下载黄p国产 | 精品一区二区免费观看| а√天堂www在线а√下载| 老鸭窝网址在线观看| 亚洲精品久久国产高清桃花| x7x7x7水蜜桃| 亚洲最大成人手机在线| 特大巨黑吊av在线直播| 亚洲精品乱码久久久v下载方式| 麻豆一二三区av精品| 激情在线观看视频在线高清| 久久热精品热| 哪里可以看免费的av片| 午夜福利视频1000在线观看| 国产亚洲精品久久久com| 国产单亲对白刺激| 亚洲国产精品久久男人天堂| 久久精品国产亚洲av涩爱 | 国产成人欧美在线观看| 成人欧美大片| 人妻夜夜爽99麻豆av| 日韩成人在线观看一区二区三区| 欧美性感艳星| 最后的刺客免费高清国语| 国产日本99.免费观看| 十八禁国产超污无遮挡网站| 国内精品一区二区在线观看| 亚洲成人精品中文字幕电影| 女生性感内裤真人,穿戴方法视频| 亚州av有码| 噜噜噜噜噜久久久久久91| 亚洲经典国产精华液单 | 99国产综合亚洲精品| 18禁裸乳无遮挡免费网站照片| 99久久精品国产亚洲精品| aaaaa片日本免费| av视频在线观看入口| 亚洲精品一卡2卡三卡4卡5卡| 欧洲精品卡2卡3卡4卡5卡区| 国产日本99.免费观看| 欧美色欧美亚洲另类二区| 免费在线观看成人毛片| 久99久视频精品免费| 男人和女人高潮做爰伦理| aaaaa片日本免费| 亚洲第一电影网av| 欧美国产日韩亚洲一区| 亚洲aⅴ乱码一区二区在线播放| eeuss影院久久| 久久久久精品国产欧美久久久| 无人区码免费观看不卡| 日韩欧美精品免费久久 | 精品乱码久久久久久99久播| 婷婷六月久久综合丁香| 久久国产乱子伦精品免费另类| 久久国产乱子免费精品| 永久网站在线| 中文字幕av成人在线电影| 757午夜福利合集在线观看| 国模一区二区三区四区视频| 99在线人妻在线中文字幕| 麻豆成人午夜福利视频| 国内揄拍国产精品人妻在线| 此物有八面人人有两片| 亚洲av成人不卡在线观看播放网| 婷婷精品国产亚洲av| 亚洲精品456在线播放app | 日韩 亚洲 欧美在线| 亚洲熟妇熟女久久| 日韩免费av在线播放| 国产极品精品免费视频能看的| 欧美精品啪啪一区二区三区| 国产成人a区在线观看| 成人高潮视频无遮挡免费网站| 国产精品99久久久久久久久| 网址你懂的国产日韩在线| 在线观看av片永久免费下载| 黄色女人牲交| 婷婷六月久久综合丁香| 欧美不卡视频在线免费观看| 欧美区成人在线视频| 国产不卡一卡二| 国产精品精品国产色婷婷| 91久久精品国产一区二区成人| 超碰av人人做人人爽久久| 日韩 亚洲 欧美在线| av专区在线播放| 国产一级毛片七仙女欲春2| 在线观看66精品国产| 热99在线观看视频| .国产精品久久| 日韩中字成人| 欧美一区二区亚洲| 国产精品av视频在线免费观看| 搡老熟女国产l中国老女人| 成熟少妇高潮喷水视频| 午夜精品在线福利| 午夜福利在线观看吧| 国产精品美女特级片免费视频播放器| 真人做人爱边吃奶动态| 在线播放国产精品三级| 午夜日韩欧美国产| 国产精品一区二区性色av| 久久中文看片网| 草草在线视频免费看| 高清日韩中文字幕在线| 国产成人影院久久av| 搡老妇女老女人老熟妇| 国产精品亚洲av一区麻豆| 91在线观看av| 欧美+亚洲+日韩+国产| 国产一区二区激情短视频| 俺也久久电影网| 好男人在线观看高清免费视频| 成人欧美大片| 又粗又爽又猛毛片免费看| 九色成人免费人妻av| 制服丝袜大香蕉在线| 欧美中文日本在线观看视频| 麻豆一二三区av精品| 我的老师免费观看完整版| 久久人人爽人人爽人人片va | 99国产精品一区二区蜜桃av| 丰满人妻熟妇乱又伦精品不卡| 看黄色毛片网站| 亚洲一区二区三区色噜噜| 少妇被粗大猛烈的视频| 欧美黄色片欧美黄色片| 国产亚洲精品久久久com| 久久久久久久亚洲中文字幕 | 免费一级毛片在线播放高清视频| 动漫黄色视频在线观看| 精品一区二区免费观看| 欧美一区二区亚洲| 99热6这里只有精品| 国产黄色小视频在线观看| 蜜桃亚洲精品一区二区三区| 波多野结衣巨乳人妻| 91av网一区二区| 舔av片在线| 夜夜夜夜夜久久久久| 国产亚洲精品综合一区在线观看| av国产免费在线观看| 女人十人毛片免费观看3o分钟| 日韩有码中文字幕| 男女那种视频在线观看| 中文字幕av成人在线电影| 国产精品三级大全| 亚洲欧美精品综合久久99| 欧美黄色淫秽网站| 99视频精品全部免费 在线| 午夜a级毛片| 亚洲精品亚洲一区二区| 在线观看66精品国产| 久久精品国产清高在天天线| a在线观看视频网站| 黄色丝袜av网址大全| 伦理电影大哥的女人| 一区二区三区激情视频| 波野结衣二区三区在线| 亚洲第一区二区三区不卡| 亚洲成人久久性| 国产精品美女特级片免费视频播放器| 俺也久久电影网| 成年女人毛片免费观看观看9| 国内精品久久久久精免费| 午夜日韩欧美国产| 国产精品一区二区性色av| 久久亚洲精品不卡| 久久精品综合一区二区三区| 18美女黄网站色大片免费观看| 51国产日韩欧美| 一级黄片播放器| 亚洲欧美精品综合久久99| 嫩草影院精品99| 欧美日韩黄片免| 亚洲电影在线观看av| 免费av不卡在线播放| 日韩免费av在线播放| 成人欧美大片| 精品免费久久久久久久清纯| 国产成人影院久久av| 最后的刺客免费高清国语| 99在线视频只有这里精品首页| 成熟少妇高潮喷水视频| 国产av麻豆久久久久久久| 欧美精品啪啪一区二区三区| 伊人久久精品亚洲午夜| 国产精品国产高清国产av| 亚洲一区二区三区色噜噜| 亚洲一区二区三区色噜噜| 精品福利观看| 神马国产精品三级电影在线观看| 国产欧美日韩一区二区三| 嫩草影院精品99| 美女cb高潮喷水在线观看| 亚洲成a人片在线一区二区| 免费看日本二区| 国产成人a区在线观看| 一个人看视频在线观看www免费| 国产极品精品免费视频能看的| 亚洲国产精品sss在线观看| 国产高清视频在线观看网站| 亚洲黑人精品在线| 亚洲精品色激情综合| 国产一区二区在线av高清观看| 精品人妻视频免费看| 国产伦在线观看视频一区| 在线观看舔阴道视频| 在线观看免费视频日本深夜| 男人舔奶头视频| 久久精品国产亚洲av香蕉五月| 久久久久久久久大av| 黄色视频,在线免费观看| 禁无遮挡网站| 亚洲第一欧美日韩一区二区三区| 欧美极品一区二区三区四区| 天堂av国产一区二区熟女人妻| 日韩国内少妇激情av| 熟妇人妻久久中文字幕3abv| 国产精品98久久久久久宅男小说| av福利片在线观看| 国产一区二区三区视频了| 男人舔奶头视频| 久久欧美精品欧美久久欧美| 搡老岳熟女国产| 久久久成人免费电影| 欧美成人a在线观看| 亚洲aⅴ乱码一区二区在线播放| 中文字幕免费在线视频6| 九色国产91popny在线| 国产欧美日韩精品一区二区| 91在线观看av| 全区人妻精品视频| 午夜免费成人在线视频| 色噜噜av男人的天堂激情| 亚洲18禁久久av| 亚洲七黄色美女视频| 露出奶头的视频| 久久久久久大精品| 搡老妇女老女人老熟妇| 亚洲精品久久国产高清桃花| 免费在线观看成人毛片| 免费黄网站久久成人精品 | 校园春色视频在线观看| 中出人妻视频一区二区| 亚洲久久久久久中文字幕| 国产精品一区二区免费欧美| 欧美日韩国产亚洲二区| 亚洲中文字幕日韩| 午夜福利18| 亚洲成人中文字幕在线播放| 久久精品91蜜桃| 欧美日韩黄片免| 日韩欧美精品v在线| 国产高潮美女av| 久久久久久久久大av| 国产一区二区亚洲精品在线观看| 欧美日本亚洲视频在线播放| 亚洲,欧美精品.| 亚洲狠狠婷婷综合久久图片| 国产精华一区二区三区| 在线看三级毛片| 99久国产av精品| 国产在视频线在精品| 久久精品国产亚洲av涩爱 | 在现免费观看毛片| 国产精品亚洲av一区麻豆| 观看免费一级毛片| 久久久久国产精品人妻aⅴ院| 国产乱人伦免费视频| 两个人的视频大全免费| 最新中文字幕久久久久| 亚洲狠狠婷婷综合久久图片| 亚洲精品影视一区二区三区av| 十八禁人妻一区二区| 在现免费观看毛片| 国产成人欧美在线观看| 欧美成人免费av一区二区三区| 欧美三级亚洲精品| 国产高清视频在线观看网站| 毛片女人毛片| 观看免费一级毛片| 国产一级毛片七仙女欲春2| 国产一区二区亚洲精品在线观看| 99久久无色码亚洲精品果冻| 首页视频小说图片口味搜索| 亚洲成av人片在线播放无| 欧美+亚洲+日韩+国产| 亚洲精品粉嫩美女一区| 国产欧美日韩一区二区精品| 在线播放国产精品三级| 亚洲精品在线美女| 高清日韩中文字幕在线| 亚洲成av人片免费观看| 最近最新中文字幕大全电影3| www.熟女人妻精品国产| 禁无遮挡网站| 国产精品永久免费网站| 国产伦精品一区二区三区四那| 91在线观看av| 亚洲成人免费电影在线观看| 悠悠久久av| 精品久久久久久久久亚洲 | av天堂在线播放| 无遮挡黄片免费观看| 精品人妻1区二区| 最新中文字幕久久久久| 亚洲av日韩精品久久久久久密| 欧美日韩乱码在线| 97热精品久久久久久| 亚洲人成伊人成综合网2020| 最新在线观看一区二区三区| 亚洲在线观看片| 成人高潮视频无遮挡免费网站| 不卡一级毛片| 999久久久精品免费观看国产| 亚洲乱码一区二区免费版| 最新中文字幕久久久久| 九色成人免费人妻av| 欧美zozozo另类| 男人和女人高潮做爰伦理| 国产真实乱freesex| 欧美黄色淫秽网站| 一本一本综合久久| 国产视频内射| 高清在线国产一区| 中文字幕高清在线视频| 又爽又黄a免费视频| 人妻久久中文字幕网| 国产视频内射| 久久人妻av系列| 久久久久久久久中文| 国产精品一区二区三区四区免费观看 | 国产精品1区2区在线观看.| 国产极品精品免费视频能看的| 亚洲国产精品999在线| 99热只有精品国产| 国产伦人伦偷精品视频| 亚洲人与动物交配视频| 麻豆av噜噜一区二区三区| 尤物成人国产欧美一区二区三区| 日日摸夜夜添夜夜添av毛片 | 一进一出好大好爽视频| 欧美+亚洲+日韩+国产| 亚洲美女视频黄频| 国产精品一区二区三区四区久久| 免费人成视频x8x8入口观看| 麻豆av噜噜一区二区三区| 久久久精品大字幕| 亚洲专区国产一区二区| 精品久久久久久,| 亚洲av免费在线观看| 亚洲国产欧洲综合997久久,| 欧美成人a在线观看| aaaaa片日本免费| 一级av片app| 精品人妻一区二区三区麻豆 | 午夜福利18| 国产av不卡久久| 日韩高清综合在线| 99久国产av精品| 精品国内亚洲2022精品成人| 一级黄色大片毛片| 久久国产乱子伦精品免费另类| 国产高清视频在线观看网站| 一区二区三区免费毛片| 麻豆国产av国片精品| 国产大屁股一区二区在线视频| av天堂中文字幕网| 亚洲精品日韩av片在线观看| 欧美性感艳星| 亚洲国产色片| 亚洲不卡免费看| 国产精品不卡视频一区二区 | 又紧又爽又黄一区二区| 免费观看精品视频网站| av天堂中文字幕网| 午夜福利在线在线| 久久中文看片网| 亚洲黑人精品在线| 精品日产1卡2卡| 久久久久免费精品人妻一区二区| 国产视频一区二区在线看| av天堂在线播放| 麻豆国产97在线/欧美| 亚洲最大成人手机在线| 成人鲁丝片一二三区免费| 9191精品国产免费久久| 欧美三级亚洲精品| 免费看a级黄色片| or卡值多少钱| 三级毛片av免费| 久久国产精品人妻蜜桃| 窝窝影院91人妻| 三级男女做爰猛烈吃奶摸视频| 免费看a级黄色片| 欧美成狂野欧美在线观看| 欧美极品一区二区三区四区| 国内精品美女久久久久久| 亚洲无线在线观看| 亚洲不卡免费看| 在线播放国产精品三级| 男人舔奶头视频| 成年人黄色毛片网站| 一个人看的www免费观看视频| 国产精品女同一区二区软件 | 精品午夜福利视频在线观看一区| 久久久久久九九精品二区国产| 97热精品久久久久久| 欧美精品啪啪一区二区三区| 国产一区二区三区在线臀色熟女| 蜜桃久久精品国产亚洲av| 亚洲电影在线观看av| 国产欧美日韩精品亚洲av| 好男人在线观看高清免费视频| 成人欧美大片| 我的女老师完整版在线观看| 欧美另类亚洲清纯唯美| 国产av在哪里看| 男女之事视频高清在线观看| 又爽又黄a免费视频| 一进一出抽搐gif免费好疼| 免费高清视频大片| 久久久国产成人免费| av专区在线播放| 亚洲成人久久性| 国产av不卡久久| 日韩欧美在线二视频| 亚洲国产精品999在线| x7x7x7水蜜桃| 成年女人看的毛片在线观看| 亚洲aⅴ乱码一区二区在线播放| 伦理电影大哥的女人| 99国产精品一区二区蜜桃av| 在线播放无遮挡| 国产精品自产拍在线观看55亚洲| 小蜜桃在线观看免费完整版高清| 有码 亚洲区| 脱女人内裤的视频| 变态另类丝袜制服| 最近最新中文字幕大全电影3| 久久99热这里只有精品18| 色综合欧美亚洲国产小说| 久久久久亚洲av毛片大全| 一个人观看的视频www高清免费观看| 亚洲一区二区三区不卡视频| 国产三级黄色录像| 观看免费一级毛片| 一区福利在线观看| 亚洲欧美精品综合久久99| 亚洲国产精品sss在线观看| 99热这里只有是精品在线观看 | 精品免费久久久久久久清纯| 国产伦在线观看视频一区| 极品教师在线免费播放| 日韩欧美在线乱码| 国产av不卡久久| 国内精品久久久久精免费| 丰满人妻一区二区三区视频av| 别揉我奶头 嗯啊视频| 色噜噜av男人的天堂激情| 性欧美人与动物交配| 无遮挡黄片免费观看| 国产亚洲欧美98| 97碰自拍视频| 欧美日韩福利视频一区二区| 亚洲 欧美 日韩 在线 免费| 午夜老司机福利剧场| .国产精品久久| 丝袜美腿在线中文| 亚洲最大成人中文| 亚洲av一区综合| netflix在线观看网站| 露出奶头的视频| 国产精品野战在线观看| 欧美黑人欧美精品刺激| 十八禁国产超污无遮挡网站| 亚洲综合色惰| 国产aⅴ精品一区二区三区波| 亚洲av.av天堂| 级片在线观看| 美女黄网站色视频| 校园春色视频在线观看| 久久精品久久久久久噜噜老黄 | 亚洲aⅴ乱码一区二区在线播放| 欧美精品啪啪一区二区三区| 欧美一区二区亚洲| 欧美精品啪啪一区二区三区| 久久6这里有精品| 黄色一级大片看看| 丝袜美腿在线中文| 麻豆成人av在线观看| 在线观看午夜福利视频| 日本黄色视频三级网站网址| 国产精华一区二区三区| 一区二区三区激情视频| 亚洲欧美日韩高清专用| 日本黄大片高清| 午夜精品在线福利| 成人av一区二区三区在线看| 亚洲国产色片| 久久精品人妻少妇| 99热这里只有是精品在线观看 | 一级av片app| 久久国产乱子免费精品| 精品99又大又爽又粗少妇毛片 | 丰满的人妻完整版| 中文字幕人妻熟人妻熟丝袜美| 欧美黄色片欧美黄色片| 天堂网av新在线| 亚洲内射少妇av| 免费观看精品视频网站| 少妇的逼好多水| 欧美潮喷喷水| 午夜福利免费观看在线| 91九色精品人成在线观看| 国产伦人伦偷精品视频| 国产毛片a区久久久久| 欧美激情在线99| 亚洲精品一卡2卡三卡4卡5卡| 在线观看午夜福利视频| 日韩精品中文字幕看吧| 欧美成人性av电影在线观看| 国产乱人伦免费视频| 亚洲精品久久国产高清桃花| 午夜福利在线观看免费完整高清在 | 国产av不卡久久| www.熟女人妻精品国产| 国产精品av视频在线免费观看| 在线播放国产精品三级| 日韩欧美三级三区| 免费一级毛片在线播放高清视频| 国产精品日韩av在线免费观看| 成人欧美大片| 欧美最黄视频在线播放免费| 日韩中字成人| 老司机福利观看| 免费av毛片视频| 男女那种视频在线观看| 国产三级黄色录像| 国产精品99久久久久久久久| 91狼人影院| 亚洲aⅴ乱码一区二区在线播放| 久久国产精品影院| 午夜免费激情av| 国内精品美女久久久久久| 中文字幕精品亚洲无线码一区| 亚洲国产精品sss在线观看| 亚洲久久久久久中文字幕| 我的女老师完整版在线观看| 夜夜躁狠狠躁天天躁| 嫩草影院精品99| 又黄又爽又免费观看的视频| 嫁个100分男人电影在线观看| 一进一出抽搐gif免费好疼| 亚洲成a人片在线一区二区| 午夜a级毛片| 天堂网av新在线| 美女xxoo啪啪120秒动态图 | 最近视频中文字幕2019在线8| 欧美绝顶高潮抽搐喷水|