第八章相关分析
协方差
命令：C = cov(X)
当X为行或列向量时，它等于 var(X) 样本标准差。
X=1:15;cov(X)
ans =
20
>> var(X)
ans =
20
当X为矩阵时，此时X的每行为一次观察值，每列为一个变量。cov（X）为协方差矩阵,它是对称矩阵。
例：x=rand(100,3);c= cov(x)
c=
0.089672 -0.012641 -0.0055434
-0.012641 0.07928 0.012326
-0.0055434 0.012326 0.082203
c的对角线为：diag(c)
ans =
0.0897
0.0793
0.0822
它等于：var(x)
ans =
0.0897 0.0793 0.0822
sqrt(diag(cov(x)))
ans =
0.2995
0.2816
0.2867
它等于：std(x)
ans =
0.2995 0.2816 0.2867
命令：c = cov(x,y)
其中x和y是等长度的列向量（不是行向量），它等于cov([x y])或cov([x,y])
例：x=[1;4;9];y=[5;8;6];
>> c=cov(x,y)
c =
16.3333 1.1667
1.1667
2.3333
>> cov([x,y])
ans =
16.3333 1.1667
1.1667
2.3333
COV(X)、 COV(X,0)[两者相等] 或COV(X,Y)、COV(X,Y,0) [两者相等]，它们都是除以n-1，而COV(X,1) or COV(X,Y,1)是除以n
x=[1;4;9];y=[5;8;6];
>> cov(x,y,1)
ans =
10.8889 0.7778
0.7778 1.5556
它的对角线与var([x y],1) 相等
ans =
10.8889 1.5556
协差阵的代数计算：
[n,p] = size(X);
X = X - ones(n,1) * mean(X);
Y = X'*X/(n-1); Y为X的协差阵
相关系数（一） 命令：r=corrcoef(x)
x 为矩阵，此时x 的每行为一次观察值，每列为一个变量。 r 为相关系数矩阵。它称为Pearson 相关系数
例：x=rand(18,3);r=corrcoef(x) r =
1.0000 0.1509 -0.2008 0.1509 1.0000 0.1142 -0.2008 0.1142 1.0000 r 为对称矩阵，主对角阵为1 命令：r=corrcoef(x,y)
其中x 和y 是等长度的列向量（不是行向量），它等于cov([x y])或cov([x,y])，或x 和y 是等长度的行向量，r=corrcoef(x,y)它则等于r=corrcoef(x ’,y ’),
r=corrcoef([x ’,y ’])
例：x=[1;4;9];y=[5;8;6]; corrcoef(x,y) ans =
1.0000 0.1890 0.1890 1.0000
()()()
()
y y
x x y y x x
2
2
∑∑∑-?
---=r
corrcoef([x,y])
ans =
1.0000 0.1890
0.1890 1.0000
C = COV(X)
Rij =C(i,j)/SQRT(C(i,i)*C(j,j))
如：X=[1 2 7 4 ;5 12 7 8;9 17 11 17]; cov(X)
ans =
16.0000 30.0000 8.0000 26.0000 30.0000 58.3333 13.3333 46.6667 8.0000 13.3333 5.3333 14.6667 26.0000 46.6667 14.6667 44.3333 corrcoef(X)
ans =
1.0000 0.9820 0.8660 0.9762 0.9820 1.0000 0.7559 0.9177 0.8660 0.7559 1.0000 0.9538 0.9762 0.9177 0.9538 1.0000 则有：30/sqrt(16*58.3333)
ans =
0.9820
命令：[r,p]=corrcoef[….]
它还将返回p值，原假设是变量之间不相关。
例：
x = [430 335 520 490 470 210 195 270 400 480];
y=[30 21 35 42 37 20 8 17 35 25];
[r,p]=corrcoef(x,y)
r =
1.0000 0.8594
0.8594 1.0000
p =
1.0000 0.0014
0.0014 1.0000
P 矩阵主对角矩阵全为1，
当总体变量X 和Y 都服从正态分布，并且总体相关系数等于0时，有：
)
2(~122
---=
n t r
n r t
P 矩阵的计算，即上例中0.0014的算法。
4.7540
8594
.012
108594.0122
2
=--?=
--=
r
n r t (1-tcdf(4.7540,8))*2 得：0.0014
在显著性水平0.05下，0.0014小于0.05，拒绝两总体不相关的原假设，即销售量与气温相关。
命令：[r,p,rlo,rup]=corrcoef(….)
rlo 与rup 是与r 矩阵大小相同的矩阵，rlo 为相关系数r 的下限，rup 为相关系数r 的上限。在缺失情况下，置信度为95%。
例：x = [430 335 520 490 470 210 195 270 400 480];
y=[30 21 35 42 37 20 8 17 35 25];
[r,p,rlo,rup]=corrcoef(x,y)
r =
1.0000 0.8594
0.8594 1.0000
p =
1.0000 0.0014
0.0014 1.0000
rlo =
1.0000 0.5006
0.5006 1.0000
rup =
1.0000 0.9662
0.9662 1.0000
因此销售量与气温相关系数95%的置信区间为
[0.5006，0.9662]
如果要求99%的销售量与气温相关系数的置信区间：[r,p,rlo,rup]=corrcoef(x,y,'alpha',0.01)
注意不是：[r,p,rlo,rup]=corrcoef(x,y,0.01)
r =
1.0000 0.8594
0.8594 1.0000
p =
1.0000 0.0014
0.0014 1.0000
rlo =
1.0000 0.3071
0.3071 1.0000
rup =
1.0000 0.9786
0.9786 1.0000
得销售量与气温相关系数99%的置信区间为
[0.3071，0.9786]比95%的置信区间[0.5006，0.9662]
更宽。
rlo和rup具体怎么算出来的还没弄明白。帮助文件是这样说的：
The confidence bounds are based on an asymptotic normal distribution of 0.5*log((1+R)/(1-R)), with an approximate variance equal to 1/(n-3).
命令：[...]=corrcoef(...,'param1',val1,'param2',val2,...) 例：[...]=corrcoef(...,'alpha',0.01) 是求在置信度为99%，求r置信区间。
[...]=corrcoef(...,'rows', 'all') 只计算两个变量（两
列）所有观察值（所有行）都存在的两个变量。缺失情况下，就是用这种方法计算。
例：
x =[ 3 9 2 34
5 NaN
6 66
NaN 4 9 7
8 7 33 9
9 14 11 22]
[r,p,rlo,rup]=corrcoef(x,'rows', 'all')
r =
NaN NaN NaN NaN
NaN NaN NaN NaN
NaN NaN 1.0000 -0.5296
NaN NaN -0.5296 1.0000
p =
NaN NaN NaN NaN
NaN NaN NaN NaN
NaN NaN 1.0000 0.3587
NaN NaN 0.3587 1.0000
rlo =
NaN NaN NaN NaN
NaN NaN NaN NaN
NaN NaN 1.0000 -0.9623
NaN NaN -0.9623 1.0000
rup =
NaN NaN NaN NaN
NaN NaN NaN NaN
NaN NaN 1.0000 0.6619
NaN NaN 0.6619 1.0000
说明：NaN表示不存在数据，在x中，只有第四列和每五列每行的数据都存在，[r,p,rlo,rup]=corrcoef(x,'rows', 'all')
只计算第三列和第四列的相关系数及相应的区间。即：
[r,p,rlo,rup]=corrcoef(x(:,3), x(:,4))
r =
1.0000 -0.5296
-0.5296 1.0000
p =
1.0000 0.3587
0.3587 1.0000
rlo =
1.0000 -0.9623
-0.9623 1.0000
rup =
1.0000 0.6619
0.6619 1.0000
命令：[r,p,rlo,rup]=corrcoef(x,'rows', 'complete') 如果某行含有N aN，则去掉所有含N aN的行，再计算相关系数和区间。
例：x =[ 3 9 2 34
5 NaN
6 66
NaN 4 9 7
8 7 33 9
9 14 11 22]
[r,p,rlo,rup]=corrcoef(x,'rows', 'complete')
r =
1.0000 0.3883 0.6080 -0.7630
0.3883 1.0000 -0.4956 0.2995
0.6080 -0.4956 1.0000 -0.9771
-0.7630 0.2995 -0.9771 1.0000
p =
1.0000 0.7462 0.5840 0.4475
0.7462 1.0000 0.6699 0.8064
0.5840 0.6699 1.0000 0.1365
0.4475 0.8064 0.1365 1.0000
rlo =
1 NaN NaN NaN
NaN 1 NaN NaN
NaN NaN 1 NaN
NaN NaN NaN 1
rup =
1 NaN NaN NaN
NaN 1 NaN NaN
NaN NaN 1 NaN
NaN NaN NaN 1
[r,p,rlo,rup]=corrcoef(x,'rows', 'complete')与去掉二、三行所得矩阵再求 corrcoef相同。即与：
[r,p,rlo,rup]=corrcoef([x(1,:);x(4,:);x(5,:)])
r =
1.0000 0.3883 0.6080 -0.7630
0.3883 1.0000 -0.4956 0.2995
0.6080 -0.4956 1.0000 -0.9771
-0.7630 0.2995 -0.9771 1.0000
p =
1.0000 0.7462 0.5840 0.4475
0.7462 1.0000 0.6699 0.8064
0.5840 0.6699 1.0000 0.1365
0.4475 0.8064 0.1365 1.0000
rlo =
1 NaN NaN NaN
NaN 1 NaN NaN
NaN NaN 1 NaN
NaN NaN NaN 1
rup =
1 NaN NaN NaN
NaN 1 NaN NaN
NaN NaN 1 NaN
NaN NaN NaN 1
这里rlo,rup中显示NaN，是因为求区间要用到n-3,这里n-3=0
命令：[r,p,rlo,rup]=corrcoef(x,'rows', 'pairwise') r(i,j)的计算是把X中i或j列中含有NaN的行去掉，再计算i与j列的相关系数等指标。
例：x =[ 3 9 2 34
5 NaN
6 66
NaN 4 9 7
8 7 33 9
9 14 11 22]
[r,p,rlo,rup]=corrcoef(x,'rows', 'pairwise') r =
1.0000 0.3883 0.6476 -0.5791
0.3883 1.0000 -0.1687 0.5868
0.6476 -0.1687 1.0000 -0.5296
-0.5791 0.5868 -0.5296 1.0000
p =
1.0000 0.7462 0.3524 0.4209
0.7462 1.0000 0.8313 0.4132
0.3524 0.8313 1.0000 0.3587
0.4209 0.4132 0.3587 1.0000
rlo =
1.0000 NaN -0.8302 -0.9895
NaN 1.0000 -0.9722 -0.8584
-0.8302 -0.9722 1.0000 -0.9623
-0.9895 -0.8584 -0.9623 1.0000
rup =
1.0000 NaN 0.9915 0.8614
NaN 1.0000 0.9457 0.9897
0.9915 0.9457 1.0000 0.6619
0.8614 0.9897 0.6619 1.0000
r矩阵中0.3883只按以下求得：
[r,p,rlo,rup]=corrcoef([3;8;9;],[9;7;14])
r =
1.0000 0.3883
0.3883 1.0000
p =
1.0000 0.7462
0.7462 1.0000
rlo =
1 NaN
NaN 1
rup =
1 NaN
NaN 1
r矩阵中0.6476只按以下求得：
[r]=corrcoef([3;5;8;9;],[2;6;33;11])
r =
1.0000 0.6476
0.6476 1.0000
相关系数（二）
corrceof只能计算pearson相关系数。corr则还可计算其
它相关系数。不加说明时都是指pearson相关系数。 corr 不能做区间估计。
命令：rho=corr(x) x为n×p矩阵，返回p×p的相关系数矩阵。它与corcoef(x)计算结果相同。
例：x=[1 2 7 4 ;5 12 7 8;9 17 11 17];
rho=corr(x)
rho =
1.0000 0.9820 0.8660 0.9762
0.9820 1.0000 0.7559 0.9177
0.8660 0.7559 1.0000 0.9538
0.9762 0.9177 0.9538 1.0000
与corcoef(x)相同
corrcoef(x)
ans =
1.0000 0.9820 0.8660 0.9762
0.9820 1.0000 0.7559 0.9177
0.8660 0.7559 1.0000 0.9538
0.9762 0.9177 0.9538 1.0000
命令：rho=corr(x,y) x为n×p1矩阵, y为n×p2矩阵,rho 为p1×p2相关矩阵。它不是对称矩阵，rho(i,j)为x的第i 列（变量）与y的第j列（变量）的相关系数。注意：corrcoef(x,y)只适合大小相同的列或行向量。
例：x=[1 2 7 4 ;5 12 7 8;9 17 11 17];y=[4 9;23 12;-2 -7];
rho=corr(x,y)
rho =
-0.2299 -0.7832
-0.0418 -0.6516
-0.6857 -0.9892
-0.4354 -0.8994
如-0.6857为x第3列与y的第1列的相关系数。它等于：rho=corr(x(:,3),y(:,1))
rho =
-0.6857
命令：[rho,pval]=corr(…)
Pval是返回原假设为两变量相关系数为0时，返回P值。例：x=[1 2 7 4 ;5 12 7 8;9 17 11 17];y=[4 9;23 12;-2 -7];
[rho,pval]=corr(x,y)
rho =
-0.2299 -0.7832
-0.0418 -0.6516
-0.6857 -0.9892
-0.4354 -0.8994
pval =
0.8523 0.4272
0.9734 0.5482
0.5190 0.0938
0.7132 0.2881
我们看pval值计算，如0.0938的计算：即x的第三列与y 的第二列相关系数所对应的值的计算。
[r,p]=corrcoef([7;7;11], [9 12 -7])
r =
1.0000 -0.9892
-0.9892 1.0000
p =
1.0000 0.0938
0.0938 1.0000
0.0938小于显著性水平0.1，即拒绝x的第三列与y第二列不相关的原假设。
命令：[...] = corr(...,'param1', val1, 'param2', val2,...)
例：《统计学原理》黄良文 p226页
通信原理课程设计报告书 课题名称 基于MATLAB 的（7，4）汉明码编 译码设计与仿真结果分析 姓 名 学 号 学 院 通信与电子工程学院 专 业 通信工程 指导教师 ※※※※※※※※※ ※ ※ ※※ ※ ※ 2009级通信工程专业 通信原理课程设计
2011年 12月 23日 一、设计任务及要求： 设计任务： 利用MATLAB编程，实现汉明码编译码设计。理解（7,4）汉明码的构造原理，掌握（7，4）汉明码的编码和译码的原理和设计步骤。并对其性能进行分析。要求： 通过MATLAB编程，设计出（7,4）汉明码的编码程序，编码后加入噪声，然后译码，画出信噪比与误比特数和信噪比与误比特率的仿真图，然后对其结果进行分析 指导教师签名： 2011年12月23日 二、指导教师评语： 指导教师签名： 年月日 三、成绩 验收盖章 年月日
基于MATLAB 的（7，4）汉明码编译码设计 与仿真结果分析 1 设计目的 （1）熟悉掌握汉明码的重要公式和基本概念。 （2）利用MATLAB 编程，实现汉明码编译码设计。 （3）理解（7,4）汉明码的构造原理，掌握（7，4）汉明码的编码和译码的原理和设计步骤。 （4）对其仿真结果进行分析。 2 设计要求 （1）通过MATLAB 编程，设计出（7,4）汉明码的编码程序。 （2）编码后加入噪声，然后译码，画出信噪比与误比特数和信噪比与误比特率的仿真图。 （3）然后对其结果进行分析。 3 设计步骤 3.1 线性分组码的一般原理 线性分组码的构造 3.1.1 H 矩阵 根据(7, 4)汉明码可知一般有 现在将上面它改写为 上式中已经将“⊕”简写成“+”。 上式可以表示成如下矩阵形式： ??? ??=⊕⊕⊕=⊕⊕⊕=⊕⊕⊕0 000346 13562456a a a a a a a a a a a a ?? ? ?? =?+?+?+?+?+?+?=?+?+?+?+?+?+?=?+?+?+?+?+?+?010011010010101100010111012345601234560123456a a a a a a a a a a a a a a a a a a a a a (1) (2)
clc; N=7; %以7级寄存器为例，并组其中的一组优选对：211，,217 connections=gfprimfd(N,'all'); f1=connections(4,:); %取一组本原多项式序列，211 f2=connections(16,:); %取另一组本原多项式序列,217 registers1=[0 0 0 0 0 0 1];%给定寄存器的初始状态 registers2=[0 0 0 0 0 0 1];%取相同的初始状态 L=2^N-1; %周期长度 sum2=0; sum1=0; for k=1:L seq1(k)=registers1(N); %第一组m序列 seq2(k)=registers2(N); %第二组序列 for j=1:N %进行模2加 sum1=sum1+f1(j+1)*registers1(j); %各级寄存器送参与模2加的值sum1=mod(sum1,2); sum2=sum2+f2(j+1)*registers2(j); %各级寄存器送参与模2加的值sum2=mod(sum2,2); end for t=N:-1:2 %寄存器移位 registers1(t)=registers1(t-1); registers2(t)=registers2(t-1); end registers1(1)=sum1; registers2(1)=sum2; sum2=0; sum1=0; end disp(f1); disp(f2); z=seq1+seq2; %m序列的相加 gold=mod(z,2); %模2运算 gold=1-2*gold; %转换为2值电平 disp(gold); R=xcorr(gold,'unbiased'); %自相关 R=R/max(R); %归一化 figure;plot(R);title('gold序列的自相关函数'); s=fftshift(abs(fft(gold,2*L)).^2); %求功率谱 s=s/max(s); figure;plot(s);title('gold序列的功率谱');
MATLAB 智能算法30个案例分析（终极版） 1 基于遗传算法的TSP算法（王辉） 2 基于遗传算法和非线性规划的函数寻优算法（史峰） 3 基于遗传算法的BP神经网络优化算法（王辉） 4 设菲尔德大学的MATLAB遗传算法工具箱（王辉） 5 基于遗传算法的LQR控制优化算法（胡斐） 6 遗传算法工具箱详解及应用（胡斐） 7 多种群遗传算法的函数优化算法（王辉） 8 基于量子遗传算法的函数寻优算法（王辉） 9 多目标Pareto最优解搜索算法（胡斐） 10 基于多目标Pareto的二维背包搜索算法（史峰） 11 基于免疫算法的柔性车间调度算法（史峰） 12 基于免疫算法的运输中心规划算法（史峰） 13 基于粒子群算法的函数寻优算法（史峰） 14 基于粒子群算法的PID控制优化算法（史峰） 15 基于混合粒子群算法的TSP寻优算法（史峰） 16 基于动态粒子群算法的动态环境寻优算法（史峰） 17 粒子群算法工具箱（史峰） 18 基于鱼群算法的函数寻优算法（王辉） 19 基于模拟退火算法的TSP算法（王辉） 20 基于遗传模拟退火算法的聚类算法（王辉） 21 基于模拟退火算法的HEV能量管理策略参数优化（胡斐）
22 蚁群算法的优化计算——旅行商问题（TSP）优化（郁磊） 23 基于蚁群算法的二维路径规划算法（史峰） 24 基于蚁群算法的三维路径规划算法（史峰） 25 有导师学习神经网络的回归拟合——基于近红外光谱的汽油辛烷值预测（郁磊） 26 有导师学习神经网络的分类——鸢尾花种类识别（郁磊） 27 无导师学习神经网络的分类——矿井突水水源判别（郁磊） 28 支持向量机的分类——基于乳腺组织电阻抗特性的乳腺癌诊断（郁磊） 29 支持向量机的回归拟合——混凝土抗压强度预测（郁磊） 30 极限学习机的回归拟合及分类——对比实验研究（郁磊） 智能算法是我们在学习中经常遇到的算法，主要包括遗传算法，免疫算法，粒子群算法，神经网络等，智能算法对于很多人来说，既爱又恨，爱是因为熟练的掌握几种智能算法，能够很方便的解决我们的论坛问题，恨是因为智能算法感觉比较“玄乎”，很难理解，更难用它来解决问题。 因此，我们组织了王辉，史峰，郁磊，胡斐四名高手共同写作MATLAB智能算法，该书包含了遗传算法，免疫算法，粒子群算法，鱼群算法，多目标pareto算法，模拟退火算法，蚁群算法，神经网络，SVM等，本书最大的特点在于以案例为导向，每个案例针对一
MATLAB程序设计实验报告
————————————————————————————————作者：————————————————————————————————日期：
MATLAB 程序设计实验报告 一、实验目的 1. 通过实验熟悉MA TLAB 仿真软件的使用方法； 2. 掌握用MATLAB 对连续信号时域分析、频域分析和s 域分析的方法，利用绘图命令绘制出典型信号的波形，了解这些信号的基本特征； 3. 掌握用MATLAB 对离散信号时域分析、频域分析和z 域分析的方法，利用绘图命令绘制出典型信号的波形，了解这些信号的基本特征； 4. 通过绘制信号运算结果的波形，了解这些信号运算对信号所起的作用。 二、实验设备 1. 计算机 2. MA TLAB R2007a 仿真软件 三、实验原理 1.MATLAB 对系统的时域分析 信号的时域运算包括信号的相加、相乘，信号的时域变换包括信号的平移、反折、倒相及信号的尺度变换。 （1）信号的相加和相乘：已知信号)(1t f 和)(2t f ，信号相加和相乘记为 )()(1t f t f =)(2t f +；)()(1 t f t f =)(2t f *。 （2）信号的微分和积分：对于连续时间信号，其微分运算是用diff 函数来完成的，其语句格式为：diff(function,’variable ’,n)，其中function 表示需要进行求导运算的信号，或者被赋值的符号表达式；variable 为求导运算的独立变量；n 为求导的阶数，默认值为求一阶导数。连续信号的积分运算用int 函数来完成，语句格式为：diff(function,’variable ’,a,b)，其中function 表示需要进行被积信号，或者被赋值的符号表达式；variable 为求导运算的独立变量；a,b 为积分上、下限，a 和b 省略时为求不定积分。 （3）信号的平移、翻转和尺度变换 信号的平移包含信号的左移与右移，信号的翻转包含信号的倒相与折叠，平移和翻转信号不会改变信号)(t f 的面积和能量。信号的尺度变换是对信号)(t f 在时间轴上的变化，可使信号压缩或扩展。)(at f 将原波形压缩a 倍，)/(a t f 将原波形扩大a 倍。 2.MATLAB 对系统频率特性的分析 （1）系统的频率响应 设线性时不变（LTI ）系统的冲激响应为)(t h ，该系统的输入（激励）信号为)(t f ,则
一、直接法 clear;clc;close all; %清除变量；清屏；关闭当前图形窗口 Fs=1000; t=0:1/Fs:1; nfft=2048; %改变nfft的值可对比不同采样值时的谱估计效果 %****************生成信号、噪声**************% x1=cos(2*pi*40*t)+3*cos(2*pi*45*t);%信号 x2=randn(size(t)); %噪声 x3=x1+x2; %信号+噪声 [Pxx,f]=periodogram(x3,window,nfft,Fs); %直接法 plot(f,10*log10(Pxx)); title('直接法 nfft=2048'); set(gca,'xlim',[1 120]); ylabel('Am/dB'); xlabel('Frequency/Hz'); 二、间接法 Fs=1000;% 采样频率 n=0:1/Fs:1;% 产生含有噪声的序列 x1=cos(2*pi*40*n)+3*cos(2*pi*45*n);%信号x2=randn(size(n)); %噪声x3=x1+x2; %信号+噪声 nfft=1024; cxn=xcorr(x3);% 计算序列的自相关函数 CXk=fft(cxn); Pxx=abs(CXk); index=0:round(nfft/2-1); f=index*Fs/nfft; plot_Pxx=10*log10(Pxx(index+1)); figure (1) plot(f,plot_Pxx); title('间接法 nfft=1024');ylabel('Am/dB'); set(gca,'xlim',[1 120]); xlabel('Frequency/Hz'); 三、Bartlett法 clear;clc;close all; %清除变量；清屏；关闭当前图形窗口 Fs=1000; t=0:1/Fs:1; nfft=1024; %****************生成信号、噪声**************% x1=cos(2*pi*40*t)+3*cos(2*pi*45*t);%信号 x2=randn(size(t)); %噪声 x3=x1+x2; %信号+噪声 window=hamming(512); %海明窗 noverlap=0; %数据无重叠 p=0.9; %置信概率 [Pxx,Pxxc]=psd(x3,nfft,Fs,window,noverlap,p); index=0:round(nfft/2-1); k=index*Fs/nfft; plot_Pxx=10*log10(Pxx(index+1)); plot_Pxxc=10*log10(Pxxc(index+1)); figure(1) plot(k,plot_Pxx);title('Bartlett法海明窗');; set(gca,'xlim',[1120]);ylabel('Am/dB'); xlabel('Frequency/Hz'); 四、Welch法
4.13 表中是某软件公司月销售额数据，其中，x为总公司的月销售额（万元）；y为某分公司的月销售额（万元）。 （1）用普通最小二乘法建立x和y的回归方程。 （2）用残差图及DW检验诊断序列的自相关性。 （3）用迭代法处理序列相关，并建立回归方程。 （4）用一阶差分法处理数据，并建立回归方程。 （5）比较以上各方法所建回归方程的优良性。 序号x y 序号x y 1 127.3 20.96 11 148.3 24.54 2 130.0 21.40 12 146.4 24.28 3 132.7 21.96 13 150.2 25.00 4 129.4 21.52 14 153.1 25.64 5 135.0 22.39 15 157.3 26.46 6 137.1 22.76 16 160. 7 26.98 7 141.1 23.48 17 164.2 27.52 8 142.8 23.66 18 165.6 27.78 9 145.5 24.10 19 168.7 28.24 10 145.3 24.01 20 172.0 28.78 (1)aa_size=size(aa,1) >> x=[ones(aa_size,1),aa(:,1)]; >> y=aa(:,2); >> b_est=inv(x'*x)*x'*y; b_est b_est = -1.4348 0.1762 (2) y_est=x*b_est; >> b1=y-y_est; >> plot(b1,'ro') p01=sum(b1(1:(aa_size-1)).*b1(2:(aa_size))); >> p02=sqrt(sum(b1(1:(aa_size-1)).^2)*sum(b1(2:aa_size).^2)); >> p=p01/p02 DW=2*(1-p) DW = 0.6793
§9. 利用Matlab 和SPSS 实现主成分分析 1.直接调用Matlab 软件实现 在软件Matlab 中实现主成分分析可以采取两种方式实现：一是通过编程来实现；二是直接调用Matlab 中自带程序实现。 通过直接调用Matlab 中的程序可以实现主成分分析： )(]2,var ,,[X princomp t iance score pc = 式中：X 为输入数据矩阵 ? ? ??? ???????=nm n n m m x x x x x x x x x X 2 1 22221 11211 （一般要求n>m ） 输出变量： ①pc 主分量f i 的系数，也叫因子系数；注意：pc T pc=单位阵 ②score 是主分量下的得分值；得分矩阵与数据矩阵X 的阶数是一致的； ③variance 是score 对应列的方差向量，即A 的特征值；容易计算方差所占的百分比 percent-v = 100*variance/sum(variance); ④t2表示检验的t2-统计量（方差分析要用） 计算过程中应用到计算模型:
ξ+????????????=??????????????m T p x x x A f f f 2121 (要求p MATLAB 智能算法30个案例分析 智能算法是我们在学习中经常遇到的算法，主要包括遗传算法，免疫算法，粒子群算法，神经网络等，智能算法对于很多人来说，既爱又恨，爱是因为熟练的掌握几种智能算法，能够很方便的解决我们的论坛问题，恨是因为智能算法感觉比较“玄乎”，很难理解，更难用它来解决问题。 因此，我们组织了王辉，史峰，郁磊，胡斐四名高手共同写作MATLAB智能算法，该书包含了遗传算法，免疫算法，粒子群算法，鱼群算法，多目标pareto算法，模拟退火算法，蚁群算法，神经网络，SVM等，本书最大的特点在于以案例为导向，每个案例针对一个实际问题，给出全部程序和求解思路，并配套相关讲解视频，使读者在读过一个案例之后能够快速掌握这种方法，并且会套用案例程序来编写自己的程序。本书作者在线，读者和会员可以向作者提问，作者做到有问必答。 本书和目录如下： 1 基于遗传算法的TSP算法（王辉） TSP (旅行商问题—Traveling Salesman Problem)，是典型的NP完全问题，即其最坏情况下的时间复杂性随着问题规模的增大按指数方式增长，到目前为止不能找到一个多项式时间的有效算法。遗传算法是一种进化算法，其基本原理是仿效生物界中的“物竞天择、适者生存”的演化法则。遗传算法的做法是把问题参数编码为染色体，再利用迭代的方式进行选择、交叉以及变异等运算来交换种群中染色体的信息，最终生成符合优化目标的染色体。实践证明，遗传算法对于解决TSP问题等组合优化问题具有较好的寻优性能。 2 基于遗传算法和非线性规划的函数寻优算法（史峰） 遗传算法提供了求解非线性规划的通用框架，它不依赖于问题的具体领域。遗传算法的优点是将问题参数编码成染色体后进行优化，而不针对参数本身，从而不受函数约束条件的限制；搜索过程从问题解的一个集合开始，而不是单个个体，具有隐含并行搜索特性，可大大减少陷入局部最小的可能性。而且优化计算时算法不依赖于梯度信息，且不要求目标函数连续及可导，使其适于求解传统搜索方法难以解决的大规模、非线性组合优化问题。 3 基于遗传算法的BP神经网络优化算法（王辉） BP模型被广泛地应用于模式分类、模式识别等方面．但BP算法收敛速度慢，且很容易陷入局部极小点，而遗传算法具有并行搜索、效率高、不存在局部收敛问题等优点而被广泛应用．遗传算法的寻优过程带有一定程度的随机性和盲从性，多数情况下只能收敛到全局次优解，且有过早收敛的现象．为了克服遗传算法寻优过程的盲从性，将有监督学习的BP算法与之结合以达到优势互补、提高算法的稳定性和全局搜索能力的目的。 4 设菲尔德大学的MATLAB遗传算法工具箱（王辉） Matlab 遗传算法(Genetic Algorithm)优化工具箱是基于基本操作及终止条件、二进制和十进制相互转换等操作的综合函数库。其实现步骤包括：通过输入及输出函数求出遗传算法主函数、初始种群的生成函数，采用选择、交叉、变异操作求得基本遗传操作函数。以函数仿真 %matlab程序：生成m序列和相关性分析 clc; N=5; connections=gfprimfd(N,'all');%生成级数为5时所有本原多项式系数序列矩阵f1=connections(2,:); %取一组本原多项式序列，此系数为45（100101）f2=connections(3,:); %取另一组本原多项式序列，此系数为75（111101）registers1=[0 0 0 0 1];%给定寄存器的初始状态 registers2=[0 0 0 0 1];%取相同的初始状态 L=2^N-1; %周期长度 sum2=0; sum1=0; for k=1:L seq1(k)=registers1(N); %第一组m序列 seq2(k)=registers2(N); %第二组序列 for j=1:N %进行模2加 sum1=sum1+f1(j+1)*registers1(j); %各级寄存器送参与模2加的值 sum1=mod(sum1,2); sum2=sum2+f2(j+1)*registers2(j); %各级寄存器送参与模2加的值 sum2=mod(sum2,2); end for t=N:-1:2 %寄存器移位 registers1(t)=registers1(t-1); registers2(t)=registers2(t-1); end registers1(1)=sum1; registers2(1)=sum2; sum1=0; sum2=0; end disp(f1); %显示反馈系数序列 disp(seq1); %显示第一组m序列 m1=xcorr(seq1，’unbiased’);%计算自相关函数 figure;stem(m1/max(m1)); title('m1的自相关函数'); %画出m1序列的自相关函数 disp(f2); disp(seq2); m2=xcorr(seq2,’unbiased’); figure;stem(m2/max(m2)); title('m2的自相关函数');%画出m2序列的自相关函数 R12= xcorr(seq1,seq2,’unbiased’) figure;plot(R12/max(R12)); title('m2,m1的互相关函数');%画出m1，m2序列的互相关特性 s=fftshift(abs(fft(seq1,2*L)).^2); s=s/max(s); figure;plot(s); §9 利用Matlab 和SPSS 实现主成分分析 1.直接调用Matlab 软件实现 在软件Matlab 中实现主成分分析可以采取两种方式实现：一是通过编程来实现；二是直接调用Matlab 中自带程序实现。 通过直接调用Matlab 中的程序可以实现主成分分析： )princomp(X t2]variance,score,[pc,= 式中：X 为输入数据矩阵 ? ? ??? ???????=nm n2n12m 2221 1m 1211 x x x x x x x x x X L M M M L L （一般要求n>m ） 输出变量： ①pc 主分量f i 的系数，也叫因子系数；注意：pc T pc=单位阵。 ②score 是主分量下的得分值；得分矩阵与数据矩阵X 的阶数是一致的。 ③variance 是score 对应列的方差向量，即A 的特征根；容易计算方差所占的百分比。 percent-v = 100*variance/sum(variance); ④t2表示检验的t-统计量（方差分析要用） 计算过程中应用到计算模型: ξx x x A f f f m 21T p 21+?????? ??????=?????? ????????M M (要求p 9.3.1 基本操作及命令 1. 访问和退出MATLAB 在大多数系统中，一旦安装了MATLAB，在调用时，应执行命令MATLAB。退出MATLAB应执行命令exit或quit。 2. 如何应用MATLAB 通常以命令驱动方式应用MATLAB。当输入单个命令时，MATLAB会立即对其进行处理，并且显示处理结果。MATLAB 也能够执行存储在文件中的命令序列。 通过键盘输入的命令，应用向上箭头键可以被存取。通过输入某个最新命令和调用特定的命令行，可以使屏幕内容向上滚动。 3. MATLAB的变量 MATLAB的一个特点是变量在应用之前不必是维数确定的。在MATLAB中，变量一旦被采用，会自动产生（如果必要，变量的维数以后还可以改变）。在命令exit或quit输入之前，这些变量将保留在存储器中。 为了得到工作空间内的变量清单，可以通过键盘输入命令who，当前存放在工作空间内的所有变量便会显示在屏幕上。 命令clear能从工作空间中清除所有非永久性变量。如果只需要从工作空间中清除某个特定变量，比如“x”，则应输入命令clear x。 4. 以“%”开始的程序行 在本书中，许多MATLAB程序在编写时附有注解和说明，这些注解和说明阐明了发生在程序中的具体进程。在MATLAB 中以“%”开始的程序行，表示注解和说明。符号“%”类似于BASIC中“REM”。以“%”开始的行，用来存储程序的注解或说明，这些注解和说明是不执行的。这就是说，在MATLAB程序行中，出现在“%”以后的一切内容都是可以忽略的。如果注解或说明需要一行以上程序行，则每一行均需以“%”为起始。 5. 应用分号操作符 分号用来取消打印。如果语句的最后一个符号是分号，则打印被取消，但是命令仍在执行，而结果不再显示。这是一个有益的特性，因为打印中间结果可能不必要。此外，在输入矩阵时，除非最后一行，分号用来指示一行的结束。 6. 应用冒号操作符 第十章典型相关分析 安庆师范学院胡云峰 习题10.2下表给出著名统计学家Rao在1952年对25个家庭的成年长子的头长（x1）、头宽（x2）、与次子头长（y1）、头宽（y2）进行调查所得数据如下: 长子次子长子次子头长头宽头长头宽头长头宽头长头宽191155179145190159195157 195149201152188151187158 181148185149163137161130 183153188149195155183158 176144171142186153173148 208157192152181145182146 189150190149175140165137 197159189152192154185152 188152197159174143178147 192150187151176139176143 179158186148197167200158 183147174147190163187150 174150185152 x=[191155179145;195149201152;181148185149;... 183153188149;176144171142;208157192152;... 189150190149;197159189152;188152197159;... 192150187151;179158186148;183147174147;... 174150185152;190159195157;188151187158;... 163137161130;195155183158;186153173148;... 181145182146;175140165137;192154185152;... 174143178147;176139176143;197167200158;... 190163187150] 第一步计算相关系数矩阵 程序R=corrcoef(x) 输出结果R= 1.00000.73460.71080.7040 0.7346 1.00000.69320.7086 0.71080.6932 1.00000.8393 0.70400.70860.8393 1.0000 计算A、B的特征值特征向量 程序 R11=R([1,2],[1,2]); §9.利用Matlab 和SPSS 实现主成分分析 1.直接调用Matlab 软件实现 在软件Matlab 中实现主成分分析可以采取两种方式实现：一是通过编程来 实现；二是直接调用Matlab 中自带程序实现。 通过直接调用Matlab 中的程序可以实现主成分分析: [pc,score, variance,t2] = princomp (X) 式中：X 为输入数据矩阵 _x i x 1m x 2m Xnm -（一般要求 n>m ） ① pc 主分量f i 的系数，也叫因子系数；注意：pc T pc=单位阵 ② score 是主分量下的得分值；得分矩阵与数据矩阵 X 的阶数是一致的； ③ v arianee 是score 对应列的方差向量，即 A 的特征值；容易计算方差所占 的百分比 1 X 12 X 21 X 22 _Xni X n2 输出变量: perce nt-v = 100*varia nce/sum(varia nee); ④t2表示检验的t2-统计量（方差分析要用） 计算过程中应用到计算模型 例：表1为某地区农业生态经济系统各区域单元相关指标数据， 运用主成分 分析方法可以用更少的指标信息较为精确地描述该地区农业生态经济的发展状 况。 表1某农业生态经济系统各区域单元的有关数据 样本 序号 X 1：人口密度 2 (人/km ) X 2：人均耕 地面积(ha) X 3：森林覆 盖率(%) X 4：农民人 均纯收入(元 /人) X 5：人均粮 食产量(kg/ 人) X 6：经济作 物占农作物 播面比例 (%) X 7：耕地占 土地面积比 率(% ) X 8：果园与 林地面积之 比(%) X 9：灌溉田 占耕地面积 之比(%) 1 363.912 0.352 16.101 192.11 295.34 26.724 18.492 2.231 26.262 2 141.503 1.684 24.301 1 752.35 452.26 32.314 14.464 1.455 27.066 3 100.695 1.067 65.601 1 181.54 270.12 18.266 0.162 7.474 12.489 4 143.739 1.336 33.205 1 436.12 354.26 17.486 11.805 1.892 17.534 5 131.412 1.623 16.607 1 405.09 586.59 40.683 14.401 0.303 22.932 6 68.337 2.032 76.204 1 540.29 216.39 8.128 4.065 0.011 4.861 7 95.416 0.801 71.106 926.35 291.52 8.135 4.063 0.012 4.862 8 62.901 1.652 73.307 1 501.24 225.25 18.352 2.645 0.034 3.201 9 86.624 0.841 68.904 897.36 196.37 16.861 5.176 0.055 6.167 10 91.394 0.812 66.502 911.24 226.51 18.279 5.643 0.076 4.477 11 76.912 0.858 50.302 103.52 217.09 19.793 4.881 0.001 6.165 12 51.274 1.041 64.609 968.33 181.38 4.005 4.066 0.015 5.402 13 68.831 0.836 62.804 957.14 194.04 9.110 4.484 0.002 5.790 14 77.301 0.623 60.102 824.37 188.09 19.409 5.721 5.055 8.413 15 76.948 1.022 68.001 1 255.42 211.55 11.102 3.133 0.010 3.425 16 99.265 0.654 60.702 1 251.03 220.91 4.383 4.615 0.011 5.593 17 118.505 0.661 63.304 1 246.47 242.16 10.706 6.053 0.154 8.701 18 141.473 0.737 54.206 814.21 193.46 11.419 6.442 0.012 12.945 19 137.761 0.598 55.901 1 124.05 228.44 9.521 7.881 0.069 12.654 20 117.612 1.245 54.503 805.67 175.23 18.106 5.789 0.048 8.461 21 122.781 0.731 49.102 1 313.11 236.29 26.724 7.162 0.092 10.078 对于上述例子，Matlab 进行主成分分析，可以得到如下结果 ①以及每一个主成分的贡献率和累计贡献率，如表 2和图1 表2.特征根及主成分贡献率 ■f j [ x j f2 a T X2 I- =A Jp 一 1 x m （要求 p §. 利用和实现主成分分析 .直接调用软件实现 在软件中实现主成分分析可以采取两种方式实现：一是通过编程来实现；二是直接调用中自带程序实现。 通过直接调用中的程序可以实现主成分分析： )(]2,var ,,[X princomp t iance score pc = 式中：为输入数据矩阵 ? ? ??? ???????=nm n n m m x x x x x x x x x X 2 1 22221 11211 （一般要求>） 输出变量： ① 主分量的系数，也叫因子系数；注意：单位阵 ②是主分量下的得分值；得分矩阵与数据矩阵的阶数是一致的； ③是对应列的方差向量，即的特征值；容易计算方差所占的百分比 *(); ④表示检验的统计量（方差分析要用） 计算过程中应用到计算模型: ξ+????????????=??????????????m T p x x x A f f f 2121 (要求<) 例：表为某地区农业生态经济系统各区域单元相关指标数据，运用主成分分析方法可以用更少的指标信息较为精确地描述该地区农业生态经济的发展状况。 表 某农业生态经济系统各区域单元的有关数据 样本序号 ：人口密度(人) ：人均耕地面积() ：森林覆盖率() ：农民人均纯收入(元人) ：人均粮食 产量 (人) ：经济作物占农作物播面比例(％) ：耕地占土地面积比率(％) ：果园与林地面积之比(％) ：灌溉田占 耕地面积之 比(％) 对于上述例子，进行主成分分析，可以得到如下结果。 ① 以及每一个主成分的贡献率和累计贡献率，如表和图。 表. 特征根及主成分贡献率 主成分 特征值 贡献率％ 累积贡献率％ 第八章相关分析 协方差 命令：C = cov(X) 当X为行或列向量时，它等于 var(X) 样本标准差。 X=1:15;cov(X) ans = 20 >> var(X) ans = 20 当X为矩阵时，此时X的每行为一次观察值，每列为一个变量。cov（X）为协方差矩阵,它是对称矩阵。 例：x=rand(100,3);c= cov(x) c= 0.089672 -0.012641 -0.0055434 -0.012641 0.07928 0.012326 -0.0055434 0.012326 0.082203 c的对角线为：diag(c) ans = 0.0897 0.0793 0.0822 它等于：var(x) ans = 0.0897 0.0793 0.0822 sqrt(diag(cov(x))) ans = 0.2995 0.2816 0.2867 它等于：std(x) ans = 0.2995 0.2816 0.2867 命令：c = cov(x,y) 其中x和y是等长度的列向量（不是行向量），它等于cov([x y])或cov([x,y]) 例：x=[1;4;9];y=[5;8;6]; >> c=cov(x,y) c = 16.3333 1.1667 1.1667 2.3333 >> cov([x,y]) ans = 16.3333 1.1667 1.1667 2.3333 COV(X)、 COV(X,0)[两者相等] 或COV(X,Y)、COV(X,Y,0) [两者相等]，它们都是除以n-1，而COV(X,1) or COV(X,Y,1)是除以n x=[1;4;9];y=[5;8;6]; >> cov(x,y,1) ans = 10.8889 0.7778 0.7778 1.5556 它的对角线与var([x y],1) 相等 ans = 10.8889 1.5556 协差阵的代数计算： [n,p] = size(X); X = X - ones(n,1) * mean(X); Y = X'*X/(n-1); Y为X的协差阵MATLAB智能算法30个案例分析
m序列的生成及其相关性的matlab分析
利用Matlab和SPSS实现主成分分析
MATLAB经典控制系统的分析和设计
matlab与典型相关分析
9利用Matlab和SPSS实现主成分分析
§9.利用Matlab和SPSS实现主成分分析
《非参数统计》与MATLAB编程 第八章 相关分析
相关主题
文本预览