华南理工大学数学实验实验七

x0
800 x0
x0
600
x0
v1
x0
400
x0
200
v2
0 0 200 400 600 800 1000 1200 1400 1600 1800 2000 || >1,|u|<1
图 3 增大老鼠的出生率
（2）将上述代码中矩阵 A 中的 1.1 改为 1.0，即减小老鼠的出生率，该系统 的演化过程如图 4 所示，从图中可以看出，由于老鼠出生率的减小，猫头鹰的食 物来源相应的减少，最终猫头鹰将灭亡。
text(X0(1,1),X0(2,1),'x0')
quiver([X(1,2),1]',[X(2,2),1]',[X(1,1)-X(1,2),0]',[X(2,1)-X
(2,2),0]',p)
set(h,'MarkerSize',6),grid,
end
end
结果及结果分析：
运行上面的代码，可以得到结果如图 1 和 2 所示，图 1 为矩阵 A 的特征值
7
xk

c1 (1.00)k
4

5

k 越大，上式的近似程度越高，所以对于足够大的 k，
xk +1

c1
(1.00)k +1

4 5

=c1
(1.00)k

4 5

=xk 上式表明，最后 xk 的每个分量（猫头鹰和老鼠的数量）几乎每个月都近似的等于
原来的数量。 xk 约为
图8
3． 实验总结和实验感悟
通过这次的实验，我重新了解并掌握了很多之前学过的线性代数的知识点， 例如特征值、特征向量、特征方程、矩阵对角化等等以及使用 MATLAB 的 eig 命令来直接求出矩阵的特征值和特征向量，理解了由差分方程所描述的动态系统 的长期行为和演化过程。在这次实验中，我还学会了利用 MATLAB 函数来调用 鼠标从而获取点击位置坐标值的方法，这是一个实用性比较强的方法，不仅仅可 以使用在本次实验中，还可以用在以后的很多方面上。
要求： （1）建立代数模型，从理论上说明最终的基因型分布。 （2）用 MATLAB 求解初始分布为 0.8，0.2，0 时，20 年后基因分布，是否
1
已经趋于稳定？
概率
后代的
AA
Aa
基因型
aa
2.1. 实验原理
1、特征值与特征向量
表 1 基因的转移 父体—母体基因型
AA-AA 1 0 0
AA-Aa 1/2 1/2 0
button = 1;
while button == 1
[xi,yi,button] = ginput(1);
%用鼠标选初始点
plot(xi,yi,'go'),hold on
X0 = [xi;yi];
X = X0;
for i=1:n
X = [A*X, X0];
%用这种方式迭代，并画图
h = plot(X(1,1),X(2,1),'R.',X(1,1:2),X(2,1:2),'r-'); hold on
x0 v2
0 0 200 400 600 800 1000 1200 1400 1600 1800 2000 || >1,|u|<1
图 2 不同初始分布的演化过程
当该系统处于图 2 所示的不稳定平衡状态时，该系统的稳定性非常脆弱，当 系统受到外界干扰或者系统自身的微小变化时，该系统将发生不同的演化并且有 截然不同的演化结果。具体分析如图 3~6 所示，分别为：图 3 老鼠出生率增大； 图 4 老鼠出生率减小；图 5 猫头鹰捕食率增大；图 6 猫头鹰捕食率减小。
2000 1800 1600 1400 1200 1000
800 600 400 200
x0 x0
x0 x0
v1
x0
x0 x0 x0 x0 x0 v2
0 0 200 400 600 800 1000 1200 1400 1600 1800 2000 || >1,|u|<1
图 4 减小老鼠的出生率
9

因为v1, v2 是 R2 的一组基，所以存在系数 c1 和 c2 ，使得
因为v1, v2 为矩阵 A 对应于特征值 =1，=0.6 的特征向量，所以
于是
一般的，
当 k → 时， (0.6)k 迅速趋于 0。假定 c1 > 0，则对于所有足够大的 k， xk 近似地
等于 c1v1 ，写为
[pc,lambda] = eig(A);
%求 A 的特征值和对应的特征向量
[Y,I] = sort(diag(abs(lambda)),'descend'); %对特征值的绝对值降序排列
temp = diag(lambda);
lambda = temp(I)
%输出按特征值的绝对值降序排列的特征值
pc = pc(:,I)
8
（1）将上述代码中矩阵 A 中的 1.1 改为 1.2，即增大老鼠的出生率，该系统 的演化过程如图 3 所示，从图中可以看出，该系统从原先的不稳定平衡状态转变 为稳定平衡状态，最终猫头鹰和老鼠的数量将达到动态平衡状态，两者以一定的 数量比例存在。
2000
1800
1600
1400
x0
1200
x0
1000
第一步：求 A 的特征值和对应的特征向量； 第二步：用 v1 和 v2 表示 x0 和 xk，k = 1, 2,... ； 第三步：解的图像表示。 代码：
%%%%%%%% 实验练习 8-1 %%%%%%%%% %捕食者-被捕食者解的图像表示 %捕食参数 p 为 0.125 close all; %关闭已经打开的所有图片 clear,clc; %清空工作区和命令行窗口 a = 0; b = 2000; c = a; d = b; p = 0.1; %确定画图范围
n→
=
P(limT n )P-1 n→
=
P


1

P
-1
0
11
图 7 为理论分析的结果，可以看出，最终 Aa 基因型和 aa 基因型都会消失。
图7
如图 8 所示为 20 年后的基因型分布，从结果可以看出，第 20 代以后基因型 趋于稳定，不但不会出现 aa 型品种，连 Aa 型品种也趋于消失。因此为了培育优 良品种，最好采用一些控制结合的手段。
AA-aa 0 1 0
2、特征值与特征向量的求法
2
3、矩阵的对角化 4、离散线性动态系统
3
4
5、eig 命令
2.2. 算法、编程与结果分析
➢ 第1题 算法：
将线性变换 x → Ax 的作用分解为易于理解的成分，其中特征值和特征向量 是分析离散动态系统的关键。根据已知信息，找出系统对应的差分方程 xk+1 = Axk ，求出 A 的特征值和对应的特征向量。再根据不同特征值的个数、绝 对值大于 1 还是小于 1、是实特征值还是复特征值的情形，分析出系统的演化过 程。主要分为三步：
pc = -pc;
z1 = pc(2,1)/pc(1,1)*x; z2 = pc(2,2)/pc(1,2)*x;
%特征向量 v1 %特征向量 v2
h = plot(x,z1),set(h,'linewidth',2), text(x(7),z1(7)-100,'v1')
h = plot(x,z2),set(h,'linewidth',2), text(x(20),z2(20)-100,'v2')

4 5

的倍数，所以
xk
的
2
个分量的比值约为
4
比
5，即每
4
只猫头鹰对应着 5000 只老鼠。图 2 为解的图像表示，从结果图中可以看出，该
系统趋向一种被称为不稳定平衡的状态。
2000 1800 1600
x0 x0
1400
x0 1200
1000 x0
x0
800
600
x0
400
v1
200
x(1) 1
,
x(1) 2
,
x(1) 3
，以后第
n
代的分
布为
x(n) 1
,
x(n) 2
,
x(n) 3
。则有，
则， 代码：
1
A
=

0
1/ 2 1/ 2
0 1 ，X
(n)
=


x(n) 1
x(n) 2
,
X
(1)
=


80%
20%

0 0 0

x(n) 3
lambda（按照特征值得绝对值降序排列）和对应的特征向量 pc；图 2 为不同初
始分布的演化过程。
6
图 1 A 的特征值以及对应的特征向量
这两个特征向量（即 pc 的第一列和第二列）是线性无关的，它们构成 R2 的 一组基。为消除小数，选取
v1
=

4 5

,
v2
=415源自n = 100;%序列迭代次数
xlabel('|\lambda| >1,|u|<1');
axis([a b c d]); %限制画图空间 grid on,hold on; %插入网格并保持显示的图像不变
x = linspace(a,b,30);