数据分析与多项式计算

5.1 数据统计分析

1.求矩阵的最大元素和最小元素

max()：求向量或矩阵的最大元素。

min()：求向量或矩阵的最小元素。

当参数为向量时，上述函数有两种调用格式：

（ 1）

y=max(X)

：返回向量X的最大值存入y，如果X中包含复数元素，则按模取最大值。

（ 2）

[y,k]=max(X)

：返回向量X的最大值存入y，最大值元素的序号存入k，如果X中包含复数元素，则按模取最大值。

例1 求向量x的最大元素，其中x=[-43,72,9,16,23,47]。

x=[-43,72,9,16,23,47]; 
>> y=max(x) 
   y= 72 
>> [y,k]=max(x) 
   y= 72 
   k= 2
           

当参数为矩阵时，函数有三种调用格式：

（ 1） max(A)：返回一个行向量，向量的第i个元素是矩阵A的第i列上的最大值。

（ 2） [Y,U]=max(A)：返回行向量Y和U， Y向量记录A中每列的最大值， U向量记录每列最大值元素的行号。

（ 3） max(A,[],dim)： dim取1或2。 dim取1时，该函数的功能和max(A)完全相同； dim取2时，该函数返回一个列向量，其第i个元素是A矩阵的第i行上的最大值。

例2 求矩阵A的每行及每列的最大元素，并求整个矩阵的最大元素。

A=[13,-56,78;25,63,-235;78,25,563;1,0,-1]; max(A) 
ans= 
     78 63 563 
max(A,[],2) 
ans= 
     78 63 563 1 
max(max(A)) 
ans=
     563
           

2.求矩阵的平均值和中值

 平均值：指算术平均值，即每项数据之和除以项数。

 中值：指在数据序列中其值的大小恰好处在中间的元素。如果数据个数为奇数，则取值为大小位于中间的元素；如果数据个数为偶数，则取中间两个元素的平均值。

在MATLAB中，求平均值和中值的函数分别为：

 mean()：求算术平均值。

 median()：求中值。

3.求和与求积

 sum():求和函数。

 prod():求积函数。

4.累加和与累乘积

在MATLAB中，求累加和与累乘积的函数分别为：

 cumsum():累加和函数。

 cumprod():累乘积函数。

5.标准差与相关系数

标准差用于计算数据偏离平均数的距离的平均值，其计算公式为

MATLAB中计算标准差的函数为std()，其调用格式为：

（ 1） std(X)：计算向量X的标准差。

（ 2） std(A)：计算矩阵A的各列的标准差。

（ 3） std(A,flag,dim)： flag取0或1，当flag=0时，按S1所列公式计算样本标准方差；当flag=1时，按S2所列公式计算总体标准方差。在默认情况下， flag=0， dim=1。

相关系数能够反映两组数据序列之间相互关系，其计算公式为

在MATLAB中，计算相关系数的函数为corrcoef()，其调用格式为：

（ 1） corrcoef(A)：返回由矩阵A所形成的一个相关系数矩阵，其中，第i行第j列的元素表示原矩阵A中第i列和第j列的相关系数。

（ 2） corrcoef(X,Y)：在这里， X,Y是向量，它们与corrcoef([X,Y])的作用一样，用于求X、 Y向量之间的相关系数。

6.排序

在MATLAB中，排序函数为sort()，其调用格式为：

（ 1） sort(X)：对向量X按升序排列。

（ 2） [Y,I]=sort(A,dim,mode)，其中dim指明对A的列还是行进行排序。 mode

指明按升序还是降序排序，若取“ ascend” ，则按升序； 若取“ descend” ，则按降序，默认为升序。输出参数中， Y是排序后的矩阵，而I记录Y中的元素在A中位置。

5.2 多项式计算

多项式的表示

在MATLAB中创建多项式向量时，要注意三点：

（1）多项式系数向量的顺序是从高到低。

（2）多项式系数向量包含0次项系数，所以其长度为多项式最高次数

加1。

（3）如果有的项没有，系数向量相应位置应用0补足。

多项式的四则运算

(1）多项式的加减运算

多项式的加减运算非常简单，即相应向量相加减。

(2）多项式乘法

conv(P1,P2)

：多项式相乘，其中， P1、 P2是两个多项式系数向量。

(3）多项式除法

[Q,r]=deconv(P1,P2)

：多项式相除，其中Q返回多项式P1除以P2的商式，r返回P1除以P2的余式。

这里，Q和r仍是多项式系数向量。deconv是conv的逆函数，因此有下式成立。

P1=conv(Q,P2)+r
           

多项式的求导

在MATLAB中，多项式求导函数为

polyder()

，其调用格式为：

（1） p=polyder( P )：求多项式P的导函数。

（2） p=polyder(P,Q)：求P× Q的导函数。

（3） [p,q]=polyder(P,Q)：求P/Q的导函数，导函数的分子存入p，分母存入q。

多项式的求值

 polyval(p,x)：代数多项式求值，其中， p为多项式系数向量， x可以是标量、向量或矩阵。若x为标量，则求多项式在该点的值；若x为向量或矩阵，则对向量或矩阵中的每个元素求多项式的值。

 polyvalm(p,x)：矩阵多项式求值，其调用格式与polyval相同，但含义不同。 polyvalm函数要求x为方阵，以方阵为自变量求多项式的值。

设A为方阵， P代表多项式x3-5x2+8，那么polyvalm(P,A)的含义是：

A*A*A-5*A*A+8*eye(size(A))
           

而polyval(P,A)的含义是：

A.*A.*A-5*A.*A+8*ones(size(A))
           

多项式的求根

在MATLAB中，多项式求根函数为roots( p )，其中p为多项式的系数向量。

例如，求多项式 x4+8x3-10的根。

>> a=[1,8,0,0,-10];
>> x=roots(a)
x =
-8.0194 + 0.0000i
1.0344 + 0.0000i
-0.5075 + 0.9736i
-0.5075 - 0.9736i
           

若已知多项式的全部根，则可以用poly函数建立起该多项式，其调用格式为：

p=poly(x)
           

例如，在以上语句之后，继续输入命令：

>> p=poly(x)
p=
1.0000 8.0000 0.0000 0.0000 -10.0000 =a
           

5.3 数据插值

引例-零件加工问题

在MATLAB中，一维插值函数为interp1()，其调用格式为：

Y1=interp1(X,Y,X1,method)
           

该语句将根据X、Y的值，计算函数在X1处的值。其中，X、Y是两个等长的已知向量，分别表示采样点和采样值。

X1是一个向量或标量，表示要插值的点。

method参数用于指定插值方法， 常用的取值有以下四种：

（ 1） linear：线性插值， 默认方法。 将与插值点靠近的两个数据点用直线连接， 然后在直线上选取对应插值点的数据。

（2） nearest：最近点插值。 选择最近样本点的值作为插值数据。

（ 3） pchip：分段3次埃尔米特插值。采用分段三次多项式，除满足插值条件，还需满足在若干节点处相邻段插值函数的一阶导数相等，使得曲线光滑的同时，还具有保形性。

（ 4）spline：3次样条插值。每个分段内构造一个三次多项式,使其插值函数除满足插值条件外， 还要求在各节点处具有连续的一阶和二阶导数。

四种方法的比较

 线性插值和最近点插值方法比较简单。 其中线性插值方法的计算量与样本点n无关。 n越大， 误差越小。

 3次埃尔米特插值和3次样条插值都能保证曲线的光滑性。相比较而言，3次埃尔米特插值具有保形性；而3次样条插值要求其二阶导数也连续，所以插值函数的性态更好。

MATLAB中的二维插值函数为interp2()，其调用格式为：

Z1=interp2(X,Y,Z,X1,Y1,method)

其中，X、Y是两个向量，表示两个参数的采样点，Z是采样点对应的函数值。X1、Y1是两个标量或向量，表示要插值的点。

应用案例-粮储仓的通风控制问题

>> x=20:10:90;
>> y=(0:5:20)';
>> z=[8.9,10.32,11.3,12.5,13.9,15.3,17.8,21.3;
8.7,10.8,11,12.1,13.2,14.8,16.55,20.8;
8.3,9.65,10.88,12,13.2,14.6,16.4,20.5;
8.1,9.4,10.7,11.9,13.1,14.5,16.2,20.3;
8.1,9.2,10.8,12,13.2,14.8,16.9,20.9];
>> xi=20:90;
>> yi=(0:20)';
>> zi=interp2(x,y,z,xi,yi,'spline');
>> surf(xi,yi,zi)
           

5.3曲线拟合

x=1790:10:2010;
y=[3.9,5.3,7.2,9.6,12.9,17.1,23.2,31.4,38.6,50.2,63.0,76.0,92.0,
105.7,122.8,131.7,150.7,179.3,203.2,226.5,248.7,281.4,308.7];
plot(x,y,'*');
p=polyfit(x,y,3);
polyval(p,2020)
plot(x,y,'*',x,polyval(p,x));
           

最小二乘法（又称最小平方法）是一种数学优化技术。它通过最小化误差的平方和来寻找数据的最佳函数匹配。

曲线拟合的实现方法

MATLAB中的多项式拟合函数为polyfit()，其功能为求得最小二乘拟合多项式系数，其调用格式为：

（ 1） P=polyfit(X,Y,m)

（ 2） [P,S]=polyfit(X,Y,m)

（ 3） [P,S,mu]=polyfit(X,Y,m)：根据样本数据X和Y，产生一个m次多项式P及其在采样点误差数据S，mu是一个二元向量，mu(1)是mean(X)，而mu(2)是std(X)。

在引例中，我们已经用polyfit()函数预测了某国2020年的人口数。这个结果是否正确呢？我们无法得到2020年的数据，但是2016年的数据已经有了。所以，不妨再预测一下该国2016年的人口数。

polyval(p,2016)

ans =

327.0964

该国2016年人口数实际为323.1

问题分析： 据研究，一个国家的人口增长具有如下特点：

（1） 发展越平稳， 人口增长越有规律。 （2） 当经济发展到一定水平时，

人口增长率反而下降。换言之，在不同的环境和经济发展水平，人口可能有不同的增长规律。

结论：在人口增长数据的拟合上，应该将二战后至今这一时期的数据与此前的数据分开处理

>> x=1950:10:2010;
>> y=[150.7,179.3,203.2,226.5,248.7,281.4,308.7];
>> p=polyfit(x,y,3)
p =
1.0e+06 *
0.0000 -0.0000 0.0034 -2.2305
>> p=polyfit(x,y,2);
>> plot(x,y,'*',x,polyval(p,x))
>> polyval(p,2016)
ans =
325.1696
>> polyval(p,2020)
ans =
336.7857
           

结论：

（1） 要对问题的背景进行详细的分析。 （2） 采样点并非越多越好，适当的时候，可以减少采样点，分段进行拟合。