function psobp
% BP neural network trained by PSO algorithm
clc
clear all
% generate training samples according to eutrophication evaluation standard
% of surface water 根据地表水富营养化的标准产生训练样本
st=cputime;
AllSamIn = [
70.00 60.00 60.00 60.00 70.00
70.00 65.00 65.00 65.00 70.00
70.00 75.00 70.00 70.00 70.00
70.00 75.00 70.00 65.00 70.00
70.00 75.00 65.00 65.00 68.00
70.00 75.00 65.00 65.00 68.00
70.00 75.00 65.00 65.00 68.00
70.00 75.00 65.00 65.00 68.00
70.00 75.00 65.00 65.00 70.00
70.00 65.00 65.00 60.00 68.00
70.00 65.00 65.00 60.00 68.00
70.00 65.00 65.00 60.00 68.00
70.00 65.00 65.00 60.00 68.00
70.00 65.00 65.00 60.00 68.00
70.00 65.00 65.00 60.00 68.00
70.00 65.00 65.00 60.00 68.00
70.00 65.00 65.00 60.00 68.00
70.00 65.00 65.00 60.00 68.00
70.00 65.00 65.00 60.00 68.00
70.00 65.00 65.00 60.00 68.00
70.00 65.00 65.00 60.00 68.00
70.00 65.00 65.00 60.00 68.00
70.00 65.00 65.00 60.00 68.00
70.00 65.00 65.00 60.00 68.00
70.00 75.00 65.00 65.00 68.00
70.00 75.00 65.00 65.00 68.00
70.00 75.00 65.00 65.00 68.00
70.00 75.00 65.00 65.00 68.00
135.00 120.00 135.00 100.00 105.00
135.00 120.00 135.00 100.00 105.00
135.00 120.00 135.00 100.00 105.00
135.00 130.00 135.00 100.00 105.00
135.00 135.00 135.00 90.00 105.00
135.00 120.00 135.00 105.00 105.00
130.00 120.00 85.00 85.00 90.00
110.00 120.00 80.00 80.00 90.00
120.00 120.00 80.00 80.00 90.00
120.00 120.00 80.00 80.00 90.00
120.00 100.00 80.00 80.00 90.00
135.00 135.00 120.00 105.00 105.00
135.00 130.00 130.00 105.00 105.00
135.00 130.00 130.00 105.00 105.00
135.00 135.00 135.00 105.00 105.00
135.00 135.00 135.00 105.00 105.00
110.00 110.00 110.00 85.00 85.00
110.00 110.00 110.00 85.00 85.00
110.00 110.00 110.00 85.00 85.00
110.00 110.00 110.00 85.00 85.00
100.00 100.00 80.00 85.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
100.00 85.00 80.00 80.00 85.00
100.00 85.00 80.00 80.00 85.00
100.00 100.00 80.00 80.00 85.00
135.00 135.00 135.00 135.00 135.00
135.00 135.00 135.00 135.00 135.00
135.00 135.00 135.00 135.00 135.00
135.00 135.00 135.00 135.00 135.00
135.00 135.00 135.00 135.00 135.00
135.00 135.00 135.00 135.00 135.00
135.00 135.00 135.00 135.00 135.00
135.00 135.00 135.00 135.00 135.00
125.00 135.00 135.00 135.00 135.00
125.00 130.00 130.00 130.00 135.00
110.00 130.00 130.00 130.00 130.00
125.00 130.00 130.00 130.00 135.00
125.00 120.00 115.00 115.00 120.00
125.00 120.00 115.00 115.00 120.00
125.00 120.00 115.00 115.00 120.00
125.00 120.00 115.00 115.00 100.00
70.00 75.00 60.00 60.00 66.00
70.00 75.00 60.00 60.00 66.00
70.00 75.00 60.00 60.00 66.00
70.00 75.00 60.00 60.00 66.00
]';
% 产生噪声,加给输出
% NoiseVar=0.1;
% rand('state',sum(100*clock));每次重启产生不同的随机数
AllSamOut=[
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
135.00
135.00
135.00
135.00
135.00
135.00
130.00
110.00
120.00
120.00
120.00
135.00
135.00
135.00
135.00
135.00
110.00
110.00
110.00
110.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
100.00
135.00
135.00
135.00
135.00
135.00
135.00
135.00
135.00
125.00
125.00
110.00
125.00
125.00
125.00
125.00
125.00
70.00
70.00
70.00
70.00
]'; %+NoiseVar*randn(1,1000);
% 预处理,压缩到[-1,+1]
global minAllSamOut;
global maxAllSamOut;
[AllSamInn,minAllSamIn,maxAllSamIn,AllSamOutn,minAllSamOut,maxAllSamOut] = premnmx(AllSamIn,AllSamOut);%将所有的输入输出值归一化
% 从总样本中抽取10%的样本作为测试样本,其余作为训练样本
TestSamIn=AllSamIn;
RealTestSamOut=[ 70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
70.00
68.00
68.00
70.00
70.00
70.00
70.00
70.00
70.00
68.00
70.00
70.00
68.00
68.00
70.00
70.00
70.00
70.00
70.00
133.00
133.00
133.00
135.00
135.00
135.00
133.00
110.00
118.00
118.00
118.00
130.00
135.00
135.00
135.00
135.00
110.00
110.00
110.00
110.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
98.00
135.00
135.00
135.00
135.00
135.00
135.00
135.00
135.00
135.00
135.00
135.00
135.00
135.00
135.00
130.00
123.00
68.00
68.00
68.00
68.00
]';
[TestSamInn,minTestSamIn,maxTestSamIn,RealTestSamOutn,minRealTestSamOut,maxRealTestSamOut] = premnmx(TestSamIn,RealTestSamOut);%对测试样本输入和测试样本输出归一化
% 训练样本
TrainSamIn=AllSamIn;%所有输入为训练样本输入
TrainSamOut=AllSamOut;%所有输出为训练样本输出
[TrainSamInn,minTrainSamIn,maxTrainSamIn,TrainSamOutn,minTrainSamOut,maxTrainSamOut] = premnmx(TrainSamIn,TrainSamOut);%将训练样本进行归一化
% Evaluation Sample
%EvaSamIn=...
% [6.5400 9.4000 8.2200 5.9100 6.1300 7.8400 8.2600 8.7900 7.2900 5.9300 4.6900 3.5100
% 0.1270 0.0390 0.0820 0.0860 0.1230 0.1370 0.0190 0.0640 0.0880 0.0500 0.0860 0.0750
% 1.2840 1.3580 1.3200 0.4500 2.2200 1.4700 1.0900 1.6600 0.8400 0.6600 0.5500 0.5200
% 3.3400 4.7200 5.0300 6.0000 7.9600 4.6600 4.0000 4.1500 4.0000 3.8100 1.6300 3.4500
%0.3500 0.4500 0.3800 0.4500 0.3000 0.4500 0.4000 0.4000 0.3000 0.3000 0.4000 0.3500];
%EvaSamInn=tramnmx(EvaSamIn,minAllSamIn,maxAllSamIn); % 预处理
%**********************************************************
% training set 训练设置
global Ptrain;
Ptrain = TrainSamInn;%Ptrain 为归一化后的训练输入
global Ttrain;
Ttrain = TrainSamOutn;%Ttrain 为归一化后的训练输出
% testing set 测试设置
Ptest = TestSamInn; %Ptest为归一化后的测试样本输入
%**********************************************************
% Initialize BPN parameters 初始化BPN因子
global indim;
indim=5;
global hiddennum;
hiddennum=3;
global outdim;
outdim=1;
%**********************************************************
% Initialize PSO parameters 初始化PSO因子
vmax=0.5; % Maximum velocity
minerr=0.0001; % Minimum error
wmax=0.90;
wmin=0.30;
global itmax; %Maximum iteration number
itmax=100;
c1=2;
c2=1.8;
%cf=c1+c2;
for iter=1:itmax
W(iter)=wmax-((wmax-wmin)/itmax)*(itmax-iter); % weight declining linearly
end
% particles are initialized between (a,b) randomly
a=-1;
b=1;
%Between (m,n), (which can also be started from zero)
m=-1;
n=1;
global N; % number of particles
N=40;
global D; % length of par