Mathematica仿真竞争性Lotka-Volterra方程(3种群)
hong在生物学上，物种（英语：species）是生物分类的基本单位。物种位于生物分类法中最底级，在属之下。笼统来讲，同一物种的任意两个性别适当的健康个体，都能够交配繁殖出具生殖能力的后代。满足以上条件的最大群体称为物种。根据演化生物学家恩斯特.麦尔的定义，物种是：能够相互配育的自然种群的类群，且这些类群与其它的类群在生殖上相互隔离。昆虫学家陈世骧（1978）对物种所下定义为：物种是繁殖单元，由又连续又间断的居群所组成；物种是进化单元，是生物系统线上的基本环节，是分类的基本单元。
在生物学上，种群（英语：population）指代在一定空间范围内同时生活的同种生物的全部个体；或者说是由个体组成，而且能够进行交配的群体。种群的个体之间一般享有同一个基因库。
直观的理解（可能不严谨），种群是物种的子集，值的是某物种在一定范围内的个体集。比如青蛙是一个物种，这片池塘里的青蛙是一个种群。
下面的论述中，并不严格区分物种和种群。
竞争性Lotka-Volterra方程是物种争夺某些共同资源的种群动态的简单模型。它们可以进一步推广到广义 Lotka-Volterra 方程以包括营养相互作用。
考虑3个具有Logistic动力学的3个种群和，在Lotka-Volterra公式的基础上添加额外的项来解释物种的相互作用。因此，竞争性Lotka-Volterra方程是：
这是一个常微分方程组。方程中并不显含时间，我们称这样的ODE是自治的。
这里代表物种对物种种群的影响。注意不必等于。因为是竞争版本，所有值都是非负的。
使用Mathematica来演示给定ODE参数和初始条件下的3物种数目随时间的演化情况。
(事实上，代码的框架是通用性质的，很容易修改来仿真其他的更加一般的ODE系统)
Manipulate[
sol=NDSolve[
{
x1'[t]==r1*x1[t]*(1-\[Alpha]11*x1[t]-\[Alpha]12 x2[t]-\[Alpha]13 x3[t]),
x2'[t]==r2*x2[t]*(1-\[Alpha]21*x1[t]-\[Alpha]22 x2[t]-\[Alpha]23 x3[t]),
x3'[t]==r3*x3[t]*(1-\[Alpha]31*x1[t]-\[Alpha]32 x2[t]-\[Alpha]33 x3[t]),
x1[0]==x10,
x2[0]==x20,
x3[0]==x30
},
{x1[t],x2[t],x3[t]},
{t,0,time}
];
Pane[
Grid[{ {
Show[
{
ParametricPlot3D[{x1[t],x2[t],x3[t]}/.sol,{t,0,time},
PlotStyle->{Thick,Blue}],
Graphics3D[
{Red,PointSize[Large],Point[{x1[t],x2[t],x3[t]}/.sol/.{t->0}],Green,Point[{x1[t],x2[t],x3[t]}/.sol/.{t->time}]}
]
},
If[range==="fixed",PlotRange->{ {0,x1max},{0,x2max},{0,x3max}},PlotRange->All],
BoxRatios->{1,1,1},
If[values,Ticks->Automatic,Ticks->None],
AxesLabel->(Style[#,Blue]&/@(Row/@Transpose[{Table[Subscript[Style["X",Italic],i],{i,3}],If[label,{"\n(Lowest\ntrophic level)","\n(Intermediate\ntrophic level)","\n(Top\npredators)"},{"","",""}]}])),
ImageSize->{400,400},ImagePadding->35
]},
{Text@Grid[{
{
Grid[{
{Row[{
Subscript[Style["X",Italic], 3],
" ("<>ToString[Round[100x3[t]/(x1[t] x2[t] x3[t])/.sol/.{t->time}][[1]]]<>"%)"
}]},
{Row[{
Subscript[Style["X",Italic], 2],
" ("<>ToString[Round[100x2[t]/(x1[t] x2[t] x3[t])/.sol/.{t->time}][[1]]]<>"%)"
}]},
{Row[{
Subscript[Style["X",Italic], 1],
" ("<>ToString[Round[100x1[t]/(x1[t] x2[t] x3[t])/.sol/.{t->time}][[1]]]<>"%)"
}]}
}],
BarChart[Round@Flatten[{x1[t],x2[t],x3[t]}/.sol/.{t->time}],
BarSpacing->0,
BarOrigin->Left,
ChartStyle->{ {RGBColor[.1,.9,.1],RGBColor[.4,.6,.4],RGBColor[.7,.3,.7]}},
Axes->None,
LabelingFunction->Left,
AspectRatio->.25,
ImageSize->{200,100}
]
}
},
Alignment->{Right,Left}]}
}],
ImageSize->400
],
Style["initial variable values",Bold],
{ {x10,2000,Subscript["X",1]},1,7000,.01,Appearance->"Labeled",ImageSize->Tiny},
{ {x20,175,Subscript["X",2]},1,7000,.01,Appearance->"Labeled",ImageSize->Tiny},
{ {x30,200,Subscript["X",3]},1,7000,.01,Appearance->"Labeled",ImageSize->Tiny},
Style["parameter values",Bold],
{ {r1,2.84,Subscript["r",1]},0,5,.01,Appearance->"Labeled",ImageSize->Tiny},
{ {r2,1.5,Subscript["r",2]},0,5,.01,Appearance->"Labeled",ImageSize->Tiny},
{ {r3,.62,Subscript["r",3]},0,5,.01,Appearance->"Labeled",ImageSize->Tiny},
{ {\[Alpha]11,.22},0,1,.01,Appearance->"Labeled",ImageSize->Tiny},
{ {\[Alpha]12,.22},0,1,.01,Appearance->"Labeled",ImageSize->Tiny},
{ {\[Alpha]13,.22},0,1,.01,Appearance->"Labeled",ImageSize->Tiny},
{ {\[Alpha]21,.22},0,1,.01,Appearance->"Labeled",ImageSize->Tiny},
{ {\[Alpha]22,.22},0,1,.01,Appearance->"Labeled",ImageSize->Tiny},
{ {\[Alpha]23,.22},0,1,.01,Appearance->"Labeled",ImageSize->Tiny},
{ {\[Alpha]31,.22},0,1,.01,Appearance->"Labeled",ImageSize->Tiny},
{ {\[Alpha]32,.22},0,1,.01,Appearance->"Labeled",ImageSize->Tiny},
{ {\[Alpha]33,.22},0,1,.01,Appearance->"Labeled",ImageSize->Tiny},
Style["time period",Bold],
{ {time,10,""},1,100,.01,Appearance->"Labeled",ImageSize->Tiny},
Style["maximum values on axes",Bold],
{ {range,"floating",""},{"fixed","floating"}},
{ {x1max,3100,Subscript["X",1]},0,20000,.01,Appearance->"Labeled",ImageSize->Tiny},
{ {x2max,300,Subscript["X",2]},0,5000,.01,Appearance->"Labeled",ImageSize->Tiny},
{ {x3max,500,Subscript["X",3]},0,5000,.01,Appearance->"Labeled",ImageSize->Tiny},
Style["label trophic levels",Bold],
{ {label,False,""},{True,False}},
Style["display values on axes",Bold],
{ {values,True,""},{True,False}},
TrackedSymbols->True,SynchronousUpdating->True,
ControlPlacement->Left,AutorunSequencing->{2,5,9,13,15}
]
从上图可以看到，从红色点标示的状态，经过10个时间单位，演化到绿色点表示的状态。随时间演化的轨迹线也在图中给了出来。
这个演示，除去必要的核心功能要点外，可以选择坐标轴是否显示数值，选择静态/动态设置显示出的坐标轴的最大值，还可以标定最后物种的数目比例。
来源：320科技工作室
