因為最近專案要使用BP神經網路來做一些飛行預測，所以今天從圖書館借來了《Matlab神經網路30個案例分析》，這本書很不錯推薦給大家，然後研究了下程式碼，使用語音分類這個例子做了原始碼實現與工具箱實現，原始碼實現過程中進行了小小的改變，工具箱用起來非常方便，但是手寫一下BP神經網路的前向後向對於理解BP神經網路還是有極大幫助的，這裡把這兩種實現方式貼出來並帶有結果截圖。顯然BP神經網路對於這種非線性擬合場合效果是非常好的。

     (1) 原始碼實現 

%% 清空環境變數
clc
clear

%% 訓練資料預測資料提取及歸一化

%下載四類語音訊號
load data1 c1
load data2 c2
load data3 c3
load data4 c4

%四個特徵訊號矩陣合成一個矩陣
data(1:500,:) = c1(1:500,:);
data(501:1000,:) = c2(1:500,:);
data(1001:1500,:) = c3(1:500,:);
data(1501:2000,:) = c4(1:500,:);

%從1到2000產生隨機數
k = rands(1,2000);
[m,n] = sort(k);

%%提取輸入、輸出資料
input= data(:,2:25);
output1= data(:,1);

%將輸出資料由一維變為四維
for i = 1:1:2000
	switch output1(i)
		case 1 
			output(i,:) = [1 0 0 0];
		case 2
			output(i,:) = [0 1 0 0];
		case 3
			output(i,:) = [0 0 1 0];
		case 4
			output(i,:) = [0 0 0 1];
	end
end

%隨機提取1500個測試資料，500個樣本為預測資料
input_train = input(n(1:1500),:)';
output_train = output(n(1:1500),:)';
input_test = input(n(1501:2000),:)';
output_test = output(n(1501:2000),:)';

%歸一化
[inputn,inputps] = mapminmax(input_train);

%變數、權值初始化
innum = 24;
midnum = 25;
outnum = 4;

w1 = rands(midnum,innum);
b1 = rands(midnum,1);
w2 = rands(outnum,midnum);
b2 = rands(outnum,1);

w1_1 = w1;
b1_1 = b1;
w2_1 = w2;
b2_1 = b2;

xite = 0.1

%%網路訓練
for ii=1:10
	E(ii)=0;
	for i = 1:1500
		x = inputn(:,i);
		for j = 1:1:midnum
			%%計算隱層值
			I(j) = inputn(:,i)'*w1(j,:)' + b1(j);
			Iout(j) =1/(1+exp(-I(j)));
		end
		
		%%計算輸出層值
		yn = w2*Iout' + b2;
		
		%%計算誤差
		e = output_train(:,i)-yn;
		E(ii) = E(ii) + sum(abs(e));
		
		%%計算權值變化率
		dw2 = e*Iout;
		db2 = e;
		
		for j = 1:1:midnum
			S = 1/(1+exp(-I(j)));
			FI(j) = S*(1-S);
		end
		
		for k = 1:innum
			for j = 1:midnum
				dw1(j,k) = FI(j)*x(k)*(w2(:,j)'*e);
				db1(j) = FI(j)*(w2(:,j)'*e);
			end
		end
		
		%%更新權值
		w1 = w1_1 + xite*dw1;
		b1 = b1_1 + xite*db1';
		w2 = w2_1 + xite*dw2;
		b2 = b2_1 + xite*db2;
		
		w1_1 = w1;
		b1_1 = b1;
		w2_1 = w2;
		b2_1 = b2;
	end
end

%%語音訊號分類
inputn_test = mapminmax('apply',input_test,inputps);

for i=1:1:500
	for j = 1:1:midnum
		I(j) = inputn_test(:,i)'*w1(j,:)' + b1(j);
		Iout(j) = 1/(1+exp(-I(j)));
	end
	
	fore(:,i) = w2*Iout' + b2;

end

%%計算誤差
for i =1:1:500
	output_fore(i) = find(fore(:,i) == max(fore(:,i)));
end
	
error =output_fore - output1(n(1501:2000))';

%畫出預測語音種類和實際語音種類的分類圖
figure(1)
plot(output_fore,'r')
hold on
plot(output1(n(1501:2000))','b')
legend('預測語音類別','實際語音類別')

%畫出誤差圖
figure(2)
plot(error)
title('BP網路分類誤差','fontsize',12)
xlabel('語音訊號','fontsize',12)
ylabel('分類誤差','fontsize',12)

%print -dtiff -r600 1-4

k=zeros(1,4);  
%找出判斷錯誤的分類屬於哪一類
for i=1:1:500
	if error(i)~=0
		[b,c]=max(output_test(:,i));
		switch c
			case 1 
                k(1) = k(1) +1;
			case 2 
                k(2) = k(2) +1;
			case 3
                k(3) = k(3) +1;
			case 4
                k(4) = k(4) +1;
		end
	end
end
				
%找出每類的個體和
kk=zeros(1,4);
for i=1:500
    [b,c]=max(output_test(:,i));
    switch c
        case 1
            kk(1)=kk(1)+1;
        case 2
            kk(2)=kk(2)+1;
        case 3
            kk(3)=kk(3)+1;
        case 4
            kk(4)=kk(4)+1;
    end
end
radio = (kk-k)./kk
			
 

  (2) 工具箱實現

%清空環境變數
clc
clear

%下載輸入輸出資料
load data1 c1
load data2 c2
load data3 c3
load data4 c4

data(1:500,:) = c1(1:500,:);
data(501:1000,:) = c2(1:500,:);
data(1001:1500,:) = c3(1:500,:);
data(1501:2000,:) = c4(1:500,:);

input = data(:,2:25);
output1 = data(:,1);

k = rands(1,2000);
[m,n] = sort(k);

input_train = input(n(1:1500),:)';
output_train = output1(n(1:1500),:)';
input_test = input(n(1501:2000),:)';
output_test = output1(n(1501:2000),:)';

[inputn,inputps] = mapminmax(input_train);
[outputn,outputps] = mapminmax(output_train);

%BP神經網路構建
net = newff(inputn,outputn,25);

%網路引數配置
net.trainParam.epochs = 100;
net.trainParam.lr = 0.1;
net.trainParam.goal = 0.00004;

%BP神經網路訓練
net = train(net,inputn,outputn);

%預測資料歸一化
inputn_test = mapminmax('apply',input_test,inputps);

%BP神經網路預測輸出
an = sim(net,inputn_test);

%輸出結果反歸一化
BPoutput = mapminmax('reverse',an,outputps);

%網路預測結果圖形
figure(1)
plot(BPoutput,':og');
hold on
plot(output_test,'-*');
legend('預測輸出','期望輸出')
title('BP網路預測輸出','fontsize',12)
ylabel('函式輸出','fontsize',12)
xlabel('樣本','fontsize',12)
%預測誤差
error = BPoutput - output_test;
figure(2)
plot(error,'-*')
title('BP網路預測誤差','fontsize',12)
ylabel('誤差','fontsize',12)
xlabel('樣本','fontsize',12)
figure(3)
plot((output_test-BPoutput)./BPoutput,'-*');
title('神經網路預測誤差百分比')
errorsum = sum(abs(error))