• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    A New Model Using Multiple Feature Clustering and Neural Networks for Forecasting Hourly PM2.5 Concentrations, and Its Applications in China

    2020-11-05 10:00:08HuiLiuZhihaoLongZhuDuanHuipengShi
    Engineering 2020年8期

    Hui Liu*, Zhihao Long, Zhu Duan, Huipeng Shi

    Institute of Artificial Intelligence and Robotics (IAIR), Key Laboratory of Traffic Safety on Track of Ministry of Education, School of Traffic and Transportation Engineering,Central South University, Changsha 410075, China

    Keywords:PM2.5 concentrations forecasting PM2.5 concentrations clustering Empirical wavelet transform Multi-step forecasting

    ABSTRACT Particulate matter with an aerodynamic diameter no greater than 2.5 μm(PM2.5)concentration forecasting is desirable for air pollution early warning. This study proposes an improved hybrid model, named multi-feature clustering decomposition (MCD)-echo state network (ESN)-particle swarm optimization(PSO), for multi-step PM2.5 concentration forecasting. The proposed model includes decomposition and optimized forecasting components. In the decomposition component, an MCD method consisting of rough sets attribute reduction (RSAR), k-means clustering (KC), and the empirical wavelet transform(EWT) is proposed for feature selection and data classification. Within the MCD, the RSAR algorithm is adopted to select significant air pollutant variables, which are then clustered by the KC algorithm. The clustered results of the PM2.5 concentration series are decomposed into several sublayers by the EWT algorithm.In the optimized forecasting component,an ESN-based predictor is built for each decomposed sublayer to complete the multi-step forecasting computation. The PSO algorithm is utilized to optimize the initial parameters of the ESN-based predictor. Real PM2.5 concentration data from four cities located in different zones in China are utilized to verify the effectiveness of the proposed model.The experimental results indicate that the proposed forecasting model is suitable for the multi-step high-precision forecasting of PM2.5 concentrations and has better performance than the benchmark models.

    1. Introduction

    With the growth of urban industry in developing countries and regions,air pollution has become a difficult issue that is attracting attention all over the world. In recent years, hazy weather has appeared in most regions of China, and air quality has become a national strategic problem. Particulate matter with an aerodynamic diameter no greater than 2.5 μm (PM2.5) contains a large number of toxic and harmful substances [1]. PM2.5is the most common air pollutant and has a negative impact on human health and air quality [2]. Previous studies have shown that PM2.5pollution has a direct impact on the respiratory and cardiovascular systems, and is closely related to the incidence and mortality of lung cancer [3]. In addition, PM2.5has a bad influence on the weather and climate. For example, PM2.5may cause abnormal rainfall and aggravate the greenhouse effect [4-7].

    Given the serious negative impact of PM2.5on people’s lives,increasing attention is being paid to PM2.5concentration forecasting. PM2.5concentration forecasting is considered to be an important and effective method for alleviating the negative effect of PM2.5[8]. This method is also important for applications of urban big data in the development of smart cities [9].

    1.1. Related works

    PM2.5concentration forecasting methods can be divided into four types: physical models, statistical models, artificial intelligence models, and hybrid models.

    Physical models focus on understanding the potentially complex emissions,transport,and conversion processes of meteorological and chemical factors[10].The physical method yields accurate prediction results. However, physical models require sufficient emission information on air pollution [11] and their calculation cost is high [12]. Statistical models overcome the disadvantage of physical methods, as they require simple samples and have a fast calculation speed [13]. However, statistical models do not sufficiently consider the covariance among various influential factors,because they are generally based on limited samples.A single artificial intelligence model can describe the rule of the nonlinear system and has great advantages in dealing with big data [14].However, the disadvantage of such a model lies in the calculation costs of the neural network,which are greater than those of a statistical model. Moreover, the training process of the neural network has a certain volatility, so its output may not be the optimal result [15].

    Considering the limitations of the above methods,hybrid models have been widely used in air pollution prediction.Hybrid models usually combine three parts: data preprocessing, feature selection,and a predictor.Data preprocessing can sort out complex data relationships in the original data and make it more stationary.Feature selection can improve the input data structure and reduce the difficulty of model training caused by a too-high dimension.Hybrid models can integrate the advantages of each algorithm to achieve better model performance. Many related works have shown that hybrid models tend to have better predictive performance [16-20]. Table 1 [16-28] lists some cutting-edge research on hybrid PM2.5concentration forecasting to better illustrate the application of hybrid methods in PM2.5concentration forecasting.

    Feature selection is rarely used in the current hybrid PM2.5concentration forecasting models listed in Table 1. However, if the input of a PM2.5concentration forecasting model includes many features, such as PM2.5, PM10, sulfur dioxide (SO2), and ozone(O3),it may cause difficulties in the PM2.5concentration forecasting model training and increase the training time.This also affects the robustness of the PM2.5concentration forecasting model [29]. At the same time, complex input data may lead to overfitting of the model and may reduce the accuracy of the model [30]. At present,common feature selection algorithms include the principal components analysis (PCA), phase space reconstruction (PSR), and gradient-boosted regression tree (GBRT). However, these methods may be unsuitable for air pollutant concentration sequences because they assume a linear system,which may lead to problems such as not achieving global optimal reduction. The rough sets attribute reduction(RSAR)algorithm,which is based on fuzzy theory, has the advantages of clear stop criteria and no parameters[31].RSAR can obtain the important attribute set of the target attribute through the dependency between different attributes. The RSAR algorithm is a hot research topic [32]. Clustering algorithms are commonly used in data mining and analysis[33].Various clustering methods exist, such as k-means clustering (KC) [34], possibilistic c-means (PCM) [35], cure clustering [36], and so forth.Compared with others, the KC algorithm has the advantages of a simple principle, fast computing speed, and excellent clustering results; thus, the KC algorithm is the most widely used clustering algorithm at present. Combining the RSAR algorithm and the KC algorithm makes it possible to use RSAR to provide reasonable clustering objects for the KC algorithm, which is a valuable research point.

    Decomposition mainly focuses on the wavelet theory method in Table 1. The decomposition algorithm can divide the original data into a series of more stable sublayers according to the different time scales. Compared with empirical mode decomposition(EMD), ensemble empirical mode decomposition (EEMD), and complex empirical mode decomposition (CEMD), the empirical wavelet transform (EWT) algorithm can adaptively divide the Fourier spectrum and select the appropriate wavelet filter banks[37]. The clustering method can also be employed for decomposition in PM2.5concentration forecasting fields. The clustering algorithm can classify the original data according to different air pollutant scenarios.The influence of sample diversity on the model training can be reduced by clustering. However, few studies use aclustering algorithm with a decomposition algorithm in hybrid PM2.5concentration forecasting.

    Table 1 Main studies on PM2.5 concentration forecasting in the past four years.

    The predictors shown in Table 1 are commonly used physical methods, machine learning,and artificial neural networks(ANNs).Although the Copernicus Atmosphere Monitoring Service (CAMS),weather research and forecasting model with chemistry (WRFChem),and Nested Air Quality Prediction Model System(NAQPMS)have accurate prediction results,they require complex data and an understanding of a variety of physical and chemical relationships.Therefore,these methods require a great deal of preparatory work and a high level of professional knowledge. Support vector machines (SVMs), support vector regression (SVR), and least-squares support vector regression (LS-SVR) are very demanding in their choice of parameters,and cannot handle problems with large data.Traditional neural networks such as the backpropagation neural network (BPNN) and evolutionary neural network (ENN) need a great deal of training to build complex neural relationships, and are easy to overfit. The echo state network (ESN) has a unique reservoir structure that consists of recurrently connected units.As a result, the training process of the ESN is simple and effective,which is suitable for nonlinear systems such as PM2.5concentration data [38]. The ESN has been used in other fields such as wind speed prediction [38]. Therefore, application of the ESN model in hybrid PM2.5concentration forecasting is very appropriate.

    1.2. The innovation of this study

    To summarize the references described above, the decomposition-based clustering algorithm, nonlinear fuzzy theory algorithm,and ESN are rarely studied in PM2.5concentration forecasting.This study aims to apply these algorithms for hybrid PM2.5concentration forecasting. The proposed hybrid PM2.5prediction model combines three methods: multi-feature clustering decomposition (MCD), ESN, and particle swarm optimization (PSO). In MCD,the RSAR algorithm is adopted to select significant air pollutant concentrations. Then the KC algorithm is used to divide the original PM2.5concentration data into several groups according to the results of the RSAR algorithm. The clustered results for the PM2.5concentration series are automatically decomposed into several sublayers by the EWT algorithm. An ESN-based predictor is built for every decomposed sublayer in each clustered group to complete the multi-step forecasting computation. The forecasted results from every sublayer are then further integrated to form the final predicting values. The PSO algorithm is utilized to optimize the initial parameters of the ESN-based predictor.The experimental results show that the proposed hybrid model can accurately predict average hourly concentrations of PM2.5.The details of the proposed model are explained in Section 2.

    2. Methodology

    2.1. Framework of the proposed model

    The construction procedure of the hybrid MCD-ESN-PSO model is as follows:

    Part A: MCD

    This part consists of the RSAR,KC,and EWT algorithms.The raw air quality data are filtered by the RSAR algorithm,and the filtered attribute data are clustered by the KC algorithm. The raw data are divided into several clusters by this step.Then the clustered data of each cluster are decomposed into sublayers by the EWT algorithm.Finally,the sublayers are utilized to establish different ESN models.Through the MCD method, the RSAR algorithm and KC algorithm can act on the raw data together to achieve the clustering of features. Through the decomposition processing of the EWT algorithm, the original time series is finally decomposed into more and better sublayers. The details of the RSAR, KC, and EWT algorithms are introduced in Sections 2.2-2.4, respectively.

    Part B: ESN

    The ESN is a basic predictor, which forecasts the decomposed PM2.5concentration data. The ESN is composed of an input layer,reserve pool, and output layer. The main idea of the ESN is to use the reserve pool to simulate a complex dynamic space that can change with the input.Referring to Ref.[38],the updated equation and output state equation of the ESN can be expressed as Eqs. (1)and (2):

    where x is input data from reserve pool to output layer;y is output;t is time; u is input data from input layer to reserve pool; f is the function of ESN; Winrepresents the connection weights of x(t - 1)to x(t); u (t + 1) is the input data; Wbackrepresents the connection weights of the input layer to the reserve pool; and Woutrepresents the connection weights of y(t - 1) to x(t).

    Part C: PSO

    Unlike the traditional ESN model, the ESN model proposed in this paper is combined with the PSO algorithm. In the ESN-PSO algorithm,the relevant parameters of the ESN model,such as input scaling,spectral radius,internal unit number,and connectivity,are optimized by the PSO algorithm.

    Finally, the forecasting results of each sublayer are combined with the corresponding original sublayer. The prediction results of each sublayer are added to obtain the final prediction results.

    2.2. Rough sets attribute reduction

    RSAR can be used to remove useless information while maintaining the quality of the sorting of the existing information [31].The obtained information is referred to as‘‘reducts.”In an information system, a set of objects are described by a set of attributes[31]. A knowledge information system is defined as follows:

    where U is a finite nonempty set of objects; V is a nonempty set of values;A is a finite nonempty set of attributes;and h is an information function that maps an object in U to exactly one value in V.

    In this study, A is a set of all attributes such as A = {PM10, CO,SO2,NO2,O3,PM2.5},and V is their value.f is the dependency function that is used to obtain γ, and γ is the dependence of the set,which is calculated in the process of RSAR.

    The reduct should maintain the quality of sorting(γ),defining a condition attributes set C ?A and an attribute set P ?C ?A.Sometimes, an information table may have more than one reduct; the intersection of all the reducts is called a ‘‘core” of a decisionmaking table and is expressed as core (P); this is the most important attributes set for an information system.

    2.3. k-means clustering

    KC is a simple iterative clustering algorithm,using distance as a similarity index [34]. Its final purpose is to find k clusters in a group of given datasets. The center of each cluster is according to the value of all clusters in which each cluster is described by the clustering center. The process of the KC algorithm is as follows:

    (1) Select the k object in the data space as the initial center;each object represents a cluster center.

    (2)Divide the data objects in the sample into the corresponding classes according to the nearest clustering center,according to the Euclidean distance between them and these cluster centers.

    where xiis the ith sample in the jth cluster;xjis the center of the jth cluster; and D represents the number of attributes of a data object.

    (3) Update the clustering center: taking the mean values of all objects in each cluster as the clustering center, calculating the value of the objective function.

    (4)Judge whether the values of the cluster centers and objective functions are equal. If they are equal, output results; otherwise,return to step (2).

    2.4. Empirical wavelet transform

    In this paper,the EWT algorithm is used for data preprocessing.The EWT,proposed by Gilles[37],is a novel signal-processing technique that builds the wavelets adaptively.The EWT is based on the theoretical framework of wavelet transform but overcomes the shortage of EMD theory and the problem of signal aliasing. The EWT adaptively divides the Fourier spectrum and selects the appropriate wavelet filter banks. The empirical scaling function and the empirical wavelets can be expressed as Eqs.(5)and(6),respectively:

    where n is divided interval;ω is frequency;β is any function in the interval[0,1]that satisfies the derivative of the order,τ is frequency coefficient; β(x) = x4(35 - 84x + 70x2- 20x3) and τ <min n [(ωn+1-ωn)/ (ωn+1+ωn)].

    2.5. Particle swarm optimization

    The PSO algorithm consists of position z,speed v,and the adaptation function. Each particle in the algorithm represents a candidate solution in the solution space. The fitness function is set according to the optimization goal. During the training of the PSO, each particle in the algorithm updates its own position by combining its current movement experience with the movement experience of the neighboring particles. The solution is realized by putting one’s own position close to the target position.The calculation formula is as follows [27]:

    3. Case study

    3.1. Study area

    Related literature studies show that the distribution of PM2.5concentrations ranges widely in China. It is mainly concentrated in North China and Central China [39,40]. In order to ensure the diversity of experimental data, the selected data should include different working conditions such as serious PM2.5pollution and weak PM2.5pollution.Beijing in the North China Plain,Guangzhou in the Pearl River Delta, Changsha in Central China, and Suzhou in the Yangtze River Delta are typical cities that are used to verify the effectiveness of the proposed model.The selected samples are spatially representative and contain PM2.5concentration data under different geographic and climatic environments, which can well verify the effectiveness of the proposed model.

    Monitoring stations record the average concentrations of six kinds of air pollutants(PM2.5, PM10, NO2, SO2, O3, and CO). Fig. 1 shows the selected datasets and related introduction.

    3.2. Data description and partitioning

    Fig. 1. Locations of the air quality monitoring stations. (a) Beijing. Beijing is the capital of China, located at the north end of the North China Plain. It has a typical warm temperate semi-humid continental monsoon climate,hot and rainy in summer,cold and dry in winter,short in spring and autumn.The average annual temperature is 10-12 °C and the average annual rainfall is more than 600 mm. (b) Changsha. Changsha is an important city in the middle reaches of the Yangtze River. It is a subtropical monsoon climate, mild climate, abundant precipitation, hot and rainy at the same time. Its average annual temperature is 17.2 °C and the average annual precipitation is 1361.6 mm.(c)Guangzhou.Guangzhou is located in the southeastern part of China,the northern edge of the Pearl River Delta,and the Pearl River passes through the urban area.It belongs to the tropical monsoon climate with high temperature,rainfall,and low wind speed.(d)Suzhou.Suzhou is located in the southeast of Jiangsu Province and the middle of the Yangtze River Delta. It is a subtropical monsoon marine climate, four distinct seasons, abundant rainfall throughout the year. Group A: models without RSAR-KC, including the ESN, LSTM, ESN-PSO, and EWT-ESN-PSO model.

    Experimental data are collected from four cities: Beijing,Guangzhou, Changsha, and Suzhou. As Shi et al. [41] have indicated, the spatial representation of surface site observation is often 0.5-16 km2, with the most frequent values being around 3 km2.The data of a single monitoring station cannot represent the air quality of the whole city.Therefore,each set of data is the mean value of all air quality monitoring stations in the corresponding city,so that the samples can represent the air quality of the whole city. For convenience, these datasets are named D1 (Beijing), D2(Guangzhou), D3 (Changsha), and D4 (Suzhou). The length of the sample data is set to one year so that the data can cover a complete set of four seasons. In this study,the sample data in 2016 are randomly selected. All experimental data includes one-hour average concentrations of PM2.5,PM10,NO2,SO2,O3,and CO collected from 1 January 2016 to 31 December 2016. All data are retrieved from the website of the China National Environmental Monitoring Center.?? http://www.cnemc.cn/

    Missing value filtering and outlier checking are implemented before data partitioning.Through the sample set analysis,it is concluded that there are 220 missing pieces of data in dataset D1.Dataset D2 is missing 158 pieces of data, dataset D3 is missing 158 pieces of data, and dataset D4 is missing 157 pieces of data.Since the number of missing samples is less than 2.5% of the total sample set, direct elimination will not have a significant impact.Through visual inspection of the outliers in Fig. 1, it is found that the outliers are mostly concentrated from January to March and from October to December. In order to ensure the training effect of the model, outliers are regarded as normal and retained.

    After removing the missing samples, D1 has 8540 samples, D2 has 8602 samples,D3 has 8602 samples,and D4 has 8603 samples.The 4001st-4600th samples of each dataset(other data are utilized for KC in group B)are used to train the models in group A(models without RSAR-KC, including the ESN, long-short term memory network(LSTM),ESN-PSO,and EWT-ESN-PSO model),which only use PM2.5concentrations.The 4601st-5000th samples are the testing set;to ensure the prediction effect,the 4601st-4900th samples are abandoned.All of the experimental data of each station is used in the RSAR and KC to preprocess the data in group B(models with RSAR-KC, including the RSAR-KC-ESN, the MCD-LSTM-PSO, and the RSAR-KC-EWT-ESN-PSO model). To ensure the effectiveness of error evaluating, each cluster is used to train an ESN model and then reconstruct the predicted results for the 4901st-5000th samples.

    In order to study the influence of different sampling processes on model accuracy, the 3001st-4000th (S1) sample points and 6001st-7000th (S2) sample points in D1 are used for comparison experiments. Fig. 2 shows the distribution of datasets S1 and S2.

    In order to further verify the effectiveness of the model,an additional dataset named D4(which contains 8603 samples)is used in the experiments. The dataset D4 selects monthly data from the spring, summer, autumn, and winter for testing; these data are named T1(1000th-1999th samples),T2(3100th-4099th samples),T3 (5000th-5999th samples), and T4 (6000th-6999th samples).They are shown in Fig. 3. Table 2 shows the descriptive statistics of the related PM2.5concentrations.

    3.3. Results and discussion

    3.3.1. Results of RSAR

    The RSAR and KC are used to preprocess the original data. The attribute decision tables of each dataset are established according to the international PM2.5classification system. According to the international PM2.5concentration index classification standard,PM2.5concentration data are classified and discretized. Like the method of classifying PM2.5concentration levels,the concentrations of the other five air pollutants are discretized according to the levels.Table 3 shows the attribute reduction table for this study.By calculating the positive region values of the other five kinds of air pollutant concentrations and PM2.5concentrations, it can be determined that the significance degrees of PM10,NO2,CO,O3,and SO2are 0.0825, 0.0948, 0.0531, 0.2189, and 0.1843, respectively.O3and SO2have great significance and are judged to be the core attributes of the established information decision system.

    It should be noted that if the correlation between the reduction attribute and the decision attribute is too strong, there is no distinction between the two.If the correlation between the reduction attributes and decision attributes is too weak, there is no correlation between them. The reduction attributes in both cases are redundant.Therefore,in order to ensure the diversity of input samples, the selection of reduction attributes needs to comprehensively consider the correlation and independence between the reduction attributes and decision attributes. For this reason, this paper uses covariance to evaluate the relationship between PM2.5concentrations and other pollutant concentrations, as shown in Table 4. The cov(PM2.5, PM10), cov(PM2.5, NO2), cov(PM2.5, CO),and cov(PM2.5, SO2) are positive. The cov(PM2.5, O3) is negative.However, the absolute value of cov(PM2.5, PM10), cov(PM2.5, NO2),and cov(PM2.5, CO) is much larger than that of cov(PM2.5, SO2)and cov(PM2.5,O3).In terms of ensuring the independence of input attributes, the result of the RSAR algorithm can be verified to be effective. In order to avoid difficulty in model training caused by dimensional disaster, these more relevant attributes are selected as the core attributes, and other data with weak correlation become the reduction attributes.

    Fig. 2. PM2.5 concentrations series of S1 and S2.

    Fig. 3. PM2.5 concentrations series of T1-T4.

    Table 2 Descriptive statistics of PM2.5 concentrations.

    Table 3 Attribute reduction table.

    Table 4 Covariance table.

    3.3.2. Results of KC

    After the attribute reduction, the original dataset becomes an N × 3 sample space; three-dimensional (3D) KC is used to divide it into several similar clusters.The sum of the squared errors(SSE)[42]and the silhouette coefficient(SC)[43]are used to choose the best value of k.Because the clustering results of the three datasets are very similar, D1 is used as an example to show the results.

    Fig.4 shows different SSE and SC when choosing different k.The k value ranges from 1 to 15, and the SSE value decreases with the increase of the k value.When k=3,the SC value is the largest,but the SSE value is larger at the same time. Considering SSE and SC together in Fig. 4, the final k value is 7. Finally, the data of station D1 are divided into seven clusters.

    When k=7,the original data D1 are divided into seven groups;the results are shown in Fig. 5. In this figure,(a) shows the results for PM2.5,while the results for SO2and O3are presented separately in(b)and(c).The results of the KC for PM2.5shown in Fig.5(a)are the key part of this paper.Intuitively,the amplitude of cluster C1 is between 0 and 200 μg·m-3,and it fluctuates gently.The amplitude of C2 is between 0 and 55 μg·m-3, and fluctuates violently. However, the short period of C2 is obvious. The amplitude of C3 is between 0 and 400 μg·m-3. The fluctuation of C3 is smooth but the periodicity is not obvious. The amplitude of C4 is between 50 and 150 μg·m-3, and its periodicity and symmetry are good. The amplitude of C5 is between 0 and 200 μg·m-3, but it fluctuates more violently than that of C1. The amplitude of C6 is between 160 and 240 μg·m-3, fluctuates violently, and has strong symmetry.The amplitude of C7 is between 0 and 100 μg·m-3,and the period is obvious but the symmetry is weak. Overall, compared with the original data in Fig. 1, each kind of data becomes more stable after clustering. The peaks and troughs of each type of data show different intervals.

    The above description is only a subjective analysis; in order to draw a more convincing conclusion, descriptive statistics are used to further analyze the clustering results of PM2.5concentration data.Table 5 shows the descriptive statistics of D1 for the different clusters.

    The mean values of the seven groups of data are 71.54, 24.00,285.74,91.47,83.90,177.00,and 34.85 μg·m-3,respectively.Combined with the maximum and minimum values of each group of data, the seven groups of data after clustering are concentrated according to the value size,thereby reducing the fluctuation range of the data within the group.This is consistent with the amplitude distribution of each group of data in Fig. 5.

    The standard deviation (SD) reflects the dispersion degree among individuals in a group. The standard deviation values of the seven groups of data after clustering are 37.02, 14.25, 47.30,20.96, 32.70, 29.81, and 19.42 μg·m-3, respectively, which are all smaller than the 71.00 μg·m-3before clustering. The data of each group after clustering are closer to their average values.This trend is reflected in Fig.5: The symmetry of the upper and lower fluctuations of the data curves of each group is stronger.

    Fig. 4. SSE and SC with different k values.

    The skewness values of the seven groups of data after clustering are 0.70, 0.72, 0.74, 0.21, 1.01, 0.22, and 0.88, respectively, which are all smaller than the 2.01 before clustering.The wave peak symmetry of the data after clustering is stronger;that is,the cycle rule is more obvious. The kurtosis values of the seven groups of data after clustering are 3.23, 2.45, 2.50, 1.98, 4.00, 1.82, and 3.12,respectively, which are all smaller than the 8.64 before clustering.The extreme distribution of data in each group of data after clustering is reduced. In Fig. 5, the fluctuation of each group of data is smooth, and there is no obvious outlier. In other words, Fig. 5 and Table 5 show that the clustered data have the above advantages.

    The MCD-ESN model is used to analyze the series length in each cluster. In order to ensure the validity of the error evaluation, the first 80% of the data in each cluster is selected for model training and the last 20%is used for model prediction performance analysis.Table 6 shows the error evaluation of each cluster.

    Fig.5. (a) Results of KC for PM2.5;(b) results of KC for PM2.5 and SO2; (c) results of KC for PM2.5 and O3.

    Table 5 Descriptive statistics of D1 for different clusters.

    Table 6 Error evaluation of each cluster for D1 of MCD-ESN.

    When the sample size is greater than 1000,the amount of data has little effect on the prediction,as in C1,C3,C5,C6,and C7.However,when the number of samples is less than 1000,the prediction effect of the model is greatly reduced; this shows that the prediction effect of the ESN network is more sensitive to the number of low samples,as in C2 and C4.When the number of samples is small after clustering,the problem can be solved by increasing the number of samples in the original sequence.This will be a problem for future research.

    3.3.3. Forecasting accuracy and analysis

    In this study, six other prediction models are provided as comparison models to investigate the prediction performance of the proposed model. In addition, to investigate the multi-step prediction performance of the proposed model, all the involved models are conducted for the step-1 to step-3 predictions. The proposed model must forget a certain amount of output results due to the characteristics of the ESN algorithm[38].Therefore,the prediction accuracy fluctuates within a certain range.In this paper,this problem is solved by averaging the results of three repeated experiments. This solution does not increase the computing time cost by too much.

    The mean absolute percentage error (MAPE), mean absolute error (MAE), root mean square error (RMSE), standard deviation of error (SDE), Pearson’s correlation coefficient (R), and index of agreement (IA) are utilized to analyze the experimental results of the prediction models; the index values of the above models for D1, D2, and D3 are given in Table 7. As can be seen from Table 7,the three datasets reflect the same model performance. In order to keep the length of the paper within a reasonable range, only D1 is selected for specific analysis. The PM2.5concentration forecasting results for D1 are shown in Fig. 6. The R and IA results of the six prediction models for S1,S2,and T1-T4 are given in Table 8.The MAPE,MAE, RMSE, and SDE results of the six prediction models for S1 and S2 are given in Fig.7.The MAPE,MAE,RMSE,and SDE results of the six prediction models for T1 and T2 are given in Fig.8.The MAPE,MAE,RMSE and SDE results of the six prediction models for T3 and T4 are given in Fig. 9. It should be noted that since the values of R and IA do not belong to the same dimension as the other four evaluation indicators, they are not shown in the form of a graph.

    In Tables 7 and 8 and Figs. 6-9, the proposed model has the smallest error evaluation values, and thus achieves accurate prediction in PM2.5concentration forecasting. Compared with the other six comparison prediction models, the proposed model has better prediction accuracy from step-1 to step-3 predictions. This shows that the methods used in the hybrid model interact positively.

    The prediction accuracy of the ESN-PSO model is better than that of the ESN model.This phenomenon indicates that the optimal parameters selected by the PSO algorithm can help to improve the prediction accuracy of the ESN model. The prediction accuracy ofthe EWT-ESN-PSO model is better than that of the ESN-PSO model, which shows that the model accuracy can be improved by adding the EWT decomposition algorithm. The sequence obtained by the EWT algorithm is more stable and less random.Therefore, it is better to take the decomposed sublayers as the model input to obtain the prediction results. The prediction accuracy of the RSAR-KC-ESN model is better than that of the ESN model, which shows that the model accuracy can be improved by adding the RSAR-KC algorithm.After clustering,there is a larger difference between different groups of data and higher similarity between the same groups of data, which can improve the prediction accuracy of the original model to a certain extent.

    Table 7 Error evaluation for the PM2.5 concentrations of D1, D2, and D3.

    Moreover, the accuracy of each prediction model decreases as the number of steps increases in Table 7 and Figs. 6-9. With the increase of the forecasting step, the error accumulation will become increasingly serious, resulting in decline of the prediction accuracy.

    The city with the best air quality is Changsha(D3),followed by Guangzhou (D2). The air quality of Beijing (D1) is relatively poor.The forecasting accuracy in Table 7 and Fig. 6 is consistent with this order.Moreover,the data in Fig.7 show that samples with different pollution levels in the same area have no effect on model accuracy. The PM2.5concentrations of S1 are smaller than those of S2, but the prediction accuracy of S2 is higher than that of S1.Therefore, it can be concluded that the prediction accuracy of the proposed model is better in cities with better air quality than in cities with heavy pollution.

    In the above analysis,Tables 7 and 8 and Figs.6 and 7 verify the validity of the data forecasts for different cities over the same time period.In order to verify the validity of the prediction for the samecity over different time periods, the experiments shown in Figs. 8 and 9 were carried out. According to the data in Figs. 8 and 9,the proposed model maintained a stable prediction effect with the change of time period. In other words, the data in Figs. 8 and 9 verify the stability and validity of the proposed model over the whole year.

    Table 8 R and IA for the PM2.5 concentrations of S1, S2, and T1-T4.

    Fig. 6. Results of the various ahead-step predictions for the PM2.5 concentrations of D1. (a) Step-1; (b) step-2; (c) step-3.

    Fig. 7. Error evaluation for the PM2.5 concentrations of (a) S1 and (b) S2.

    Fig. 8. Error evaluation for the PM2.5 concentrations of (a) T1 and (b) T2.

    Fig. 9. Error evaluation for the PM2.5 concentrations of (a) T3 and (b) T4.

    In this study, the following computing is implemented in the simulation environments: Intel i5-6500 CPU 3.2 GHz, RAM 8 GB.Table 9 gives the computation time of the comparison models in D1. Because both RSAR-KC and PSO are offline processing, the computation time is not compare with them.

    The computation speed of the ESN is much faster than that of the LSTM, owing to the advantages of the ESN network itself.Because of the existence of the reserve pool, it is only necessary to train the output weight in the training process of the ESN network, which greatly improves the computing speed.

    After adding the EWT decomposition algorithm, the computational speed of the model decreases to a certain extent. Because each decomposition layer needs to be trained and predicted, the computational speed of the original model plays a vital role here,which further reflects the superiority of the ESN.

    The change in prediction steps has little effect on the calculation speed of the model. This may be because the computational capacity of the algorithm model is relatively large.

    4. Conclusions

    This study establishes an improved hybrid ESN forecasting model to predict and analyze the hourly average concentrationsof PM2.5based on the MCD method and the PSO. The proposed hybrid model is compared with several benchmark models to prove its effectiveness. The attribute reduction results show that the concentrations of SO2and O3play an important role in predicting the concentrations of PM2.5. Moreover, research on the influence of relevant meteorological parameters will be carried out in future studies. The clustering results show that PM2.5concentration data become stationary and regular after the clustering processing. These features are useful and conducive for ESN-based deep training. The prediction results show that: ① The MCD method can improve the accuracy of models; ②the proposed hybrid model has better prediction accuracy than other relevant deep learning or single models; ③the proposed hybrid model has achieved good experimental results with PM2.5pollutant concentration data from four cities in China; and ④the proposed hybrid PM2.5forecasting framework can also be applied in other air pollution time series multi-step predictions. The forecasted results can be embedded in relevant early warning systems for urban air pollution management.

    Table 9 Computation time of the comparison models in D1.

    The main contributions of this study can be summarized as follows:

    (1)A novel PM2.5concentrations multi-step prediction model is developed based on the MCD,ESN,and PSO,which yields accurate forecasting performance for hourly average PM2.5concentrations.The multi-step forecasting results can be used for the development of PM2.5pollution warning systems.

    (2)A novel decomposition method in hybrid PM2.5concentration forecasting named MCD is developed. This method integrates the feature extraction into decomposition.Multi-dimensional KC clustering is carried out using the feature extraction results of the RSAR algorithm,which not only guarantees the effectiveness of the clustering results, but also considers the influence of multidimensional features. The EWT algorithm-based KC algorithm is then employed for data preprocessing. The clustering algorithm is used to group the original PM2.5concentrations according to different PM2.5concentration scenarios. Next, combined with the EWT decomposition algorithm,raw PM2.5concentrations data are distinguished according to different characteristics in the timescale.Finally,the optimization function of the decomposition is realized.

    (3) In the proposed hybrid PM2.5concentration forecasting model, the ESN is employed as a predictor. The sparse connection of neurons in the reservoir of the ESN not only improves the convergence of the neural network model, but also enhances the model generalization. This characteristic can reduce the probability of the overfitting problem in the process of model training.Moreover, the ESN has good real-time performance in computing.

    Acknowledgements

    The study is fully supported by the National Natural Science Foundation of China (61873283), the Changsha Science &Technology Project (KQ1707017) and the Innovation Driven Project of the Central South University (2019CX005).

    Compliance with ethics guidelines

    Hui Liu, Zhihao Long, Zhu Duan, and Huipeng Shi declare that they have no conflict of interest or financial conflicts to disclose.

    黑人欧美特级aaaaaa片| av免费在线观看网站| 久久精品亚洲精品国产色婷小说| 热99国产精品久久久久久7| 亚洲午夜理论影院| www日本在线高清视频| 曰老女人黄片| 亚洲av成人一区二区三| 在线十欧美十亚洲十日本专区| 国产精品亚洲一级av第二区| 悠悠久久av| av天堂在线播放| 精品一品国产午夜福利视频| 欧美激情 高清一区二区三区| 老汉色av国产亚洲站长工具| 日韩视频一区二区在线观看| 欧美精品人与动牲交sv欧美| 波多野结衣一区麻豆| 夜夜夜夜夜久久久久| 成人影院久久| 窝窝影院91人妻| 一本大道久久a久久精品| 欧美日韩精品网址| 成年人午夜在线观看视频| 亚洲午夜理论影院| 精品国产国语对白av| 黄色毛片三级朝国网站| 久久青草综合色| 欧美 亚洲 国产 日韩一| 女性被躁到高潮视频| 999久久久精品免费观看国产| 久久人人97超碰香蕉20202| 日韩精品免费视频一区二区三区| 窝窝影院91人妻| 王馨瑶露胸无遮挡在线观看| 曰老女人黄片| 999久久久国产精品视频| 欧美在线黄色| 深夜精品福利| 国产av精品麻豆| 嫁个100分男人电影在线观看| 色婷婷av一区二区三区视频| 亚洲avbb在线观看| 久久青草综合色| 下体分泌物呈黄色| 天天躁夜夜躁狠狠躁躁| 日本欧美视频一区| 久久久久精品人妻al黑| 国产伦人伦偷精品视频| 视频在线观看一区二区三区| 黄色片一级片一级黄色片| 成年女人毛片免费观看观看9 | 美女高潮到喷水免费观看| 最近最新中文字幕大全免费视频| 后天国语完整版免费观看| 成人国语在线视频| 亚洲精品久久午夜乱码| 一本大道久久a久久精品| 超碰97精品在线观看| tube8黄色片| 男女床上黄色一级片免费看| 中文字幕另类日韩欧美亚洲嫩草| 国产黄频视频在线观看| 国产一区二区在线观看av| 国产精品欧美亚洲77777| 精品一区二区三区av网在线观看 | 精品一品国产午夜福利视频| 久久久国产精品麻豆| 最近最新免费中文字幕在线| 日本wwww免费看| 国产亚洲欧美精品永久| 啦啦啦视频在线资源免费观看| 成年版毛片免费区| 黑人巨大精品欧美一区二区蜜桃| 两个人免费观看高清视频| 女人久久www免费人成看片| 久久久久久久精品吃奶| 别揉我奶头~嗯~啊~动态视频| 色尼玛亚洲综合影院| 美女福利国产在线| 9热在线视频观看99| 男人舔女人的私密视频| 精品人妻在线不人妻| 国产一区二区激情短视频| 亚洲国产欧美一区二区综合| 菩萨蛮人人尽说江南好唐韦庄| 成人18禁高潮啪啪吃奶动态图| 女人爽到高潮嗷嗷叫在线视频| 午夜两性在线视频| 亚洲成人免费av在线播放| 国产伦人伦偷精品视频| 国产色视频综合| 国产激情久久老熟女| 久久 成人 亚洲| 少妇 在线观看| 日韩欧美国产一区二区入口| 人人妻人人爽人人添夜夜欢视频| 性色av乱码一区二区三区2| 男人操女人黄网站| 欧美成人免费av一区二区三区 | 一本—道久久a久久精品蜜桃钙片| 91老司机精品| 在线 av 中文字幕| 日本五十路高清| 日韩欧美三级三区| 亚洲国产精品一区二区三区在线| 精品国产乱码久久久久久小说| 欧美中文综合在线视频| netflix在线观看网站| 蜜桃国产av成人99| 考比视频在线观看| 欧美日韩亚洲综合一区二区三区_| 中文欧美无线码| 午夜久久久在线观看| 午夜福利影视在线免费观看| 妹子高潮喷水视频| 久久狼人影院| 俄罗斯特黄特色一大片| 99久久国产精品久久久| 成人18禁高潮啪啪吃奶动态图| 在线 av 中文字幕| 97在线人人人人妻| 男女床上黄色一级片免费看| 色精品久久人妻99蜜桃| 欧美黄色淫秽网站| 啦啦啦中文免费视频观看日本| 一二三四在线观看免费中文在| 亚洲成人免费av在线播放| 国产午夜精品久久久久久| 中国美女看黄片| 91麻豆av在线| 亚洲成人手机| 婷婷成人精品国产| 久热这里只有精品99| 久久人妻福利社区极品人妻图片| 99精品久久久久人妻精品| 中文字幕最新亚洲高清| 激情在线观看视频在线高清 | 在线永久观看黄色视频| 亚洲性夜色夜夜综合| 国产国语露脸激情在线看| 久久午夜综合久久蜜桃| 欧美激情久久久久久爽电影 | cao死你这个sao货| 亚洲国产av影院在线观看| 2018国产大陆天天弄谢| 9191精品国产免费久久| 少妇粗大呻吟视频| 十八禁网站免费在线| 十八禁人妻一区二区| av又黄又爽大尺度在线免费看| 青青草视频在线视频观看| 黄色视频在线播放观看不卡| 无遮挡黄片免费观看| 99在线人妻在线中文字幕 | 国产成人精品无人区| 亚洲国产av新网站| 怎么达到女性高潮| 日韩一区二区三区影片| 久热这里只有精品99| 老鸭窝网址在线观看| 国产成人欧美在线观看 | 国产精品久久久久成人av| 国产在线一区二区三区精| 国产精品国产av在线观看| 天堂中文最新版在线下载| 一进一出抽搐动态| 黄色怎么调成土黄色| 不卡一级毛片| 欧美老熟妇乱子伦牲交| 久久久久国内视频| 大片电影免费在线观看免费| 成人永久免费在线观看视频 | 色婷婷av一区二区三区视频| 最近最新中文字幕大全电影3 | 麻豆国产av国片精品| 久久人人爽av亚洲精品天堂| 国产精品免费视频内射| 午夜福利一区二区在线看| 亚洲av欧美aⅴ国产| 婷婷丁香在线五月| 无限看片的www在线观看| 亚洲少妇的诱惑av| 午夜福利免费观看在线| 69精品国产乱码久久久| 最新美女视频免费是黄的| 一边摸一边做爽爽视频免费| 巨乳人妻的诱惑在线观看| 一夜夜www| 中文字幕色久视频| 亚洲自偷自拍图片 自拍| 国产精品影院久久| 在线天堂中文资源库| 后天国语完整版免费观看| 99精品在免费线老司机午夜| 日韩一区二区三区影片| 欧美激情 高清一区二区三区| 国产精品偷伦视频观看了| www.熟女人妻精品国产| 曰老女人黄片| av欧美777| 日本a在线网址| 女人精品久久久久毛片| 99国产精品99久久久久| 91大片在线观看| av免费在线观看网站| 欧美变态另类bdsm刘玥| 国产精品久久久久久人妻精品电影 | 黄色怎么调成土黄色| 老司机靠b影院| 捣出白浆h1v1| 久久香蕉激情| 久久 成人 亚洲| av线在线观看网站| 免费高清在线观看日韩| 久久久精品国产亚洲av高清涩受| 久久国产亚洲av麻豆专区| 欧美日韩亚洲高清精品| 欧美激情极品国产一区二区三区| 久久久精品国产亚洲av高清涩受| 丁香六月欧美| 亚洲国产成人一精品久久久| 五月开心婷婷网| 日本wwww免费看| av电影中文网址| 在线观看一区二区三区激情| 亚洲av第一区精品v没综合| av网站在线播放免费| 久久精品亚洲熟妇少妇任你| 极品少妇高潮喷水抽搐| 亚洲性夜色夜夜综合| 亚洲国产欧美在线一区| 曰老女人黄片| 少妇被粗大的猛进出69影院| 免费人妻精品一区二区三区视频| 久久午夜综合久久蜜桃| 每晚都被弄得嗷嗷叫到高潮| 成年人免费黄色播放视频| 精品国产一区二区久久| netflix在线观看网站| 国产成人系列免费观看| 满18在线观看网站| 国产亚洲欧美精品永久| 一级黄色大片毛片| 80岁老熟妇乱子伦牲交| 国产色视频综合| 啦啦啦中文免费视频观看日本| 国产欧美日韩一区二区三区在线| 99国产综合亚洲精品| 桃花免费在线播放| 国产成+人综合+亚洲专区| 日韩视频在线欧美| 一级,二级,三级黄色视频| av天堂久久9| 在线观看www视频免费| 99re6热这里在线精品视频| 欧美日韩一级在线毛片| 色播在线永久视频| 色在线成人网| 午夜福利在线观看吧| 午夜免费鲁丝| 亚洲全国av大片| 天天添夜夜摸| 午夜激情久久久久久久| 欧美日韩精品网址| 国产伦理片在线播放av一区| 日本av免费视频播放| 热re99久久国产66热| videos熟女内射| 51午夜福利影视在线观看| 国产精品 欧美亚洲| 亚洲精品在线美女| 国产免费现黄频在线看| 日韩一卡2卡3卡4卡2021年| av电影中文网址| 国产免费av片在线观看野外av| 麻豆国产av国片精品| 色播在线永久视频| 啦啦啦免费观看视频1| 夜夜爽天天搞| 我的亚洲天堂| 亚洲全国av大片| 免费女性裸体啪啪无遮挡网站| 久久久久精品国产欧美久久久| 在线观看www视频免费| 久久精品人人爽人人爽视色| 99久久精品国产亚洲精品| 超碰成人久久| 91麻豆精品激情在线观看国产 | 亚洲精品久久成人aⅴ小说| 亚洲五月婷婷丁香| 丁香六月欧美| 国产成人系列免费观看| 1024视频免费在线观看| 国产在线精品亚洲第一网站| 999久久久国产精品视频| 亚洲av日韩在线播放| 久久久久久久久免费视频了| 91成人精品电影| 欧美精品人与动牲交sv欧美| 久久 成人 亚洲| 午夜福利欧美成人| 丝袜人妻中文字幕| 国产免费视频播放在线视频| tube8黄色片| 99riav亚洲国产免费| 久久 成人 亚洲| 一个人免费看片子| 免费看十八禁软件| 亚洲色图综合在线观看| 新久久久久国产一级毛片| 亚洲精品中文字幕一二三四区 | 国产免费av片在线观看野外av| 亚洲熟女精品中文字幕| 亚洲精品av麻豆狂野| 国产在线观看jvid| 汤姆久久久久久久影院中文字幕| 在线观看免费高清a一片| 天天影视国产精品| 视频在线观看一区二区三区| 亚洲精品av麻豆狂野| 精品一区二区三卡| 69精品国产乱码久久久| 一边摸一边抽搐一进一出视频| 18禁国产床啪视频网站| 三级毛片av免费| 两性午夜刺激爽爽歪歪视频在线观看 | 成在线人永久免费视频| 美女扒开内裤让男人捅视频| 国产无遮挡羞羞视频在线观看| 视频区欧美日本亚洲| 亚洲精品乱久久久久久| 亚洲欧美精品综合一区二区三区| 国产午夜精品久久久久久| 午夜激情久久久久久久| 美女视频免费永久观看网站| 免费观看av网站的网址| 国产一区二区 视频在线| 自线自在国产av| 久久天躁狠狠躁夜夜2o2o| 我的亚洲天堂| 欧美日韩成人在线一区二区| 在线观看免费高清a一片| 精品人妻在线不人妻| 视频区图区小说| 下体分泌物呈黄色| 日本黄色日本黄色录像| 亚洲国产欧美在线一区| 欧美亚洲日本最大视频资源| 黄频高清免费视频| 亚洲情色 制服丝袜| 极品少妇高潮喷水抽搐| 国产精品偷伦视频观看了| 夜夜骑夜夜射夜夜干| 亚洲成人国产一区在线观看| 亚洲成a人片在线一区二区| 日韩中文字幕欧美一区二区| 午夜福利免费观看在线| 搡老乐熟女国产| 国产精品影院久久| 视频在线观看一区二区三区| 丰满少妇做爰视频| 亚洲精品一二三| 女人精品久久久久毛片| 一区二区三区激情视频| 欧美日韩亚洲综合一区二区三区_| 人人妻人人添人人爽欧美一区卜| 久久国产精品男人的天堂亚洲| 国产欧美日韩一区二区精品| av不卡在线播放| 男女床上黄色一级片免费看| 午夜久久久在线观看| 亚洲 国产 在线| 久久亚洲真实| 亚洲av成人一区二区三| 18禁国产床啪视频网站| 在线av久久热| 999精品在线视频| 性色av乱码一区二区三区2| av国产精品久久久久影院| 日韩 欧美 亚洲 中文字幕| 亚洲人成电影免费在线| 国产精品1区2区在线观看. | 黄片小视频在线播放| 国产区一区二久久| 欧美亚洲日本最大视频资源| 纵有疾风起免费观看全集完整版| 日日夜夜操网爽| 国产精品免费一区二区三区在线 | 国产黄色免费在线视频| 99国产精品一区二区三区| 天堂动漫精品| 亚洲精品国产一区二区精华液| 欧美+亚洲+日韩+国产| 丝袜人妻中文字幕| 成年人午夜在线观看视频| 一个人免费在线观看的高清视频| 人人澡人人妻人| 天天躁狠狠躁夜夜躁狠狠躁| 亚洲天堂av无毛| 亚洲精品一卡2卡三卡4卡5卡| 大陆偷拍与自拍| 青青草视频在线视频观看| 欧美在线黄色| 久久这里只有精品19| 国产无遮挡羞羞视频在线观看| 人成视频在线观看免费观看| 精品高清国产在线一区| av视频免费观看在线观看| 中文字幕色久视频| 精品人妻1区二区| 欧美精品一区二区免费开放| 在线十欧美十亚洲十日本专区| 国产免费现黄频在线看| 黄色怎么调成土黄色| 久久毛片免费看一区二区三区| 夜夜夜夜夜久久久久| 亚洲国产成人一精品久久久| 亚洲人成电影观看| 亚洲va日本ⅴa欧美va伊人久久| 色94色欧美一区二区| 久久热在线av| 国产又色又爽无遮挡免费看| 久久久精品国产亚洲av高清涩受| 超碰97精品在线观看| av超薄肉色丝袜交足视频| 99久久99久久久精品蜜桃| 国产欧美亚洲国产| 50天的宝宝边吃奶边哭怎么回事| 免费少妇av软件| 性少妇av在线| 悠悠久久av| 欧美精品一区二区免费开放| 久热爱精品视频在线9| 肉色欧美久久久久久久蜜桃| 丁香六月欧美| www.999成人在线观看| 高清在线国产一区| 日韩视频在线欧美| 国产麻豆69| 亚洲国产成人一精品久久久| 99riav亚洲国产免费| 欧美国产精品va在线观看不卡| 亚洲av日韩精品久久久久久密| 国产一区二区 视频在线| 久热这里只有精品99| 欧美成人免费av一区二区三区 | 一边摸一边抽搐一进一出视频| 精品国产乱码久久久久久男人| 亚洲欧美一区二区三区久久| 亚洲第一青青草原| 亚洲欧洲精品一区二区精品久久久| 日韩制服丝袜自拍偷拍| 免费一级毛片在线播放高清视频 | av超薄肉色丝袜交足视频| 欧美 亚洲 国产 日韩一| 人人妻人人澡人人爽人人夜夜| 久久午夜综合久久蜜桃| 国产91精品成人一区二区三区 | 99re在线观看精品视频| 久久久久视频综合| 涩涩av久久男人的天堂| 美女主播在线视频| 久久婷婷成人综合色麻豆| 久久性视频一级片| 变态另类成人亚洲欧美熟女 | 久久久久视频综合| 嫩草影视91久久| 日韩有码中文字幕| 激情视频va一区二区三区| 欧美在线一区亚洲| 一本综合久久免费| 午夜福利免费观看在线| 久久国产亚洲av麻豆专区| 日本黄色日本黄色录像| 久久午夜综合久久蜜桃| 考比视频在线观看| 黄色丝袜av网址大全| 亚洲精品久久成人aⅴ小说| 亚洲国产欧美在线一区| 国产精品久久久av美女十八| 久久久久久亚洲精品国产蜜桃av| 精品人妻1区二区| 国产区一区二久久| 欧美精品高潮呻吟av久久| 欧美黄色淫秽网站| 97在线人人人人妻| 老司机午夜十八禁免费视频| 天天躁夜夜躁狠狠躁躁| bbb黄色大片| 大型黄色视频在线免费观看| 欧美日韩成人在线一区二区| 999久久久精品免费观看国产| 国产熟女午夜一区二区三区| 又紧又爽又黄一区二区| 久久免费观看电影| 一边摸一边抽搐一进一出视频| 男女免费视频国产| 五月天丁香电影| 亚洲精品中文字幕一二三四区 | 51午夜福利影视在线观看| www日本在线高清视频| 久久人妻熟女aⅴ| 欧美性长视频在线观看| 久久精品国产亚洲av高清一级| 日韩免费高清中文字幕av| 日本a在线网址| 女性被躁到高潮视频| 欧美精品人与动牲交sv欧美| 欧美精品亚洲一区二区| 亚洲av欧美aⅴ国产| 一级a爱视频在线免费观看| a级毛片在线看网站| 欧美成人午夜精品| 日韩大码丰满熟妇| 999精品在线视频| 美女主播在线视频| 久久久久久久大尺度免费视频| 汤姆久久久久久久影院中文字幕| 黄网站色视频无遮挡免费观看| 欧美日韩成人在线一区二区| 午夜两性在线视频| 久久九九热精品免费| 999久久久精品免费观看国产| 99精品久久久久人妻精品| 在线亚洲精品国产二区图片欧美| 99精品久久久久人妻精品| 日韩制服丝袜自拍偷拍| 最新在线观看一区二区三区| 精品一区二区三卡| 亚洲avbb在线观看| 国产av一区二区精品久久| 免费av中文字幕在线| 亚洲精品一二三| 亚洲色图综合在线观看| 国产在线视频一区二区| 丁香六月欧美| 丝袜美腿诱惑在线| 亚洲一区中文字幕在线| av视频免费观看在线观看| 久久ye,这里只有精品| a级片在线免费高清观看视频| 制服人妻中文乱码| 国产在线免费精品| 中文字幕高清在线视频| 啦啦啦 在线观看视频| 大片免费播放器 马上看| 成在线人永久免费视频| 两性午夜刺激爽爽歪歪视频在线观看 | 色婷婷av一区二区三区视频| 国产野战对白在线观看| 性少妇av在线| 国产精品麻豆人妻色哟哟久久| 十分钟在线观看高清视频www| 国产真人三级小视频在线观看| 免费久久久久久久精品成人欧美视频| 亚洲视频免费观看视频| 50天的宝宝边吃奶边哭怎么回事| 国产三级黄色录像| 一本久久精品| 天天操日日干夜夜撸| 亚洲五月色婷婷综合| 欧美精品一区二区免费开放| 黑人巨大精品欧美一区二区蜜桃| 人人妻人人澡人人看| 免费看a级黄色片| 可以免费在线观看a视频的电影网站| 亚洲国产看品久久| 熟女少妇亚洲综合色aaa.| 蜜桃在线观看..| 久久国产亚洲av麻豆专区| 亚洲七黄色美女视频| 水蜜桃什么品种好| 青草久久国产| 欧美日韩黄片免| 国产欧美日韩一区二区精品| 十八禁高潮呻吟视频| 一本色道久久久久久精品综合| 欧美在线黄色| www.熟女人妻精品国产| 一二三四社区在线视频社区8| 91精品国产国语对白视频| 黑丝袜美女国产一区| 曰老女人黄片| 桃红色精品国产亚洲av| 日韩欧美国产一区二区入口| 91麻豆av在线| 国产一卡二卡三卡精品| 精品高清国产在线一区| 狠狠狠狠99中文字幕| 国产精品一区二区精品视频观看| 国产精品成人在线| 亚洲av成人不卡在线观看播放网| 18禁裸乳无遮挡动漫免费视频| 久久久欧美国产精品| av国产精品久久久久影院| 免费观看av网站的网址| 侵犯人妻中文字幕一二三四区| 后天国语完整版免费观看| 免费看a级黄色片| 免费久久久久久久精品成人欧美视频| 国产伦理片在线播放av一区| 脱女人内裤的视频| 欧美激情高清一区二区三区| 黄色成人免费大全| 亚洲一卡2卡3卡4卡5卡精品中文| 国产日韩欧美亚洲二区| 国产高清国产精品国产三级| 美女高潮到喷水免费观看| 激情在线观看视频在线高清 | 99热网站在线观看| avwww免费| 中文字幕最新亚洲高清| 90打野战视频偷拍视频| 一二三四在线观看免费中文在|