• <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.

    亚洲精品粉嫩美女一区| 久久久久久久久久久久大奶| 女性生殖器流出的白浆| 欧美日韩国产mv在线观看视频| 王馨瑶露胸无遮挡在线观看| 日韩中文字幕视频在线看片| 汤姆久久久久久久影院中文字幕| 欧美激情 高清一区二区三区| 在线永久观看黄色视频| 欧美激情 高清一区二区三区| 精品人妻熟女毛片av久久网站| 亚洲欧美激情在线| 日韩欧美国产一区二区入口| 久久天堂一区二区三区四区| 日韩熟女老妇一区二区性免费视频| 满18在线观看网站| 女性生殖器流出的白浆| 久久精品国产亚洲av高清一级| 亚洲色图 男人天堂 中文字幕| 午夜激情久久久久久久| 十八禁人妻一区二区| 一本色道久久久久久精品综合| 国产一区二区 视频在线| 亚洲av成人一区二区三| 亚洲午夜精品一区,二区,三区| 中文字幕精品免费在线观看视频| 国产aⅴ精品一区二区三区波| 国产成人免费无遮挡视频| 91精品国产国语对白视频| 国产99久久九九免费精品| 淫妇啪啪啪对白视频| 91精品国产国语对白视频| 免费一级毛片在线播放高清视频 | 日韩一区二区三区影片| 成人av一区二区三区在线看| 999久久久国产精品视频| 好男人电影高清在线观看| 日本五十路高清| 80岁老熟妇乱子伦牲交| 人妻一区二区av| 成人亚洲精品一区在线观看| 极品教师在线免费播放| 99国产精品一区二区蜜桃av | 精品国产乱码久久久久久小说| 一本色道久久久久久精品综合| 美女扒开内裤让男人捅视频| a级毛片黄视频| 国产在线观看jvid| 纵有疾风起免费观看全集完整版| 自拍欧美九色日韩亚洲蝌蚪91| 91老司机精品| 美女国产高潮福利片在线看| 日韩精品免费视频一区二区三区| 最新的欧美精品一区二区| 亚洲欧洲精品一区二区精品久久久| av免费在线观看网站| 国产淫语在线视频| av一本久久久久| 一级黄色大片毛片| 18禁国产床啪视频网站| 中文欧美无线码| 欧美大码av| 超碰成人久久| 美女午夜性视频免费| 亚洲av第一区精品v没综合| 亚洲免费av在线视频| 色播在线永久视频| 搡老岳熟女国产| 亚洲精品中文字幕一二三四区 | 欧美午夜高清在线| 后天国语完整版免费观看| 怎么达到女性高潮| 午夜福利视频在线观看免费| 欧美精品一区二区大全| www.精华液| 国产老妇伦熟女老妇高清| av在线播放免费不卡| 亚洲一码二码三码区别大吗| 91大片在线观看| 亚洲欧美精品综合一区二区三区| www.999成人在线观看| 91国产中文字幕| 色视频在线一区二区三区| 欧美亚洲日本最大视频资源| 19禁男女啪啪无遮挡网站| bbb黄色大片| 亚洲中文日韩欧美视频| 波多野结衣av一区二区av| 肉色欧美久久久久久久蜜桃| 在线播放国产精品三级| 国产成人欧美| 性高湖久久久久久久久免费观看| 日韩中文字幕视频在线看片| 在线观看66精品国产| 精品乱码久久久久久99久播| 99久久国产精品久久久| 在线十欧美十亚洲十日本专区| 精品国产一区二区三区四区第35| a在线观看视频网站| 久久久久久久久免费视频了| 嫩草影视91久久| 99精国产麻豆久久婷婷| 免费久久久久久久精品成人欧美视频| 黑人巨大精品欧美一区二区mp4| 亚洲精品一二三| 成人手机av| 成人18禁高潮啪啪吃奶动态图| 啦啦啦在线免费观看视频4| 中文字幕av电影在线播放| 久久国产精品人妻蜜桃| 亚洲avbb在线观看| 亚洲av成人一区二区三| 黄色 视频免费看| 欧美黑人精品巨大| 日韩视频一区二区在线观看| 中文欧美无线码| svipshipincom国产片| 中文字幕精品免费在线观看视频| 老熟妇仑乱视频hdxx| 欧美变态另类bdsm刘玥| 大陆偷拍与自拍| 久久午夜亚洲精品久久| 男女高潮啪啪啪动态图| 成人永久免费在线观看视频 | 美女主播在线视频| 国产一区二区 视频在线| 777久久人妻少妇嫩草av网站| 亚洲精品在线观看二区| av有码第一页| 久久这里只有精品19| 一本久久精品| 久久中文看片网| 在线观看免费午夜福利视频| 老熟妇乱子伦视频在线观看| 中文字幕人妻丝袜一区二区| 久久久久视频综合| 不卡一级毛片| 亚洲av第一区精品v没综合| 丰满人妻熟妇乱又伦精品不卡| 不卡一级毛片| 亚洲精品中文字幕一二三四区 | 女性生殖器流出的白浆| 国产精品欧美亚洲77777| 久久久精品区二区三区| 欧美日本中文国产一区发布| 欧美另类亚洲清纯唯美| 丝袜在线中文字幕| 久久久久久久精品吃奶| 国产精品电影一区二区三区 | 国产伦人伦偷精品视频| 伦理电影免费视频| 人人妻,人人澡人人爽秒播| 国产免费福利视频在线观看| 高潮久久久久久久久久久不卡| 黑人巨大精品欧美一区二区蜜桃| 精品免费久久久久久久清纯 | 精品少妇久久久久久888优播| 757午夜福利合集在线观看| 免费看a级黄色片| 日韩视频一区二区在线观看| 可以免费在线观看a视频的电影网站| 色综合婷婷激情| 国产精品一区二区精品视频观看| 一进一出好大好爽视频| 久久精品成人免费网站| 黄色丝袜av网址大全| 色老头精品视频在线观看| 精品亚洲成国产av| 国产欧美亚洲国产| 1024视频免费在线观看| 高清av免费在线| 亚洲va日本ⅴa欧美va伊人久久| 女人被躁到高潮嗷嗷叫费观| 亚洲中文字幕日韩| 天天影视国产精品| 亚洲九九香蕉| 女人精品久久久久毛片| 欧美乱码精品一区二区三区| 午夜免费成人在线视频| 露出奶头的视频| 国产麻豆69| 亚洲国产av影院在线观看| av片东京热男人的天堂| 久久热在线av| 国产精品免费视频内射| 欧美午夜高清在线| av有码第一页| 人人妻人人澡人人看| 搡老熟女国产l中国老女人| 久久久国产欧美日韩av| 国产免费现黄频在线看| 日本一区二区免费在线视频| 亚洲综合色网址| 伦理电影免费视频| 日日夜夜操网爽| 亚洲熟女毛片儿| 久久精品aⅴ一区二区三区四区| 久久99一区二区三区| 午夜激情久久久久久久| 久久久久网色| 一区福利在线观看| 久久毛片免费看一区二区三区| 国产男靠女视频免费网站| 99国产精品免费福利视频| 久久亚洲真实| av又黄又爽大尺度在线免费看| 免费在线观看视频国产中文字幕亚洲| 日本黄色视频三级网站网址 | 老熟妇仑乱视频hdxx| www.999成人在线观看| 日本av免费视频播放| 精品国内亚洲2022精品成人 | 精品久久久精品久久久| tube8黄色片| 国产有黄有色有爽视频| 国产av精品麻豆| 免费日韩欧美在线观看| 人妻 亚洲 视频| 黄色视频,在线免费观看| 夜夜夜夜夜久久久久| 夜夜骑夜夜射夜夜干| 欧美在线黄色| 男女免费视频国产| 午夜久久久在线观看| 青青草视频在线视频观看| 超碰97精品在线观看| 午夜福利一区二区在线看| 精品少妇黑人巨大在线播放| 视频在线观看一区二区三区| 欧美中文综合在线视频| 欧美乱码精品一区二区三区| 久久久精品区二区三区| 国产主播在线观看一区二区| 久久天躁狠狠躁夜夜2o2o| 在线播放国产精品三级| 亚洲专区国产一区二区| 久久久精品免费免费高清| 高清黄色对白视频在线免费看| 国产精品偷伦视频观看了| 高潮久久久久久久久久久不卡| 久久久久久久大尺度免费视频| 99精品欧美一区二区三区四区| 午夜福利一区二区在线看| av视频免费观看在线观看| 如日韩欧美国产精品一区二区三区| 夜夜夜夜夜久久久久| 亚洲成人手机| 在线观看www视频免费| 精品亚洲乱码少妇综合久久| 大陆偷拍与自拍| 国产欧美日韩综合在线一区二区| 国产97色在线日韩免费| 亚洲成人手机| av福利片在线| 国产黄色免费在线视频| 精品欧美一区二区三区在线| 丝袜美足系列| 久久精品亚洲熟妇少妇任你| 久久久久久免费高清国产稀缺| 巨乳人妻的诱惑在线观看| 欧美黄色片欧美黄色片| 久久精品成人免费网站| 国产亚洲午夜精品一区二区久久| 成人永久免费在线观看视频 | 亚洲国产欧美网| 午夜激情av网站| 亚洲 国产 在线| 亚洲av日韩精品久久久久久密| 在线av久久热| 亚洲专区字幕在线| 少妇精品久久久久久久| 亚洲精品美女久久av网站| 肉色欧美久久久久久久蜜桃| 亚洲人成伊人成综合网2020| 淫妇啪啪啪对白视频| 无人区码免费观看不卡 | 两性午夜刺激爽爽歪歪视频在线观看 | av一本久久久久| 亚洲成人免费电影在线观看| 日韩视频在线欧美| 国内毛片毛片毛片毛片毛片| av网站免费在线观看视频| 男女无遮挡免费网站观看| 国产精品99久久99久久久不卡| 欧美日韩av久久| 午夜福利视频在线观看免费| 久久毛片免费看一区二区三区| 757午夜福利合集在线观看| 亚洲成人免费av在线播放| 丝袜美腿诱惑在线| 女人精品久久久久毛片| 国产精品熟女久久久久浪| 一区福利在线观看| 国产熟女午夜一区二区三区| 一本大道久久a久久精品| 日韩精品免费视频一区二区三区| 这个男人来自地球电影免费观看| 搡老乐熟女国产| 老汉色∧v一级毛片| 黄色视频不卡| 人人妻人人澡人人看| av电影中文网址| 精品视频人人做人人爽| 视频在线观看一区二区三区| 9色porny在线观看| 日韩制服丝袜自拍偷拍| 国产成人精品久久二区二区91| 99久久99久久久精品蜜桃| 91精品三级在线观看| 美女扒开内裤让男人捅视频| 精品视频人人做人人爽| 高清av免费在线| 好男人电影高清在线观看| 亚洲熟女精品中文字幕| 国产一区二区三区视频了| 亚洲av国产av综合av卡| 精品福利观看| 欧美精品一区二区大全| 精品福利观看| 久久免费观看电影| 下体分泌物呈黄色| 男女下面插进去视频免费观看| 日韩免费高清中文字幕av| 国产成人精品在线电影| 亚洲成人手机| 免费高清在线观看日韩| 欧美日韩国产mv在线观看视频| 在线天堂中文资源库| 久久人妻av系列| 在线观看免费视频网站a站| 中文字幕色久视频| 欧美精品高潮呻吟av久久| 99久久精品国产亚洲精品| 好男人电影高清在线观看| 久久久国产成人免费| 成人特级黄色片久久久久久久 | 搡老乐熟女国产| 久久久精品免费免费高清| 人妻久久中文字幕网| 人人妻,人人澡人人爽秒播| 高清黄色对白视频在线免费看| av天堂久久9| 嫩草影视91久久| 亚洲五月婷婷丁香| 国产亚洲精品一区二区www | 欧美日韩亚洲综合一区二区三区_| 熟女少妇亚洲综合色aaa.| 女人被躁到高潮嗷嗷叫费观| 久久99一区二区三区| 一级毛片女人18水好多| 99国产极品粉嫩在线观看| 久久性视频一级片| 高清欧美精品videossex| 99国产精品免费福利视频| 午夜激情久久久久久久| 国产亚洲精品一区二区www | 美女高潮喷水抽搐中文字幕| 日韩中文字幕视频在线看片| 黄色 视频免费看| 亚洲国产精品一区二区三区在线| 免费一级毛片在线播放高清视频 | 高清黄色对白视频在线免费看| 高潮久久久久久久久久久不卡| 在线观看一区二区三区激情| 一级,二级,三级黄色视频| 亚洲av成人一区二区三| 日韩欧美一区二区三区在线观看 | 久久香蕉激情| 亚洲午夜精品一区,二区,三区| 亚洲熟妇熟女久久| 欧美日韩亚洲国产一区二区在线观看 | 无人区码免费观看不卡 | 这个男人来自地球电影免费观看| 王馨瑶露胸无遮挡在线观看| 亚洲国产欧美日韩在线播放| 精品福利永久在线观看| 亚洲成人手机| 老汉色av国产亚洲站长工具| 欧美国产精品一级二级三级| 老熟妇仑乱视频hdxx| 国产成人影院久久av| 国产视频一区二区在线看| 1024视频免费在线观看| 日本黄色视频三级网站网址 | 人人妻人人添人人爽欧美一区卜| 十八禁网站网址无遮挡| 一区二区三区国产精品乱码| 天天添夜夜摸| 久久精品亚洲精品国产色婷小说| 亚洲国产精品一区二区三区在线| 国产精品免费视频内射| 在线观看免费日韩欧美大片| 亚洲国产成人一精品久久久| 久久久精品国产亚洲av高清涩受| 18禁观看日本| 夜夜夜夜夜久久久久| 午夜福利视频在线观看免费| 美女扒开内裤让男人捅视频| 午夜福利一区二区在线看| 成人18禁在线播放| 欧美黑人精品巨大| 9色porny在线观看| 国产一区二区激情短视频| 国产成人精品久久二区二区免费| 91av网站免费观看| 亚洲一区中文字幕在线| 视频在线观看一区二区三区| 婷婷成人精品国产| 午夜老司机福利片| 免费黄频网站在线观看国产| 欧美中文综合在线视频| 久久av网站| 黄色a级毛片大全视频| 老汉色av国产亚洲站长工具| 男男h啪啪无遮挡| 黑丝袜美女国产一区| 亚洲自偷自拍图片 自拍| 高清黄色对白视频在线免费看| 两性午夜刺激爽爽歪歪视频在线观看 | 久久精品国产99精品国产亚洲性色 | av又黄又爽大尺度在线免费看| 国产精品.久久久| 777米奇影视久久| 丝瓜视频免费看黄片| 欧美乱妇无乱码| 国产精品一区二区精品视频观看| 天天添夜夜摸| 99精品在免费线老司机午夜| 国产91精品成人一区二区三区 | 国产精品久久久久久精品古装| 最近最新免费中文字幕在线| 香蕉久久夜色| 女人久久www免费人成看片| 成人18禁高潮啪啪吃奶动态图| 成年动漫av网址| 亚洲成国产人片在线观看| www.自偷自拍.com| 欧美日韩一级在线毛片| 亚洲精品成人av观看孕妇| 亚洲七黄色美女视频| 中文字幕精品免费在线观看视频| 十八禁网站免费在线| 高清在线国产一区| 欧美黑人精品巨大| 首页视频小说图片口味搜索| av在线播放免费不卡| 最新美女视频免费是黄的| 91精品三级在线观看| 国产成人精品久久二区二区免费| 最近最新中文字幕大全免费视频| 成人特级黄色片久久久久久久 | 亚洲中文字幕日韩| 日韩有码中文字幕| 国产精品久久久av美女十八| 男女高潮啪啪啪动态图| 一二三四社区在线视频社区8| 亚洲国产欧美网| 久久精品国产亚洲av高清一级| 亚洲欧美日韩另类电影网站| 成年人黄色毛片网站| 日韩免费av在线播放| 亚洲精品乱久久久久久| 1024香蕉在线观看| 美女福利国产在线| 99精国产麻豆久久婷婷| 国产精品免费视频内射| 日韩有码中文字幕| 欧美在线黄色| 国产在线精品亚洲第一网站| 国产一区二区三区在线臀色熟女 | 国产人伦9x9x在线观看| 性高湖久久久久久久久免费观看| 欧美人与性动交α欧美软件| 欧美日韩成人在线一区二区| 99国产极品粉嫩在线观看| 熟女少妇亚洲综合色aaa.| 一本一本久久a久久精品综合妖精| 真人做人爱边吃奶动态| 欧美日韩亚洲综合一区二区三区_| 老熟妇乱子伦视频在线观看| 久久久国产成人免费| 亚洲美女黄片视频| av天堂久久9| 视频区欧美日本亚洲| 极品少妇高潮喷水抽搐| svipshipincom国产片| 国产真人三级小视频在线观看| 少妇的丰满在线观看| 巨乳人妻的诱惑在线观看| 美女主播在线视频| e午夜精品久久久久久久| 韩国精品一区二区三区| 欧美日韩一级在线毛片| 国产在线一区二区三区精| 男人操女人黄网站| 欧美日韩黄片免| 最新的欧美精品一区二区| 久久国产精品影院| 叶爱在线成人免费视频播放| 亚洲欧美一区二区三区久久| 欧美激情 高清一区二区三区| 最近最新中文字幕大全免费视频| 在线观看舔阴道视频| videosex国产| 午夜免费鲁丝| 久久久水蜜桃国产精品网| 久久热在线av| 欧美变态另类bdsm刘玥| 成人18禁高潮啪啪吃奶动态图| 精品乱码久久久久久99久播| 免费在线观看日本一区| 女性生殖器流出的白浆| 丁香六月欧美| 精品久久久久久久毛片微露脸| 99九九在线精品视频| 老司机影院毛片| 老司机深夜福利视频在线观看| 欧美变态另类bdsm刘玥| e午夜精品久久久久久久| 天天影视国产精品| 午夜成年电影在线免费观看| kizo精华| 丰满迷人的少妇在线观看| 男女床上黄色一级片免费看| 黄色怎么调成土黄色| 日本黄色日本黄色录像| 久久青草综合色| a级片在线免费高清观看视频| 99热网站在线观看| 国产精品免费视频内射| 99国产精品一区二区蜜桃av | 午夜日韩欧美国产| 色94色欧美一区二区| 亚洲中文日韩欧美视频| 色94色欧美一区二区| 久久人人爽av亚洲精品天堂| 日韩欧美一区二区三区在线观看 | bbb黄色大片| 精品国产亚洲在线| 在线观看免费视频网站a站| 多毛熟女@视频| 美女午夜性视频免费| 法律面前人人平等表现在哪些方面| 日韩欧美一区视频在线观看| 欧美黄色淫秽网站| 涩涩av久久男人的天堂| 大码成人一级视频| 亚洲精品久久成人aⅴ小说| 91麻豆av在线| 91成人精品电影| 国产男女内射视频| 国产精品麻豆人妻色哟哟久久| 久久久国产成人免费| av免费在线观看网站| 高清视频免费观看一区二区| 一级毛片女人18水好多| 十八禁网站免费在线| 淫妇啪啪啪对白视频| av在线播放免费不卡| 欧美大码av| 美女高潮到喷水免费观看| 一区二区三区乱码不卡18| 大型av网站在线播放| 成年版毛片免费区| 老鸭窝网址在线观看| 建设人人有责人人尽责人人享有的| 中文字幕人妻熟女乱码| 91成人精品电影| 国产亚洲午夜精品一区二区久久| 欧美激情 高清一区二区三区| 国产成人av激情在线播放| 久久久水蜜桃国产精品网| netflix在线观看网站| 一夜夜www| 国产精品免费视频内射| 亚洲熟女精品中文字幕| 成人国产一区最新在线观看| 久久久久久亚洲精品国产蜜桃av| 啦啦啦免费观看视频1| 在线播放国产精品三级| 天堂动漫精品| 9191精品国产免费久久| 99国产精品一区二区蜜桃av | 日韩制服丝袜自拍偷拍| 久久久久视频综合| 两个人免费观看高清视频| 日本av手机在线免费观看| 亚洲人成77777在线视频| 国产精品熟女久久久久浪| 丝袜人妻中文字幕| 婷婷成人精品国产| 亚洲av欧美aⅴ国产| 蜜桃国产av成人99| 性少妇av在线| h视频一区二区三区| 黄片播放在线免费| 五月开心婷婷网| 免费在线观看视频国产中文字幕亚洲| 欧美在线一区亚洲| 老汉色∧v一级毛片| 亚洲欧洲日产国产| 丝袜美腿诱惑在线| 悠悠久久av| av一本久久久久| 精品免费久久久久久久清纯 | 捣出白浆h1v1| 欧美日韩av久久| 国产97色在线日韩免费| 啦啦啦 在线观看视频| 精品少妇一区二区三区视频日本电影| 国产单亲对白刺激| 女人久久www免费人成看片| 建设人人有责人人尽责人人享有的| 人妻久久中文字幕网| 黑丝袜美女国产一区|