## 3 The LBOP-ARIMA Model
### The Savitzky-Golay Filtering  
Compared with the standard ARIMA modeling process, the Savitzky-Golay smoothing filtering is proposed to be incorporated in the ARIMA model to take care of the data volatility problem. The main principle is to use a \(k^{\text{th}}\) order polynomial to fit the data points collected within a window of customized length. The S-G filtering is essentially a weighted average algorithm combined with a sliding window. The coefficient is obtained by fitting the data points in the sliding window with the least square method of a given high-order polynomial.  
The multi-order weighting coefficient calculation is given here. The window length is set to be \(n=2m+1\), and \(m\) is the integer constant (\(m\geq 1\)). The sampling time is defined as \(x=(x_{1},x_{2},\ldots,x_{n})\), where \((x_{1},x_{2},\ldots,x_{n})=(-m,-m+1,\ldots,0,\ldots,m)\). The objective of the \((k-1)^{\text{th}}\) degree polynomial is to fit the data points in the window. The fitting function is shown in Equation (5).  
\[y=\partial_{0}+\partial_{1}x+\partial_{2}x^{2}+\ldots+\partial_{k-1}x^{(k-1)} \tag{5}\]  
where \(\partial_{k-1}\) is the \((k-1)^{\text{th}}\) coefficient and \(y\) represents the time series value at time \(x\).  
According to Equation (5), each measurement point is fitted to obtain \(n\)-related equations. To generate the solutions for the equations, it is necessary to set \(n\) as no less than \(k\). Usually \(n>k\) is selected, and the fitting parameter \(A\) is calculated according to the least square method, \(A=[\partial_{0},\partial_{1},\ldots,\partial_{k-1}]^{T}\),  
\[\begin{bmatrix}y_{-m}\\ y_{-m+1}\\ \vdots\\ y_{m}\end{bmatrix}=\begin{bmatrix}1&-m&\cdots&(-m)^{k-1}\\ 1&-m+1&\cdots&(-m+1)^{k-1}\\ \vdots&\vdots&\vdots&\vdots\\ 1&m&\cdots&m^{k-1}\end{bmatrix}\begin{bmatrix}\partial_{0}\\ \partial_{1}\\ \vdots\\ \partial_{k-1}\end{bmatrix}+\begin{bmatrix}e_{-m+1}\\ e_{-m+1}\\ \vdots\\ e_{m}\end{bmatrix} \tag{6}\]  
where \(e_{m}\) represents the constant terms in the linear equations.  
Figure 4: The ARIMA modeling step.  
The least square solution \(\widehat{A}\) of \(A\) is:  
\[\widetilde{A}=\left(X^{T}\cdot X\right)^{-1}\cdot X^{T}\cdot Y \tag{7}\]  
\[\text{where }X=\begin{bmatrix}1&-m&\cdots&(-m)^{k-1}\\ 1&-m+1&\cdots&(-m+1)^{k-1}\\ \vdots&\vdots&\vdots&\vdots\\ 1&m&\cdots&m^{k-1}\end{bmatrix},Y=[y_{-m},\ldots,y_{m}]^{T}.\]  
The model predicted or filtered value \(\widetilde{Y}\) of \(Y\) is:  
\[\widetilde{Y}=X\cdot\widetilde{A}=X\cdot\left(X^{T}\cdot X\right)^{-1}\cdot X ^{T}\cdot Y=B\cdot Y \tag{8}\]  
\[B=X\cdot\left(X^{T}\cdot X\right)^{-1}\cdot X^{T} \tag{9}\]  
To run verifications, here the collected 12-h data set of network traffic status is trained and the time resolution is 10 s. The length of the sliding window is set when \(m=5\). The polynomial order is set as 3. The amount of downlink data in the specified window is calculated to obtain the sequence data of the bandwidth occupancy rate of link \(l_{j}\). The original sequence is fitted and reconstructed. The above S-G smoothing filtering method is introduced, and the effect after processing is shown in Figure 5. It can be seen that the proposed strategy is effective at reducing the data fluctuation and thus improving the data stability.