徑向基網絡（RBF network）

來源：http://blog.csdn.net/zouxy09/article/details/13297881

1、徑向基函數

徑向基函數（Radical Basis Function，RBF）方法是Powell在1985年提出的。所謂徑向基函數，其實就是某種沿徑向對稱的標量函數。通常定義為空間中任一點x到某一中心c之間歐氏距離的單調函數，可記作k(||x-c||)，其作用往往是局部的，即當x遠離c時函數取值很小。例如高斯徑向基函數：

       當年徑向基函數的誕生主要是為了解決多變量插值的問題。可以看下面的圖。具體的話是先在每個樣本上面放一個基函數，圖中每個藍色的點是一個樣本，然后中間那個圖中綠色虛線對應的，就表示的是每個訓練樣本對應一個高斯函數（高斯函數中心就是樣本點）。然后假設真實的擬合這些訓練數據的曲線是藍色的那根（最右邊的圖），如果我們有一個新的數據x1，我們想知道它對應的f(x1)是多少，也就是a點的縱坐標是多少。那么由圖可以看到，a點的縱坐標等於b點的縱坐標加上c點的縱坐標。而b的縱坐標是第一個樣本點的高斯函數的值乘以一個大點權值得到的，c的縱坐標是第二個樣本點的高斯函數的值乘以另一個小點的權值得到。而其他樣本點的權值全是0，因為我們要插值的點x1在第一和第二個樣本點之間，遠離其他的樣本點，那么插值影響最大的就是離得近的點，離的遠的就沒什么貢獻了。所以x1點的函數值由附近的b和c兩個點就可以確定了。拓展到任意的新的x，這些紅色的高斯函數乘以一個權值后再在對應的x地方加起來，就可以完美的擬合真實的函數曲線了。

二、徑向基網絡

       到了1988年， Moody和 Darken提出了一種神經網絡結構，即RBF神經網絡，屬於前向神經網絡類型，它能夠以任意精度逼近任意連續函數，特別適合於解決分類問題。

       RBF網絡的結構與多層前向網絡類似，它是一種三層前向網絡。輸入層由信號源結點組成；第二層為隱含層，隱單元數視所描述問題的需要而定，隱單元的變換函數是RBF徑向基函數，它是對中心點徑向對稱且衰減的非負非線性函數；第三層為輸出層，它對輸入模式的作用作出響應。從輸人空間到隱含層空間的變換是非線性的，而從隱含層空間到輸出層空間變換是線性的。

 

       RBF網絡的基本思想是：用RBF作為隱單元的“基”構成隱含層空間，這樣就可將輸入矢量直接（即不需要通過權連接）映射到隱空間。根據Cover定理，低維空間不可分的數據到了高維空間會更有可能變得可分。換句話來說，RBF網絡的隱層的功能就是將低維空間的輸入通過非線性函數映射到一個高維空間。然后再在這個高維空間進行曲線的擬合。它等價於在一個隱含的高維空間尋找一個能最佳擬合訓練數據的表面。這點與普通的多層感知機MLP是不同的。

       當RBF的中心點確定以后，這種映射關系也就確定了。而隱含層空間到輸出空間的映射是線性的，即網絡的輸出是隱單元輸出的線性加權和，此處的權即為網絡可調參數。由此可見，從總體上看，網絡由輸人到輸出的映射是非線性的，而網絡輸出對可調參數而言卻又是線性的。這樣網絡的權就可由線性方程組直接解出，從而大大加快學習速度並避免局部極小問題。

       從另一個方面也可以這樣理解，多層感知器（包括BP神經網絡）的隱節點基函數采用線性函數，激活函數則采用Sigmoid函數或硬極限函數。而RBF網絡的隱節點的基函數采用距離函數（如歐氏距離），並使用徑向基函數（如Gaussian函數）作為激活函數。徑向基函數關於n維空間的一個中心點具有徑向對稱性，而且神經元的輸入離該中心點越遠，神經元的激活程度就越低。隱節點的這一特性常被稱為“局部特性”。

 

三、RBF網絡的設計與求解

       RBF的設計主要包括兩個方面，一個是結構設計，也就是說隱藏層含有幾個節點合適。另一個就是參數設計，也就是對網絡各參數進行求解。由上面的輸入到輸出的網絡映射函數公式可以看到，網絡的參數主要包括三種：徑向基函數的中心、方差和隱含層到輸出層的權值。到目前為止，出現了很多求解這三種參數的方法，主要可以分為以下兩大類：

1、方法一：

       通過非監督方法得到徑向基函數的中心和方差，通過監督方法（最小均方誤差）得到隱含層到輸出層的權值。具體如下：

（1）在訓練樣本集中隨機選擇h個樣本作為h個徑向基函數的中心。更好的方法是通過聚類，例如K-means聚類得到h個聚類中心，將這些聚類中心當成徑向基函數的h個中心。

（2）RBF神經網絡的基函數為高斯函數時，方差可由下式求解：

       式中c_max 為所選取中心之間的最大距離，h是隱層節點的個數。擴展常數這么計算是為了避免徑向基函數太尖或太平。

（3）隱含層至輸出層之間神經元的連接權值可以用最小均方誤差LMS直接計算得到，計算公式如下：（計算偽逆）（d是我們期待的輸出值）

 

2、方法二：

         采用監督學習算法對網絡所有的參數（徑向基函數的中心、方差和隱含層到輸出層的權值）進行訓練。主要是對代價函數（均方誤差）進行梯度下降，然后修正每個參數。具體如下：

（1）隨機初始化徑向基函數的中心、方差和隱含層到輸出層的權值。當然了，也可以選用方法一中的（1）來初始化徑向基函數的中心。

（2）通過梯度下降來對網絡中的三種參數都進行監督訓練優化。代價函數是網絡輸出和期望輸出的均方誤差：

然后每次迭代，在誤差梯度的負方向已一定的學習率調整參數。

 

四、代碼實現：

1、第一種方法

       第一種方法在zhangchaoyang的博客上面有C++的實現，只是上面針對的是標量的數據（輸入和輸出都是一維的）。而在Matlab中也提供了第一種方法的改進版（呵呵，個人覺得，大家可以在Matlab中運行open newrb查看下源代碼）。

       Matlab提供的一個函數是newrb()。它有個技能就是可以自動增加網絡的隱層神經元數目直到均方差滿足我們要求的精度或者神經元數數目達到最大（也就是我們提供的樣本數目，當神經元數目和我們的樣本數目一致時，rbf網絡此時的均方誤差為0）為止。它使用方法也能簡單：

rbf = newrb(train_x, train_y);

output = rbf(test_x);

      直接把訓練樣本給它就可以得到一個rbf網絡了。然后我們把輸入給它就可以得到網絡的輸出了。

 

2、第二種方法

       第二種方法在zhangchaoyang的博客上面也有C++的實現，只是上面針對的還是標量的數據（輸入和輸出都是一維的）。但我是做圖像的，網絡需要接受高維的輸入，而且在Matlab中，向量的運算要比for訓練的運算要快很多。所以我就自己寫了個可以接受向量輸入和向量輸出的通過BP算法監督訓練的版本。BP算法可以參考這里：BackpropagationAlgorithm ，主要是計算每層每個節點的殘差就可以了。另外，我的代碼是可以通過梯度檢查的，但在某些訓練集上面，代價函數值卻會隨着迭代次數上升，這就很奇怪了，然后降低了學習率還是一樣。但在某些簡單點的訓練集上面還是可以工作的，雖然訓練誤差也挺大的（沒有完全擬合訓練樣本）。所以大家如果發現代碼里面有錯誤的部分，還望大家告知下。

  主要代碼見下面：

learnRBF.m

 

[cpp] view plain copy

1. %// This is a RBF network trained by BP algorithm    

2. %// Author : zouxy    

3. %// Date   : 2013-10-28    

4. %// HomePage : http://blog.csdn.net/zouxy09    

5. %// Email  : [e-mail]    

6.   

7. close all; clear; clc;  

8.   

9. %%% ************************************************  

10. %%% ************ step 0: load data ****************  

11. display('step 0: load data...');  

12. % train_x = [1 2 3 4 5 6 7 8]; % each sample arranged as a column of train_x  

13. % train_y = 2 * train_x;  

14. train_x = rand(5, 10);  

15. train_y = 2 * train_x;  

16. test_x = train_x;  

17. test_y = train_y;  

18.   

19. %% from matlab  

20. % rbf = newrb(train_x, train_y);  

21. % output = rbf(test_x);  

22.   

23.   

24. %%% ************************************************  

25. %%% ******** step 1: initialize parameters ********   

26. display('step 1: initialize parameters...');  

27. numSamples = size(train_x, 2);  

28. rbf.inputSize = size(train_x, 1);  

29. rbf.hiddenSize = numSamples;        % num of Radial Basis function  

30. rbf.outputSize = size(train_y, 1);  

31. rbf.alpha = 0.1;  % learning rate (should not be large!)  

32.   

33. %% centre of RBF  

34. for i = 1 : rbf.hiddenSize  

35.     % randomly pick up some samples to initialize centres of RBF  

36.     index = randi([1, numSamples]);   

37.     rbf.center(:, i) =  train_x(:, index);  

38. end  

39.   

40. %% delta of RBF  

41. rbf.delta = rand(1, rbf.hiddenSize);  

42.   

43. %% weight of RBF  

44. r = 1.0; % random number between [-r, r]  

45. rbf.weight = rand(rbf.outputSize, rbf.hiddenSize) * 2 * r - r;  

46.   

47.   

48. %%% ************************************************  

49. %%% ************ step 2: start training ************  

50. display('step 2: start training...');  

51. maxIter = 400;  

52. preCost = 0;  

53. for i = 1 : maxIter  

54.     fprintf(1, 'Iteration %d ,', i);  

55.     rbf = trainRBF(rbf, train_x, train_y);  

56.     fprintf(1, 'the cost is %d \n', rbf.cost);  

57.       

58.     curCost = rbf.cost;  

59.     if abs(curCost - preCost) < 1e-8  

60.         disp('Reached iteration termination condition and Termination now!');  

61.         break;  

62.     end  

63.     preCost = curCost;  

64. end  

65.   

66.   

67. %%% ************************************************  

68. %%% ************ step 3: start testing ************   

69. display('step 3: start testing...');  

70. Green = zeros(rbf.hiddenSize, 1);  

71. for i = 1 : size(test_x, 2)  

72.     for j = 1 : rbf.hiddenSize  

73.         Green(j, 1) = green(test_x(:, i), rbf.center(:, j), rbf.delta(j));  

74.     end   

75.     output(:, i) = rbf.weight * Green;  

76. end  

77. disp(test_y);  

78. disp(output);  

 

 

trainRBF.m

 

[cpp] view plain copy

1. function [rbf] = trainRBF(rbf, train_x, train_y)  

2.     %%% step 1: calculate gradient  

3.     numSamples = size(train_x, 2);  

4.     Green = zeros(rbf.hiddenSize, 1);  

5.     output = zeros(rbf.outputSize, 1);  

6.     delta_weight = zeros(rbf.outputSize, rbf.hiddenSize);  

7.     delta_center = zeros(rbf.inputSize, rbf.hiddenSize);  

8.     delta_delta =  zeros(1, rbf.hiddenSize);  

9.     rbf.cost = 0;  

10.     for i = 1 : numSamples  

11.         %% Feed forward  

12.         for j = 1 : rbf.hiddenSize  

13.             Green(j, 1) = green(train_x(:, i), rbf.center(:, j), rbf.delta(j));  

14.         end   

15.         output = rbf.weight * Green;      

16.           

17.         %% Back propagation  

18.         delta3 = -(train_y(:, i) - output);  

19.         rbf.cost = rbf.cost + sum(delta3.^2);  

20.         delta_weight = delta_weight + delta3 * Green';  

21.         delta2 = rbf.weight' * delta3 .* Green;  

22.         for j = 1 : rbf.hiddenSize  

23.             delta_center(:, j) = delta_center(:, j) + delta2(j) .* (train_x(:, i) - rbf.center(:, j)) ./ rbf.delta(j)^2;  

24.             delta_delta(j) = delta_delta(j)+ delta2(j) * sum((train_x(:, i) - rbf.center(:, j)).^2) ./ rbf.delta(j)^3;  

25.         end  

26.     end  

27.   

28.     %%% step 2: update parameters  

29.     rbf.cost = 0.5 * rbf.cost ./ numSamples;  

30.     rbf.weight = rbf.weight - rbf.alpha .* delta_weight ./ numSamples;  

31.     rbf.center = rbf.center - rbf.alpha .* delta_center ./ numSamples;  

32.     rbf.delta = rbf.delta - rbf.alpha .* delta_delta ./ numSamples;  

33. end  

 

 

green.m

 

[plain] view plain copy

1. function greenValue = green(x, c, delta)  

2.     greenValue = exp(-1.0 * sum((x - c).^2) / (2 * delta^2));  

3. end  

 

 

五、代碼測試

      首先，我測試了一維的輸入，需要擬合的函數很簡單，就是y=2x。

train_x = [1 2 3 4 5 6 7 8];

train_y = 2 * train_x;

     所以期待的輸出就是：

2    4     6     8   10    12    14   16

     我代碼訓練迭代200次后的網絡輸出是：

2.0042   4.0239    5.9250    8.0214  10.0692   11.9351   14.0179  15.9958

      Matlab的newrb的輸出是：

2.0000   4.0000    6.0000    8.0000  10.0000   12.0000   14.0000  16.0000

      可以看到，Matlab的是完美擬合啊。我的那個還是均方誤差還是挺大的。

      然后，我測試了高維的輸入，訓練樣本是通過Matlab的rand(5, 10)來得到的，它生成的是5行10列[0 1]之間的隨機數。也就是說我們的樣本是10個，每個樣本的維度是5維。我們測試的也是很簡單的函數y=2x。結果如下：

 

 

      關於這個結果，我也不說什么了。期待大家發現代碼里面錯誤的地方，然后告知下，非常感謝。

 

以上均為拷貝過來的，博主提供的matlab代碼還可以改善的，但是不失為一篇學習的好文。