# Using Genetic Programming with Prior Formula Knowledge to Solve Symbolic Regression Problem

Authors: Qiang Lu; Jun Ren; Zhiguang Wang
Journal: Computational Intelligence and Neuroscience
 (2016)

---

## Abstract

A researcher can infer mathematical expressions of
functions quickly by using his professional knowledge (called
Prior Knowledge). But the results he finds may be biased and
restricted to his research field due to limitation of his knowledge.
In contrast, Genetic Programming method can discover fitted
mathematical expressions from the huge search space through
running evolutionary algorithms. And its results can be generalized
to accommodate different fields of knowledge. However,
sinceGP has to search a huge space, its speed of finding the
results is rather slow. Therefore, in this paper, a framework
of connection between Prior Formula Knowledge and GP (PFK-GP)
is proposed to reduce the space of GP searching. The PFK
is built based on the Deep Belief Network (DBN) which can
identify candidate formulas that are consistent with the features
of experimental data. By using these candidate formulas as the
seed of a randomly generated population, PFK-GP finds the right
formulas quickly by exploring the search space of data features.
We have compared PFK-GP with Pareto GP on regression of
eight benchmark problems. The experimental results confirm
that the PFK-GP can reduce the search space and obtain the
significant improvement in the quality of SR.


---

## Body

## 1. Introduction

Symbolic regression (SR) is used to discover mathematical expressions of functions that can fit the given data based on the rules of accuracy, simplicity, and generalization. As distinct from linear or nonlinear regression that efficiently optimizes the parameters in the prespecified model, SR tries to seek appropriate models and their parameters simultaneously for a purpose of getting better insights into the dataset. Without any prior knowledge of physics, kinematics, and geometry, some natural laws described by mathematical expressions, such as Hamiltonians, Lagrangians, and other laws of geometric and momentum conservation, can be distilled from experimental data by the Genetic Programming (GP) method on SR [1].Since SR is an NP-hard problem, some evolutionary algorithms were proposed to find approximate solutions to the problem, such as Genetic Programming (GP) [2], Gene Expression Programming (GEP) [3], Grammatical Evolution (GE) [4, 5], Analytic Programming (AP) [6], and Fast Evolutionary Programming (FEP) [7]. Moreover, recent researches in SR problem have taken into account machine learning (ML) algorithms [8–10]. All of the above algorithms randomly generate candidate population. But none of them can use various features of known functions to construct mathematical expressions adapted for describing the features of given data. Therefore, these algorithms may exploit huge search space that consists of all possible combinations of functions and its parameters.Nevertheless, a researcher always analyzes data, infers mathematical expressions, and obtains results according to his professional knowledge. After getting experimental data, he observes the data distribution and their features and analyzes them with his knowledge. Then, he tries to create some mathematical models based on natural laws. He can obtain the values of coefficients in these models through regression analysis methods or other mathematical methods. And he evaluates the formulas which are mathematical models with the values by using various fitness functions. If the researcher finds some of the formulas that fit the experimental data, he can transform and simplify these formulas and then obtain the final formula that can represent the data. Furthermore, his rich experience and knowledge can help him to reduce the searching space complexity so that he can find the best fit mathematical expression rapidly. As the researchers use their knowledge to discover the best fitted formulas, the methods that inject domain knowledge into the process of SR problem solving have been proposed to improve performance and scalability in complex problem [11–13]. The domain knowledge, which is manually created by the researcher’s intuition and experience, is of various formulas which are prior solutions to special problems. If the domain knowledge automatically generates some fitted formulas that are used in evolutionary search without the researcher involvement, the speed of solving SR problem will be quickened. A key challenge is how to build and utilize the domain knowledge just like the researcher does.In this paper, we present a framework of connection between Prior Formula Knowledge andGP (PFK-GP) to address the challenge:(i)
We classify a researcher’s domain knowledge into PFK Base (PFKB) and inference ability after analyzing the process that a research discovered formulas from experimental data (Section 2).PFKB contains two primary functions: classification and recognition. The aim of two functions is to generate feature functions which can represent the feature of experimental data.(ii)
In order to implement classification, we use the deep learning methodDBN [14, 15] which, compared with other shallow learning methods (Section 3.1), can classify experimental data into a special mathematical model that is consistent with data features. However, the classification method may lead to overfitting because the method can only categorize experimental data into known formula models which come from the set of training formula models.(iii)
Therefore, recognition is used to overcome the overfitting. It can extract mathematical models of functions that can show full or partial features of experimental data. Three algorithmsGenerateFs,CountSamePartF, andCountSpecU (see Algorithms 2, 3, and 4) are designed to implement recognition. For example, from the dataset generated by f(x)=exp⁡(sin⁡(x)+x3/(8∗105)), the basic functions sin, exp, and cube can be found by the above three algorithms. In Figure 1, the function sin shows the periodicity of data, and exp or cube shows the growth rate of data. Therefore, these basic functions (called feature functions) can describe some features of the dataset.(iv)
The inference ability is concluded to the searching ability of evolutionary algorithm. As researches infer mathematical models,GP is used to combine, transform, and verify these models. These feature functions that are generated byPFKB are selected to be combined into the candidate population in the light of algorithmrandomGenP (see Algorithm 5). With the candidate population,GP can get convergent result quickly because it searches answers in a limit space which is composed of various feature functions. Through experiment on eight benchmark problems (Table 5  E1–E8), the results demonstrate thatPFK-GP, compared with Pareto optimizationGP [16, 17], shows the significant improvement in accuracy and convergence.Figure 1
The functionf(x)=exp⁡(sin⁡(x)+x3/(8∗105)).


## 2. Background

### 2.1. Definition and Representation of Mathematical Expression

In this section, we will define concepts about SR problem and show how to represent these concepts by applyingBNF expression. For SR problem, the word “formula” is the general term which describes mathematical expression that fits the given data. We define a formula model is a special mathematical model in which formulas have the same relationships and variables except for different coefficient values. Relationships can be described by operators, such as algebraic operators, functions, and differential operators (http://en.wikipedia.org/wiki/Mathematical_model). Therefore, a formula model is a set where each element is a formula. For example, the two formulas 0.1∗sin⁡x+0.7∗log⁡(x) and 0.3∗sin⁡x+0.9∗log⁡(x) belong to the formula model a1∗sin⁡x+a2∗log⁡(x). Data that are represented by different formulas in one formula model may have similar features which are data distributions, data relationships between different variables, data change laws, and so on, because these formulas have the same relationships.In order to represent a formula model and its corresponding formulas, we define the followingBNF expressions:F≔C∣S,C≔S"("C")"∣S,S≔B”("AX,AX”)"∣U”("AX”)",B≔"+"∣"-"∣"∗"∣"/",U≔”sqrt”∣”log”∣”tanh”∣”sin”∣”cos”∣”exp”∣”tan”∣”abs”∣”quart” ∣ ”cube”∣”square”⋯,A≔a1∣a2∣a3∣⋯∣an,X≔x1∣x2∣x3∣⋯∣xm,whereF is a formula model. X is a parameter set. A indicates a coefficient set. B is a set of binary functions, while U is a set of unary functions. S is a set of atomic functions which does not contain any subfunctions. C is a set of complex functions which contains complex functions in C and atomic functions in S. With the above definitions, any formulas and its corresponding model can be shown by theseBNF expressions. For instance, the formula exp⁡(sin⁡(x)+x3/(8∗105)) is represented by F and C, and its subfunction sin⁡(x) is represented by U. The constants 8, 3, and 105 are shown by elements in A. With these BNF expressions, a formula model can be transformed into one tree. And the tree is a candidate individual in population ofGP solving SR problem. Every subtree in the tree is a subformula which can show special data features. A subtree that shows features of experimental data is called feature-subtree. If a tree has more feature-subtrees, the tree is more likely to fit the data. How to construct the tree consisting of feature-subtrees is a key step in our method which is implemented by the algorithmrandomGenP (see Algorithm 5).


### 2.2. The Process of Researcher Analyzing Data

The process that a researcher tries to solve SR problems is shown in Figure2. He depends heavily on his experience which is obtained through a long-term accumulation of study and research. After a researcher collected experimental data, he discovers regular patterns from data by using the methods of data analysis and visualization. He then constructs formula models which were consistent with these regular patterns according to his experiences. After that, he computes the coefficient values in formula models by using appropriate regression methods and obtains some formulas from different formula models. According to results of evaluating these formulas, he chooses the formula that is most fitted to the data. If the formula cannot represent data features, he needs to reselect a new formula model and do the above steps until one fitting formula is found.Figure 2
The process that researchers study the SR problem.We think the researcher’s experience and knowledge have two roles in processing SR problem. One role is Prior Formula Knowledge (PFK) which can help a researcher to quickly find fitted formulas that match experimental data features. Through study and work, the researcher accumulates his domain knowledge of various characteristics of formula model. When the researcher observes experimental data, he can apply his domain knowledge to recognize and classify the data. The other is the ability of inference and deduction which can help the researcher to combine, transform, and verify mathematical expression. We conclude that the PFK contains two primary functions: classification and recognition.Classification. when experimental data features are in accord with characteristics of one formula model in PFK, the dataset can be categorized into the model. The prerequisite of classification is that different formula models have different characteristics in PFK Base. As shown in Figure 3, six families of curves are generated by six formula models taking different coefficient values. The curves in the same family show similar data features while the curves in different families show different data features. Therefore, we can infer that the curves (including surfaces and hypersurfaces) generated by different formula models can be classified according to their data features.Figure 3
Curves are generated by six formula models with different coefficient value.Although many machine learning algorithms such as linear regression [18], SVM [19], Boosting [20], and PCVMs [21] can be used to identify and classify data, it is difficult for these algorithms to classify these curves. That is because these algorithms depend on features that are extracted manually from data, while these features from different complex curves are difficult to be represented by a feature vector which is built based on the researcher’s experiences. In contrast to these algorithms, DL can automatically extract features and have a good performance for the recognition of complex curves, such as image [15], speech [22], and natural language [23]. TheGenerateFs algorithm (see Algorithm 2) based onDBN is shown to classify the data.Recognition. Some formulas can represent remark features of curves generated by formula model. For example, after observing the curve in Figure 1, a researcher can easily infer that the formula sin or cos is one of formulas that constitute the curve because data in curve show periodicity. Therefore, these formulas are called feature functions that can be recognized or extracted by PFK. AlgorithmsCountSamePartF andCountSpecU (see Algorithms 3 and 4) are built to recognize the feature functions.Recognition can help the researcher overcome overfitting of results that are generated by classification because classification can help researcher to only identify formula models from training set while recognition can help the researcher identify subformula models that are consistent with local data features.The ability of inference and deduction is one of main measurements for evaluating performance of artificial intelligence methods. In the SR problem,GP, compared with other methods such as logical reasoning, statistical learning, and genetic algorithm, is a revolutionary method of searching fitting formulas because it can seek the appropriate formula models and their coefficient values simultaneously by evolving the population of formula individuals. Therefore, in the paper, we useGP as the method of inferring and deducing formulas.To optimizeGP, researchers have proposed various approaches, such as optimal parsimony pressure [24]; Pareto front optimization [17] and its age-fitness method [25] are used to control bloat and premature convergence inGP. In order to reduce the space complexity of searching formulas, the methods of arithmetic [26] and machine learning are injected intoGP. In the paper, with the algorithmrandomGenP (see Algorithm 5) about generating population and the method of Pareto front optimization,PFK-GP can research the formula model in the appropriate space and can find right formulas quickly.


## 2.1. Definition and Representation of Mathematical Expression

In this section, we will define concepts about SR problem and show how to represent these concepts by applyingBNF expression. For SR problem, the word “formula” is the general term which describes mathematical expression that fits the given data. We define a formula model is a special mathematical model in which formulas have the same relationships and variables except for different coefficient values. Relationships can be described by operators, such as algebraic operators, functions, and differential operators (http://en.wikipedia.org/wiki/Mathematical_model). Therefore, a formula model is a set where each element is a formula. For example, the two formulas 0.1∗sin⁡x+0.7∗log⁡(x) and 0.3∗sin⁡x+0.9∗log⁡(x) belong to the formula model a1∗sin⁡x+a2∗log⁡(x). Data that are represented by different formulas in one formula model may have similar features which are data distributions, data relationships between different variables, data change laws, and so on, because these formulas have the same relationships.In order to represent a formula model and its corresponding formulas, we define the followingBNF expressions:F≔C∣S,C≔S"("C")"∣S,S≔B”("AX,AX”)"∣U”("AX”)",B≔"+"∣"-"∣"∗"∣"/",U≔”sqrt”∣”log”∣”tanh”∣”sin”∣”cos”∣”exp”∣”tan”∣”abs”∣”quart” ∣ ”cube”∣”square”⋯,A≔a1∣a2∣a3∣⋯∣an,X≔x1∣x2∣x3∣⋯∣xm,whereF is a formula model. X is a parameter set. A indicates a coefficient set. B is a set of binary functions, while U is a set of unary functions. S is a set of atomic functions which does not contain any subfunctions. C is a set of complex functions which contains complex functions in C and atomic functions in S. With the above definitions, any formulas and its corresponding model can be shown by theseBNF expressions. For instance, the formula exp⁡(sin⁡(x)+x3/(8∗105)) is represented by F and C, and its subfunction sin⁡(x) is represented by U. The constants 8, 3, and 105 are shown by elements in A. With these BNF expressions, a formula model can be transformed into one tree. And the tree is a candidate individual in population ofGP solving SR problem. Every subtree in the tree is a subformula which can show special data features. A subtree that shows features of experimental data is called feature-subtree. If a tree has more feature-subtrees, the tree is more likely to fit the data. How to construct the tree consisting of feature-subtrees is a key step in our method which is implemented by the algorithmrandomGenP (see Algorithm 5).


## 2.2. The Process of Researcher Analyzing Data

The process that a researcher tries to solve SR problems is shown in Figure2. He depends heavily on his experience which is obtained through a long-term accumulation of study and research. After a researcher collected experimental data, he discovers regular patterns from data by using the methods of data analysis and visualization. He then constructs formula models which were consistent with these regular patterns according to his experiences. After that, he computes the coefficient values in formula models by using appropriate regression methods and obtains some formulas from different formula models. According to results of evaluating these formulas, he chooses the formula that is most fitted to the data. If the formula cannot represent data features, he needs to reselect a new formula model and do the above steps until one fitting formula is found.Figure 2
The process that researchers study the SR problem.We think the researcher’s experience and knowledge have two roles in processing SR problem. One role is Prior Formula Knowledge (PFK) which can help a researcher to quickly find fitted formulas that match experimental data features. Through study and work, the researcher accumulates his domain knowledge of various characteristics of formula model. When the researcher observes experimental data, he can apply his domain knowledge to recognize and classify the data. The other is the ability of inference and deduction which can help the researcher to combine, transform, and verify mathematical expression. We conclude that the PFK contains two primary functions: classification and recognition.Classification. when experimental data features are in accord with characteristics of one formula model in PFK, the dataset can be categorized into the model. The prerequisite of classification is that different formula models have different characteristics in PFK Base. As shown in Figure 3, six families of curves are generated by six formula models taking different coefficient values. The curves in the same family show similar data features while the curves in different families show different data features. Therefore, we can infer that the curves (including surfaces and hypersurfaces) generated by different formula models can be classified according to their data features.Figure 3
Curves are generated by six formula models with different coefficient value.Although many machine learning algorithms such as linear regression [18], SVM [19], Boosting [20], and PCVMs [21] can be used to identify and classify data, it is difficult for these algorithms to classify these curves. That is because these algorithms depend on features that are extracted manually from data, while these features from different complex curves are difficult to be represented by a feature vector which is built based on the researcher’s experiences. In contrast to these algorithms, DL can automatically extract features and have a good performance for the recognition of complex curves, such as image [15], speech [22], and natural language [23]. TheGenerateFs algorithm (see Algorithm 2) based onDBN is shown to classify the data.Recognition. Some formulas can represent remark features of curves generated by formula model. For example, after observing the curve in Figure 1, a researcher can easily infer that the formula sin or cos is one of formulas that constitute the curve because data in curve show periodicity. Therefore, these formulas are called feature functions that can be recognized or extracted by PFK. AlgorithmsCountSamePartF andCountSpecU (see Algorithms 3 and 4) are built to recognize the feature functions.Recognition can help the researcher overcome overfitting of results that are generated by classification because classification can help researcher to only identify formula models from training set while recognition can help the researcher identify subformula models that are consistent with local data features.The ability of inference and deduction is one of main measurements for evaluating performance of artificial intelligence methods. In the SR problem,GP, compared with other methods such as logical reasoning, statistical learning, and genetic algorithm, is a revolutionary method of searching fitting formulas because it can seek the appropriate formula models and their coefficient values simultaneously by evolving the population of formula individuals. Therefore, in the paper, we useGP as the method of inferring and deducing formulas.To optimizeGP, researchers have proposed various approaches, such as optimal parsimony pressure [24]; Pareto front optimization [17] and its age-fitness method [25] are used to control bloat and premature convergence inGP. In order to reduce the space complexity of searching formulas, the methods of arithmetic [26] and machine learning are injected intoGP. In the paper, with the algorithmrandomGenP (see Algorithm 5) about generating population and the method of Pareto front optimization,PFK-GP can research the formula model in the appropriate space and can find right formulas quickly.


## 3. Genetic Programming with Prior Formula Knowledge Base

### 3.1. Formula Prior Knowledge Base

The FPK needs to have the ability of identifying and classifying the formula modelF based on data features. Although the features between formula models are different, it is difficult to extract features from data which are generated by these models because different formula models represent seemly similar but different features. Based on the above definitions in formula model, the features among functions in set S are different. The features between the function s∈S and the function c∈C may be similar if c is the parameter of s. As shown in Figure 4, these functions sin⁡(x), cube(x), and exp⁡(x) belonging to S constitute C, such as exp⁡(sin⁡(x)) and sin⁡(x)+cube(x)/(8∗105), and are shaped into the final function exp⁡(sin⁡(x)+cube(x)/(8∗105)). sin⁡(x) shows periodicity in results; log⁡(x) shows slow variation trend in results; cube(x) shows high variation trend in results. So, these features in the three functions are different. However, there are similar features between cube(x) and sin⁡(x)+cube(x)/(8∗105), because the two functions have high variation trend in result. The shallow learning method such as SVM can hardly identify complex features of functions shown in Figure 5.Figure 4
Illustration of theDBN framework.Figure 5
Accuracy result of SVM andDBN classifying P1–P39.In this paper,DBN is used to classify data into a special formula model according to features that are automatically extracted from data. Generally,DBN is made up of multiple layers that are represented by Restricted Boltzmann Machine (RBM). As RBM is an universal approximate of discrete models [27], we can speculate thatDBN which has many layers of RBM can recognize the features of data generated by complex function just like Convolutional Neural Network (CNN) classifying image [28]. The process ofDBN recognizing formulas is illustrated in Figure 4.DBN can extract features from data, layer by layer. So, the lower RBM layer can represent the simple functions and the higher RBM layer can transform the simple functions into more complex functions.We use the data generated by the different formula modelsF (see Table 7) as training samples.DBN is trained by these samples. The model that is finally gained byDBN training methods isPFKB, which is aimed at identifying the formula model that can represent features of the data. The process ofDBN training is outlined in algorithmTrainDBN(Algorithm 1) which uses the same steps as mentioned in [14, 15].Algorithm 1:
TrainDBN: training DBN to generate the PFKB.
Input: X1, Y1, X2, Y2 (X1, Y1 are training data; X2, Y2 are testing data)
Output: PFKB
(1) Initial(DL, opts) // initial the structure of DL and the parameters opts
(2) DL = DLsetup(DBN, X1, opts) // layer-wise pre-training DL
(3) DL = DLtrain(DBN,X1, opts) // build up each layer of DL to train
(4) nn=DLunfoldtonn(DBN, opts) // after training each layer, passing the parameters to nn
(5) PFKB = nntrain(nn,X1,Y1, opts) // fine-tune the whole deep architecture
(6) accuracy = nntest(PFKB,X2, Y2) // accuracy is the criterion of the quality of PFKM, if it is too
small, then re-training after adjusting the model architecture or parameters
(7) return PFKBAlgorithm 2:
G
e
n
e
r
a
t
e
F
s.
Input: X, PFKB, s (X is the dataset; s is the number of formula models whose features can show the dataset)
Output:   Fs (formula model vector which are used in generating the initial population of GP)
(1)  Ftemp = predictModel(PFKB, X) // as the intermediate data, it is the original result that exploit the PFKB
to predict without sorting by the fitness
(2)  Ftemp = sortModelByFit(Ftemp) // sort the models in order of decreasing fitness
(3) fori = 1 : s
(4)  Fs(i) = Ftemp(i)
(5) end
(6) returnFsAlgorithm 3:
C
o
u
n
t
S
a
m
e
P
a
r
t
F.
Input: t, Fs
Output:  C (ordered local expressions set that are sorted according to frequency of occurrence which is larger than t)
(1)    C=F=∅
(2)    for each  pairfi,fjinFs
(3)     Fij=fi∩fj
(4)     for eachcminFij
(5)      if  cm∈F
(6)       changeF’s elementcmvtocmv+1 // v indicates the number of times that cm appears
(7)      else
(8)       addcm1intoF
(9)     end
(10) end
(11) for eachcmvinF
(12)  if  v≥t
(13)   add  cmvintoC
(14) end
(15) sort(C)
(16) return CAlgorithm 4:
C
o
u
n
t
S
p
e
c
U.
Input: t, Fs
Output:  specU (ordered specU function set that are sorted according to frequency of occurrence which is larger than t)
(1)    U=F=∅
(2)    for each  fi  in  Fs
(3)     for eachuminfi
(4)      if  um∈F
(5)       changeF’s elementumvtoumv+1 // v indicates the number of times that um appears
(6)      else
(7)       addum1intoF
(8)     end
(9)    end
(10) for eachumvinF
(11)  if  v≥t
(12)   add  umvintospecU
(13) end
(14) returnspecUAlgorithm 5:
r
a
n
d
o
m
G
e
n
P.
Input:  Fs, C, specU, B, U,n where, n represents the number of individuals which are generated, B and U is the candidate
function library)
Output:  P (population)
(1) I=Fs∪C
(2) for eachiinI
(3)  addiintoP
(4) end
(5) k=specU
(6) Q=specU+B+U // add the elements of specU, B and U into the queue Q successively, umv∈U, add umv into Q
(7) P_temp = traditionalRandomIndividual(Q, k)
(8) addP_temp intoP
(9) k=n-I-k
(10) P_temp = traditionalRandomIndividual(B+U,k)
(11) addP_temp intoP
(12) returnPPFKB is only changed with formula models. If there are no new trained formula models in an application, the algorithmTrainDBN will not be executed. When the number of trained formula models is large enough, little new formula model will appear, andPFKB will seldom be changed. In the paper,TrainDBN is performed exactly once in order to generatePFKB.


### 3.2. Classification and Recognition withPFKB

In order to deal with the problem of how to classify and recognize formula model from data, we should consider the problem from two aspects. One situation is that data can be represented by a special formula model fromPFKB, while the other one is that data cannot be represented by a formula model fromPFKB. In the first case, we exploitPFKB to identify formula models of data byDBN classification. Based on ordered results ofDBN classification, we gain a set of formula models (Fs=f1,…,fs) which are most similar to features of the data. The process that deals with the first case is outlined in algorithmGenerateFs. The algorithm is fast becausePFKB has been built byTrainDBN, and s is small integer value.In the second case, when a researcher observes laws that are hidden in experimental data, he often tries to find some formulasC which are consistent with partial features of the data. Therefore, we propose the two assumptions as follows.Assumption 1.
More formula modelsfs have the same subformula modelpf in the setFs which is the result ofGenerateFs running, more strongly that thepf can express features of data.
In order to compute the samepf inFs, we express the formula model as the string of expression and seek the same part of them by using intersection between the two strings (without considering the elements in setsX andA). Define the intersection between two expressions as follows:(1)fi∩fj=c1,…,ck,cm∩cn=∅,m≠n,cn∈C,cn∉S,1≤m,n≤k,1≤i,j≤k.For example, f1=z+a∗cos⁡x+tan⁡x/(exp⁡x+log⁡x), f2=z+a∗cos⁡x+abs(x)/(exp⁡x+log⁡(x)), f1∩f2={z+a∗cos⁡x, exp⁡(x)+log⁡(x)}. The method, which obtainspf whose frequency of occurrence inf is larger than thresholdt, is described as the algorithmCountSamePartF.For verifying Assumption1, we apply the dataset from E0 (see Table 5) as the testing set and get the identifying results Fs={P18,P15,P34,P11,P3} from Table 7 through the algorithmGenerateFs. The intersections between the top two formulas P18∩P15 are (sqrt(x0))/(tan⁡(x1)) and (tan⁡x0)/(exp⁡(x1)), which are partial mathematical expressions in E0. And we use the dataset from Table 6 to test E0 and get the identifying results Fs={T7,T6,T8,T9,T1} throughGenerateFs algorithm. The intersection between two expressions is as follows: T7∩T6={sqrt(x0)/tan⁡(x1)}, T8∩T9={tan⁡x0/exp⁡(x1), cos⁡(x0)∗sin⁡(x1)}. We find that the elements which have more frequency of occurrence in the intersections set are more likely to express some data features. The above two experiments illustrate that Assumption 1 is rational.Assumption 2.
If functionu∈U exists inFs obtained byGenerateFs and the number of the same u is larger than threshold t, we can conclude that u can show some local data features.
The functionb∈B except xy is common function, which has a high probability of occurrence in mathematical expressions. Therefore, it is difficult to express special data features. Compared with B, the function u∈U can show obvious features of data. For instance, sin⁡(x) presents the periodicity of data and log⁡(x) represents data features about extreme increase or decrease. The method, which obtains the special function u that can show the local data features, is outlined as the algorithmCountSpecU.For verifying Assumption2, we also choose the dataset which are generated from E0 (see Table 5) as the testing data and apply theCountSpecU algorithm to calculate the special u among Fs={P18,P15,P34,P11,P3}. The result of the algorithm is shown in Table 1. We find the result specU={tan,cos,sqrt,exp,sin} (sin and cos are the operators of the same kind) is part of E0. Hence, we can discover that theu set, which is gained by the algorithmCountSpecU, can show local features of the dataset.Table 1
The result ofU in Fs computed by algorithm CountSpecU  (see Algorithm 4).
u in Fs
Frequency of occurrence
tan
3
cos
2
sqrt
1
exp
1
log
1


### 3.3.GP with Prior Formula Knowledge Base

In order to deal with SR problem,GP is executed to automatically composite and evolve mathematical expression. The process ofGP is similar to the process that a researcher transforms formula models and obtains fitting formulas based on his knowledge. Since those algorithms inPFKB, which is created based on analyzing the process of how a research infers fitted formulas, can recognize formula models that are consistent with data features, we combine these formula models ofPFKB recognizing into the process ofGP in order to reduce the searching space and increase the speed of discovering right solutions.When initializingGP algorithm, we can select candidate formulas fromFs,C, andspecU as individuals in population ofGP. The setsFs,C, andspecU are gained by the above algorithms inPFKB. Therefore, thePFKB is injected into the process of population generating. And this population can contribute to reserving data features as much as possible and reducing the searching space because these individuals commonly have good fitness value. With the population,GP algorithm can speed up the convergence and improve the accuracy of SR results. However, it may lead to the bias results. To overcome the problem, some random individuals must be imported into the population. The process of population creating is as follows.Firstly, the elements in setsFs and C are inserted into the population. Then, the setspecU and the candidate function sets B andU are merged into the new candidate function queueQ. And the number of elements inspecU is twice as much as the other elements inQ because B∪C⊆specU. Those elements inspecU are more likely to be part of individuals in the population after applying the methodtraditionalRandomIndividual [16] which is designed to generate randomly k individuals from the special function set. At last, the rest of individuals of population are created bytraditionalRandomIndividual with setsB andU. The process of population generating is described as the algorithmrandomGenP.Generally,|Fs|+|C|+|specU|<n/2, where n is the number of individuals in population. Furthermore, in order to enhance the affection ofPFKB in the process ofGP evolution, the methodrandomGenP is used to create new individuals in every few generations of evolutionary computation. Meanwhile, the method of Pareto front [17] is introduced into the algorithmPFK-GP to balance the accuracy against the complexity of model. The detail of algorithmPFK-GP is shown in Algorithm 6.Algorithm 6:
PFK-GP.
Input:  data, PFKB, t1, t2, B, U, n, k, g, interval
Output: F (candidate formulas set)
(1)    Fs=GenerateFs(data, PFKB)
(2)    C = CountSamePartF(t1, Fs)
(3)    U=CountSpecU(t2, Fs)
(4)    P=randomGenP(Fs, C, specU, B+U, n)
(5)    while (bestFitness  <=  threshold&&  i<g)
(6)     P = crossover(P)
(7)     P = mutate(P)
(8)     Pt = ParetoOptimise(P) // prevent the formula model too complex
(9)     Pt_fitnees = EvaluatePopulation(Pt)
(10)  bestFitness,F = Selectbest(Pt, Pt_fitness, k) // choose the best k individuals and get the best fitness value from the individuals
(11)  ifimod interval
(12)    P1=randomGenP(F, C, specU, B+U,n/2)
(13)    P2 = traditionalRandomIndividual(B+U,n/2)
(14)    P=P1∪P2
(15)  else
(16)   P1 = traditionalRandomIndividual(B+U,n-k)
(17)   P=P1∪P
(18)  end
(19)  i++
(20) end
(21) returnF


## 3.1. Formula Prior Knowledge Base

The FPK needs to have the ability of identifying and classifying the formula modelF based on data features. Although the features between formula models are different, it is difficult to extract features from data which are generated by these models because different formula models represent seemly similar but different features. Based on the above definitions in formula model, the features among functions in set S are different. The features between the function s∈S and the function c∈C may be similar if c is the parameter of s. As shown in Figure 4, these functions sin⁡(x), cube(x), and exp⁡(x) belonging to S constitute C, such as exp⁡(sin⁡(x)) and sin⁡(x)+cube(x)/(8∗105), and are shaped into the final function exp⁡(sin⁡(x)+cube(x)/(8∗105)). sin⁡(x) shows periodicity in results; log⁡(x) shows slow variation trend in results; cube(x) shows high variation trend in results. So, these features in the three functions are different. However, there are similar features between cube(x) and sin⁡(x)+cube(x)/(8∗105), because the two functions have high variation trend in result. The shallow learning method such as SVM can hardly identify complex features of functions shown in Figure 5.Figure 4
Illustration of theDBN framework.Figure 5
Accuracy result of SVM andDBN classifying P1–P39.In this paper,DBN is used to classify data into a special formula model according to features that are automatically extracted from data. Generally,DBN is made up of multiple layers that are represented by Restricted Boltzmann Machine (RBM). As RBM is an universal approximate of discrete models [27], we can speculate thatDBN which has many layers of RBM can recognize the features of data generated by complex function just like Convolutional Neural Network (CNN) classifying image [28]. The process ofDBN recognizing formulas is illustrated in Figure 4.DBN can extract features from data, layer by layer. So, the lower RBM layer can represent the simple functions and the higher RBM layer can transform the simple functions into more complex functions.We use the data generated by the different formula modelsF (see Table 7) as training samples.DBN is trained by these samples. The model that is finally gained byDBN training methods isPFKB, which is aimed at identifying the formula model that can represent features of the data. The process ofDBN training is outlined in algorithmTrainDBN(Algorithm 1) which uses the same steps as mentioned in [14, 15].Algorithm 1:
TrainDBN: training DBN to generate the PFKB.
Input: X1, Y1, X2, Y2 (X1, Y1 are training data; X2, Y2 are testing data)
Output: PFKB
(1) Initial(DL, opts) // initial the structure of DL and the parameters opts
(2) DL = DLsetup(DBN, X1, opts) // layer-wise pre-training DL
(3) DL = DLtrain(DBN,X1, opts) // build up each layer of DL to train
(4) nn=DLunfoldtonn(DBN, opts) // after training each layer, passing the parameters to nn
(5) PFKB = nntrain(nn,X1,Y1, opts) // fine-tune the whole deep architecture
(6) accuracy = nntest(PFKB,X2, Y2) // accuracy is the criterion of the quality of PFKM, if it is too
small, then re-training after adjusting the model architecture or parameters
(7) return PFKBAlgorithm 2:
G
e
n
e
r
a
t
e
F
s.
Input: X, PFKB, s (X is the dataset; s is the number of formula models whose features can show the dataset)
Output:   Fs (formula model vector which are used in generating the initial population of GP)
(1)  Ftemp = predictModel(PFKB, X) // as the intermediate data, it is the original result that exploit the PFKB
to predict without sorting by the fitness
(2)  Ftemp = sortModelByFit(Ftemp) // sort the models in order of decreasing fitness
(3) fori = 1 : s
(4)  Fs(i) = Ftemp(i)
(5) end
(6) returnFsAlgorithm 3:
C
o
u
n
t
S
a
m
e
P
a
r
t
F.
Input: t, Fs
Output:  C (ordered local expressions set that are sorted according to frequency of occurrence which is larger than t)
(1)    C=F=∅
(2)    for each  pairfi,fjinFs
(3)     Fij=fi∩fj
(4)     for eachcminFij
(5)      if  cm∈F
(6)       changeF’s elementcmvtocmv+1 // v indicates the number of times that cm appears
(7)      else
(8)       addcm1intoF
(9)     end
(10) end
(11) for eachcmvinF
(12)  if  v≥t
(13)   add  cmvintoC
(14) end
(15) sort(C)
(16) return CAlgorithm 4:
C
o
u
n
t
S
p
e
c
U.
Input: t, Fs
Output:  specU (ordered specU function set that are sorted according to frequency of occurrence which is larger than t)
(1)    U=F=∅
(2)    for each  fi  in  Fs
(3)     for eachuminfi
(4)      if  um∈F
(5)       changeF’s elementumvtoumv+1 // v indicates the number of times that um appears
(6)      else
(7)       addum1intoF
(8)     end
(9)    end
(10) for eachumvinF
(11)  if  v≥t
(12)   add  umvintospecU
(13) end
(14) returnspecUAlgorithm 5:
r
a
n
d
o
m
G
e
n
P.
Input:  Fs, C, specU, B, U,n where, n represents the number of individuals which are generated, B and U is the candidate
function library)
Output:  P (population)
(1) I=Fs∪C
(2) for eachiinI
(3)  addiintoP
(4) end
(5) k=specU
(6) Q=specU+B+U // add the elements of specU, B and U into the queue Q successively, umv∈U, add umv into Q
(7) P_temp = traditionalRandomIndividual(Q, k)
(8) addP_temp intoP
(9) k=n-I-k
(10) P_temp = traditionalRandomIndividual(B+U,k)
(11) addP_temp intoP
(12) returnPPFKB is only changed with formula models. If there are no new trained formula models in an application, the algorithmTrainDBN will not be executed. When the number of trained formula models is large enough, little new formula model will appear, andPFKB will seldom be changed. In the paper,TrainDBN is performed exactly once in order to generatePFKB.


## 3.2. Classification and Recognition withPFKB

In order to deal with the problem of how to classify and recognize formula model from data, we should consider the problem from two aspects. One situation is that data can be represented by a special formula model fromPFKB, while the other one is that data cannot be represented by a formula model fromPFKB. In the first case, we exploitPFKB to identify formula models of data byDBN classification. Based on ordered results ofDBN classification, we gain a set of formula models (Fs=f1,…,fs) which are most similar to features of the data. The process that deals with the first case is outlined in algorithmGenerateFs. The algorithm is fast becausePFKB has been built byTrainDBN, and s is small integer value.In the second case, when a researcher observes laws that are hidden in experimental data, he often tries to find some formulasC which are consistent with partial features of the data. Therefore, we propose the two assumptions as follows.Assumption 1.
More formula modelsfs have the same subformula modelpf in the setFs which is the result ofGenerateFs running, more strongly that thepf can express features of data.
In order to compute the samepf inFs, we express the formula model as the string of expression and seek the same part of them by using intersection between the two strings (without considering the elements in setsX andA). Define the intersection between two expressions as follows:(1)fi∩fj=c1,…,ck,cm∩cn=∅,m≠n,cn∈C,cn∉S,1≤m,n≤k,1≤i,j≤k.For example, f1=z+a∗cos⁡x+tan⁡x/(exp⁡x+log⁡x), f2=z+a∗cos⁡x+abs(x)/(exp⁡x+log⁡(x)), f1∩f2={z+a∗cos⁡x, exp⁡(x)+log⁡(x)}. The method, which obtainspf whose frequency of occurrence inf is larger than thresholdt, is described as the algorithmCountSamePartF.For verifying Assumption1, we apply the dataset from E0 (see Table 5) as the testing set and get the identifying results Fs={P18,P15,P34,P11,P3} from Table 7 through the algorithmGenerateFs. The intersections between the top two formulas P18∩P15 are (sqrt(x0))/(tan⁡(x1)) and (tan⁡x0)/(exp⁡(x1)), which are partial mathematical expressions in E0. And we use the dataset from Table 6 to test E0 and get the identifying results Fs={T7,T6,T8,T9,T1} throughGenerateFs algorithm. The intersection between two expressions is as follows: T7∩T6={sqrt(x0)/tan⁡(x1)}, T8∩T9={tan⁡x0/exp⁡(x1), cos⁡(x0)∗sin⁡(x1)}. We find that the elements which have more frequency of occurrence in the intersections set are more likely to express some data features. The above two experiments illustrate that Assumption 1 is rational.Assumption 2.
If functionu∈U exists inFs obtained byGenerateFs and the number of the same u is larger than threshold t, we can conclude that u can show some local data features.
The functionb∈B except xy is common function, which has a high probability of occurrence in mathematical expressions. Therefore, it is difficult to express special data features. Compared with B, the function u∈U can show obvious features of data. For instance, sin⁡(x) presents the periodicity of data and log⁡(x) represents data features about extreme increase or decrease. The method, which obtains the special function u that can show the local data features, is outlined as the algorithmCountSpecU.For verifying Assumption2, we also choose the dataset which are generated from E0 (see Table 5) as the testing data and apply theCountSpecU algorithm to calculate the special u among Fs={P18,P15,P34,P11,P3}. The result of the algorithm is shown in Table 1. We find the result specU={tan,cos,sqrt,exp,sin} (sin and cos are the operators of the same kind) is part of E0. Hence, we can discover that theu set, which is gained by the algorithmCountSpecU, can show local features of the dataset.Table 1
The result ofU in Fs computed by algorithm CountSpecU  (see Algorithm 4).
u in Fs
Frequency of occurrence
tan
3
cos
2
sqrt
1
exp
1
log
1


## 3.3.GP with Prior Formula Knowledge Base

In order to deal with SR problem,GP is executed to automatically composite and evolve mathematical expression. The process ofGP is similar to the process that a researcher transforms formula models and obtains fitting formulas based on his knowledge. Since those algorithms inPFKB, which is created based on analyzing the process of how a research infers fitted formulas, can recognize formula models that are consistent with data features, we combine these formula models ofPFKB recognizing into the process ofGP in order to reduce the searching space and increase the speed of discovering right solutions.When initializingGP algorithm, we can select candidate formulas fromFs,C, andspecU as individuals in population ofGP. The setsFs,C, andspecU are gained by the above algorithms inPFKB. Therefore, thePFKB is injected into the process of population generating. And this population can contribute to reserving data features as much as possible and reducing the searching space because these individuals commonly have good fitness value. With the population,GP algorithm can speed up the convergence and improve the accuracy of SR results. However, it may lead to the bias results. To overcome the problem, some random individuals must be imported into the population. The process of population creating is as follows.Firstly, the elements in setsFs and C are inserted into the population. Then, the setspecU and the candidate function sets B andU are merged into the new candidate function queueQ. And the number of elements inspecU is twice as much as the other elements inQ because B∪C⊆specU. Those elements inspecU are more likely to be part of individuals in the population after applying the methodtraditionalRandomIndividual [16] which is designed to generate randomly k individuals from the special function set. At last, the rest of individuals of population are created bytraditionalRandomIndividual with setsB andU. The process of population generating is described as the algorithmrandomGenP.Generally,|Fs|+|C|+|specU|<n/2, where n is the number of individuals in population. Furthermore, in order to enhance the affection ofPFKB in the process ofGP evolution, the methodrandomGenP is used to create new individuals in every few generations of evolutionary computation. Meanwhile, the method of Pareto front [17] is introduced into the algorithmPFK-GP to balance the accuracy against the complexity of model. The detail of algorithmPFK-GP is shown in Algorithm 6.Algorithm 6:
PFK-GP.
Input:  data, PFKB, t1, t2, B, U, n, k, g, interval
Output: F (candidate formulas set)
(1)    Fs=GenerateFs(data, PFKB)
(2)    C = CountSamePartF(t1, Fs)
(3)    U=CountSpecU(t2, Fs)
(4)    P=randomGenP(Fs, C, specU, B+U, n)
(5)    while (bestFitness  <=  threshold&&  i<g)
(6)     P = crossover(P)
(7)     P = mutate(P)
(8)     Pt = ParetoOptimise(P) // prevent the formula model too complex
(9)     Pt_fitnees = EvaluatePopulation(Pt)
(10)  bestFitness,F = Selectbest(Pt, Pt_fitness, k) // choose the best k individuals and get the best fitness value from the individuals
(11)  ifimod interval
(12)    P1=randomGenP(F, C, specU, B+U,n/2)
(13)    P2 = traditionalRandomIndividual(B+U,n/2)
(14)    P=P1∪P2
(15)  else
(16)   P1 = traditionalRandomIndividual(B+U,n-k)
(17)   P=P1∪P
(18)  end
(19)  i++
(20) end
(21) returnF


## 4. Experiments

In the experiments, we employDBN in the DeepLearnToolbox [30] to classify formula models and build the algorithmPFK-GP based on GPTIPS [29]. The 39 formula models in Table 7 are composed of formulas from [31, 32] and some formulas are created by ourselves. The data generated by these 39 formula models is used as training data of algorithmDBN to createPFKB. The formula models in Table 5 are used to generate the testing data for verifying accuracy of algorithmsGenerateFs andPFK-GP. The formula models in Table 6 are devoted to validating the two algorithmsCountSamePartF andCountSpecU (see Algorithms 3 and 4).For most formula models from Tables5, 6, and 7, we sampled them by equal step taking their parameter values from the range [−49, 50]. For some particular formulas, we also sample them with a special equal step from special numerical scope. For example, the value x in sqrt(x) is in the range [0, 99], the value x in log⁡(x) ranges between 1 and 100. We create 500 groups of different parameters value in each formula model. The coefficients in these formula models are fetched with equal step from the range [−2.996  3.0]. When all coefficients of a formula model take special values, the formula model generates a formula, namely, a sample of the formula model. We create 7500 groups of different coefficients in each formula model. So, each formula model has 7500 samples where each sample has 500 groups of different parameters value. We take 6000 samples of these samples as training data and the others as test data.We adoptDBN as the classification model and compare it with SVM that is implemented by the tool libsvm [33]. The training and testing data for the two algorithms are originated from formula models P1–P39. The parameter values inDBN and SVM are illustrated as Table 2. We take the first five formulas fromFs generated byGenerateFs as a result set of recognition. If the test formula is included in the set, we think that the recognition result is correct. The accuracy of recognition results ofDBN and SVM is showed as Figure 5. TheDBN method can help to classify all kinds of test data into its fitted formula models. However, the SVM method can only correctly classify several kinds of test data. The overall average accuracy ofDBN classification is 99.65%, while the accuracy of SVM is 26.72%. The result demonstrates thatDBN is more suitable for recognizing data generated by mathematical expression, becauseDBN can automatically extract features from the data, layer by layer, and is similar to composition of formula which is constituted by its subformulas.Table 2
The parameter values in algorithm ofDBN and SVM.
D
B
N parameters
Value
SVM parameters
Value
The number ofDBN layers
4
svm_type
c-svc
The size ofDBN hide nodes
50
Kernel
Gaussian
The number of epochs
200
Gamma
0.07
Batch size
40
Coef
0
Momentum
0
Cost
1.0
Alpha
1
Degree
3.0
activation_function
sigm
Shrinking
1We set parameters inGP with Pareto optimization (PO-GP) [29] andPFK-GP as shown in Table 3. For data generated by P13 (see Table 7, coefficient z is −2.098, a is −2.998),PO-GP andPFK-GP deal with the SR problem, respectively. The result was illustrated as Figure 6, where (2)Trainerror=meanytrain-ypredtrain2.Table 3
Parameter values inGP and PFK-GP.
Parameter
Value
GPTIPS [29] multigene syntax
Representation
Number of genes: 1
Maximum tree depth: 5
Population size
50
Number of generations
1000
Selection
Lexicographic tournament selection
Tournament size
3
Crossover operator
Subtree crossover
Crosser probability
0.85
Mutation operator
Subtree mutation
Mutation probability
0.1
Reproduction probability
0.05
Fitness
1
N
∑
y
-
y
^
2
Elitism
Keep 1 best individualFigure 6
The evolutionary result ofP13 withPO-GP andPFK-GP.We could find that after processing test data of formula modelP13,PFK-GP found its best model at the first generation and its fitness is higher, whilePO-GP found its best model until 718th generation and its fitness is much lower than that inPFK-GP. ThePFK-GP can get the right formulas quickly because the model P13 recognized by the algorithmGenerateFs is inserted into the initialized population of evolutionary computation. For the formula models whose characteristics are consistent with data features inPFKB, they can be recognized with high probability and can be combined into population ofPFK-GP. ThePFK-GP can firstly search the coefficients in these formula models and get the mathematical expression with good fitness value. Therefore, the algorithmGenerateFs can speed up the process ofPFK-GP dealing with SR and can improve the accuracy of SR results.In order to test whetherPFK-GP can overcome overfitting or not, a dataset is created by E1 which has not existed in the training models ofPFKB. The two algorithmsPO-GP andPFK-GP are, respectively, applied to process the dataset. The two algorithms, which run, respectively, 100 and 1000 generations, have similar convergence curves in Figure 8. However,PFK-GP can find better fitness results compared withPO-GP, becausePFK-GP searches fitted solution in the space includes more functions whose data features are in accord with E1. Since the initial population, which is generated by the algorithms (CountSamePartF andCountSpecU) inPFKB, contains subformulas in formula models which are recognized byPFKB and represents data features of these subformulas,PFK-GP can find the right formulas which are more fitted to the raw dataset.In order to observe overall performance of thePFK-GP, we select six datasets as testing set. Three of them generated by formula models (P9,P13,P19) from Table 7 are involved in the process of trainingDBN, while the other three generated by formula models (E1,E4,E6) from Table 5 are not involved in that process. The two algorithmsPFK-GP andPO-GP are executed, respectively, ten times in order to gain the right formulas from the six different datasets. The six results of mean training error gained by the two algorithms are shown in Figure 9. And the average results from six groups of mean training errors are listed in Figure 7. ThePFK-GP(E) andGP(E) are the average results of E1,E4, and E6, whilePFK-GP(P) andPO-GP(P) are the average results of P9,P13, and P19. We can conclude that the comprehensive performance of thePFK-GP is better than that of thePO-GP based on the results in Figures 7 and 9, because the algorithmPFK-GP utilizes the methodGenerateFs to find the fitted formula model directly and the methodsCountSamePartF andCountSpecU to identify subformula models which have data features consistent with test set. The best mathematical expressionsPFK-GP andPO-GP found are listed in Table 4.Table 4
The best mathematical expression ofPFK-GP finding.
Number
The best mathematical expression
E
1
y
=
0.7001
∗
tan
⁡
x
1
∗
x
4
-
5.049
-
0.7001     ∗x1∗cube⁡2.575∗cube⁡x4-1.001
E
2
y
=
7.214
∗
sin
⁡
x
1
+
1.001
∗
tan
⁡
x
2
E
3
y
=
0.25
∗
square
⁡
x
1
+
x
2
-
6
-
0.2179
E
4
y
=
1.001
-
1
.
332
∗
4
∗
x
1
+
log
⁡
square
⁡
x
2
2
∗
square
⁡
x
2
+
5.585
E
5
y
=
x
2
-
2.092
tanh
⁡
square
⁡
sin
⁡
x
1
+
0.8795
E
6
y
=
6
∗
cos
⁡
x
2
∗
sin
⁡
x
1
-
0.00444
E
7
y
=
sin
⁡
x
1
-
6
∗
x
1
+
square
⁡
x
1
+
14
E
8
y
=
log
⁡
x
2
+
sqrt
⁡
x
1
+
sin
⁡
x
1
+
0.1823Table 5
Test data used inPFK-GP.
Number
Formula
E
0
y
=
-
1.97
+
1.25
∗
sqrt
⁡
x
0
tan
⁡
x
1
+
tan
⁡
x
0
exp
⁡
x
1
+
cos
⁡
x
0
∗sin⁡x1
E
1
y
=
exp
⁡
2
∗
x
1
∗
sin
⁡
p
i
∗
x
4
+
sin
⁡
x
2
∗
x
3
E
2
y
=
3.56
+
7.23
∗
sin
⁡
x
0
+
tan
⁡
x
3
E
3
y
=
x
0
-
3
∗
x
3
-
3
+
2
∗
sin
⁡
x
0
-
4
∗
x
3
-
4
E
4
y
=
quart
⁡
x
1
-
3
+
cube
⁡
x
2
-
3
-
x
2
-
3
quart
⁡
x
2
-
4
+
10
E
5
y
=
tanh
⁡
cos
⁡
2
∗
x
0
+
x
3
E
6
y
=
6
∗
sin
⁡
x
0
∗
cos
⁡
x
3
E
7
y
=
sin
⁡
x
0
+
square
⁡
x
0
+
5
E
8
y
=
sqrt
⁡
x
0
+
log
⁡
1.2
∗
x
3
+
sin
⁡
x
0Table 6
Test data of two algorithmsCountSamePartF and CountSpecU  (see Algorithms 3 and 4).
Number
Formula
T
1
P
1
T
2
P
2
T
3
P
3
T
4
P
4
T
5
P
5
T
6
y
=
z
+
a
∗
sqrt
⁡
x
0
tan
⁡
x
1
+
sin
⁡
x
1
T
7
y
=
z
+
a
∗
sqrt
⁡
x
0
tan
⁡
x
1
+
x
1
T
8
y
=
z
+
a
∗
tan
⁡
x
0
exp
⁡
x
1
+
cos
⁡
x
0
∗
sin
⁡
x
1
+
sin
⁡
x
1
T
9
y
=
z
+
a
∗
tan
⁡
x
0
exp
⁡
x
1
+
cos
⁡
x
0
∗
sin
⁡
x
1
+
x
1
T
10
y
=
z
+
a
∗
cos
⁡
x
0
∗
sin
⁡
x
1
+
square
⁡
x
1Table 7
Training data of algorithmTrainDBN  (see Algorithm 1).
Number
Formula
P
1
y
=
z
+
a
∗
x
0
P
2
y
=
z
+
a
1
∗
x
3
+
x
1
a
2
∗
x
4
P
3
y
=
z
+
a
1
∗
x
3
-
x
0
+
x
1
/
x
4
a
2
∗
x
4
P
4
y
=
z
+
a
∗
sin
⁡
x
P
5
y
=
z
+
a
∗
log
⁡
x
P
6
y
=
z
+
a
∗
sqrt
⁡
x
P
7
y
=
z
-
a
1
∗
exp
⁡
a
2
∗
x
0
P
8
y
=
z
+
a
1
∗
sqrt
⁡
a
2
∗
x
0
∗
x
3
∗
x
4
P
9
y
=
a
1
∗
sqrt
⁡
x
0
a
2
∗
log
⁡
x
1
∗
a
3
∗
exp
⁡
x
2
a
4
∗
square
⁡
x
3
P
10
y
=
z
+
a
1
∗
a
2
∗
x
1
+
a
3
∗
square
⁡
x
2
a
4
∗
cube
⁡
x
3
+
a
5
∗
quart
⁡
x
4
P
11
y
=
z
+
a
1
∗
cos
⁡
a
2
∗
x
0
∗
x
0
∗
x
0
P
12
y
=
z
-
a
1
∗
cos
⁡
a
2
∗
x
0
∗
sin
⁡
a
3
∗
x
4
P
13
y
=
z
-
a
∗
tan
⁡
x
0
tan
⁡
x
1
∗
tan
⁡
x
2
tan
⁡
x
3
P
14
y
=
z
-
a
∗
cos
⁡
x
0
-
tan
⁡
x
1
∗
tanh
⁡
x
2
sin
⁡
x
3
P
15
y
=
z
-
a
∗
tan
⁡
x
0
exp
⁡
x
1
∗
log
⁡
x
2
-
tan
⁡
x
3
P
16
y
=
a
∗
x
3
P
17
y
=
a
1
∗
x
1
+
a
2
∗
x
4
P
18
y
=
sqrt
⁡
x
2
tan
⁡
x
5
/
a
P
19
y
=
cos
⁡
x
2
cube
⁡
x
5
/
a
P
20
y
=
tanh
⁡
x
2
∗
a
∗
cube
⁡
x
5
+
abs
⁡
x
1
P
21
y
=
tanh
⁡
abs
⁡
x
2
∗
a
+
x
5
∗
cube
⁡
x
5
+
abs
⁡
x
1
P
22
y
=
tanh
⁡
tan
⁡
x
5
a
∗
cube
⁡
x
5
+
abs
⁡
x
1
P
23
y
=
tanh
⁡
cos
⁡
x
2
∗
a
∗
cube
⁡
sqrt
⁡
x
2
P
24
y
=
tanh
⁡
cos
⁡
x
2
∗
a
∗
cube
⁡
x
5
+
abe
⁡
x
1
P
25
y
=
z
P
26
y
=
z
+
x
2
P
27
y
=
z
+
x
2
x
0
∗
x
2
z
+
x
2
x
0
∗
x
2
X
X
X
P
28
y
=
x
0
-
z
1
x
0
+
x
2
∗
x
5
-
z
2
x
0
∗
a
1
P
29
y
=
a
∗
sqrt
⁡
x
P
30
y
=
a
∗
log
⁡
x
P
31
y
=
a
∗
square
⁡
x
P
32
y
=
a
∗
tanh
⁡
x
P
33
y
=
a
∗
sin
⁡
x
P
34
y
=
a
∗
cos
⁡
x
P
35
y
=
a
∗
exp
⁡
x
P
36
y
=
a
∗
cube
⁡
x
P
37
y
=
a
∗
quart
⁡
x
P
38
y
=
a
∗
tan
⁡
x
P
39
y
=
a
∗
abs
⁡
xFigure 7
Average results from six groups of means training errors inPO-GP andPFK-GP.Figure 8
The SR evolutionary process ofE1 withPO-GP orPFK-GP under different generations.Figure 9
Training error results in which six datasets generated byE1, E4, E6, P1, P13, and P19 dealt withPO-GP andPFK-GP, respectively.In order to measure relativity between experimental data and predictive data, the formula Training Variation Explained (TVE) is defined as follows:(3)TVE=1-sumytrain-ypredtrain2sumytrain-meanytrain2.The higher the TVE value, the more valid the predictive data.PO-GP andPFK-GP are run ten times, respectively, in the dataset generated from eight prediction models (see Table 6  E1–E8). The eight results of different dataset processed by the above two algorithms are listed in Figure 10. And the maximum, minimum, and average results of TVE are listed in Figure 11. From the results in the two figures, the formulas thatPFK-GP finds are more relative to the experimental formula models than thosePO-GP finds.Figure 10
TVE for modelsE1–E8 contrast the random population in traditionalPO-GP with thePFK-GP.Figure 11
TVE results ofPO-GP compared withPFK-GP in eight formula models.


## 5. Related Work

The search space of SR is huge even for rather simple basis functions [31]. In order to avoid search space that is too far from the desired output range determined by the training dataset, the interval arithmetic [34] and the affine arithmetic [26], which can compute the bounds ofGP tree expression, are imported into SR. Although the method based on affine arithmetic can generate the tighter bounds of the expression in comparison with the interval arithmetic method, its accuracy often leads to high computational complexity [35]. Moreover, the size of search space is still huge because there are plentiful candidate expressions which fit to the data bound computed by the above two arithmetic methods.In addition to the above arithmetic method, machine learning methods are used to compact or reduce the search space of SR. FFX technology uses pathwise regularized learning algorithm to rapidly prune a huge set of candidate basis functions down to compact model based on the generalization linearly model (GLM); hence the technology outperforms GP-SR in speed and scalability due to its simplicity and deterministic nature [8]. However, it may abandon correct expressions and make them not in the space of GLM. A hybrid deterministic GP-SR algorithm [36] is proposed to overcome the problem of missing correct expression. The hybrid algorithm extracts candidate basis expressions by using FFX and inputs the expressions into the GP-SR. The hybrid algorithm utilizes the candidate expression generated by the linear regression method (pathwise regulation), while our algorithm utilizes the candidate expression by applying the algorithmsCountSamePartF, GenerateFs, andCountSpecU.By applying expectation-maximization (EM) framework to SR, the clustered SR (CSR) can identify and infer symbolic repression of piecewise function from unlabelled, time-series data [9]. The CSR can reduce the space of searching piecewise function owing to the fact that the EM can search simultaneously the subfunction parameters and latent variables that represent the information of function segment. The abstract expression grammar (AEG) SR is proposed to perform the process of genetic algorithm (GA), allowing user control of the search space and the final output formulas [37]. On understanding the given application, users can specify the goal expression of SR and limit the size of search space by using abstract expression grammars. Compared with manually assigning expression and limiting the search space with AEGSR, in the paper, the methods about PFK can automatically extract the candidate expression from dataset by using statistical method and dynamically adjust the search space by usingGP.The methods that inject prior or expert knowledge in evolutionary search [12, 13] are introduced to find effective solutions that can show mathematical expression more compactable and interpretable. In these papers, the prior and expert knowledge are the solutions which are mathematical expressions in some applications. The knowledge is merged intoGP by inserting randomized pieces of the approximate solution into population. One of the major differences between these methods and our method is how prior or expert knowledge is created. The knowledge in [12, 13] is the existing formula model that comes from the previous solutions and can be called static knowledge. However, the knowledge in our method is the formula model which is consistent with data features that are originated from the algorithmsGenerateFs, CountSamePartF, andCountSpecU and can be called dynamical knowledge that is changed with the features of test dataset. Therefore, our methods can insert more suitable knowledge into theGP.


## 6. Conclusion

In this paper, aPFK-GP method is proposed to deal with the problem of symbolic regression based on analyzing the process of how a researcher constructs a mathematical model. The method can understand experimental data features and can extract some formulas consistent with experimental data features. In order to implement the function of understanding data features,PFK-GP, through theDBN method, firstly createsPFKB that can extract features from test dataset generated by training formula models. The experiment results confirm, compared with SVM, thatDBN can produce better results that extract features from formula models and classify test data into its corresponding formula model. Then, the methods of classification and recognition are implemented to find some formula models that are similar or related to experimental data features as much as possible. For the classification, we exploit the algorithmGenerateFs based onDBN to match the experimental data with formula models inPFKB. With regard to recognition, we propose the algorithms ofCountSamePartF andCountSpecU to obtain some subformula models which have local features consistent with experimental data. The classification can helpPFKB to find formula models that are consistent with whole data features while the recognition can helpPFKB to find subformula models consistent with local data features. At last, the algorithmrandomGenP is used to generate individuals of evolutionary population according to the result of the above three algorithms. Through combining and transforming these individuals,GP can automatically obtain approximate formulas that are best fitting to the experimental data.Compared with ParetoGP,PFK-GP, which is built on thePFKB with the functions of classification and recognition, can explore formulas in the search space of data features. So, it can accelerate the speed of convergence and improve the accuracy of formula obtained.Obviously, the high efficiency ofPFK-GP depends on the powerful methods of classification and recognition based onPFKB. Therefore, it is an important part of the future work to improve the accuracy of the above two methods. The two methods depend on the representation of data features of formula model. In the paper, the two assumptions based on statistics and counts are used to obtain the formulas which can show the data features. The features of formula model are not defined explicitly. And the two assumption are not proved by formal proofs. There are some uncertainties in those assumptions. Therefore, the new representation which can show whole or local features of formula models will be researched to find formulas which can better fit to experiment data. In addition, the rules of formulas transforming and inferring that are similar to researchers’ methods will be explored in the evolution ofGP.



---

*Source: 1021378-2015-12-24.xml*