国产日韩欧美一区二区三区三州_亚洲少妇熟女av_久久久久亚洲av国产精品_波多野结衣网站一区二区_亚洲欧美色片在线91_国产亚洲精品精品国产优播av_日本一区二区三区波多野结衣 _久久国产av不卡

?

基于高光譜特征和光合參數(shù)監(jiān)測松小蠹不同危害時期的相關(guān)分析

2020-04-21 04:00:10劉夢盈馬云強(qiáng)張忠和劉雪蓮
林業(yè)科學(xué)研究 2020年2期
關(guān)鍵詞:云南松針葉反射率

劉夢盈,石 雷,馬云強(qiáng),張忠和,劉雪蓮,姚 俊,杜 婷

(1. 中國林業(yè)科學(xué)研究院資源昆蟲研究所,云南 昆明 650224;2. 南京林業(yè)大學(xué),江蘇 南京 210037;3. 西南林業(yè)大學(xué),云南 昆明 650224)

云南松(Pinus yunnanensisFranch.)是我國西南地區(qū)主要樹種,云南省主要森林植被類型之一,占林地面積52%,在林業(yè)生產(chǎn)中占有舉足輕重的地位[1]。松小蠹(Tomicusspp.)是云南松的主要害蟲之一,云南至今已有超過600 km2的云南松林死亡[2]。松小蠹蛀梢和蛀干習(xí)性使其長時間隱藏在樹體內(nèi)部,不易發(fā)現(xiàn)和防治。5—11月為蛀梢為害期,成蟲蛀食當(dāng)年生枝梢;11月至次年4月為蛀干為害期,成蟲在樹干韌皮部和木質(zhì)部鉆蛀、產(chǎn)卵,蛀梢和蛀干交替出現(xiàn),歷時長,危害嚴(yán)重。松小蠹還具備遷飛習(xí)性,擴(kuò)散快,導(dǎo)致松林成片死亡[3]。

高光譜遙感技術(shù)主要用于研究植被受病蟲害危害前后變化及受害程度和植被光譜曲線、植被指數(shù)變化之間的關(guān)系,確定不同樹種的危害敏感波段和時期,建立病蟲害監(jiān)測模型。劉暢[4]、Ju[5]等利用非成像高光譜儀測量銳齒櫟、馬尾松葉片光譜曲線,測定其光合色素含量,建立了光譜和光合的聯(lián)系。王震[6]、伍南[7]等基于高光譜特征建立了馬尾松松材線蟲病、赤枯病的病情程度反演模型。地面成像光譜能獲取地物光譜數(shù)據(jù)及實(shí)現(xiàn)“圖譜合一”,利于混合光譜分析和病蟲害監(jiān)測模型反演。徐爽[8]、Shafri[9]等利用地面成像高光譜儀實(shí)現(xiàn)了長棗表面蟲傷、靈芝基莖腐病監(jiān)測。于泉州[10]、宋小寧[11]基于高光譜數(shù)據(jù)進(jìn)行森林資源調(diào)查,提取森林冠層光譜特征、冠層含水量等信息。劉良云等[12]利用多時相高光譜航空影像,實(shí)現(xiàn)了冬小麥銹病監(jiān)測。

目前,利用遙感技術(shù)監(jiān)測松小蠹危害的研究多基于多光譜遙感進(jìn)行,而對地面高光譜成像儀監(jiān)測松小蠹危害情況的關(guān)注較少。因此,本研究通過實(shí)地調(diào)查松小蠹危害時期和發(fā)生情況,室內(nèi)采集云南松枝梢針葉光譜反射率,分析光譜特征,建立松小蠹危害時期監(jiān)測模型,實(shí)現(xiàn)危害早期監(jiān)測。并結(jié)合野外測定的不同危害時期云南松針葉光合參數(shù),建立光譜與光合的聯(lián)系,利用光譜數(shù)據(jù)實(shí)現(xiàn)光合參數(shù)預(yù)測,估算云南松的生長健康狀況。

1 試驗(yàn)地概況與研究方法

1.1 試驗(yàn)地概況

試驗(yàn)地位于石林彝族自治縣西街口鎮(zhèn),海拔1 700~1 950 m,地處 103°29′20″ E,24°49′25″ N。亞熱帶低緯度高原山地季風(fēng)氣候,干濕分明,5—10月為雨季,11月到次年4月為旱季,年均氣溫14.7 ℃。試驗(yàn)地以云南松為優(yōu)勢種,屬針葉林,大多為中幼齡林,約80%為純林,林分單一,且云南松大多生長在土壤貧瘠、水分條件缺乏,立地條件差的地段,易受病蟲害侵?jǐn)_,具體采樣地點(diǎn)如圖1所示。

1.2 研究方法

1.2.1 材料 試驗(yàn)材料為受松小蠹危害的云南松枝梢,根據(jù)松小蠹對云南松的實(shí)際危害情況,以單梢為單位,將松小蠹危害時期分為4個時期,健康時期、危害早期、危害中期與危害末期(表1)[13]。通過前期調(diào)查,根據(jù)外觀特征分別選擇不同危害時期的云南松枝梢,每個時期各選擇30個枝梢。

1.2.2 云南松枝梢光譜反射率測定 本研究采用的高光譜成像系統(tǒng)由SOC710VP光譜儀、2個500 W的鎢燈、計(jì)算機(jī)、暗箱等組成。針葉光譜測定在可控光照條件的實(shí)驗(yàn)室進(jìn)行,光譜儀架設(shè)在暗箱中,用鎢燈進(jìn)行照明。在圖像采集時,設(shè)置以下參數(shù):曝光時間35 ms,光圈大小5.6,物鏡高度50 cm;鎢燈光源入射點(diǎn)與鏡頭垂直于載物臺的焦點(diǎn)重合。將針葉平鋪在載物臺上,打開鎢燈,用標(biāo)準(zhǔn)參考板校正。啟動光譜儀,進(jìn)行光譜采集。

圖 1 采樣點(diǎn)示意圖Fig. 1 Schematic diagram of sampling points

表 1 松小蠹危害時期劃分Table 1 Hazard time division of Tomicus spp.

1.2.3 氣體交換實(shí)驗(yàn) 2018年10—12月,于晴朗天氣下午15:00—17:00,在野外使用LI-6400便攜式光合作用測定系統(tǒng)測定健康時期、危害早期、危害中期3個時期的云南松枝梢針葉凈光合速率(Photo)、氣孔導(dǎo)度(Cond)、胞間CO2濃度(Ci)及蒸騰速率(Trmmol)等參數(shù)(危害末期云南松枝梢接近枯死,故不測定其氣體交換參數(shù))。在自然光源條件下,選取陽面枝梢針葉進(jìn)行測定,每組處理中每個枝梢重復(fù)測定5次,取平均值。利用Photo和Trmmol計(jì)算水分利用效率[14](WUE=Photo/Trmmol);利用 Ci和空氣中的CO2濃度(Ca)計(jì)算氣孔限制值[15](Ls=1-Ci/Ca)。

1.3 數(shù)據(jù)分析

1.3.1 高光譜數(shù)據(jù)處理 利用SRAnal710軟件,對高光譜圖像數(shù)據(jù)進(jìn)行光譜、黑場、空間光譜輻射標(biāo)定等預(yù)處理,將DN值轉(zhuǎn)換為反射率。采用Savitzky-Golay濾波法對針葉光譜反射率曲線進(jìn)行平滑,以消除原始數(shù)據(jù)中的噪聲干擾[16],在OriginPro8.0中繪制針葉光譜曲線。

采用Pearson相關(guān)分析與逐步回歸分析篩選出特征變量,用以建模。本研究基于前人研究基礎(chǔ),結(jié)合云南松光譜反射率數(shù)據(jù)定義了8個植被指數(shù)、23個高光譜特征參數(shù)[17-23](表2)。

表 2 特征變量表Table 2 Characteristic variable scale

1.3.2 監(jiān)測模型及精度評價 將數(shù)據(jù)量的70%作為建立模型的數(shù)據(jù),30%作為精度驗(yàn)證的數(shù)據(jù)。主要涉及以下模型:

多元線性回歸模型:y=β0+β1x1+β2x2+β3x3+...+βnxn

線性函數(shù)模型:y=β0+β1x

二次函數(shù)模型:y=β0+β1x+β2x2

三次函數(shù)模型:y=β0+β1x+β2x2+β3x3

式中:y代表擬合值;x,x1,x2,x3, ...,xn代表主要特征參數(shù);β0,β1,β2,β3, ...,βn代表擬合系數(shù)。

本研究使用決定系數(shù)R2和均方根誤差RMSE進(jìn)行精度驗(yàn)證。R2值在0~1之間,R2越接近1,說明預(yù)測模型參考價值越高。RMSE值越小,說明預(yù)測模型精度更高。

2 結(jié)果與分析

2.1 云南松受松小蠹攻擊后不同危害時期針葉的光譜變化

2.1.1 不同危害時期云南松光譜反射率分析 不同危害時期針葉光譜反射率見圖2。針葉光譜反射率變化可反映松小蠹危害時期變化,健康時期、危害早期和危害中期針葉光譜反射率曲線具有一致性:在740~1 036 nm,隨著危害時期變化,反射率逐漸降低;在680~740 nm(“紅邊”),健康時期斜率最大,然后依次為危害早期和危害中期;在510~560 nm( “綠 峰 ”) 和 640~680 nm( “紅谷”),健康時期光譜反射率低于危害中期和危害早期。危害末期光譜反射率曲線“綠峰”、“紅谷”消失,大致為一條緩慢上升的曲線。從整體上看,反射率大小依次為健康時期、危害早期、危害中期、危害末期。

2.1.2 不同危害時期的云南松光譜微分曲線分析不同危害時期云南松光譜一階微分曲線見圖3,在500~540 nm,健康時期、危害早期和危害中期存在明顯波峰,且值大小依次為危害早期、危害中期、健康時期,而危害末期曲線不存在明顯波峰;在660~740 nm,存在明顯波峰,且峰值達(dá)到最高,其值大小依次為健康時期、危害早期、危害中期。隨著危害時期變化,出現(xiàn)“紅邊”“藍(lán)移 ”現(xiàn)象,峰值向短波方向移動且逐漸降低;在760~1 036 nm,各危害時期的波峰、波谷出現(xiàn)波段范圍大致相似。

2.1.3 不同危害時期的云南松光譜相關(guān)性分析 松小蠹危害時期與光譜反射率、光譜一階微分值相關(guān)性分析見圖4、圖5,在399、409、503~539、553~574、610~677、687~697、708~791、807~838、849~870、881~988、1 004、1 025~1 036 nm,松小蠹危害時期與針葉光譜反射率及光譜一階微分值極顯著相關(guān),可有效反映松小蠹危害時期,利用以上波段構(gòu)建的高光譜參數(shù)建立松小蠹危害時期監(jiān)測模型是可行的。

圖 2 不同危害時期云南松枝梢針葉光譜反射率曲線Fig. 2 Spectral reflectance curve of Pinus yunnanensis in different damage periods

圖 3 不同危害時期針葉光譜反射率一階微分曲線Fig. 3 First order differential curve of spectral reflectance of Pinus yunnanensis in different damage periods

圖 4 松小蠹不同危害時期與光譜反射率的相關(guān)性分析Fig. 4 The correlation between hazard time division of Tomicus spp. and spectral reflectance

圖 5 不同危害時期與光譜一階微分值的相關(guān)性分析Fig. 5 The correlation between hazard time division of Tomicus spp. and the first derivative of the spectrum

2.1.4 松小蠹危害時期診斷模型建立及驗(yàn)證 通過相關(guān)分析和逐步回歸分析(a)分別對植被指數(shù)和高光譜特征參數(shù)篩選用于建模的變量,獲得顯著或極 顯著變 量 :NDVI2、PRI;D、Dr、SDb、SDy(表3)。將植被指數(shù)和高光譜特征變量共同輸入SPSS中,進(jìn)行逐步回歸分析(b),獲得顯著或極顯著變量:D、Dr、NDVI705(表3)。

基于篩選的變量分別建立松小蠹危害時期的預(yù)測模型,并進(jìn)行精度驗(yàn)證,根據(jù)擬合精度和預(yù)測精度的R2、RMSE選擇精度最高的預(yù)測模型(表4),其中以D、Dr、NDVI705為變量的多元線性回歸模型精度最高。植被指數(shù)類參數(shù)中,以NDVI2為變量的模型精度更高;高光譜特征位置參數(shù)中,以D為變量的模型精度更高。幾個模型中,多元線性回歸模型精度較高。

表 3 植被指數(shù)變量和高光譜特征參數(shù)變量篩選結(jié)果Table 3 The screening results of vegetation index and hyperspectral characteristic parameter

表 4 危害時期監(jiān)測模型(n=84)及檢驗(yàn)(n=36)Table 4 Monitoring models(n = 84) and validation(n = 36) of damage time

基于逐步回歸分析(b)得到的4個參數(shù)(D、Dr、NDVI705、(SDr-SDb)/(SDr+SDb))建立危害時期的多元線性回歸模型(表4),其模型精度大于以D、Dr、NDVI705為變量的多元線性回歸模型。對上述4個參數(shù)進(jìn)行建模重要性分析,重要性依次為 0.129、 0.114、 0.575、 0.183, 結(jié) 果 表 明 :(SDr-SDb)/(SDr+SDb)參數(shù)雖未與松小蠹的危害時期顯著相關(guān),但仍具有一定重要性,將其納入模型變量中可提高模型精度。

2.2 松小蠹不同危害時期光合變化分析

2.2.1 松小蠹不同危害時期針葉光合參數(shù)變化 主要光合指標(biāo)見表5,健康時期、危害早期、危害中期針葉Photo 分別為4.621、2.089、0.643 μmol·m-2s-1,健康時期是危害中期的7倍,受松小蠹危害后,云南松針葉Photo逐漸下降。方差分析及Duncan’s多重比較顯示:松小蠹不同危害時期的Photo、Cond、Ci及Ls存在顯著差異(P<0.05),健康時期、危害早期、危害中期的值依次下降;相較于健康時期針葉,危害早期和危害中期Trmmol明顯下降,存在顯著差異;危害中期針葉WUE顯著低于健康時期和危害早期。結(jié)合松小蠹不同時期針葉Ci和Ls分析,表明:針葉Photo降低是松小蠹危害云南松導(dǎo)致葉肉細(xì)胞中與光合作用相關(guān)的酶活性下降引起,氣孔因素并非主要因素,其下降是為了減少針葉蒸騰引起的水分損失。

2.2.2 松小蠹不同危害時期針葉主要光合參數(shù)的相關(guān)性分析 不同危害時期針葉Photo與其他光合參數(shù)的相關(guān)性發(fā)生了不同程度變化(表6)。簡單相關(guān)分析表明:3個不同危害時期針葉Photo與Cond、Ci、Trmmol均呈不同相關(guān)性,Photo與Ci呈顯著負(fù)相關(guān)(r健康時期=-0.989**,r危害早期=-0.476**,r危害中期=-0.508**),而 Photo與其他因子在不同危害時期呈不一致相關(guān)性和顯著性。偏相關(guān)分析表明:健康時期、危害早期針葉Photo與Ci均呈極顯著負(fù)相關(guān)(r健康=-0.732**,r危害早期=-0.730**),而危害中期呈正相關(guān)(r危害中期=0.069);健康時期針葉Photo與Cond呈顯著負(fù)相關(guān)(r健康=-0.410*),危害早期呈極顯著負(fù)相關(guān)(r危害早期=-0.520**),而危害中期呈正相關(guān)(r危害中期=0.205);針葉Photo與Trmmol在3個不同危害時期呈不一致相關(guān)性和顯著性。說明云南松受到松小蠹危害后針葉凈光合速率變化是多個因子共同作用的結(jié)果。

表 5 云南松受松小蠹危害的不同時期的針葉氣體交換參數(shù)的影響Table 5 The effects of gas exchange parameters in different periods of Pinus yunnanensis damaged by Tomicus spp.

表 6 針葉主要光合參數(shù)的相關(guān)性分析Table 6 The correlation analysis of main photosynthetic parameters of needles

對松小蠹危害時期和針葉主要光合參數(shù)進(jìn)行簡單相關(guān)分析、偏相關(guān)分析及復(fù)相關(guān)分析(表7)。簡單相關(guān)分析表明:危害時期與Photo、Cond、Trmmol呈極顯著正相關(guān)(rPhoto=0.943**,rCond=-0.797**,rTrmmol=0.688**),與Ci呈極顯著負(fù)正相關(guān)(rCi=-0.764**)。偏相關(guān)分析表明:危害時期與Photo呈極顯著正相關(guān)(r=0.588**),與其他光合參數(shù)相關(guān)性均未達(dá)到顯著性水平;復(fù)相關(guān)分析表明危害時期與光合參數(shù)的相關(guān)性達(dá)到極顯著水平(r=0.909**)。

2.2.3 監(jiān)測模型的建立與精度驗(yàn)證 分別建立危害時期與Photo、Cond、Ci及Trmmol的預(yù)測模型,擬合和預(yù)測精度的R2、RMSE選擇精度最高的預(yù)測模型(表8):危害時期與Photo、Cond及Trmmol的三次函數(shù)模型,以及與Ci的線性函數(shù)模型。其中危害時期與Photo的三次函數(shù)模型的精度最高(擬合:R2=0.931,RMSE=0.204;預(yù)測:R2=0.961,RMSE=0.240)。

另根據(jù)危害時期與Photo、Cond、Ci及Trmmol的復(fù)相關(guān)分析結(jié)果,建立多元線性回歸模型(擬合:R2=0.918,RMSE=0.236;預(yù) 測:R2=0.930,RMSE=0.272),能較好的模擬預(yù)測松小蠹的危害時期。此外,還對Photo、Cond、Ci、Trmmol在建模分析中進(jìn)行了重要性分析,前三個參數(shù)重要依次為0.707、0.247、0.047,而Trmmol在模型中的可以忽略,因此建立基于Photo、Cond、Ci的多元線性回歸模型(擬合:R2=0.916,RMSE=0.225;預(yù)測:R2=0.911,RMSE=0.261)。

表 7 松小蠹危害時期和針葉主要光合參數(shù)的相關(guān)性分析Table 7 The correlation between damage time of Tomicus spp.andmain photosynthetic parameters

表 8 危害時期監(jiān)測模型(n=63)及檢驗(yàn)(n=27)Table 8 Monitoring models(n = 63) and validation(n = 27) of damage time

2.3 云南松針葉光譜指數(shù)與光合參數(shù)的相關(guān)分析

Hernandez[18]在針葉樹木上研究發(fā)現(xiàn)PRI與氣孔導(dǎo)度顯著相關(guān);NDVI與SR均能與Photo建立一定的相關(guān)性。結(jié)合本研究,對植被指數(shù)、高光譜特征參數(shù)與Photo、Cond、Ci及Trmmol進(jìn)行相關(guān)性分析和逐步回歸分析,篩選顯著或極顯著變量(表9),分別建立植被指數(shù)高光譜特征參數(shù)的擬合方程,并進(jìn)行精度驗(yàn)證,根據(jù)擬合和預(yù)測精度的R2、RMSE選擇精度最高的擬合方程(表10)。其中以參數(shù)NDVI705為變量的Photo擬合方程精度最高,擬合和預(yù)測R2均大于0.8。對Ci植被指數(shù)擬合方程的NDVI705、HI及高光譜特征參數(shù)擬合方程的λr、SDb、SDnir/SDb進(jìn)行重要性分析,結(jié)果表明:NDVI705和HI的重要性分別為0.825、0.125;λr、SDb、SDnir/SDb的重要性依次為 0.370、0.419、0.211。HI、SDnir/SDb雖與Ci未呈現(xiàn)顯著相關(guān),但該變量對其預(yù)測具有一定貢獻(xiàn),將該變量納入方程中,提高方程擬合效果。

表 9 特征變量篩選結(jié)果Table 9 The screening results of characteristic parameter

表 10 光合參數(shù)估算方程(n=63)及驗(yàn)證(n=27)Table 10 Estimation equation(n = 63) and validation(n = 27) of photosynthetic parameters

3 討論

傳統(tǒng)森林病蟲害監(jiān)測費(fèi)時費(fèi)力,誤差大。近年來,高光譜遙感在農(nóng)林業(yè)上得到廣泛應(yīng)用,可快速獲取植被、作物的有效信息,植被光譜及一階微分信息能有效反映病蟲害危害情況,光譜各波段組合可以估測植被生理特性。伍南[23]等在對油茶炭蛆病研究時發(fā)現(xiàn)光譜一階微分信息可反映炭疽病病情指數(shù),在496 nm和680 nm附近存在波峰,547 nm附近存在波谷,且隨著病情指數(shù)增大,波谷逐漸消失,峰值逐漸降低,“紅邊”出現(xiàn)“藍(lán)移”。本研究發(fā)現(xiàn),隨著松小蠹危害時期變化,云南松枝梢針葉光譜反射率和光譜一階微分曲線總體降低,光譜曲線在 510~560 nm、640~680 nm 和 740~1 036 nm 區(qū)域變化明顯;在660~740 nm,光譜一階微分曲線峰值逐漸向短波方向移動并逐漸降低,出現(xiàn)“紅邊”“藍(lán)移”現(xiàn)象。表明云南松枝梢針葉光譜和光譜一階微分信息變化可以反映松小蠹危害時期變化。隨著危害時期變化,云南松枝梢針葉顏色、光譜信息變化明顯。

以往研究表明:植被光譜指數(shù)與光合參數(shù)存在一定關(guān)系,可用線性或非線性表示,能反映植被健康狀況。本研究對8個植被指數(shù)和23個高光譜特征參數(shù)與Photo、Cond、Ci、Trmmol進(jìn)行相關(guān)分析,Photo與NDVI705和λr、Cond與SR2和SDr/SDb、Trmmol與SR1和SDr/SDb均呈極顯著正相關(guān),Ci與NDVI705和λr呈極顯著負(fù)相關(guān);以NDVI705為參數(shù)的Photo三次函數(shù)方程、以SR2為參數(shù)的Cond三次函數(shù)方程、以NDVI705和HI為參數(shù)的Ci多元線性回歸方程及以SDr、SDr/SDb為參數(shù)的Trmmol多元線性回歸方程精度最高。

4 結(jié)論

研究松小蠹不同危害時期的針葉光譜與光合特征,基于光譜特征及光合參數(shù)建立松小蠹危害時期監(jiān)測模型。針葉光譜曲線在510~560 nm、640~680 nm和740~1 036 nm區(qū)域變化明顯;危害時間越久,光譜一階微分曲線在660~740 nm出現(xiàn)“紅邊”“藍(lán)移”,且峰值逐漸降低;以D、Dr、NDVI705、(SDr-SDb)/(SDr+SDb)為變量的松小蠹危害時期多元線性回歸模型擬合效果最優(yōu)(擬合R2=0.969,預(yù)測R2=0.946、RMSE=0.212)。松小蠹危害時期與Photo密切相關(guān)(r簡單相關(guān)分析=0.943**,r偏相關(guān)分析=0.588**),其三次函數(shù)模型擬合效果較好(擬合R2=0.931,預(yù)測R2=0.961,RMSE=0.240)。基于光譜指數(shù)建立了光合參數(shù)的擬合方程,Photo的擬合效果最優(yōu)。研究結(jié)果對于利用成像高光譜技術(shù)監(jiān)測森林病蟲害具有重要意義,為今后“天空地”一體化監(jiān)測奠定了一定基礎(chǔ)。光合參數(shù)能反應(yīng)植被的生長狀況,可作為森林健康監(jiān)測的一個指標(biāo),建立高光譜特征參數(shù)與光合參數(shù)之間的模型簡單可行,可預(yù)估森林植被健康狀況,有利于森林健康監(jiān)測。

致謝:感謝石林縣林業(yè)局森防站在實(shí)地調(diào)查中給予努力支持,作者在此表示衷心的感謝。

猜你喜歡
云南松針葉反射率
影響Mini LED板油墨層反射率的因素
近岸水體異源遙感反射率產(chǎn)品的融合方法研究
海洋通報(2022年4期)2022-10-10 07:40:26
具有顏色恒常性的光譜反射率重建
更 正
云南松正交膠合木膠合工藝探索?
風(fēng)雨海棠
云南松天然更新幼苗和幼樹期的生長動態(tài)分析
云南松及其變種栽培與適應(yīng)性研究進(jìn)展
化學(xué)腐蝕硅表面結(jié)構(gòu)反射率影響因素的研究*
電子器件(2017年2期)2017-04-25 08:58:37
廣西細(xì)葉云南松森林資源動態(tài)分析
石阡县| 宁阳县| 五华县| 始兴县| 桦南县| 保德县| 芒康县| 龙岩市| 墨脱县| 无为县| 佛学| 平潭县| 田东县| 白朗县| 五台县| 咸宁市| 施甸县| 兴城市| 康平县| 大庆市| 肥乡县| 织金县| 灵山县| 武宣县| 石首市| 宁波市| 抚宁县| 修文县| 谷城县| 乌鲁木齐县| 北京市| 南宁市| 彰武县| 社旗县| 鸡泽县| 江西省| 新和县| 如东县| 本溪| 大厂| 双柏县|