mlp_bp1.m

基本bp算法,实现函数的拟合,提供给大家分享

%BP算法作业
%题目一：训练函数y=sinx，x∈（0，2π）
%%%%%%%%%输入输出变换采用sigmoid函数，采用9个训练样本，361个测试样本
%%%%%%%%%性能指标精度ε取0.01
%%%%%%%%%BP采用单隐层
function MLP_BP1
global N Nt IN HN ON miu eita ep;
[w_H,w_Hold,w_O,w_Oold,theta_H,theta_Hold,theta_O,theta_Oold]=ini_BP;
[P,T]=get_Trainsample();
judge=0;%判断标志
epoch=1;%训练周期数
iN=1;%标志训练样本
delta_O=zeros(ON,N);%训练误差
delta_H=zeros(HN,N);%训练误差
while judge==0  
    for iN=1:N
        [input_H,output_H,input_O,output_O]=get_Trainoutput(w_H,P,theta_H,w_O,theta_O,iN);
        [delta_O,delta_H]=train_Error(T,output_O,w_O,output_H,iN);
        [w_H,w_Hold,w_O,w_Oold,theta_H,theta_Hold,theta_O,theta_Oold]=modify_Para(delta_H,P,delta_O,output_H,w_H,w_Hold,w_O,w_Oold,theta_H,theta_Hold,theta_O,theta_Oold,iN);
    end   %完成一次训练（一个训练周期）
    [judge,error]=E_Judge(T,output_O);
    error_Store(epoch)=error;
    epoch=epoch+1;
end
test_Output;


%初始化BP网络
function [w_H,w_Hold,w_O,w_Oold,theta_H,theta_Hold,theta_O,theta_Oold]=ini_BP
N=9;%训练样本个数
Nt=361;%测试样本个数

IN=1;%输入层神经元数目
HN=4;%隐层神经元个数
ON=1;%输出层神经元个数

%alpha_O=0.01;%输出层至隐层的学习效率
%alpha_H=0.01;%隐层至输入层学习效率

ep=0.01;%设定精度

miu=0.1;%μ
eita=0.05;%η

%以下初始化权值和阈值（Nguyen和Widrow的权值初始化算法）
lamda=0.7^IN*sqrt(HN);
w_H=rand(HN,IN)-0.5*ones(HN,IN);%在-0.5至0.5范围内
w_Hold=zeros(HN,IN);
w_O=rand(ON,HN)-0.5*ones(HN,IN);%在-0.5至0.5范围内
w_Oold=zeros(ON,HN);
%输入层到隐层的权值初始化
for i=1:HN    
    for j=1:IN
        sum=0;
        for k=1:HN
            sum=w_H(k,j)^2+sum;
        end
        w_H(i,j)=lamda.*w_H(i,j)./sqrt(sum);
    end
end
%隐层到输出层的权值初始化
for i=1:ON    
    for j=1:HN
        sum=0;
        for k=1:ON
            sum=w_O(k,j)^2+sum;
        end
        w_O(i,j)=lamda.*w_O(i,j)./sqrt(sum);
    end
end
%隐层阈值初始化
theta_H=zeros(HN,1);
theta_Hold=theta_H;
for i=1:HN
    j=ceil(IN*rand);
    theta_H(i,1)=w_H(i,j)*(2*rand-1);
end
%输出层阈值初始化
theta_O=zeros(ON,1);
theta_Oold=theta_O;
for i=1:ON
    j=ceil(HN*rand);
    theta_O(i,1)=w_O(i,j)*(2*rand-1);
end
%训练样本初始化
function [P,T]=get_Trainsample()
global IN ON N;
P=zeros(IN,N);%单个样本输入数据
T=zeros(ON,N);%单个样本教师数据,期望输出
for i=1:IN    
    P(i,:)=0:(2*pi/(N-1)):(2*pi);  
end
for i=1:ON
    T(i,:)=sin(P(i,:));
end
plor(P,T,'+');
%计算各层输出
function  [input_H,output_H,input_O,output_O]=get_Trainoutput(w_H,P,theta_H,w_O,theta_O,iN)
global IN HN ON N;
input_H=zeros(HN,N);%隐层的输入
output_H=zeros(HN,N);%隐层的输出
input_O=zeros(ON,N);%输出层的输入
output_O=zeros(ON,N);%输出层的输出
%计算隐层输入输出
input_H(1)=0;
for i=1:HN
    for j=1:IN
        input_H(i,iN)=w_H(i,j)*(j,iN)+input_H(i,iN);
    end
    output_H(i,iN)=input_H(i,iN)-theta_H(i);
    output_H(i,iN)=1/(1+exp(-output_H(i,iN)));
end
%计算输出层输入输出
input_O(1)=0;
for i=1:ON
    for j=1:HN
        input_O(i,iN)=w_O(i,j)*input_H(j,iN)+input_O(i,iN);
    end
    output_O(i,iN)=input_O(i,iN)-theta_O(i);
    output_O(i,iN)=1/(1+exp(-output_O(i,iN)));
end
return;
%计算各层训练误差
function [delta_O,delta_H]=train_Error(T,output_O,w_O,output_H,iN)
global ON HN;
%计算输出层训练误差
for i=1:ON    
    delta_O(i,iN)=(T(i,iN)-output_O(i,iN))*output_O(i,iN)*(1-output_O(i,iN));
end
%计算隐层训练误差
for i=1:HN
    for j=1:ON
    delta_H(i,iN)=delta_H(i,iN)+delta_O(j,iN)*w_O(j,i);
    end
    delta_H(i,iN)=output_H(i,iN)*(1-output_H(i,iN))*delta_H(i,iN);
end

%修正权值和阈值
function [w_H,w_Hold,w_O,w_Oold,theta_H,theta_Hold,theta_O,theta_Oold]=modify_Para(delta_H,P,delta_O,output_H,w_H,w_Hold,w_O,w_Oold,theta_H,theta_Hold,theta_O,theta_Oold,iN)
global IN HN ON miu eita;
%修正权值
%输入层到隐层
w_Hold=w_H;
for i=1:HN
    for j=1:IN
        w_H(i,j)=w_H(i,j)+miu*delta_H(i,iN)*P(j,iN)+eita*(w_H(i,j)-w_Hold(i,j));
    end
end
%隐层到输出层
w_Oold=w_O;
for i=1:ON
    for j=1:HN
        w_O(i,j)=w_O(i,j)+miu*delta_O(i,iN)*output_H(j,iN)+eita*(w_O(i,j)-w_Oold(i,j));
    end
end
%修正阈值
%隐层
theta_Hold=theta_H;
for i=1:HN
    theta_H(i)=theta_H(i)+miu*delta_H(i,iN)+eita*(theta_H(i)-theta_Hold(i));
end
%输出层
theta_Oold=theta_O;
for i=1:ON
    theta_O(i)=theta_O(i)+miu*delta_O(i,iN)+eita*(theta_O(i)-theta_Oold(i));
end
%判断精度有无满足要求
function [judge,error]=E_Judge(T,output_O)
global N ON ep;
error=0;%样本误差(代价函数的值)
for j=1:N
    for i=1:ON
        error=0.5*(T(i,j)-output_O(i,j))^2+error;
    end
end
if error <= ep
    judge=1;
else
    judge=0;
end
return;

function test_Output(epoch,error_Store,)
global IN ON Nt;
Pt=zeros(IN,Nt);%测试单个样本输入数据
Tt=zeros(ON,Nt);
Tt=zeros(ON,Nt);%测试单个样本输出数据
for i=1:IN    
    P(i,:)=0:(2*pi/(Nt-1)):(2*pi);  
end
for i=1:ON
    T(i,:)=sin(P(i,:));
end
plor(P,T,'+');