三．举例
例1：求解线性规划问题：
max f=2x1+5x2
s.t
先将目标函数转化成最小值问题：min(-f)=- 2x1-5x2程序：
f=[-2 -5];
A=[1 0;0 1;1 2];
b=[4;3;8];
[x,fval]=linprog(f,A,b)
f=fval*(-1)
结果：x = 2
3
fval = -19.0000
maxf = 19
例2：minf=5x1-x2+2x3+3x4-8x5
s.t –2x1+x2-x3+x4-3x5≤6
2x1+x2-x3+4x4+x5≤7
0≤x j≤15 j=1,2,3,4,5
程序：
f=[5 -1 2 3 -8];
A=[-2 1 -1 1 -3;2 1 -1 4 1];
b=[6;7];
lb=[0 0 0 0 0];
ub=[15 15 15 15 15];
[x,fval]=linprog(f,A,b,[],[],lb,ub)
结果：x =
0.0000
0.0000
8.0000
0.0000
15.0000
minf =
-104
例3：求解线性规划问题：
minf=5x1+x2+2x3+3x4+x5
s.t –2x1+x2-x3+x4-3x5≤1
2x1+3x2-x3+2x4+x5≤-2
0≤x j≤1 j=1,2,3,4,5
程序：
f=[5 1 2 3 1];
A=[-2 1 -1 1 -3;2 3 -1 2 1];
b=[1;-2];
lb=[0 0 0 0 0];
ub=[1 1 1 1 1];
[x,fval,exitflag,output,lambda]=linprog(f,A,b,[],[],lb,ub) 运行结果：Exiting: One or more of the residuals, duality gap, or total relative error
has grown 100000 times greater than its minimum value so far: the primal appears to be infeasible (and the dual unbounded).
(The dual residual < TolFun=1.00e-008.)
x = 0.0000
0.0000
1.1987
0.0000
0.0000
fval =
2.3975
exitflag =
-1
output =
iterations: 7
cgiterations: 0
algorithm: 'lipsol'
lambda =
ineqlin: [2x1 double]
eqlin: [0x1 double]
upper: [5x1 double]
lower: [5x1 double]
显示的信息表明该问题无可行解。所给出的是对约束破坏最小的解。
例4：求解实例1的生产计划问题
建立数学模型：
设x1、x2分别为生产甲、乙产品的件数。f为该厂所获总润。
max f=70x1+120x2
s.t 9x1+4x2≤3600
4x1+5x2≤2000
3x1+10x2≤3000
x1,x2≥0
将其转换为标准形式：
min f=-70x1-120x2
s.t 9x1+4x2≤3600
4x1+5x2≤2000
3x1+10x2≤3000
x1,x2≥0
程序：f=[-70 -120];
A=[9 4 ;4 5;3 10 ];
b=[3600;2000;3000];
lb=[0 0];
ub=[];
[x,fval,exitflag]=linprog(f,A,b,[],[],lb,ub) maxf=-fval
结果：x =
200.0000
240.0000
fval =
-4.2800e+004
exitflag =
1
maxf =
4.2800e+004
例5：求解实例2
建立数学模型：
max f=0.15x1+0.1x2+0.08 x3+0.12 x4
s.t x1-x2- x3- x4≤0
x2+ x3- x4≥0
x1+x2+x3+ x4=1
x j≥0 j=1,2,3,4
将其转换为标准形式：
min z=-0.15x1-0.1x2-0.08 x3-0.12 x4
s.t x1-x2- x3- x4≤0
-x2- x3+ x4≤0
x1+x2+x3+ x4=1
x j≥0 j=1,2,3,4
程序： f = [-0.15;-0.1;-0.08;-0.12];
A = [1 -1 -1 -1
0 -1 -1 1];
b = [0; 0];
Aeq=[1 1 1 1];
beq=[1];
lb = zeros(4,1);
[x,fval,exitflag] = linprog(f,A,b,Aeq,beq,lb)
f=-fval
结果：x =
0.5000
0.2500
0.0000
0.2500
fval =
-0.1300
exitflag =
1
f =
0.1300
即4个项目的投资百分数分别为50%，25%，0, 25%时可使该公司获得最大的收益，其最大收益可到达13%。过程正常收敛。
例6：求解实例3
建立数学模型：
设a i j为由工厂i运到市场j的费用，x i j 是由工厂i运到市场j的箱数。b i是工厂i的产量，d j是市场j的需求量。
b= ( 60 40 50 )T d= ( 20 35 33 34 )T
s.t
x i j≥0
程序：A=[2 1 3 2;1 3 2 1;3 4 1 1];
f=A(:);
B=[ 1 0 0 1 0 0 1 0 0 1 0 0
0 1 0 0 1 0 0 1 0 0 1 0
0 0 1 0 0 1 0 0 1 0 0 1];
D=[1 1 1 0 0 0 0 0 0 0 0 0
0 0 0 1 1 1 0 0 0 0 0 0
0 0 0 0 0 0 1 1 1 0 0 0
0 0 0 0 0 0 0 0 0 1 1 1];
b=[60;40;50];
d=[20;35;33;34];
lb=zeros(12,1);
[x,fval,exitflag]=linprog(f,B,b,D,d,lb) 结果：x =
0.0000
20.0000
0.0000
35.0000
0.0000
0.0000
0.0000
0.0000
33.0000
0.0000
18.4682
15.5318
fval =
122.0000
exitflag =
1
即运输方案为：甲市场的货由B厂送20箱；乙市场的货由A厂送35箱；丙商场的货由C厂送33箱；丁市场的货由B厂送18箱，再由C厂送16箱。
最低总运费为：122元。
§2 非线性规划模型
一．非线性规划课题
实例1 表面积为36平方米的最大长方体体积。
建立数学模型：
设x、y、z分别为长方体的三个棱长，f为长方体体积。
max f = x y (36-2 x y)/2 (x+y)
实例2 投资决策问题
某公司准备用5000万元用于A、B两个项目的投资，设x1、x2分别表示配给项目A、B的投资。预计项目A、B的年收益分别为20%和16%。同时，投资后总的风险损失将随着总投资和单位投资的增加而增加，已知总的风险损失为2x12+x22+(x1+x2)2.问应如何分配资金，才能使期望的收益最大，同时使风险损失为最小。
建立数学模型：
max f=20x1+16x2-λ[2x12+x22+(x1+x2)2]
s.t x1+x2≤5000
x1≥0,x2≥0
目标函数中的λ≥0是权重系数。
由以上实例去掉实际背景，其目标函数与约束条件至少有一处是非线性的，称其为非线性问题。
非线性规划问题可分为无约束问题和有约束问题。实例1为无约束问题，实例2为有约束问题。二．无约束非线性规划问题：
求解无约束最优化问题的方法主要有两类：直接搜索法(Search method)和梯度法(Gradient method).
1．fminunc函数
调用格式：x=fminunc(fun,x0)
x=fminunc(fun,x0,options)
x=fminunc(fun,x0,options,P1,P2)
[x,fval]=fminunc(…)
[x,fval, exitflag]=fminunc(…)
[x,fval, exitflag,output]=fminunc(…)
[x,fval, exitflag,output,grad]=fminunc(…)
[x,fval, exitflag,output,grad,hessian]=fminunc(…)
说明：fun为需最小化的目标函数，x0为给定的搜索的初始点。options指定优化参数。
返回的x为最优解向量；fval为x处的目标函数值；exitflag描述函数的输出条件；output返回优化信息；grad返回目标函数在x处的梯度。Hessian返回在x处目标函数的Hessian矩阵信息。
例1 ：求
程序：编辑ff1.m文件
function f=ff1(x)
f=8*x(1)-4*x(2) +x(1)^2+3*x(2)^2;
通过绘图确定一个初始点：
[x,y]=meshgrid(-10:.5:10);
z= 8*x-4*y +x.^2+3*y.^2;
surf(x,y,z)
选初始点：x0=(0,0)
x0=[0,0];
[x,fval,exitflag]=fminunc(@ff1,x0) 结果：x =
-4.0000 0.6667
fval =
-17.3333
exitflag =
1
例2：
程序：编辑ff2.m文件：
function f=ff2(x)
f=4*x(1)^2+5*x(1)*x(2)+2*x(2)^2; 取初始点：x0=(1,1)
x0=[1,1];
[x,fval,exitflag]=fminunc(@ff2,x0)
结果： x =
1.0e-007 *
-0.1721 0.1896 fval =
2.7239e-016 exitflag =
1
例3：将上例用提供的梯度g最小化函数进行优化计算。
修改M文件为：
function [f,g]=ff3(x)
f=4*x(1)^2+5*x(1)*x(2)+2*x(2)^2;
if nargut >1
g(1)=8*x(1)+5*x(2);
g(2)=5*x(1)+4*x(2);
end
通过下面将优化选项结构options.GradObj设置为’on’来得到梯度值。
options=optimset(‘Gradobj’,’on’);
x0=[1,1];
[x,fval,exitflag]=fminunc(@ff3,x0,options)
结果：x =
1.0e-015 *
-0.2220 -0.2220
fval =
5.4234e-031
exitflag =
1
2．minsearch函数
调用格式：x=fminsearch(fun,x0)
x=fminsearch(fun,x0,options)
x=fminsearch(fun,x0,options,P1,P2)
[x,fval]=fminsearch(…)
[x,fval, exitflag]=fminsearch(…)
[x,fval, exitflag,output]=fminsearch(…)
[x,fval, exitflag,out put,grad]=fminsearch(…)
[x,fval, exitflag,output,grad,hessian]=fminsearch(…)
说明：参数及返回变量同上一函数。对求解二次以上的问题，fminsearch函数比fminunc函数有效。
3．多元非线性最小二乘问题：
非线线性最小二乘问题的数学模型为：
其中L为常数。
调用格式： x=lsqnonlin(fun,x0)
x=lsqnonlin(fun,x0,lb,ub)
x=lsqnonlin(fun,x0,options)
x=lsqnonlin(fun,x0,options,P1,P2)
[x,resnorm]=lsqnonlin(…)
[x,resnorm, residual,exitflag]=lsqnonlin(…)
[x,resnorm, residual , exitflag,output]=lsqnonlin(…)
[x,resnorm, residual,exitflag, output,lambda]=lsqnonlin(…)
[x,resnorm, r esidual,exitflag, output,lambda,jacobian]=lsqnonlin(…)
说明：x返回解向量；resnorm返回x处残差的平方范数值：sum(fun(x).^2)；residual返回x处的残差值fun(x)；lambda返回包含x处拉格朗日乘子的结构参数；jacobian返回解x处的fun函数的雅可比矩阵。
lsqnonlin默认时选择大型优化算法。Lsqnonlin通过将https://www.doczj.com/doc/e412550215.html,rgeScale设置为’off’来作中型优化算法。其采用一维搜索法。
例4．求minf=4(x2-x1)2+(x2-4)2，选择初始点x0(1,1)
程序：
f ='4*(x(2)-x(1))^2+(x(2)-4)^2'
[x,reshorm]=lsqnonlin(f,[1,1])
结果：x =
3.9896 3.9912
reshorm =
5.0037e-009
例5：求，选择初始点x0(0.2,0.3)
求解：先编辑ff5.m文件：
function f=ff5(x)
k=1:10;
f=2+2*k-exp(k*x(1))-exp(k*x(2)); 然后作程序：x0=[0.2,0.3];
[x,resnorm]=lsqnonlin(@ff5,x0) 结果：x =
0.2578 0.2578
resnorm =
124.3622
二．有约束非线性规划问题：
数学模型： min F(x)
s.t G i (x) ≤0i=1,…,m
G j (x) =0 j=m+1,…,n
x l≤x≤x u
其中：F(x)为多元实值函数，G(x)为向量值函数，
在有约束非线性规划问题中，通常要将该问题转换为更简单的子问题，这些子问题可以求并作为迭代过程的基础。其基于K-T方程解的方法。它的K-T方程可表达为：
方程第一行描述了目标函数和约束条件在解处梯度的取消。由于梯度取消，需要用拉格朗日乘子λi 来平衡目标函数与约束梯度间大小的差异。
调用格式: x=fmincon(f,x0,A,b)
x=fmincon(f,x0,A,b,Aeq,beq)
x=fmincon(f,x0,A,b,Aeq,beq,lb,ub)
x=fmincon(f,x0,A,b,Aeq,beq,lb,ub,nonlcon)
x=fmincon(f,x0,A,b,Aeq,beq,lb,ub,nonlcon,options)
[x,fval]=fmincon(…)
[x, fval, exitflag]=fmincon(…)
[x, fval, exitflag, output]=fmincon(…)
[x, fval, exitflag, output, lambda]=fmincon(…)
说明：x=fmincon(f,x0,A,b)返回值x为最优解向量。其中：x0为初始点。A,b为不等式约束的系数矩阵和右端列向量。
x=fmincon(f,x0,A,b,Aeq,beq) 作有等式约束的问题。若没有不等式约束，则令A=[ ]、b=[ ] 。
x=fmincon(f, x0,A,b,Aeq,beq,lb,ub, nonlcon ,options) 中lb ,ub为变量x的下界和上界；nonlcon=@fun,由M文件fun.m给定非线性不等式约束c (x) ≤0和等式约束g(x)=0；options为指定优化参数进行最小化。
例6：求解：min 100(x2-x12 )2+(1-x1)2
s.t x1≤2;
x2≤2
程序：首先建立ff6.m文件：
function f=ff6(x)
f=100*(x(2)-x(2)^2)^2+(1-x(1))^2;
然后在工作空间键入程序：
x0=[1.1,1.1];
A=[1 0;0 1];
b=[2;2];
[x,fval]=fmincon(@ff6,x0,A,b)
结果： x =
1.0000 1.0000
fval =
3.1936e-011
例7：求解：
首先建立目标函数文件ff7.m文件： function f=ff7(x)
f=-x(1)*x(2)*x(3)
然后将约束条件改写成如下不等式： -x1-2x2-2x3≤0
x1+2x2+2x3≤72
在工作空间键入程序：
A=[-1 –2 –2;1 2 2];
b=[0;72];
x0=[10;10;10];
[x,fval]=fmincon(@ff71,x0,A,b) 结果： x =
24.0000
12.0000
12.0000
fval =
-3456
例8求解：minf=e x1(6x12+3x22+2x1x2+4x2+1)
s.t x1x2-x1-x2+1≤0
-2x1x2-5≤0
程序：首先建立目标函数文件ff8.m文件：
function f=ff8(x)
f=exp(x(1))*(6*x(1)^2+3*x(2)^2+2*x(1)*x(2)+4*x(2)+1);
再建立非线性的约束条件文件：ff8g.m
function [c,g]=ff8g(x)
c(1)=x(1)*x(2)-x(1)-x(2)+1；
c(2)=-2*x(1)*x(2)-5；
g=[];
然后在工作空间键入程序：
x0=[1,1];
nonlcon=@ff8g
[x,fval] =fmincon(@ff8,x0,[],[],[],[],[],[],nonlcon)
结果： x =
-2.5000 1.0000
fval =
3.3244
exitflag =
1
当有等式约束时，要放在矩阵g的位置，如上例中加等式约束：
x(1)+2*x(1)=0
程序：首先建立fun1.m文件：
function[c,g]=ff8g1(x)
c(1)=x(1)*x(2)-x(1)-x(2)+1；
c(2)=-2*x(1)*x(2)-5；
g(1)=x(1)+2*x(2);
然后在工作空间键入程序：
x0=[-1,1];
nonlcon=@ff8g1;
[x,fval,exitflag] =fmincon(@ff8,x0,[],[],[],[],[],[], nonlcon)
结果： x =
-2.2361 1.1180
fval =
3.6576
exitflag =
1
§3 二次规划模型
数学模型：
其中H为二次型矩阵，A、Aeq分别为不等式约束与等式约束系数矩阵，f,b,beq,lb,ub,x为向量。
求解二次规划问题函数为quadprog( )
调用格式：X= quadprog(H,f,A,b)
X= quadprog(H,f,A,b,Aeq,beq)
X= quadprog(H,f,A,b,Aeq,beq,lb,ub)
X= quadprog(H,f,A,b,Aeq,beq,lb,ub,x0)
X= quadprog(H,f,A,b,Aeq,beq,lb,ub,x0,options)
[x,fval]= quadprog(…)
[x,fval,exitflag]= quadprog(…)
[x,fval,exitflag,output]= quadprog(…)
[x,fval,exitflag,output,lambda]= quadprog(…)
说明：输入参数中，x0为初始点；若无等式约束或无不等式约束，就将相应的矩阵和向量设置为空；opti ons为指定优化参数。输出参数中，x是返回最优解；fval是返回解所对应的目标函数值；exitflag是描述搜索是否收敛；output是返回包含优化信息的结构。Lambda是返回解x入包含拉格朗日乘子的参数。
例1：求解：二次规划问题
min f(x)= x1-3x2+3x12+4x22-2x1x2
s.t 2x1+x2≤2
-x1+4x2≤3
程序：f=[1;-3]
H=[6 -2;-2 8]
A=[2 1;-1 4]
b=[2;3]
[X,fval,exitflag]=quadprog(H,f,A,b)
结果：X =
-0.0455
0.3636
fval =
-0.5682
exitflag =
1
例2：求解：二次规划问题
min +x12+2x22-2x1x2-4x1-12x2
s.t x1+x2≤2
-x1+2x2≤2
2x1+x2≤3
0≤x1, 0≤x2
程序：H=[2 -2;-2 4];
f=[-4;-12];
A=[1 1;-1 2;2 1];
b=[2;2;3];
lb=zeros(2,1);
[x,fval,exitflag]=quadprog(H,f,A,b,[],[],lb)
结果： x =
0.6667
1.3333
fval =
-16.4444
exitflag =
1
§4 多目标规划模型
多目标规划定义为在一组约束下，多个不同的目标函数进行优化设计。数学模型：
s.t g j(x) ≤0j=1, 2, … ,k
其中x=(x1 ,x2 , … ,x n)为一个n维向量；f i(x)为目标函数，i=1, 2, … ,m; g j(x)为系统约束, j=1, 2, …,k。
当目标函数处于冲突状态时，不存在最优解使所有目标函数同时达到最优。于是我们寻求有效解(又称非劣解或非支配解或帕累托解)
定义：若(∈Ω)的邻域内不存在Δx，使得(+Δx∈Ω)，且
则称为有效解。
多目标规划问题的几种常用解法：
(1) 主要目标法
其基本思想是：在多目标问题中，根据问题的实际情况，确定一个目标为主要目标，而把其余目标作为次要目标，并且根据经验，选取一定的界限值。这样就可以把次要目标作为约束来处理，于是就将原来的多目标问题转化为一个在新的约束下的单目标最优化问题。
(2) 线性加权和法
其基本思想是：按照多目标f i(x) (i=1, 2, … ,m)的重要程度，分别乘以一组权系数λj(j=1, 2, … ,m)然后相加作为目标函数而构成单目标规划问题。即，其中
例1：某钢铁厂准备用5000万用于A、B两个项目的技术改造投资。设x1、x2分别表示分配给项目A、B 的投资。据专家预估计，投资项目A、B的年收益分别为70%和66%。同时，投资后总的风险损失将随着总投资和单项投资的增加而增加，已知总的风险损失为0.02x12+0.01x22+0.04(x1+x2)2，问应如何分配资金才能使期望的收益最大，同时使风险损失为最小。
建立数学模型
max f1(x)=70x1+66x2
min f2(x)= 0.02x12+0.01x22+0.04(x1+x2)2
s.t x1+x2≤5000
0≤x1, 0≤x2
线性加权构造目标函数： max f=0.5f1(x) –0.5f2(x)
化最小值问题：min (-f)=- 0.5f1(x) +0.5f2(x)
优化工具箱的使用 MATLAB的优化工具箱提供了各种优化函数，这些优化函数可以通过在命令行输入相应的函数名加以调用；此外为了使用方便，MA TLAB还提供了图形界面的优化工具（GUI Optimization tool）。 1 GUI优化工具 GUI优化工具的启动 有两种启动方法： （1）在命令行输入optimtool； （2）在MA TLAB主界面单击左下角的“Start”按钮，然后依次选择“Toolboxes→Optimization→Optimization tool” GUI优化工具的界面 界面分为三大块： 左边（Problem Setup and Results）为优化问题的描述及计算结果显示； 中间（Options）为优化选项的设置； 右边（Quick Reference）为帮助。为了界面的简洁，可以单击右上角“<<”、“>>”的按钮将帮助隐藏或显示。 1、优化问题的描述及计算结果显示 此板块主要包括选择求解器、目标函数描述、约束条件描述等部分。 选择合适的求解器以及恰当的优化算法，是进行优化问题求解的首要工作。 ?Solver：选择优化问题的种类，每类优化问题对应不同的求解函数。 ?Algorithm：选择算法，对于不同的求解函数，可用的算法也不同。 Problem框组用于描述优化问题，包括以下内容： ?Objective function: 输入目标函数。 ?Derivatives: 选择目标函数微分（或梯度）的计算方式。 ?Start point: 初始点。 Constraints框组用于描述约束条件，包括以下内容： ?Linear inequalities: 线性不等式约束，其中A为约束系数矩阵，b代表约束向量。 ?Linear equalities: 线性等式约束，其中Aeq为约束系数矩阵，beq代表约束向量。 ?Bounds: 自变量上下界约束。 ?Nonlinear Constraints function; 非线性约束函数。 ?Derivatives: 非线性约束函数的微分（或梯度）的计算方式。 Run solver and view results框组用于显示求解过程和结果。 （对于不同的优化问题类型，此板块可能会不同，这是因为各个求解函数需要的参数个数不一样，如Fminunc 函数就没有Constraints框组。） 2、优化选项（Options） ?Stopping criteria: 停止准则。
Matlab频谱分析程序
Matlab 信号处理工具箱 谱估计专题 频谱分析 Spectral estimation （谱估计）的目标是基于一个有限的数据集合描述一个信号的功率（在频率上的）分布。功率谱估计在很多场合下都是有用的，包括对宽带噪声湮没下的信号的检测。 从数学上看，一个平稳随机过程n x 的power spectrum （功率谱）和correlation sequence （相关序列）通过discrete-time Fourier transform （离散时间傅立叶变换）构成联系。从normalized frequency （归一化角频率）角度看，有下式 ()()j m xx xx m S R m e ωω∞ -=-∞ = ∑ 注：()() 2 xx S X ωω=，其中 ()/2 /2 lim N j n n N N X x e N ωω=-=∑ πωπ -<≤。 其matlab 近似为X=fft(x,N)/sqrt(N)，在下文中()L X f 就是指matlab fft 函数的计算结果了 使用关系2/s f f ωπ=可以写成物理频率f 的函数，
其中s f 是采样频率 ()()2/s jfm f xx xx m S f R m e π∞ -=-∞ = ∑ 相关序列可以从功率谱用IDFT 变换求得： ()()()/2 2//2 2s s s f jfm f j m xx xx xx s f S e S f e R m d df f πωππ ωωπ- -= =?? 序列n x 在整个Nyquist 间隔上的平均功率可以 表示为 ()()() /2 /2 02s s f xx xx xx s f S S f R d df f ππ ωωπ- -= =?? 上式中的 ()()2xx xx S P ωωπ = 以及()()xx xx s S f P f f = 被定义为平稳随机信号n x 的power spectral density (PSD)（功率谱密度） 一个信号在频带[]1 2 1 2 ,,0ωωωω π ≤<≤上的平均功率 可以通过对PSD 在频带上积分求出 []()()2 1 121 2 ,xx xx P P d P d ωωωωωω ωωωω-- = +?? 从上式中可以看出()xx P ω是一个信号在一个无 穷小频带上的功率浓度，这也是为什么它叫做功率谱密度。
第九章最优化方法的MatIab实现 在生活和工作中，人们对于同一个问题往往会提出多个解决方案，并通过各方面的论证从中提取最佳方案。最优化方法就是专门研究如何从多个方案中科学合理地提取出最佳方案的科学。由于优化问题无所不在，目前最优化方法的应用和研究已经深入到了生产和科研的各个领域，如土木工程、机械工程、化学工程、运输调度、生产控制、经济规划、经济管理等，并取得了显著的经济效益和社会效益。 用最优化方法解决最优化问题的技术称为最优化技术，它包含两个方面的内容： 1）建立数学模型即用数学语言来描述最优化问题。模型中的数学关系式反映了最优化问题所要达到的目标和各种约束条件。 2）数学求解数学模型建好以后，选择合理的最优化方法进行求解。 最优化方法的发展很快，现在已经包含有多个分支，如线性规划、整数规划、非线性规划、动态规划、多目标规划等。 9.1 概述 利用Matlab的优化工具箱，可以求解线性规划、非线性规划和多目标规划问题。 具体而言，包括线性、非线性最小化，最大最小化，二次规划，半无限问题，线性、非线性方程（组）的求解，线性、非线性的最小二乘问题。另外，该工具箱还提供了线性、非线性最小化，方程求解，曲线拟合，二次规划等问题中大型课题的求解方法，为优化方法在工程中的实际应用提供了更方便快捷的途径。 9.1.1优化工具箱中的函数 优化工具箱中的函数包括下面几类： 1 ?最小化函数
2.方程求解函数 3.最小—乘（曲线拟合）函数
4?实用函数 5 ?大型方法的演示函数 6.中型方法的演示函数 9.1.3参数设置 利用OPtimSet函数，可以创建和编辑参数结构；利用OPtimget函数，可以获得o PtiOns优化参数。 ? OPtimget 函数 功能：获得OPtiOns优化参数。 语法：
信号检测与估值 仿真报告 题目信号检测与估值的MATLAB仿真学院通信工程学院 专业通信与信息系统 学生姓名 学号 导师姓名
作业1 试编写程序，画出相干移频键控、非相干移频键控（无衰落）和瑞利衰落信道下非相干移频键控的性能曲线。 （1）根据理论分析公式画性能曲线； （2）信噪比范围（0dB-10dB），间隔是1dB； （3）信噪比计算SNR=10lg（Es/N0） 一、脚本文件 1、主程序 %******************************************************** %二元移频信号检测性能曲线(理论分析) %FSK_theo.m %******************************************************** clear all; clc; SNRindB=0:1:20; Pe_CFSK=zeros(1,length(SNRindB)); Pe_NCFSK=zeros(1,length(SNRindB)); Pe_NCFSK_Rayleigh=zeros(1,length(SNRindB)); for i=1:length(SNRindB) EsN0=exp(SNRindB(i)*log(10)/10); Es_aveN0=exp(SNRindB(i)*log(10)/10); Pe_CFSK(i)=Qfunct(sqrt(EsN0));%相干移频键控系统 Pe_NCFSK(i)=0.5*exp(-EsN0/2);%非相干移频键控系统（无衰落） Pe_NCFSK_Rayleigh(i)=1/(2+Es_aveN0);%非相干移频键控系统（瑞利衰落）end semilogy(SNRindB,Pe_CFSK,'-o',SNRindB,Pe_NCFSK,'-*',SNRindB,Pe_NCFSK_Rayleigh ,'-'); xlabel('Es/No或平均Es/No(dB)'); ylabel('最小平均错误概率Pe'); legend('相干移频','非相干移频(无衰落)','非相干移频(瑞利衰落)'); title('二元移频信号检测性能曲线'); axis([0 20 10^-7 1]); grid on; 2、调用子函数 %******************************************************** %Q函数 %Qfunct.m %********************************************************
核心函数： (1)function [pop]=initializega(num,bounds,eevalFN,eevalOps,options)--初始种群的生成函数 【输出参数】 pop--生成的初始种群 【输入参数】 num--种群中的个体数目 bounds--代表变量的上下界的矩阵 eevalFN--适应度函数 eevalOps--传递给适应度函数的参数 options--选择编码形式(浮点编码或是二进制编码)[precision F_or_B],如 precision--变量进行二进制编码时指定的精度 F_or_B--为1时选择浮点编码，否则为二进制编码,由precision指定精度) (2)function [x,endPop,bPop,traceInfo] = ga(bounds,evalFN,evalOps,startPop,opts,... termFN,termOps,selectFN,selectOps,xOverFNs,xOverOps,mutFNs,mutOps)--遗传算法函数 【输出参数】 x--求得的最优解 endPop--最终得到的种群 bPop--最优种群的一个搜索轨迹 【输入参数】 bounds--代表变量上下界的矩阵 evalFN--适应度函数 evalOps--传递给适应度函数的参数 startPop-初始种群 opts[epsilon prob_ops display]--opts(1:2)等同于initializega的options参数，第三个参数控制是否输出，一般为0。如[1e-6 1 0] termFN--终止函数的名称,如[maxGenTerm] termOps--传递个终止函数的参数,如[100] selectFN--选择函数的名称,如[normGeomSelect] selectOps--传递个选择函数的参数,如[0.08] xOverFNs--交*函数名称表，以空格分开，如[arithXover heuristicXover simpleXover] xOverOps--传递给交*函数的参数表，如[2 02 32 0] mutFNs--变异函数表，如[boundaryMutation multiNonUnifMutation nonUnifMutation unifMutation] mutOps--传递给交*函数的参数表,如[4 0 06 100 34 100 34 0 0] 注意】matlab工具箱函数必须放在工作目录下 【问题】求f(x)=x+10*sin(5x)+7*cos(4x)的最大值，其中0<=x<=9 【分析】选择二进制编码，种群中的个体数目为10，二进制编码长度为20，交*概率为0.95,变异概率为0.08 【程序清单】 编写目标函数 function[sol,eval]=fitness(sol,options) x=sol(1) eval=x+10*sin(5*x)+7*cos(4*x) 把上述函数存储为fitness.m文件并放在工作目录下
实验报告三 实验目的： 实现常规波束形成及基于MUSIC 方法的方位估计。 实验内容： 1）若干阵元的接收阵，信号频率为10KHz ，波束主轴12度，仿真给出常规波束形成的波束图。 2）16个阵元的均匀线列阵，信号频率为10KHz ，信号方位为12度，用MUSIC 方法完成目标定向，信噪比-5dB ，0dB ，5dB 。 i) 波束形成时的阵型设计为两种，一种是均匀线列阵，阵元16个；一种是均匀圆阵，阵元数为16个，比较这两种阵型的波束图。 ii ）比较不同信噪比下MUSIC 方法估计的性能（统计100次）。 实验原理： i)常规波束形成： 如图所示，基阵的输出),(θt v 。 ∑∑=*=* ==M m i i M m i i w t x t x w t v 1 1 ) ()()()(),(θθθ 采用向量符号则有， )()()()(),(H H θθθw x x w t t t v == 式中，x(t)和w(q )分别为观测数据向量和加权系数向量， ) ,(θt v 图 1 波束形成器基本原理图
T M 21])()()([)(t x t x t x t Λ=x T M 21])()()([)(θθθθw w w Λ =w 基阵输出端的空间功率谱表示为： ) ()( )()]()([)( )]()()()([ )],(),([ ] ),([)(H H H H H *2 θθθθθθθθθθRw w w x x w w x x w =====t t E t t E t v t v E t v E P 式中，R 为观测数据的协方差矩阵。 ii ）基于MUSIC 方法的方位估计： )()()()(1 t n t s a t x i d i +=∑=θ T M 21])()()([)(t x t x t x t Λ =x )()()()(t n t s A t x +=θ 假设: (1 ) 信号源的数目d 是已知的, 且d < M ; (2 ) 各信号的方向矢量是相互独立的, 即)(θA 是一个列满秩矩阵; (3 ) 噪声)(t n 是空间平稳随机过程, 为具有各态历经性的均值为零、方差为σ2n 的高斯过程; (4 ) 噪声各取样间是统计独立的。 在上述假设条件下, 基阵输出的协方差矩阵可表示为： I A AR t x t x E R H s H 2])()([α+== 其中, R s 为信号的协方差矩阵；I 为单位矩阵。对R 进行特征分解, 并以特 征值降值排列可得 H m m M d m m H m m d m m e e e e R ∑∑+==+ =1 1λ λ 信号子空间与噪声子空间正交。 若噪声子空间记为E N , 即 ∑+== M d m H m m N e e E 1
最优化方法(Matlab)实验报告 ——Fibonacci 法 一、实验目的： 用MATLAB 程序实现一维搜索中用Fibonacc 法求解一元单峰函数的极小值问题。二、实验原理： （一）、构造Fibonacci 数列：设数列{}k F ，满足条件： 1、011F F == 2、11 k k k F F F +-=+则称数列{}k F 为Fibonacci 数列。（二）、迭代过程： 首先由下面的迭代公式确定出迭代点： 1 1 1 (),1,...,1(),1,...,1n k k k k k n k n k k k k k n k F a b a k n F F u a b a k n F λ---+--+=+ -=-=+ -=-易验证，用上述迭代公式进行迭代时，第k 次迭代的区间长度缩短比率恰好为 1 n k n k F F --+。故可设迭代次数为n ，因此有11121211221111223231 ()()......()()n n n n n n n n n F F F F F F b a b a b a b a b a F F F F F F F ------= -=?-==?-=-若设精度为L ，则有第n 次迭代得区间长度111 ()n n n b a L b a L F -≤-≤，即 就是 111 ()n b a L F -≤，由此便可确定出迭代次数n 。
假设第k 次迭代时已确定出区间[,]k k a b 以及试探点,[,]k k k k u a b λ∈并且k k u λ<。计算试探点处的函数值，有以下两种可能：（1）若()()k k f f u λ>，则令 111111111,,()() () k k k k k k k k n k k k k k n k a b b f f F a b a F λλμλμμ++++--++++-=====+-计算1()k f μ+的值。（2）()()k k f f u λ≤，则令 111121111,,()() () k k k k k k k k n k k k k k n k a a b f f F a b a F μμλμλλ++++--++++-=====+-计算1()k f λ+的值。 又因为第一次迭代确定出了两个迭代点，以后每迭代一次，新增加一个迭代点，这样在迭代n-1后便计算完了n 个迭代点。因此第n 次迭代中，选用第n-1次的迭代点以及辨别常数δ构造n λ和n μ： 1 1n n n n λλμλδ --==+再用同样的方法进行判断：（1）、若()n f λ>()n f μ则令 1 n n n n a b b λ-==（2）、若()n f λ<=()n f μ则令 1n n n n a a b μ-==这样便可确定出最优解的存在区间[,]n n a b 。
Matlab 信号处理工具箱 谱估计专题 频谱分析 Spectral estimation （谱估计）的目标是基于一个有限的数据集合描述一个信号的功率（在频率上的）分布。功率谱估计在很多场合下都是有用的，包括对宽带噪声湮没下的信号的检测。 从数学上看，一个平稳随机过程n x 的power spectrum （功率谱）和correlation sequence （相关序列）通过discrete-time Fourier transform （离散时间傅立叶变换）构成联系。从normalized frequency （归一化角频率）角度看，有下式 ()()j m xx xx m S R m e ωω∞ -=-∞ = ∑ 注：()() 2 xx S X ωω=，其中( )/2 /2 lim N j n n N n N X x e ωω=-=∑ πωπ-<≤。其matlab 近似为X=fft(x,N)/sqrt(N)，在下文中()L X f 就是指matlab fft 函数的计算结果了 使用关系2/s f f ωπ=可以写成物理频率f 的函数，其中s f 是采样频率 ()()2/s jfm f xx xx m S f R m e π∞ -=-∞ = ∑ 相关序列可以从功率谱用IDFT 变换求得： ()()()/2 2//2 2s s s f jfm f j m xx xx xx s f S e S f e R m d df f πωπ π ωωπ--= =? ? 序列n x 在整个Nyquist 间隔上的平均功率可以表示为 ()()() /2 /2 02s s f xx xx xx s f S S f R d df f π π ωωπ--= =? ?
%模拟2ASK % Pe=zeros(1,26); jishu=1; for snr=-10:0.5:15 max = 10000; s=round(rand(1,max));%长度为max的随机二进制序列 f=100;%载波频率 nsamp = 1000;每个载波的取样点数 tc=0:2*pi/999:2*pi;tc的个数应与nsamp相同 cm=zeros(1,nsamp*max); cp=zeros(1,nsamp*max); mod=zeros(1,nsamp*max); for n=1:max; if s(n)==0; m=zeros(1,nsamp); b=zeros(1,nsamp); else if s(n)==1; m=ones(1,nsamp); b=ones(1,nsamp); end end c = sin(f*tc); cm((n-1)*nsamp+1:n*nsamp)=m; cp((n-1)*nsamp+1:n*nsamp)=b; mod((n-1)*nsamp+1:n*nsamp)=c; end tiaoz=cm.*mod;%2ASK调制 t = linspace(0,length(s),length(s)*nsamp); tz=awgn(tiaoz,snr);%信号tiaoz中加入白噪声，信噪比为SNR=10dB jiet = 2*mod.*tz; %相干解调 [N,Wn]=buttord(0.2,0.3,1,15); [b,a]=butter(N,Wn); dpsk=filter(b,a,jiet);%低通滤波 % 抽样判决，判决门限为0.5 depsk = zeros(1,nsamp*max); for m = nsamp/2:nsamp:nsamp*max; if dpsk(m) < 0.5; for i = 1:nsamp depsk((m-500)+i) = 0; end
优化方法上机大作业 学院： 姓名： 学号： 指导老师：肖现涛
第一题 源程序如下： function zy_x = di1ti(x) %di1ti是用来求解优化作业第一题的函数。 x0=x; yimuxulong=0.000001; g0=g(x0);s0=-g0; A=2*ones(100,100); k=0; while k<100 lanmed=-(g0)'*s0/(s0'*A*s0); x=x0+lanmed*s0; g=g(x); k=k+1; if norm(g) break; end miu=norm(g)^2/norm(g0)^2; s=-g+miu*s0; g0=g; s0=s;x0=x; end function f=f(x) f=(x'*ones(100,1))^2-x'*ones(100,1); function g=g(x) g=(2*x'*ones(100,1))*ones(100,1)-ones(100,1); 代入x0,运行结果如下： >> x=zeros(100,1); >> di1ti(x) After 1 iterations,obtain the optimal solution. The optimal solution is -0.250000. The optimal "x" is "ans". ans =0.005*ones(100,1). 实用最优化方法 ——matlab编程作业 题一、 初值为[-1;1] 其中g0、g1分别为不同x值下得导数，f0、f1为函数值 MATLAB程序： x0=[-1;1]; s0=[1;1]; c1=0.1;c2=0.5;a=0;b=inf;d=1;n=0; x1=x0+d*s0; g0=[-400*(x0(2)-x0(1)^2)*x0(1)-2*(1-x0(1));200*(x0(2)-x0(1) ^2)]; g1=[-400*(x1(2)-x1(1)^2)*x1(1)-2*(1-x1(1));200*(x1(2)-x1(1) ^2)]; f1=100*(x1(2)-x1(1)^2)^2+(1-x1(1))^2; f0=100*(x0(2)-x0(1)^2)^2+(1-x0(1))^2; while((f0-f1<-c1*d*g0'*s0)||(g1'*s0 题目：分别用最速下降法、FR 共轭梯度法、DFP 法和BFGS 法求解问题： 22112212minf(x)x 2x x 4x x 3x =-++- 取初始点(1)T x (1,1)=，通过Matlab 编程实现求解过程。 公用函数如下： 1、function f= fun( X ) %所求问题目标函数 f=X(1)^2-2*X(1)*X(2)+4*X(2)^2+X(1)-3*X(2); end 2、function g= gfun( X ) %所求问题目标函数梯度 g=[2*X(1)-2*X(2)+1,-2*X(1)+8*X(2)-3]; end 3、function He = Hess( X ) %所求问题目标函数Hesse 矩阵 n=length(X); He=zeros(n,n); He=[2,-2; -2,4]; End 解法一：最速下降法 function [ x,val,k ] = grad( fun,gfun,x0 ) %功能：用最速下降法求无约束问题最小值 %输入：x0是初始点，fun 和gfun 分别是目标函数和梯度 %输出：x 、val 分别是最优点和最优值，k 是迭代次数 maxk=5000;%最大迭代次数 rho=0.5;sigma=0.4; k=0;eps=10e-6; while (k 第九章最优化方法的Matlab实现 在生活和工作中，人们对于同一个问题往往会提出多个解决方案，并通过各方面的论证从中提取最佳方案。最优化方法就是专门研究如何从多个方案中科学合理地提取出最佳方案的科学。由于优化问题无所不在，目前最优化方法的应用和研究已经深入到了生产和科研的各个领域，如土木工程、机械工程、化学工程、运输调度、生产控制、经济规划、经济管理等，并取得了显著的经济效益和社会效益。 用最优化方法解决最优化问题的技术称为最优化技术，它包含两个方面的内容：1）建立数学模型即用数学语言来描述最优化问题。模型中的数学关系式反映了最优化问题所要达到的目标和各种约束条件。 2）数学求解数学模型建好以后，选择合理的最优化方法进行求解。 最优化方法的发展很快，现在已经包含有多个分支，如线性规划、整数规划、非线性规划、动态规划、多目标规划等。 9.1 概述 利用Matlab的优化工具箱，可以求解线性规划、非线性规划和多目标规划问题。具体而言，包括线性、非线性最小化，最大最小化，二次规划，半无限问题，线性、非线性方程（组）的求解，线性、非线性的最小二乘问题。另外，该工具箱还提供了线性、 非线性最小化，方程求解，曲线拟合，二次规划等问题中大型课题的求解方法，为优化方法在工程中的实际应用提供了更方便快捷的途径。 9.1.1 优化工具箱中的函数 优化工具箱中的函数包括下面几类： 1．最小化函数 表9-1 最小化函数表 2．方程求解函数 表9-2 方程求解函数表 3．最小二乘（曲线拟合）函数 表9-3 最小二乘函数表 4．实用函数 表9-4 实用函数表 5．大型方法的演示函数 表9-5 大型方法的演示函数表 NSGA-II 算法实例 目前的多目标优化算法有很多, Kalyanmoy Deb 的带精英策略的快速非支配排序遗传算法(NSGA-II) 无疑是其中应用最为广泛也是最为成功的一种。本文用的算法是MATLAB 自带的函数gamultiobj ，该函数是基于NSGA-II 改进的一种多目标优化算法。 一、 数值例子 多目标优化问题 424221********* 4224212212112 12min (,)10min (,)55..55 f x x x x x x x x x f x x x x x x x x x s t x =-++-=-++-≤≤??-≤≤? 二、 Matlab 文件 1． 适应值函数m 文件： function y=f(x) y(1)=x(1)^4-10*x(1)^2+x(1)*x(2)+x(2)^4-x(1)^2*x(2)^2; y(2)=x(2)^4-x(1)^2*x(2)^2+x(1)^4+x(1)*x(2); 2． 调用gamultiobj 函数，及参数设置： clear clc fitnessfcn=@f; %适应度函数句柄 nvars=2; %变量个数 lb=[-5,-5]; %下限 ub=[5,5]; %上限 A=[];b=[]; %线性不等式约束 Aeq=[];beq=[]; %线性等式约束 options=gaoptimset('paretoFraction',0.3,'populationsize',100,'generations', 200,'stallGenLimit',200,'TolFun',1e-100,'PlotFcns',@gaplotpareto); % 最优个体系数paretoFraction 为0.3；种群大小populationsize 为100，最大进化代数generations 为200， % 停止代数stallGenLimit 为200， 适应度函数偏差TolFun 设为1e-100，函数gaplotpareto ：绘制Pareto 前端 [x,fval]=gamultiobj(fitnessfcn,nvars,A,b,Aeq,beq,lb,ub,options) methmodeltjut@https://www.doczj.com/doc/e412550215.html, mathematica MATLAB优化应用 §1 线性规划模型 一、线性规划问题： 实例1：生产计划问题 假设某厂计划生产甲、乙两种产品，现库存主要材料有A类3600公斤，B 类2000公斤，C类3000公斤。每件甲产品需用材料A类9公斤，B类4公斤，C类3公斤。每件乙产品，需用材料A类4公斤，B类5公斤，C类10公斤。甲单位产品的利润70元，乙单位产品的利润120元。问如何安排生产，才能使该厂所获的利润最大。 建立数学模型： 设x1、x2分别为生产甲、乙产品的件数。f为该厂所获总润。 max f=70x1+120x2 s.t 9x1+4x2≤3600 4x1+5x2≤2000 3x1+10x2≤3000 x1,x2≥0 实例2：投资问题 某公司有一批资金用于4个工程项目的投资，其投资各项目时所得的净收益(投入资金百分比)如下表：工程项目收益表 由于某种原因，决定用于项目A的投资不大于其他各项投资之和而用于项目B和C的投资要大于项目D的投资。试确定该公司收益最大的投资分配方案。 建立数学模型： 设x1、x2 、x3 、x4分别代表用于项目A、B、C、D的投资百分数。 max f=0.15x1+0.1x2+0.08 x3+0.12 x4 s.t x1-x2- x3- x4≤0 x2+ x3- x4≥0 x1+x2+x3+ x4=1 x j≥0 j=1,2,3,4 实例3：运输问题 有A、B、C三个食品加工厂，负责供给甲、乙、丙、丁四个市场。三个厂每天生产食品箱数上限如下表： 四个市场每天的需求量如下表： 从各厂运到各市场的运输费(元/每箱)由下表给出： 求在基本满足供需平衡的约束条件下使总运输费用最小。 建立数学模型： 设a i j为由工厂i运到市场j的费用，x i j 是由工厂i运到市场j的箱数。b i 是工厂i的产量，d j是市场j的需求量。 -20-15-10-50510152000.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 Kalman滤波器在matlab仿真程序下的不同信噪比时的误码率： %multiuser_dectect.m clc; clear all; hold on BER_sum=zeros(1,13);%设定求和误码率的零矩阵； BER_ave=zeros(1,13); %设定平均误码率的零矩阵； for m=1:10；%m的长度为1到10 间隔为1； snr_in_db=-20:3:16;%定义信噪比的长度为-20到16 间隔为3；snr_in_db是信噪比用db表示 for i=1:length(snr_in_db);%i的长度为1到信噪比的长度 BER(i)= Kalman_S1(snr_in_db(i));%卡尔曼的误码率函数； end BER_sum=BER_sum+BER;%误码率求和的算法 end; BER_ave=0.1* BER_sum ; %误码率平均值的算法 semilogy( snr_in_db,BER_ave,'rd-')；%y轴维数坐标图定义横坐标为信噪比，纵坐标为误码率； %Kalman_S1.m %Kalman algorithm %synchronous CDMA同步cdma %channel: White Gaussis Noise function [p] = Kalman_S1(snr_in_dB) SNR=10^(snr_in_dB/10); %信噪比由dB形式转化 sgma=1; % noise standard deviation is fixed 定义方差 Eb=sgma^2*SNR; A=[sqrt(Eb),sqrt(Eb),sqrt(Eb),sqrt(Eb),sqrt(Eb),sqrt(Eb),sqrt(Eb),sqrt(Eb),sqrt(Eb),sqr t(Eb),sqrt(Eb),sqrt(Eb),sqrt(Eb),sqrt(Eb),sqrt(Eb)]; K=length(A); 功率谱密度估计方法的MATLAB实现 在应用数学和物理学中，谱密度、功率谱密度和能量谱密度是一个用于信号的通用概念，它表示每赫兹的功率、每赫兹的能量这样的物理量纲。在物理学中，信号通常是波的形式，例如电磁波、随机振动或者声波。当波的频谱密度乘以一个适当的系数后将得到每单位频率波携带的功率，这被称为信号的功率谱密度（power spectral density, PSD）或者谱功率分布（spectral power distribution, SPD）。功率谱密度的单位通常用每赫兹的瓦特数（W/Hz）表示，或者使用波长而不是频率，即每纳米的瓦特数（W/nm）来表示。信号的功率谱密度当且仅当信号是广义的平稳过程的时候才存在。如果信号不是平稳过程，那么自相关函数一定是两个变量的函数，这样就不存在功率谱密度，但是可以使用类似的技术估计时变谱密度。信号功率谱的概念和应用是电子工程的基础，尤其是在电子通信系统中，例如无线电和微波通信、雷达以及相关系统。因此学习如何进行功率谱密度估计十分重要，借助于Matlab工具可以实现各种谱估计方法的模拟仿真并输出结果。下面对周期图法、修正周期图法、最大熵法、Levinson递推法和Burg法的功率谱密度估计方法进行程序设计及仿真并给出仿真结果。 以下程序运行平台：Matlab R2015a（8.5.0.197613） 一、周期图法谱估计程序 1、源程序 Fs=100000; %采样频率100kHz N=1024; %数据长度N=1024 n=0:N-1; t=n/Fs; xn=sin(2000*2*pi*t); %正弦波,f=2000Hz Y=awgn(xn,10); %加入信噪比为10db的高斯白噪声 subplot(2,1,1); plot(n,Y) title('信号') xlabel('时间');ylabel('幅度'); 用遗传算法优化BP神经网络的 Matlab编程实例 由于BP网络的权值优化是一个无约束优化问题，而且权值要采用实数编码，所以直接利用Matlab遗传算法工具箱。以下贴出的代码是为一个19输入变量，1个输出变量情况下的非线性回归而设计的，如果要应用于其它情况，只需改动编解码函数即可。 程序一：GA训练BP权值的主函数 function net=GABPNET(XX,YY) %-------------------------------------------------------------------------- % GABPNET.m % 使用遗传算法对BP网络权值阈值进行优化，再用BP 算法训练网络 %-------------------------------------------------------------------------- %数据归一化预处理 nntwarn off XX=premnmx(XX); YY=premnmx(YY); %创建网络 net=newff(minmax(XX),[19,25,1],{'tansig','tansig','purelin'},' trainlm'); %下面使用遗传算法对网络进行优化 P=XX; T=YY; R=size(P,1); S2=size(T,1); S1=25;%隐含层节点数 S=R*S1+S1*S2+S1+S2;%遗传算法编码长度 aa=ones(S,1)*[-1,1]; popu=50;%种群规模 initPpp=initializega(popu,aa,'gabpEval');%初始化种群 gen=100;%遗传代数 %下面调用gaot工具箱，其中目标函数定义为gabpEval [x,endPop,bPop,trace]=ga(aa,'gabpEval',[],initPpp,[1e-6 1 1],'maxGenTerm',gen,... 'normGeomSelect',[0.09],['arithXover'],[2],'nonUnifMutatio n',[2 gen 3]); %绘收敛曲线图 figure(1) plot(trace(:,1),1./trace(:,3),'r-'); hold on plot(trace(:,1),1./trace(:,2),'b-'); xlabel('Generation'); ylabel('Sum-Squared Error'); figure(2) plot(trace(:,1),trace(:,3),'r-'); hold on plot(trace(:,1),trace(:,2),'b-'); xlabel('Generation'); ylabel('Fittness'); %下面将初步得到的权值矩阵赋给尚未开始训练的BP网络 [W1,B1,W2,B2,P,T,A1,A2,SE,val]=gadecod(x); net.LW{2,1}=W1; net.LW{3,2}=W2; net.b{2,1}=B1; net.b{3,1}=B2; XX=P; YY=T; %设置训练参数 net.trainParam.show=1; net.trainParam.lr=1; net.trainParam.epochs=50; net.trainParam.goal=0.001; %训练网络 net=train(net,XX,YY); 程序二：适应值函数 function [sol, val] = gabpEval(sol,options) % val - the fittness of this individual % sol - the individual, returned to allow for Lamarckian evolution % options - [current_generation] load data2 nntwarn off XX=premnmx(XX); YY=premnmx(YY); P=XX; T=YY; R=size(P,1); S2=size(T,1); S1=25;%隐含层节点数 S=R*S1+S1*S2+S1+S2;%遗传算法编码长度 for i=1:S, x(i)=sol(i); end; [W1, B1, W2, B2, P, T, A1, A2, SE, val]=gadecod(x); 实验二 语音信号的频域处理 一、 实验目的、要求 （1）掌握语音信号频域分析方法 （2）了解语音信号频域的特点 （3）了解谱减法作为频域语音增强的原理与编程实现 （3）了解谱减法的缺点，并分析产生该缺点的原因 二、实验原理 语音虽然是一个时变、非平稳的随机过程。但在短时间内可近似看作是平稳的。因此如果能从带噪语音的短时谱中估计出“纯净”语音的短时谱，即可达到语音增强的目的。由于噪声也是随机过程，因此这种估计只能建立在统计模型基础上。利用人耳感知对语音频谱分量的相位不敏感的特性，这类语音增强算法主要针对短时谱的幅度估计。 短时话幅度估计概述 设一帧加窗后的带噪语音为 ()()() 01y n s n d n n N =+≤≤- (2.1) 其中()s n 为纯净语音，()d n 假设为平稳加性高斯噪声。 将()y n 在一组基{()}k n φ上展开，使展对系数为各不相关的随机变量。设()y n 的相关函数为(,)y R n m ，由K －L 展开得知{()}k n φ满足 1 ()()(,)()N k y k m K n R n m m λφφ-==∑ (2.2) 则()y n 的展开式为 1 1 0()()()()N k k K N k k n y n Y n Y y n n φφ-=-=? =????=?? ∑∑ (2.3) 如果()y n 的相关长度小于帧长N ，则()k n φ的近似函数为 2()k nk n j N π??? = ?? ? (2.4) 可见()y n 的展开过程实际上相当于离散博里叶交换，其展开系数（为傅里叶变换系数。由()()()y n s n d n =+，则有：k k k Y S N =+。 其中[]||exp k k k Y Y j θ=、[]||exp k k k S S j α=、k N 分别为()y n 、()s n 及()d n 的傅里叶交换系数。由于假设噪声是高斯分布的，其傅里叶系数k N 相当于多个高最优化方法matlab作业
Matlab最优化编程例子
最优化方法的Matlab实现(公式完整版)
多目标优化实例和matlab程序
MATLAB优化应用
kalman滤波在不同信噪比时的误码率matlab仿真程序
功率谱密度估计方法的MATLAB实现
用遗传算法优化BP神经网络的Matlab编程实例
谱相减MATLAB代码以及信噪比计算
相关主题
文本预览