当前位置：文档之家› 基于遗传算法BP网络优化

基于遗传算法BP网络优化

班级：控制5班学号：2111504213 姓名：张睿

设计一个BP 神经网络监督控制系统，被控对象为：

1000

()s 87.35s 10470G s s

=++ 采样时间1ms ，输入信号为方波信号，幅值0.5，频率2hz 。设计一个BP 神经网络监督控制系统，并采用遗传算法进行BP 神经网络参数及权值的优化设计，并进行matlab 仿真。需要说明控制系统结构，遗传算法优化BP 网络的具体步骤，并对仿真结果做出分析。

解决过程及思路如下：

1 BP 网络算法

以第p 个样本为例，用于训练的BP 网络结构如图1所示。

图1 具有一个隐含层和输出层的BP 神经网络结构

网络的学习算法如下： (1)信息的正向传播

隐含层神经元的输入为所有输入加权之和，

j ij i i

x w x =∑

隐层神经元的输出'j x 采用S 函数激发j x ，则

()1j j j x x f x e

-==

+ '''(1)

j j j j

x x x x ?=-?

…...

...... (i)

输入层

隐含层

输出层

… …

…

ij w

jk w

x k

j x

输出层的神经元输出为

'k jk j j

x w x =∑

网络输出与理想输出误差为

()k k e k x x =-

误差性能指标函数为

1()2N p k E e k ==∑

上式的N 表示网络输出层的个数。

（2）利用梯度下降法调整各层间权值的反向传播对从第j 个输入到第k 个输出的权值有：

()N

p k jk k j k k jk

jk E x w e e k x w w η

ηη==???=-==??∑

∑ 其中，η为学习速率，[]0,1η∈。

K+1时刻网络权值为

(1)()jk jk jk w k w k w +=+?

对从第i 个输入到第j 个输出的权值有：

p k

ij k

k ij

E x w e w w η

η=???=-=??∑ 式中，'

'''

(1)j j

k k jk j j i ij j j ij

x x x x w x x x w x x w ????=??=??-????? t+1时刻网络权值为

(1)()ij ij ij w k w k w +=+?

2.BP 网络的监督控制系统结构

设计的BP 网络监督控制系统结构如图2所示。

图2 BP 神经网络监督控制

在BP 网络结构中，取网络的输入为r （k ），实际输出为y （k ），PID 控制输出为up （k ），隐层神经元的输出采用S 函数激发，

()1j j j x x f x e

-==

+ 网络的权向量为W1, W2 BP 的网络输出为

控制律为

u(k)=up(k)+yn(k)

采用梯度下降法调整网络的权值为

神经网络权值的调整过程为

3.遗传算法对BP 网络权值的优化过程（1）取逼近总步骤为:N=100 （2）终止代数：G=80 （3）样本个数：Size=30 （4）交叉概率：Pc=0.70

（5）二进制编码长度：Godel=10

（6）变异概率：Pm=0.001-[1:1:Size]*0.001/Size

（7）用于优化的BP 网络结构为：1-4-1 i=1 j=1,2,3 ,4 （8）网络权值W1的取值范围为：[-1,+1] （9）网络权值W2的取值范围为：[-0.5,+0.5]

(10)取BP 网络误差绝对值为参数选择的最小目标函数：

∑=j

j x w '2k yn ）（'222)()(j

j k j j x k e w x k e w E w ??-=????-=??-=?ηηηij n ij ij w y k e w E w ???

?-=??-=?)(ηη))1()(()()1(22222--+?+=+k w k w w k w k w j j j j j α))

1()(()()1(--+?+=+k w k w w k w k w ij ij ij ij ij α

式中，N 为逼近的总步骤，ee （i ）为第i 步BP 网络逼近误差。

（11）需要优化参数为：],,,,,,,,,,,[p 232221201716151413121110j w w w w w w w w w w w w 4 遗传算法优化BP 网络权值的步骤 (1)初始化种群；

(2)计算其适应值，保留最优个体，判断是否达到最优解； (3)交叉、变异产生新个体；

(4)重新计算种群中每个个体的适应值并保留最优个体；

(5)交叉、变异前后的种群放在一起进行二人竞赛选择法，直到填满新的种群；

(6)转2)直到找到最优解BestJ 。 5 MATLAB 仿真结果

BP 监督网络遗传算法优化程序包括3部分，即遗传算法优化程序ga_bp.m ，BP 网络逼近函数程序bp_a 和BP 网络逼近测试程序bp_test 。

输入信号为r （t ）=0.5*sign(sin(2*2*pi*k*ts)采样时间ts=0.001s ，η=0.30，a=0.05,kp=1,kd=1.

经遗传算法优化后，对象p 的值为

P=[-0.2160,0.7576,0.5230,0.9863,-0.0714,0.2551,0.6911,-0.3627，0.2146，0.3338 ，-0.0875，-0.0582]

仿真结果图：

4050

790

800810820830840850

860870880890Times

B e s t J

图3 代价函数J 的优化过程

00.10.20.30.4

0.50.60.70.80.91

-0.8

-0.6

-0.4

-0.2

0.2

0.4

0.6

time(s)

r a n d y

图4 方波位置跟踪

0.1

0.2

0.3

0.4

0.50.6

0.7

0.8

0.9

-10010

time(s)

y n

0.1

0.2

0.3

0.4

0.50.6

0.7

0.8

0.9

-50050

time(s)

u p

0.1

0.2

0.3

0.4

0.50.6

0.7

0.8

0.9

-10010time(s)

u 0

0.1

0.2

0.3

0.4

0.50.6

0.7

0.8

0.9

-101time(s)

e r r o r

图5 BP 网络,PD 及总控制器输出的比较以及误差曲线

结论：采用遗专算法可以实现BP 网络参数初始值的优化，节约计算量。并由仿真结果可知，其误差大部分趋于0，但局部有三个地区的误差比较大，产生原因可能与遗传算法的运行参数有关。代码：

1、bp_a.m

function[p,BsJ]=rbf_gaf(p,BsJ)

ts=0.001;

alfa=0.05;

xite=0.30;

sys=tf(1000,[1,87.35,10470]);

dsys=c2d(sys,ts,'z');

[num,den]=tfdata(dsys,'v');

y_1=0;y_2=0;

u_1=0;u_2=0;

e_1=0;

xi=0;

x=[0,0]';

I=[0,0,0,0]';

Iout=[0,0,0,0]';

FI=[0,0,0,0]';

kp=25;

kd=0.3;

w1=[p(1),p(2),p(3),p(4);p(5),p(6),p(7),p(8)];

w1_1=w1;w1_2=w1;

w2=[p(9);p(10);p(11);p(12)];

w2_1=w2;w2_2=w2_1;

for k=1:1:1000

timef(k)=k*ts;

Y=1;

if Y==1

r(k)=0.5*sign(sin(2*2*pi*k*ts)); %Square Signal elseif Y==2

r(k)=0.5*(sin(3*2*pi*k*ts)); %Square Signal end

y(k)=-den(2)*y_1-den(3)*y_2+num(2)*u_1+num(3)*u_2; e(k)=r(k)-y(k);

xi=r(k);

for j=1:1:4

I(j)=x'*w1(:,j);

Iout(j)=1/(1+exp(-I(j)));

end

yn(k)=w2'*Iout; % Output of NNI networks

%PD Controller

up(k)=kp*x(1)+kd*x(2);

M=2;

if M==1 %Only Using PID Control

u(k)=up(k);

elseif M==2 %Total control output

u(k)=up(k)+yn(k);

end

if u(k)>=10

u(k)=10;

end

if u(k)<=-10

u(k)=-10;

end

if k==400

u(k)=u(k)+5.0;

end

%Update NN Weight

ee(k)=u(k)-yn(k);

w2=w2_1+(xite*ee(k))*Iout+alfa*(w2_1-w2_2);

for j=1:1:4

FI(j)=exp(-I(j))/(1+exp(-I(j)))^2;

end

dw1=0*w1;

for i=1:1:2

for j=1:1:4

dw1(i,j)=ee(k)*xite*FI(j)*w2(j)*x(i);

end

w1=w1_1+dw1+alfa*(w1_1-w1_2);

w1_2=w1_1;w1_1=w1;

w2_2=w2_1;w2_1=w2;

u_2=u_1;

u_1=u(k);

y_2=y_1;

y_1=y(k);

x(1)=e(k); %Calculating P

x(2)=(e(k)-e_1)/0.05; %Calculating D

e_1=e(k);

end

B=0;

for i=1:1:100

Ji(i)=abs(ee(i));

B=B+50*Ji(i);

end

BsJ=B;

2、bp_test.m

clear all;

close all;

load pfile;

ts=0.001;

alfa=0.05;

xite=0.30;

kp=25;

kd=0.3;

sys=tf(1000,[1,87.35,10470]);

dsys=c2d(sys,ts,'z');

[num,den]=tfdata(dsys,'v');

N=1;

if N==1

w1=[p(1),p(2),p(3),p(4);p(5),p(6),p(7),p(8)]; w2=[p(9);p(10);p(11);p(12)];

else if N==2

w1=rand(2,4);

w2=rand(1,4);

end

w1_1=w1;w1_2=w1;

w2_1=w2;w2_2=w2_1;

y_1=0;y_2=0;

u_1=0;u_2=0;

e_1=0;

x=[0,0]';

I=[0,0,0,0]';

Iout=[0,0,0,0]';

FI=[0,0,0,0]';

ts=0.001;

for k=1:1:1000

time(k)=k*ts;

Y=1;

if Y==1

r(k)=0.5*sign(sin(2*2*pi*k*ts)); %Square Signal elseif Y==2

r(k)=0.5*(sin(3*2*pi*k*ts)); %Square Signal end

y(k)=-den(2)*y_1-den(3)*y_2+num(2)*u_1+num(3)*u_2; e(k)=r(k)-y(k);

for j=1:1:4

I(j)=x'*w1(:,j);

Iout(j)=1/(1+exp(-I(j)));

end

yn(k)=w2'*Iout; % Output of NNI networks

%PD Controller

up(k)=kp*x(1)+kd*x(2);

M=2;

if M==1 %Only Using PID Control

u(k)=up(k);

elseif M==2 %Total control output

u(k)=up(k)+yn(k);

end

if u(k)>=10

u(k)=10;

end

if u(k)<=-10

u(k)=-10;

end

if k==400

u(k)=u(k)+6.0;

end

%Update NN Weight

ee(k)=u(k)-yn(k);

w2=w2_1+(xite*ee(k))*Iout+alfa*(w2_1-w2_2);

for j=1:1:4

FI(j)=exp(-I(j))/(1+exp(-I(j)))^2;

end

dw1=0*w1;

for i=1:1:2

for j=1:1:4

dw1(i,j)=ee(k)*xite*FI(j)*w2(j)*x(i);

end

w1=w1_1+dw1+alfa*(w1_1-w1_2);

w1_2=w1_1;w1_1=w1;

w2_2=w2_1;w2_1=w2;

u_2=u_1;

u_1=u(k);

y_2=y_1;

y_1=y(k);

x(1)=e(k); %Calculating P

x(2)=(e(k)-e_1)/0.05; %Calculating D

e_1=e(k);

end

figure(1);

plot(time,r,'r',time,y,'b');

xlabel('time(s)');ylabel('r and y');

figure(2);

subplot(411);

plot(time,yn,'b');

xlabel('time(s)');ylabel('yn');

subplot(412);

plot(time,up,'k');

xlabel('time(s)');ylabel('up');

subplot(413);

plot(time,u,'r');

xlabel('time(s)');ylabel('u');

subplot(414);

plot(time,r-y,'b');

xlabel('time(s)');ylabel('error');

3、ga_bp.m

%GA(Generic Algorithm) to Optimize Initial Parameters of RBF Approaching clear all;

close all;

Size=30; % èoì?×ü??ì?êyá?

G=80; % ?êDí×?′óμü′ú′?êy

CodeL=10; % ????±?á?óμóDμ?è?é?ì?êyá?

for i=1:1:8

umax(i)=ones(1); % á??óè¨?μμ?×?′ó?μ

umin(i)=-ones(1); % è¨?μμ?×?D??μ

end

for i=9:1:12

umax(i)=0.5*ones(1); % ?D?μμ?×?′ó?μ

umin(i)=-0.5*ones(1); % ?D?μμ?×?D??μ

end

% ???ù?ó2?êy??DD3?ê??t????±à??

E=round(rand(Size,12*CodeL)); % 3?ê??ˉ??ì?′ú??￡¨è?é?ì?￡?

% ò?′?é??-í????μá·?aê?

BsJ=0;

for kg=1:1:G

time(kg)=kg;

for s=1:1:Size

m=E(s,:);

for j=1:1:12

yoy(j)=0;

mj=m((j-1)*CodeL+1:1:j*CodeL);

for i=1:1:CodeL

yoy(j)=yoy(j)+mj(i)*2^(i-1);

end

f(s,j)=(umax(j)-umin(j))*yoy(j)/1023+umin(j);

end

clear yoy;

p=f(s,:);

[p,BsJ]=bp_a(p,BsJ);

BsJi(s)=BsJ;

end

[OderJi,IndexJi]=sort(BsJi);

BestJ(kg)=OderJi(1);

BJ=BestJ(kg);

Ji=BsJi+1e-005;

fi=1./Ji;

[Oderfi,Indexfi]=sort(fi); %Arranging fi small to bigger

Bestfi=Oderfi(Size); % Let Bestfi=max(fi)

BestS=E(Indexfi(Size),:); % Let BestS=E(m), m is the Indexfi belong to max(fi)

%****** Step 2 : Select and Reproduct Operation******

fi_sum=sum(fi);

fi_Size=(Oderfi/fi_sum)*Size;

fi_S=floor(fi_Size);

kk=1;

for i=1:1:Size

for j=1:1:fi_S(i)

TempE(kk,:)=E(Indexfi(i),:);

kk=kk+1;

end

%************ Step 3 : Crossover Operation ************

pc=0.70;

n=ceil(20*rand);

for i=1:2:(Size-1)

temp=rand;

if pc>temp %Crossover Condition

for j=n:1:20

TempE(i,j)=E(i+1,j);

TempE(i+1,j)=E(i,j);

end

TempE(Size,:)=BestS;

E=TempE;

%************ Step 4: Mutation Operation ************** pm=0.001-[1:1:Size]*(0.001)/Size; %Bigger fi, smaller pm for i=1:1:Size

for j=1:1:12*CodeL

temp=rand;

if pm>temp %Mutation Condition

if TempE(i,j)==0

TempE(i,j)=1;

else

TempE(i,j)=0;

end

%Guarantee TempE(Size,:) belong to the best individual TempE(Size,:)=BestS;

E=TempE;

%******************************************************* end

Bestfi

BestS

Best_J=BestJ(G)

figure(1);

plot(time,BestJ);

xlabel('Times');ylabel('Best J');

save pfile p;