matlab编写程序.doc

mathematicsBasic Matrix Operations a=1 2 3 4 5生成矩阵； b=a+2矩阵加上数字 plot(b)画三点图 grid on生成网格 bar(b)生成条状图 xlabel(sample#) 给X轴加标注 ylabel(pound) 给Y轴加标注 title(bar plot)加标题 plot(b,*)用*表示点 axis(0 10 10 20 0 20)各个轴的范围 A = 1 2 0; 2 5 -1; 4 10 -1 B=A转置 C=A*B矩阵相乘 C=A.*B数组相乘 X=inv(A)逆 I=inv(A)*A单位矩阵 eig(A)特征值 svd(A) the singular value decomposition. 奇异值分解 p = round(poly(A)生成特征多项式的系数 roots(p) 特征多项式的根，即矩阵的特征值 q = conv(p,p) 向量的卷积 r = conv(p,q) 再向量的卷积 plot(r) who 变量列表 whos 变量的详细列表 sqrt(-1)； sqrt(2) 开方Matrix Manipulation A = magic(3) 生成魔方矩阵graphics2-D PlotsLine Plot of a Chirp（线状图） doc plot 在帮助中查看plot更多的信息 x=0:0.05:5; y=sin(x.2); plot(x,y); plot(x,y,o); xlabel(Time) ylabel(Amplitude)Bar Plot of a Bell Shaped Curve（条状图） x = -2.9:0.2:2.9; bar(x,exp(-x.*x);函数中包括了计算，plot直接完成Stairstep Plot of a Sine Wave（阶梯图） x=0:0.25:10; stairs(x,sin(x);Errorbar Plot x=-2:0.1:2; y=erf(x)； y=erf(x);e = rand(size(x)/10 errorbar(x,y,e);Polar Plot t=0:0.01:2*pi; polar(t,abs(sin(2*t).*cos(2*t); abs绝对值Stem Plot x = 0:0.1:4; y = sin(x.2).*exp(-x); stem(x,y)Scatter Plot load count.dat scatter(count(:,1),count(:,2),r*) xlabel(Number of Cars on Street A); ylabel(Number of Cars on Street B);更多plot的类型可以在命令行输入doc praph2d查看3-D PlotsMesh Plot of Peaks z=peaks(25);peaks is a function of two variables, obtained by translating and scaling Gaussian distributions, which is useful for demonstrating mesh, surf, pcolor, contour, and so on. mesh(z); meshc(z)；带有轮廓线 colormap(hsv)；对mesh的图象使用不同的色彩，如gray, hot, bone. EtcSurface Plot of Peaks z=peaks(25); surf(z); 3-D shaded surface plot surfc(z);带有轮廓线 colormap jet;同 colormap(jet); 对mesh的图象使用不同的色彩Surface Plot (with Shading) of Peaks z=peaks(25); surfl(z); Surface plot with colormap-based lighting shading interp; Set color shading properties colormap(pink);Contour Plot of Peaks z=peaks(25); contour(z,16) 轮廓线 colormap(hsv)Quiver （Quiver or velocity plot） x = -2:.2:2; y = -1:.2:1; xx,yy = meshgrid(x,y); zz = xx.*exp(-xx.2-yy.2); px,py = gradient(zz,.2,.2); quiver(x,y,px,py,2);Slice（Volumetric slice plot） x,y,z = meshgrid(-2:.2:2,-2:.25:2,-2:.16:2); v = x.*exp(-x.2-y.2-z.2); xslice = -1.2,.8,2; yslice = 2; zslice = -2,0; x,y,z = meshgrid(-2:.2:2,-2:.25:2,-2:.16:2); slice(x,y,z,v,xslice,yslice,zslice) colormap hsvProgrammingManipulating Multidimensional ArraysCreating Multi-Dimensional ArraysThe CAT function is a useful tool for building multidimensional arrays. B = cat(DIM,A1,A2,.) builds a multidimensional array by concatenating A1, A2 . along the dimension DIM. Calls to CAT can be nested. B = cat( 3, 2 8; 0 5, 1 3; 7 9, 2 3; 4 6);B(:,:,1) = 2 8 0 5B(:,:,2) = 1 3 7 9B(:,:,3) = 2 3 4 6 A = cat(3,9 2; 6 5, 7 1; 8 4); B = cat(3,3 5; 0 1, 5 6; 2 1); C = cat(4,A,B,cat(3,1 2; 3 4, 4 3; 2 1);C(:,:,1,1) = 9 2 6 5C(:,:,2,1) = 7 1 8 4C(:,:,1,2) = 3 5 0 1C(:,:,2,2) = 5 6 2 1C(:,:,1,3) = 1 2 3 4C(:,:,2,3) = 4 3 2 1 A = cat(2,9 2; 6 5, 7 1; 8 4); B = cat(2,3 5; 0 1, 5 6; 2 1); C = cat(4,A,B,cat(2,1 2; 3 4, 4 3; 2 1)C(:,:,1,1) = 9 2 7 1 6 5 8 4C(:,:,1,2) = 3 5 5 6 0 1 2 1C(:,:,1,3) = 1 2 4 3 3 4 2 1 C = cat(2,A,B,cat(2,1 2; 3 4, 4 3; 2 1)C =9 2 7 1 3 5 5 6 1 2 4 3;6 5 8 4 0 1 2 1 3 4 2 1 C = cat(1,A,B,cat(2,1 2; 3 4, 4 3; 2 1)C = 9 2 7 1 6 5 8 4 3 5 5 6 0 1 2 1 1 2 4 3 3 4 2 1Finding the DimensionsSzA = size(A)DimsA = ndims(A)SzC = size(C)DimsC = ndims(C)Accessing Elements K = C(:,:,1,1 3)K(:,:,1,1) = 9 2 6 5K(:,:,1,2) = 1 2 3 4Manipulating Multi-Dimensional Arrays A = rand(3,3,2); B = permute(A, 2 1 3); C = permute(A, 3 2 1);Selecting 2D Matrices From Multi-Dimensional Arrays A = cat( 3, 1 2 3; 9 8 7; 4 6 5, 0 3 2; 8 8 4; 5 3 5, . 6 4 7; 6 8 5; 5 4 3);% The EIG function is applied to each of the horizontal slices of A.for i = 1:3 eig(squeeze(A(i,:,:)end x1 = -2*pi:pi/10:0; x2 = 2*pi:pi/10:4*pi; x3 = 0:pi/10:2*pi; x1,x2,x3 = ndgrid(x1,x2,x3); z = x1 + exp(cos(2*x2.2) + sin(x3.3); slice(z,5 10 15, 10, 5 12); axis tight;Structures% Draw a visualization of a structure. strucdem_helper(1) = John Doe; patient.billing = 127.00; patient.test = 79 75 73; 180 178 177.5; 172 170 169; patientpatient = name: John Doe billing: 127 test: 3x3 double patient(2).name = Ann Lane; patient(2).billing = 28.50; patient(2).test = 68 70 68; 118 118 119; 172 170 169;% Update the visualization. strucdem_helper(2); fnames1 = fieldnames(patient) patient2 = rmfield(patient,test); fnames2 = fieldnames(patient2)fnames1 = name billing testfnames2 = name billing A = struct( data, 3 4 7; 8 0 1, 9 3 2; 7 6 5, . nest, . struct( testnum, Test 1, . xdata, 4 2 8, ydata, 7 1 6), . struct( testnum, Test 2, . xdata, 3 4 2, ydata, 5 0 9)% Update the visualization. strucdem_helper(3)聚类分析K(C)均值算法：x=-7.82 -4.58 -3.97;-6.68 3.16 2.17;4.36 -2.19 2.09;6.72 0.88 2.80;. -8.64 3.06 3.5;-6.87 0.57 -5.45;4.47 -2.62 5.76;6.73 -2.01 4.18;. -7.71 2.34 -6.33;-6.91 -0.49 -5.68;6.18 2.81 5.82;6.72 -0.93 -4.04;. -6.25 -0.26 0.56;-6.94 -1.22 1.13;8.09 0.20 2.25;6.81 0.17 -4.15;. -5.19 4.24 4.04;-6.38 -1.74 1.43;4.08 1.3 5.33;6.27 0.93 -2.78；idx,ctrs=kmeans(x,3)如何编写求K-均值聚类算法的Matlab程序在聚类分析中，K-均值聚类算法（k-means algorithm）是无监督分类中的一种基本方法，其也称为C-均值算法，其基本思想是：通过迭代的方法，逐次更新各聚类中心的值，直至得到最好的聚类结果。假设要把样本集分为c个类别，算法如下：（1）适当选择c个类的初始中心；（2）在第k次迭代中，对任意一个样本，求其到c个中心的距离，将该样本归到距离最短的中心所在的类，（3）利用均值等方法更新该类的中心值；（4）对于所有的c个聚类中心，如果利用（2）（3）的迭代法更新后，值保持不变，则迭代结束，否则继续迭代。下面介绍作者编写的一个分两类的程序，可以把其作为函数调用。% function samp1,samp2=kmeans(samp); 作为调用函数时去掉注释符samp=11.1506 6.7222 2.3139 5.9018 11.0827 5.7459 13.2174 13.8243 4.8005 0.9370 12.3576; %样本集l0 l=size(samp);%利用均值把样本分为两类，再将每类的均值作为聚类中心th0=mean(samp);n1=0;n2=0;c1=0.0;c1=double(c1);c2=c1;for i=1:l if samp(i)th0 c1=c1+samp(i);n1=n1+1; else c2=c2+samp(i);n2=n2+1; endendc1=c1/n1;c2=c2/n2; %初始聚类中心t=0;cl1=c1;cl2=c2;c11=c1;c22=c2; %聚类中心while t=0 samp1=zeros(1,l); samp2=samp1; n1=1;n2=1; for i=1:l if abs(samp(i)-c11) load hald covx=cov(ingredients) pc, latent, explained=pcacov(covx) pc, score, latent, tsquare= princomp (X)pc the principal component coefficientsThe scores are the data formed by transforming the original data into the space of theprincipal components.The values of the vector latent are the variance of the columns of SCORE.Tsquare contains Hotellings T2 statistic for each data point. clear load cities whoYour variables are:categories names ratings boxplot(ratings,0,+,0) set(gca,yticklabel,categories) stdr=std(ratings); %计算各列的标准差： sr=ratings./stdr(ones(329,1),:);%当原始数据的量级和量纲存在较大差异时，需要先对数据进行标准化，然后进行主成分析。标准化的方法是将原始数据的各列除以各列的标准差： pcs,newdata,var,ts=princomp(sr);% 现在寻找主要成分：第1个输出主成分pcs，是pp矩阵；第2个输出主成分得分(newdata)，主成分得分是原始数据在主成分所定义的新坐标系中的确定的数据，其大小与输入数据矩阵大小相同。下面我们看看newdata的前两列数据作为前两个主成分时的结果： x=newdata(:,1); y=newdata(:,2); plot(x,y,+) xlabel(第一主成分) ylabel(第二主成分)此时将在图像中生成一个十字交叉线，交点跟随鼠标移动。在散点附近单击，将标注该点的字符串。标注结束后，敲回车。结果显示如下：对于这些异常值我们可以直接删除，也即是将那几行的元素直接置空，比如New York对应第213行： rsubset=ratings; rsubset(213,:)=;第3个输出主成分方差(var)主成分方差var是有newdata的对应列所解释的包含方程的向量： percent_explained=100*var/sum(var)；% 计算每个主成分所解释的总方差的百分比： pareto(percent_explained)；%用帕累托图描述的百分数： xlabel(主成分) ylabel(方差解释)第4个输出Hotelling的T2检验（ts）：Hotelling的T2检验统计量是描述每一测量值与数据中心距离的统计量，用它可以找到数据中的极值点： st,ind=sort(ts); st=flipud(st); ind=flipud(ind); ex=ind(1); ex=ind(1)ex = 213 names(ex,:)ans =New York, NY 说明New York, NY是离数据中心最远的城市。用pls toolbox来做主成分分析可能会简单一些。XL,YL,XS,YS,BETA,PCTVAR,MSE = plsregress(X,Y,ncomp)邻近支持向量机（Proximal Support Vector Machine，PSVM）（数据挖掘基础教程，印度）范明，孟小峰译，P194）functionw,gamma=psvm(A,D,nu)% psvm linear & nonlinear classification% input:A,D,nu. output:w,gammaA=-1.7984 -1.6730 1;-1.0791 -0.5937 1;-0.5995 0.7556 1;1.0791 -1.4032 1;. 0.1199 0.2159 1;0.3597 0.4857 -1;-0.3597 1.5651 -1;0.5995 0.4857 -1;. 0.1199 -0.3238 -1;1.55876 0.4857 -1;d=eye(10,10);d1=d(1:5,:);d2=-1*d(6:10,:);D=d1;d2;nu=0.1;% w,gamma=psvm(A,D,nu)m,n=size(A);e=ones(m,1);H=D*A-e;r=sum(H); %r=H*e;r=(speye(n+1)/nu+H*H)r; %solve(I/nu+H*H)r=H*eu=nu*(1-(H*r);s=D*u;w=(s*A); %w=A*D*ugamma=-sum(s); %gamma=-e*D*ux=A;z=sign(w*x-gamma)Matlab的BP网络实例建立一个网络训练（已知的样本点）验证Load 装载数据创建网络net = newfit(houseInputs,houseTargets,20); 网络训练net=train(net,houseInputs,houseTargets);仿真网络sim注意： net = newcf(houseInputs,houseTargets,13);Warning: feval on script names will not work, or may work differently, in a future version of MATLAB. To make your code insensitive to any change and to suppress this warning message: - Either change the script to a function. - Or use eval instead of feval. The script file in question is tansig. In network.subsasgnsetLayerTransferFcn at 1272 In network.subsasgn at 189 In newcfnew_5p1 at 175 In newcf at 108? Error using = fevalAttempt to execute SCRIPT tansig as a function:C:Documents and SettingszhangMy DocumentsMATLABtansig.mError in = network.subsasgnsetLayerTransferFcn at 1272fpv = feval(transferFcn,fpdefaults);Error in = network.subsasgn at 189 net,err = setLayerTransferFcn(net,i,transferFcn);Error in = newcfnew_5p1 at 175 net.layersi.transferFcn = tansig;Error in = newcf at 108 net = new_5p1(varargin:); 发生该错误的原因是C:Documents and SettingszhangMy DocumentsMATLABtansig.m，tansig.m文件时MATLAB内部定义函数，在这里是自己编写的程序命名为tansig，自己编写的程序不能再以他命名。同样有下面的错误：net = newff(houseInputs,houseTargets,20);？Warning: feval on script names will not work, or may work differently, in a future version of MATLAB. To make your code insensitive to any change and to suppress this warning message: - Either change the script to a function. - Or use eval instead of feval. The script file in question is tansig. In network.subsasgnsetLayerTransferFcn at 1272 In network.subsasgn at 189 In newffnew_5p1 at 159 In newff at 89? Error using = fevalAttempt to execute SCRIPT tansig as a function:C:Documents and SettingszhangMy DocumentsMATLABtansig.mError in = network.subsasgnsetLayerTransferFcn at 1272fpv = feval(transferFcn,fpdefaults);Error in = network.subsasgn at 189 net,err = setLayerTransferFcn(net,i,transferFcn);Error in = newffnew_5p1 at 159 net.layersi.transferFcn = tansig;Error in = newff at 89 net = new_5p1(varargin:); Fitting a Function（训练前馈网络的第一步是建立网络对象。函数newff建立一个可训练的前馈网络。生成一个感知器神经网络net=newp(pr,s,tf,lf)；网络仿真sim: Y,Pf,Af=sim(net,P,Pi,Ai) 神经网络仿真函数；train神经网络训练函数BP神经网络的重要函数：咯logsig对数S型传递函数；purelin线性传递函数；simuff前向网络仿真；tansig正切S型传递函数；trainbp利用BP算法训练前向网络；trainbpx利用快速BP算法训练前向网络；trainlm利用Levenberg-Marquardt规则训练前向网络错误：net = newff(houseInputs,houseTargets,20);？Warning: feval on script names will not work, or may work differently, in a future version of MATLAB. To make your code insensitive to any change and to suppress this warning message: - Either change the script to a function. - Or use eval instead of feval. The script file in question is tansig. In network.subsasgnsetLayerTransferFcn at 1272 In network.subsasgn at 189 In newffnew_5p1 at 159 In newff at 89? Error using = fevalAttempt to execute SCRIPT tansig as a function:C:Documents and SettingszhangMy DocumentsMATLABtansig.mError in = network.subsasgnsetLayerTransferFcn at 1272fpv = feval(transferFcn,fpdefaults);Error in = network.subsasgn at 189 net,err = setLayerTransferFcn(net,i,transferFcn);Error in = newffnew_5p1 at 159 net.layersi.transferFcn = tansig;Error in = newff at 89 net = new_5p1(varargin:);clearclcclose allwarning off%数据输入huanghe_p=370 503 434 575 490 420 560 640 558 343 326 405 446 423 422 697 598 377 435 472 451 667 601 689 541 485 425 389 382 707 422;huanghe_t=515 713 586 753 720 567 717 987 810 489 453 589 639 568 595 982 849 519 615 652 599 941 893 999 758 701 630 561 520 1040 535;%归一化处理p=(huanghe_p-min(huanghe_p)/(max(huanghe_p)-min(huanghe_p);t=(huanghe_t-min(huanghe_t)/(max(huanghe_t)-min(huanghe_t);%网络有关参数EPOCHS=10000;GOAL=0.000005;LR=0.01;MAX_FAIL=100;%建立bp神经网络，并训练，仿真。其中输入为p，输出为t%-隐层神经元确定-s=3:15;%s 为常向量，表示神经元的个数res=zeros(size(s);%res将要存储误差向量，这里先置零pn=p(1:5);p(6:10);p(11:15);p(16:20);tn=t(1:5);t(6:10);t(11:15);t(16:20);for i=1:length(s)%7.0版本使用下面代码%输出层的神经元个数必须等于tn的行数%每层的传递函数可以自己指定%创建了2层网路，隐层s(i)(任意)；输出层4(由tn决定不可改，但需输入)%net=newff(minmax(pn),s(i),4,tansig,purelin,trainlm);%7.6版本使用下面的代码%输出层神经元由tn决定，不需要输入，我们只需要确定输入层和隐层的神经元个数%所有层的传递函数都可以自己定义%创建3层网络，输入层8(任意)；隐层s(i)(任意)；输出层4(由tn决定，不可改，不需要输入)%其实两个版本的区别在于7.6的bp网络创建函数直接给出了t，故最后一层不需要指定了，故相同的参数时，7.6的会多一层%net=newff(pn,tn,4,s(i),tansig,tansig,purelin,trainlm); net.iw1,1=zeros(size(net.iw1,1)+0.5;net.lw2,1=zeros(size(net.lw2,1)+0.75;net.b1,1=zeros(size(net.b1,1)+0.5;net.b2,1=zeros(size(net.b2,1);net.trainParam.epochs=EPOCHS;net.trainParam.goal =GOAL;net.trainParam.lr=LR;net.trainParam.max_fail=MAX_FAIL;net=train(net,pn,tn);y=sim(net,pn);e=tn-y;error=mse(e,net);res(i)=norm(error);end%选取最优神经元数，number为使得误差最小的隐层神经元个数ttmp2,ind=min(res);no=s(ind);%选定隐层神经元数目后，建立网络，训练仿真。%7.0版本%net=newff(minmax(pn),no,4,tansig,purelin,trainlm);%7.6版本net=newff(pn,tn,4,no,tansig,tansig,purelin,trainlm); net.iw1,1=zeros(size(net.iw1,1)+0.5;net.lw2,1=zeros(size(net.lw2,1)+0.75;net.b1,1=zeros(size(net.b1,1)+0.5;net.b2,1=zeros(size(net.b2,1);net.trainParam.epochs=EPOCHS;net.trainParam.goal =GOAL;net.trainParam.lr=LR;net.trainParam.max_fail=MAX_FAIL;net=train(net,pn,tn);y=sim(net,pn);e=tn-y;error=mse(e,net)%error为网络的误差向量r=norm(error);%r为网络的整体误差save net %保存最好的网络%预测input=p(11:15);p(16:20);p(21:25);p(26:30);yuce=sim(net,input);%结果反归一化y_norm=y(1, y(2, y(3, y(4,;yuce_norm=yuce(1, yuce(2, yuce(3, yuce(4,;%训练数据的仿真结果t_val=y_norm*(max(huanghe_t(1:20)-min(huanghe_t(1:20)+min(huanghe_t(1:20);%预测数据的仿真结果yuce=yuce_norm*(max(huanghe_t(1:20)-min(huanghe_t(1:20)+min(huanghe_t(1:20);%计算误差%训练数据仿真的相对误差wucha=abs(t_val-huanghe_t(1:20)./huanghe_t(1:20);b=minmax(wucha);average_wucha=mean(wucha);%作图figure(1)plot(1:20,huanghe_t(1:20),*-,1:20,t_val,o:)title(训练数据仿真结果)legend(原始数据,仿真结果)figure(2)plot(1:20,huanghe_t(11:30),*-,1:20,yuce,o:)title(预测数据仿真结果)legend(原始数据,仿真结果)对于BP网络，存在一个重要的结论，即单隐层的BP网络可以逼近任意的非线性映射，前提是隐含层的神经元个数可以随意调整。本文一个单隐层的BP神经网络设计为例，介绍里利用神经网络工具箱进行BP网络设计及分析的过程。此处，利用一个单隐层的BP网络来逼近一个函数。1. 问题描述通过对函数进行采样得到了网络的输入变量P和目标变量T，分别为： CODE:P = -1 : 0.1 : 1;T = -0.9602 -0.5770 -0.0729 0.3771 0.6405 0.6600 0.4609 0.1336 -0.2013 -0.4344 -0.5000 -0.3930 -0.1647 0.0988 0.3072 0.3960 0.3449 0.1816 -0.0312 -0.2189 -0.3201;每组向量都有21维数据，可以将输入向量和目标向量绘制在一起，如下图：P1.JPG2. 网络设计网络的输入层和输出层的神经元个数均为1，根据隐含层设计经验，解决该问题的网络的隐层神经元个数应该在38之间。下面设计一个隐含层神经元数目可变的BP网络，通过误差对比，确定最佳的隐含层神经元个数，并检查隐含层神经元个数对网络性能的影响。网络设计及训练代码如下： CODE:P = -1 : 0.1 : 1;T = -0.9602 -0.5770 -0.0729 0.3771 0.6405 0.6600 0.4609 0.1336 -0.2013 -0.4344 -0.5000 -0.3930