顏懋瑤,彭曉昶,代光輝,張志明
(云南大學 生態學與環境學院,云南 昆明 650500)
達摩麝鳳蝶(Byasadaemonius)是鳳蝶科(Epicopeiidae)麝鳳蝶屬(Byasa)的一種昆蟲,主要分布于云南、四川等地區。其整個生活史嚴重依賴一種植物——貫葉馬兜鈴(Aristolochiadelayayi),它在幼蟲階段僅僅取食貫葉馬兜鈴的葉片,吸收其中的馬兜鈴酸富集到體內,以達到抵抗螞蟻以及鳥類這些天敵的捕食和破壞的作用。貫葉馬兜鈴是達摩麝鳳蝶的特異性寄主[1]。貫葉馬兜鈴因具有和草果(Amomumtsaoko)相似的獨特氣味又名山草果。其主要產于海拔1 640~2 250 m的云南省西北部及四川省西南部的金沙江干旱河谷中,是一種極具地域特色的野生食用香料,但人工培育極難存活,生境范圍狹窄,一般僅生長于山坡裸巖或半風化的石灰巖區域,直到現在當地村民仍延續著采食野生貫葉馬兜鈴的習慣[2],這也使得貫葉馬兜鈴的生存難度增加。貫葉馬兜鈴也是金沙江干旱河谷區域的特有植物[3]。目前該種已被《世界自然保護聯盟紅色名錄》(IUCN2014)列為瀕危物種。而一直依賴貫葉馬兜鈴生存的達摩麝鳳蝶也隨著其特異性寄主植物的瀕危,被《中國紅色物種名錄》評為易危物種[4]。目前,人們對于這2個亟待保護的物種的自然資源分布狀況卻不甚了解,僅有零星的野外分布記錄,這對保護工作的開展極為不利。
物種分布模型(Species distribution models)被廣泛應用于預測物種的潛在分布,其中最大熵模型(MaxEnt)利用物種的已知存在點數據和相關的氣候因子、地形因子等環境變量作為約束條件,根據一定的算法構建生態位模型,相較于其他的物種分布模型,最大熵模型對于在野外分布點極少、實地調查困難的物種也能獲得很好的預測結果,已經被許多研究所使用[5-7]。但在之前的研究中,大多數的潛在分布預測都針對獨立的物種,在大尺度上(>200 km)僅需要考慮氣候、地形等非生物因子,然而在最近的如Wisz M S等[8]的研究發現中了解到即使大尺度的物種潛在分布預測,如果存在種間關系也是需要作為約束條件加入模型的。如李國慶等[8]研究中發現當物種間存在強相互作用時,僅考慮非生物因素進行物種分布預測時,會高估物種的潛在分布區[9]。達摩麝鳳蝶和貫葉馬兜鈴之間存在這種不可忽略的特異性寄生的關系,所以在預測時還應該考慮種間關系的影響。在僅考慮非生物因子或種間關系情況下得到的預測結果可能存在很大差異。此外,不同因子對預測結果的影響程度不同,在分析種間關系對達摩麝鳳蝶潛在分布影響時,還需計算不同因子對預測結果影響程度的大小。通徑分析(path analysis)可以探討不同因素對于物種潛在分布的影響。通徑分析是數量遺傳學家Sewall Wright于1921年提出的一種多元統計技術,當許多自變量共同影響一個因變量時,每個自變量對因變量的重要性是不同的,通徑分析可用于分析各個自變量對因變量的直接重要性和間接重要性[10]。
目前對貫葉馬兜鈴和達摩麝鳳蝶兩個保護物種的基礎研究十分匱乏,主要關注貫葉馬兜鈴的揮發物成分[11]、快速離體繁殖[12]或是其與達摩麝鳳蝶的共生機制[1]方面,尚無對2個物種自然資源分布和共生對其潛在分布影響方面的研究。在本研究中,擬解決以下關鍵科學問題:(1)貫葉馬兜鈴作為達摩麝鳳蝶的特異性寄主,它們的種間關系是否對達摩麝鳳蝶的潛在適生區預測有顯著影響?(2)非生物因子、種間關系因子對達摩麝鳳蝶潛在適生區的影響程度有何差異?通過MaxEnt在不同因子條件下設置模型對達摩麝鳳蝶的適生區進行預測,并檢驗非生物因子、種間關系因子和對達摩麝鳳蝶潛在分布的影響程度。旨在研究達摩麝鳳蝶的潛在自然資源分布情況,并探討寄主植物與環境因子及其相互作用對其潛在分布的影響。
本文2個研究物種為達摩麝鳳蝶和貫葉馬兜鈴(圖1),根據2個物種的標本記錄,確定研究區位于中國西南25°00′~31°00′N,99°00′~104°00′E的區域內。實際采樣點主要分布于麗江、迪慶地區的金沙江干旱河谷地段。這一地區由于焚風效應、局地環流等原因形成了與其水平地帶、垂直地帶完全不同的小氣候生境。年均溫21~23 ℃,年降水量600~800 mm,年均蒸發量2 700~3 800 mm,存在著很典型的稀樹灌木草叢植被群落類型,是生態學研究的熱點區域[13]。

圖1 達摩麝鳳蝶在貫葉馬兜鈴葉片上產卵Fig.1 B.daemonius were laying eggs on A.delayayi
通過實地調查,記錄2個物種所在地的GPS點作為存在數據。數據采集于大理市鶴慶縣、大理市賓川縣、香格里拉市虎跳峽鎮、德欽縣奔子欄鎮等地剔除掉距離小于1 km的存在點,最終獲得貫葉馬兜鈴的有效數據12個,達摩麝鳳蝶的有效數據9個(表1)。

表1 貫葉馬兜鈴和達摩麝鳳蝶存在點GPS數據Tab.1 GPS points of Aristolochia delayayi &Byasa daemonius
19個氣候因子從worldclim(http://www.world-clim.org)下載,海拔、坡度數據來源NASA的SRTM(http://www.jpl.nasa.gov/srtm/),分辨率為1 km×1 km,首先用 Erdas 2013 進行裁剪出所需范圍,然后再用 Arcgis 10.5的表面分析功能提取坡度、坡向的數據,最后將海拔、坡度、坡向圖層都轉換為MaxEnt所需的asc格式。
考慮到環境因子間的空間自相關,先用SPSS對19個氣候因子、3個高程數據進行相關性分析,選取|r|>0.9的環境因子,僅選擇一個與物種分布關聯緊密或便于模型解釋的用于模型預測[14]。最終保留年平均溫度(Bio 01)、年溫度季節變化(Bio 04)、最干季平均溫度(Bio 09)、最濕月降水量(Bio 13)、最冷季降水量(Bio 19 )、海拔(Alt)這6個非生物因素。
2.3.1 MaxEnt模型的構建與參數選擇
本文針對存在特異性寄生關系的貫葉馬兜鈴和達摩麝鳳蝶,設置了(1)Ma模型(僅考慮非生物因子):由達摩麝鳳蝶存在點數據和非生物因子創建的模型,即僅使用達摩麝鳳蝶存在點與氣候因子及地形因子作為約束條件進行預測的傳統潛在分布預測模型。(2)Mb模型(僅考慮種間關系):達摩麝鳳蝶存在點數據和生物因子模型,這其實需要兩個步驟,先創建Mc模型(寄主植物的潛在分布),即貫葉馬兜鈴存在點和非生物因子模型,即貫葉馬兜鈴的存在點與氣候因子及地形因子作為約束條件進行潛在分布預測,得到貫葉馬兜鈴潛在分布預測結果(即為生物因子);然后用達摩麝鳳蝶的存在點和僅使用生物因子作為約束條件進行潛在分布預測。(3)Md全模型(同時考慮非生物因子和種間關系):是使用達摩麝鳳蝶的存在點與非生物因子及生物因子一起作為約束條件的潛在分布預測模型。
為提高對于小樣本數據模型的準確性,將最大熵模型中的特征組合(feature combination,FC)選項中的“linear feature”(線型,用于約束均值)、“quadratic feature”(二次型,用于約束方差)、“product feature”(乘積型,用于約束協方差)3個參數加入模型,同時又引入“hinge features”(片段化)參數以平滑預測結果,減少訓練數據的過擬合[15]。
2.3.2 MaxEnt模型的精度驗證
(1)AUC曲線 模型精度的評估指標有受試者工作特征曲線下面積(area under the receiver operating characteristic curve,AUC)AUC值是指分類方案成功將預測正確的點與預測錯誤的點分開的概率,區間在0.5~1.0之間,0.5代表一個完全隨機的預測,1.0代表一個完美的預測,AUC值在0.5~0.7之間,表示預測精度較差,0.7~0.9之間,表示預測精度中等,>0.9則表型的預測精度非常高[16]。AUC值由MaxEnt模型輸出。
(2)赤池信息準則 赤池信息量準則(akaike information criterion)即AIC,是衡量統計模型擬合優良性的一種標準,是由日本統計學家赤池弘次創立和發展的。赤池信息量準則建立在熵的概念基礎上,可以權衡所估計模型的復雜度和此模型擬合數據的優良性[17]。AIC值越小,說明模型的復雜度和擬合度越高。AIC值通過Arcgis 10.5輸出。
2.3.3 適生區的劃分
使用Arcgis 10.5的重分類,自然段點法將潛在分布區劃分為高適生區(≥0.9)、次適生區(0.7~0.9)、低適生區(0.5~0.7)和非適生區(<0.5)[18]。
2.3.4 不同因子在達摩麝鳳蝶潛在適生區預測中的重要性——主導因子選取
刀切法Jackknife,每次從樣本集中刪除一個或幾個樣本,剩余的樣本成為“刀切”樣本,由這一系列“刀切”樣本測算變量重要性可以得到每個因子對模型構建的貢獻度[19],根據每個因子貢獻度的不同可以得到影響預測結果的主導因子。結果由MaxEnt模型輸出。
參考Hélida等[20]的方法,Md全模型可拆分為Ma、Mb和Mc等3個模型,將符合正態分布的Md全模型作為因變量,Ma、Mb、Mc作為代表因子的因變量使用逐步回歸分析法通徑分析,得到各個自變量對因變量的直接通徑系數和間接通徑系數以及決策系數(R2)來量化各自變量對因變量的綜合影響大小。結果通過SPSS 24軟件輸出。
各模型的AUC值(受試者工作特征曲線下面積)均達到0.9以上,說明各模型預測的精度都較高。AIC值(赤池信息量準則)Ma>Mc>Md>Mb,即考慮了種間關系的2個模型擬合結果表現較好,其中Mb模型的擬合度最高(表 2)。

表2 模型精度Tab.2 Model accuracy
得到的不同模型的潛在分布(圖2 a-d)。

圖2 不同模型潛在分布預測結果注:a為Ma模型,b為Mb模型,c為Mc模型,d為Md全模型。Fig.2 The results of potential distribution prediction graphs for different models
通過對預測結果的對比可知,高適宜區面積:Ma>Mc>Mb>Md;中適宜區面積:Ma>Mc>Mb>Md(圖2,表3)。考慮了種間關系的模型得到的適生區面積遠小于僅考慮非生物因子的模型。

表3 中高適宜潛在分布區面積Tab.3 Distribution area of medium and high potential
刀切法結果顯示(圖3),寄主植物貫葉馬兜鈴的種間關系是貢獻率最高的影響因子,且其貢獻率顯著高于其他非生物因子。

圖3 基于刀切法的不同環境因子重要性注:深藍色為因子單獨建模時的貢獻率,其長度越長貢獻率越高;淺綠色為除該因子外其他環境因子貢獻率的總和。Fig.3 Importance analysis of environmental factors by Jackknife method
Md模型的正態性檢驗結果表明,Md全模型服從正態分布(P<0.01)(表4),可以進行通徑分析。

表4 正態性檢驗結果Tab.4 Normality test results
結果表明,Ma、Mb、Mc對Md模型的直接通徑系數分別為0.215、0.04、0.578,決策系數R2分別為0.229、0.053、0.546。Mc模型的直接通徑系數和決策系數均為最高,其通過Ma和Mb的間接通徑系數也最高。Mc模型,即種間關系因子對Md全模型的影響最大,種間關系是對于達摩麝鳳蝶潛在適生區分布最重要的決定性因素(表5)。

表5 Md全模型與各自變量通徑分析Tab.5 Md model and path analysis of independent variables
通過對達摩麝鳳蝶的潛在分布預測研究可知,4種模型的AUC值都很高,預測結果都顯示高適生區在滇西北、川西和藏東南的干熱和干暖河谷地區,與實際分布情況[21]吻合,預測結果都是具有可信度的,但在反映物種實際分布的能力上存在差異。考慮了種間關系的Mb和Md模型擬合度更高,說明種間關系的加入提高了模型的識別能力和精確度,減少了遺漏誤差。這一結果在適生區的實際地理分布上也有反映。貫葉馬兜鈴是金沙江干旱河谷的特有種。Ma模型輸出的潛在適生區除了金沙江干旱河谷外,還包含了峽谷之上的高山地帶。而Mb模型輸出的潛在適生區范圍明顯遵循金沙江干旱河谷分布范圍,更加符合達摩麝鳳蝶的實際生境。比較Mc與Mb模型的預測結果,發現貫葉馬兜鈴的潛在分布區范圍大于且包含達摩麝鳳蝶的潛在分布區,符合達摩麝鳳蝶與貫葉馬兜鈴特異性寄生的種間關系,也表明僅考慮非生物因子的Ma模型對存在特異性寄生關系這種種間關系的物種潛在分布范圍的預測結果可能是過大的,與此前的研究結果吻合[8]。結合通徑分析結果,證明種間關系不僅是影響達摩麝鳳蝶潛在分布最為重要的影響因素,還提高了模型預測的準確性。由于無法獲取當地更精確的氣候數據和更多的物種分布數據,本文的研究空間尺度無法再縮小。因此,種間關系在更小尺度上對預測結果的影響表現以及在影響預測結果的主要因子的全面性等問題上還需要進一步的探討。
Zurell等[22]的研究認為,物種分布模型在對存在種間競爭關系和互利共生關系物種潛在適生區預測時可以達到較好的模擬效果,但是對捕食-被捕食關系的預測結果可能是不理想的。但本文的研究結果與其相反,考慮了寄主植物的種間關系后,達摩麝鳳蝶的潛在適生區更加符合其真實分布。這可能是由于,與一般的捕食-被捕食關系不同,達摩麝鳳蝶食源單一,高度特化,生活史嚴重依賴于貫葉馬兜鈴的存在,種間關系不可避免地成為限制達摩麝鳳蝶存在的關鍵因子。種間關系是多種多樣的,本文只討論了其中一種。不同的種間關系對物種的影響模式不同,將種間關系作為生物因子加入預測時,應考慮種間關系是否作為影響物種分布的限制因子。
一些研究認為,加入種間關系等生物因子會和氣候因子等非生物因子產生相關性上的混淆[23],本文采用通徑分析的方法,把生物因子對物種潛在分布的直接和間接影響計算出來,避免了生物因子與非生物因子之間的共線性混淆。
通過實地調查的情況了解到,達摩麝鳳蝶生存地區氣候環境惡劣,其生境破碎化程度嚴重。寄主植物貫葉馬兜鈴作為當地傳統的香料植物,長期遭受到人為過度采摘,其生存空間持續縮小,這也直接影響到達摩麝鳳蝶的生存,本文研究結果可作為就地保護的依據,再結合指導當地群眾適度利用、科學管理,以達到可持續發展的目標。