田 昊,張葛祥*,榮海娜,Mario J.P REZ-JIM NEZ,Luis VALENCIA-CABRERA,陳 鵬,侯 蓉,齊敦武
(1.西南交通大學(xué)電氣工程學(xué)院,成都610031; 2.塞維利亞大學(xué)計(jì)算機(jī)科學(xué)與人工智能系,塞維利亞41012,西班牙;3.成都大熊貓繁育研究基地,成都610081; 4.四川省瀕危野生動(dòng)物保護(hù)生物學(xué)重點(diǎn)實(shí)驗(yàn)室,成都610081)
(*通信作者電子郵箱zhgxdylan@126.com)
大熊貓是世界自然基金會(huì)的形象大使和當(dāng)之無(wú)愧的世界生物多樣性保護(hù)的旗艦物種[1-6]。在研究人員從事與大熊貓活動(dòng)相關(guān)的工作中,需要對(duì)大熊貓的生活習(xí)性、種群數(shù)量、種群年齡結(jié)構(gòu)組成以及引起種群數(shù)量變化的各種因素進(jìn)行分析和預(yù)測(cè)[7-9]。這些研究對(duì)大熊貓的保護(hù)工作具有重要意義[10-11]。
在生態(tài)系統(tǒng)種群動(dòng)態(tài)性建模中最廣泛使用的方法是常微分方程法和偏微分方程法[12-13],相關(guān)科研人員曾廣泛使用萊斯利矩陣對(duì)不同大熊貓種群進(jìn)行動(dòng)態(tài)預(yù)測(cè)[14-15]。常微分?jǐn)?shù)學(xué)方程通過(guò)簡(jiǎn)化的數(shù)學(xué)分析很容易用簡(jiǎn)單的公式描述復(fù)雜系統(tǒng),但該方法無(wú)法捕捉到空間動(dòng)態(tài)性和隨機(jī)效應(yīng)。偏微分方程建模方法雖克服了常微分方程法建模方法中不能捕獲生態(tài)系統(tǒng)空間動(dòng)態(tài)性的缺點(diǎn),但建模過(guò)程中存在非常復(fù)雜的數(shù)學(xué)分析。隨著生態(tài)系統(tǒng)復(fù)雜性不斷增加,采用常微分方程與偏微分方程方法建立的模型會(huì)變得更加難以解決問(wèn)題,以至于這種建模框架在許多情況下不適用。
目前,關(guān)于大熊貓種群動(dòng)態(tài)的研究文獻(xiàn)還比較少。文獻(xiàn)[16-17]利用計(jì)算機(jī)仿真分析了大熊貓種群生存力;文獻(xiàn)[18-19]分析討論了大熊貓種群動(dòng)態(tài)與它們主要食物來(lái)源竹子之間的相互影響作用;文獻(xiàn)[20-21]通過(guò)計(jì)算機(jī)仿真研究了棲息地對(duì)大熊貓種群動(dòng)態(tài)的影響;文獻(xiàn)[22-23]研究了汶川地震對(duì)大熊貓種群的影響;在文獻(xiàn)[24-25]中,相關(guān)研究人員分析了氣候變化對(duì)大熊貓種群及其保護(hù)政策的威脅和挑戰(zhàn)。這些研究都無(wú)法很好地反映大熊貓種群的空間動(dòng)態(tài)性和隨機(jī)效應(yīng),國(guó)外研究人員已經(jīng)將膜計(jì)算數(shù)據(jù)建模應(yīng)用于生態(tài)系統(tǒng)種群動(dòng)態(tài)性的研究中[26-31]。在國(guó)內(nèi),本文的實(shí)驗(yàn)室團(tuán)隊(duì)已經(jīng)使用概率膜系統(tǒng)對(duì)成都大熊貓繁育研究基地的圈養(yǎng)大熊貓生態(tài)系統(tǒng)進(jìn)行建模研究[32-33],該模型粗略地考慮了大熊貓生命周期中的繁殖、死亡、進(jìn)食、救護(hù)行為,繁殖模塊設(shè)為每年固定數(shù)量的新生大熊貓,并且假定大熊貓為單胞胎繁殖,救護(hù)模塊中一年實(shí)際救護(hù)數(shù)量最多為1只個(gè)體。根據(jù)對(duì)譜系數(shù)據(jù)更深入的研究發(fā)現(xiàn)大熊貓雙胞胎現(xiàn)象較為普遍,并且2005—2008年大熊貓的出生率較之后的年份明顯高出很多,查證大熊貓繁殖期后設(shè)置合理的大熊貓繁殖年齡,根據(jù)處于繁殖期的雌性大熊貓數(shù)量作為主要因素之一設(shè)計(jì)進(jìn)化規(guī)則。救護(hù)模塊中加入同一年內(nèi)出現(xiàn)兩只個(gè)體的可能性,最為重要的是隨著野放現(xiàn)象的出現(xiàn),在計(jì)算模型中需要考慮加入野放模塊。而且根據(jù)概率膜系統(tǒng)的定義,其無(wú)法對(duì)兩個(gè)大熊貓種群之間的相互作用與交流進(jìn)行數(shù)據(jù)建模,如若后期需要對(duì)多個(gè)種群進(jìn)行建模研究,此時(shí)概率膜系統(tǒng)則無(wú)法實(shí)現(xiàn)。
在此提出使用種群動(dòng)態(tài)膜系統(tǒng)對(duì)中國(guó)大熊貓保護(hù)研究中心的圈養(yǎng)大熊貓進(jìn)行數(shù)據(jù)建模,新模型從個(gè)體行為上來(lái)模擬中國(guó)大熊貓保護(hù)研究中心圈養(yǎng)大熊貓的生態(tài)系統(tǒng),更詳細(xì)地分析討論了2005—2016年的大熊貓譜系數(shù)據(jù),充分考慮了可能會(huì)影響該種群動(dòng)態(tài)性的相關(guān)因素。該模型更完整地體現(xiàn)出中國(guó)大熊貓保護(hù)研究中心圈養(yǎng)大熊貓生態(tài)系統(tǒng)的特點(diǎn),更加符合該生態(tài)系統(tǒng)的實(shí)際情況,同時(shí)適用于接下來(lái)對(duì)中國(guó)大熊貓保護(hù)研究中心圈養(yǎng)大熊貓生態(tài)系統(tǒng)的研究。種群動(dòng)態(tài)膜系統(tǒng)(Population Dynamics P system,PDP system)比微分方程法更復(fù)雜,但是使用的語(yǔ)言更接近專(zhuān)家語(yǔ)言(生態(tài)學(xué)家)。通俗地說(shuō),用PDP系統(tǒng)建模方法通過(guò)對(duì)象集和一些規(guī)則來(lái)控制它們的進(jìn)化和相互影響來(lái)模擬大熊貓種群現(xiàn)實(shí)生活中的進(jìn)化過(guò)程。在復(fù)雜系統(tǒng)建模中對(duì)象是最重要的組成部分(大熊貓、空間資源、食物等),它們同時(shí)可以根據(jù)一系列的重寫(xiě)規(guī)則來(lái)描述系統(tǒng)的動(dòng)態(tài)性。每一個(gè)相關(guān)規(guī)則采用概率函數(shù)試圖捕捉自然中固有的隨機(jī)性。
本文以中國(guó)動(dòng)物園協(xié)會(huì)已發(fā)布的全球圈養(yǎng)大熊貓譜系數(shù)據(jù)為依據(jù),從中搜集統(tǒng)計(jì)中國(guó)大熊貓保護(hù)研究中心大熊貓種群數(shù)據(jù),并以其作為研究對(duì)象。根據(jù)對(duì)大熊貓譜系數(shù)據(jù)中年齡結(jié)構(gòu)、繁殖率和死亡率的分析考量,以及與大熊貓專(zhuān)家研討得到的相關(guān)意見(jiàn),在大熊貓從出生到死亡的整個(gè)生命過(guò)程中,依照其不同時(shí)期的行為特點(diǎn)將其分為幼年期、亞成年期、青年期、中年期、中老年期和老年期共6個(gè)年齡階段,并確定繁殖期所處的年齡范圍。以6個(gè)年齡階段為基礎(chǔ),追蹤每只大熊貓的行為,統(tǒng)計(jì)確定大熊貓?jiān)诟麟A段的食物種類(lèi)、進(jìn)食數(shù)量、死亡率以及在繁殖期內(nèi)的大熊貓繁殖率。
同時(shí),由于人為因素的干擾,在此需要考慮野外大熊貓救護(hù)行為對(duì)大熊貓種群動(dòng)態(tài)性產(chǎn)生影響的因素,根據(jù)往年救護(hù)大熊貓的數(shù)量、年齡和性別的統(tǒng)計(jì)和分析,確定每年可能救護(hù)的大熊貓數(shù)量和年齡分布。野外的雌性大熊貓由于生存環(huán)境等原因無(wú)法同時(shí)撫養(yǎng)雙胞胎幼體,之前都是考慮大熊貓是單胞胎的哺乳動(dòng)物,但在圈養(yǎng)大熊貓生態(tài)系統(tǒng)中隨著生存環(huán)境的改善和人工繁殖技術(shù)的發(fā)展,雙胞胎個(gè)體都能得到良好的成長(zhǎng)環(huán)境,故在此需要把雙胞胎的可能性納入模型。通過(guò)與中國(guó)大熊貓保護(hù)研究中心專(zhuān)家討論研究發(fā)現(xiàn),不同時(shí)間段的大熊貓出生率有較大差別,2005—2008年出生率較高,之后隨著大熊貓種群數(shù)量的增長(zhǎng)、年齡結(jié)構(gòu)的穩(wěn)定,出生率也隨之下降并穩(wěn)定下來(lái),故本文首次在種群動(dòng)態(tài)膜系統(tǒng)的計(jì)算模型中加入分段式的進(jìn)化規(guī)則,對(duì)2008年前后兩個(gè)時(shí)期的大熊貓種群參數(shù)進(jìn)行不同的參數(shù)設(shè)置,使其更符合中國(guó)大熊貓保護(hù)研究中心圈養(yǎng)大熊貓生態(tài)系統(tǒng)的種群特點(diǎn)。圈養(yǎng)大熊貓除了保護(hù)種群這個(gè)重要原因,更是希望通過(guò)野放讓大熊貓回歸自然,復(fù)壯野外小種群。故隨著野化技術(shù)的成熟,從2012年開(kāi)始選擇合適的大熊貓進(jìn)行野化放歸。
模型中所需的參數(shù)如下所示:
i表示圈養(yǎng)大熊貓的性別,其中i=1表示雄性大熊貓,i=2表示雌性大熊貓。
ki,1表示圈養(yǎng)大熊貓年齡達(dá)到亞成年階段;ki,2表示圈養(yǎng)大熊貓年齡達(dá)到成年階段;ki,3表示圈養(yǎng)大熊貓年齡達(dá)到中年階段;ki,4,1表示圈養(yǎng)大熊貓年齡達(dá)到中老年階段;ki,4,2表示圈養(yǎng)大熊貓年齡達(dá)到老年階段;ki,5表示圈養(yǎng)大熊貓的最大年齡;ki,6表示圈養(yǎng)大熊貓?zhí)幱谟啄觌A段的死亡率;ki,7表示圈養(yǎng)大熊貓?zhí)幱趤喅赡觌A段的死亡率;ki,8表示圈養(yǎng)大熊貓?zhí)幱诔赡觌A段的死亡率;ki,9表示圈養(yǎng)大熊貓?zhí)幱谥心觌A段的死亡率;ki,10表示圈養(yǎng)大熊貓?zhí)幱谥欣夏觌A段的死亡率;ki,11表示圈養(yǎng)大熊貓?zhí)幱诶夏觌A段的死亡率;ki,12表示圈養(yǎng)大熊貓進(jìn)入繁殖期的年齡;ki,13表示圈養(yǎng)大熊貓結(jié)束繁殖期的年齡。
g1表示繁殖期大熊貓生育單胞胎的概率;g2表示繁殖期大熊貓生育雙胞胎的概率;g3表示一年中供給圈養(yǎng)大熊貓的竹筍量(kg);g4表示一年中供給圈養(yǎng)大熊貓的竹子量(kg);g5表示一年中供給圈養(yǎng)大熊貓的其他食物(例如水果、牛奶等)量(kg)。
fi,1表示每只幼年大熊貓一年中需要消耗的竹筍總量;fi,2表示每只幼年大熊貓一年中需要消耗的竹子總量;fi,3表示每只幼年大熊貓一年中需要消耗的其他食物總量;fi,4表示每只亞成年大熊貓一年中需要消耗的竹筍總量;fi,5表示每只亞成年大熊貓一年中需要消耗的竹子總量;fi,6表示每只亞成年大熊貓一年中需要消耗的其他食物總量;fi,7表示每只成年大熊貓和老年大熊貓一年中需要消耗的竹筍總量;fi,8表示每只成年大熊貓和老年大熊貓一年中需要消耗的竹子總量;fi,9表示每只成年大熊貓和老年大熊貓一年中需要消耗的其他食物總量。
cmin表示每年救護(hù)大熊貓的最少數(shù)量;cmax表示每年救護(hù)大熊貓的最大數(shù)量;cmaxage表示救護(hù)大熊貓的最大年齡;pcc表示一年中救護(hù)c只大熊貓的概率;pgi表示救護(hù)大熊貓性別為i的概率;paj表示救護(hù)大熊貓年齡為j的概率;pdd表示一年中野放d只大熊貓的概率;psi表示野放大熊貓性別為i的概率;pyj表示野放大熊貓年齡為j的概率。
膜計(jì)算是自然計(jì)算的新分支,由歐洲科學(xué)院院士、羅馬尼亞科學(xué)院院士Paun提出[34]。膜計(jì)算受生物細(xì)胞的啟發(fā),旨在從細(xì)胞組織的結(jié)構(gòu)和功能中抽象出計(jì)算模型[35-36]。膜系統(tǒng)(計(jì)算模型)由膜結(jié)構(gòu)、初始對(duì)象集和進(jìn)化規(guī)則集三要素構(gòu)成,膜系統(tǒng)中膜結(jié)構(gòu)分隔了對(duì)象和規(guī)則所處的區(qū)域,對(duì)象按照所在區(qū)域的規(guī)則不斷進(jìn)化,進(jìn)化規(guī)則以極大并行的方式執(zhí)行。膜系統(tǒng)運(yùn)行開(kāi)始前的狀態(tài)稱(chēng)為初始格局,膜系統(tǒng)自初始格局開(kāi)始根據(jù)不同的進(jìn)化規(guī)則從一個(gè)格局進(jìn)化到下一個(gè)格局,直到?jīng)]有規(guī)則可以執(zhí)行時(shí)系統(tǒng)終止,此時(shí)的膜系統(tǒng)狀態(tài)稱(chēng)為終止格局,同時(shí)終止格局時(shí)可以指定膜系統(tǒng)里不同區(qū)域中的對(duì)象數(shù)目作為膜系統(tǒng)的輸出結(jié)果。多環(huán)境概率膜系統(tǒng)是在膜系統(tǒng)的基礎(chǔ)上引入概率函數(shù)而得到的一種新的膜系統(tǒng)概念,種群動(dòng)態(tài)膜系統(tǒng)是多環(huán)境膜系統(tǒng)的一個(gè)變種[37]。
一個(gè)度為(q,m),同時(shí)滿(mǎn)足q,m≥1且T≥1種群動(dòng)態(tài)膜系統(tǒng)可形式化定義為:
其中:
G=(V,S) 是一個(gè) V={e1,e2,…,em} 的有向圖。
Γ和Σ是字母表。Γ是包含有限對(duì)象的非空字母表,由種群動(dòng)態(tài)膜系統(tǒng)所有區(qū)域中的所有對(duì)象組成;Σ代表字母表Γ中可以出現(xiàn)在所有環(huán)境中的對(duì)象組成的字母表。
Τ表示生態(tài)系統(tǒng)進(jìn)化過(guò)程被仿真的年數(shù)。
RE表示有限規(guī)則集,形如
其中:x,y1,y2,…,yh∈ Σ,(ej,ejl) ∈ S,1 ≤l≤ h,1 ≤ k≤ n,且p,p'是定義域?yàn)閧1,2,…,T}的可計(jì)算的函數(shù)。
u是集合{1,2,…,q}×{0,+,-}中元素單射標(biāo)號(hào)的有q個(gè)節(jié)點(diǎn)的根樹(shù),如果膜標(biāo)號(hào)為(i,α)那么這個(gè)膜將表示為同時(shí)該膜的標(biāo)號(hào)為i控制電荷為α,樹(shù)的根節(jié)點(diǎn)設(shè)為1作為關(guān)聯(lián)標(biāo)簽。
對(duì)于每一個(gè)r∈R和1≤j≤m,fr,j是一個(gè)定義域?yàn)閧1,2,…,T}的可計(jì)算函數(shù)。
對(duì)于每一個(gè) i,j(1 ≤ i≤ q,1 ≤ j≤ m),Mi,j是字母表 Γ中的有限多重集。
對(duì)于每一個(gè)j,1≤j≤m,Ej是字母表Σ中的有限多重集。
整個(gè)模型由初始模塊、繁殖模塊、救護(hù)模塊、進(jìn)食模塊、死亡模塊、野放模塊和更新模塊6個(gè)部分組成。
計(jì)算模型考慮使用度為(2,1)時(shí)間單元為T(mén)的種群動(dòng)態(tài)膜系統(tǒng)模型
定義如下:
G是一個(gè)空表。
工作字母表Γ為:
其中:Xi,j,y表示繁殖模塊之前性別為i年齡為j的大熊貓個(gè)體,Yi,j,y表示進(jìn)入死亡模塊性別為i年齡為j的大熊貓個(gè)體,Zi,j,y表示存活性別為 i年齡為 j的大熊貓個(gè)體,Wi,j,y表示進(jìn)食模塊之后性別為i年齡為j的大熊貓個(gè)體,Ci,j,y表示救護(hù)性別為i年齡為j的大熊貓個(gè)體,y表該個(gè)體目前所處的仿真周期。對(duì)象S,B,O表示不同類(lèi)型的食物:分別為竹筍、竹子和其他食物。對(duì)象F,A1,G1是輔助符號(hào)分別表示如下:F是一個(gè)符號(hào)用來(lái)作為在每個(gè)周期初始產(chǎn)生新的食物,A1是一個(gè)符號(hào)用來(lái)作為觸發(fā)救護(hù)個(gè)體的行為;G1是一個(gè)符號(hào)用來(lái)作為觸發(fā)野放個(gè)體的行為,下標(biāo)1表示目前為第一個(gè)仿真周期。
T表示在進(jìn)化生態(tài)系統(tǒng)中仿真的年數(shù)。
u=[[]1]2表示采用兩層膜結(jié)構(gòu)。
M1,1和M2,1是Γ中的有限多重集,描述區(qū)域u中的初始對(duì)象:
qi,j表示性別為i年齡為j的初始大熊貓數(shù)量
M2,1={F A1G1}
1)初始規(guī)則。
為大熊貓?zhí)峁┦澄?/p>
2)繁殖規(guī)則。
未達(dá)到繁殖期年齡的個(gè)體
為了驗(yàn)證上述設(shè)計(jì)的中國(guó)大熊貓保護(hù)研究中心圈養(yǎng)大熊貓種群動(dòng)態(tài)膜系統(tǒng)計(jì)算模型的有效性和正確性,以2005—2016年全球大熊貓譜系數(shù)據(jù)為依據(jù),分別統(tǒng)計(jì)了各年份中國(guó)大熊貓保護(hù)研究中心圈養(yǎng)大熊貓的數(shù)量,用實(shí)際數(shù)據(jù)和實(shí)驗(yàn)仿真結(jié)果進(jìn)行對(duì)比。把2005年各年齡段的大熊貓數(shù)量作為初始參數(shù)輸入,2006—2016年各年齡段的大熊貓數(shù)量作為結(jié)果輸出。仿真實(shí)驗(yàn)所用計(jì)算機(jī)為戴爾 Inspiron N5110,2.10 GHz,3.16 GB,使用 MeCoSim 仿真軟件作為仿真平臺(tái)。
實(shí)驗(yàn)中的參數(shù)設(shè)置為:仿真周期cycles設(shè)置為11,實(shí)驗(yàn)中一個(gè)仿真周期代表大熊貓生態(tài)系統(tǒng)中的一個(gè)自然年,周期運(yùn)行步數(shù) steps為5,為保證實(shí)驗(yàn)結(jié)果的可靠性仿真次數(shù)simulations設(shè)置為1000次。
大熊貓生態(tài)系統(tǒng)通常涉及參數(shù)由于自然變化以及許多因素影響大熊貓個(gè)體的健康、死亡率和生物學(xué)現(xiàn)象,故所研究的過(guò)程本質(zhì)上是隨機(jī)的。在種群定投膜系統(tǒng)中,進(jìn)化規(guī)則以概率方式執(zhí)行,可以有效處理含隨機(jī)性、不確定性的數(shù)據(jù)。部分進(jìn)化規(guī)則的關(guān)鍵參數(shù)設(shè)置如表1~3所示。
表2 大熊貓出生率及其繁殖期年齡Tab.2 Birth ratio and reproductive age of giant panda
表3 野放大熊貓數(shù)量及年齡Tab.3 Quantity and age of released giant panda
不同年份雄性大熊貓和雌性大熊貓數(shù)量變化的實(shí)驗(yàn)仿真結(jié)果與實(shí)際數(shù)據(jù)對(duì)比分別如圖1所示。
圖1 雌雄大熊貓數(shù)量實(shí)驗(yàn)結(jié)果與實(shí)際數(shù)據(jù)對(duì)比Fig.1 Comparison of experimental and real results with respect to the numbers of female and male giant pandas
不同年份大熊貓總體數(shù)量變化實(shí)驗(yàn)仿真結(jié)果與實(shí)際數(shù)據(jù)進(jìn)行對(duì)比,結(jié)果如圖2所示。
圖2 大熊貓總體數(shù)量變化Fig.2 Changes in the overall number of giant pandas
從圖1~2的仿真結(jié)果可以看出,本文設(shè)計(jì)的種群動(dòng)態(tài)膜系統(tǒng)模型仿真結(jié)果與真實(shí)數(shù)據(jù)相比較,雄性大熊貓仿真結(jié)果與實(shí)際數(shù)據(jù)的相對(duì)誤差最大不超過(guò)±8.95%,大部分控制在±2.5%以?xún)?nèi);雌性大熊貓仿真結(jié)果與實(shí)際數(shù)據(jù)的相對(duì)誤差最大不超過(guò)±7.81%,大部分控制在±2.4%以?xún)?nèi);全體大熊貓仿真結(jié)果與實(shí)際數(shù)據(jù)的相對(duì)誤差最大不超過(guò)±4.13%,大部分控制在±2.7%以?xún)?nèi)。同時(shí)由于種群基數(shù)較小導(dǎo)致相對(duì)誤差較為明顯,但依然遠(yuǎn)小于±10%的相對(duì)誤差標(biāo)準(zhǔn)。同時(shí)與之前采用概率膜系統(tǒng)設(shè)計(jì)的簡(jiǎn)單模型相比,該模型實(shí)驗(yàn)結(jié)果誤差更小,擬合度更高,故本文設(shè)計(jì)的大熊貓種群動(dòng)態(tài)膜系統(tǒng)較好地模擬了中國(guó)大熊貓保護(hù)研究中心大熊貓生態(tài)系統(tǒng)種群動(dòng)態(tài)變化趨勢(shì)。
本文提出了一種研究大熊貓種群動(dòng)態(tài)的新方法,種群動(dòng)態(tài)膜系統(tǒng)模型能夠有效預(yù)測(cè)種群未來(lái)的變化趨勢(shì),可以及時(shí)評(píng)估種群發(fā)展的合理性,設(shè)計(jì)出一個(gè)基于大熊貓譜系數(shù)據(jù)的具有雙層膜結(jié)構(gòu)、對(duì)象集和一系列進(jìn)化規(guī)則的種群動(dòng)態(tài)膜系統(tǒng)模擬實(shí)際生態(tài)系統(tǒng)中大熊貓固有的隨機(jī)進(jìn)化過(guò)程。這是首次搜集和分析如此完整的譜系數(shù)據(jù),并以此為基礎(chǔ)研究大熊貓的種群動(dòng)態(tài)。同時(shí)有專(zhuān)門(mén)開(kāi)發(fā)的仿真軟件MeCoSim來(lái)協(xié)助設(shè)計(jì)該計(jì)算模型并仿真驗(yàn)證所提出的計(jì)算模型。它還可以進(jìn)行基于仿真實(shí)驗(yàn)結(jié)果的決策,診斷種群異常的原因,分析外界因素對(duì)種群動(dòng)態(tài)性的影響,為大熊貓基地的工作人員提供決策參考,避免盲目決策。