作者简介:孙明娟(1981—),女,延安大学讲师.研究方向:应用概率统计. E-mail:sunmingjuan 0535@163.com
中文责编:方 圆; 英文责编:木 南
延安大学数学与计算机科学学院, 陕西延安 716000
Sun Mingjuan and Dong QinglaiSchool of Mathematics and Computer Science, Yan'an University, Yan'an 716000, Shaanxi Province, P.R.China
stochastic differential equation; ratio-dependent functional response; stochastic Chemostat; stochastically asymptotic stable; Itô; formula; dynamical behavior
DOI: 10.3724/SP.J.1249.2016.04425
研究随机噪声对微生物连续培养的影响,引入白噪声描述营养消耗率受随机噪声的干扰,建立一类具有比率型功能反应函数的随机恒化器模型. 通过构造Lyapunov函数,利用Itô公式证明系统正解的全局存在唯一性,并论证了系统的绝灭平衡点是全局随机渐近稳定的; 探讨了噪声强度大小对随机模型的解围绕相应确定性模型的平衡点振荡行为的影响.通过数值模拟验证所得理论结果的正确性.
To explore the influence of random noise on microbial continuous culture, we establish a stochastic Chemostat model with ratio-dependent functional response in which white noise is introduced to describe nutrition conversion rate influenced by random noise. By constructing stochastic Lyapunov function and using the Itô formula, we prove that there is a unique positive solution in the system with positive initial value, and the washout equilibrium is stochastically asymptotic stable. Furthermore, we investigate the impact of white noise on the behavior of the solution spirals around the positive equilibrium of deterministic system. The numerical simulations support the proposed theoretical results.
恒化器(Chemostat)是实验室中连续培养微生物的一种实验装置,长期以来恒化器模型也是生物与数学工作者研究的热点[1-3]. 随着随机微分方程广泛用于生物数学, 研究注意到微生物培养过程中如果忽略环境中的随机因素,可能导致对微生物增长规律的描述不够准确,如文献[4- 6]指出,噪声强度较大时可能会导致恒化器内微生物的整体溢出. 文献[7- 8]考虑到微生物和营养之间存在的捕食和被捕食关系,根据比率依赖理论分别研究流出率和半饱和常数受随机噪声的干扰,且具有比率型功能反应函数的随机恒化器模型,研究结果一致表明,环境噪声对微生物培养有重要影响,环境噪声太大会使随机恒化器模型的解围绕确定性模型解的振荡幅度增大,一定条件下会导致培养失败.因此,建立随机模型时考虑随机噪声的影响来刻画恒化器中微生物的培养过程更符合客观实际,更能准确把握种群增长规律.
文献[1-2,9-12]对具有常数消耗率的恒化器模型研究方法及相关结论做了介绍. 而文献[13]通过研究利用啤酒酵母进行乙醇发酵的动力学方程后,指出细胞得率系数(消耗率)不应是常数,改进后的模型出现振荡现象. 因此,在基本模型的基础上考察消耗率受随机噪声干扰的情形更符合实际生物背景. 文献[14]考虑一类营养的转化率受到随机噪声干扰,具有Holling II 型功能反应函数的随机恒化器模型,揭示了在不同条件下模型的解围绕其相应的确定性模型各类平衡点的振荡行为,但营养消耗率受随机噪声干扰的具有比率型功能反应函数的恒化器动力学模型的研究尚未见报道.本研究讨论随机噪声的存在对微生物连续培养动力行为的影响,研究如下的随机恒化器模型:
{dS=(D(S0-S)-1/δ(mSx)/(S+x))dt-α(mSx)/(S+x)dBt
dx=x((mS)/(S+x)-D)dt(1)
其中, S(t)表示t时刻培养皿中营养的浓度; x(t)表示t时刻微生物的密度; S0为输入的初始营养浓度; D表示流出率; 1/δ表示营养吸收转化率; m是半饱和常数, 且S0、 δ和m均为正常数; Bt和α分别为定义在完备概率空间上的标准布朗运动和噪声强度, 且环境噪声对系统的影响表现为对营养消耗率的影响上,即受环境影响1/δ变为1/δ+αB·t.
1 预备知识[15]
对于n维随机系统
dX(t)=f(t,X(t))dt+g(t,X(t))dBt(2)
其中,
f(t,X(t)):[t0,+∞)×Rn→Rn;
g(t,X(t)):[t0,+∞)×Rn→Rn×m.
f(t,X(t))和g(t,X(t))满足局部Lipschitz条件.将系统(2)线性化得
dX(t)=F(t)X(t)dt+G(t)X(t)dBt(3)
引理 设F(t)=F, G(t)=G, 如果系统(3)的零解X(t)=0是渐近稳定的,且在X(t)=0的足够小邻域内,存在一个足够小的常数l使系统(2)和系统(3)的系数满足
|f(t,X)-F·X|+|g(t,X)-G·X|<
l|X|(4)
则系统(2)的零解X(t)=0是全局随机渐近稳定的.
本研究假设(Ω,{Ft}t≥0,P)为完备的概率空间,其中, {Ft}t≥0是Ω上的一个σ代数流且满足通常条件(即右连续, F0含所有零测集),并记
R2+={x∈R2:xi>0,i=1,2}
R-2+={x∈R2:xi≥0,i=1,2}
对于系统(1)的正解的全局存在唯一性有如下定理.
定理1 对于任意给定的初值(S(0),x(0))∈R2+, 模型(1)存在唯一解(S(t),x(t))(t≥0), 且该解以概率1位于R2+中, 即对于所有的t≥0,(S(x),x(t))∈R2+a.s..
【证】显然,模型(1)对于任意的t∈[0,τe)存在唯一的局部解(S(t),x(t)). 其中, τe为爆破时间, 下面证明τe=∞a.s..
取m0>0足够大,使
S(0)∈[1/(m0), m0], x(0)∈[1/(m0), m0]
又对于每一个整数m≥m0, 定义停时τm=inf {t∈[0,τe):S(t)(1/m, m)或x(t)(1/m, m)},对于空集, 定义inf=∞. 显然, 当m→∞时, τm是递增的.定义τ∞=limm→∞τm, 由于τ∞≤τea.s., 因此要证明τe=∞a.s., 只需证明τ∞=∞a.s.. 为此定义如下的C2函数V: R2+→R-+
V=A(S-1-lnS)+(x-1-lnx)
其中, A为待定的正常数.由于对任意的u>0, u-1-lnu>0, 因此函数V是正定的.
由Itô公式得
dV=AdS-A/SdS+A/(2S2)(dS)2+
dx-(dx)/x+((dx)2)/(2x2)=
LVdt-[(AαmSx)/(δ(S+x))-(Aαmx)/(δ(S+x))]dBt(5)
其中,
LV=ADS0-ADS-A/δ·(mSx)/(S+x)-(ADS0)/S+AD+(Amx)/(δ(S+x))+(Aα2m2x2)/(2(S+x)2)+(mSx)/(S+x)-xD-(mS)/(S+x)+D≤ADS0-A/δ·(mSx)/(S+x)+AD+(Am)/δ+(Aα2m2)/2+(mSx)/(S+x)+D
令A=δ,则
LV≤δDS0+δD+m+(δα2m2)/2+DK
由式(5)得
dV≤Kdt+((αmx)/(S+x)-(αmSx)/(S+x))dBt
上式两端从0到τm∧t积分,再取期望得
E(V(S(τm∧t),x(τm∧t)))≤
V(S(0),x(0))+Kt(6)
又因
E(V(S(τm∧t),x(τm∧t)))≥
E(I<sub>χ(τm≤t)V(S(τm∧t),x(τm∧t)))=
V(S(τm),x(τm))P(τm≤t),(7)
且
V(S(τm),x(τm))≥
A(m-1-lnm)∧(1/m-1-ln1/m)
其中, I<sub>χ(τ</sub>m≤t)是χ(τm≤t)的示性函数. 显然,当m→∞时, V(S(τm),x(τm))→∞. 由式(6)和式(7)有
P(τm≤t)≤(V(S(0),x(0))+Kt)/(V(S(τm),x(τm))),
当m→∞时, P(τ∞≤t)=0, 所以τ∞=∞,故τe=∞. 定理证毕.
当m<D时, E0=(S0,0)是随机模型(1)相应确定性模型的唯一平衡点且全局渐近稳定,显然它仍然是随机模型(1)的平衡点.以下定理说明,当m<D时,在噪声强度满足一定条件时,随机模型(1)的平衡点E0是全局随机渐近稳定的.
定理2 当m<D时,若α2<(2(D-m))/(m2), 则模型(1)的绝灭平衡点E0=(S0,0)是随机全局渐近稳定的.
【证】设(S(t),x(t))是模型(1)满足任意给定初值(S(0),x(0))∈R2+的解.
令x1=S-S0, x2=x, 则随机模型(1)转化为
{Zdx1=(-Dx1-1/δ(m(x1+S0)x2)/(x1+S0+x2))dt-
α(m(x1+S0)x2)/(x1+S0+x2)dBt
dx2=x2((m(x1+S0))/(x1+S0+x2)-D)dt(8)
将式(8)线性化得
{Zdx1=(-Dx1-m/δx2)dt-αmx2dBt
dx2=(mx2-Dx2)dt(9)
下面证明模型(9)的零解全局随机渐近稳定.
考虑Lyapunov函数
V=x21+x22+Mx2
其中, M为待定的正常数. 显然, V是正定函数. 由Itô公式得
dV=2x1dx1+(dx1)2+2x2dx2+(dx2)2+
Mdx2=LVdt-2αmx1x2dBt
其中,LV=-2Dx21-(2D-m2α2-2m)x22+M(m-D)x2-(2m)/δx1x2, 可以证明当α2<(2(D-m))/(m2)时, LV为负定函数.事实上:
① 当x1≥0时,若α2<(2(D-m))/(m2),则
LV≤-2Dx21-(2D-m2α2-2m)x22≤0
② 当x1<0时,由S=x1+S0>0, 所以,-x1<S0, 则
LV<-2Dx21-(2D-m2α2-2m)x22+M(m-D)x2+(2m)/δS0x2
令M=(2mS0)/(δ(D-m)),则有
LV<-2Dx21-(2D-m2α2-2m)x22≤0.
由①和②可知, LV≤0且等号成立的充要条件为x1=x2=0, 所以LV负定.根据文献[15], 模型(9)的零解随机渐近稳定.
将模型(8)和模型(9)写成式(2)和式(3)的形式,对应式(4)的左端
|f(t,X)-F·X|+|g(t,X)-G·X|=
((1/δ(m(x1+S0)x2)/(x1+S0+x2)-m/δx2)2+((mx2(x1+S0))/(x1+S0+x2)-mx2)2)1/2+((αmx2-(αmx2(x1+S0))/(x1+S0+x2))2)1/2=
(m2x22(x22)/((x1+S0+x2)2))1/2((α2)1/2+(1/(δ2)+1)1/2)
对任意小的正数ε<0, 当|x1|<ε, |x2|<ε时,存在常数
l=m((α2)1/2+(1/(δ2)+1)1/2)
使得|f(t,X)-F·X|+|g(t,X)-G·X|<l|X|.
由引理,模型(8)的零解全局随机渐近稳定.从而模型(1)的绝灭平衡点E0=(S0,0)全局随机渐近稳定.
平衡点的振荡行为
当m>D时, E*=(S*, x*)=((D δ S0)/(m-D+δD), δ(S0-(D δ S0)/(m-D+δD)))是模型(1)所对应的确定性系统的唯一全局渐近稳定的正平衡点. 此时, E*不再是随机模型(1)的平衡点.但可通过研究模型(1)的解在E*=(S*, x*)附近的振荡行为来了解其解的渐近性态.
定理3 当m>D时,若α2<D/(m2δ2), 则对任意给定的初值((S(0),x(0))∈R2+, 模型(1)的解(S(t),x(t))具有以下性质:
limt→∞sup1/tE∫0t[(S(u)-S*)2+1/(δ2)(x(u)-x*)2]du≤(α2m2δ2S20)/(D-α2m2δ2).
【证】考虑Lyapunov函数
V=(S+x/δ-S0)2+Λ(x-x*-x*lnx/(x*)),
其中, Λ是待定正数,显然V是正定函数,由Itô公式得
dV=2(S+x/δ-S0)d(S+x/δ)+(d(S+x/δ))2+
Λ[(1-(x*)/x)dx+(x*)/(2x2)(dx)2]=
LVdt-2(S+x/δ-S0)(amSx)/(S+x)dBt(10)
其中,
LV=-2D(S-S*+x/δ-(x*)/δ)2+((αmSx)/(S+x))2+Λ(x-x*)((mS)/(S+x)-D).
注意到S0=S*+(x*)/δ, D=(mS*)/(S*+x*), 所以,
((αmSx)/(S+x))2≤α2m2x2≤α2m2δ2(S+x/δ)2≤
2α2m2δ2[(S+x/δ-S0)2+S20]=
2α2m2δ2[(S-S*+x/δ-(x*)/δ)2+S20],
其次,对于Λ(x-x*)((mS)/(S+x)-D):
1)若(S-S*)(x-x*)<0, 则
Λ(x-x*)((mS)/(S+x)-D)=Λ(S-S*)(x-x*)
[(mS)/((S+x)(S-S*))-(mS*)/((S*+x*)(S-S*))]≤
Λ(S-S*)(x-x*)
[(mS)/((S+x)S)-(mS*)/(-(S*+x*)S*)]≤
(Λm)/(S*+x*)(S-S*)(x-x*);
2)若(S-S*)(x-x*)≥0, 则
Λ(x-x*)((mS)/(S+x)-D)=
Λm(x-x*)[S/(S+x)-(S*)/(S*+x*)]
① 若x≥x*, 则S≥S*, 那么
S/(S+x)-(S*)/(S*+x*)≤S/(S*+x*)-(S*)/(S*+x*)
得
Λ(x-x*)((mS)/(S+x)-D)≤
(Λm)/(S*+x*)(S-S*)(x-x*)
② 若x<x*, 则S≤S*,则
S/(S+x)-(S*)/(S*+x*)≥S/(S*+x*)-(S*)/(S*+x*)
同样得
Λ(x-x*)((mS)/(S+x)-D)≤
(Λm)/(S*+x*)(S-S*)(x-x*)
所以
LV≤-(2D-2α2m2δ2)
[(S-S*)2+1/(δ2)(x-x*)2]+2α2m2δ2S20-
[4/δ(D-α2m2δ2)-(Λm)/(S*+x*)]
(S-S*)(x-x*)
当D-α2m2δ2>0, 即α2<D/(m2δ2)时,取
Λ=(4(D-α2m2δ2)(S*+x2))/(mδ)
则
LV≤-2(D-α2m2δ2)
[(S-S*)2+1/(δ2)(x-x*)2]+
2α2m2δ2S20(11)
对式(10)两端从0到t积分,再取期望得
E∫0tdV=E∫0tLVdt
因此,由式(11)得
limt→∞ sup1/t E∫0t[(S(u)-S*)2+1/(δ2)(x(u)-x*)2]du≤(α2m2δ2S20)/(D-α2m2δ2)
研究一类营养消耗率受随机噪声干扰的具有比率型功能反应函数的随机恒化器模型. 结果表明,随机系统(1)正解的全局存在唯一性不受噪声影响[定理1]. 当噪声强度不大时, 随机系统(1)的绝灭平衡点是全局随机渐近稳定性[定理2]. 尽管随机系统(1)不存在正平衡点, 但当噪声强度不大时, 随机系统(1)的解围绕确定系统正平衡点做随机振荡, 且振幅是有限的[定理3].
从生物学意义来看,正平衡点的存在及全局渐近稳定在一定程度上反映了微生物能够培养成功,而灭绝平衡点的存在及全局渐近稳定反映了微生物的培养是失败的.本研究结果表明,环境噪声的确对微生物的培养有一定影响,但是控制各种环境噪声足够小时,不会改变原有确定性模型的性质.而在确定性模型中,控制溢流量的大小是培养成功的关键,当m<D时,培养失败; 当m>D时,培养成功.本研究得出控制环境噪声是保证培养成功的重要因素, 噪声越小培养成功的几率越大.
为验证本研究所得理论结果,根据文献[16]的方法进行一系列的数值模拟,模拟过程采用Euler-Maruyama方法.
在数值模拟中, 初始值为(S(0), x(0))=(0.5,0.5), 选取参数m=0.4、 D=0.5、 S0=1及δ=0.9, 则系统(1)仅存在边界平衡点. 选取α=0.5, 则定理3的条件满足,图1的数值模拟结果表明随机模型(1)的解是随机全局渐近稳定的. 选取α=5, 则定理3的条件不满足,图2的数值模拟结果显示振荡幅度明显大于图1的情形.
在数值模拟中,初始值为(S(0), x(0))=(0.5,0.5), 选取参数m=0.5、 D=0.15、 S0=1及δ=0.9, 则系统(1)所对应的确定性系统存在正平衡点且S*=0.649 5, x*=0.278 4. 选取α=0.1, 则定理3条件满足,图3数值模拟结果说明,随机模型(1)的解围绕确定性模型的平衡点振荡. 选取α=0.9, 则定理3的条件不满足, 图4数值模拟结果显示振荡幅度明显大于图3的情形.
图1 m<D, α=0.5, 确定性模型与随机模型的解曲线
Fig.1 Solution curves of deterministic model and stochastic model while m<D, α=0.5
图2 m<D, α=5, 确定性模型与随机模型的解曲线
Fig.2 Solution curves of deterministic model and stochastic model while m<D, α=5
图3 m>D, α=0.1, 确定性模型与随机模型的解曲线
Fig.3 Solution curves of deterministic model and stochastic model while m>D, α=0.1
深圳大学学报理工版
JOURNAL OF SHENZHEN UNIVERSITY SCIENCE AND ENGINEERING
(1984年创刊 双月刊)
主 管 深圳大学
主 办 深圳大学
编辑出版 深圳大学学报理工版编辑部
主 编 阮双琛
国内发行 深圳市邮电局
国外发行 中国国际图书贸易集团有限公司(北京399信箱)
地 址 北京东黄城根北街16号
邮 编 100717
电 话 0755-26732266
0755-26538306
Email journal@szu.edu.cn
标准刊号 ISSN 1000-2618
CN 44-1401/N