氣候環境變化對冬蟲夏草產量的影響與蟲草產量預測模型研究

朱壽東　黃璐琦　郭蘭萍　馬興田　郝慶秀　樂智勇　張小波　楊光　張燕　陳美蘭

[摘要] 冬蟲夏草是我國特有的名貴中藥材，其道地產區是青藏高原那曲、昌都，青海玉樹等地區。最近10年，蟲草資源在盲目和過度的采挖中已經到達臨界點，如何加快保護、發展和合理利用冬蟲夏草這一寶貴資源，已提到了冬蟲夏草研究和利用的重要議事日程。青藏高原生態環境和氣候變化比其他地區具有超前變化趨勢，這意味著冬蟲夏草種群的分布和演化與其他種群比起來更明顯和劇烈，該文基于RS（remote sensing，遙感）/GIS（geographic information system，地理信息系統）技術，利用雪線高程，采收期（4，5月）的平均氣溫，降水量和日照時長等環境因子和蟲草產量之間的關系，建立加權幾何平均模型，在每年的6月上旬預測當年的蟲草產量，經過檢驗，該模型的精度可達82.16%以上。該研究可為蟲草產業的健康、可持續發展提供基礎數據和信息。

[關鍵詞] 冬蟲夏草； 遙感； 環境因子； 加權幾何平均； 產量預測

Climate change impacts on yield of Cordyceps sinensis and

research on yield prediction model of C. sinensis

ZHU Shoudong1， HUANG Luqi1*， GUO Lanping1， MA Xingtian2， HAO Qingxiu1，

LE Zhiyong2， ZHANG Xiaobo1， YANG Guang1， ZHANG Yan1， CHEN Meilan1

（1.State Key Laboratory Breeding Base of Daodi Herbs， National Resource Center for Chinese

Materia Medica， China Academy of Chinese Medical Sciences， Beijng 100700， China；

2.Kangmei Pharmaceutical Co.， Ltd.， Shenzhen 518000， China）

[Abstract] Cordyceps sinensis is a Chinese unique precious herbal material， its genuine producing areas covering Naqu， Changdu in Qinghai Tibet Plateau， Yushu in Qinghai province and other regions. In recent 10 years， C. sinensis resources is decreasing as a result of the blindly and excessively perennial dug. How to rationally protect， develop and utilize of the valuable resources of C. sinensis has been referred to an important field of research on C. sinensis. The ecological environment and climate change trend of Qinghai Tibet plateau happens prior to other regions， which means that the distribution and evolution of C. sinensis are more obvious and intense than those of the other populations. Based on RS （remote sensing）/GIS（geographic information system） technology， this paper utilized the relationship between the snowline elevation， the average temperature， precipitation and sunshine hours in harvest period （April and may） of C. sinensis and the actual production of C. sinensis to establish a weighted geometric mean model. The model′s prediction accuracy can reach 82.16% at least in forecasting C. sinensis year yield in Naqu area in every early June. This study can provide basic datum and information for supporting the C. sinensis industry healthful， sustainable development.

[Key words] Cordyceps sinensis； remote sensing； environmental factor； weighted geometric mean； yield prediction

目前世界上已知的廣義蟲草屬有500多種，其中冬蟲夏草是我國特有的名貴中藥資源，是由冬蟲夏草菌Cordyceps sinensis（Berk.） Sacc.侵染鱗翅目蝙蝠蛾幼蟲后發育而成的由子座和菌核組成的菌體。冬草夏草一般生長在海拔3 000 m以上的山地陰坡、半陰坡的灌叢和草甸中，在我國主要分布于青藏高原及邊緣地區[1]，零星分布于境外的尼泊爾、不丹的喜馬拉雅山脈高寒草甸，我國冬蟲夏草的采挖量占世界總采挖量的98%以上[2]。

由于天然冬蟲夏草的生長有嚴格的寄生性及要求特殊的生態地理環境，野生資源十分有限[3]。另一方面，目前冬蟲夏草不能進行規模化人工培育[3]，只能自然生長，于是其蘊藏量在日益趨高的價格逼迫下迅速下降。據統計，20世紀50年代以前全國產量曾達到1 000 t以上，60年代初為5 080 t，自70年代以來，因全球氣候變化、生態環境破壞，以及掠奪式的采挖，加之其自身生長緩慢，自然資源更新能力低下，更使天然蘊藏量迅速減少[4]。據80年代中期第三次全國中草藥資源普查結果顯示其野生蘊藏量僅約400 t，90年代就僅有300 t[5]，總體還呈逐年下降趨勢。目前冬蟲夏草天然資源已瀕于滅絕，被列為國家二級重點保護野生植物。

由于各種原因從第三次中藥資源普查至今，沒有對蟲草進行過全面的調查，特別是蟲草在西藏的具體分布和蘊藏量，有關部門只能從各鄉群眾采挖范圍上獲取一些基本情況，基本上是經驗數據。因此，準確掌握蟲草的產量信息是可持續開發冬蟲夏草資源的基礎，也是一個亟待解決的問題。

楊大榮等[6]通過研究表明，冬蟲夏草的分布受到地形、地貌、海拔、氣候等眾多因子的影響。吳慶貴等[7]通過樣方法指出植被蓋度、植物數量、坡度、天敵因子、溫度和濕度是影響冬蟲夏草種群分布的主導因子。李芬等[8]通過文獻調研和野外調查相結合的方法選擇海拔高度、植被和土壤作為影響冬蟲夏草生長的環境因子。李暉等[9]指出西藏地區冬蟲夏草分布于年均溫在[-1.9 ℃， 4.7 ℃]，年降水量在[370 mm， 850 mm]。綜合分析發現，冬蟲夏草的分布和氣候環境變化關系較大，因此本文選擇雪線高程，采收期（4月和5月）的平均氣溫，降水量和日照時長等環境因子，分別建立這些因子與蟲草產量的線性關系，利用獨立性權數法確定系數（相關系數比值作為權重），建立加權幾何平均模型，來預測當年的蟲草產量。

1 材料

1.1 實地調查與采樣

本文在2013—2014年通過實地考察、問卷調查等方法對眾多環境因子進行篩選，結合冬蟲夏草的生長周期，根據經驗初步選出蟲草所在地區的5月份雪線高程，4，5，6月的降水量，4，5，6月的平均氣溫，4，5，6月的日照時長共10個因子作為影響蟲草產量的因子，供后期建立預測模型使用。

1.2 降水量、平均氣溫、日照時長等相關因子

根據蟲草所在地區的經緯度，從中國氣象科學數據共享服務網（http：//cdc.nmic.cn）找到該區域內的氣象臺站，利用這些氣象臺站的平均值作為環境因子的值。按照這樣的方法得到4，5，6月的降水量，平均氣溫，日照時長，以及由此衍生的4月5月的降水量和，平均氣溫和，日照時長和。

1.3 研究區介紹

那曲地區地處西藏北部，位于青藏高原腹地，是長江、怒江、拉薩河、易貢河等大江大河的源頭。與新疆維吾爾自治區和青海省交界，東鄰昌都市，南接拉薩市、日喀則市、林芝市，西與阿里地區相連。境域面積45萬km2，地理位置約E83°55′—E95°5′， N29°55′—N36°30′。那曲地區除東部少量的半農半牧區外，基本上是一個純牧區，擁有草地面積5億多畝。那曲地區是冬蟲夏草的道地產區，并且是蟲草分布核心帶。

2 方法

2.1 構建雪線高程

本文利用NDSI指數（normalized difference snow index， 歸一化雪被指數）提取雪線[10]。NDSI是觀測冰雪定量指標，積雪有很強的可見光反射和很強的短波紅外吸收特性。積雪在0.5 μm附近有高反射率，在1.6，2.1 μm附近反射率較低。

NDSI=（Ref0.555 μm-Ref1.640 μm）/（Ref0.555 μm+Ref1.640 μm） （1）

公式（1）中Ref0.555 μm，Ref1.640 μm分別為0.555，1.64 μm處的反射率，NDSI≥0.4的區域即為冰雪覆蓋區。

2.1.1 把DN值圖像轉化成反射率圖像 根據遙感影像的來源，有2種情況，一種是針對Landsat TM5/Landsat TM7影像，另一種是針對Landsat TM8影像。

對于Landsat TM5和Landsat TM7影像：0.555，1.64 μm波段分別對應著影像的第2和第5波段。

第一步，分別計算各個波段每個像元的輻射亮度L值。

L=Lmax-LminQCALmax-QCALmin×（QCAL-QCALmax）+Lmin （2）

公式（2）中，QCAL為某一像元的DN值（digital number， 像元亮度值），即QCAL=DN。QCALmax為像元可以取的最大值，QCALmin為像元可以取的最小值。Lmax為該波段的最大輻射亮度，Lmin為該波段的最小輻射亮度，可以從遙感影像的頭文件中獲取以上參數的信息。

第二步，計算各波段反射率（反照率、反射率）ρ。

ρi=π×L×D2ESUN×cosθ（i為第i波段） （3）

公式（3）中，D為日地距離，可從影像的頭文件中獲取該參數信息。ESUN（solar exoatmospheric spectral irradiances，或ESUNI）為大氣頂層太陽輻照度可從遙感權威單位定期測定并公布的信息中獲取，見表1。

2.1.2 在生成的反射率圖像上計算NDSI 對于Landsat TM5/7，NDSI=（B2-B5）/（B2+B5）；對于Landsat TM8，NDSI=（B3-B6）/（B3+B6）。式中B2， B5， B3， B6分別代表各波段反射率圖像。

2.1.3 提取NDSI上大于等于0.4的區域 根據MODIS產品說明，認為NDSI大于等于0.4的區域被積雪覆蓋。利用ENVI軟件里的band math工具來提取積雪覆蓋的區域。

2.1.4 提取積雪區域的邊界 由于提取的積雪覆蓋的區域是柵格圖形，利用ENVI軟件轉化為矢量后是多邊形，接著利用Arctools工具把多邊形轉化為線類型的矢量圖形，至此雪線已提取完成。

2.1.5 利用空間分析提取雪線處最低高程 利用Arcmap 10軟件，加載DEM（digital elevation model，數字高程模型）數據，與雪線分布數據進行空間分析，提取雪線處的高程信息，最后去除圖像邊緣的異常值。

2.2 氣候環境因子與蟲草產量的關系

利用相關分析工具，分別確定各因子與蟲草產量之間的相關系數，按照相關系數大小排序來確定入選模型的因子，各因子與蟲草產量的相關系數，就是該因子的權重。確定了因子后，基于時間序列的蟲草產量數據和各因子數據建立蟲草產量的單因子線性回歸模型。

2.3 蟲草產量預測模型的建立

目前主流的賦權方法主要分為2類，一類是主觀賦權法，主要是以德爾菲法（Delphi），層次分析法（analytic hierarchy process， AHP），專家會議法為代表。而另一類客觀賦權法則是以主成分分析法，因子分析法，熵值法，離差、均差法等純數學分析為代表[11]。本文通過對比試驗，采用客觀賦權法中的獨立性權系數法求出各因子的權重系數，然后采用加權幾何平均模型來預測蟲草產產量。

2.3.1 獨立性權系數法原理 獨立性權系數法[3]是根據各指標與其他指標之間的共線性強弱來確定指標權重的。

設有指標項，若指標與其他指標的復相關系數越大，則說明與其他指標之間的共線性關系越強，越容易由其他指標的線性組合表示，重復信息越多，因此該指標的權重也就應該越小。

即若指標與其他指標的復相關系數R越大，該指標的權重越小。

R=（y-）（-）（y-）2（-）2

由于本文中強調各相關因子與蟲草產量的關系，因此用各因子與蟲草產量的相關系數替代了各因子間的復相關系數。

2.3.2 加權幾何平均模型 如果預測因子與蟲草產量的線性函數為fi，i=1， 2， …， N，N為最終入選模型的因子數目。各個因子權重比值為wi，則蟲草產量Y為公式（6）。

Y=（∏Nf wii）1Ni=1wi（6）

2.4 蟲草產量數據調研和采集

在2013，2014年6月，通過實地走訪那曲地區聶榮縣、比如縣、索縣和那曲蟲草交易市場，和采挖蟲草的藏民進行交流，獲取了蟲草的采挖情況。結合查閱資料獲取了近6年的那曲蟲草產量[1521]，見表2。

2.5 生態因子與蟲草產量的關系

2.5.1 雪線高程 根據蟲草的采挖時間，選擇收集那曲地區的5月份的landsat衛星影像數據，如果5月份數據含云量特別大，可選擇4月份的數據。那曲地區所在的影像條帶號是：WRS_PATH=137，WRS_ROW=038。最終選擇的數據情況見表3。

按照2.1的方法對遙感進行處理，得到雪線處的高程信息見表4，提取雪線的過程中雖然刪除了影像邊緣的異常值，但是影像零散的冰雪覆蓋區對雪線的邊界有影響，根據正態分布原理[12]，把提取出來的雪線高程的“平均值-2×標準差”作為雪線的最低高程。按照2012年、2014年、2011年、2009年、2010年的順序，雪線的最低高程依次上升，見圖1。

由于2013年的影像含云量較大，建模過程中舍去這一年的數據。利用剩余的5期數據，計算得雪線高程和蟲草產量的相關系數r=0.40，線性回歸關函數為y=0.030 1x-99.367。

2.5.2 4月和5月降水量之和 選擇那曲地區附近的5個氣象臺站（編號為55294，55299，55493，56106，56202）的降水量數據的平均值作為該因子的值。通過比較單月降水量、各月降水量之和與蟲草

產量的相關系數發現，4月和5月降水量之和比單獨月降水量與蟲草產量的相關性更大，于是選擇4月和5月降水量之和作為預測因子。那曲地區各年4月和5月降水量之和見表5。

計算得到4月和5月降水量的和與蟲草產量的相關系數r=0.95，線性回歸關函數為y=0.13x-49.44。

2.5.3 5月平均氣溫 同降水量一樣，選擇5個氣象臺站的數據平均值作為該因子的值。通過比較單月平均氣溫、各月平均氣溫之和與蟲草產量的相關系數發現，5月份的平均氣溫比其他月、月份之和與蟲草產量的相關性更大，于是選擇5月平均氣溫作為預測因子，那曲地區各年5月平均氣溫見表5。計算得到5月平均氣溫與蟲草產量的相關系數r=-0.69，線性回歸關函數為y=-1.477 1x+114.09。

2.5.4 4月和5月日照時長之和 通過比較單月日照時長、各月日照時長之和與蟲草產量的相關系數發現，4月和5月日照時長之和比其他月、月份的和與蟲草產量的相關性更大，于是選擇4月和5月日照時長之和作為預測因子，那曲地區各年4月和5月日照時長之和見表5。計算得到5月平均氣溫與蟲草產量的相關系數r=-0.50，線性回歸關函數為y=-0.022 1x+144.4。

2.6 構建蟲草產量預測模型

根據各個因子與蟲草產量的相關系數，得到相關系數的比值是r1∶r2∶r3∶r4=0.4∶0.95∶|-0.69|∶|-0.50|=1∶3∶2∶1.于是各因子權重比值w1∶w

上式中，X1是雪線高程，X2是4月和5月降水量的和，X3是5月平均氣溫，X4是4月和5月日照時長的和，Y是預測的蟲草產量（單位t）。

3 結果

3.1 模型精度驗證

把2009—2014年的環境因子數據帶入得到的模型中，得到的預測值和實際值見表6。特別地，2013年的遙感影像含云量較大，無法提取雪線高程，于是對預測模型進行修正，去掉f1，只利用剩余的3個因子進行預測，得到2013年的蟲草產量預測值。綜合過去6年的預測結果和真實值，計算得到模型預測的平均偏差是7.98%，最大偏差是17.84%。

3.2 2015年蟲草產量預測

在2015年6月上旬，根據生態因子的構建方法，得到X1雪線高程是4 574 m，X2 4月和5月降水量之和是80.72 mm，X3 5月平均氣溫是5.08 ℃，X4 4月和5月日照時長之和由于數據暫缺無法獲取，所以去掉f4，只利用前3個因子進行預測，得到2015年那曲地區的蟲草產量是46.42 t，比去年增產約60%，是2013年產量的85%左右。綜合分析來看，今年高產與今年蟲草生長季降雨量豐沛有較大關系。

4 討論

本文通過RS/GIS技術來構建氣候環境因子，把氣候環境因子與蟲草產量進行線性擬合，配合幾何加權平均模型，得到了蟲草產量的預測模型。經過檢驗，發現模型預測偏差較小。而且，如果有非重要的生態因子因為實際情況（天氣或者時間關系）暫缺時，模型的精度不會受到較大影響，具有較高的穩定性和靈活性。

[參考文獻]

[1] 章力建，李兵，胡育驕，等.中國冬蟲夏草資源管理概況[J]. 中國草地學報， 2010， 32（增刊）： 1.

[2] 劉天平，卓嘎，郭健斌. 西藏冬蟲夏草資源可持續利用探討[J]. 西藏研究， 2010（3）： 114.

[3] 蔡佩云，孫世英. 青海省冬蟲夏草資源持續開發利用淺析[J]. 中國草地學報， 2010， 32（增刊）： 6.

[4] 張古忍，吳光國，古德祥，等. 冬蟲夏草資源可持續利用面臨的問題及其對策[J]. 中國草地學報[J]， 2010， 32（增刊）： 51.

[5] 陳仕江，鐘國躍，馬開森. 冬蟲夏草資源可持續利用的思考[J]. 中國草地學報， 2010， 32（增刊）： 44.

[6] 楊大榮，彭艷瓊，陳吉岳，等. 中國冬蟲夏草分布格局與環境變化對其分布的影響[J]. 中國草地學報， 2010，32（增刊）： 22.

[7] 吳慶貴，蘇智先，蘇瑞軍，等. 冬蟲夏草生境選擇主導因子[J]. 廣西植物， 2009，29（3）： 331.

[8] 李芬，吳志豐，徐翠，等. 三江源區冬蟲夏草資源適宜性空間分布[J]. 生態學報， 2014，34（5）： 1318.

[9] 李暉，王立輝. 西藏冬蟲夏草生境的氣候[J]. 西藏科技， 2011（6）： 68.

[10] 趙英時. 遙感應用分析原理與方法[M]. 北京： 科學出版社， 2003.

[11] 李亮. 評價中權系數理論與方法比較[M]. 上海： 上海交通大學， 2009.

[12] 劉自遠，劉成福. 綜合評價中指標權重系數確定方法探討[J]. 中國衛生質量管理， 2006， 13（2）： 44.

[13] 蔡輝，沈婉蘭. 醫院綜合效益評價中的權數（二）——獨立性權數[J]. 中國醫院統計， 1997， 4（3）：131.

[14] 杜子芳. 抽樣技術及其應用[M]. 北京：清華大學出版社， 2005.

[15] 李芬，吳志豐，徐翠，等. 三江源區冬蟲夏草資源適宜性空間分布[J]. 生態學報， 2014，34（5）： 1318.

[16] 新華網. 西藏冬蟲夏草價格進入下行通道[EB/OL]. （20100705） [20160209].http：//news.xinhuanet.com/fortune/201007/05/c_12300368.htm .

[17] 西藏商城藏善堂. 新產冬蟲夏草價格2010年5月走勢分析[EB/OL]. （20100423） [20160209 ].http：//www.xizangshop.com/article138.html.

[18] 冬蟲夏草屋. 同仁堂蟲草價格調價漲幅達三成[EB/OL]. （20120602） [20160209].http：//www.dcxc5.com/article488.html.

[19] 冬蟲夏草屋. 蟲草產量連年下降，冬蟲夏草資源日趨稀缺[EB/OL]. （20120602）[20160209].http：//www.dcxc5.com/article585.html.

[20] 中國經濟網. 青海計劃今年采集冬蟲夏草35噸[EB/OL]. （20140529） [20160209].http：//city.ce.cn/ztpd/ysgjccj/ysxwlb/yshdxc/201405/29/t20140529_1653483.shtml.

[21] 中國經濟網. 2014年新草動態[EB/OL]. （20140529） [20160209].http：//city.ce.cn/ztpd/ysgjccj/ysxwlb/yshdxc/201405/29/t20140529_1653566.shtml.

[責任編輯 呂冬梅]