计算机在生命科学中应用MATLAB

目录
实验一：数据采集 (2)
实验二：绘图 (8)
实验三：函数 (14)
实验四：解非线性方程 (21)
实验五：解线性方程组 (25)
实验六：插值计算 (29)
实验七：计算积分 (35)
实验八：数据拟合 (46)
实验一：数据采集
实验学时：2
实验类型：（验证）
实验要求：（必修）
一、实验目的
通过本实验的学习，使学生掌握MATLAB软件的操作界面，系统帮助信息的获取方法。

掌握矩阵的操作方法，命令的输入方法。

掌握M文件的编辑、操作方法。

为以后的操作打下基础。

二、实验内容
（1）掌握MA TLAB软件的操作界面
（2）掌握M文件的编辑、操作方法
（3）数据的输入方法
三、实验原理、方法和手段
根据MATLAB命令的输入要求进行操作。

四、实验组织运行要求
采用集中授课形式。

五、实验条件
PⅣ计算机40台，MATLAB软件。

六、实验步骤
（一）基本功能
双击MATLAB图标，打开MATLAB Command Window，它是用户输入命令的地方， MATLAB将计算结果也显示在此。

共有，view，web, Windows, Help五个主要功能。

1 简易数学
>> 1+2+3
ans =
6
>> 1*10 + 2*20 + 3*30
ans =
140
>> x=1+2+3
x =
6
如果在上述的例子结尾加上;，则计算结果不会显示在命令窗口中，要得知计算值只须键入该变量名即可。

>> x=1+2+3;
>> x
x =
6
MATLAB提供基本的算术运算有：
加 (+)、减 (-)、乘 (*)、除 (/)、幂次方 (^)，例如：5+3, 5-3, 5*3, 5/3, 5^3
要计算面积Area = ,半径r = 2，则可键入
>> r=2;
>> area=pi*r^2;
>> area
12.5664
我们也可以将上述命令打在同一行，以, 或是; 分开，例如
>> r=2, area=pi*r^2
>> r=2; area=pi*r^2;
请注意上述二式的差异，前者有计算值显示，而后者无。

如果一个命令过长可以在结尾加上... ，例如
>> r=2;
>> area = pi ...
*r^2
另外一个符号注解是由%起头，也就是说在%之后的任何文字都被视为程序的注解。

例如>> r=2; % 键入半径
>> area=pi*r^2; % 计算面积
MATLAB可以将计算结果以不同的精确度的数字格式显示，在命令窗口键入以下显示格式的命令，以π值为例
命令数字值说明
format short 3.1416 预设的 4 位有效小数位数
format long 3.149 15 位有效小数位数
format short e 3.1416e+000
4 位有效小数位数加上指数表格式
观察下列命令后pi结果的变化：
>>format long
>>pi
>>format short
>>pi
2 变量
MATLAB对使用变量名称的规定：
1.变量名称的英文大小写是有区别的（apple, Apple, AppLe，三个变量不同）。

2.变量的长度上限为 19 个字符。

3.变量名的第一个字符必须是英文字符，随后可以英文字符、数字或下划线。

3其它功能
MATLAB利用了↑↓二个光标移动键将所执行过的命令重复使用。

按下↑前一次命令重新出现，之后再按Enter键，即再执行前一次的命令。

键入who可以查看所有定义过的变量名称。

而键入clear则是清除所有定义过的变量名称；如果只是要去除x及y 二个变量，则可以键入clear x y。

Ctrl-C（即同时按Ctrl及C二个键）可以用来中止执行中的MATLAB的工作。

4帮助
利用help命令，如果你要找题材 (topic)，直接键入help <topic>。

利用命令窗口的功能菜单中的Help，从中选取Table of Contents（目录）或是Index（索引）。

例如
>> help sqrt
SQRT Square root.
SQRT(X) is the square root of the elements of X. Complex
results are produced if X is not positive.
（二）数组与矩阵输入
1 数组与矩阵的定义
MATLAB的运算是以数组及矩阵方式，而这二者在MATLAB的基本运算性质不同，数组强调元素对元素的运算，而矩阵则采用线性代数的运算方式。

定义一变量为数组或是矩阵时，须用中括号[ ] 将元素置于其中。

数组为一维元素所构成，而矩阵为多维元素所组成，例如
>> x = [1 2 3] % 一维 1x3 数组
>> x = [1 2 3; 4 5 6] % 二维 2x3 矩阵
假设要计算y = sin (x), 0至π而x = 0, 0.2π, 0.4π,...,π，即可用数组方式运算，例如
>> x = [0 0.2*pi 0.4*pi 0.6*pi 0.8*pi pi] % 注意数组内也可作运算
x =
0 0.6283 1.2566 1.8850 2.5133 3.1416
>> y=sin(x)
y =
0 0.5878 0.9511 0.9511 0.5878 0.0000
要找出数组的某个元素或数个元素：
>> x(3) % 第三个x的元素
ans =
>> y(5) % 第五个y的元素
ans =
0.5878
>> x(1:5) % 列出第一到第五个x的元素
ans =
0 0.6283 1.2566 1.8850 2.5133
>> y(3:-1:1) % 列出第三到第一个y的元素，3为起始值，1为终止值，-1为增量
ans =
0.9511 0.5878 0
如果要建立的数组的元素多达数百个，则须采用以下的方式。

>> x=(0:0.2:1) % 以:区隔起始值=0、增量值=0.2、终止值=1
>> x=linspace(0,1,51) % 利用linspace，以区隔起始值=0终止值=1之间的元素数目=51 >> x=(0:0.01:1)*pi % 注意数组外也可作运算
>> a=1:5, b=1:2:9 % 这二种方式更直接
a =
1 2 3 4 5
b =
1 3 5 7 9
>> c=[b a] % 可利用先前建立的数组 a 及数组 b ，组成新数组
c =
1 3 5 7 9 1
2
3
4 5
2 数组运算符
以下将数组的运算符号及其意义列出，除了加减符号外其余的数组运算符号均须多加 . 符号。

数组运算功能
+ 加
- 减
.* 乘
./ 左除
.^ 次方
.' 转置
>> a=1:5; a-2 % 从数组a减2
ans =
-1 0 1 2 3
>> 2*a-1 % 以2乘数组a再减1
ans =
1 3 5 7 9
>> b=1:2:9; a+b % 数组a加数组b
ans =
2 5 8 11 14
>> a.*b % 数组a及b中的元素与元素相乘
1 6 15 28 45
>> a.^2 % 数组中的各个元素作二次方
ans =
1 4 9 16 25
3 特殊矩阵
zeros函数是形成元素皆为0 的矩阵；ones函数是形成元素皆为 1 的矩阵； eye则是产生一个单位矩阵，如zeros(m)可以产生一个m×m的方阵，而zeros(m,n)产生的是m×n的矩阵。

>> B=zeros(2,3)
B =
0 0 0
0 0 0
>> C=[1 2; 3 4; 5 6];
>> size(C) % 使用 size 命令得到C矩阵的大小
ans =
3 2
>> A=ones(2), B=ones(2,3) % 1 的矩阵
A =
1 1
1 1
B =
1 1 1
1 1 1
>> A=eye(2), B=eye(2,3) % 单位矩阵
A =
1 0
0 1
B =
1 0 0
0 1 0
（四）编写M-file
MATLAB提供了 M-file 的方式，可让使用者自行将命令及算式写成程序然后储存，其扩展各为m，如 test.m，其中的test就是文件名称。

在命令窗口中选择File再选择New，当程序写完后要存档时，必须以.m 名称储存。

以下的tutex1.m是一个简易绘图程序做为示范使用M-file
x=linspace(0,2*pi,20); y=sin(x);
plot(x,y,'r+')
xlabel('x-value')
ylabel('y-value')
title('2D plot')
写好上述程序后即可在命令窗口下键入tutex1，即可执行已建立的tutex1.m程序。

（五）设置工作目录
当在执行M-file时，我们最好是将自己的M-file储存在自己的工作目录下，而不要放在MATLAB内建的目录下，要在自己的工作目录执行程序可分为二个步骤：(1)建立搜寻路径，(2) 切换目录。

建立搜寻路径
MATLAB 将许多内建函数分门别类放在不同的子目录下，因此它在工作时须依次的搜寻这些目录，这个过程称为「搜寻路径」。

MATLAB的命令path可以让我们将自己的工作目录加在原来MATLAB 的搜寻路径之前或之后，如
先在D盘中创建文件夹“stu01”,然后输入下列命令：
>> path(path,'d:\stu01') % 将自己的目录 \stu01加在MATLAB的搜寻路径之后
>> path
>> path('d:\stu01',path) % 将自己的目录 \stu01加在MATLAB的搜寻路径之前
>> path
七、实验报告
要求实验前做好实验预习工作，熟悉本次实验的基本内容，重点和难点，实验过程中认真分析实验结果，写出详细的实验报告。

实验二：绘图
实验学时：2
实验类型：（验证）
实验要求：（必修）
一、实验目的
通过本实验的学习，使学生掌握熟悉绘图命令的操作方法。

了解三维的曲线绘图。

为以后的MATLAB绘图打下基础。

二、实验内容
（1）掌握二维绘图命令的操作方法
（2）掌握三维绘图命令的操作方法
（3）掌握图形控制方法
三、实验原理、方法和手段
根据绘图命令的要求进行操作。

四、实验组织运行要求
采用集中授课形式。

五、实验条件
PⅣ计算机40台，MATLAB软件。

六、实验步骤
（一）二维绘图
plot是用来绘函数x对函数y的二维图，例如要绘出 y = sin (x), 0至2π。

plot可以在一个图上绘数条曲线，且以不同的符号及颜色来表示曲线，如要在x、y轴及图上加注说明，则可利用命令xlabel, ylabel, title，使用命令grid 加上网格线。

______________________________________________
字母颜色标点线型
y 黄色·点线
m 粉红○圈线
c 亮蓝××线
r 大红＋＋字线
g 绿色－实线
b 蓝色*星形线
w 白色：虚线
k 黑色－.(--) 点划线
___________________________________________________
>> v1=linspace(0,2*pi,20); v2=sin(v1); % 建立 v1 及 v2 数组
>> plot(v1,v2) % 利用 plot，输入的变量为 x 轴， y 轴
>> v3=cos(v1); % 建立 v3 数组
>> plot(v1,v2,v1,v3) % 绘二条曲线
>> plot(v1,v2,v1,v2,'+') % 一样绘二条曲线，不过第二条曲线以符号 + 标示
>> plot(v1,v2,v1,v2.*v3,'--') % 绘二条曲线，一条代表 v1-v2 函数关系，一条
% 代表 v1-(v2.*v3) 函数关系
>> xlabel('x-axis') % 加上x轴的说明，在二个单引号'之间键入文字的说明
>> ylabel('y-axis') % 加上y轴的说明
>> title('2D plot') % 加上图的说明
gtext则是依据鼠标或上下左右游标键来放置文字说明，其语法为gtext('string')。

>> x=linspace(0,2*pi,30); y=sin(x); z=cos(x);
>> plot(x,y,x,z) % 绘二条曲线 y=sin(x), z=cos(x)
>> text(2.5,0.7,'sin(x)') % (2.5,0.7)是依据绘图大小的座标值
>> gtext('cos(x)') % 将鼠标移至适当位置再按鼠标键
一般的 x-y 图在横轴及纵轴皆是以线性尺度来绘图，如果要绘图的数据的 x 或 y 值变化范围太大，就须要改用对数 (log) 尺度来绘图才可得到合理的图。

MATLAB 提供三种对数尺度的绘图命令：semilogx,semilogy, loglog，它们的作用分别是x轴以对数尺度绘图，y轴以对数尺度绘图，x和y轴以对数尺度绘图。

>> y=0:0.1:10; x=10.^y
>> plot(x,y) %
>> semilogx(x,y) % 改以对数尺度绘图
>> x=[0 2 5 7 10 12 15 17 20 21];
>> y=[0.1 0.2 0.5 0.6 0.9 1 1.2 1.26 1.22 1.2];
>> plot(x,y) % 先以线性尺度绘图，再分别以三种对数尺度绘
>> semilogx(x,y) % 图，注意各个图像会改变
>> semilogy(x,y)
>> loglog(x,y)
（二）绘图选项
1 MATLAB 有许多的绘图选项，如将二条曲线划在同一个图中。

>> x=linspace(0,2*pi,30);
>> y=sin(x); z=cos(x);
>> plot(x,y,x,z) % 将 y=sin(x) 及 z=cos(x) 二函数分布绘图
>> plot(x,y,'g:',x,z,'r--') % 加上不同的颜色及符号来区别二条曲线title ——给图形加标题
xlable ——给x轴加标注
ylable ——给y轴加标注
text ——在图形指定位置加标注
gtext ——将标注加到图形任意位置
grid on(off) ——打开、关闭坐标网格线
legend ——添加图例
axis ——控制坐标轴的刻度
2横轴和纵轴的控制
要控制绘图的横轴及纵轴比例，可以用 axis配合下列的相关的选项：
axis([xmin xmax ymin ymax]) 以 xmin xmax 设定横轴的下限及上限，以
ymin ymax 设定纵轴的下限及上限
axis auto
横轴及纵轴依照数据大小的上下限来订定，横轴及纵轴比例是4:3
axis square 横轴及纵轴比例是 1:1
axis equal 将横轴纵轴的尺度比例设成相同值
axis normal 以预设值画纵轴及横轴
axis off 将纵轴及横轴取消
axis on 恢复纵轴及横轴
上述的各个命令的语法也可以将关键字放在括弧内的单引号之间，如axis(' ')。

>> x=linspace(0,2*pi,30); y=sin(x); z=cos(x);
>> plot(x,y,x,z)
>> axis off
>> axis on
>> axis('square','equal')
>> axis('xy','normal')
3 子图
要将数个相关的图画在同一页时，可以用subplot这个命令。

其语法为 subplot(m,n,p)，其中 m, n代表绘图成 m x n 个子图，m表示在 y方向有 m 个图， n表示在 x 方向有n个图，p 是代表第几个子图。

下例是以 subplot分别画出线性及对数尺度的四个子图：>> x=[0 2 5 7 10 12 15 17 20 21];
>> y=[0.1 0.2 0.5 0.6 0.9 1 1.2 1.26 1.22 1.2];
>> subplot(2,2,1), plot(x,y) % 画左上角的图
>> subplot(2,2,2), semilogx(x,y) % 画右上角的图
>> subplot(2,2,3), semilogy(x,y) % 画左下角的图
>> subplot(2,2,4), loglog(x,y) % 画右下角的图
4图形放大及缩小
zoom 命令可以将图形放大或缩小，若要将图形放大时用 zoom on，zoom out，当不再要放大或缩小图形时用 zoom off。

>> M=peaks(25); % peaks 是MATLAB内建的一个像山峰的特别函数，25是这个
>> plot(M) % 函数矩阵的大小，如果数值愈大则画出的山峰图愈平滑
>> zoom on % 开始放大图形，每按一次鼠标左键图形就放大一次
>> zoom out % 开始缩小图形，每按一次鼠标右键图形就缩小一次
>> zoom off % 停止图形放大或缩小功能
5函数分布的快速绘图
fplot的命令可以用来自动的画一个已定义的函数分布图，而无须产生绘图所须要的一组数据做为变量。

其语法为fplot('fun',[xmin xmax ymin ymax])，其中 fun为一已定义的函数名称，例如 sin, cos等等；而 xmin, xmax, ymin, ymax 则是设定绘图横轴及纵轴的下限及上限。

以下的例子是将一函数 f(x)=sin(x)/x在-20 x 20,-0.4 y 1.2之间画出：>> fplot('sin(x)./x',[-20 20 -0.4 1.2])
>> title('Fplot of f(x)=sin(x)/x')
>> xlabel('x'), ylabel('f(x)')
（三）三维绘图
1三维的曲线绘图
plot3 可以用来画一个三维的曲线，它的格式类似 plot ，不过多了 z方向的数据。

其与法可以是 plot3(X,Y,Z) 或是 plot3(X,Y,Z,'linetype')，其中的 linetype是设定画线的符号和颜色。

下面的例子说明一个三维的曲线图：
>> t=0:pi/50:10*pi;
>> plot3(sin(t),cos(t),t)
>> title('Helix'), xlabel('sin(t)', ylabel('cos(t)'), zlabel('t')
2曲面及等值线绘图
如果要画一个三维的曲面，MATLAB是以meshgrid配合与mesh或surf命令来绘图。

先要以meshgrid产生在x-y平面的二维的网格数据，再以一组z轴的数据对应到这个二维的网格，即可画出三维的曲面。

>> x=-7.5:0.5:7.5; y=x; % 先产生x及y二个数组
>> [X,Y]=meshgrid(x,y); % 再以meshgrid形成二维的网格数据
>> R=sqrt(X.^2+Y.^2)+eps; % 加上eps可避免当R在分母时趋近零时会无法定义
>> Z=sin(R)./R; % 产生z轴的数据
>> mesh(X,Y,Z) % 将z轴的变化值以网格方式画出
>> surf(X,Y,Z) % 将z轴的变化值以曲面方式画出
七、实验报告
要求实验前做好实验预习工作，熟悉本次实验的基本内容，重点和难点，实验过程中认真分析实验结果，写出详细的实验报告。

实验三：函数
实验学时：2
实验类型：（验证）
实验要求：（必修）
一、实验目的
通过本实验的学习，使学生掌握掌握MATLAB常用函数的功能。

熟悉多项式的表示方法，多项式的运算，数据分析函数，关系及逻辑运算，用户自定义函数，矩阵运算函数，for 循环, while 循环的使用方法。

二、实验内容
（1）多项式的表示方法，多项式的运算函数
（2）数据分析函数
（3）矩阵运算函数
三、实验原理、方法和手段
根据MATLAB函数的使用方法进行操作。

四、实验组织运行要求
采用集中授课形式。

五、实验条件
PⅣ计算机40台，MATLAB软件。

六、实验步骤
（一）多项式函数
1多项式的表示方法
在生命科学实验分析中，多项式常被用来模拟一个现象的函数。

令p(x) 代表一个多项式：
MATLAB 以一最简便方式代表上述的多项式 p=[1 4 -7 -10]，其中的数值是多项式的各阶项（从高到低）的各个系数，其实p 也是一个数组，用以代表这个多项式。

2 多项式的加减乘除运算
有了多项式的表示式后，即可来计算其函数值。

可以用函数 polyval直接做运算，语法为polyval(p,x)，其中p 代表多项式各阶系数的数组。

因此：
>> x=linspace(-1,3);
>> p=[1 4 7 -10];
>> v=polyval(p,x);
函数conv做乘法运算以及函数deconv做除法运算。

当二多项式相乘，其语法为conv(a,b)，其中a, b代表二个多项式的数组。

而二多项式相除有deconv函数，其语法稍有不同
[q,r]=deconv(a,b)，其中q,r分别代表整除多项式及余数多项式。

>> a=[1 2 3 4]; b=[1 4 9 16];
>> c=a+b
c =
2 6 12 20
>> d=a-b
d =
0 -2 -6 -12
>> e=conv(a,b)
e =
1 6 20 50 75 84 64
>> g=e+[0 0 0 c]
g =
1 6 20 5
2 81 96 84
>> [f,r]=deconv(e,b)
f =
1 2 3 4
r =
0 0 0 0 0 0 0 % 因为是整除所以余数多项式的各系数皆为零
>> [h,r]=deconv(g,a)
f =
1 4 9 18
r =
0 0 0 0 2 6 12 % 余数多项式为 2*x^2 + 6*x + 12
（二）数据分析函数
1极值、平均、总和、连乘及排序
最大值max，最小值min，平均值 mean，一组数据的中位数median，总和值sum，连乘值prod，累积总和值cumsum，累积连乘值cumprod，排序函数sort。

max(x) 找出x数组的最大值
max(x,y) 找出x及y数组的最大值，会有二个极值分属x及y数组
[y,i]=max(x) 找出x数组的最大值以y显示，其在x数组的位置以i显示
min(x) 找出x数组的最小值
min(x,y) 找出x及y数组的最小值，会有二个极值分属x及y数组
[y,i]=min(x) 找出x数组的最小值以y显示，其在x数组的位置以i显示
mean(x) 找出x数组的平均值
median(x) 找出x数组的中位数
sum(x) 计算x数组的总和值
prod(x) 计算x数组的连乘值
cumsum(x) 计算x数组的累积总和值
cumprod(x) 计算x数组的累积连乘值
>> rains=1000*rand(2,6) % rains为一个2x6的数组
rains =
126.8 148.5 173.0 148.4 194.7 208.9
328.8 300.7 268.3 210.5 278.4 321.5
>> avg_rain=mean(rains) % 将rains数组中的每一行的平均值列出
avg_rain =
227.8000 224.6000 220.6500 179.4500 236.5500 265.2000
>> avg_rain=mean(avg_rain) % 将上述数组中的平均值列出
avg_rain =
225.7083
>> max_rain=max(rains) % 将rains数组中的每一行的最大值列出
max_rain =
328.8000 300.7000 268.3000 210.5000 278.4000 321.5000
>> [max_rain,x]=max(rains) % 将rains数组中的每一行的最大值及其位置列出
max_rain =
328.8000 300.7000 268.3000 210.5000 278.4000 321.5000
x =
2 2 2 2 2 2
>> min_rain=min(rains) % 将rains数组中的每一行的最小值列出
min_rain =
126.8000 148.5000 173.0000 148.4000 194.7000 208.9000
>> s_sort=sort(rains) % 将rains数组的值由小到大做排序
s_sort =
126.8000 148.5000 173.0000 148.4000 194.7000 208.9000
328.8000 300.7000 268.3000 210.5000 278.4000 321.5000
>> x=[1 2 3 4 5];
>> sum(x) % 将x数组的值做总和
ans =
15
>> prod(x) % 将x数组的值做连乘
ans =
120
>> cumsum(x) % 将x数组的值累积后做总和
ans =
1 3 6 10 15
>> cumprod(x) % 将x数组的值累积后做连乘
ans =
1 2 6 24 120
（三）选择命令及函数
关系及逻辑运算
在执行关系及逻辑运算时，MATLAB 将输入的不为零的数值都视为真 (True)而为零的数值则视为否 (False)。

运算的输出值将判断为真者以 1 表示，而判断为否者以 0 表示。

MATLAB 提供以下的关系判断及逻辑的运算符：
符号关系的意义
< 小于
<= 小于等于
> 大于
>= 大于等于
== 等于
~= 不等于
上述的各个运算符需用在二个大小相同的数组或是矩阵的比较。

>> a=1:5, b=5-a,
a =
1 2 3 4 5
b =
4 3 2 1 0
>> tf= a>4
tf =
0 0 0 0 1
>> tf= a==b
tf =
0 0 0 0 0
>> tf= b-(a>2)
tf =
4 3 1 0 -1
>> tf= ~(a>4)
tf =
1 1 1 1 0
>> tf= (a>2)&(a<6)
tf =
0 0 1 1 1
以下是算式利用关系及逻辑运算产生一不连续的讯号。

>> x=linspace(0,10,100); % 产生数据
>> y=sin(x); % 产生 sine 函数
>> z=(y>=0).*y; % 将 sin(x) 的负值设为零
>> z=z + 0.5*(y<0); % 再将上式的值加上0.5
>> z=(x<8).*z; % 将大于 x=8 以后的值设为零
>> hold on
>> plot(x,z)
>> xlabel('x'),ylabel('z=f(x)')
>> title('A discontinuous signal')
>> hold off
（四）用户自定义函数
我们可以定义一函数cirarea计算圆的面积，以下的 M-file: cirarea.m就是定义这个函数
function c=cirarea(r)
c=pi*r.^2;
M-file定义的函数有其语法的一些规定：
第一行命令以function开始，接着是输出的变量，等号，函数名称，输入的变量是接着函数名称放在括号之内。

如：function out1=userfun(in1)，这行的out1是输出的变量，userfun
是函数名称，in1是输入的变量。

function [out1, out2]= serfun(in1, in2) 如果输出变量[out1,out2] 和输入变量 (in1, in2)不只一个时，则在输出变量部分需加上[ ]。

（五）随机数
随机数依其统计分布特性分为：均匀(uniform) 随机数，常态 (normal) 随机数。

均匀随机数是指其值平均分布于一区间，而常态随机数的值则呈现高斯分布，形状像一个中间高二头低的山丘。

1均匀随机数
用 MATLAB 函数 rand产生在区间 [0, 1] 的均匀随机数，它是平均分布在 [0, 1]之间。

均匀随机数函数的语法为rand(n), rand(m,n)，其结果分别产生一矩阵含nxn个随机数和一矩阵含mxn的随机数。

注意每次产生随机数的值都不会一样，这些值代表的是随机且不可预期的，这正是我们用随机数的目的。

>> rand(1,6) % 第一次使用随机数产生器
ans =
0.2190 0.0470 0.6789 0.6793 0.9347 0.3835
>>hist(ans)
>>plot(ans)
>> rand(1,6) % 第二次使用随机数产生器，注意每次产生的随机数值皆不同
ans =
0.5194 0.8310 0.0346 0.0535 0.5297 0.6711
2常态随机数
用 MATLAB 函数 randn产生常态随机数，它是以高斯分布在随机数出现的上下限区间。

randn(n)和randn(n,m) 是分别产生一矩阵含 nxn个随机数和一矩阵含mxn的常态随机数，其平均值为0变异数为1。

>> x=-2.9:0.2:2.9; % 这个例子用到 hist 函数来画出二种随机数的分布图
>> y=randn(1,5000);
>> hist(y,x)
>> title('Histogram of Normal Random Data')
>> y1=rand(1,5000);
>> hist(y1,x)
>> title('Histogram of Uniform Random Data')
七、实验报告
要求实验前做好实验预习工作，熟悉本次实验的基本内容，重点和难点，实验过程中认真分析实验结果，写出详细的实验报告。

实验四：解非线性方程
实验学时：2
实验类型：（验证）
实验要求：（必修）
一、实验目的
通过实验，使学生ezplot函数、fzero函数、roots函数．fminbnd的使用方法，求解非线性方程和非线性方程组。

要求学生能独立操作完成全部实验。

二、实验内容
（1）ezplot函数
（2）fzero函数
（3）fminbnd函数
三、实验原理、方法和手段
应用MATLAB函数进行矩阵的操作。

四、实验组织运行要求
采用集中授课形式。

五、实验条件
PⅣ计算机40台，MATLAB软件。

六、实验步骤
（一）多项式计算
1 多项式的表示
>>p =[1, -6,-72, -27]
p是多项式p(x)=x3-6x2-72x-27的matlab描述方法，我们可用：
>>p1=poly2str(p,‘x’) —函数文件，显示数学多项式的形式
p1 =x^3 - 6 x^2 - 72 x – 27
2多项式的计算
>> p =[1, -6,-72, -27]
>>r=roots(p)
ans =
12.1229
-5.7345
-0.3884
>>roots([1 0 -15 14])
ans =
-4.2749
3.2749
1.0000
（二）一般求解
求解非线性方程的函数fzero
fzero是求单变量函数的零点（即方程之根）。

格式：options=optimset(‘display’,’iter’,’tolx’,delta)
z=fzero(‘fun’,x,options),
delta是一个正数，收敛容差，默认值为eps。

disply 可选取iter显示每一步的迭代结果,final只显示最终结果,off不打印信息。

fun为单变量实函数f(x)的文件名，x为初始估计值。

如果函数不存在零点，则返回NaN。

X是长度为2的向量，用于指定间隔区间，并且要求f(x(1))和f(x(2))不同号，当x为标量时，则x作为起点。

例：计算sin(x)从x=3开始的零点，命令：x=fzero(‘sin(x)’,3)。

tol指定了容许的误差容限。

trace显示每次的迭代情况。

例如：计算方程x^3-15*x+14=0在3附近的根。

命令：
>>fzero('x^3-15*x+14',3,optimset('display','iter'))
ans =
3.2749
>>fplot(‘x^3-15*x+14',[-5,5])
（三）极值计算
1单变量函数最小值
x=fminbnd(‘fun’,x1,x2),求出指定函数fun在区间[x1,x2]中的局部最小点。

>>f=fminbnd('x^2+3*x+2',-5,5)
f =
-1.5000
>>[x,fval,exitflag]=fminbnd('-(3-2*x)^2*x',0,1.5)
x =
0.5000
fval =
-2.0000
exitflag =
1
exitflag>0 收敛
exitflag=0 已达设定的次数
exitflag<0 不收敛
>>fplot('(3-2*x)^2*x',[0,1.5])
2 多变量函数最小值
计算函数f=-(sin(x)+sin(y))-x/10在区间[x(-2，2),y(-2,2)]的最小值。

编写M文件
function f=myfun(x)
f=-(sin(x(1))+sin(x(2)))-x(1)/10
>> [x,fval,exitflag]=fminsearch(@myfun,[1,1])
x =
1.6710 1.5708
fval =
-2.1621
exitflag =
1
绘制图形：
>>x=-2:0.1:2, y=x
>> [X,Y]=meshgrid(x,y),
>>Z=(sin(X)+sin(Y))-X/10
>>mesh(X,Y,Z)
七、实验报告
要求实验前做好实验预习工作，熟悉本次实验的基本内容，重点和难点，实验过程中认真分析实验结果，写出详细的实验报告。

实验五：解线性方程组
实验学时：2
实验类型：（验证）
实验要求：（必修）
一、实验目的
生命科学数学模型中有大量的线性方程组的求解问题，通过本实验的学习，使学生掌握解线性方程组的MATLAB 运算方法。

掌握线性方程组的seidel 迭代法的各程序设计方法。

用来解决生命科学中的实际问题。

二、实验内容
（1）利用矩阵左除 \ 做运算求解
（2）以逆矩阵运算
（3）seidel 迭代法求解
三、实验原理、方法和手段
利用高斯消去法原理进行线性方程组的直接计算，利用seidel 迭代法设计线性方程组的迭代算法。

四、实验组织运行要求
采用集中授课形式。

五、实验条件
P Ⅳ计算机40台，MATLAB 软件。

六、实验步骤
对于方程ax ＝b ，a 为a n ×m 矩阵，有三种情况：
当n=m 时，此方程成为“恰定”方程
当n>m 时，此方程成为“超定”方程
当n<m 时，此方程成为“欠定”方程
恰定方程组的解
方程ax ＝b(a 为非奇异)
x=a -1 b
两种解:
x=inv(a)*b — 采用求逆运算解方程
x=a\b — 采用左除运算解方程 （一）利用矩阵解法
(1)1()1()()k k +--=-+-X D L UX D L d
假设一线性方程组为
我们将以上方程组以矩阵方式表示：AX=B
其中 A 为等式左边各方程的系数项，X为要求解的未知项，B代表等式右边之已知项。

要解上述的线性方程组，我们可以利用矩阵左除运算符 \ 做运算，即：X=A\B。

若以逆矩阵运算求解AX=B, X=B，即是X=inv(A)*B。

>> A=[3 2 -1; -1 3 2; 1 -1 -1]; % 将等式的左边系数键入
>> B=[10 5 -1]'; % 将等式右边之已知项键入，B要做转置
>> X=A\B % 左除运算符求解
X = % 注意X为行向量
-2
5
-6
>> C=A*X % 验算解是否正确
C = % C=B
10
5
-1
>> X=inv(A)*B; %逆矩阵运算求解
(二)迭代法求解
利用seidel迭代法编写MATLAB程序seidel.m
function s=seidel(a,d,x0)
D=diag(diag(a))
u=-triu(a,1)
l=-tril(a,-1)
c=inv(D-l)
b=c*u
g=c*d
s=b*x0+g
n=1
while max(s-x0)>=1.0e-5
x0=s
s=b*x0+g
n=n+1
end
n,s
输入以下命令：
>>a=[10,-1,-2;-1,10,-2;-1,-1,5]
>> d=[7.2 8.3 4.2]'
>>x0=[0,0,0]’
>>seidel(a,d,x0)
>>n =
8
s =
1.1000
1.2000
1.3000
七、实验报告
要求实验前做好实验预习工作，熟悉本次实验的基本内容，重点和难点，实验过程中认真分析实验结果，写出详细的实验报告。

实验六：插值计算
实验学时：2
实验类型：（验证）
实验要求：（必修）
一、实验目的
生命科学数学模型中有大量的插值和拟合问题，通过本实验的学习，使学生掌握掌值中interp1一维插值和spline 样条插值的操作方法。

用来解决生命科学中的实际问题。

二、实验内容
（1）interp1一维插值
（2）interp2二维插值
（3）spline 样条插值
三、实验原理、方法和手段
应用MATLAB函数interp1，spline，interp2进行插值操作，根据函数参数的要求进行操作。

四、实验组织运行要求
采用集中授课形式。

五、实验条件
PⅣ计算机40台，MATLAB软件。

六、实验步骤
假设一组已知的数据为
假设某些点并不属于上述的
但是
，如果我们要估计这些点的函数值
，就要做内插 (interpolation)。

（一）一维内插
MATLAB的一维内插函数是interp1，其语法为
interp1(x,y,xi),interp1(x,y,xi,'method')；其中的x,y是原已知的数据的x,y值，而xi 则是要内插的数据点，另外method可以设定内插方法有 linear,cubic,spline，分别是一次、三次方程和spline函数，其中预设方法是linear。

如果数据的变化较大，以 spline函数内插所形成的曲线最平滑，所以效果最好。

而三次方程序得到的内插曲线平滑度介于线性与spline函数之间。

假设某生物工程反应器使用时间（单位为sec）与反应器内部温度的测量值如下
时间 0 1 2 3 4 5
温度 0 20 60 68 77 110
其中温度的数据从 20oC变化到 110oC，如果要估计在t=2.6, 4.9 sec的反应器内部温度，可用下列命令：
>> x=[0 1 2 3 4 5]'; % 键入时间
>> y=[0 20 60 68 77 110]'; % 键入温度
>> y1=interp1(x,y,2.6) % 要内插的数据点为2.6
y1 = % 对应 2.6 的函数值为 64.8
64.8
>> y1=interp1(x,y,[2.6 4.9]) % 内插数据点为 2.6, 4.9，注意用[ ]将多个内插点放在其中
y1 =
64.8
106.7
>> y1=interp1(x,y,2.6,'cubic') % 以三次方程对数据点 2.6 作内插
y1 = % 对应2.6 的函数值为 66.264
66.264
>> y1=interp1(x,y,2.6,'spline') % 以spline函数对数据点2.6 作内插
y1 = % 对应2.6 的函数值为66.368
66.368
以下的例子还配合绘图功能，用以比较不同内插方法的差异。

>> h=1:12;
>> temp=[5 8 9 15 25 29 31 30 22 25 27 24]; % 这组温度数据变化较大
>> plot(h,temp,'--',h,temp,'+') % 将线性内插结果绘图
>> h_3=1:0.1:12 % 要每0.1小时估计一次温度值
>> t_3=interp1(h,temp,h_3,'cubic') % 以三次方程做内插
>> t_s=interp1(h,temp,h_3,'spline') % 以spline函数做内插
>> hold on
>> subplot(1,2,1)
>> plot(h,temp,'--',h,temp,'+',h_3,t_3) % 将线性及三次方程内插绘图
>> subplot(1,2,2)
>> plot(h,temp,'--',h,temp,'+',h_3,t_s) % 将线性方程及spline内插绘图
>> hold off
（二）二维内插
二维数据内插interp2
格式：zi=interp2(x,y,z,xi,yi,method),给定数据对（x,y,z）,可求出对应于(xi,yi)对应点处的zi值。

（三）Spline 内插
spline内插可以用interp1指定内插方式为spline来做。

另一种方式也可以用spline(x,y,xi)来做，其中的x,y,xi的用法与interp1中的语法相同。

>> x=[0 1 2 3 4 5]';
>> y=[0 20 60 68 77 110]';
>> y1=spline(x,y,2.6)
y1 =
67.3
>> y1=spline(x,y,[2.6,4.9])
y1 =
67.3 105.2
区别不同插值计算
>>x=1:7
>>y=[1 3 0 20 20 4 18]
>>x1=1:0.2:7
>>y1=interp1(x,y,x1,’cubic’)
>>y2=spline(x,y,x1)
>>plot(x,y,’*’,x1,y1,’r’,x1,y2,’-.’)。