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

        ?

        基于泛克里格方法的樹種多樣性空間估計(jì)研究

        2013-12-29 02:59:19張會(huì)儒雷相東李曉冬
        關(guān)鍵詞:模型

        賀 鵬,張會(huì)儒,雷相東,李曉冬

        (中國林業(yè)科學(xué)研究院 資源信息研究所,北京 100091)

        基于泛克里格方法的樹種多樣性空間估計(jì)研究

        賀 鵬,張會(huì)儒,雷相東,李曉冬

        (中國林業(yè)科學(xué)研究院 資源信息研究所,北京 100091)

        生物多樣性保護(hù)已成為森林可持續(xù)經(jīng)營的一個(gè)重要目標(biāo)。生物多樣性的空間分布估計(jì)對于森林經(jīng)營者了解生物多樣性的空間格局和變化趨勢具有重要意義。本研究主要基于汪清林業(yè)局二類調(diào)查數(shù)據(jù)和局級(jí)固定樣地?cái)?shù)據(jù),以樹種多樣性為例,利用泛克里格方法對全局樹種多樣性進(jìn)行空間分布估計(jì)。結(jié)果表明,泛克里格方法的預(yù)測精度較高(R2=0.697 7),樹種多樣性在固定樣地的抽樣尺度下,存在一定的空間自相關(guān)性,預(yù)測的樹種多樣性最高值為1.82,最低值為0.221,空間分布具有明顯的空間異質(zhì)性,呈現(xiàn)出北部低,東部,南部和西部高的空間分布格局。研究結(jié)果為區(qū)域尺度基于固定樣地和地統(tǒng)計(jì)學(xué)的樹種多樣性空間分布估計(jì)提供了方法和參照

        樹種多樣性;空間分布格局;地統(tǒng)計(jì)學(xué);泛克里格方法

        生物多樣性保護(hù)已成為森林可持續(xù)經(jīng)營的一個(gè)重要目標(biāo)。森林經(jīng)營措施與生物多樣性密切相關(guān),所以森林經(jīng)營者需要了解森林經(jīng)營活動(dòng)后生物多樣性的變化, 確保物種和種群不會(huì)因森林采伐和更新而處于危險(xiǎn)境地, 從而對經(jīng)營和保護(hù)做出正確決策[1]。因此準(zhǔn)確地估計(jì)生物多樣性和分析生物多樣性的空間格局對于森林經(jīng)營者從整體上了解生物多樣性的分布格局及變化趨勢和制定適應(yīng)性經(jīng)營策略具有重要意義。

        地統(tǒng)計(jì)學(xué)的理論和方法可以定量地描述和解釋空間異質(zhì)性或者空間相關(guān)性;也可以建立各種空間預(yù)測模型,并進(jìn)行空間數(shù)據(jù)插值和估計(jì);同時(shí)能夠?qū)臻g格局的尺度、幾何形狀、變異方向進(jìn)行定量地分析和有效的估計(jì);有利于我們深刻地了解生命有機(jī)體(個(gè)體、種群和群落)的空間分布情況和空間異質(zhì)的機(jī)制[2]。因此地統(tǒng)計(jì)學(xué)在林業(yè)中得到了比較廣泛地應(yīng)用,許多學(xué)者經(jīng)常運(yùn)用其對森林結(jié)構(gòu)進(jìn)行分析,例如胸徑地理變異規(guī)律分析[3],胸高斷面積空間插值和空間格局分析[4-6],蓄積空間插值和空間格局分析[7-8],和生物量空間插值和空間格局分析[9-10]等。

        生物多樣性受到各種生物因素和非生物因素的影響,所以生物多樣性不是純隨機(jī)變量,而是區(qū)域化變量,具有隨機(jī)性和結(jié)構(gòu)性雙重屬性。因此本研究以樹種多樣性為例,基于吉林省汪清林業(yè)局局級(jí)固定樣地和二類調(diào)查數(shù)據(jù),利用地統(tǒng)計(jì)學(xué)的原理和方法對整個(gè)林業(yè)局13個(gè)林場的樹種多樣性進(jìn)行空間估計(jì)。實(shí)現(xiàn)由樣地“點(diǎn)”上樹種多樣性信息擴(kuò)展到整個(gè)區(qū)域“面”上樹種多樣性信息。同時(shí)利用插值得到整個(gè)林場的樹種多樣性空間格局分布圖,為在區(qū)域尺度內(nèi),基于固定樣地和地統(tǒng)計(jì)學(xué)的樹種多樣性空間分布估計(jì)提供方法和參照。

        1 研究區(qū)概況

        汪清林業(yè)局位于吉林省延邊自治州的東部,所處的地理坐標(biāo)為東經(jīng) 123°56′~ 131°04′,北緯43°05′~ 43°40′,東西跨度約為 85 km, 南北跨度約為60 km。該局下設(shè)13個(gè)林場,全局經(jīng)營面積為304 173 hm2。汪清林業(yè)局屬于長白山系的中低丘陵地區(qū), 海拔為360~1 477 m,氣候?qū)儆跍貛Т箨懶约撅L(fēng)氣候,平均氣溫為3.9 ℃, 最低溫度為-37.5 ℃,平均降水量為547 mm,且主要集中在夏季,占全年總降水量的80%。該區(qū)森林覆蓋率極高且植被種類繁多,主要以天然林為主,占有林地面積的93.7%。而天然林中以闊葉混交林、蒙古櫟林和針闊混交林為主,這三林分占全局有林地面積的73.8%、蓄積的77.6%。主要針葉樹種有長白落葉松Larix olgensis、冷杉Abies nephrolepis、紅松Pinus koraiensis、云杉Picea jezoensis等。主要闊葉樹種蒙古櫟Quercus mongolica、白樺Betula platyphylla、楊樹Populus ussuriensis、楓樺Betula costata、色木Acer mono、水曲柳Fraxinus mandshurica、黃菠蘿Phellodendron amurense等。

        2 數(shù)據(jù)與方法

        2.1 數(shù)據(jù)來源

        研究數(shù)據(jù)來源于吉林省汪清林業(yè)局2007年二類調(diào)查數(shù)據(jù),包括局級(jí)固定樣地?cái)?shù)據(jù)。其中固定樣地共1 391個(gè)(剔除固定樣地?cái)?shù)據(jù)中蓄積為0的樣地后),為系統(tǒng)布設(shè)樣地,間距為1 km×2 km,樣地面積為0.06 hm2,樣地空間分布圖如圖1a(圖中樣地顏色越深表示樹種多樣性越高)。固定樣地主要調(diào)查因子為樹種、胸徑、坡向、坡度、海拔等,起測徑階為5 cm。

        汪清林業(yè)局林場圖件:林相圖及相關(guān)的資源數(shù)據(jù)庫。圖1b表示全局的純林和混交林的空間分布圖。

        圖1 樣地分布(a)和純林及混交林分布(b)Fig. 1 Spatial distribution of plots (a) and distribution of forest styles (b)

        2.2 樣地樹種多樣性指標(biāo)計(jì)算及統(tǒng)計(jì)

        根據(jù)固定樣地每木檢尺的數(shù)據(jù),選取Shanonwiener指數(shù)(公式1)計(jì)算每個(gè)樣地的樹種多樣性指數(shù)。計(jì)算后樣地樹種多樣性指數(shù)的統(tǒng)計(jì)描述如表1,同時(shí)隨機(jī)預(yù)留出10%的樣地用來檢驗(yàn)方法的預(yù)測精度。

        表1 固定樣地樹種多樣性指數(shù)統(tǒng)計(jì)Table 1 Statistical characteristics of tree species diversity indexes for studied area

        式(1)中HSpecies表示樹種多樣性指數(shù),pi表示各樹種胸高斷面積占林分胸高斷面積的比例。

        2.3 泛克里格方法

        在泛克里格方法(Universal Kriging)中, 預(yù)測值z(s0)由兩部分組成。一部分是m(s0)估計(jì)值,通過不同環(huán)境解釋變量對主變量的回歸估計(jì),s0表示待估測點(diǎn)位置。第二部分是通過克里格插值殘差e(s0)。泛克里格估計(jì)的公式如下[11]:

        式(3)中,βk表示趨勢模型的估計(jì)參數(shù),本研究中選擇普通最小二乘法進(jìn)行估計(jì)。wi表示不同采樣點(diǎn)位置殘差的權(quán)重,此權(quán)重由殘差半變異函數(shù)確定。e表示回歸模型的殘差。本文中所有回歸估計(jì)和泛克里格插值均在R統(tǒng)計(jì)軟件gstat包中進(jìn)行[12]。

        地統(tǒng)計(jì)建模中的第一步是對樣本計(jì)算試驗(yàn)半變異函數(shù)值。采用以下公式4:

        式(4)中:r(h)為回歸殘差的變異函數(shù);N(h)是在某一空間方位上相距距離為h的固定樣地對數(shù);z(xi)和z(xi+h)分別是固定樣地xi和與xi相隔距離h處的回歸殘差值。采用不同的模型,如指數(shù)模型,球狀模型等對試驗(yàn)半變異函數(shù)值進(jìn)行擬合,得到理論半變異函數(shù)模型。一般半變異函數(shù)r(h)是一個(gè)單調(diào)遞增函數(shù),隨距離h增加而增加,當(dāng)超出一定距離(變程a)時(shí),r(h)不再增加,此時(shí)的值稱為“基臺(tái)值(C)”。而理論變異函數(shù)在原點(diǎn)的取值稱為“塊金值(C0)”,塊金值表示存在小于抽樣尺度的變異或者存在測量誤差。基臺(tái)值和塊金值之差稱為“偏基臺(tái)值(C1)”[13]。

        本研究中對樹種多樣性趨勢模型采用的解釋變量為每個(gè)樣地的坐標(biāo)值和森林類型,森林類型變量為一個(gè)分類變量,取0表示純林,取1表示混交林。趨勢模型的具體表達(dá)式如下:

        m(s)=β0+β1x2(s)+ β2y2(s)+ β3xy(s)+β4x(s)+β5y(s)+ β6f(s) (5)

        式(5)中:β表示各個(gè)解釋變量的系數(shù),f(s)的值取決于s點(diǎn)的森林類型。

        2.4 模型評(píng)價(jià)和檢驗(yàn)

        采用交叉檢驗(yàn)(cross-validation)來評(píng)價(jià)變異函數(shù)模型擬合效果。交叉檢驗(yàn)分析重復(fù)從已知數(shù)據(jù)集中刪除一個(gè)采樣點(diǎn),用剩余的采樣點(diǎn)估測刪除點(diǎn)的數(shù)值[14-16]。對于泛克里格方法預(yù)測的精度, 采用隨機(jī)預(yù)留的140個(gè)檢驗(yàn)樣本進(jìn)行檢驗(yàn)。本研究中通過計(jì)算平均誤差(ME), 平均絕對誤差(MAE),均方根誤差(RMSE)和決定系數(shù)(R2)來對模型進(jìn)行檢驗(yàn)和預(yù)測精度進(jìn)行評(píng)估。

        各式中Xf和X0分別表示樹種多樣性的預(yù)測值和實(shí)測值,表示實(shí)測值的平均值。

        3 結(jié)果與分析

        3.1 回歸模型和殘差半變異函數(shù)分析

        表2 回歸模型各變量參數(shù)的估計(jì)值Table 2 Estimation of parameters of regression models

        表2中的結(jié)果為通過普通最小二乘法估計(jì)回歸模型式(4)中各個(gè)變量的系數(shù)。β0截距包含了變量f取0(純林)時(shí)的信息。表明樹種多樣性與純林存在著負(fù)相關(guān)關(guān)系,β6為變量f取1(混交林)的系數(shù),表明樹種多樣性與混交林存在正相關(guān)關(guān)系。

        通過趨勢回歸模型對每個(gè)樣地的樹種多樣性成分進(jìn)行回歸估計(jì),從而得到每個(gè)樣地的樹種多樣性殘差,在此基礎(chǔ)上對殘差進(jìn)行半變異函數(shù)結(jié)構(gòu)分析。本研究采用了不同類型了半變異函數(shù)理論模型對殘差的試驗(yàn)半變異函數(shù)進(jìn)行擬合,包括指數(shù)模型,高斯模型,球狀模型,圓形模型。其中指數(shù)模型擬合效果最好,塊金值(C0)為0.070 31,偏基臺(tái)值(C1)為 0.022 87, 變程(a)達(dá) 3 000 m。這表明在變程范圍內(nèi),當(dāng)前固定樣地抽樣尺度下各個(gè)樣地之間存在著一定的空間自相關(guān)性。表3為理論半變異模型的參數(shù)及交叉檢驗(yàn)統(tǒng)計(jì)量。圖3表示通過試驗(yàn)半變異函數(shù)值擬合的理論半變異函數(shù)模型。 從表3和圖2中均可以得出趨勢回歸模型殘差存在著空間自相關(guān)性。并且交叉檢驗(yàn)表明理論半變異函數(shù)擬合效果比較理想。

        表3 理論半變異函數(shù)參數(shù)和交叉檢驗(yàn)統(tǒng)計(jì)量Table 3 Parameters of theoretical semivariance model and statistical quantity of cross-validation

        圖2 趨勢回歸模型殘差的半變異函數(shù)Fig.2 Semivariances of residuals of regression model

        3.2 檢驗(yàn)

        泛克里格空間預(yù)測的精度,采用隨機(jī)預(yù)留的140個(gè)樣地的實(shí)測值與預(yù)測進(jìn)行統(tǒng)計(jì)比較,結(jié)果如圖 3, y=0.408 78+0.699 72x, R2=0.697 7,表明泛克里格方法預(yù)測精度較高,是一種可行的方法。

        圖3 檢驗(yàn)樣本的實(shí)測值與預(yù)測值統(tǒng)計(jì)比較Fig. 3 Statistical comparison between measured and predicted diversity in validation dataset

        3.3 空間分布估計(jì)及分析

        通過泛克里格方法空間插值得到整個(gè)全局的樹種多樣性空間分布格局圖(圖4a)和標(biāo)準(zhǔn)誤差預(yù)測圖(圖4b)。從圖4a可看出,樹種多樣性空間分布具有明顯的異質(zhì)性,預(yù)測的空間分布圖中樹種多樣性最高值為1.82,最低值為0.221。在北部,西北部和南部邊緣地區(qū),樹種多樣性較低。這些的地區(qū)主要以蒙古櫟純林為主,這是導(dǎo)致這些地區(qū)物種多樣性普遍偏低的主要原因。而且主要道路和耕地兩邊的樹種多樣性也普遍偏低,這些地區(qū)可能受到人為干擾較大。而東部,西部,中部和南部的樹種多樣比較高,主要因?yàn)檫@些地區(qū)的森林類型以闊葉混交林,針闊混交林和針葉混交林為主。從整體上看,汪清林業(yè)局地區(qū)森林的樹種多樣性比較高,這與林業(yè)局主要以天然林為主,占有林地面積的93.7%,天然林保護(hù)類型主要以限伐區(qū)和禁伐區(qū)為主, 因此受到人為干擾比較小。從圖4b可以看出,在林場邊緣地帶預(yù)測標(biāo)準(zhǔn)誤差較大,這與樣地?cái)?shù)量較少有關(guān),影響插值精度。同時(shí)標(biāo)準(zhǔn)誤差預(yù)測圖,也給出來對于每個(gè)點(diǎn)的預(yù)測的可信度。標(biāo)準(zhǔn)誤差越高的地方表明需要加密樣地,來提高預(yù)測精度。

        圖4 樹種生物多樣性空間分布(a)和標(biāo)準(zhǔn)預(yù)測誤差分布(b)Fig. 4 Prediction spatial distribution (a) and standard error distribution (b)

        4 結(jié)論與討論

        樹種多樣性的空間分布估計(jì)可以為森林經(jīng)營者了解森林經(jīng)營活動(dòng)對森林生態(tài)系統(tǒng)生物多樣性的影響和制定合理的可持續(xù)經(jīng)營的方案提供依據(jù)。本文基于汪清林業(yè)局二類調(diào)查數(shù)據(jù)和局級(jí)固定樣地?cái)?shù)據(jù),采用泛克里格方法對整個(gè)汪清林業(yè)局13個(gè)林場的樹種多樣性指數(shù)進(jìn)行估計(jì),并得到樹種多樣性的空間格局分布圖。研究結(jié)果表明泛克里格方法空間插值方法是一種可行的方法,它可以解決尺度推移的問題,根據(jù)二類調(diào)查固定樣地的數(shù)據(jù),就可以得到整個(gè)林業(yè)局甚至是更大區(qū)域水平上的樹種多樣性估計(jì),現(xiàn)實(shí)了由“點(diǎn)”上的樹種多樣性信息擴(kuò)展到“面”上樹種多樣性信息。

        從空間估計(jì)的結(jié)果來看,汪清林業(yè)局地區(qū)森林的樹種多樣性具有明顯的空間格局和空間異質(zhì)性,這與當(dāng)?shù)氐膶?shí)際情況相符合,蒙古櫟林中樹種多樣性較小,而闊葉混交林,針葉混交林和針闊混交林中樹種多樣性較高。

        采用泛克里格方法對樹種多樣性進(jìn)行估計(jì),首先對趨勢進(jìn)行回歸分析,然后對殘差進(jìn)行半變異函數(shù)結(jié)構(gòu)分析, 進(jìn)而對樣地實(shí)驗(yàn)半變異函數(shù)值擬合理論半變異函數(shù)模型,最后進(jìn)行局部空間估計(jì)。這種方法不僅考慮了解釋變量(空間位置和森林類型)對樹種多樣性的影響,同時(shí)也考慮了樹種多樣性的空間自相性。今后還可以結(jié)合不同的環(huán)境變量(如海拔,氣候因子等)和遙感信息進(jìn)行研究,可以有效地提高預(yù)測精度。

        [1] 雷相東, 唐守正. 林分結(jié)構(gòu)多樣性指標(biāo)研究綜述[J]. 林業(yè)科學(xué),2002,38(3)∶140-146.

        [2] 馮益明. 空間統(tǒng)計(jì)學(xué)理論及其在林業(yè)中的應(yīng)用[M]. 北京∶ 中國林業(yè)出版社, 2008,1-123.

        [3] 洪 偉, 吳承幀. 杉木種源胸徑生長地理變異規(guī)律的研究[J]. 植物生態(tài)學(xué)報(bào) , 1998,22(2)∶ 186-192.

        [4]Malhi Y, Wood D, Baker T, et al. The regional variation of aboveground live biomass in old-growth Amazonian forests [J]. Global Change Biology, 2006,12∶1107-1138.

        [5]Meng Q, Cieszewski C, Madden M. Large area forest inventory using Landsat ETM+∶ A geostatistical approach[J]. ΙSPRS Journal of Photogrammetry and Remote Sensing,2009, 64∶27-36.

        [6]Akhavan R, Kia-Daliri H. Spatial variability and estimation of tree attributes in a plantation forest in the Caspian region of Ιran using geostatistical analysis[J]. Caspian Journal of Environmental Sciences,2010, 8(2)∶ 163-172.

        [7]李明陽 , 劉 敏 , 劉米蘭 . 基于 GΙS 的森林調(diào)查因子地統(tǒng)計(jì)學(xué)分析 [J]. 南京林業(yè)大學(xué)學(xué)報(bào):自然科學(xué)版 , 2010, (6)∶ 66-70.

        [8]Sales M H, Souza Jr C M, Kyriakidis P C. Ιmproving spatial distribution estimation of forest biomass with geostatistics∶ A case study for Rond?nia, Brazil[J]. Ecological Modelling, 2007,205(1/2)∶ 221-230.

        [9]Du H Q, Zhou G M, Fan W Y, et al. Spatial heterogeneity and carbon contribution of aboveground biomass of moso bamboo by using geostatistical theory[J]. Plant Ecol, 2010,207∶131-139.

        [10] 劉曉梅 , 布仁倉 , 鄧華衛(wèi) , 等 . 基于地統(tǒng)計(jì)學(xué)豐林自然保護(hù)區(qū)森林生物量估測及空間格局分析[J]. 生態(tài)學(xué)報(bào), 2011,31(16)∶ 4783-4790.

        [11] Cressie N A C. Statistic for spatial data[M]. John Wiley&Sons, New York,1993. 151-172.

        [12] Pebesma E J. Multivariable Geostatistics in S∶ the Gstat Package[J]. Computer & Geoderma, 2004,67∶ 215-226.

        [13] 王政權(quán). 地統(tǒng)計(jì)學(xué)及在生態(tài)學(xué)中的應(yīng)用[M]. 北京∶ 科學(xué)出版社 , 1999.

        [14] 湯國安 , 楊 昕 . ArcGΙS 地理信息系統(tǒng)空間分析實(shí)驗(yàn)教程[M]. 北京:科學(xué)出版社 , 2006.

        [15] 羅 梅 , 鄭小賢 . 八達(dá)嶺遼東櫟 -油松混交林空間結(jié)構(gòu)及其多樣性 [J]. 中南林業(yè)科技大學(xué)學(xué)報(bào) , 2012, 32(9)∶ 55-58.

        [16] 李際平 , 趙春燕 , 袁曉紅 , 等 . 杉闊間物種多樣性影響因子灰色關(guān)聯(lián)分析 [J]. 中南林業(yè)科技大學(xué)學(xué)報(bào) , 2011, 31(10)∶ 10-14.

        Estimation of spatial distribution of tree species diversity based on Universal Krige Model

        HE Peng, ZHANG Hui-ru, LEΙ Xiang-dong, LΙ Xiao-dong
        (Ιnstitute of Forest Resource Ιnformation Techniques, Chinese Academy of Forestry, Beijing 100091, China)

        Biodiversity conservation is one of major goals of forest sustainable management. Estimation of spatial distribution is extremely signif i cant for forest manager to get better sense of the spatial pattern and dynamical changes of biodiversity. By taking the tree species diversity as an example and using Universal Krige (UK) , the spatial distribution of tree species diversity based on bureaulevel permanent plot data of forest inventory in Wangqing forest enterprise, Jilin Province were estimated. The results demonstrate that the UK for the studied area was a feasible method with a higher prediction accuracy (R2=0.6977); Spatially, tree species diversity varied obviously in this region from 1.82 to 0.221; Ιn addition, the spatial distribution map shows that there was an obvious trend in the studied area∶ lower tree species diversity in north and higher in east, west and south. Meanwhile, the results provide a method and reference for the region-scale tree species diversity estimation based on the permanent plot and geo-statistics.

        tree species diversity; spatial distribution pattern; geo-statistics; universal kriging

        S718

        A

        1673-923X(2013)12-0067-05

        2013-03-25

        林業(yè)公益性行業(yè)專項(xiàng)“我國典型森林類型健康經(jīng)營關(guān)鍵技術(shù)研究”(201004002)

        賀 鵬(1988-),男, 碩士,主要研究方向:森林資源經(jīng)營管理與決策

        張會(huì)儒,博士,研究員,主要研究方向:森林可持續(xù)經(jīng)營研究; E-mail:huiru@caf.ac.cn

        [本文編校:吳 彬 ]

        猜你喜歡
        模型
        一半模型
        一種去中心化的域名服務(wù)本地化模型
        適用于BDS-3 PPP的隨機(jī)模型
        提煉模型 突破難點(diǎn)
        函數(shù)模型及應(yīng)用
        p150Glued在帕金森病模型中的表達(dá)及分布
        函數(shù)模型及應(yīng)用
        重要模型『一線三等角』
        重尾非線性自回歸模型自加權(quán)M-估計(jì)的漸近分布
        3D打印中的模型分割與打包
        中文字幕精品一区二区日本| 97se亚洲精品一区| 97se在线| 亚洲一级无码AV毛片久久| 久久综合五月天啪网亚洲精品| 新婚人妻不戴套国产精品| 久久久人妻精品一区bav| 蜜桃视频一区二区在线观看| 成av免费大片黄在线观看| 国产偷国产偷亚洲欧美高清| 亚洲第一女人天堂av| 我和丰满妇女激情视频| 久激情内射婷内射蜜桃人妖| 人妻熟妇乱系列| 中文字幕国产精品专区| 青青草精品视频在线播放| 99精品国产一区二区三区a片| 亚洲AV无码久久久一区二不卡 | 日产精品一区二区免费| 亚洲av老熟女一区二区三区| 国产精品一卡二卡三卡| 麻豆精品久久久久久久99蜜桃| 亚洲国产不卡av一区二区三区| 亚洲国产天堂久久综合网| 成 人 免费 在线电影| 久久久精品3d动漫一区二区三区| 精品视频一区二区杨幂| 日韩一区二区中文字幕视频| 亚洲精品国产一二三无码AV| 亚洲饱满人妻视频| 精品亚洲不卡一区二区| 人成在线免费视频网站| 国产精品免费一区二区三区四区| 欧美jizzhd精品欧美| 精品久久久久久午夜| 国产亚洲综合另类色专区 | 一本一道av无码中文字幕﹣百度| 无码 制服 丝袜 国产 另类 | 国产在线无码免费视频2021| 伊人影院成人在线观看| 日本爽快片100色毛片|