数学建模常用模型matlab程序总结
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2005年A题 长江水污染
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
A
%灰色预测
clc,clear;
%X=[6.7 7.3 7.8 5.1 4.6 7.8];%0
%
%X=[9.6 13.3 12.8 25.4 26.5 22.5 23.5];
X=[174 179 183 207 234 220.5 256 270 285];
%X=log(x+1);
n=length(X);%计算矩阵长度
x1=cumsum(X);%累加后的矩阵
for i=2:n
rou(i)=X(i)/x1(i-1);%光滑性检验
end
for i=2:n
sigma(i)=x1(i)/x1(i-1);%准指数检验
end
i=2:n;
y(i-1)=X(i);
y1=y';%将原行向量转化成对应列向量
i=1:n-1; %计算数据矩阵B的第一列数据
C(i)=-0.5*(x1(i)+x1(i+1));
B=[C' ones(n-1,1)];%构造矩阵B
au=((B'*B)\B')*y1; %计算参数a,u矩阵
au;
a=au(1);
u=au(2);
i=1:n+10; %计算预测累加数列的值
ago(i)=(X(1)-u/a)*exp(-a*(i-1))+u/a;%累加数列值
yc(1)=ago(1);
i=1:n+10-1; %还原数列的值
yc(i+1)=ago(i+1)-ago(i);
i=2:n;
error(i)=yc(i)-X(i);%计算 残差值
yc(1)=ago(1);
i=1:n;%修正的还原数列值 !!!!!
yc(i+1)=ago(i+1)-ago(i);%预测值
c=std(error)/std(X);%计算后验差比
p=0;
for i=1:2n
if(abs(error(i)-mean(error))<0.6745*std(X))
p=p+1;
end
end
p=p/(n-1);%小误差概率
w1=min(abs(error));
w2=max(abs(error));
i=1:n;
w(i)=(w1+0.5*w2)./(abs(error(i))+0.5*w2);
w=sum(w)/(n);%关联度
disp('光滑性检验')
rou%光滑性检验
disp('准指数检验')
sigma%准指数检验
disp('参数a,u')
au %参数a,u
disp('累加数列ago的值')
%ago %累加数列ago的值
x1
disp('原始序列')
X %原始序列
disp('预测数值')
yc %预测数值
disp('残差值')
error%残差值
disp('后验差比值')
c %后验差比值
disp('小误差概率的值')
p %小误差概率的值
disp('关联度')
w %关联度w
QQ=abs((error./X)*100)%相对误差百分比
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
B TOPSIS 评价模型
clc,clear;
%读取表格中的相关数据
num2=xlsread('changjiang.xlsx','A2:A18');% DO
num3=xlsread('changjiang.xlsx','B2:B18');% CODMn
num4=xlsread('changjiang.xlsx','C2:C18');% NH3-N
b=max(num2)-num2;%1./num2; %对DO对应数据进行去倒数处理
r2=(b-min(b))./(max(b)-min(b));
r3=(num3-min(num3))./(max(num3)-min(num3));
r4=(num4-min(num4))./(max(num4)-min(num4));
R=[r2 r3 r4]; %将所有计算结果返回矩阵R
%构建加权规范决策矩阵Z
Q=[0.02 0.15 0.40 0.72 0.91 1];
%对第2列数据处理
b2=size(num2);
b2(1);
for i=1:1:b2(1);
if num2(i)>=7.5;
num2(i)=Q(1);
elseif num2(i)>=6&&num2(i)<7;
num2(i)=Q(2);
elseif num2(i)>=5&&num2(i)<6;
num2(i)=Q(3);
elseif num2(i)>=3&&num2(i)<5;
num2(i)=Q(4);
elseif num2(i)>=2&&num2(i)<3;
num2(i)=Q(5);
elseif num2(i)>=0&&num2(i)<2;
num2(i)=Q(6);
end
end
num2;
%对第3列数据处理
b3=size(num3);
b3(1);
for i=1:1:b3(1);
if n
um3(i)<=2;
num3(i)=Q(1);
elseif num3(i)>2&&num3(i)<=4;
num3(i)=Q(2);
elseif num3(i)>4&&num3(i)<=6;
num3(i)=Q(3);
elseif num3(i)>6&&num3(i)<=10;
num3(i)=Q(4);
elseif num3(i)>10&&num3(i)<=15;
num3(i)=Q(5);
elseif num3(i)>15;
num3(i)=Q(6);
end
end
num3;
%对第4列数据进行处理
b4=size(num4);
b4(1);
for i=1:1:b4(1);
if num4(i)<=0.15;
num4(i)=Q(1);
elseif num4(i)>0.15&&num4(i)<=0.5;
num4(i)=Q(2);
elseif num4(i)>0.5&&num4(i)<=1;
num4(i)=Q(3);
elseif num4(i)>1&&num4(i)<=1.5;
num4(i)=Q(4);
elseif num4(i)>1.5&&num4(i)<=2;
num4(i)=Q(5);
elseif num4(i)>2;
num4(i)=Q(6);
end
end
num4;
W=[num2 num3 num4];%[num1 num2 num3 num4];%权重矩阵
%确定正、负理想解
r11=min(R(:,1));%正理想解
r12=min(R(:,2));
r13=min(R(:,3));
%步骤五：计算距离
d11=(R(:,1)-r11).^2;
d12=(R(:,2)-r12).^2;
d13=(R(:,3)-r13).^2;
d01=[d11 d12 d13];
P1=W.*d01;%得到评价矩阵
Pingjia=(sum(P1'))';%最终评价矩阵17*1
Pingjia
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
C 时间序列模型
clc,clear;
A=[...];%待处理数据
a=A';
figure
subplot(1,2,1)
autocorr(A)
title('自相关函数图')
subplot(1,2,2)
parcorr(A)
title('偏相关函数图')
da=a;
n=length(da);
for i=0:2
for j=0:3
spec=garchset('R',i,'M',j,'Display','off');
[coeffX,errorsX,LLFX]=garchfit(spec,da);
num=garchcount(coeffX);
[aic,bic]=aicbic(LLFX,num,n);
fprintf('R=%d,M=%d,AIC=%f,BIC=%f\n',i,j,aic,bic);
end
end
%计算AIC的循环可以省去不要
r=input('输入阶数R=');
m=input('输入阶数M=');
spec2=garchset('R',r,'M',m,'Display','off');
[coeffX,errorsX,LLFX]=garchfit(spec2,da);
[sigmaForecast,w_Forecast]=garchpred(coeffX,da,5);
w_Forecast
x_pred=a(end)+cumsum(w_Forecast)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2006年B题 艾滋病
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
投入产出模型
A DEA 1
clc,clear;
x=xlsread('111.xls','C13:F32');%多目标输入
y=xlsread('111.xls','G13:G32');%多目标产出
X=x';
Y=y';
n=size(X',1);
m=size(X,1);
s=size(Y,1);
A=[-X' Y'];
b=zeros(n,1);
LB=zeros(m+s,1);
UB=[];
for i=1:n;
f=[zeros(1,m) -Y(:,i)'];
Aeq=[X(:,i)' zeros(1,s)];
beq=1;
w(:,i)=linprog(f,A,b,Aeq,beq,LB,UB);%解线性规划，得到DMUi的最佳权向量wi
E(i,i)=Y(:,i)'*w(m+1:m+s,i);%求出DMUi的相对效率值Eii
end
for i=1:20
xiaolv(i)=E(i,i);
end
w;%输出最佳权向量
E;%输出相对效率值
xiaolv
omega=w(1:m,:);%输出投入权向量
mu=w(m+1:m+s,:);%输出产出权向量
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
B DEA 2
clc,clear;
X=[...];%多指标输入矩阵X
Y=[...];%多指标输出矩阵Y
n=size(X',1);
m=size(X,1);
s=size(Y,1);
epsilon=10^-10;%定义非阿基米德无
穷小
f=[zeros(1,n) -epsilon*ones(1,m+s) 1];
A=zeros(1,n+m+s+1);
b=0;
LB=zeros(n+m+s+1,1);
UB=[];
LB(n+m+s+1)=-Inf;
for i=1:n;
Aeq=[X eye(m) zeros(m,s) -X(:,i);
Y zeros(s,m) -eye(s) zeros(s,1)];
beq=[zeros(m,1);
Y(:,i)];
w(:,i)=linprog(f,A,b,Aeq,beq,LB,UB);%解线性规划，的DMUi的最佳权向量wi
end
w%输出最佳权向量
lambda=w(1:n,:)%输出 .*
s_minus=w(n+1:n+m,:)%输出S*-
s_plus=w(n+m+1:n+m+s,:)%输出S*+
theta=w(n+m+s+1,:)%输出 .*
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2011年A题 重金属
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
A 计算灰色关联度
clc,clear;
A=xlsread('三问.xlsx','N2:O315');
B=[3.6 35];%母序列
[m,n]=size(A);
derta=zeros(m,n);
for i=1:n
derta(:,i)=abs(A(:,i)-B(i));
end
derta
for i=1:n
a(:,i)=min(derta(:,i));
b(:,i)=max(derta(:,i));
a=min(a);
b=max(b);
end
a
b
v=(a+0.5*b)./(derta+0.5*b)
Q=[0.0042 0.3341];%权重
q=Q./sum(Q);
for i=1:n
U(:,i)=v(:,1)*q(i);
end
u=sum(U');
u'
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
B 绘制重金属污染物的空间分布图（三维立体图）
clc
clear
close all
D=xlsread('2011A附件_数据.xls','A2:F314');
x=D(:,2);
y=D(:,3);
z=D(:,4);
C=D(:,6);
c=D(:,5);
xi=linspace(min(x),max(x),100); %进行插值
yi=linspace(min(y),max(y),100);
[xi,yi]=meshgrid(xi,yi);
zi=griddata(x,y,z,xi,yi); %不规则数据的曲面插值
ci=griddata(x,y,c,xi,yi);
marker={'d','o','s','^','h'};
color={'k','r','y','c','w'};
mat={'As','Cd','Cr','Cu','Hg','Ni','Pb','Zn'};
str={'等高线','生活区','工业区','山林区','交通区','绿地区'};
% 三维立体图
figure
h=surf(xi,yi,zi);
set(h,'cdata',ci);
colormap cool
title('三维立体图（颜色条表示分类）')
xlabel('X')
ylabel('Y')
colorbar
hidden off
hold on
for i=1:5
loc=c==i;
hold on
for j=1:5
Loc=C==j;
plot3(x(Loc),y(Loc),z(Loc),marker{j},'markerfacecolor',color{j});
end
str{1}='三维图';
legend(str,'location','best')
end
%%%%%%%%%%%%%%%%%%%%%%%%
C 模糊综合评价程序（程序可以运行，但是太过冗长亟待改进）
%模糊综合评价程序（降半梯形法）
clc,clear;
DATA=[1.73 2.17 1.72 3.6 2.62 1.46 2.11 3.17
2 2.88 1.51 4.46 7.86 1.56 2.74 3.56
1.12 1.17 1.26 1.31 1.17 1.26 1.18 1.06
1.58 2.67 1.67 3.94 3.59 1.35 1.99 3.26
1.74 2.16 1.41 2.29 3.29 1.24 1.96 2.24];
[m,n]=size(DATA);
for i=5
C1=DATA(i,1);
C2=DATA(i,2);
C3=DATA(i,3);
C4=DATA(i,4);
C5=DATA(i,5);
C6=DATA(i,6);
C7=DATA(i,7);
C8=DATA(i,8);
%第一个地区 As1
if C1<=1;
r11=1;
elseif C1>1&C1<2;
r11=(2-C1);
else C1>2;
r11=0;
end;
if C1<=1|C1>=3;
r12=0;
elseif C1>1&C1<2;
r12=C1-1;
else C1>2&C1<3;
r12=3-C1;
end;
if C1>3;
r13=0;
elseif C1>2&C1<3;
r13=(C1-2);
else C1<=2;
r13=0;
end;
if C1>=3;
r14=1;
else C1<3;
r14=0;
end;
R1=[r11 r12 r13 r14];
%第一个地区 Cd2
if C2<=1;
r21=1;
elseif C2>1&C2<2;
r21=(2-C2);
else C3>2;
r21=0;
end;
if C2<=1|C2>=3;
r22=0;
elseif C2>1&C2<2;
r22=C2-1;
else C2>2&C2<3;
r22=3-C2;
end;
if C2>3;
r23=0;
elseif C2>2&C2<3;
r23=(C2-2);
else C2<=2;
r23=0;
end;
if C2>=3;
r24=1;
else C2<3;
r24=0;
end;
R2=[r21 r22 r23 r24];
%第一个地区 Cr3
if C3<=1;
r31=1;
elseif C3>1&C3<2;
r31=(2-C3);
else C3>2;
r31=0;
end;
if C3<=1|C3>=3;
r32=0;
elseif C3>1&C3<2;
r32=C3-1;
else C3>2&C3<3;
r32=3-C3;
end;
if C3>3;
r33=0;
elseif C3>2&C3<3;
r33=(C3-2);
else C3<=2;
r33=0;
end;
if C3>=3;
r34=1;
else C3<3;
r34=0;
end;
R3=[r31 r32 r33 r34];
%第一个地区 Cu4
if C4<=1;
r41=1;
elseif C4>1&C4<2;
r41=(2-C4);
else C4>2;
r41=0;
end;
if C4<=1|C4>=3;
r42=0;
elseif C4>1&C4<2;
r42=C4-1;
else C4>2&C4<3;
r42=3-C4;
end;
if C4>3;
r43=0;
elseif C4>2&C4<3;
r43=(C4-2);
else C4<=2;
r43=0;
end;
if C4>=3;
r44=1;
else C4<3;
r44=0;
end;
R4=[r41 r42 r43 r44];
%第一个地区 Hg5
if C5<=1;
r51=1;
elseif C5>1&C5<2;
r51=(2-C5);
else C5>2;
r51=0;
end;
if C5<=1|C5>=3;
r52=0;
elseif C5>1&C5<2;
r52=C5-1;
else C5>2&C5<3;
r52=3-C5;
end;
if C5>3;
r53=0;
elseif C5>2&C5<3;
r53=(C5-2);
else C5<=2;
r53=0;
end;
if C5>=3;
r54=1;
else C5<3;
r54=0;
end;
R5=[r51 r52 r53 r54];
%第一个地区 Ni6
if C6<=1;
r61=1;
elseif C6>1&C6<2;
r61=(2-C6);
else C6>2;
r61=0;
end;
if C6<=1|C6>=3;
r62=0;
elseif C6>1&C6<2;
r62=C6-1;
else C6>2&C6<3;
r62=3-C6;
end;
if C6>3;
r63=0;
elseif C6>2&C6<3;
r63=(C6-2);
else C6<=2;
r63=0;
end;
if C6>=3;
r64=1;
else C6<3;
r64=0;
end;
R6=[r61 r62 r63 r64];
%第一个地区 Pb7
if C7<=1;
r71=1;
elseif C7>1&C7<2;
r71=(2-C7);
else C7>2;
r71=0;
end;
if C7<=1|C7>=3;
r72=0;
elseif C7>1&C7<2;
r72=C7-1;
else C7>2&C7<3;
r72=3-C7;
end;
if C7>3;
r73=0;
elseif C7>2&C7<3;
r73=(C7-2);
else C7<=2;
r73=0;
end;
if C7>=3;
r74=1;
else C7<3;
r74=0;
end;
R7=[r71 r72 r73 r74];
%第一个地区 Zn8
if C8<=1;
r81=1;
elseif C8>1&C8<2;
r81=(2-C8);
else C8>2;
r81=0;
end;
if C8<=1|C8>=3;
r82=0;
elseif C8>1&C8<2;
r82=C8-1;
else C8>2&C8<3;
r82=3-C8;
end;
if C8>3;
r83=0;
elseif C8>2&C8<3;
r83=(C8-2);
else C8<=2;
r83=0;
end;
if C8>=3;
r84=1;
else C8<3;
r84=0;
end;
R8=[r81 r82 r83 r84];
R=[R1;R2;R3;R4;R5;R6;R7;R8];
end
R
%%%%%%%%%%%%%%%%%%%%%%%%%%
D 佛洛依德 求最短路问题
clc,clear;
A=xlsread('first.xlsx',1);
B=xlsread('first.xlsx',2,'B1:C92');
n=size(A);%计算矩阵的行列数
L=[];
for i=1:n(1)
a=A(i);%左侧结点
b=A(i+140);%右侧结点
ax=B(a)%左侧点的横坐标
ay=B(a+92)%左侧点的纵坐标
bx=B(b)%右侧点的横坐标
by=B(b+92)%右侧点的纵坐标
L(a,b)=sqrt((bx-ax)^2+(by-ay)^2);%计算两点之间的距离L
L(92,:)=0;%因为剔除一点之后对应A中没有92行对应的数据
end
Lout=L+L';%对L进行转置相加
for x1=1:92 %找出输出矩阵中的0，用inf将其代替
for x2=1:92
if Lout(x1,x2)==0
Lout(x1,x2)=inf;
end
end
end
for x1=1:92%令输出矩阵中的对角线上的元素全为0
Lout(x1,x1)=0;
end
L
Lout;
size(L);
%佛洛依德算法
n=92;
A=Lout;% MATLAB 中, Inf 表示∞
D=A; %赋初值
for(i=1:n)
for(j=1:n)
R(i,j)=j;
end;
end %赋路径初值
for(k=1:n)
for(i=1:n)
for(j=1:n)
if(D(i,k)+D(k,j)
R(i,j)=k;
end;
end;
end %更新 rij
k; %显示迭代步数
D %显示每步迭代后的路长
R; %显示每步迭代后的路径
pd=0;
for i=1:n %含有负权时
if(D(i,i)<0)pd=1;
break;
end;
end %存在一条含有顶点 vi 的负回路
if(pd)break;
end %存在一条负回路, 终止程序
end %程序结束
D
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
E 0_1整数规划模型（模型还存在一些问题，还需要调试）
%0_1规划
function [xmin,f]=L01n(c,A,b,N,pre)%求0_1整数规划
if nargin<5%% nargin
pre=0;
if nargin<4
N=0;
end
end
c=c(:);%将矩阵的所有元素按列首尾相接依次排列
b=b(:);
[m,n]=size(A);%计算矩阵A的阶数
f=sum(abs(c));%对矩阵c中的元素的绝对值进行求和
x=zeros(1,n);%生成1行n列全0矩阵，与矩阵A的列数相同
ft=0;%给ft赋初值0
while 1
jj=0;
ft=dot(c,x);%dot(c,x)=c'*x 得到ft为1行n列全0矩阵 %%计算结果
为0
t1=ft-f;%%c中所有元素绝对值得和的负数
while(t1<=0)&(jj
t1=A(jj,:).*x-b(jj);%t1为矩阵A的第jj行所有列元素乘以矩阵x减去矩阵b的第jj个元素
%应该得到1行n列矩阵元素全部为 -b（jj）
if jj>N%N若为0 ，则应该肯定
if t1>0%不可能
jj=0;
end
else
if abs(t1)>pre%pre若为0，则应该肯定成立
jj=0;break;
end
end
end
end
if jj==m%当jj的值等于A的行数时
f=ft;%把ft的值0赋给f
xmin=x;%x为1行全0矩阵
jj=0;
end
k=1;
while x(k)==1
x(k)=0;
if k==n
return;
end
k=k+1;
end
x(k)=1;
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
F 计算最近围堵点
clc,clear;
A=xlsread('围堵.xlsx','A2:VJ81');
D=xlsread('围堵.xlsx','A84:VJ84');
n=size(A);
for i=1:n(1)
for j=1:n(2)
a(1,j)=min(A(:,j));
if A(i,j)==min(A(:,j))
b(i,j)=min(A(:,j));
elseif A(i,j)~=min(A(:,j))
b(i,j)=0;
end
end
end
b;
xlswrite('围堵.xlsx',b,2)
a;
xlswrite('围堵.xlsx',a,3)
N=D-30;
M=a-(D-30);
B=M';
xlswrite('围堵.xlsx',B,4)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
典型相关性分析
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%多元统计分析
clc,clear;
%load data;
%disp ('对数据进行标准化处理,并计算样本的协方差矩阵')
X=xlsread('G:\学习\数学建模\重要资料\典型相关分析\附件2-指标总表.xls','J2:AJ28');
Y=xlsread('G:\学习\数学建模\重要资料\典型相关分析\附件2-指标总表.xls','A2:I28');
%X=[x1 x2 x3 x4];
%Y=[y1 y2 y3 y4];
N=size(X);
n1=N(1);%X矩阵的行数
n2=N(2);%X矩阵的列数
M=size(Y);
m1=M(1);%Y矩阵的行数
m2=M(2);%Y矩阵的列数
All=[zscore(X) zscore(Y)];%zscore公式(X-MEAN(X)) ./ STD(X)
%Q=size(All);
%q1=Q(1);%All矩阵的行数
%q2=Q(2);%All矩阵的列数
DXY=cov(All) %协方差矩阵
%2）、计算矩阵A和B：
%disp ('计算矩阵A和B')
V11=DXY(1:n1,1:n2);
V12=DXY(1:n1,n2+1:n2+m2);
V21=DXY(n1+1:n1+m1,1:n2);
V22=DXY(n1+1:n1+m1,n2+1:n2+m2);
A=inv(V11)*V12*inv(V22)*V21;%
B=inv(V22)*V21*inv(V11)*V12;%
%3）、计算矩阵A和B的特征值和特征向量，进行排序并且判断其是否相等：
%disp ('矩阵A的特征向量和特征值')
[U_A,D_A]=eig(A);%
%disp ('矩阵B的特征向量和特征值')
[U_B,D_B]=eig(B);%
%对特征值进行排序：
[mA1,nA1]=sort(D_A,1);% 数字1表示从上到下升序排列
[mA2,nA2]=sort(mA1,2);% 数字2表示从上到下降序排列
[mB1,nB1]=sort(D_B,1);
[mB2,nB2]=sort(mB1,2);
%[Y,I]=sort（X,a）
%返回矩阵Y是对每列元素进行升序或降序排列
%返回矩阵I是得出Y中各个元素在X中属于哪一行
%判断A与B的特征值是否相等：
if mA2~=mB2
disp
('错误：A与B的特征值不相等')
else
disp ('经验证，A与B的特征值相等')
end
%4）、进行规格化处理，并求u和v：
%规格化处理,并求u和v：
for i=1:n2
a=(U_A(:,nA2(n2,n2-i+1)))'*V11*U_A(:,nA2(n2,n2-i+1));
u(:,i)=U_A(:,nA2(n2,n2-i+1))./sqrtm(a);
v(:,i)=inv(V22)*V21*u(:,i)./sqrtm(mA2(n2,n2-i+1));
end
%5）、结果输出：
%disp ('典型变量系数及相应典型相关系数依次为：')
for i=1:n2
disp ('-----------分割线-----------')
disp (i)
disp ('u=')
disp (u(:,i))
disp ('v=')
disp (v(:,i))
disp ('典型相关系数为：')
disp (mA2(4,4-i+1))
end
disp ('原变量在典型变量上的负荷矩阵为')
inv(u)
inv(v)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
夹角余弦法
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%多目标规划中用夹角余弦法确定权重
clc,clear;
A=[5.20,10.08,5.25,9.72,6.60;5.20,6.70,4.20,5.25,3.75;
4.73,5.71,3.82,5.54,3.30;0.473,1.599,0.473,1.313,0.803];
a1=A(1,:);
a2=A(2,:);
a3=A(3,:);
a4=A(4,:);
u11=min(a1);
u12=max(a2);
u13=max(a3);
u14=min(a4);
u1=[u11;u12;u13;u14];
u21=max(a1);
u22=min(a2);
u23=min(a3);
u24=max(a4);
u2=[u21;u22;u23;u24];
r1=abs(u11-a1)./(max(a1)-min(a1));
r2=abs(u12-a2)./(max(a2)-min(a2));
r3=abs(u13-a3)./(max(a3)-min(a3));
r4=abs(u14-a4)./(max(a4)-min(a4));
r=[r1;r2;r3;r4]
derta1=abs(u21-a1)./(max(a1)-min(a1));
derta2=abs(u22-a2)./(max(a2)-min(a2));
derta3=abs(u23-a3)./(max(a3)-min(a3));
derta4=abs(u24-a4)./(max(a4)-min(a4));
derta=[derta1;derta2;derta3;derta4]
dj=r.*derta;
Cs=sum(dj')';
Cx1=sqrt(sum((r.^2)')');
Cx2=sqrt(sum((derta.^2)')');
c=Cs./(Cx1.*Cx2)
cq=sum((r.*derta)')'./(sqrt((sum((r.^2)'))').*sqrt((sum((derta.^2)'))'))
w=c./sum(c)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
熵权法
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%商权法确定权重
%确定指标层权重
clear;
x=[6.22 282.43 53.31 47.52 91.58 18.01 65.40 218.95
7.18 373.88 46.78 58.94 274.99 19.19 84.80 245.36
4.04 152.32 38.96 17.32 40.96 15.45 36.56 73.29
5.69 346.81 51.65 52.07 125.71 16.63 61.54 224.99
6.26 280.54 43.64 30.19 114.99 15.29 60.71 154.24 ];
x = x';
y=[];
[m,n]=size(x);
for i=1:n
y(:,i)=x(:,i)/sum([x(:,i)]); %原始矩阵归一化
end
for l=1:n
s(1,l)=0;
for j=1:m
p(1,l)=y(j,l)*log(y(j,l));
s(1,l)=s(1,l)+p(1,l);
end
end
k=(log(m))^(-1);
e=-k*s;
h=ones(1,n)-e;
w=h/sum(h);%指标权重值
sum(w);
%计算综合评价值
g=y*w'
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
BP 神经网络模型
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
close all ;
clear ;
echo on ;
clc ;
% NEWFF——生成一个新的前向神经网络
% TRAIN——对 BP 神经网络进行训练
% SIM——对 BP 神经网络进行仿真
pause
% 敲任意键开始
clc
% 定义训练样本
% P 为输入矢量
p=
[ ];
% T 为目标矢量
t= [];
%训练样本的归一化
for i=1:（训练样本的指标数）
P(i,:)=(p(i,:)-min(p(i,:)))/(max(p(i,:))-min(p(i,:)));
end
pause
clc
% 创建一个新的前向神经网络
net=newff(minmax(P),[15,2],{'tansig','purelin'},'traingda'); %这些参数要自己设置
% 当前输入层权值和阈值
inputWeights=net.IW{1,1};
inputbias=net.b{1} ;
% 当前网络层权值和阈值
layerWeights=net.LW{2,1} ;
layerbias=net.b{2} ;
pause
clc
% 设置训练参数 也要自己设置
net.trainParam.show = 50;
net.trainParam.lr = 0.05;
%net.trainParam.mc = 0.9; % 附加动量因子
net.trainParam.epochs =5000;
net.trainParam.goal = 1e-4;
pause
clc
% 调用 TRAINGDM 算法训练 BP 网络
[net,tr]=train(net,P,t);
pause
clc
% 对 BP 网络进行仿真
p_test=[];
for i=1:6
P_test(i,:)=(p_test(i,:)-min(p_test(i,:)))/(max(p_test(i,:))-min(p_test(i,:)));
end
t_test=[];
A = sim(net,P_test)
% 计算仿真误差
E = t - A
error=mse(E)
pause
clc
echo off
%反归一化
for i=1:2
predict(i,:)=A(i,:)*(max(t(i,:))-min(t(i,:)))+ min(t(i,:));
end
predict %即仿真结果
pause
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
层次分析法
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%矩阵归一化处理
clc，clear；
%若已知矩阵为A，则[V,D]=eig(A)，
%其中对角阵D的对角元素为A的特征值，
%而矩阵D的每一列向量为其所对应的特征向量.
clc,clear;
A=input('请输入待处理矩阵：A=');
[V, D]=eig(A);
disp('最大特征值 e')
e=max(eig(A))%最大特征值
disp('最大特征值e对应的特征向量 v')
v=null(A-e*eye(length(A)))% e 对应特征向量
disp('对特征向量进行归一化处理')
vv=v./sum(v)%归一化处理
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
！！！！娱乐！！！
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 中华人民共和国国歌
clc,clear all
fs = 44100; dt = 1/fs; f0 = 320; T16 = 0.125;
t16 = [0:dt:T16]; [temp m] = size(t16);
t2 = linspace(0,8*T16,8*m); [temp i] = size(t2);
t4 = linspace(0,4*T16,4*m); [temp j] = size(t4);
t8 = linspace(0,2*T16,2*m); [temp k] = size(t8);
t12 = linspace(0,4/3*T16,4/3*m); [temp l] = size(t12);
% Modification functions
mod2 = sin(pi*t2/t2(end));
mod4 = sin(pi*t4/t4(end));
mod8 = sin(pi*t8/t8(end));
mod12 = sin(pi*t12/t12(end));
mod16 = sin(pi*t16/t16(end));
ScaleTable = [2/3 3/4 5/6 15/16 ...
1 9/8 5/4 4/3 3/2 5/3 9/5 15/8 ...
2 9/4 5/2 8/3 3 10/3 15/4 4 ...
1/2 9/16 5/8];
do0g = mod2.*cos(2*pi*ScaleTable(21)*f0*t2);re0g = mod2.*cos(2*pi*ScaleTable(22)*f0*t2);mi0g = mod2.*cos(2*pi*ScaleTable(23)*f0*t2);fa0g = mod2.*cos(2*pi*ScaleTable(1)*f0*t2);so0g = mod2.*cos(2*pi*ScaleTable(2)*f0*t2);la0g = mod2.*cos(2*pi*ScaleTable(3)*f0*t2);ti0g = mod2.*cos(2*pi*ScaleTable(4)*f0*t2);
do1g =
mod2.*cos(2*pi*ScaleTable(5)*f0*t2);re1g = mod2.*cos(2*pi*ScaleTable(6)*f0*t2);mi1g = mod2.*cos(2*pi*ScaleTable(7)*f0*t2);fa1g = mod2.*cos(2*pi*ScaleTable(8)*f0*t2);so1g = mod2.*cos(2*pi*ScaleTable(9)*f0*t2);la1g = mod2.*cos(2*pi*ScaleTable(10)*f0*t2);% tb1g = mod2.*cos(2*pi*ScaleTable(11)*f0*t2);ti1g = mod2.*cos(2*pi*ScaleTable(12)*f0*t2);
do2g = mod2.*cos(2*pi*ScaleTable(13)*f0*t2);re2g = mod2.*cos(2*pi*ScaleTable(14)*f0*t2);mi2g = mod2.*cos(2*pi*ScaleTable(15)*f0*t2);fa2g = mod2.*cos(2*pi*ScaleTable(16)*f0*t2);so2g = mod2.*cos(2*pi*ScaleTable(17)*f0*t2);la2g = mod2.*cos(2*pi*ScaleTable(18)*f0*t2);ti2g = mod2.*cos(2*pi*ScaleTable(19)*f0*t2);
do3g = mod2.*cos(2*pi*ScaleTable(20)*f0*t2);
I = zeros(1,i);
do0f = mod4.*cos(2*pi*ScaleTable(21)*f0*t4);re0f = mod4.*cos(2*pi*ScaleTable(22)*f0*t4);mi0f = mod4.*cos(2*pi*ScaleTable(23)*f0*t4);fa0f = mod4.*cos(2*pi*ScaleTable(1)*f0*t4);so0f = mod4.*cos(2*pi*ScaleTable(2)*f0*t4);la0f = mod4.*cos(2*pi*ScaleTable(3)*f0*t4);ti0f = mod4.*cos(2*pi*ScaleTable(4)*f0*t4);
do1f = mod4.*cos(2*pi*ScaleTable(5)*f0*t4);re1f = mod4.*cos(2*pi*ScaleTable(6)*f0*t4);mi1f = mod4.*cos(2*pi*ScaleTable(7)*f0*t4);fa1f = mod4.*cos(2*pi*ScaleTable(8)*f0*t4);so1f = mod4.*cos(2*pi*ScaleTable(9)*f0*t4);la1f = mod4.*cos(2*pi*ScaleTable(10)*f0*t4);% tb1f = mod4.*cos(2*pi*ScaleTable(11)*f0*t4);ti1f = mod4.*cos(2*pi*ScaleTable(12)*f0*t4);
do2f = mod4.*cos(2*pi*ScaleTable(13)*f0*t4);re2f = mod4.*cos(2*pi*ScaleTable(14)*f0*t4);mi2f = mod4.*cos(2*pi*ScaleTable(15)*f0*t4);fa2f = mod4.*cos(2*pi*ScaleTable(16)*f0*t4);so2f = mod4.*cos(2*pi*ScaleTable(17)*f0*t4);la2f = mod4.*cos(2*pi*ScaleTable(18)*f0*t4);ti2f = mod4.*cos(2*pi*ScaleTable(19)*f0*t4);
do3f = mod4.*cos(2*pi*ScaleTable(20)*f0*t4);
J = zeros(1,j);
do0e = mod8.*cos(2*pi*ScaleTable(21)*f0*t8);re0e = mod8.*cos(2*pi*ScaleTable(22)*f0*t8);mi0e = mod8.*cos(2*pi*ScaleTable(23)*f0*t8);fa0e = mod8.*cos(2*pi*ScaleTable(1)*f0*t8);so0e = mod8.*cos(2*pi*ScaleTable(2)*f0*t8);la0e = mod8.*cos(2*pi*ScaleTable(3)*f0*t8);ti0e = mod8.*cos(2*pi*ScaleTable(4)*f0*t8);
do1e = mod8.*cos(2*pi*ScaleTable(5)*f0*t8);re1e = mod8.*cos(2*pi*ScaleTable(6)*f0*t8);mi1e = mod8.*cos(2*pi*ScaleTable(7)*f0*t8);fa1e = mod8.*cos(2*pi*ScaleTable(8)*f0*t8);so1e = mod8.*cos(2*pi*ScaleTable(9)*f0*t8);la1e = mod8.*cos(2*pi*ScaleTable(10)*f0*t8);% tb1e = mod8.*cos(2*pi*ScaleTable(11)*f0*t8);ti1e = mod8.*cos(2*pi*ScaleTable(12)*f0*t8);
do2e = mod8.*cos(2*pi*ScaleTable(13)*f0*t8);re2e = mod8.*cos(2*pi*ScaleTable(14)*f0*t8);mi2e = mod8.*cos(2*pi*ScaleTable(15)*f0*t8);fa2e = mod8.*cos(2*pi*ScaleTable(16)*f0*t8);so2e = mod8.*cos(2*pi*ScaleTable(17)*f0*t8);la2e = mod8.*cos(2*pi*ScaleTable(18)*f0*t8);ti2e = mod8.*cos(2*pi*ScaleTable(19)*f0*t8);
do3e = mod8.*cos(2*pi*ScaleTable(20)*f0*t8);
K = zeros(1,k);
do0d = mod12.*cos(2*pi*ScaleTable(21)*f0*t12);re0d = mod12.*cos(2*pi*ScaleTable(22)*f0*t12);mi0d = mod12.*cos(2*pi*ScaleTable(23)*f0*t12);fa0d = mod12.*cos(2*pi*ScaleTable(1)*f0*t12);so0d = mod12.*cos(2*pi*ScaleTable(2)
*f0*t12);la0d = mod12.*cos(2*pi*ScaleTable(3)*f0*t12);ti0d = mod12.*cos(2*pi*ScaleTable(4)*f0*t12);
do1d = mod12.*cos(2*pi*ScaleTable(5)*f0*t12);re1d = mod12.*cos(2*pi*ScaleTable(6)*f0*t12);mi1d = mod12.*cos(2*pi*ScaleTable(7)*f0*t12);fa1d = mod12.*cos(2*pi*ScaleTable(8)*f0*t12);so1d = mod12.*cos(2*pi*ScaleTable(9)*f0*t12);la1d = mod12.*cos(2*pi*ScaleTable(10)*f0*t12);% tb1d = mod12.*cos(2*pi*ScaleTable(11)*f0*t12);ti1d = mod12.*cos(2*pi*ScaleTable(12)*f0*t12);
do2d = mod12.*cos(2*pi*ScaleTable(13)*f0*t12);re2d = mod12.*cos(2*pi*ScaleTable(14)*f0*t12);mi2d = mod12.*cos(2*pi*ScaleTable(15)*f0*t12);fa2d = mod12.*cos(2*pi*ScaleTable(16)*f0*t12);so2d = mod12.*cos(2*pi*ScaleTable(17)*f0*t12);la2d = mod12.*cos(2*pi*ScaleTable(18)*f0*t12);ti2d = mod12.*cos(2*pi*ScaleTable(19)*f0*t12);
do3d = mod12.*cos(2*pi*ScaleTable(20)*f0*t12);
L = zeros(1,l);
do0s = mod16.*cos(2*pi*ScaleTable(21)*f0*t16);re0s = mod16.*cos(2*pi*ScaleTable(22)*f0*t16);mi0s = mod16.*cos(2*pi*ScaleTable(23)*f0*t16);fa0s = mod16.*cos(2*pi*ScaleTable(1)*f0*t16);so0s = mod16.*cos(2*pi*ScaleTable(2)*f0*t16);la0s = mod16.*cos(2*pi*ScaleTable(3)*f0*t16);ti0s = mod16.*cos(2*pi*ScaleTable(4)*f0*t16);
do1s = mod16.*cos(2*pi*ScaleTable(5)*f0*t16);re1s = mod16.*cos(2*pi*ScaleTable(6)*f0*t16);mi1s = mod16.*cos(2*pi*ScaleTable(7)*f0*t16);fa1s = mod16.*cos(2*pi*ScaleTable(8)*f0*t16);so1s = mod16.*cos(2*pi*ScaleTable(9)*f0*t16);la1s = mod16.*cos(2*pi*ScaleTable(10)*f0*t16);% tb1s = mod16.*cos(2*pi*ScaleTable(11)*f0*t16);ti1s = mod16.*cos(2*pi*ScaleTable(12)*f0*t16);
do2s = mod16.*cos(2*pi*ScaleTable(13)*f0*t16);re2s = mod16.*cos(2*pi*ScaleTable(14)*f0*t16);mi2s = mod16.*cos(2*pi*ScaleTable(15)*f0*t16);fa2s = mod16.*cos(2*pi*ScaleTable(16)*f0*t16);so2s = mod16.*cos(2*pi*ScaleTable(17)*f0*t16);la2s = mod16.*cos(2*pi*ScaleTable(18)*f0*t16);ti2s = mod16.*cos(2*pi*ScaleTable(19)*f0*t16);
do3s = mod16.*cos(2*pi*ScaleTable(20)*f0*t16);
M = zeros(1,m);
v1 = [do0e L mi0e so0e so0e M,la0f M so0f M,mi0e M do0e so0d so0d so0d M,mi0f do0f L,so0d so0d so0d so0d so0d so0d,do0f K M so0e,do1f J do1e,do1e M do1e so0d M la0d ti0d M,do1f do1f,J mi1e M do1d re1d mi1d M,so1f L so1f M,mi1e M mi1e do1e M mi1e,so1e M mi1e re1f M,re1f J M la1f M so1f M,re1f M mi1f M,so1e M mi1e K so1e,M mi1e re1d mi1d do1e K,mi1f J M so1e M la1e do1e do1e M,mi1e M mi1e so1e so1e M,re1e re1e re1e la1f ,re1f J so1e,do1f J do1e,mi1f J mi1e,so1f J do1e M mi1e so1e M so1e M,la1f M so1f,L mi1e M do1e so1e so1e so1e,mi1e K do1e J,so1f do1f K mi1e L do1e so1e so1e so1e,mi1e J do1e J,so1f do1f,K so1f do1f,K so1f do1f,do1f];
s = v1;s = s/max(s);sound(s,fs);clc
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%玲珑剔透球
%玲珑剔透球
[X0,Y0,Z0]=sphere(30); %产生单位球面的三维坐标
X=2*X0;Y=2*Y0;Z=2*Z0; %产生半径为2的球面的三维坐标
surf(X0,Y0,Z0); %画单位球面
shading interp %采用插补明暗处理
hold on;
mesh(X,Y,Z);
hold off %画外球面
hidden off %产生透视效果
axis off %不显示坐标轴
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
以上即为本人近段时间准备数学建模过程中总结的一些程序，还望对参加数学建模的同学有一定的帮助！！！
祝你成功！！
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%