# A Novel Efficient Feature Dimensionality Reduction Method and Its Application in Engineering

Authors: Zhun Cheng; Zhixiong Lu
Journal: Complexity
 (2018)

---

## Abstract

In the engineering field, excessive data dimensions affect the efficiency of machine learning and analysis of the relationships between data or features. To render feature dimensionality reduction more effective and faster, this paper proposes a new feature dimensionality reduction approach combining a sampling survey method with a heuristic intelligent optimization algorithm. Drawing on feature selection, this method builds a feature-scoring system and a reduced-dimension length-scoring system based on the sampling survey method. According to feature scores and reduced-dimension lengths, the method selects a number of features and reduced-dimension lengths that are ranked in the front with high scores. This feature dimensionality reduction method allows for in-depth optimal selection of features and reduced-dimension lengths with high scores using an improved heuristic intelligent optimization algorithm. To verify the effectiveness of the dimensionality reduction method, this paper applies it to road roughness time-domain estimation based on vehicle dynamic response and gene-selection research in bioengineering. Results in the first case show that the proposed method can improve the accuracy of road roughness time-domain estimation to above 0.99 and reduce measured data of the vehicle dynamic response, reducing the experimental workload significantly. Results in the second case show that the method can select a set of genes quickly and effectively with high disease recognition accuracy.


---

## Body

## 1. Introduction

The curse of dimensionality is a common problem in engineering. Ineffective and unreasonable feature dimensionality reduction will compromise machine-learning efficiency, pattern recognition accuracy, and data-mining efficiency while increasing the workload of measured data experiments to some extent. A set of high-dimension features possesses the following problems: too many features but few samples, too many features with few or no relations to the mining task, and excessive redundancy among features [1–4].Feature dimensionality reduction methods can be classified into two types: feature selection and feature projection. Feature projection is also calledfeature reparameterization and transforms the original feature space into one with a smaller dimensionality and more independent dimensions [5–7]. Feature projection generally uses principal component analysis (PCA), kernel principal component analysis, or other improved principal component analysis methods for linear or nonlinear transformations. Although these methods can realize feature dimensionality reduction to some extent, they require extensive calculations and high computational complexity. These dimensionality reduction methods also cannot reduce necessary experimental data and are inconvenient when identifying features heavily correlated to mining tasks. By contrast, feature selection obtains low-dimension features more effectively compared to machine learning or data mining and carries low computational complexity [8–10]. However, existing feature selection methods are still time-intensive with respect to calculations, and some data dimensionality reduction problems in engineering remain difficult to solve using these methods.The road roughness time-domain estimation of vehicle dynamic response is taken as an example in this paper, capturing feature dimensionality reduction in multivariate time series. Acquiring road roughness information based on vehicle dynamic response is an economical and practical method, but most researchers have only improved the precision of road roughness estimation by using different neural networks or enhanced estimation methods [11–14]; they have thus far failed to evaluate the types and quantities of dynamic responses to be assessed and the positions of points to be measured in the vehicle. Measurement points can change positions freely in the vehicle, and each position can offer three dynamic responses (vertical vibration acceleration, speed, and displacement); as such, the example has theoretically infinite feature dimensions. In this case, to improve the accuracy and speed of road roughness estimation based on neural networks, key features must be selected. Additionally, the subsequent experiment needs to test and measure the selected features, which guides the experiment and reduces the workload substantially.For instance, in gene-selection research on cancer, the mining of gene expression profile data can identify cancer types. The gene expression profile has many data dimensions but small samples. For instance, leukemia has 7129 gene dimensions but only 72 samples for analysis, which increases the challenges of feature dimensionality reduction. Various methods exist for gene selection, but they are computationally costly and complex [15–18]. Many researchers have introduced information to facilitate gene selection, but these methods cannot properly reflect the length of a set of gene features, namely, the number of dimensions after dimensionality reduction. Overall, engineering projects in road roughness time-domain estimation based on vehicle dynamic response and biological gene selection suffer from a curse of dimensionality with an extensive number of features. In addition to their computational cost and complexity, common dimensionality reduction methods (e.g., PCA and cluster analysis) fail to effectively reduce dimensionality (i.e., experimental data) and thus cannot improve estimation accuracy. If researchers use feature selection, uncertainty persists regarding the number of features selected such that the length of the final feature set is unclear. Several variables must be measured to determine which are uncertain.To solve the problem of excessive dimensions, this paper proposes a new feature selection method based on a sampling survey method and heuristic intelligent optimization algorithm. The general process is as follows: first, we build a feature-scoring system and a reduced-dimension length-scoring system based on the sampling survey method; second, we sort features and reduced-dimension lengths according to their scores and select the features and reduced-dimension lengths that are ranked in the front; third, to apply selection optimization to features and reduced-dimension lengths simultaneously, we redefine the meaning of the information carried by individuals in the heuristic intelligent optimization algorithm population and improve the algorithm accordingly. To better explain and verify the effect of the dimensionality reduction method, this paper conducts feature dimensionality reduction using road roughness time-domain estimation of vehicle dynamic response and gene-selection research in bioengineering. Results show that the new feature dimensionality reduction method proposed can select features quickly and effectively: the accuracy of road roughness time-domain estimation was generally higher than 0.99 and effectively reduced vehicle dynamic response measurement data, and a reasonable number of genes was selected.


## 2. New Feature Dimensionality Reduction Method Based on Sampling Survey Method and Improved Heuristic Intelligent Optimization Algorithm

### 2.1. Scoring System Based on Sampling Survey Method

The sampling survey method takes part of the population as samples for analysis [19]. In a typical dimensionality problem, given a vast number of features, the final number remains unknown after dimensionality reduction; therefore, using enumeration to analyze all samples requires an intensive workload and may fail. By employing the sampling survey method and a certain sample scale, this paper analyzes the features and number of them after dimensionality reduction. All samples must have the same chance of being selected. To obtain useful information for feature dimensionality reduction, we propose a feature-scoring system and a scoring system for reduced-dimension length.The feature-scoring system operates as follows:Step 1
Drawn1 sets of features with any number of dimensions taken randomly from features with equal probability. Calculate the fitness of each sample based on classification recognition accuracy, data-mining efficiency, or neural network estimation accuracy.
Step 2
Set a threshold valuep and select all samples with fitness larger than p; these samples comprise the set of valid feature scores.
Step 3
Score each feature based on the set of valid feature scores. The score of each feature is composed of the following two parts:
(1)Ffeature_score_1=∑i=1N1N1−Ri−Ffeature_score_1_minFfeature_score_1_max−Ffeature_score_1_min,where Ffeature_score_1 is the first part of the score (Score 1) of a feature, N1 is the number of sets of valid feature scores, Ri is the fitness ranking by the ith set of valid feature scores, Ffeature_score_1_min is the minimum of all Score 1s, and Ffeature_score_1_max is the maximum of all Score 1s.
(2)Ffeature_score_2=∑i=1N1fscore2i−Ffeature_score_2_minFfeature_score_2_max−Ffeature_score_2_min,where Ffeature_score_2 is the second part of the score (Score 2) of a feature; fscore2i is a piecewise function such that if the ith set of valid feature scores contains the feature, then fscore2i=1 and fscore2i=0 otherwise; Ffeature_score_2_min is the minimum of all Score 2s; and Ffeature_score_2_max is the maximum of all Score 2s.The total score of a feature is represented by the following formula:(3)Ffeature_score=Ffeature_score_1+Ffeature_score_2,where Ffeature_score is the total feature score.According to Formulas (1) and (2), Ffeature_score_1 and Ffeature_score_2, strictly speaking, each belong to the range of [0, 1]. Ffeature_score_1 reflects the degree of contribution of a feature to fitness. Ffeature_score_2 reflects the frequency of a feature in a set of valid feature scores. A feature with a higher Ffeature_score_1 or Ffeature_score_2 can realize a better fitness value.
Step 4
Rank features according to their total scores and choose the firstn2 features for the set of valid features.The scoring system for reduced-dimension length is as follows:Steps 1 and 2
These steps are the same as in the feature-scoring system
Step 3
Score the length of each reduced dimension according to the set of valid feature scores. The score of each dimension length is as follows:
(4)Fdimension_score=∑i=1N2ffitnessi,where Fdimension_score is the total reduced-dimension length score, N2 is the number of sets with the same dimension length in a set of valid feature scores, and ffitness is the fitness of each sample.
Step 4
Rank each dimension length according toFdimension_score, select those with a proportion of cumulative gross scores greater than 80% out of the total score of all dimension lengths, and identify reasonable dimension lengths after dimensionality reduction.To verify that a scoring system based on the sampling survey method can be used to effectively select features with the greatest influence on results, this paper uses a calculation example for simple mathematical analysis: the classical problem of Hald cement. Hald cement contains four chemical components: 3CaO·Al2O3 (x1), 3CaO·SiO2 (x2), 4CaO·Al2O3·Fe2O3 (x3), and 2CaO·SiO2 (x4). Table 1 lists data on relations between heat y (calorie) released by each gram of cement and the content of these four components.Table 1
Data on relations between contentx and heat y of four components of Hald cement.
Sample no.
x1
x2
x3
x4
y
1
7
26
6
60
78.5
2
1
29
15
52
74.3
3
11
56
8
20
104.3
4
11
31
8
47
87.6
5
7
52
6
33
95.9
6
11
55
9
22
109.2
7
3
71
17
6
102.7
8
1
31
22
44
72.5
9
2
54
18
22
93.1
10
21
47
4
26
115.9
11
1
40
23
34
83.8
12
11
66
9
12
113.3
13
10
68
8
12
109.4In statistics, variance reflects the degree to which two or more data points are different. Then, using the idea of calculating the degree of contribution in PCA, this paper takes the proportion of variance of one variable from the total variance of all variables to represent the degree of influence of that variable:(5)di=Varxi∑Varxi,where di is the degree of influence of the ith feature, Varxi is the variance of the ith feature, and ∑Varxi is the total variance of all features.Formula (5) can be used to calculate the degrees of influence of four components of Hald cement (x1, x2, x3, and x4), which are 0.0579, 0.4050, 0.0686, and 0.4686, respectively. x2 and x4 have high degrees of influence, x3 has less influence, and x1 has the least influence. To further reflect the degree of influence of these variables on heat y, we employed the radial basis function neural network (RBFNN) and considered several parts of the four components of Hald cement as input features to estimate y, the heat released by each gram of each component. If a feature does not have typical representativeness, then the estimation accuracy of RBFNN will be low. Therefore, we can rank the four features in terms of their influence degree according to the estimation accuracy after RBFNN machine learning. Considering a single feature as input into the neural network and the coefficient of determination as the evaluation index for the degree of closeness of heat y’s predicted value to the actual value, we obtained coefficients of determination for x1, x2, x3, and x4 as 0.7765, 0.9580, 0.8454, and 0.9495, respectively. The results are similar to the variance analysis; therefore, we can conclude that x2 and x4 exert great influences on heat y, x3 has less influence, and x1 has the least influence.According to Step 1 of the scoring system based on the sampling survey method proposed in this paper, we must also calculate the fitness of each sample by considering the accuracy of classification recognition, data-mining efficiency, or neural network estimation accuracy as the fitness to analyze the four features. To better compare their degrees of influence, we considered the output approximation degree of RBFNN when using the scoring system; Table2 lists the results.Table 2
Scores of four components of Hald cement.
x1
x2
x3
x4
0.0920
1.8962
0.6642
1.1559When ranking the four features by degree of influence in descending order, we obtainedx2, x4, x3, and x1. The results of dimensionality length show that when the dimensionality is 2, we can estimate heat y using the RBFNN with zero error. In this case, two features input into the neural network must contain either x2 or x4. The above analysis confirms that the proposed scoring system based on the sampling survey method can be used to evaluate the degree of feature importance to select an ideal feature-set length.


### 2.2. Improving the Heuristic Intelligent Optimization Algorithm to Enhance Algorithm Execution Efficiency

This paper uses an artificial fish swarm algorithm (AFSA) and particle swarm optimization (PSO) to explain and verify the new feature dimensionality reduction method [20–24]. To improve the rate of convergence and precision of the algorithm, we refine AFSA and PSO.AFSA is improved as follows:(1)
Add the memory function to all AFSA behaviors and store historical optimal positions and corresponding optimal values of individuals and the entire fish swarm
(2)
Introduce the speed and position calculation formulas of PSO into swarm foraging behavior and update and calculate the movement speed and position of the next-generation fish swarm
(3)
According to the idea of variation in a genetic algorithm, if the global optimum remains the same in multiple records, we shall disturb multiple fish in the swarm with a certain probabilityThe PSO is improved as follows:(1)
To avoid the prematurity phenomenon in PSO, introduce the probability judgment Metropolis rule of a simulated annealing (SA) algorithm
(6)P=exp−ΔET,where P is the probability of accepting a poor solution, ΔE is the difference in fitness between two generations, and T is the temperature variable in SA, which declines as the number of iterations increases.
(2)
Change the update criterion of temperature variableT to
(7)Ti=T0−TendN−iN+Tend,where Ti is the value of T of the ith iteration, T0 is the initial value of T, Tend is the final value of T, and N is the total number of iterations.
(3)
Taking the nonuniform mutation idea of a genetic algorithm as reference, disturb learning parameterC of flight speed in PSO to intensify the local search ability. The following formula calculates learning parameter C:
(8)Ci+1=Ci⋅rand⋅a1⋅N−iN,where Ci+1 is the value of C of the i+1th iteration, Ci is the value of C of the ith iteration, and a1 is a constant, such that a1=1 in this paper.To theoretically confirm the excellent optimum searching performance of the improved heuristic intelligent optimization algorithm, this paper analyzes convergence of the improved PSO and improved AFSA. Solis and Wets [25] studied the stochastic optimization algorithm and provided convergence criteria for stochastic optimization. For a type of stochastic optimization algorithm to converge to the global optimal solution, the following two conditions must be satisfied:
(1)
Condition H1:fDxk,ξ≤fxk, and if ξ∈A, then fDxk,ξ≤fξ. In this condition, A is the feasible solution space, D is the stochastic optimization algorithm, xk is the result of the kth iteration, xk+1=Dxk,ξ, and ξ is the solution searched by the algorithm in this iteration.Condition H1 ensures that the fitnessfxk in each iteration of the stochastic optimization algorithm is nonincremental. The stochastic optimization algorithm with convergence ability can get sequence fxkk=1∞ with a probability of 1, and the sequence converges to a certain value. A stochastic optimization algorithm satisfying condition H1 has convergence ability; although the final convergence value is not always a globally optimal solution, it may be locally optimal.
(2)
Condition H2: ForA’s arbitrary Borel subset B, if its measure >0, then ∏k=0∞1−μkB=0. μkB is the probability measure of the algorithm’s kth iteration result in set B.Condition H2 indicates that if a stochastic optimization algorithm can search and find the globally optimal solution, the algorithm goes through the globally optimal solution at least once when the number of iterations approaches infinity; that is, the algorithm has ergodicity. The convergence of the improved PSO is confirmed below.In the PSO, we define each particle by updating its speed and position based on the following formulas:(9)vk+1=c0vk+c1pbest−xk+c2gbest−xk,(10)xk+1=xk+vk+1,where vk+1 is the speed of a particle in the k+1th iteration, c0 is the inertia weight factor, c1 and c2 are learning factors, pbest is the historical best position of an individual, gbest is the historical best position of the swarm, and xk+1 is the position of a particle in the k+1th iteration.The PSO stores the historical optimal position of a swarm,gbest, thus ensuring that the optimization result obtained in each step of the algorithm is nonincremental; hence, the PSO satisfies condition H1.Then, we judge the convergence of particle speed and position according to convergence in probability. With an inertia weightc0 and the proposed learning factor calculation formula, when the number of iterations approaches infinity, we have limk→∞Pvk−0<ε=1, where ε is an arbitrarily small positive value. Suppose x∗ is a point of convergence in the PSO; then, we have limk→∞Pxk−x∗<ε=1.If the number of iterations is sufficiently large, then the particle’s position and update speed will converge tox∗ and 0, respectively, with a probability of 1.Additionally, we analyze and confirm the convergence of PSO using the difference equation [26]. Although vk and xk are multidimensional variables, the dimensions are mutually independent; therefore, we simplify algorithm analysis to be a 1D case. To simplify the calculation, we suppose that pbest, the position of optimal solution of a particle, and gbest, the position of the optimal solution of the whole swarm, are constant.Then, according to Formulas (9) and (10), we get
(11)vk+2=c0vk+1+c1pbest−xk+1+c2gbest−xk+1,(12)xk+2=xk+1+vk+2.Substituting Formulas (9)–(11) into Formula (12), we get
(13)xk+2+−c0+c1+c2−1xk+1+c0xk=c1pbest+c2gbest.It is a second-order constant-coefficient nonhomogeneous difference equation, and this paper uses the characteristic equation method to solve it.Equation (13)’s characteristic equation is λ2+−c0+c1+c2−1λ+c0=0. A quadratic equation with one unknown quantity can be divided into the following three cases: Δ=0, Δ>0, and Δ<0 (Δ=−c02+c1+c2−12−4c0).IfΔ=0, then we have λ1=λ2=−−c0+c1+c2−1/2. In this case, we obtain the following equation:
(14)xk=x0+1−c1−c2x0+c0v0+c1pbest+c2gbestλ1kλ1k,where x0 is the initial value of the algorithm iteration and v0 is the initial speed of the iteration.IfΔ>0, then we have λ1=c0−c1−c2+1+Δ/2 and λ2=c0−c1−c2+1−Δ/2. In this case, we get the following equations:
(15)xk=c1pbest+c2gbestc1+c2+λ2b1−b2λ2−λ1λ1k+b2−λ1b1λ2−λ1λ2k,b1=x0−c1pbest+c2gbestc1+c2,b2=1−c1−c2x0+c0v0+c1pbest+c2gbest−c1pbest+c2gbestc1+c2.IfΔ<0, then we have λ1=c0−c1−c2+1+i−Δ/2 and λ2=c0−c1−c2+1−i−Δ/2. xk’s calculation formula is the same as the solution in the case of Δ>0.Ifxk has a limit when k approaches infinity, then xk requires convergence in all three cases. We get λ1<1 and λ2<1. Essentially, the PSO’s convergence conditions are c0<1, c1+c2>0, and 2c0−c1−c2+2>0.We propose an improved PSO to improve learning factorsc1 and c2. With an increase in the number of iterations, c1 and c2 approach 0. Therefore, the improved PSO surely satisfies the above convergence conditions only if the inertia factor c0 has a proper value.We prove the convergence of the improved PSO proposed using the three methods above. However, if the algorithm is required to converge to the global optimal solution, it must satisfy condition H2. In this case, the algorithm goes through the global optimal solution at least once when its number of iterations approaches infinity; that is, the particle in the PSO should demonstrate ergodicity. Apparently, standard PSO cannot guarantee identification of the global optimal value, resulting in the aforementioned prematurity phenomenon. We introduce SA into the evolutionary criteria of PSO because SA has been theoretically proven to converge to the global optimal solution with a probability of 1 in certain conditions [27].Next, this paper proposes modifying the learning factor to provide the algorithm with a fast flying speed in early iterations and then obtaining a large search scope in the search space to search and find the global optimal solution when the search scope declines as the search proceeds. According to the derivation of reference [28], the standard AFSA is convergent. Proposed improvements introduce the speed and position calculation formulas of PSO into the foraging behaviors of the fish swarm. This way, in addition to clustering and tailgating behaviors, the entire fish swarm can share information about the historical optimal position of the swarm and the optimal positions of individuals during foraging behavior. Moreover, with the first improvement to AFSA (i.e., storing the historical optimal positions of individuals and the whole fish swarm and corresponding optimal values) and the convergence analysis of PSO, we can ensure that the improved AFSA satisfies the convergence judging condition H1 in the stochastic optimization algorithm.Because we use the speed and position update formulas of PSO in the fish swarm’s foraging behaviors where the fish swarm may also exhibit clustering and tailgating behaviors, all of which move in the direction with higher fitness, we suppose thatx∗ is the convergence point of AFSA. Then, we have limk→∞Pxk−x∗<ε=1 where ε is an arbitrarily small positive value; therefore, the AFSA possesses convergence.Similarly, to satisfy condition H2, the position of the fish swarm should have ergodicity in AFSA. Like PSO, the standard AFSA also suffers from the prematurity phenomenon. If a super individual appears in the early iterations of the algorithm, then the whole fish swarm will converge quickly to the value. Therefore, we propose a third AFSA improvement: if the historical optimal solution does not change in multiple iterations, then the algorithm randomly changes some individuals in the fish swarm. When the iteration approaches infinity, the fish swarm position will take all possible values in the definition domain with a theoretic probability of 1. Essentially, the improved AFSA is globally convergent with a certain probability.To verify the computational complexity of the improved PSO and improved AFSA, we choose the following test function [29]:
(16)maxfx=xsin10πx+2,where x is an independent variable within the domain of [−1, 2].For analysis, we used a Lenovo computer with an Intel(R) Core(TM) i5-4590 processor, CPU at 3.30 GHz, and 8 GB memory. The PSO and AFSA used 20 particles over 30 iterations. Figure1 presents the algorithm test results.Figure 1
Iterative evolution curves of test function.The standard AFSA and standard PSO were trapped in the area of the local optimal value in early iterations and slowly identified the global optimal value as the number of iterations increased, whereas the improved AFSA and improved PSO quickly identified the area of global optimal value in early iterations and finally converged to the global optimal value. The improved PSO and improved AFSA spent only 0.0904 s and 0.0719 s, respectively, converging to the global optimal value, whereas the standard PSO and standard AFSA took 0.0990 s and 0.5321 s, respectively.


### 2.3. Improving the Heuristic Intelligent Optimization Algorithm Based on the Need for Feature Dimensionality Reduction

To apply the heuristic intelligent optimization algorithm to optimal feature selection, we improved the algorithm based on classification recognition accuracy, data-mining efficiency, or neural network estimation precision as the target value of algorithm optimization. We redefined the individuals of each population in the intelligent optimization algorithm (i.e., redefined information carried by either the fish swarm in AFSA or each particle in PSO). By imitating the transcription and translation process of genetic information in biology, information carried by individuals is defined as shown in Figure2.Figure 2
Information carried by individuals in the algorithm.Here,S1 is the first information of an individual within a value range of [0, 1]. M is the total number of elements in the valid feature set. Table 3 lists reasonable dimension lengths after dimensionality reduction.Table 3
Selected intervals of reasonable dimension length after dimensionality reduction.
Lower limit of interval
Upper limit of interval
Reasonable dimension length
b1=0
c1
L1
b2
c2
L2
…
…
…
bn
cn=1
LnIn this case,L1, L2, and Ln represent reasonable dimension lengths after the first, second, and nth dimensionality reduction (i.e., the number of reasonable features after dimensionality reduction), respectively. Parameters are b1=0, c1=b2, c2=b3, …, cn−1=bn, and cn=1. c1−b1 is equal to the proportion of the total score of L1 out of the total score of all reasonable dimension lengths. Similarly, cn−bn is equal to the proportion of the total score of Ln out of the total score of all reasonable dimension lengths. The following is the calculation formula of ci−bi:
(17)ci−bi=Fdimension_score_Li∑j=1nFdimension_score_Lj,where ci and bi are the upper and lower limits Li, respectively, of the ith reasonable dimension length; Fdimension_score_Li is the total score of Li, and n is the total number of reasonable dimension lengths.If the value ofs1 is within the interval of the ith line in Table 1, then the corresponding reasonable dimension length Li is selected. Moreover, the heuristic intelligent algorithm will extract information from s2 to sLi+1 according to the value of Li, and the information will be used to calculate the target value of the algorithm.According to continuous iterative calculation of the heuristic intelligent optimization algorithm, the information carried by each individual will evolve in the direction of the optimal target value to obtain the most reasonable dimension length after dimensionality reduction and feature selection.


## 2.1. Scoring System Based on Sampling Survey Method

The sampling survey method takes part of the population as samples for analysis [19]. In a typical dimensionality problem, given a vast number of features, the final number remains unknown after dimensionality reduction; therefore, using enumeration to analyze all samples requires an intensive workload and may fail. By employing the sampling survey method and a certain sample scale, this paper analyzes the features and number of them after dimensionality reduction. All samples must have the same chance of being selected. To obtain useful information for feature dimensionality reduction, we propose a feature-scoring system and a scoring system for reduced-dimension length.The feature-scoring system operates as follows:Step 1
Drawn1 sets of features with any number of dimensions taken randomly from features with equal probability. Calculate the fitness of each sample based on classification recognition accuracy, data-mining efficiency, or neural network estimation accuracy.
Step 2
Set a threshold valuep and select all samples with fitness larger than p; these samples comprise the set of valid feature scores.
Step 3
Score each feature based on the set of valid feature scores. The score of each feature is composed of the following two parts:
(1)Ffeature_score_1=∑i=1N1N1−Ri−Ffeature_score_1_minFfeature_score_1_max−Ffeature_score_1_min,where Ffeature_score_1 is the first part of the score (Score 1) of a feature, N1 is the number of sets of valid feature scores, Ri is the fitness ranking by the ith set of valid feature scores, Ffeature_score_1_min is the minimum of all Score 1s, and Ffeature_score_1_max is the maximum of all Score 1s.
(2)Ffeature_score_2=∑i=1N1fscore2i−Ffeature_score_2_minFfeature_score_2_max−Ffeature_score_2_min,where Ffeature_score_2 is the second part of the score (Score 2) of a feature; fscore2i is a piecewise function such that if the ith set of valid feature scores contains the feature, then fscore2i=1 and fscore2i=0 otherwise; Ffeature_score_2_min is the minimum of all Score 2s; and Ffeature_score_2_max is the maximum of all Score 2s.The total score of a feature is represented by the following formula:(3)Ffeature_score=Ffeature_score_1+Ffeature_score_2,where Ffeature_score is the total feature score.According to Formulas (1) and (2), Ffeature_score_1 and Ffeature_score_2, strictly speaking, each belong to the range of [0, 1]. Ffeature_score_1 reflects the degree of contribution of a feature to fitness. Ffeature_score_2 reflects the frequency of a feature in a set of valid feature scores. A feature with a higher Ffeature_score_1 or Ffeature_score_2 can realize a better fitness value.
Step 4
Rank features according to their total scores and choose the firstn2 features for the set of valid features.The scoring system for reduced-dimension length is as follows:Steps 1 and 2
These steps are the same as in the feature-scoring system
Step 3
Score the length of each reduced dimension according to the set of valid feature scores. The score of each dimension length is as follows:
(4)Fdimension_score=∑i=1N2ffitnessi,where Fdimension_score is the total reduced-dimension length score, N2 is the number of sets with the same dimension length in a set of valid feature scores, and ffitness is the fitness of each sample.
Step 4
Rank each dimension length according toFdimension_score, select those with a proportion of cumulative gross scores greater than 80% out of the total score of all dimension lengths, and identify reasonable dimension lengths after dimensionality reduction.To verify that a scoring system based on the sampling survey method can be used to effectively select features with the greatest influence on results, this paper uses a calculation example for simple mathematical analysis: the classical problem of Hald cement. Hald cement contains four chemical components: 3CaO·Al2O3 (x1), 3CaO·SiO2 (x2), 4CaO·Al2O3·Fe2O3 (x3), and 2CaO·SiO2 (x4). Table 1 lists data on relations between heat y (calorie) released by each gram of cement and the content of these four components.Table 1
Data on relations between contentx and heat y of four components of Hald cement.
Sample no.
x1
x2
x3
x4
y
1
7
26
6
60
78.5
2
1
29
15
52
74.3
3
11
56
8
20
104.3
4
11
31
8
47
87.6
5
7
52
6
33
95.9
6
11
55
9
22
109.2
7
3
71
17
6
102.7
8
1
31
22
44
72.5
9
2
54
18
22
93.1
10
21
47
4
26
115.9
11
1
40
23
34
83.8
12
11
66
9
12
113.3
13
10
68
8
12
109.4In statistics, variance reflects the degree to which two or more data points are different. Then, using the idea of calculating the degree of contribution in PCA, this paper takes the proportion of variance of one variable from the total variance of all variables to represent the degree of influence of that variable:(5)di=Varxi∑Varxi,where di is the degree of influence of the ith feature, Varxi is the variance of the ith feature, and ∑Varxi is the total variance of all features.Formula (5) can be used to calculate the degrees of influence of four components of Hald cement (x1, x2, x3, and x4), which are 0.0579, 0.4050, 0.0686, and 0.4686, respectively. x2 and x4 have high degrees of influence, x3 has less influence, and x1 has the least influence. To further reflect the degree of influence of these variables on heat y, we employed the radial basis function neural network (RBFNN) and considered several parts of the four components of Hald cement as input features to estimate y, the heat released by each gram of each component. If a feature does not have typical representativeness, then the estimation accuracy of RBFNN will be low. Therefore, we can rank the four features in terms of their influence degree according to the estimation accuracy after RBFNN machine learning. Considering a single feature as input into the neural network and the coefficient of determination as the evaluation index for the degree of closeness of heat y’s predicted value to the actual value, we obtained coefficients of determination for x1, x2, x3, and x4 as 0.7765, 0.9580, 0.8454, and 0.9495, respectively. The results are similar to the variance analysis; therefore, we can conclude that x2 and x4 exert great influences on heat y, x3 has less influence, and x1 has the least influence.According to Step 1 of the scoring system based on the sampling survey method proposed in this paper, we must also calculate the fitness of each sample by considering the accuracy of classification recognition, data-mining efficiency, or neural network estimation accuracy as the fitness to analyze the four features. To better compare their degrees of influence, we considered the output approximation degree of RBFNN when using the scoring system; Table2 lists the results.Table 2
Scores of four components of Hald cement.
x1
x2
x3
x4
0.0920
1.8962
0.6642
1.1559When ranking the four features by degree of influence in descending order, we obtainedx2, x4, x3, and x1. The results of dimensionality length show that when the dimensionality is 2, we can estimate heat y using the RBFNN with zero error. In this case, two features input into the neural network must contain either x2 or x4. The above analysis confirms that the proposed scoring system based on the sampling survey method can be used to evaluate the degree of feature importance to select an ideal feature-set length.


## 2.2. Improving the Heuristic Intelligent Optimization Algorithm to Enhance Algorithm Execution Efficiency

This paper uses an artificial fish swarm algorithm (AFSA) and particle swarm optimization (PSO) to explain and verify the new feature dimensionality reduction method [20–24]. To improve the rate of convergence and precision of the algorithm, we refine AFSA and PSO.AFSA is improved as follows:(1)
Add the memory function to all AFSA behaviors and store historical optimal positions and corresponding optimal values of individuals and the entire fish swarm
(2)
Introduce the speed and position calculation formulas of PSO into swarm foraging behavior and update and calculate the movement speed and position of the next-generation fish swarm
(3)
According to the idea of variation in a genetic algorithm, if the global optimum remains the same in multiple records, we shall disturb multiple fish in the swarm with a certain probabilityThe PSO is improved as follows:(1)
To avoid the prematurity phenomenon in PSO, introduce the probability judgment Metropolis rule of a simulated annealing (SA) algorithm
(6)P=exp−ΔET,where P is the probability of accepting a poor solution, ΔE is the difference in fitness between two generations, and T is the temperature variable in SA, which declines as the number of iterations increases.
(2)
Change the update criterion of temperature variableT to
(7)Ti=T0−TendN−iN+Tend,where Ti is the value of T of the ith iteration, T0 is the initial value of T, Tend is the final value of T, and N is the total number of iterations.
(3)
Taking the nonuniform mutation idea of a genetic algorithm as reference, disturb learning parameterC of flight speed in PSO to intensify the local search ability. The following formula calculates learning parameter C:
(8)Ci+1=Ci⋅rand⋅a1⋅N−iN,where Ci+1 is the value of C of the i+1th iteration, Ci is the value of C of the ith iteration, and a1 is a constant, such that a1=1 in this paper.To theoretically confirm the excellent optimum searching performance of the improved heuristic intelligent optimization algorithm, this paper analyzes convergence of the improved PSO and improved AFSA. Solis and Wets [25] studied the stochastic optimization algorithm and provided convergence criteria for stochastic optimization. For a type of stochastic optimization algorithm to converge to the global optimal solution, the following two conditions must be satisfied:
(1)
Condition H1:fDxk,ξ≤fxk, and if ξ∈A, then fDxk,ξ≤fξ. In this condition, A is the feasible solution space, D is the stochastic optimization algorithm, xk is the result of the kth iteration, xk+1=Dxk,ξ, and ξ is the solution searched by the algorithm in this iteration.Condition H1 ensures that the fitnessfxk in each iteration of the stochastic optimization algorithm is nonincremental. The stochastic optimization algorithm with convergence ability can get sequence fxkk=1∞ with a probability of 1, and the sequence converges to a certain value. A stochastic optimization algorithm satisfying condition H1 has convergence ability; although the final convergence value is not always a globally optimal solution, it may be locally optimal.
(2)
Condition H2: ForA’s arbitrary Borel subset B, if its measure >0, then ∏k=0∞1−μkB=0. μkB is the probability measure of the algorithm’s kth iteration result in set B.Condition H2 indicates that if a stochastic optimization algorithm can search and find the globally optimal solution, the algorithm goes through the globally optimal solution at least once when the number of iterations approaches infinity; that is, the algorithm has ergodicity. The convergence of the improved PSO is confirmed below.In the PSO, we define each particle by updating its speed and position based on the following formulas:(9)vk+1=c0vk+c1pbest−xk+c2gbest−xk,(10)xk+1=xk+vk+1,where vk+1 is the speed of a particle in the k+1th iteration, c0 is the inertia weight factor, c1 and c2 are learning factors, pbest is the historical best position of an individual, gbest is the historical best position of the swarm, and xk+1 is the position of a particle in the k+1th iteration.The PSO stores the historical optimal position of a swarm,gbest, thus ensuring that the optimization result obtained in each step of the algorithm is nonincremental; hence, the PSO satisfies condition H1.Then, we judge the convergence of particle speed and position according to convergence in probability. With an inertia weightc0 and the proposed learning factor calculation formula, when the number of iterations approaches infinity, we have limk→∞Pvk−0<ε=1, where ε is an arbitrarily small positive value. Suppose x∗ is a point of convergence in the PSO; then, we have limk→∞Pxk−x∗<ε=1.If the number of iterations is sufficiently large, then the particle’s position and update speed will converge tox∗ and 0, respectively, with a probability of 1.Additionally, we analyze and confirm the convergence of PSO using the difference equation [26]. Although vk and xk are multidimensional variables, the dimensions are mutually independent; therefore, we simplify algorithm analysis to be a 1D case. To simplify the calculation, we suppose that pbest, the position of optimal solution of a particle, and gbest, the position of the optimal solution of the whole swarm, are constant.Then, according to Formulas (9) and (10), we get
(11)vk+2=c0vk+1+c1pbest−xk+1+c2gbest−xk+1,(12)xk+2=xk+1+vk+2.Substituting Formulas (9)–(11) into Formula (12), we get
(13)xk+2+−c0+c1+c2−1xk+1+c0xk=c1pbest+c2gbest.It is a second-order constant-coefficient nonhomogeneous difference equation, and this paper uses the characteristic equation method to solve it.Equation (13)’s characteristic equation is λ2+−c0+c1+c2−1λ+c0=0. A quadratic equation with one unknown quantity can be divided into the following three cases: Δ=0, Δ>0, and Δ<0 (Δ=−c02+c1+c2−12−4c0).IfΔ=0, then we have λ1=λ2=−−c0+c1+c2−1/2. In this case, we obtain the following equation:
(14)xk=x0+1−c1−c2x0+c0v0+c1pbest+c2gbestλ1kλ1k,where x0 is the initial value of the algorithm iteration and v0 is the initial speed of the iteration.IfΔ>0, then we have λ1=c0−c1−c2+1+Δ/2 and λ2=c0−c1−c2+1−Δ/2. In this case, we get the following equations:
(15)xk=c1pbest+c2gbestc1+c2+λ2b1−b2λ2−λ1λ1k+b2−λ1b1λ2−λ1λ2k,b1=x0−c1pbest+c2gbestc1+c2,b2=1−c1−c2x0+c0v0+c1pbest+c2gbest−c1pbest+c2gbestc1+c2.IfΔ<0, then we have λ1=c0−c1−c2+1+i−Δ/2 and λ2=c0−c1−c2+1−i−Δ/2. xk’s calculation formula is the same as the solution in the case of Δ>0.Ifxk has a limit when k approaches infinity, then xk requires convergence in all three cases. We get λ1<1 and λ2<1. Essentially, the PSO’s convergence conditions are c0<1, c1+c2>0, and 2c0−c1−c2+2>0.We propose an improved PSO to improve learning factorsc1 and c2. With an increase in the number of iterations, c1 and c2 approach 0. Therefore, the improved PSO surely satisfies the above convergence conditions only if the inertia factor c0 has a proper value.We prove the convergence of the improved PSO proposed using the three methods above. However, if the algorithm is required to converge to the global optimal solution, it must satisfy condition H2. In this case, the algorithm goes through the global optimal solution at least once when its number of iterations approaches infinity; that is, the particle in the PSO should demonstrate ergodicity. Apparently, standard PSO cannot guarantee identification of the global optimal value, resulting in the aforementioned prematurity phenomenon. We introduce SA into the evolutionary criteria of PSO because SA has been theoretically proven to converge to the global optimal solution with a probability of 1 in certain conditions [27].Next, this paper proposes modifying the learning factor to provide the algorithm with a fast flying speed in early iterations and then obtaining a large search scope in the search space to search and find the global optimal solution when the search scope declines as the search proceeds. According to the derivation of reference [28], the standard AFSA is convergent. Proposed improvements introduce the speed and position calculation formulas of PSO into the foraging behaviors of the fish swarm. This way, in addition to clustering and tailgating behaviors, the entire fish swarm can share information about the historical optimal position of the swarm and the optimal positions of individuals during foraging behavior. Moreover, with the first improvement to AFSA (i.e., storing the historical optimal positions of individuals and the whole fish swarm and corresponding optimal values) and the convergence analysis of PSO, we can ensure that the improved AFSA satisfies the convergence judging condition H1 in the stochastic optimization algorithm.Because we use the speed and position update formulas of PSO in the fish swarm’s foraging behaviors where the fish swarm may also exhibit clustering and tailgating behaviors, all of which move in the direction with higher fitness, we suppose thatx∗ is the convergence point of AFSA. Then, we have limk→∞Pxk−x∗<ε=1 where ε is an arbitrarily small positive value; therefore, the AFSA possesses convergence.Similarly, to satisfy condition H2, the position of the fish swarm should have ergodicity in AFSA. Like PSO, the standard AFSA also suffers from the prematurity phenomenon. If a super individual appears in the early iterations of the algorithm, then the whole fish swarm will converge quickly to the value. Therefore, we propose a third AFSA improvement: if the historical optimal solution does not change in multiple iterations, then the algorithm randomly changes some individuals in the fish swarm. When the iteration approaches infinity, the fish swarm position will take all possible values in the definition domain with a theoretic probability of 1. Essentially, the improved AFSA is globally convergent with a certain probability.To verify the computational complexity of the improved PSO and improved AFSA, we choose the following test function [29]:
(16)maxfx=xsin10πx+2,where x is an independent variable within the domain of [−1, 2].For analysis, we used a Lenovo computer with an Intel(R) Core(TM) i5-4590 processor, CPU at 3.30 GHz, and 8 GB memory. The PSO and AFSA used 20 particles over 30 iterations. Figure1 presents the algorithm test results.Figure 1
Iterative evolution curves of test function.The standard AFSA and standard PSO were trapped in the area of the local optimal value in early iterations and slowly identified the global optimal value as the number of iterations increased, whereas the improved AFSA and improved PSO quickly identified the area of global optimal value in early iterations and finally converged to the global optimal value. The improved PSO and improved AFSA spent only 0.0904 s and 0.0719 s, respectively, converging to the global optimal value, whereas the standard PSO and standard AFSA took 0.0990 s and 0.5321 s, respectively.


## 2.3. Improving the Heuristic Intelligent Optimization Algorithm Based on the Need for Feature Dimensionality Reduction

To apply the heuristic intelligent optimization algorithm to optimal feature selection, we improved the algorithm based on classification recognition accuracy, data-mining efficiency, or neural network estimation precision as the target value of algorithm optimization. We redefined the individuals of each population in the intelligent optimization algorithm (i.e., redefined information carried by either the fish swarm in AFSA or each particle in PSO). By imitating the transcription and translation process of genetic information in biology, information carried by individuals is defined as shown in Figure2.Figure 2
Information carried by individuals in the algorithm.Here,S1 is the first information of an individual within a value range of [0, 1]. M is the total number of elements in the valid feature set. Table 3 lists reasonable dimension lengths after dimensionality reduction.Table 3
Selected intervals of reasonable dimension length after dimensionality reduction.
Lower limit of interval
Upper limit of interval
Reasonable dimension length
b1=0
c1
L1
b2
c2
L2
…
…
…
bn
cn=1
LnIn this case,L1, L2, and Ln represent reasonable dimension lengths after the first, second, and nth dimensionality reduction (i.e., the number of reasonable features after dimensionality reduction), respectively. Parameters are b1=0, c1=b2, c2=b3, …, cn−1=bn, and cn=1. c1−b1 is equal to the proportion of the total score of L1 out of the total score of all reasonable dimension lengths. Similarly, cn−bn is equal to the proportion of the total score of Ln out of the total score of all reasonable dimension lengths. The following is the calculation formula of ci−bi:
(17)ci−bi=Fdimension_score_Li∑j=1nFdimension_score_Lj,where ci and bi are the upper and lower limits Li, respectively, of the ith reasonable dimension length; Fdimension_score_Li is the total score of Li, and n is the total number of reasonable dimension lengths.If the value ofs1 is within the interval of the ith line in Table 1, then the corresponding reasonable dimension length Li is selected. Moreover, the heuristic intelligent algorithm will extract information from s2 to sLi+1 according to the value of Li, and the information will be used to calculate the target value of the algorithm.According to continuous iterative calculation of the heuristic intelligent optimization algorithm, the information carried by each individual will evolve in the direction of the optimal target value to obtain the most reasonable dimension length after dimensionality reduction and feature selection.


## 3. Road Roughness Time-Domain Estimation Based on Proposed Dimensionality Reduction Method

Researchers generally measure contacting or noncontacting road roughness to determine roughness; however, this method returns measurement results slowly and requires expensive equipment. A method combining a dynamic response and neural network can quickly estimate road roughness. Yet little research exists about the types and quantity of dynamic response and the points to be measured in the vehicle. Because measurement points can change positions freely in the vehicle, and each position can offer three dynamic responses (i.e., vertical vibration acceleration, speed, and displacement), this example theoretically has infinite feature dimensions.


### 3.1. Building a Vehicle Vibration Model with Freely Changing Measurement Points in the Vehicle

Assume the vehicle is symmetric along the longitudinal axis, and the left and right wheels are under the same road conditions while driving. Given these assumptions, we simplify the vehicle model into a half-vehicle model and the vehicle body into a rigid rod. The vibration model in this paper only considers degrees of freedom: the vehicle body’s vertical vibration, pitch-angle vibration, vertical vibration of the seat system, vertical vibration of front wheel, and vibration of the rear wheel (Figure3). We consider the points in the vehicle, of which the distance to the center of mass is e, as the dynamic response measurement points and build the vibration equation based on Lagrange’s second equation.Figure 3
Vehicle vibration model with five degrees of freedom.In Figure3, Z1, Z2, Z3, Zc, and α denote the vertical vibration displacements (m) of the front wheel, rear wheel, seat system, center of mass of the vehicle body, and pitch-angle displacement of the vehicle body (rad), respectively; m1, m2, m3, mc, and J are the mass of the front wheel, rear wheel, seat system, vehicle body (kg), and rotational inertia of the vehicle body’s mass center around the horizontal axis (kg·m2), respectively; a, b, l, and e represent the distance from the front axle, back axle, seat system, and measurement point to the center of mass (m), respectively; C1, C2, C3, Cf, and Cr are the damping of the front suspension, rear suspension, seat system, front wheel, and rear wheel (N·s/m), respectively; K1, K2, K3, Kf, and Kr are the rigidity of the front suspension, rear suspension, seat system, front wheel, and rear wheel (kN/m), respectively; and Qf and Qr are the vertical displacement of road roughness to the front wheel and rear wheel (m), respectively.The following is Lagrange’s second equation:(18)ddt∂L∂q̇i−∂L∂qi+∂U∂q̇i=Qi,where q denotes the generalized coordinates; L=T−V represents the difference between the expression of kinetic energy T and potential energy V; U is the dissipated energy, referring to energy consumed by the damping element; and Qi is the generalized force.Because the pitch angle has a smaller vibration, the paper takessinα≈α. We assume the vertical vibration displacement of the measurement point to be ZD and the distances from the measurement point to the front axle and the rear axle to be l1 and l2, respectively.The kinetic energy of the system is(19)T=12mcŻD+eθ̇D2+12Jα2+12m1Ż12+12m2Ż22+12m3Ż32.The potential energy of the system is(20)V=12K1ZD−l1α−Z12+12K2ZD+l2α−Z22+12K3ZD+lα−Z32+12KfZ1−Qf2+12KrZ2−Qr2.The dissipated energy of the system is(21)U=12C1ŻD−l1α̇−Ż12+12C2ŻD+l2α̇−Ż22+12C3ŻD+lα̇−Ż32+12CfŻ1−Q̇f+12CrŻ2−Q̇r.Substituting Equations (19)–(21) into Equation (18), we obtain the vehicle vibration differential equation in the example. Figure 4 illustrates the simulation model built in Matlab/Simulink.Figure 4
Vibration simulation model.According to the vibration simulation model, we can obtain dynamic response data conveniently. The vibration model uses the filtered white noise method to produce the road roughness signal, meeting experiment requirements, and considers a level-B road as the road to be estimated.


### 3.2. Research on Feature Dimensionality Reduction in Road Roughness Estimation Using the Proposed Method

This example uses the improved AFSA for optimal feature selection. We number the vertical accelerated speed, vertical speed, and vertical displacement of the positions of measurement point, center of the front wheel, center of the rear wheel, and center of the driver’s seat as shown in Table4.Table 4
Numbers of dynamic responses.
Measurement point
Front wheel
Rear wheel
Driver’s seat
Accelerated speed
1
3
4
5
Speed
6
8
9
10
Displacement
11
13
14
15Numbers 2, 7, and 12 are the pitch-angle accelerated speed, pitch-angle speed, and pitch angle of the measurement point, respectively.The example uses the improved AFSA and assumes information carried by theith fish to be Xi.
(22)Xi=xi1,xi2,xi3,…,xi17,where xi1 is the value of the first information, within the range of [0, 1]; xi2,…,xi16 are positive integers within [1, 15], the values of which represent the corresponding dynamic response numbers; and xi17 is the distance from the measurement point to the center of mass, e, and e∈−1.25,1.51 in the example.In the iteration process, the improved AFSA finds the corresponding reasonable dimension lengthLi from Table 3 according to the value of xi1 and selects the first Li elements from xi2,…,xi16 accordingly.The example uses the RBFNN for time-domain estimation of road roughness. We used the improved AFSA and conducted 12 total optimization experiments; results appear in Table5.Table 5
Optimization selection results of dynamic response features.
Test number
Features selected
R2
Distancee
1
[1, 3, 4, 7, 8, 12, 13]
0.9959
−0.5574
2
[1, 3, 4, 12, 13, 15]
0.9934
0.3557
3
[2, 3, 4, 13]
0.9939
0.9049
4
[1, 4, 8, 9]
0.9948
−1.4350
5
[2, 3, 4, 13]
0.9944
−0.0074
6
[1, 3, 4, 8, 9, 14]
0.9890
−1.2500
7
[2, 3, 4, 13]
0.9944
−0.5321
8
[2, 3, 4, 6, 10, 12]
0.9904
−1.0486
9
[1, 4, 8, 9]
0.9932
−1.4712
10
[1, 3, 4, 9, 11, 13]
0.9823
−0.3755
11
[2, 4, 7, 13, 14]
0.9892
−0.9510
12
[3, 4, 8, 13, 14, 15]
0.9649
0.4839In Table5, R2 is the coefficient of determination. The feature dimensionality reduction method proposed can choose features with great influence on the precision of road roughness estimation. The sets of features selected in Experiments 3 and 4 were applied in Experiments 7 and 9; thus, we used the features selected in Experiments 3 and 4 for time-domain estimation of road roughness. Figures 5 and 6 display the estimation results.Figure 5
Road roughness estimation results with features in Experiment 3.Figure 6
Road roughness estimation results with features in Experiment 4.The actual value was close to the estimated value. The coefficient of determination exceeded 0.99 in both tests, indicating that the proposed feature dimensionality reduction method can select a feature set with high contributions to road roughness using only four sets of features; thus, we can accurately estimate road roughness using only four dynamic response values. To further confirm this result, we used the random forest method to estimate road roughness based on the results of Experiments 3 and 4. Estimation results are presented in Figures7 and 8.Figure 7
Road roughness estimation results based on the random forest method and features in Experiment 3.Figure 8
Road roughness estimation results based on the random forest method and features in Experiment 4.The road roughness estimation results obtained with feature sets from Experiments 3 and 4 based on the random forest method were 0.9489 and 0.9708, respectively. Although the estimation accuracy of the random forest method was slightly lower than that of RBFNN, it was still quite high, suggesting that the set of dynamic response features selected with the proposed method can match other estimation data methods perfectly. When the proposed feature selection method was combined with a prediction method, only a few dynamic response parameters were needed to estimate road roughness accurately.


## 3.1. Building a Vehicle Vibration Model with Freely Changing Measurement Points in the Vehicle

Assume the vehicle is symmetric along the longitudinal axis, and the left and right wheels are under the same road conditions while driving. Given these assumptions, we simplify the vehicle model into a half-vehicle model and the vehicle body into a rigid rod. The vibration model in this paper only considers degrees of freedom: the vehicle body’s vertical vibration, pitch-angle vibration, vertical vibration of the seat system, vertical vibration of front wheel, and vibration of the rear wheel (Figure3). We consider the points in the vehicle, of which the distance to the center of mass is e, as the dynamic response measurement points and build the vibration equation based on Lagrange’s second equation.Figure 3
Vehicle vibration model with five degrees of freedom.In Figure3, Z1, Z2, Z3, Zc, and α denote the vertical vibration displacements (m) of the front wheel, rear wheel, seat system, center of mass of the vehicle body, and pitch-angle displacement of the vehicle body (rad), respectively; m1, m2, m3, mc, and J are the mass of the front wheel, rear wheel, seat system, vehicle body (kg), and rotational inertia of the vehicle body’s mass center around the horizontal axis (kg·m2), respectively; a, b, l, and e represent the distance from the front axle, back axle, seat system, and measurement point to the center of mass (m), respectively; C1, C2, C3, Cf, and Cr are the damping of the front suspension, rear suspension, seat system, front wheel, and rear wheel (N·s/m), respectively; K1, K2, K3, Kf, and Kr are the rigidity of the front suspension, rear suspension, seat system, front wheel, and rear wheel (kN/m), respectively; and Qf and Qr are the vertical displacement of road roughness to the front wheel and rear wheel (m), respectively.The following is Lagrange’s second equation:(18)ddt∂L∂q̇i−∂L∂qi+∂U∂q̇i=Qi,where q denotes the generalized coordinates; L=T−V represents the difference between the expression of kinetic energy T and potential energy V; U is the dissipated energy, referring to energy consumed by the damping element; and Qi is the generalized force.Because the pitch angle has a smaller vibration, the paper takessinα≈α. We assume the vertical vibration displacement of the measurement point to be ZD and the distances from the measurement point to the front axle and the rear axle to be l1 and l2, respectively.The kinetic energy of the system is(19)T=12mcŻD+eθ̇D2+12Jα2+12m1Ż12+12m2Ż22+12m3Ż32.The potential energy of the system is(20)V=12K1ZD−l1α−Z12+12K2ZD+l2α−Z22+12K3ZD+lα−Z32+12KfZ1−Qf2+12KrZ2−Qr2.The dissipated energy of the system is(21)U=12C1ŻD−l1α̇−Ż12+12C2ŻD+l2α̇−Ż22+12C3ŻD+lα̇−Ż32+12CfŻ1−Q̇f+12CrŻ2−Q̇r.Substituting Equations (19)–(21) into Equation (18), we obtain the vehicle vibration differential equation in the example. Figure 4 illustrates the simulation model built in Matlab/Simulink.Figure 4
Vibration simulation model.According to the vibration simulation model, we can obtain dynamic response data conveniently. The vibration model uses the filtered white noise method to produce the road roughness signal, meeting experiment requirements, and considers a level-B road as the road to be estimated.


## 3.2. Research on Feature Dimensionality Reduction in Road Roughness Estimation Using the Proposed Method

This example uses the improved AFSA for optimal feature selection. We number the vertical accelerated speed, vertical speed, and vertical displacement of the positions of measurement point, center of the front wheel, center of the rear wheel, and center of the driver’s seat as shown in Table4.Table 4
Numbers of dynamic responses.
Measurement point
Front wheel
Rear wheel
Driver’s seat
Accelerated speed
1
3
4
5
Speed
6
8
9
10
Displacement
11
13
14
15Numbers 2, 7, and 12 are the pitch-angle accelerated speed, pitch-angle speed, and pitch angle of the measurement point, respectively.The example uses the improved AFSA and assumes information carried by theith fish to be Xi.
(22)Xi=xi1,xi2,xi3,…,xi17,where xi1 is the value of the first information, within the range of [0, 1]; xi2,…,xi16 are positive integers within [1, 15], the values of which represent the corresponding dynamic response numbers; and xi17 is the distance from the measurement point to the center of mass, e, and e∈−1.25,1.51 in the example.In the iteration process, the improved AFSA finds the corresponding reasonable dimension lengthLi from Table 3 according to the value of xi1 and selects the first Li elements from xi2,…,xi16 accordingly.The example uses the RBFNN for time-domain estimation of road roughness. We used the improved AFSA and conducted 12 total optimization experiments; results appear in Table5.Table 5
Optimization selection results of dynamic response features.
Test number
Features selected
R2
Distancee
1
[1, 3, 4, 7, 8, 12, 13]
0.9959
−0.5574
2
[1, 3, 4, 12, 13, 15]
0.9934
0.3557
3
[2, 3, 4, 13]
0.9939
0.9049
4
[1, 4, 8, 9]
0.9948
−1.4350
5
[2, 3, 4, 13]
0.9944
−0.0074
6
[1, 3, 4, 8, 9, 14]
0.9890
−1.2500
7
[2, 3, 4, 13]
0.9944
−0.5321
8
[2, 3, 4, 6, 10, 12]
0.9904
−1.0486
9
[1, 4, 8, 9]
0.9932
−1.4712
10
[1, 3, 4, 9, 11, 13]
0.9823
−0.3755
11
[2, 4, 7, 13, 14]
0.9892
−0.9510
12
[3, 4, 8, 13, 14, 15]
0.9649
0.4839In Table5, R2 is the coefficient of determination. The feature dimensionality reduction method proposed can choose features with great influence on the precision of road roughness estimation. The sets of features selected in Experiments 3 and 4 were applied in Experiments 7 and 9; thus, we used the features selected in Experiments 3 and 4 for time-domain estimation of road roughness. Figures 5 and 6 display the estimation results.Figure 5
Road roughness estimation results with features in Experiment 3.Figure 6
Road roughness estimation results with features in Experiment 4.The actual value was close to the estimated value. The coefficient of determination exceeded 0.99 in both tests, indicating that the proposed feature dimensionality reduction method can select a feature set with high contributions to road roughness using only four sets of features; thus, we can accurately estimate road roughness using only four dynamic response values. To further confirm this result, we used the random forest method to estimate road roughness based on the results of Experiments 3 and 4. Estimation results are presented in Figures7 and 8.Figure 7
Road roughness estimation results based on the random forest method and features in Experiment 3.Figure 8
Road roughness estimation results based on the random forest method and features in Experiment 4.The road roughness estimation results obtained with feature sets from Experiments 3 and 4 based on the random forest method were 0.9489 and 0.9708, respectively. Although the estimation accuracy of the random forest method was slightly lower than that of RBFNN, it was still quite high, suggesting that the set of dynamic response features selected with the proposed method can match other estimation data methods perfectly. When the proposed feature selection method was combined with a prediction method, only a few dynamic response parameters were needed to estimate road roughness accurately.


## 4. Gene-Selection Research in Bioengineering Using the Proposed Dimensionality Reduction Method

The experimental subjects in the example consisted of four open microarray datasets on leukemia, the colon, SRBCT, and brain cancer. Leukemia [30] contained 72 samples, each with 7129 genes. Acute myeloid leukemia (AML) contained 25 samples, and acute lymphoblastic leukemia contained 47. The colon [31] contained 62 samples (40 of tumor tissue and 22 of normal tissue), each with 2000 genes. SRBCT [32] included 83 samples, each with 2308 genes: 29 samples of the Ewing family of tumors, 18 samples of neuroblastoma, 25 samples of rhabdomyosarcoma, and 11 samples of Burkitt lymphomas.Brain cancer [33] contained 60 samples, each with 7219 genes: 46 samples of patients with classic brain cancer and 14 samples with desmoplastic brain cancer. Although brain cancer only had two classes, the dataset has proven difficult to classify in current research with classification accuracy lower than 85% in much of the literature. Lung cancer [34] included 203 samples, each with 3312 genes. Samples were classified into five subtypes: 139 adenocarcinoma, 20 pulmonary carcinoids, 21 squamous-cell lung carcinomas, six small-cell lung carcinomas, and 17 normal lungs. Lymphoma had 58 samples, each with 7129 genes. The samples included two subtypes: 32 cured patients and 26 not cured.Because each type of tumor disease is related to many genes, we roughly extracted genes using the information index to classification (IIC) [34]. This example uses the improved PSO and assumes information carried by the ith particle to be Yi.
(23)Yi=yi1,yi2,yi3,…,yin1,where yi1 is the value of the first information, within the range of [0, 1]; yi2,…,yin1 are positive integers within the range of [1, n1], and their values represent corresponding gene numbers; and n1 is the number of roughly extracted genes based on IIC.The example uses the extreme learning machine (ELM) to recognize tumor subtypes and takes the mean and standard deviations of the 5-fold CV accuracy of 100 experiments to verify the accuracy of extracted genes.To demonstrate the advantage of the proposed method in sample classification accuracy, we chose five methods with high gene-selection accuracy, including binary particle swarm optimization BPSO-ELM, K-means-PSO-ELM, K-means-BPSO-ELM, BPSO-GCSI-ELM, and K-means-GCSI-MBPSO-ELM, to compare with the proposed method. Results of 100 independent repeated experiments were taken from [16]; Table 6 shows the comparison results.Table 6
Comparisons between proposed feature dimensionality reduction method and six other methods.
Mean accuracy (%) + std (%) (number of genes)
Colon
Brain cancer
Lung cancer
Lymphoma
BPSO-ELM
93.34 ± 1.99 (9)
85.45 ± 2.33 (7)
94.80 ± 0.72 (11)
83.50 ± 2.72 (8)
K-Means-PSO-ELM
93.94 ± 1.17 (5)
86.55 ± 2.35 (5)
95.67 ± 0.72 (11)
83.72 ± 2.33 (6)
K-Means-BPSO-ELM
93.50 ± 2.02 (9)
87.23 ± 2.34 (8)
95.64 ± 0.56 (12)
85.14 ± 2.87 (6)
BPSO-GCSI-ELM
97.83 ± 1.33 (9)
89.89 ± 2.23 (7)
96.13 ± 0.42 (12)
—
K-Means-GCSI-MBPSO-ELM
97.61 ± 1.37 (6)
88.63 ± 2.16 (3)
97.10 ± 0.63 (11)
86.97 ± 2.44 (8)
Proposed method
99.08 ± 1.23 (13)
93.73 ± 2.64 (4)
98.67 ± 0.43 (11)
93.60 ± 2.33 (7)BPSO is a PSO algorithm that solves discrete optimization problems using particles formed by the binary system [35]. BPSO-ELM optimally selects disease-related genes using BPSO and evaluates the classification results of genes selected each time using ELM [36]. K-Means-PSO-ELM and K-means-BPSO-ELM each complete clustering analysis on all disease-related genes using K-means first followed by optimal selection with standard PSO and BPSO combined with the ELM classifier [37]. GCSI refers to gene-to-class sensitivity information; the BPSO-GCSI-ELM method generates an initial particle swarm using the GCSI of each gene and influences the updating and evolutionary process of particles [38]. Combining the above methods, K-means-GCSI-MBPSO-ELM first completes K-means clustering analysis on all disease-related genes and then introduces the GCSI into the modified BPSO to optimally select genes in combination with the ELM classifier [39].Because other gene extraction methods have demonstrated good effects on gene extraction of leukemia and SRBCT, this paper does not compare the results of these two types of diseases. Tables6 and 7 indicate that we can select genes with greater influences on tumor subtype classification using the proposed feature dimensionality reduction method. Gene sets selected with the proposed method exhibited high classification accuracy, substantially higher than other methods. The standard deviation of classification accuracy of the proposed method was smaller than that of the other methods except brain cancer. For the four diseases, the average accuracy of the proposed method was higher than the highest average accuracy of the other methods by 1.25%, 3.84%, 1.57%, and 6.63%, respectively (from left to right in Table 6). For brain cancer and lymphoma, the average classification accuracy improved to approximately 93%; for colon and lung cancer, the average classification accuracy reached approximately 99%. The classification accuracy was fairly high, greatly improving the practicability of gene-based disease detection.Table 7
Gene sets selected using the proposed feature dimensionality reduction method with classification accuracy.
Dataset
Selected gene sets
Accuracy of test set (%) ± std
5-fold CV accuracy (%) ± std
Leukemia
4050, 2642, 2121
100 ± 0.00
100 ± 0.00
4050, 2642, 1882
100 ± 0.00
100 ± 0.00
4050, 2642, 3258
100 ± 0.00
100 ± 0.00
2335, 2642, 1843, 4050
100 ± 0.00
100 ± 0.00
Brain cancer
3041, 3052, 3692
94.36 ± 0.043
92.47 ± 0.021
3052, 973, 3041, 3692, 4796
92.00 ± 0.023
91.78 ± 0.046
3041, 3052, 1014, 3692
92.42 ± 0.042
92.73 ± 0.026
7129, 2881, 3052, 865, 1970, 2935, 4871
92.78 ± 0.012
91.88 ± 0.019
Colon
493, 792, 739, 513, 280, 1346, 1334, 1954, 1110, 1563, 223, 175, 1668
97.13 ± 0.028
99.08 ± 0.012
493, 792, 377, 1976, 280, 1039, 1334, 698, 1110, 1563, 175, 1668, 1617
98.18 ± 0.024
98.34 ± 0.009
792, 1423, 14, 1976, 1909, 1110, 1589, 102, 107, 1916, 175, 1151, 251
93.73 ± 0.031
98.71 ± 0.013
792, 14, 1976, 765, 1909, 1524, 1110, 175, 43, 53, 1293, 1740, 251
96.86 ± 0.033
99.05 ± 0.011
SRBCT
742, 1003, 1954, 430, 2050, 123
100 ± 0.00
100 ± 0.00
545, 1955, 1434, 509, 971, 255
100 ± 0.00
100 ± 0.00
1003, 545, 1911, 153, 123, 1489, 2161
100 ± 0.00
100 ± 0.00
1955, 2050, 545, 2144, 2045, 123, 1489
100 ± 0.00
100 ± 0.00
Lung cancer
1765, 2779, 2841, 1474, 2045, 3191, 2763, 2817, 525, 1552, 1630
100 ± 0.00
98.27 ± 0.004
525, 1493, 607, 2763, 792, 580, 867, 368, 3279, 2158, 1225
100 ± 0.00
98.39 ± 0.003
1765, 883, 2763, 792, 580, 867, 985, 3279, 2988, 2045, 814
100 ± 0.00
98.67 ± 0.004
1765, 525, 2763, 2841, 1474, 2583, 867, 985, 2045, 814, 918
100 ± 0.00
98.67 ± 0.00
Lymphoma
5279, 4862, 6965, 2374, 1855, 2060
93.45 ± 0.025
93.39 ± 0.018
2828, 2437, 152, 1855, 4998, 2416
97.55 ± 0.028
92.62 ± 0.019
1855, 2828, 152, 2437, 80, 530, 1102
92.30 ± 0.041
92.36 ± 0.029
5279, 4687, 4940, 5449, 1133, 1855, 4519
90.50 ± 0.031
93.60 ± 0.023To compare the ELM classifier with other classifiers based on disease classification results, we selected the random forest classifier and support vector machine (SVM) classifier to identify tumor diseases based on the gene sets in Table7. We used the average and standard deviation of 5-fold CV accuracy from 100 tests to evaluate the classification results of different classifiers. Table 8 shows the results of the random forest classifier and SVM classifier.Table 8
Tumor disease recognition results of random forest classifier and SVM classifier.
Dataset
Gene sets selected
Random forest classifier (%) ± std
SVM classifier (%) ± std
Leukemia
4050, 2642, 2121
98.06 ± 0.72
100 ± 0.00
4050, 2642, 1882
95.56 ± 1.10
99.44 ± 0.72
4050, 2642, 3258
96.39 ± 1.34
98.75 ± 0.79
2335, 2642, 1843, 4050
96.53 ± 1.18
97.50 ± 0.88
Brain cancer
3041, 3052, 3692
82.65 ± 3.15
65.48 ± 5.18
3052, 973, 3041, 3692, 4796
81.75 ± 3.52
72.75 ± 3.88
3041, 3052, 1014, 3692
82.68 ± 3.28
75.83 ± 3.76
7129, 2881, 3052, 865, 1970, 2935, 4871
82.20 ± 3.42
84.95 ± 3.28
Colon
493, 792, 739, 513, 280, 1346, 1334, 1954, 1110, 1563, 223, 175, 1668
85.32 ± 2.99
94.52 ± 2.18
493, 792, 377, 1976, 280, 1039, 1334, 698, 1110, 1563, 175, 1668, 1617
79.52 ± 0.038
90.16 ± 2.46
792, 1423, 14, 1976, 1909, 1110, 1589, 102, 107, 1916, 175, 1151, 251
83.55 ± 5.03
93.55 ± 2.74
792, 14, 1976, 765, 1909, 1524, 1110, 175, 43, 53, 1293, 1740, 251
85.65 ± 3.52
94.35 ± 2.77
SRBCT
742, 1003, 1954, 430, 2050, 123
94.12 ± 2.17
99.26 ± 0.52
545, 1955, 1434, 509, 971, 255
95.41 ± 1.78
98.94 ± 1.04
1003, 545, 1911, 153, 123, 1489, 2161
94.72 ± 1.72
97.63 ± 1.15
1955, 2050, 545, 2144, 2045, 123, 1489
94.78 ± 1.64
98.01 ± 1.10
Lung cancer
1765, 2779, 2841, 1474, 2045, 3191, 2763, 2817, 525, 1552, 1630
94.91 ± 0.91
98.74 ± 0.33
525, 1493, 607, 2763, 792, 580, 867, 368, 3279, 2158, 1225
94.26 ± 1.12
98.55 ± 0.48
1765, 883, 2763, 792, 580, 867, 985, 3279, 2988, 2045, 814
94.83 ± 1.03
99.16 ± 0.31
1765, 525, 2763, 2841, 1474, 2583, 867, 985, 2045, 814, 918
94.52 ± 1.00
99.18 ± 0.43
Lymphoma
5279, 4862, 6965, 2374, 1855, 2060
79.22 ± 3.61
83.71 ± 2.67
2828, 2437, 152, 1855, 4998, 2416
73.19 ± 4.06
78.98 ± 3.03
1855, 2828, 152, 2437, 80, 530, 1102
72.76 ± 4.65
71.57 ± 3.63
5279, 4687, 4940, 5449, 1133, 1855, 4519
73.24 ± 4.50
75.47 ± 3.56In terms of tumor disease recognition accuracy based on the selected gene sets using the proposed method, the random forest classifier and SVM classifier were generally poorer than the ELM classifier. For lung cancer, however, the SVM classifier’s average accuracy was slightly higher than that of the ELM classifier and its standard deviation was smaller.


## 5. Conclusions

(1)
This paper proposes a new feature dimensionality reduction method based on a sampling survey method and heuristic intelligent optimization algorithm. First, we build a feature-scoring system and a reduced-dimension length-scoring system based on the sampling survey method and select features ranking in the front with reasonable dimension lengths. Then, we select features according to the improved heuristic intelligent optimization algorithm. Our method has only two steps with simple operation and low computational complexity.
(2)
We successfully apply the proposed feature dimensionality reduction method to road roughness time-domain estimation research. Results show that we only need as few as four dynamic responses to estimate road roughness accurately. The experimental results of 12 feature dimensionality reductions reveal the mean accuracy of road roughness time-domain estimation to be 0.9897.
(3)
We successfully apply the proposed feature dimensionality reduction method to gene-selection research in bioengineering. Results show that the gene sets selected with the proposed method have high classification accuracy with mean classification accuracy substantially higher than that of other gene-selection methods.



---

*Source: 2879640-2018-10-08.xml*