# Short-Term Wind Speed Forecast Based on B-Spline Neural Network Optimized by PSO

Authors: Zhongqiang Wu; Wenjing Jia; Liru Zhao; Changhan Wu
Journal: Mathematical Problems in Engineering
 (2015)

---

## Abstract

Considering the randomness and volatility of wind, a method based on B-spline neural network optimized by particle swarm optimization is proposed to predict the short-term wind speed. The B-spline neural network can change the division of input space and the definition of basis function flexibly. For any input, only a few outputs of hidden layers are nonzero, the outputs are simple, and the convergence speed is fast, but it is easy to fall into local minimum. The traditional method to divide the input space is thoughtless and it will influence the final prediction accuracy. Particle swarm optimization is adopted to solve the problem by optimizing the nodes. Simulated results show that it has higher prediction accuracy than traditional B-spline neural network and BP neural network.


---

## Body

## 1. Introduction

Wind power is a kind of clean, free, and renewable natural resource. In wind power generation systems, the randomness and volatility of wind will affect the quality of power and the reliability of power system. Improving the prediction accuracy of short-term wind speed is of great significance for the real-time scheduling of power, reliability improving of power supply, and cost lowering of the wind power generation [1–3].The common methods of wind speed forecasting include time series method, support vector machine (SVM) method, and wavelet analysis method. Time series method is adopted in [4, 5], making full use of the sample data and the prediction model is simple which can reduce the calculation time, but it does not take the correlation among the sample data into account, and the estimation of model orders is quite complicated. In [6, 7] SVM is adopted, SVM can make indivisible vector in low dimension space map into divided function in high dimension space by choosing appropriate kernel function, the generalizing ability is improved and the local minimization problem is solved, and it has higher prediction accuracy than time series method. The multiresolution analysis character of wavelet transform and the generalizing ability of SVM are combined to realize the short-term wind speed forecasting in [8]. The method can change the nonlinear phenomenon very well and improve the prediction accuracy, but the design process is rather complex.In view of the various defects of prediction methods above in practical application, a prediction method based on B-spline neural network optimized by particle swarm optimization (PSO-BSNN) is proposed. Firstly, the wind speed data are used to calculate correlation integral function [9, 10] and prove that the discrete time series is chaotic. By phase space reconstruction [11, 12], the high dimension chaotic attractors are recovered. BSNN can change the division of input space and the definition of basis function flexibly [13, 14]. For any input, only a few outputs of the hidden layers are nonzero, so the outputs are simple and convergence speed is fast, but traditional method which divides the input space into linear space is thoughtless and will influence the final prediction accuracy. PSO is an intelligent search method. It has strong global search ability and is also easy to be realized [15–20]. It is used to optimize the nodes of BSNN. Simulated results show that PSO-BSNN has higher prediction accuracy than BSNN and BPNN.


## 2. Phase-Space Reconstruction

A series of continuous sample data are collected from a field of wind farm in Colorado State. Set the sample interval of wind speed as 10 min; there will be 6 sample data in each hour. Let the 6 data be averaged; then the mean value of each hour is got. The discrete sample data of 10 days are recorded asx
-
=
x
-
1
,
x
-
2
,
…
,
x
-
n
n
=
240. For better usage of nonlinear nature of neural network transfer function, the sample data are normalized into the interval [
0,1
], and the normalization formula can be expressed as follows:(1)
x
i
=
x
-
i
-
x
-
min
⁡
x
-
max
⁡
-
x
-
min
⁡
,
0
≤
i
≤
n
,where x
-
i is sample data, x
-
min
⁡ is minimum of all sample data, x
-
max
⁡ is the maximum of all sample data, x
i is normalization data in the interval [
0,1
], and n is the length of time sequence. The normalized wind speed v is displayed in Figure 1.Figure 1
Wind speed time series.


### 2.1. Takens Theorem

Takens has proved that an appropriate embedding dimensionm can be found which satisfies m
≥
2
d
+
1, and d is the dimension of dynamic system. The principle of the dynamic system can be recovered by phase-space reconstruction. If τ is the delay time and m is embedding dimension, then x
=
x
1
,
x
2
,
…
,
x
n can be described as X
i, after phase-space reconstruction. X
i is phase point of the restructured time series, which can be described as follows:(2)
X
i
=
x
i
,
x
i
+
τ
,
x
i
+
2
τ
,
…
,
x
i
+
m
-
1
τ
∈
R
m
,
i
=
1,2
,
…
,
M
,where M
=
n
-
m
-
1
τ is the number of phase point.


### 2.2. Parameter Selection

The selection of embedding dimensionm and delay time τ is very important for the time series x to be reconstructed, which will influence the final prediction accuracy. The saturation correlation dimension method is given to calculate embedding dimension m, and the correlation integral function can be stated as follows:(3)
C
r
=
1
M
2
∑
i
,
j
=
1
M
θ
r
-
a
i
,
j
,
r
>
0
,where r is searching radius, a
i
,
j
=
X
i
-
X
j, and θ
· is Heaviside unit function:(4)
θ
r
-
a
i
,
j
=
0
,
r
-
a
i
,
j
≤
0
,
1
,
r
-
a
i
,
j
>
0
.For a proper scope ofr, the attractor dimension is e
=
ln
⁡
C
r
/
ln
⁡
r; the relation curve between ln
⁡
r and ln
⁡
C
r is obtained as in Figure 2. The relation curve between m and e is obtained as in Figure 3.Figure 2
The relation curve betweenln
⁡
r and ln
⁡
C
r.Figure 3
The relation curve betweene and m.It can be seen from Figure3 that when m
=
5, the attractor dimension e reaches the maximum, so m
=
5 is the saturation correlation dimension. The logarithm value of correlation integral function is less than zero; namely, C
r
∈
0,1 and then the change trend of C
r is decay exponentially, which indicates that the time series x is chaotic rather than noise series.Mutual information method is used to determine the delay timeτ. The information analysis of time series X
i can be defined as the entropy function:(5)
H
X
i
=
-
∑
l
=
0
m
-
1
P
x
i
+
l
τ
log
⁡
P
x
i
+
l
τ
,where H
· is entropy function and P
· is probability density function.IfX
i and X
j are discrete random variables, its joint entropy function can be expressed as follows:(6)
H
X
i
,
X
j
=
-
∑
l
=
0
m
-
1
∑
μ
=
0
m
-
1
P
x
i
+
l
τ
,
x
j
+
μ
τ
log
⁡
P
x
i
+
l
τ
,
x
j
+
μ
τ
.Mutual information function is(7)
I
τ
=
H
X
i
+
H
X
j
-
H
X
i
,
X
j
.The relation curve betweenτ and I
τ is displayed in Figure 4.Figure 4
The relation curve betweenτ and I
τ.Mutual information functionI
τ signifies the correlation between x
t and x
t
+
τ. If I
τ
=
0, they are completely irrelevant, and when I
τ reaches the minimum for the first time, the value τ is the delay time required. It can be seen from Figure 4, and τ
=
4.


## 2.1. Takens Theorem

Takens has proved that an appropriate embedding dimensionm can be found which satisfies m
≥
2
d
+
1, and d is the dimension of dynamic system. The principle of the dynamic system can be recovered by phase-space reconstruction. If τ is the delay time and m is embedding dimension, then x
=
x
1
,
x
2
,
…
,
x
n can be described as X
i, after phase-space reconstruction. X
i is phase point of the restructured time series, which can be described as follows:(2)
X
i
=
x
i
,
x
i
+
τ
,
x
i
+
2
τ
,
…
,
x
i
+
m
-
1
τ
∈
R
m
,
i
=
1,2
,
…
,
M
,where M
=
n
-
m
-
1
τ is the number of phase point.


## 2.2. Parameter Selection

The selection of embedding dimensionm and delay time τ is very important for the time series x to be reconstructed, which will influence the final prediction accuracy. The saturation correlation dimension method is given to calculate embedding dimension m, and the correlation integral function can be stated as follows:(3)
C
r
=
1
M
2
∑
i
,
j
=
1
M
θ
r
-
a
i
,
j
,
r
>
0
,where r is searching radius, a
i
,
j
=
X
i
-
X
j, and θ
· is Heaviside unit function:(4)
θ
r
-
a
i
,
j
=
0
,
r
-
a
i
,
j
≤
0
,
1
,
r
-
a
i
,
j
>
0
.For a proper scope ofr, the attractor dimension is e
=
ln
⁡
C
r
/
ln
⁡
r; the relation curve between ln
⁡
r and ln
⁡
C
r is obtained as in Figure 2. The relation curve between m and e is obtained as in Figure 3.Figure 2
The relation curve betweenln
⁡
r and ln
⁡
C
r.Figure 3
The relation curve betweene and m.It can be seen from Figure3 that when m
=
5, the attractor dimension e reaches the maximum, so m
=
5 is the saturation correlation dimension. The logarithm value of correlation integral function is less than zero; namely, C
r
∈
0,1 and then the change trend of C
r is decay exponentially, which indicates that the time series x is chaotic rather than noise series.Mutual information method is used to determine the delay timeτ. The information analysis of time series X
i can be defined as the entropy function:(5)
H
X
i
=
-
∑
l
=
0
m
-
1
P
x
i
+
l
τ
log
⁡
P
x
i
+
l
τ
,where H
· is entropy function and P
· is probability density function.IfX
i and X
j are discrete random variables, its joint entropy function can be expressed as follows:(6)
H
X
i
,
X
j
=
-
∑
l
=
0
m
-
1
∑
μ
=
0
m
-
1
P
x
i
+
l
τ
,
x
j
+
μ
τ
log
⁡
P
x
i
+
l
τ
,
x
j
+
μ
τ
.Mutual information function is(7)
I
τ
=
H
X
i
+
H
X
j
-
H
X
i
,
X
j
.The relation curve betweenτ and I
τ is displayed in Figure 4.Figure 4
The relation curve betweenτ and I
τ.Mutual information functionI
τ signifies the correlation between x
t and x
t
+
τ. If I
τ
=
0, they are completely irrelevant, and when I
τ reaches the minimum for the first time, the value τ is the delay time required. It can be seen from Figure 4, and τ
=
4.


## 3. B-Spline Neural Network

BSNN is designed based on spline interpolation theory and it can change the division of input space and definition of the basis function flexibly. For any input, only a few outputs of the hidden layers are nonzero so that the outputs are simple and convergence speed is faster. It can be used to realize the nonlinear function approximation with high precision.


### 3.1. Basis Function of BSNN

According to (2), BSNN prediction model is structured as shown in Figure 5, where X
=
x
0
,
…
,
x
l
,
…
,
x
m
-
1
∈
R
m is the input vector and y is the output of the network and corresponds to the mean value of wind speed in the next hour.Figure 5
Prediction model of BSNN.(
1) The Division of Input Space. The input space has to be divided before the basis function being defined. U
l is a finite interval for x
l:(8)
U
l
=
x
l
:
x
l
min
⁡
≤
x
l
≤
x
l
max
⁡
,
0
≤
l
≤
m
-
1
.Then intervalU
l can be divided as follows:(9)
x
l
min
⁡
<
λ
l
,
1
≤
⋯
λ
l
,
ρ
⋯
≤
λ
l
,
q
l
<
x
l
max
⁡
,where λ
l
,
ρ is the ρth internal node, and the external node can be defined as follows:(10)
⋯
λ
l
,
-
1
≤
λ
l
,
0
≤
x
l
min
⁡
,
x
l
max
⁡
≤
λ
l
,
q
l
+
1
≤
λ
l
,
q
l
+
2
⋯
.The intervalU
l can be divided into q
l
+
1 subinterval by all of the nodes as follows:(11)
U
l
,
ρ
=
x
l
:
x
l
∈
λ
l
,
ρ
,
λ
l
,
ρ
+
1
,
ρ
=
0,1
,
…
,
q
l
-
1
,
x
l
∈
λ
l
,
ρ
,
λ
l
,
ρ
+
1
,
ρ
=
q
l
.(
2) The Multivariable Basis Function of BSNN. Let α
x
l indicate k-order basis function for x
l, its ρth basis function is N
ρ
k
x
l which is defined by λ
l
,
(
ρ
-
k
η
)
,
λ
l
,
(
ρ
-
k
η
+
1
)
,
…
,
λ
l
,
ρ, η is expansion factor, and the recurrence relations can be expressed as follows:(12)
N
ρ
k
x
l
=
x
l
-
λ
l
,
ρ
-
k
η
λ
l
,
ρ
-
η
-
λ
l
,
ρ
-
k
η
N
ρ
-
η
k
-
1
x
l
+
λ
l
,
ρ
-
x
l
λ
l
,
ρ
-
λ
l
,
ρ
-
k
-
1
η
N
ρ
k
-
1
x
l
,
N
ρ
1
x
l
=
0
,
x
l
∈
λ
l
,
ρ
-
η
,
λ
l
,
ρ
,
1
,
x
l
∉
λ
l
,
ρ
-
η
,
λ
l
,
ρ
.If the number of internal nodes is 2 and expansion factorη
=
1, the basis function graph of 3 degrees BSNN for single variable x
l is shown in Figure 6.Figure 6
The basis function of BSNN.The number of basis function forX is(13)
q
=
∏
l
=
0
m
-
1
q
l
.The basis function ofX is got as the result that single variable basis functions combine with each other:(14)
α
0,1
,
…
,
m
-
1
X
=
∏
l
=
0
m
-
1
α
x
l
=
∏
l
=
0
m
-
1
∏
ρ
=
1
q
l
N
ρ
k
x
l
.


### 3.2. The Mapping of BSNN

The network from input to output can be divided into two parts.(1)
X
→
α
X.
In this section, it realizes the nonlinear mapping from input layer to hidden layer and setα
X
=
α
1
X
,
…
,
α
r
X
,
…
,
α
q
X
T
∈
R
q, where α
r
X is output value of the rth basis function for X. According to the positivity and tightness of basis function, the number of output nonzero values is b
=
(
k
η
)
m, α
X belongs to the interval [
0,1
], and the other output values are zero. If we let η
=
1, b
=
k
m is got.
(2)
α
X
→
y.
In this section, it realizes linear mapping from hidden layer to output layer:(15)
y
=
ω
t
α
X
,where ω
(
t
)
=
ω
0
(
t
)
,
…
,
ω
r
(
t
)
,
…
,
ω
q
(
t
) is the weight value of the network. The weight value can be updated as follows:(16)
ω
t
+
1
=
ω
t
+
β
y
d
-
y
α
T
X
α
T
X
α
X
,where t is the number of iteration, y
d is the expected value, y is the actual value, and β
0
<
β
<
2 is learning rate and ensures that the iterative algorithm can converge. The verification is as follows.Setε
t
=
y
d
-
y
=
y
d
-
ω
t
α
X; then(17)
Δ
ε
t
=
ε
t
+
1
-
ε
t
=
-
ω
t
+
1
-
ω
t
α
X
=
-
Δ
ω
t
α
X
.Set the Lyapunov function asV
t
=
ε
2
t; then(18)
Δ
V
t
=
V
t
+
1
-
V
t
=
ε
2
t
+
1
-
ε
2
t
=
ε
t
+
1
-
ε
t
2
ε
t
+
ε
t
+
1
-
ε
t
=
Δ
ε
t
2
ε
t
+
Δ
ε
t
=
-
Δ
ω
t
α
X
2
ε
t
-
Δ
ω
t
α
X
=
-
β
ε
t
2
ε
t
+
β
ε
t
=
-
β
2
-
β
ε
2
t
.Since0
<
β
<
2, then Δ
V
t
<
0, and the parameters can converge.


## 3.1. Basis Function of BSNN

According to (2), BSNN prediction model is structured as shown in Figure 5, where X
=
x
0
,
…
,
x
l
,
…
,
x
m
-
1
∈
R
m is the input vector and y is the output of the network and corresponds to the mean value of wind speed in the next hour.Figure 5
Prediction model of BSNN.(
1) The Division of Input Space. The input space has to be divided before the basis function being defined. U
l is a finite interval for x
l:(8)
U
l
=
x
l
:
x
l
min
⁡
≤
x
l
≤
x
l
max
⁡
,
0
≤
l
≤
m
-
1
.Then intervalU
l can be divided as follows:(9)
x
l
min
⁡
<
λ
l
,
1
≤
⋯
λ
l
,
ρ
⋯
≤
λ
l
,
q
l
<
x
l
max
⁡
,where λ
l
,
ρ is the ρth internal node, and the external node can be defined as follows:(10)
⋯
λ
l
,
-
1
≤
λ
l
,
0
≤
x
l
min
⁡
,
x
l
max
⁡
≤
λ
l
,
q
l
+
1
≤
λ
l
,
q
l
+
2
⋯
.The intervalU
l can be divided into q
l
+
1 subinterval by all of the nodes as follows:(11)
U
l
,
ρ
=
x
l
:
x
l
∈
λ
l
,
ρ
,
λ
l
,
ρ
+
1
,
ρ
=
0,1
,
…
,
q
l
-
1
,
x
l
∈
λ
l
,
ρ
,
λ
l
,
ρ
+
1
,
ρ
=
q
l
.(
2) The Multivariable Basis Function of BSNN. Let α
x
l indicate k-order basis function for x
l, its ρth basis function is N
ρ
k
x
l which is defined by λ
l
,
(
ρ
-
k
η
)
,
λ
l
,
(
ρ
-
k
η
+
1
)
,
…
,
λ
l
,
ρ, η is expansion factor, and the recurrence relations can be expressed as follows:(12)
N
ρ
k
x
l
=
x
l
-
λ
l
,
ρ
-
k
η
λ
l
,
ρ
-
η
-
λ
l
,
ρ
-
k
η
N
ρ
-
η
k
-
1
x
l
+
λ
l
,
ρ
-
x
l
λ
l
,
ρ
-
λ
l
,
ρ
-
k
-
1
η
N
ρ
k
-
1
x
l
,
N
ρ
1
x
l
=
0
,
x
l
∈
λ
l
,
ρ
-
η
,
λ
l
,
ρ
,
1
,
x
l
∉
λ
l
,
ρ
-
η
,
λ
l
,
ρ
.If the number of internal nodes is 2 and expansion factorη
=
1, the basis function graph of 3 degrees BSNN for single variable x
l is shown in Figure 6.Figure 6
The basis function of BSNN.The number of basis function forX is(13)
q
=
∏
l
=
0
m
-
1
q
l
.The basis function ofX is got as the result that single variable basis functions combine with each other:(14)
α
0,1
,
…
,
m
-
1
X
=
∏
l
=
0
m
-
1
α
x
l
=
∏
l
=
0
m
-
1
∏
ρ
=
1
q
l
N
ρ
k
x
l
.


## 3.2. The Mapping of BSNN

The network from input to output can be divided into two parts.(1)
X
→
α
X.
In this section, it realizes the nonlinear mapping from input layer to hidden layer and setα
X
=
α
1
X
,
…
,
α
r
X
,
…
,
α
q
X
T
∈
R
q, where α
r
X is output value of the rth basis function for X. According to the positivity and tightness of basis function, the number of output nonzero values is b
=
(
k
η
)
m, α
X belongs to the interval [
0,1
], and the other output values are zero. If we let η
=
1, b
=
k
m is got.
(2)
α
X
→
y.
In this section, it realizes linear mapping from hidden layer to output layer:(15)
y
=
ω
t
α
X
,where ω
(
t
)
=
ω
0
(
t
)
,
…
,
ω
r
(
t
)
,
…
,
ω
q
(
t
) is the weight value of the network. The weight value can be updated as follows:(16)
ω
t
+
1
=
ω
t
+
β
y
d
-
y
α
T
X
α
T
X
α
X
,where t is the number of iteration, y
d is the expected value, y is the actual value, and β
0
<
β
<
2 is learning rate and ensures that the iterative algorithm can converge. The verification is as follows.Setε
t
=
y
d
-
y
=
y
d
-
ω
t
α
X; then(17)
Δ
ε
t
=
ε
t
+
1
-
ε
t
=
-
ω
t
+
1
-
ω
t
α
X
=
-
Δ
ω
t
α
X
.Set the Lyapunov function asV
t
=
ε
2
t; then(18)
Δ
V
t
=
V
t
+
1
-
V
t
=
ε
2
t
+
1
-
ε
2
t
=
ε
t
+
1
-
ε
t
2
ε
t
+
ε
t
+
1
-
ε
t
=
Δ
ε
t
2
ε
t
+
Δ
ε
t
=
-
Δ
ω
t
α
X
2
ε
t
-
Δ
ω
t
α
X
=
-
β
ε
t
2
ε
t
+
β
ε
t
=
-
β
2
-
β
ε
2
t
.Since0
<
β
<
2, then Δ
V
t
<
0, and the parameters can converge.


## 4. B-Spline Neural Network Optimized by PSO

The internal nodesλ
l
,
ρ
ρ
=
1
,
…
,
q
l of BSNN have a great relationship with the prediction accuracy. The traditional method which divides the input space into linear space is thoughtless and it is difficult to get higher approximation accuracy. PSO is used to optimize the parameters, which can improve the accuracy and avoid falling into local minimum.


### 4.1. The Basic Principle of PSO

PSO stems from the enlightenment of birds foraging behavior, which is mainly used to solve the optimization problem. Each particle is a point with a certain speed in the solution space and it also has the related individual fitness which corresponds to the objective function; the effectiveness of the solution is determined by the fitness function that is defined on the basis of optimized goals. Particles are following the current optimum particles in solution space and finding the optimal solution by iteration.Set the current position asp
s
=
p
s
1
,
…
,
p
s
d
,
…
,
p
s
D, the current speed as v
s
=
v
s
1
,
…
,
v
s
d
,
…
,
v
s
D, and the best position that the particle has experienced as p
best
=
p
best
1
,
…
,
p
best
d
,
…
,
p
best
D in D dimension search space. Set the number of all particles as N; the best position of all the particles experienced is g
best
s
=
g
best
s
1
,
…
,
g
best
d
,
…
,
g
best
D. The particles can update themselves by tracking individual extremum p
best and global extremum g
best during each iteration. In the process of finding the two extremums, the particles update their own speed and position as follows:(19)
v
s
d
t
+
1
=
w
v
s
d
t
+
c
1
r
1
p
best
d
t
-
p
s
d
t
+
c
2
r
2
g
best
d
t
-
p
s
d
t
,
p
s
d
t
+
1
=
p
s
d
t
+
v
s
d
t
+
1
,where c
1
,
c
2 are acceleration constant; v
s
d
t is the speed of the dth particle; p
s
d
t is the position of the dth particle; r
1
,
r
2 are random number and they belong to the interval 0,1; w is weight value and is used to balance global searching and local searching, and w can update itself as follows:(20)
w
=
w
max
⁡
-
t
t
max
⁡
w
max
⁡
-
w
min
⁡
,where w
max
⁡ is the maximum of w and w
min
⁡ is the minimum of w; t is the current number of iteration and t
max
⁡ is the maximum of iteration.In order to prevent the particle far from the search space, the speed of each particle will be limited in a maximumv
max
⁡
v
max
⁡
>
0; if the updated speed exceeds the specified speed, then set the speed of the particle as v
max
⁡, v
max
⁡
=
δ
p
max
⁡, where δ belongs to the interval 0,1 and p
max
⁡ is the maximum of p
s
d.


### 4.2. PSO-BSNN Algorithm Implementation

The traditional method for BSNN to divide the internal nodes of input space is thoughtless. If the internal nodes cannot be determined accurately, the BSNN is easy to fall into local minimum and then the prediction accuracy is influenced. PSO is a kind of intelligent search method; the internal nodes can be optimized by using PSO algorithm and then get the optimal location. The dimension of each particle vector isD, and each vector is composed of internal and external nodes λ
l
,
ρ of BSNN. For the single input variable x
i
l, the number of internal nodes is 2 and the external nodes is 4; then the vector dimension D
=
6. The optimization process is as follows.(1) Initialization: the basic initial parameters include the positionp
s
=
p
s
1
,
…
,
p
s
6
=
λ
i
l
,
-
1
,
λ
i
l
,
0
,
λ
i
l
,
1
,
λ
i
l
,
2
,
λ
i
l
,
3
,
λ
i
l
,
4, speed of each particle v
s, the population of particles N, the largest iteration number t
max
⁡, and the acceleration constants c
1
,
c
2.(2) Fitness Evaluation: the error between expected value and actual value is selected as objective function which is used to measure the effectiveness of optimization. The mean square errore
MSE is adopted as the fitness function:(21)
f
=
e
MSE
=
1
N
∑
i
=
1
N
y
i
d
-
y
i
2
0
≤
i
≤
N
,where N is the number of sample date, y
i
d is the expected value for the ith particle, and y
i is the actual output value for the ith particle.(3) Update the parameters according to formula (17) and limit the speed according to v
max
⁡.(4) Calculate the fitness evaluation value, comparing it with the optimal fitness evaluation value; if better, update the particle location; otherwise do not change.(5) Judge whether the result satisfies the termination condition. The termination condition is whether it reaches the maximum number of iterations or satisfies mean square error. If it does, then it is over, or return to step (2).


## 4.1. The Basic Principle of PSO

PSO stems from the enlightenment of birds foraging behavior, which is mainly used to solve the optimization problem. Each particle is a point with a certain speed in the solution space and it also has the related individual fitness which corresponds to the objective function; the effectiveness of the solution is determined by the fitness function that is defined on the basis of optimized goals. Particles are following the current optimum particles in solution space and finding the optimal solution by iteration.Set the current position asp
s
=
p
s
1
,
…
,
p
s
d
,
…
,
p
s
D, the current speed as v
s
=
v
s
1
,
…
,
v
s
d
,
…
,
v
s
D, and the best position that the particle has experienced as p
best
=
p
best
1
,
…
,
p
best
d
,
…
,
p
best
D in D dimension search space. Set the number of all particles as N; the best position of all the particles experienced is g
best
s
=
g
best
s
1
,
…
,
g
best
d
,
…
,
g
best
D. The particles can update themselves by tracking individual extremum p
best and global extremum g
best during each iteration. In the process of finding the two extremums, the particles update their own speed and position as follows:(19)
v
s
d
t
+
1
=
w
v
s
d
t
+
c
1
r
1
p
best
d
t
-
p
s
d
t
+
c
2
r
2
g
best
d
t
-
p
s
d
t
,
p
s
d
t
+
1
=
p
s
d
t
+
v
s
d
t
+
1
,where c
1
,
c
2 are acceleration constant; v
s
d
t is the speed of the dth particle; p
s
d
t is the position of the dth particle; r
1
,
r
2 are random number and they belong to the interval 0,1; w is weight value and is used to balance global searching and local searching, and w can update itself as follows:(20)
w
=
w
max
⁡
-
t
t
max
⁡
w
max
⁡
-
w
min
⁡
,where w
max
⁡ is the maximum of w and w
min
⁡ is the minimum of w; t is the current number of iteration and t
max
⁡ is the maximum of iteration.In order to prevent the particle far from the search space, the speed of each particle will be limited in a maximumv
max
⁡
v
max
⁡
>
0; if the updated speed exceeds the specified speed, then set the speed of the particle as v
max
⁡, v
max
⁡
=
δ
p
max
⁡, where δ belongs to the interval 0,1 and p
max
⁡ is the maximum of p
s
d.


## 4.2. PSO-BSNN Algorithm Implementation

The traditional method for BSNN to divide the internal nodes of input space is thoughtless. If the internal nodes cannot be determined accurately, the BSNN is easy to fall into local minimum and then the prediction accuracy is influenced. PSO is a kind of intelligent search method; the internal nodes can be optimized by using PSO algorithm and then get the optimal location. The dimension of each particle vector isD, and each vector is composed of internal and external nodes λ
l
,
ρ of BSNN. For the single input variable x
i
l, the number of internal nodes is 2 and the external nodes is 4; then the vector dimension D
=
6. The optimization process is as follows.(1) Initialization: the basic initial parameters include the positionp
s
=
p
s
1
,
…
,
p
s
6
=
λ
i
l
,
-
1
,
λ
i
l
,
0
,
λ
i
l
,
1
,
λ
i
l
,
2
,
λ
i
l
,
3
,
λ
i
l
,
4, speed of each particle v
s, the population of particles N, the largest iteration number t
max
⁡, and the acceleration constants c
1
,
c
2.(2) Fitness Evaluation: the error between expected value and actual value is selected as objective function which is used to measure the effectiveness of optimization. The mean square errore
MSE is adopted as the fitness function:(21)
f
=
e
MSE
=
1
N
∑
i
=
1
N
y
i
d
-
y
i
2
0
≤
i
≤
N
,where N is the number of sample date, y
i
d is the expected value for the ith particle, and y
i is the actual output value for the ith particle.(3) Update the parameters according to formula (17) and limit the speed according to v
max
⁡.(4) Calculate the fitness evaluation value, comparing it with the optimal fitness evaluation value; if better, update the particle location; otherwise do not change.(5) Judge whether the result satisfies the termination condition. The termination condition is whether it reaches the maximum number of iterations or satisfies mean square error. If it does, then it is over, or return to step (2).


## 5. Simulation Results and Analysis

Let the largest iteration numbert
max
⁡
=
100, the acceleration constant c
1
=
c
2
=
3, and the population of particles N
=
30, and let r
1
=
0.3, r
2
=
0.6, and the initial position and speed of each particle are got randomly in the defined interval. Adopting 3 orders B-spline neural network, the number of neural network input is 5. The internal node number of each input variable is 2 and expansion factor is 1; then the nodes number of hidden layer units is 3
5 and the nonzero output nodes number of hidden layer is 2
5. There is only one output unit.Three kinds of prediction error are adopted as criteria for evaluation:(22)
e
RMSE
=
1
n
∑
i
=
1
n
y
i
-
y
i
∗
2
,
e
MAE
=
1
n
∑
i
=
1
n
y
i
-
y
i
∗
,
e
MRE
=
1
n
∑
i
=
1
n
y
i
-
y
i
∗
y
i
,where e
RMSE is root mean square error and it is used to measure the deviation between prediction value and real value, e
MAE is mean absolute error that reflects the error of prediction value deviating from real value, and e
MRE is mean relative error that reflects the reliability of mean absolute error value.According to the wind speed date of 10 days in Figure1 that is provided from the wind form in Colorado State, the first data of 7 days is used as training sample to train neural network and the last date of 3 days is used as test sample to examine the results of prediction. Simulation results based on PSO-BSNN, BSNN, and BPNN are displayed in Figure 7. They all adopt 3-layer network structures.Figure 7
Prediction result of 3 days.The prediction error is listed in Table1.Table 1
Forecasting error of each prediction model.
Prediction model
e
RMSE/(m/s)
e
MAE/(m/s)
e
MRE/%
PSO-BSNN
0.0319
0.0348
0.0826
BSNN
0.0700
0.0664
0.1555
BPNN
0.1060
0.1017
0.2267It can be seen from Table1 that the prediction model based on BSNN has better performance than the prediction model based on BPNN, e
RMSE and e
MAE have fallen by 0.2455 m/s and 0.2614 m/s, respectively, e
MRE has fallen by 14.54%, and the prediction effects have some improvements. Comparing with the prediction model based on BSNN, the prediction model based on PSO-BSNN enhances the prediction performance largely, because the input space of BSNN has been optimized by particle swarm. It gets higher prediction accuracy than the prediction model based on BPNN, e
RMSE and e
MAE have fallen by 0.4198 m/s and 0.4233 m/s, respectively, and e
MRE has fallen by 36.50%.Now BPNN is set 4-layer network structures, and the simulation result is shown in Figure8.Figure 8
The prediction result of 3 days with 4-layer BPNN.The prediction error is listed in Table2.Table 2
Forecasting error of BPNN prediction model.
Prediction model
e
RMSE/(m/s)
e
MAE/(m/s)
e
MRE/%
BPNN
0.0810
0.0697
0.1603It can be seen from Table2 that the simulation result based on 4-layer BPNN structure is very similar to BSNN, but still poorer than the prediction result based on PSO-BSNN.The prediction of 3-day wind speed date belongs to short-term wind speed prediction. Now, the wind speed date of 10 hours from Sinkiang in China is used in prediction. The first data of 7 hours is used as training sample and last wind speed date of 3 hours is used as test sample. It belongs to super short-term wind speed prediction. The sample interval of wind speed is still set as 10 min, there will be 6 sample data in each hour, and the number of sample dates of 3 hours is 18. The simulation result is shown in Figure9.Figure 9
Prediction result of 3 hours.The prediction error is listed in Table3.Table 3
Forecasting error of each prediction model.
Prediction model
e
RMSE/(m/s)
e
MAE/(m/s)
e
MRE/%
PSO-BSNN
0.0419
0.0448
0.0926
BSNN
0.0912
0.1064
0.1555
BPNN
0.1377
0.1410
0.2431It can be seen from Table3 that comparing with BPNN, in the simulation results based on BSNN, e
RMSE and e
MAE have fallen by 0.0465 m/s and 0.0346 m/s, respectively, and e
MRE has fallen by 8.76%. Comparing with BSNN, in the simulation results based on PSO-BSNN, e
RMSE and e
MAE have fallen by 0.0493 m/s and 0.0616 m/s, respectively, and e
MRE has fallen by 6.29%.Very short-term wind speed prediction results based on 4-layer BP neural network structure are shown in Figure10.Figure 10
Prediction result of 3 hours with 4-layer BPNN.The prediction error is listed in Table4.Table 4
Forecasting error of BPNN prediction model.
Prediction model
e
RMSE/(m/s)
e
MAE/(m/s)
e
MRE/%
BPNN
0.1023
0.1152
0.1825It can be seen from Table4 that the super short-term wind speed prediction based on 4-layer BPNN is very similar to BSNN, but still poorer than the prediction result based on PSO-BSNN.According to the simulation results of short-term and super short-term wind speed prediction, the prediction model based on PSO-BSNN always gets better results than BSNN and BPNN.


## 6. Conclusions

A prediction model based on BSNN optimized by PSO is present in this paper. BSNN is designed based on B-spline interpolation theory, it can change the input space division flexibly, and the traditional method to divide the input space into linear space is thoughtless which will influence the final prediction accuracy. PSO is adopted to solve the problem. Simulation results show that the PSO-BSNN has better performance than BSNN and BPNN, and the forecasting accuracy is improved effectively.



---

*Source: 278635-2015-04-28.xml*