第一篇：?jiǎn)渭冃畏╩atlab程序
算法實(shí)現(xiàn)與分析
算法1．單純形法 具體算例：
minz=?3x1+x2+2x3 3x1+2x2?3x3=6 x1?2x2+x3+x5=4
x1,x2,x3≥0標(biāo)準(zhǔn)化后：
min z=?3x1+x2+2x3+Mx4+Mx5
3x1+2x2?3x3+x4=6 x1?2x2+x3+x5=4
x1,x2,x3,x4,x5≥0用單純形法求解，程序如下： clear clc
M=1000000;
A=[3,2,-3,1,0;1,-2,1,0,1];%系數(shù)矩陣 C=[-3,1,2,M,M,0];%價(jià)值矩陣 B=[6;4];Xt=[4 5];
for i=1:length(C)-1 D=0;
for j=1:length(Xt)
D=D+A(j,i)*C(Xt(j));
end
xi(i)=C(i)-D;end s=[];
for i=1:length(xi)
if xi(i)<0 s=[s,i];
end end
f=length(s);h=1;
while(f)
for k=1:length(s)j=1;A x=[];
for i=1:length(Xt)
if A(i,s(k))>0 x(j)=i;j=j+1;
end end x
if(length(x)+1==1)break;end y=1 x
for i=1:length(x)
if B(x(i))/A(x(i),s(k)) end end y=x(y);end y1=Xt(y);%??3?±?á? s k aa=A(y,s(k))%s(k)?a??è?±?á? A(y,:)=A(y,:)./aa;B(y,:)=B(y,:)./aa;z=[]; for i=1:length(Xt)z=[z,i];end z z(y)=[];z Xt for i=1:length(z);yz=-A(z(i),s(k)) A(z(i),:)=A(z(i),:)+A(y,:).*yz B(z(i))B(y)yz B(z(i))=B(z(i))+B(y).*yz end for i=1:length(Xt) if Xt(i)==y1 Xt(i)=s(k);break end end Xt disp('×a??oó')A=A B=B AB=[A,B]; for i=1:length(C)D=0; for j=1:length(Xt)D=D+AB(j,i)*C(Xt(j)); end xi(i)=C(i)-D; end xi s=[]; for i=1:length(xi)-1 if xi(i)<0 s=[s,i]; end end s vpa([A,B;C]);f=length(s);h=h+1; if h==5 break end end -xi(length(xi)) clear clc X=[1 2 3 4 5];A=[ 1 2 1 0 0;4 0 0 1 0;0 4 0 0 1];C=[2 3 0 0 0 ];b=[8;16;12];t=[3 4 5];B0=A(:,t);while 1 CB0=C(:,t);XN01=X; for i=1:length(t); for j=1:length(X); if XN01(j)==t(i) XN01(j)=0; end end end j=1; for i=1:length(X); if XN01(i)~=0 XN0(j)=XN01(i); j=j+1; end end for j=1:length(XN0); CN0(j)=C(XN0(j)); end N0=[]; for i=1:length(XN0); N0=[N0,A(:,XN0(i))]; end xiN0=CN0-CB0*B0*N0;j=1;z=[]; for i=1:length(xiN0) if xiN0(i)>0 z(j)=i; j=j+1; end end if length(z)+1==1; break; end n=1; for i=1:length(z) if z(i)>z(n) n=i; end end k=XN0(z(n));%換入變量 B=B0*b; P=B0*A(:,k);j=1; for i=1:length(P) if P(i)>0 x(j)=i; j=j+1; end end y=1; for i=1:length(x) if B(x(y))/P(x(y))>B(x(i))/P(x(i)) y=i; end end y1=x(y); y=t(y1);%換出變量 for i=1:length(t) if t(i)==y m=i; break; end end t(m)=k; P2=B0*A(:,k);q=P2(y1);P2(y1)=-1;P2=-P2./q; E=[1 0 0;0 1 0;0 0 1];E(:,m)=P2;B0=E*B0;end CB0*B0*b function [xx,fm]=myprgmh(m,n,A,b,c)B0=A(:,1:m);cb=c(:,1:m);xx=1:n;sgm=c-cb*B0^-1*A;h=-1;sta=ones(m,1);for i=m+1:n if sgm(i)>0 h=1; end end while h>0 [msg,mk]=max(sgm); for i=1:m sta(i)=b(i)/A(i,mk); end [mst,mr]=min(sta); zy=A(mr,mk); for i=1:m if i==mr for j=1:n A(i,j)=A(i,j)/zy; end b(i)=b(i)/zy; end end for i=1:m if i~=mr for j=1:n A(i,j)=A(i,j)-A(i,mk)*A(mr,j); end b(i)=b(i)-A(i,mk)*b(mr); end end B1=A(:,1:m); cb(mr)=c(mk); xx(mr)=mk; sgm=c-cb*B1*A; for i=m+1:n if sgm(i)>0 h=1; end end end fm=c*xx; 線性規(guī)劃單純形法matlab解法 %單純形法matlab程序-ssimplex % 求解標(biāo)準(zhǔn)型線性規(guī)劃:max c*x;s.t.A*x=b;x>=0 % 本函數(shù)中的A是單純初始表，包括:最后一行是初始的檢驗(yàn)數(shù)，最后一列是資源向量b % N是初始的基變量的下標(biāo) % 輸出變量sol是最優(yōu)解, 其中松弛變量（或剩余變量）可能不為0 % 輸出變量val是最優(yōu)目標(biāo)值，kk是迭代次數(shù) % 例：max 2*x1+3*x2 % s.t.x1+2*x2<=8 % 4*x1<=16 % 4*x2<=12 % x1,x2>=0 % 加入松馳變量，化為標(biāo)準(zhǔn)型，得到 % A=[1 2 1 0 0 8;% 4 0 0 1 0 16;% 0 4 0 0 1 12;% 2 3 0 0 0 0];% N=[3 4 5];% [sol,val,kk]=ssimplex(A,N)% 然后執(zhí)行 [sol,val,kk]=ssimplex(A,N)就可以了。function [sol,val,kk]=ssimplex(A,N)[mA,nA]=size(A);kk=0;% 迭代次數(shù) flag=1;while flag kk=kk+1;if A(mA,:)<=0 % 已找到最優(yōu)解 flag=0;sol=zeros(1,nA-1);for i=1:mA-1 sol(N(i))=A(i,nA);end val=-A(mA,nA);else for i=1:nA-1 if A(mA,i)>0&A(1:mA-1,i)<=0 % 問(wèn)題有無(wú)界解 disp('have infinite solution!');flag=0;break;end end if flag % 還不是最優(yōu)表，進(jìn)行轉(zhuǎn)軸運(yùn)算 temp=0;for i=1:nA-1 if A(mA,i)>temp temp=A(mA,i);inb=i;% 進(jìn)基變量的下標(biāo) end end sita=zeros(1,mA-1);for i=1:mA-1 if A(i,inb)>0 sita(i)=A(i,nA)/A(i,inb);end end temp=inf;for i=1:mA-1 if sita(i)>0&sita(i) A(outb,:)=A(outb,:)/A(outb,inb);for i=1:mA if i~=outb A(i,:)=A(i,:)-A(outb,:)*A(i,inb);End End End End end %******************************讀取數(shù)據(jù)************************************* MATLAB讀取數(shù)據(jù) Xlsread(‘lujing’,’mingcheng’)%******************************輸出數(shù)據(jù)************************************* MATLAB輸出Excel數(shù)據(jù) Xlswrite(‘eular.xla’,[A],’sheet 1’)%******************************三維曲面************************************* 三維曲面 x1=[0.05 0.05 0.05 0.1 0.05 0.1 0.2 0 0 0 0 0 0 0 0];x2=[0.1 0.25 0.12 0.12 0 0 0 0 0.05 0.08 0.1 0.12 0.15 0.2 0.25];y1=[153 130 146 133 160 154 140 133 165 169 171 186 175 156 152];y2=[7.22 5.57 6.66 6 6.6 8.21 5.1 4.66 7.58 8.26 8.99 8.73 7.71 7 6.49];y3=[112 100 131 88 80 72 72 80 116 83 80 80 64 68 74];scatter3(x1,x2,y2)figure [X,Y,Z]=griddata(x1,x2,y2,linspace(min(x1),max(x1))',linspace(min(x2),max(x2)),'v4');%插值 pcolor(X,Y,Z);shading interp%偽彩色圖 figure,contourf(X,Y,Z)%等高線圖 figure,surf(X,Y,Z);%三維曲面 %******************************各種距離************************************* pdist(X)— 樣本X中各n維向量的歐氏距離 pdist(X, 'cityblock')— 各n維向量的絕對(duì)距離 pdist(X, ‘minkowski'，r)— 閔可夫斯基距離 pdist(X, 'mahal')— 各n維向量的馬氏距離 mean(A)— A中各列向量的均值 Var(A)— A中各列向量的方差 Std(A)— A中各列向量的標(biāo)準(zhǔn)差 Cov(A)— A中各列向量的協(xié)方差矩陣 Corrcoef(A)— A中各列向量的相關(guān)矩陣 命令：[x,fval] = fminbnd(F,a,b)%******************************泰勒公式************************************* 計(jì)算F在a,b之間取極小值時(shí)的x與y(即fval).格式：taylor(F,a)功能：F在x=a處的泰勒級(jí)數(shù)前5項(xiàng) 格式：taylor(F,v)功能：F對(duì)變量v的泰勒展式前5項(xiàng) 格式：taylor(F,v,n)功能：求F的n 階泰勒展式,且 (n缺省時(shí)默認(rèn)為 5） %******************************一維插值************************************* x=[-2*pi:0.5*pi:2*pi];y=cos(x);xi=-2*pi:0.3*pi:2*pi;y_nearest=interp1(x,y, xi,‘nearest’)%最近鄰插值 y_linear= interp1(x,y,xi)%雙線性差值 y_spline= interp1(x,y,xi, ‘spline’)%三次樣條插值 y_cubic= interp1(x,y,xi, ‘cubic’)%雙三次插值 plot(x,y,'o',xi,y_nearest,'-',xi,y_linear, 'r* ', xi,y_spline,'k:',xi,y_cubic,'k-');legend('original data','nearest','linear','spline','cubic')%******************************曲線擬合************************************* x1=[2:16];y1=[6.42,8.2,9.58,9.5,9.7,10,9.93,9.99,10.49,10.59,10.6,10.8,10.6,10.9,10.76];b01=[0.1435,0.084]; %初始參數(shù)值 fun1=inline('x./(b(1)+b(2)*x)','b','x');% 定義函數(shù) [b1,r1,j1]=nlinfit(x1,y1,fun1,b01);y=x1./(0.1152+0.0845*x1);%根據(jù)b1寫(xiě)出具體函數(shù) plot(x1,y1,'*',x1,y,'-or');R2=1-sum((y-y1).^2)/sum((y1-mean(y1)).^2)%*****************************殘差與置信區(qū)間******************************** X=[ones(30,1), x1‘,x2’,x3‘]; %輸入自變量（注意1與自變量組成的矩陣）[b,bint,r,rint,s]=regress(y',X);%多元線性回歸 b,bint,s rcoplot(r,rint)%殘差及其置信區(qū)間作圖 a=[1,3:9,11:30];X1=X(a,:);y1=y(a);%剔除異常點(diǎn) [b1,bint1,r1,rint1,s1]=regress(y1',X1);b1,bint1,s1 %再次線性回歸 rcoplot(r1,rint1) %------------% Copula理論及其在matlab中的實(shí)現(xiàn)程序應(yīng)用實(shí)例 %------------ %******************************讀取數(shù)據(jù)************************************* % 從文件hushi.xls中讀取數(shù)據(jù) hushi = xlsread('hushi.xls');% 提取矩陣hushi的第5列數(shù)據(jù)，即滬市的日收益率數(shù)據(jù) X = hushi(:,5);% 從文件shenshi.xls中讀取數(shù)據(jù) shenshi = xlsread('shenshi.xls');% 提取矩陣shenshi的第5列數(shù)據(jù)，即深市的日收益率數(shù)據(jù) Y = shenshi(:,5); %****************************繪制頻率直方圖********************************* % 調(diào)用ecdf函數(shù)和ecdfhist函數(shù)繪制滬、深兩市日收益率的頻率直方圖 [fx, xc] = ecdf(X);figure;ecdfhist(fx, xc, 30);xlabel('滬市日收益率');% 為X軸加標(biāo)簽 ylabel('f(x)');% 為Y軸加標(biāo)簽 [fy, yc] = ecdf(Y);figure;ecdfhist(fy, yc, 30);xlabel('深市日收益率');% 為X軸加標(biāo)簽 ylabel('f(y)');% 為Y軸加標(biāo)簽 %****************************計(jì)算偏度和峰度********************************* % 計(jì)算X和Y的偏度 xs = skewness(X)ys = skewness(Y) % 計(jì)算X和Y的峰度 kx = kurtosis(X)ky = kurtosis(Y) %******************************正態(tài)性檢驗(yàn)*********************************** % 分別調(diào)用jbtest、kstest和lillietest函數(shù)對(duì)X進(jìn)行正態(tài)性檢驗(yàn) [h,p] = jbtest(X)% Jarque-Bera檢驗(yàn) [h,p] = kstest(X,[X,normcdf(X,mean(X),std(X))])% Kolmogorov-Smirnov檢驗(yàn) [h, p] = lillietest(X)% Lilliefors檢驗(yàn) % 分別調(diào)用jbtest、kstest和lillietest函數(shù)對(duì)Y進(jìn)行正態(tài)性檢驗(yàn) [h,p] = jbtest(Y)% Jarque-Bera檢驗(yàn) [h,p] = kstest(Y,[Y,normcdf(Y,mean(Y),std(Y))])% Kolmogorov-Smirnov檢驗(yàn) [h, p] = lillietest(Y)% Lilliefors檢驗(yàn) %****************************求經(jīng)驗(yàn)分布函數(shù)值******************************* % 調(diào)用ecdf函數(shù)求X和Y的經(jīng)驗(yàn)分布函數(shù) [fx, Xsort] = ecdf(X);[fy, Ysort] = ecdf(Y);% 調(diào)用spline函數(shù)，利用樣條插值法求原始樣本點(diǎn)處的經(jīng)驗(yàn)分布函數(shù)值 U1 = spline(Xsort(2:end),fx(2:end),X);V1 = spline(Ysort(2:end),fy(2:end),Y); % 調(diào)用ecdf函數(shù)求X和Y的經(jīng)驗(yàn)分布函數(shù) [fx, Xsort] = ecdf(X);[fy, Ysort] = ecdf(Y);% 提取fx和fy的第2個(gè)至最后一個(gè)元素，即排序后樣本點(diǎn)處的經(jīng)驗(yàn)分布函數(shù)值 fx = fx(2:end);fy = fy(2:end); % 通過(guò)排序和反排序恢復(fù)原始樣本點(diǎn)處的經(jīng)驗(yàn)分布函數(shù)值U1和V1 [Xsort,id] = sort(X);[idsort,id] = sort(id);U1 = fx(id);[Ysort,id] = sort(Y);[idsort,id] = sort(id);V1 = fy(id); %*******************************核分布估計(jì)********************************** % 調(diào)用ksdensity函數(shù)分別計(jì)算原始樣本X和Y處的核分布估計(jì)值 U2 = ksdensity(X,X,'function','cdf');V2 = ksdensity(Y,Y,'function','cdf'); % **********************繪制經(jīng)驗(yàn)分布函數(shù)圖和核分布估計(jì)圖********************** [Xsort,id] = sort(X);% 為了作圖的需要，對(duì)X進(jìn)行排序 figure;% 新建一個(gè)圖形窗口 plot(Xsort,U1(id),'c','LineWidth',5);% 繪制滬市日收益率的經(jīng)驗(yàn)分布函數(shù)圖 hold on plot(Xsort,U2(id),'k-.','LineWidth',2);% 繪制滬市日收益率的核分布估計(jì)圖 legend('經(jīng)驗(yàn)分布函數(shù)','核分布估計(jì)', 'Location','NorthWest');% 加標(biāo)注框 xlabel('滬市日收益率');% 為X軸加標(biāo)簽 ylabel('F(x)');% 為Y軸加標(biāo)簽 [Ysort,id] = sort(Y);% 為了作圖的需要，對(duì)Y進(jìn)行排序 figure;% 新建一個(gè)圖形窗口 plot(Ysort,V1(id),'c','LineWidth',5);% 繪制深市日收益率的經(jīng)驗(yàn)分布函數(shù)圖 hold on plot(Ysort,V2(id),'k-.','LineWidth',2);% 繪制深市日收益率的核分布估計(jì)圖 legend('經(jīng)驗(yàn)分布函數(shù)','核分布估計(jì)', 'Location','NorthWest');% 加標(biāo)注框 xlabel('深市日收益率');% 為X軸加標(biāo)簽 ylabel('F(x)');% 為Y軸加標(biāo)簽 %****************************繪制二元頻數(shù)直方圖***************************** % 調(diào)用ksdensity函數(shù)分別計(jì)算原始樣本X和Y處的核分布估計(jì)值 U = ksdensity(X,X,'function','cdf');V = ksdensity(Y,Y,'function','cdf');figure;% 新建一個(gè)圖形窗口 % 繪制邊緣分布的二元頻數(shù)直方圖，hist3([U(:)V(:)],[30,30])xlabel('U（滬市）');% 為X軸加標(biāo)簽 ylabel('V（深市）');% 為Y軸加標(biāo)簽 zlabel('頻數(shù)');% 為z軸加標(biāo)簽 %****************************繪制二元頻率直方圖***************************** figure;% 新建一個(gè)圖形窗口 % 繪制邊緣分布的二元頻數(shù)直方圖，hist3([U(:)V(:)],[30,30])h = get(gca, 'Children');% 獲取頻數(shù)直方圖的句柄值 cuv = get(h, 'ZData');% 獲取頻數(shù)直方圖的Z軸坐標(biāo) set(h,'ZData',cuv*30*30/length(X));% 對(duì)頻數(shù)直方圖的Z軸坐標(biāo)作變換 xlabel('U（滬市）');% 為X軸加標(biāo)簽 ylabel('V（深市）');% 為Y軸加標(biāo)簽 zlabel('c(u,v)');% 為z軸加標(biāo)簽 %***********************求Copula中參數(shù)的估計(jì)值****************************** % 調(diào)用copulafit函數(shù)估計(jì)二元正態(tài)Copula中的線性相關(guān)參數(shù) rho_norm = copulafit('Gaussian',[U(:), V(:)])% 調(diào)用copulafit函數(shù)估計(jì)二元t-Copula中的線性相關(guān)參數(shù)和自由度 [rho_t,nuhat,nuci] = copulafit('t',[U(:), V(:)]) %********************繪制Copula的密度函數(shù)和分布函數(shù)圖************************ [Udata,Vdata] = meshgrid(linspace(0,1,31));% 為繪圖需要，產(chǎn)生新的網(wǎng)格數(shù)據(jù) % 調(diào)用copulapdf函數(shù)計(jì)算網(wǎng)格點(diǎn)上的二元正態(tài)Copula密度函數(shù)值 Cpdf_norm = copulapdf('Gaussian',[Udata(:), Vdata(:)],rho_norm);% 調(diào)用copulacdf函數(shù)計(jì)算網(wǎng)格點(diǎn)上的二元正態(tài)Copula分布函數(shù)值 Ccdf_norm = copulacdf('Gaussian',[Udata(:), Vdata(:)],rho_norm);% 調(diào)用copulapdf函數(shù)計(jì)算網(wǎng)格點(diǎn)上的二元t-Copula密度函數(shù)值 Cpdf_t = copulapdf('t',[Udata(:), Vdata(:)],rho_t,nuhat);% 調(diào)用copulacdf函數(shù)計(jì)算網(wǎng)格點(diǎn)上的二元t-Copula分布函數(shù)值 Ccdf_t = copulacdf('t',[Udata(:), Vdata(:)],rho_t,nuhat);% 繪制二元正態(tài)Copula的密度函數(shù)和分布函數(shù)圖 figure;% 新建圖形窗口 surf(Udata,Vdata,reshape(Cpdf_norm,size(Udata)));% 繪制二元正態(tài)Copula密度函數(shù)圖 xlabel('U');% 為X軸加標(biāo)簽 ylabel('V');% 為Y軸加標(biāo)簽 zlabel('c(u,v)');% 為z軸加標(biāo)簽 figure;% 新建圖形窗口 surf(Udata,Vdata,reshape(Ccdf_norm,size(Udata)));% 繪制二元正態(tài)Copula分布函數(shù)圖 xlabel('U');% 為X軸加標(biāo)簽 ylabel('V');% 為Y軸加標(biāo)簽 zlabel('C(u,v)');% 為z軸加標(biāo)簽 % 繪制二元t-Copula的密度函數(shù)和分布函數(shù)圖 figure;% 新建圖形窗口 surf(Udata,Vdata,reshape(Cpdf_t,size(Udata)));% 繪制二元t-Copula密度函數(shù)圖 xlabel('U');% 為X軸加標(biāo)簽 ylabel('V');% 為Y軸加標(biāo)簽 zlabel('c(u,v)');% 為z軸加標(biāo)簽 figure;% 新建圖形窗口 surf(Udata,Vdata,reshape(Ccdf_t,size(Udata)));% 繪制二元t-Copula分布函數(shù)圖 xlabel('U');% 為X軸加標(biāo)簽 ylabel('V');% 為Y軸加標(biāo)簽 zlabel('C(u,v)');% 為z軸加標(biāo)簽 %**************求Kendall秩相關(guān)系數(shù)和Spearman秩相關(guān)系數(shù)*********************** % 調(diào)用copulastat函數(shù)求二元正態(tài)Copula對(duì)應(yīng)的Kendall秩相關(guān)系數(shù) Kendall_norm = copulastat('Gaussian',rho_norm)% 調(diào)用copulastat函數(shù)求二元正態(tài)Copula對(duì)應(yīng)的Spearman秩相關(guān)系數(shù) Spearman_norm = copulastat('Gaussian',rho_norm,'type','Spearman')% 調(diào)用copulastat函數(shù)求二元t-Copula對(duì)應(yīng)的Kendall秩相關(guān)系數(shù) Kendall_t = copulastat('t',rho_t)% 調(diào)用copulastat函數(shù)求二元t-Copula對(duì)應(yīng)的Spearman秩相關(guān)系數(shù) Spearman_t = copulastat('t',rho_t,'type','Spearman') % 直接根據(jù)滬、深兩市日收益率的原始觀測(cè)數(shù)據(jù)，調(diào)用corr函數(shù)求Kendall秩相關(guān)系數(shù) Kendall = corr([X,Y],'type','Kendall')% 直接根據(jù)滬、深兩市日收益率的原始觀測(cè)數(shù)據(jù)，調(diào)用corr函數(shù)求Spearman秩相關(guān)系數(shù) Spearman = corr([X,Y],'type','Spearman') %******************************模型評(píng)價(jià)************************************* % 調(diào)用ecdf函數(shù)求X和Y的經(jīng)驗(yàn)分布函數(shù) [fx, Xsort] = ecdf(X);[fy, Ysort] = ecdf(Y);% 調(diào)用spline函數(shù)，利用樣條插值法求原始樣本點(diǎn)處的經(jīng)驗(yàn)分布函數(shù)值 U = spline(Xsort(2:end),fx(2:end),X);V = spline(Ysort(2:end),fy(2:end),Y);% 定義經(jīng)驗(yàn)Copula函數(shù)C(u,v)C = @(u,v)mean((U <= u).*(V <= v));% 為作圖的需要，產(chǎn)生新的網(wǎng)格數(shù)據(jù) [Udata,Vdata] = meshgrid(linspace(0,1,31));% 通過(guò)循環(huán)計(jì)算經(jīng)驗(yàn)Copula函數(shù)在新產(chǎn)生的網(wǎng)格點(diǎn)處的函數(shù)值 for i=1:numel(Udata) CopulaEmpirical(i)= C(Udata(i),Vdata(i));end figure;% 新建圖形窗口 % 繪制經(jīng)驗(yàn)Copula分布函數(shù)圖像 surf(Udata,Vdata,reshape(CopulaEmpirical,size(Udata)))xlabel('U');% 為X軸加標(biāo)簽 ylabel('V');% 為Y軸加標(biāo)簽 zlabel('Empirical Copula C(u,v)');% 為z軸加標(biāo)簽 % 通過(guò)循環(huán)計(jì)算經(jīng)驗(yàn)Copula函數(shù)在原始樣本點(diǎn)處的函數(shù)值 CUV = zeros(size(U(:)));for i=1:numel(U) CUV(i)= C(U(i),V(i));end % 計(jì)算線性相關(guān)參數(shù)為0.9264的二元正態(tài)Copula函數(shù)在原始樣本點(diǎn)處的函數(shù)值 rho_norm = 0.9264;Cgau = copulacdf('Gaussian',[U(:), V(:)],rho_norm);% 計(jì)算線性相關(guān)參數(shù)為0.9325，自由度為4的二元t-Copula函數(shù)在原始樣本點(diǎn)處的函數(shù)值 rho_t = 0.9325;k = 4.0089;Ct = copulacdf('t',[U(:), V(:)],rho_t,k);% 計(jì)算平方歐氏距離 dgau2 =(CUV-Cgau)'*(CUV-Cgau)dt2 =(CUV-Ct)'*(CUV-Ct) %******************************灰色預(yù)測(cè)************************************* A=[ 96 144 194 276 383 466 554 652 747 832 972 169 235 400 459 565 695 805 881 1011 1139 151 238 335 425 541 641 739 866 975 1087 1238 164 263 376 531 600 711 913 1038 1173 1296 1497 182 318 445 576 708 856 1000 1145 1292 1435 1667 216 361 504 642 818 979 1142 1305 1479 1644 1920] m=size(A,2);x0=mean(A,2);x1=cumsum(x0);alpha=0.4;n=length(x0);z1=alpha*x1(2:n)+(1-alpha)*x1(1:n-1);Y=x0(2:n);B=[-z1,ones(n-1,1)];ab=BY x_A=(x0(1)-ab(2)/ab(1))*(exp(-ab(1)*n)-exp(-ab(1)*(n-1))) 圖標(biāo)參數(shù)************************************* y 黃-實(shí)線.點(diǎn)< 小于號(hào) m 紫： 點(diǎn)線o 圓s 正方形 c 青-.點(diǎn)劃線x 叉號(hào)d 菱形 r 紅--虛線+ 加號(hào)h 六角星 g 綠* 星號(hào)p 五角星 b 藍(lán)v 向下三角形 w 白^ 向上三角形 k 黑> 大于號(hào) %******************************圖形修飾************************************* 函數(shù) 含義 grid on(/off)給當(dāng)前圖形標(biāo)記添加（取消）網(wǎng)絡(luò) xlable(‘string’)標(biāo)記橫坐標(biāo) ylabel(‘string’)標(biāo)記縱坐標(biāo) title(‘string’)給圖形添加標(biāo)題 text(x,y,’string’)在圖形的任意位置增加說(shuō)明性文本信息 gtext(‘string’)利用鼠標(biāo)添加說(shuō)明性文本信息 axis([xmin xmax ymin ymax])設(shè)置坐標(biāo)軸的最小最大值 例5.1.2 給例5.1.1 的圖形中加入網(wǎng)絡(luò)和標(biāo)記。（見(jiàn)圖5.1.3 和5.1.4）>> x=0:pi/10:2*pi;>> y1=sin(x);>> y2=cos(x);>> plot(x,y1,x,y2)>> grid on >> xlabel('independent variable X')>> ylabel('Dependent Variable Y1 & Y2')>> title('Sine and Cosine Curve')>> text(1.5,0.3,'cos(x)')>> gtext('sin(x)')>> axis([0 2*pi-0.9 0.9]) 圖5.1.3 使用了圖形修飾的plot 函數(shù)繪制的正弦曲線 例5.2.1 繪制方程x=t y=sin(t)z=cos(t) 在t=[0,2*pi]上的空間方程。（見(jiàn)圖5.2.1） >> clf >> x=0:pi/10:2*pi;>> y1=sin(x);>> y2=cos(x);>> plot3(y1,y2,x,'m:p')>> grid on >> xlabel('Dependent Variable Y1')>> ylabel('Dependent Variable Y2')>> zlabel('Independent Variable X')>> title('Sine and Cosine Curve') 圖5.2.1 函數(shù)plot 繪制的三維曲線圖 %******************************三維圖形************************************* 例5.2.3 繪制方程 sin((x^2+y^2)^(1/2))z =---------------------(x^2+y^2)^(1/2) 在x∈[-7.5,7.5];y∈[-7.5,7.5] 的圖形。 >> x=-7.5:0.5:7.5;y=x;>> [X,Y]=meshgrid(x,y);>> R=sqrt(X.^2+Y.^2)+eps;>> Z=sin(R)./R;>> surf(X,Y,Z)>> xlabel('X 軸方向')>> ylabel('Y 軸方向')>> zlabel('Z 軸方向')(見(jiàn)圖5.2.4) 圖5.2.4 例5.2.4 繪制由方程形成的立體圖。（見(jiàn)圖5.2.5）z=x*exp(-(x^2+y^2)) >> clear >> x=-2:0.1:2;y=x;>> [X,Y]=meshgrid(x,y);>> Z=X.*exp(-X.^2-Y.^2);>> surf(X,Y,Z) 圖5.2.5 %******************************三維多角度觀察********************************* MTALAB 允許用戶設(shè)置觀察點(diǎn)，其指令是： view(azimuth,elevation)其中方位角azimuth 是觀察點(diǎn)和坐標(biāo)原點(diǎn)連線在x-y平面的投影和y 軸負(fù)方向的夾角，仰 角elevation 是觀察點(diǎn)與坐標(biāo)原點(diǎn)的連線和x-y平面的夾角。對(duì)于這兩個(gè)角度，三維圖形的 默認(rèn)值分別是-37.5 和30，二維圖形的默認(rèn)值是0 和90。 例5.2.5 從不同的角度觀察高斯矩陣的曲面。 >> z=peaks(40);>> subplot(2,2,1);>> mesh(z);>> subplot(2,2,2);>> mesh(z);>> view(-37.5,-30);>> subplot(2,2,3);>> mesh(z);>> view(180,0);>> subplot(2,2,4);>> mesh(z);>> view(0,90); 圖5.2.6 對(duì)應(yīng)不同觀察點(diǎn)的三維曲面圖 %******************************其他圖形************************************* 除了plot 繪圖函數(shù)以外，在有些場(chǎng)合對(duì)繪制的曲線會(huì)有一些特殊要求，這就要其他函 數(shù)來(lái)實(shí)現(xiàn)，常用的幾種函數(shù)如下（見(jiàn)表5.3.1） 表5.3.1 其他圖形函數(shù)表 函數(shù) 含義 loglog 使用對(duì)數(shù)坐標(biāo)系繪圖 semilogx 橫坐標(biāo)為對(duì)數(shù)坐標(biāo)軸，縱坐標(biāo)為線性坐標(biāo)軸 semilogy 橫坐標(biāo)為線性坐標(biāo)軸，縱坐標(biāo)為對(duì)數(shù)坐標(biāo)軸 polar 繪制極坐標(biāo)圖 fill 繪制實(shí)心圖 bar 繪制直方圖 pie 繪制餅圖 area 繪制面積圖 quiver 繪制向量場(chǎng)圖 stairs 繪制階梯圖 sterm 繪制火柴桿圖 例5.3.1 >> x=0:pi/10:2*pi;>> y1=sin(x);>> subplot(2,2,1);>> plot(x,y1);>> subplot(2,2,2);>> bar(x,y1);>> subplot(2,2,3);>> fill(x,y1,'g');>> subplot(2,2,4);>> stairs(x,y1,'k'); %******************************直方圖************************************* 函數(shù)bar(x)可以繪制直方圖，這對(duì)統(tǒng)計(jì)或者數(shù)據(jù)采集非常直觀實(shí)用。它共有四種形式： bar,bar3,barh 和bar3h，其中bar 和bar3 分別用來(lái)繪制二維和三維豎直方圖，barh 和b ar3h 分別用來(lái)繪制二維和三維水平直方圖，調(diào)用格式是： bar(x,y)其中x 必須單調(diào)遞增或遞減，y 為n m× 矩陣，可視化結(jié)果為m 組，每 組n 個(gè)垂直柱，也就是把y 的行畫(huà)在一起，同一列的數(shù)據(jù)用相同的顏色表示； bar(x,y,width)（或bar(y,width)）指定每個(gè)直方條的寬度，如width>1，則直方條會(huì)重 疊，默認(rèn)值為width=0.8； bar(…,’grouped’)使同一組直方條緊緊靠在一起； bar(…,’stack’)把同一組數(shù)據(jù)描述在一個(gè)直方條上。 例5.3.2 >> y=[5 3 2 9;4 7 2 7;1 5 7 3];>> subplot(2,2,1),bar(y)>> x=[5 9 11];>> subplot(2,2,2),bar3(x,y)>> subplot(2,2,3),bar(x,y,'grouped')>> subplot(2,2,4),bar(rand(2,3),.75,'stack') %******************************面積圖************************************* 5.3.2 面積圖 函數(shù)area 用來(lái)繪制面積圖，面積圖在plot 的基礎(chǔ)上填充x 軸和曲線之間的面積，該圖 用于查看某個(gè)數(shù)在該列所有數(shù)的總和中所占的比例。 例5.3.3 >> x=-3:3;>> y=[3 2 5;6 1 8;7 4 9;6 3 7;8 2 9;4 2 9;3 1 7];>> area(x,y) %******************************餅圖************************************* 5.3.3 餅圖 函數(shù)pie 用來(lái)繪制餅圖，它可以形象地表示出向量中各元素所占比例。其調(diào)用格式是： pie(x)x 中的元素通過(guò)x/sum(x)進(jìn)行歸一化，以確定餅圖中的份額； pie(x,explode)向量explode 和x 元素?cái)?shù)相同，用來(lái)指出需要分開(kāi)的餅片，explode 中 不為零的部分會(huì)被分開(kāi)。 例5.3.4 設(shè)某班的某課程的考試成績(jī)?nèi)缦拢?0 分以上有32 人，81 至90 有58 人，71 至80 分有27 人，60 至70 分為21 人，60 分以下有16 人，畫(huà)出餅圖。（見(jiàn)圖5.3.4） >> x=[32 58 27 21 16];>> explode0=[1 0 0 0 0];>> subplot(1,2,1)>> pie(x,explode0)>> explode1=[0 0 0 0 1];>> subplot(1,2,2)>> pie(x,explode1) 5.3.4 不同坐標(biāo)系中的繪圖 Semilogx,semilogy,loglo,polar(theta,rho)的使用方法和plot 完全類(lèi)似，不同的只是繪 制到不同的圖形坐標(biāo)上。函數(shù)semilogx 繪制x 軸為對(duì)數(shù)標(biāo)度的圖形，在半對(duì)數(shù)坐標(biāo)系中繪圖； 函數(shù)semilogy 繪制y 軸為對(duì)數(shù)標(biāo)度的圖形；函數(shù)loglog 繪制兩個(gè)軸都為對(duì)數(shù)間隔的圖形； 函數(shù)polar(theta,rho)繪制極坐標(biāo)圖形，其中theta 為相角，rho 為其對(duì)應(yīng)的半徑。 例5.3.5 繪制ρ=acos(3θ),a=2 的圖形。（見(jiàn)圖5.3.5） >> theta=-pi:pi/80:pi;>> polar(theta,2*cos(3*theta)) 圖5.3.5 極坐標(biāo)圖 5.4 符號(hào)表達(dá)式繪圖 MATLAB 軟件提供了將表達(dá)式進(jìn)行圖形顯示的功能。完成此功能需調(diào)用fplot 函數(shù)和 ezplot 函數(shù)。 函數(shù)fplot 用來(lái)繪制數(shù)學(xué)函數(shù)，其調(diào)用格式為： fplot(fun,lims)其中fun 就是所要繪制的函數(shù)，可以是定義函數(shù)的M 文件名，也可以是以x 為變量的可計(jì) 算字符串。例如’diric(x,10)’或’[sin(x),cos(x)]’，對(duì)于向量x 的每個(gè)元素，函數(shù) fun(x)必須返回一個(gè)行向量。如果fun 返回[f1(x),f2(x),f3(x)]，輸入[x1;x2]，就會(huì)返回矩陣 f1(x1)f2(x1)f3(x1)f1(x2)f2(x2)f3(x2)lims=[XMIN XMAX YMIN YMAX]限定了x,y 軸上的繪圖空間。 例5.4.1 >> subplot(2,2,1),fplot('humps',[0 1])>> subplot(2,2,2),fplot('abs(exp(-j*x*(0:9))*ones(10,1))',[0 2*pi])>> subplot(2,2,3),fplot('[tan(x),sin(x),cos(x)]',2*pi*[-1 1-1 1])>> subplot(2,2,4),fplot('sin(1./x)',[0.01 0.1],1e-3) 圖5.4.1 fplot 函數(shù)繪制表達(dá)式圖形 ezplot 函數(shù)是簡(jiǎn)捷繪圖指令之一，它無(wú)需數(shù)據(jù)準(zhǔn)備，直接畫(huà)出函數(shù)圖形，基本調(diào)用格式 為ezplot(f)其中f 是字符串或代表數(shù)學(xué)函數(shù)的符號(hào)表達(dá)式，只有一個(gè)符號(hào)變量，可以是x，缺省情況下 x 軸的繪圖區(qū)域?yàn)閇-π, π ]，但我們可以用ezplot(f,xmin,xmax)或ezplot(f,[xmin,xmax])來(lái)指定x 的范圍。 例5.4.2 >> y='x^2';>> subplot(1,2,1)>> ezplot(y)>> subplot(1,2,2)>> y='sin(x)';>> ezplot(y,[-pi,pi]) 圖5.4.2 ezplot 函數(shù)繪制表達(dá)式圖形 5.5 plot 函數(shù) MATLAB 對(duì)數(shù)據(jù)是按列存儲(chǔ)和計(jì)算的，運(yùn)用plot(x)時(shí)，當(dāng)x 為一個(gè)向量時(shí)，以其元 素為縱坐標(biāo)，其序號(hào)為橫坐標(biāo)值繪制曲線。當(dāng)x 為實(shí)矩陣時(shí)，則以其序號(hào)為橫坐標(biāo)，按列 繪制每列元素相對(duì)于序號(hào)的曲線，當(dāng)x 為n m× 矩陣時(shí)，就有n 條曲線。 如果x,y 是同維向量，plot(x,y)指令以x 元素為橫坐標(biāo)值，y 元素為縱坐標(biāo)值繪制曲線。 如x 是向量，y 是有一維與x 元素?cái)?shù)量相等的矩陣，則以x 為共同橫坐標(biāo)，按列繪制y 每 列元素值，曲線數(shù)為y 的另一維的元素?cái)?shù)。如果x,y 是同維矩陣，則以x,y 對(duì)應(yīng)列元素為、縱坐標(biāo)分別繪制曲線，數(shù)目等于矩陣的列數(shù)。 例5.5.1 >> x=[3 5 10 8];>> subplot(2,2,1)>> plot(x)>> x=[3 5 10 8;7 2 9 4;2 7 2 7]';>> subplot(2,2,2)>> plot(x)>> x=[3 5 6 8];>> y=[1 5 10 4];>> subplot(2,2,3)>> plot(x,y)>> x=[1 3 5 7;2 4 6 8]';>> y=[6 2 5 10;3 5 2 6]';>> subplot(2,2,4)>> plot(x,y,'k:*')第二篇：改進(jìn)單純形法matlab程序
第三篇：運(yùn)籌學(xué)單純形法matlab程序
第四篇：線性規(guī)劃單純形法matlab解法
第五篇：MATLAB程序總結(jié)
點(diǎn)擊咨詢 