# A Novel Plant Root Foraging Algorithm for Image Segmentation Problems

Authors: Lianbo Ma; Kunyuan Hu; Yunlong Zhu; Hanning Chen; Maowei He
Journal: Mathematical Problems in Engineering
 (2014)

---

## Abstract

This paper presents a new type of biologically-inspired global optimization methodology for image segmentation based on plant root foraging behavior, namely, artificial root foraging algorithm (ARFO). The essential motive of ARFO is to imitate the significant characteristics of plant root foraging behavior including branching, regrowing, and tropisms for constructing a heuristic algorithm for multidimensional and multimodal problems. A mathematical model is firstly designed to abstract various plant root foraging patterns. Then, the basic process of ARFO algorithm derived in the model is described in details. When tested against ten benchmark functions, ARFO shows the superiority to other state-of-the-art algorithms on several benchmark functions. Further, we employed the ARFO algorithm to deal with multilevel threshold image segmentation problem. Experimental results of the new algorithm on a variety of images demonstrated the suitability of the proposed method for solving such problem.


---

## Body

## 1. Introduction

In recent years, biologically-inspired optimization algorithms have been widely applied for solving complex optimization problems. These algorithms have achieved great development and have been proven to be more powerful than the conventional methods based on formal logics or mathematical programming [1–4]. For instances, genetic algorithm (GA), originally proposed by Holland [5], represents a fairly abstract model of Darwinian evolution and biological genetics; ant colony optimization (ACO), conceived by Dorigo et al. [6, 7], is inspired by the ants’ social behaviors in finding shortest paths; particle swarm optimization (PSO), proposed by Eberhart and Kennedy [8], gleans ideas from social behavior of bird flocking and fish schooling. However, the existing bionic optimization models usually simulate animal behaviors, and the research of some natural plant-related evolution mechanisms does not obtain enough attention, such as plant root foraging behavior [9]. Currently, capturing the novel properties of potential plants to construct new intelligent computation models has attracted more and more attention.In plant ecosystem, plant root foraging is viewed as one significant characteristic of plant behavior, which utilizes the morphological and physiological plastic responses affecting individual plant growth [10]. Plastic root foraging responses have been widely recognized as an important strategy for plants to explore heterogeneously distributed resources [11–14]. Root morphological plasticity generates higher elongate-length and branching biomass per unit soil volume in nutrient-rich patches compared to nutrient-poor patches. If a randomly growing root encounters an area high in moisture, it is likely to branch profusely and grow, which is called root tropism. After many generations, the spatial configuration of root system in the soil is shaped with higher root biomass in the area higher in moisture. Logically, such evolutionary principles can be modeled as an optimization process. Here, a novel optimization model based on plant root foraging behaviors is established, called ARFO. This method provides an open framework to utilize research in plant behavioral ecology to tackle complex problems, which set up a new information section mechanism based on auxin to determine how to select new growing points and branching number of roots.Image segmentation is regarded as a basic operation for meaningful analysis and interpretation of image acquired, responsible for discriminating an object from other objects having distinct gray levels or different colors. As an important and effective tool among the existing segmentation techniques, multilevel threshold segmentation technique can be roughly classified into two types: optimal threshold methods [15–18] and property based threshold methods [19–21]. Those multilevel threshold methods are fast and suitable for the case that the number of thresholds is hard to determine and should be specified in advance. Several algorithms have so far been proposed in literatures that have addressed the issue of optimal threshold [22, 23]. Many of these algorithms attempt to achieve optimal threshold such that the thresholded classes achieve some desired characteristic. References [18, 22, 23] proposed some novel methods derived from optimizing an objective function. However, all of these methods have a common problem, that the computational complexity rises exponentially when extended to multilevel threshold due to the exhaustive search employed, which limits the multilevel threshold applications.Recently, to address this concerned issue, several swarm intelligence (SI) algorithms as the powerful optimization tools have been introduced to the field of image segmentation owing to their predominant abilities of coping with complex nonlinear optimizations [24, 25]. In this paper, the ARFO algorithm is employed for multilevel threshold image segmentation to prove its effectiveness and robustness, which is based on entropy criterion based fitness measure. Several other mentioned algorithms are tested as a comparison.The rest of the paper is organized as follows. In Section2, the root foraging model and the implementation details of ARFO algorithm will be given. Section 3 experimentally compares ARFO with other well-known algorithms on ten benchmark functions. And its application to image segmentation has been presented in Section 4. Finally, Section 5 outlines the conclusion.


## 2. Root Foraging Model for Optimization

### 2.1. Root Foraging Model

The root foraging behavior concept is proposed by McNickle et al. based on MVT theory [9, 26], in which plant roots exhibit morphological responses to environmental heterogeneity through the selective allocation of root biomass per unit time to the nitrogen-rich patch. It is worthy noted that there already exist some effective approaches to mimic such root foraging behaviors, such asL-system [27–29]. But the key issue for optimization application is not resolved, which is that how to select new growing points to grow and how to ensure the roots towards the optimal position.In the plant biology, the tropism behaviors of roots system, especially the gravitropism and hydrotropism, are the significantly important characteristics of root behaviors [30]. Gravitropism induces plant roots grow along a specific vector relative to gravity [31] while hydrotropism is the root morphological responses is guided by a stimuli or gradient in nutrition concentration [32]. Note that a kind of essential hormone information-auxin plays an important role in conducting these tropism behaviors, influencing simultaneously the spatial configuration of the roots system [33]. In order to idealize above plant root foraging behaviors, some criteria are presented as follows.(
1) Criterion-1: Auxin Regulation. Auxin conducts the root’s adaptive growth and the information exchange among root tips. The roots spatial structure is regulated by the auxin concentration, which is not static, but dynamically reallocated after new roots germinate and grow.(
2) Criterion-2: Root Growth Behavior. Regrowing means roots can elongate with a certain angle, and branching is the main progress for multiplying the number of roots, in which new roots originate from the branching of the parent root tip by a certain elongate-unit in new direction.(
3) Criterion-3: Classification of the Root System. The whole plant root system is divided into three categories sorted by the auxin distribution: mainroots, the lateral-roots, and the dead-roots. Each root tip of mainroots can implement both regrowing and branching, in which the branching number is determined by auxin concentration. The root tips of lateral-roots can only implement regrowing operation.(
4) Criterion-4: Root Tropisms. Root affected by various tropisms evolves towards the most nutrient rich position. Hydrotropism makes the growing direction of the root tips towards the optimal individual position. Meantime gravitropism imposes that roots grow in the direction of gravitational pull (i.e., downward), which means that the growing angle of each root tip is limited to a certain range.The root foraging model assumes that the optimization problem, whose feasible region is used as soil environment, is a minimization problem withD-dimensional vectors (d
1
,
d
2
,
…
,
d
D) and a growing point (i.e., root tip) corresponds to a potential solution instance to the problem. Through the regrowing and branching operations, the root tips of the proposed model will quickly spread to the global optimum in the feasible region of optimization problems, until no new branches are generated. Note that, for more precise modeling of the auxin information mechanism, the nutrient concentration strategy is incorporated.


### 2.2. Auxin Information

The auxin information is conducted to regulate the growth points and branching number of roots. We assumed that the sum of auxin concentrationF
i is considered as 1 in the mathematical model of root systems. In order to calculate the F
i value of every root the next equation is used:
(1)
f
i
=
fi
t
i
-
f
worst
f
best
-
f
worst
,
(2)
N
i
=
Nutitio
n
i
-
Nutitio
n
worst
Nutitio
n
best
-
Nutitio
n
worst
,
(3)
F
i
=
f
i
∑
j
=
1
N
f
i
*
ξ
+
N
i
∑
j
=
1
N
N
i
(
1
-
ξ
)
,
ξ
⊂
(
0,1
)
,
where i is the position of the growing point, ξ is a uniform random quantity, N is the total number of the points, fit (·) is the fitness value of the point. and f
worst and f
best are the maximum and minimum of the current points, respectively. Nutritio
n
i is the current nutrient concentration of individual i.In addition to the relative position of each point (objective function value), the auxin concentration is also determined by the gradient of function fitness. We establish the concept of nutrient concentration to simulate the moisture gradient of soil environment, expressed as(4)
Nutitio
n
i
(
t
+
1
)
=
{
Nutitio
n
i
(
t
)
+
1
,
if
f
i
(
t
+
1
)
<
f
i
(
t
)
,
Nutitio
n
i
(
t
)
-
1
,
if
f
i
(
t
+
1
)
>
f
i
(
t
)
.Nutrient concentration reveals that even if the amount of moisture is not best, the higher moisture gradient of the soil is also conducive to the root cell division to absorb water and nutrition.


### 2.3. Mainroots Growth Operations

According toCriterion-3 of Section 2, at each iteration, the root population consists of three kinds of members: mainroots, lateral-roots, and dead-roots. The root located in the most promising area, conferring the highest auxin concentration value, is chosen as mainroot. It then implements the regrowing operator and branching operator to explore and exploit the whole space. The procedures of simulating foraging behavior of mainroots are abstracted as follows.


#### 2.3.1. Regrowing Operator

Since the mainroots with higher auxin concentration values tend to mature, their growth uses larger elongate-length unit, taking on significant hydrotropism and gravitropism characteristics. In order to simulate the influences of these tropisms in the mainroots growth process, the mathematic expression of regrowing operator can be described as follows.Step 1.
In each cycle of roots growth process, all root tips are sorted by auxin concentration values. Some excellent root tips with higher auxin concentration values have higher probability to be selected as mainroots for branching operator. The number of the roots chosen as mainroots is defined by(5)
S
=
N
*
C
r
,
where S is the size of mainroots group, N is the total number of the roots, and C
r is the selection probability.Step 2.
Under hydrotropism, select some of current mainroots to search towards the optimal position of individuals:(6)
X
i
(
t
+
1
)
=
X
i
(
t
)
+
R
3
(
X
best
-
X
i
(
t
)
)
.Step 3.
Under gravitropism, new growing points of rest mainroots can be formulated as(7)
X
i
(
t
+
1
)
=
X
i
(
t
)
+
R
4
*
Gro
w
max
⁡
*
D
i
(
φ
)
.
And the growth angle representation of the mainrooti is defined by
(8)
φ
i
(
t
+
1
)
=
φ
i
(
t
)
+
R
5
*
α
max
⁡
+
R
6
*
π
,
α
max
⁡
<
π
,
where R
3, R
5, and R
6 are random coefficients varying from 0 to 1, R
4 is a normally distributed random number with mean 0 and standard deviation 1, and Gro
w
max
⁡ is the maximum of elongate-length (i.e., objective function boundary range). To mimic the phenomenon that mainroot tip in the next generation tends to grow downwards under gravitropism, the entire growing angle space α
max
⁡ is limited to π.From (6), the mainroot can rapidly approach a near-optimal solution and stay for some growth cycles in several domains containing local optima. From (7), it can escape from these local optima to explore other domains perhaps containing the global optimum, but it was not able to stop there.


#### 2.3.2. Branching Operator

In the root foraging model, the branching operator is the local search strategy to exploit around the optimal area, in which the mainroots branch new roots using a small elongate-length unit, forming a certain density distribution slowly in the promising region’s immediate vicinity. This procedure can be described as follows.Step 1.
The nutrient concentration (Nutrition
i) of mainroot i is compared with the threshold value BranchG (satisfied 0
<
Branch
G
<
1) to determine whether it performs branching operator by
(9)
branching
,
if
Nutitio
n
i
>
Branch
G
,
nobranching
,
otherelse
.Step 2.
Determine the branching number according to the auxin concentration distribution. In principle, the position with higher moisture content and gradient in soil is more conducive to roots branching. The branching numberw
i is calculated as
(10)
w
i
=
{
S
min
⁡
F
i
*
r
*
(
S
max
⁡
-
S
min
⁡
)
+
S
min
⁡
,
if
Nutitio
n
i
>
Branch
G
,
0
,
else
,
where S
max
⁡ and S
min
⁡ are the maximum and minimum of the new growing points, respectively, and r is random distribution coefficient. F
i is auxin concentrationof mainroot i.Step 3.
Calculate the fitness of new branching points. Regarding the growth direction of the initial mainroot as zero degree, the searching angle space is divided intoS
max
⁡ subzones and the angle of new growing point is randomly falling within one of S
max
⁡ angle subzones. The mathematic equation of new growing point yields
(11)
X
i
(
t
+
1
)
=
X
i
(
t
)
+
R
1
*
D
i
(
φ
i
)
,
where X
i
(
t
+
1
) is the new growing point from X
i
(
t
), R
1 is the elongate-length unit, which is a random varying from 0 to 1, φ
i is the growth angle (
φ
i
1
,
φ
i
2
,
…
,
φ
i
(
D
-
1
)
), and the growth direction D
i
(
φ
i
)
=
(
d
i
1
,
d
i
2
,
…
,
d
i
D
) which can be calculated from φ
i via a Polar to Cartesian coordinates transformation [34]:
(12)
d
i
1
=
∏
k
=
1
D
-
1
cos
⁡
(
φ
i
k
)
,
d
i
j
=
sin
(
φ
i
j
)
∏
k
=
1
D
-
1
cos
⁡
⁡
(
φ
i
k
)
,
d
i
D
=
sin
(
φ
i
D
-
1
)
.
The growth angleφ
i is calculated as follows:
(13)
φ
(
t
+
1
)
=
φ
(
t
)
+
R
2
α
init
+
K
*
β
max
⁡
S
max
⁡
,
where R
2 is a random coefficient varying from 0 to 1, α
init is original growth angle value of the initial mainroot as zero degree, K is randomly parameter selecting the subzone, S
max
⁡ is subzones number, and β
max
⁡ is the maximum growing turning angle. Similar to (8), considering gravitropism, the growing angle space β
max
⁡ is limited to π to make mainroots tend to grow downwards.


### 2.4. Lateral-Roots Growth Operation

Due to the fact that the lateral-roots group with lower auxin concentration stays in the initial growth phase, its elongation step is relatively small and the influence of the hydrotropism can be ignored. The lateral-roots with the smaller elongation step are exploiting a local optimal domain, which closed to its original point.The mathematic expression, which simulates above behavior, can be defined as (14):
(14)
X
i
(
t
+
1
)
=
X
i
(
t
)
+
R
5
*
Gro
w
max
⁡
*
D
i
(
φ
i
)
,
(15)
φ
i
(
t
+
1
)
=
φ
i
(
t
)
+
R
6
α
max
⁡
,
where R
5 is normal distribution function (mean is 0, deviation is 1) and R
6 is a random coefficient.


### 2.5. Dead-Roots Elimination

In the proposed root foraging model, it is assumed thatN
i is the current population size, N
i will increase by one if a root tip splits and reduce by one if a root dies determined by auxin distribution, and it will vary in the searching process [35]. The branching criterion and dead-roots eliminating criterion are listed in
(16)
N
i
=
N
i
+
w
i
,
if
F
i
>
Branch
G
,
N
i
=
N
i
-
1
,
if
F
i
<
N
mority
,
where w
i is the branching number defined by (10), Branch
G is the branching threshold value defined by (9), and N
mority is the death threshold value parameter.


### 2.6. ARFO Algorithm

We implement an instantiation algorithm of the proposed root foraging model, namely, artificial root foraging optimizer (ARFO), for numerical function optimization. The pseudocode of the detailed procedure is given as follows.Step 1 (initialization).
Initialize the positions of root tips population and evaluate the nutrient (fitness) of the population. Set maximum number of cycles (LimitC).Step 2.
Iteration = 0.Step 3 (auxin concentration calculation).
(
1
) Calculate the nutrient concentration and auxin concentration values of the population by (2) and (3).
(
2
) Divide the population into two groups: mainroots group and lateral-roots group by (5).Step 4 (mainroots operation).
Loop over each mainroot tip(1)
Select a random portion of mainroots group to take regrowing operator by (6), (7), and (8).
(2)
Evaluate the fitness of the renewal mainroots and apply greedy selection.
(3)
Adjust the corresponding nutrient concentration value by (4).
(4)
If the condition of branching determined by (9) is met, continue; otherwise, go to Step 5.
(5)
Calculate the branching number by (10) and the branching new roots by (11).
(6)
The population size increases by the corresponding number calculated by (10).Step 5 (lateral-roots operation).
Loop over each root tip of lateral-roots(1)
Lateral-root takes regrowing operator by (14) and (15).
(2)
Evaluate the fitness of the renewal lateral-root and apply greedy selection.
(3)
Adjust the corresponding nutrient concentration value by (4).Step 6 (dead roots eliminating).
Remove the dead individuals from the population according to their auxin concentration values by (16).Step 7.
Iteration = Iteration + 1;Step 8.
If the iteration is greater than LimitC, stop the procedure; otherwise, go to Step3.Step 9.
Output the best solution achieved.By using the regrowing and branching operators, ARFO not only alleviates the premature convergence in the entire solution space, but also improves the accuracy of optimum in the local solution space, at least to some extent.(1)
Exploring strategy. The regrowing operator implemented by mainroots employs a large elongate-length unit to explore the previously unscanned regions in the search space.
(2)
Exploiting strategy. To perform fine-tuning exploitation of the global optimum, both the branching operator and the lateral-root’s regrowing operator employ a small elongate-length unit and a random angle.


## 2.1. Root Foraging Model

The root foraging behavior concept is proposed by McNickle et al. based on MVT theory [9, 26], in which plant roots exhibit morphological responses to environmental heterogeneity through the selective allocation of root biomass per unit time to the nitrogen-rich patch. It is worthy noted that there already exist some effective approaches to mimic such root foraging behaviors, such asL-system [27–29]. But the key issue for optimization application is not resolved, which is that how to select new growing points to grow and how to ensure the roots towards the optimal position.In the plant biology, the tropism behaviors of roots system, especially the gravitropism and hydrotropism, are the significantly important characteristics of root behaviors [30]. Gravitropism induces plant roots grow along a specific vector relative to gravity [31] while hydrotropism is the root morphological responses is guided by a stimuli or gradient in nutrition concentration [32]. Note that a kind of essential hormone information-auxin plays an important role in conducting these tropism behaviors, influencing simultaneously the spatial configuration of the roots system [33]. In order to idealize above plant root foraging behaviors, some criteria are presented as follows.(
1) Criterion-1: Auxin Regulation. Auxin conducts the root’s adaptive growth and the information exchange among root tips. The roots spatial structure is regulated by the auxin concentration, which is not static, but dynamically reallocated after new roots germinate and grow.(
2) Criterion-2: Root Growth Behavior. Regrowing means roots can elongate with a certain angle, and branching is the main progress for multiplying the number of roots, in which new roots originate from the branching of the parent root tip by a certain elongate-unit in new direction.(
3) Criterion-3: Classification of the Root System. The whole plant root system is divided into three categories sorted by the auxin distribution: mainroots, the lateral-roots, and the dead-roots. Each root tip of mainroots can implement both regrowing and branching, in which the branching number is determined by auxin concentration. The root tips of lateral-roots can only implement regrowing operation.(
4) Criterion-4: Root Tropisms. Root affected by various tropisms evolves towards the most nutrient rich position. Hydrotropism makes the growing direction of the root tips towards the optimal individual position. Meantime gravitropism imposes that roots grow in the direction of gravitational pull (i.e., downward), which means that the growing angle of each root tip is limited to a certain range.The root foraging model assumes that the optimization problem, whose feasible region is used as soil environment, is a minimization problem withD-dimensional vectors (d
1
,
d
2
,
…
,
d
D) and a growing point (i.e., root tip) corresponds to a potential solution instance to the problem. Through the regrowing and branching operations, the root tips of the proposed model will quickly spread to the global optimum in the feasible region of optimization problems, until no new branches are generated. Note that, for more precise modeling of the auxin information mechanism, the nutrient concentration strategy is incorporated.


## 2.2. Auxin Information

The auxin information is conducted to regulate the growth points and branching number of roots. We assumed that the sum of auxin concentrationF
i is considered as 1 in the mathematical model of root systems. In order to calculate the F
i value of every root the next equation is used:
(1)
f
i
=
fi
t
i
-
f
worst
f
best
-
f
worst
,
(2)
N
i
=
Nutitio
n
i
-
Nutitio
n
worst
Nutitio
n
best
-
Nutitio
n
worst
,
(3)
F
i
=
f
i
∑
j
=
1
N
f
i
*
ξ
+
N
i
∑
j
=
1
N
N
i
(
1
-
ξ
)
,
ξ
⊂
(
0,1
)
,
where i is the position of the growing point, ξ is a uniform random quantity, N is the total number of the points, fit (·) is the fitness value of the point. and f
worst and f
best are the maximum and minimum of the current points, respectively. Nutritio
n
i is the current nutrient concentration of individual i.In addition to the relative position of each point (objective function value), the auxin concentration is also determined by the gradient of function fitness. We establish the concept of nutrient concentration to simulate the moisture gradient of soil environment, expressed as(4)
Nutitio
n
i
(
t
+
1
)
=
{
Nutitio
n
i
(
t
)
+
1
,
if
f
i
(
t
+
1
)
<
f
i
(
t
)
,
Nutitio
n
i
(
t
)
-
1
,
if
f
i
(
t
+
1
)
>
f
i
(
t
)
.Nutrient concentration reveals that even if the amount of moisture is not best, the higher moisture gradient of the soil is also conducive to the root cell division to absorb water and nutrition.


## 2.3. Mainroots Growth Operations

According toCriterion-3 of Section 2, at each iteration, the root population consists of three kinds of members: mainroots, lateral-roots, and dead-roots. The root located in the most promising area, conferring the highest auxin concentration value, is chosen as mainroot. It then implements the regrowing operator and branching operator to explore and exploit the whole space. The procedures of simulating foraging behavior of mainroots are abstracted as follows.


### 2.3.1. Regrowing Operator

Since the mainroots with higher auxin concentration values tend to mature, their growth uses larger elongate-length unit, taking on significant hydrotropism and gravitropism characteristics. In order to simulate the influences of these tropisms in the mainroots growth process, the mathematic expression of regrowing operator can be described as follows.Step 1.
In each cycle of roots growth process, all root tips are sorted by auxin concentration values. Some excellent root tips with higher auxin concentration values have higher probability to be selected as mainroots for branching operator. The number of the roots chosen as mainroots is defined by(5)
S
=
N
*
C
r
,
where S is the size of mainroots group, N is the total number of the roots, and C
r is the selection probability.Step 2.
Under hydrotropism, select some of current mainroots to search towards the optimal position of individuals:(6)
X
i
(
t
+
1
)
=
X
i
(
t
)
+
R
3
(
X
best
-
X
i
(
t
)
)
.Step 3.
Under gravitropism, new growing points of rest mainroots can be formulated as(7)
X
i
(
t
+
1
)
=
X
i
(
t
)
+
R
4
*
Gro
w
max
⁡
*
D
i
(
φ
)
.
And the growth angle representation of the mainrooti is defined by
(8)
φ
i
(
t
+
1
)
=
φ
i
(
t
)
+
R
5
*
α
max
⁡
+
R
6
*
π
,
α
max
⁡
<
π
,
where R
3, R
5, and R
6 are random coefficients varying from 0 to 1, R
4 is a normally distributed random number with mean 0 and standard deviation 1, and Gro
w
max
⁡ is the maximum of elongate-length (i.e., objective function boundary range). To mimic the phenomenon that mainroot tip in the next generation tends to grow downwards under gravitropism, the entire growing angle space α
max
⁡ is limited to π.From (6), the mainroot can rapidly approach a near-optimal solution and stay for some growth cycles in several domains containing local optima. From (7), it can escape from these local optima to explore other domains perhaps containing the global optimum, but it was not able to stop there.


### 2.3.2. Branching Operator

In the root foraging model, the branching operator is the local search strategy to exploit around the optimal area, in which the mainroots branch new roots using a small elongate-length unit, forming a certain density distribution slowly in the promising region’s immediate vicinity. This procedure can be described as follows.Step 1.
The nutrient concentration (Nutrition
i) of mainroot i is compared with the threshold value BranchG (satisfied 0
<
Branch
G
<
1) to determine whether it performs branching operator by
(9)
branching
,
if
Nutitio
n
i
>
Branch
G
,
nobranching
,
otherelse
.Step 2.
Determine the branching number according to the auxin concentration distribution. In principle, the position with higher moisture content and gradient in soil is more conducive to roots branching. The branching numberw
i is calculated as
(10)
w
i
=
{
S
min
⁡
F
i
*
r
*
(
S
max
⁡
-
S
min
⁡
)
+
S
min
⁡
,
if
Nutitio
n
i
>
Branch
G
,
0
,
else
,
where S
max
⁡ and S
min
⁡ are the maximum and minimum of the new growing points, respectively, and r is random distribution coefficient. F
i is auxin concentrationof mainroot i.Step 3.
Calculate the fitness of new branching points. Regarding the growth direction of the initial mainroot as zero degree, the searching angle space is divided intoS
max
⁡ subzones and the angle of new growing point is randomly falling within one of S
max
⁡ angle subzones. The mathematic equation of new growing point yields
(11)
X
i
(
t
+
1
)
=
X
i
(
t
)
+
R
1
*
D
i
(
φ
i
)
,
where X
i
(
t
+
1
) is the new growing point from X
i
(
t
), R
1 is the elongate-length unit, which is a random varying from 0 to 1, φ
i is the growth angle (
φ
i
1
,
φ
i
2
,
…
,
φ
i
(
D
-
1
)
), and the growth direction D
i
(
φ
i
)
=
(
d
i
1
,
d
i
2
,
…
,
d
i
D
) which can be calculated from φ
i via a Polar to Cartesian coordinates transformation [34]:
(12)
d
i
1
=
∏
k
=
1
D
-
1
cos
⁡
(
φ
i
k
)
,
d
i
j
=
sin
(
φ
i
j
)
∏
k
=
1
D
-
1
cos
⁡
⁡
(
φ
i
k
)
,
d
i
D
=
sin
(
φ
i
D
-
1
)
.
The growth angleφ
i is calculated as follows:
(13)
φ
(
t
+
1
)
=
φ
(
t
)
+
R
2
α
init
+
K
*
β
max
⁡
S
max
⁡
,
where R
2 is a random coefficient varying from 0 to 1, α
init is original growth angle value of the initial mainroot as zero degree, K is randomly parameter selecting the subzone, S
max
⁡ is subzones number, and β
max
⁡ is the maximum growing turning angle. Similar to (8), considering gravitropism, the growing angle space β
max
⁡ is limited to π to make mainroots tend to grow downwards.


## 2.3.1. Regrowing Operator

Since the mainroots with higher auxin concentration values tend to mature, their growth uses larger elongate-length unit, taking on significant hydrotropism and gravitropism characteristics. In order to simulate the influences of these tropisms in the mainroots growth process, the mathematic expression of regrowing operator can be described as follows.Step 1.
In each cycle of roots growth process, all root tips are sorted by auxin concentration values. Some excellent root tips with higher auxin concentration values have higher probability to be selected as mainroots for branching operator. The number of the roots chosen as mainroots is defined by(5)
S
=
N
*
C
r
,
where S is the size of mainroots group, N is the total number of the roots, and C
r is the selection probability.Step 2.
Under hydrotropism, select some of current mainroots to search towards the optimal position of individuals:(6)
X
i
(
t
+
1
)
=
X
i
(
t
)
+
R
3
(
X
best
-
X
i
(
t
)
)
.Step 3.
Under gravitropism, new growing points of rest mainroots can be formulated as(7)
X
i
(
t
+
1
)
=
X
i
(
t
)
+
R
4
*
Gro
w
max
⁡
*
D
i
(
φ
)
.
And the growth angle representation of the mainrooti is defined by
(8)
φ
i
(
t
+
1
)
=
φ
i
(
t
)
+
R
5
*
α
max
⁡
+
R
6
*
π
,
α
max
⁡
<
π
,
where R
3, R
5, and R
6 are random coefficients varying from 0 to 1, R
4 is a normally distributed random number with mean 0 and standard deviation 1, and Gro
w
max
⁡ is the maximum of elongate-length (i.e., objective function boundary range). To mimic the phenomenon that mainroot tip in the next generation tends to grow downwards under gravitropism, the entire growing angle space α
max
⁡ is limited to π.From (6), the mainroot can rapidly approach a near-optimal solution and stay for some growth cycles in several domains containing local optima. From (7), it can escape from these local optima to explore other domains perhaps containing the global optimum, but it was not able to stop there.


## 2.3.2. Branching Operator

In the root foraging model, the branching operator is the local search strategy to exploit around the optimal area, in which the mainroots branch new roots using a small elongate-length unit, forming a certain density distribution slowly in the promising region’s immediate vicinity. This procedure can be described as follows.Step 1.
The nutrient concentration (Nutrition
i) of mainroot i is compared with the threshold value BranchG (satisfied 0
<
Branch
G
<
1) to determine whether it performs branching operator by
(9)
branching
,
if
Nutitio
n
i
>
Branch
G
,
nobranching
,
otherelse
.Step 2.
Determine the branching number according to the auxin concentration distribution. In principle, the position with higher moisture content and gradient in soil is more conducive to roots branching. The branching numberw
i is calculated as
(10)
w
i
=
{
S
min
⁡
F
i
*
r
*
(
S
max
⁡
-
S
min
⁡
)
+
S
min
⁡
,
if
Nutitio
n
i
>
Branch
G
,
0
,
else
,
where S
max
⁡ and S
min
⁡ are the maximum and minimum of the new growing points, respectively, and r is random distribution coefficient. F
i is auxin concentrationof mainroot i.Step 3.
Calculate the fitness of new branching points. Regarding the growth direction of the initial mainroot as zero degree, the searching angle space is divided intoS
max
⁡ subzones and the angle of new growing point is randomly falling within one of S
max
⁡ angle subzones. The mathematic equation of new growing point yields
(11)
X
i
(
t
+
1
)
=
X
i
(
t
)
+
R
1
*
D
i
(
φ
i
)
,
where X
i
(
t
+
1
) is the new growing point from X
i
(
t
), R
1 is the elongate-length unit, which is a random varying from 0 to 1, φ
i is the growth angle (
φ
i
1
,
φ
i
2
,
…
,
φ
i
(
D
-
1
)
), and the growth direction D
i
(
φ
i
)
=
(
d
i
1
,
d
i
2
,
…
,
d
i
D
) which can be calculated from φ
i via a Polar to Cartesian coordinates transformation [34]:
(12)
d
i
1
=
∏
k
=
1
D
-
1
cos
⁡
(
φ
i
k
)
,
d
i
j
=
sin
(
φ
i
j
)
∏
k
=
1
D
-
1
cos
⁡
⁡
(
φ
i
k
)
,
d
i
D
=
sin
(
φ
i
D
-
1
)
.
The growth angleφ
i is calculated as follows:
(13)
φ
(
t
+
1
)
=
φ
(
t
)
+
R
2
α
init
+
K
*
β
max
⁡
S
max
⁡
,
where R
2 is a random coefficient varying from 0 to 1, α
init is original growth angle value of the initial mainroot as zero degree, K is randomly parameter selecting the subzone, S
max
⁡ is subzones number, and β
max
⁡ is the maximum growing turning angle. Similar to (8), considering gravitropism, the growing angle space β
max
⁡ is limited to π to make mainroots tend to grow downwards.


## 2.4. Lateral-Roots Growth Operation

Due to the fact that the lateral-roots group with lower auxin concentration stays in the initial growth phase, its elongation step is relatively small and the influence of the hydrotropism can be ignored. The lateral-roots with the smaller elongation step are exploiting a local optimal domain, which closed to its original point.The mathematic expression, which simulates above behavior, can be defined as (14):
(14)
X
i
(
t
+
1
)
=
X
i
(
t
)
+
R
5
*
Gro
w
max
⁡
*
D
i
(
φ
i
)
,
(15)
φ
i
(
t
+
1
)
=
φ
i
(
t
)
+
R
6
α
max
⁡
,
where R
5 is normal distribution function (mean is 0, deviation is 1) and R
6 is a random coefficient.


## 2.5. Dead-Roots Elimination

In the proposed root foraging model, it is assumed thatN
i is the current population size, N
i will increase by one if a root tip splits and reduce by one if a root dies determined by auxin distribution, and it will vary in the searching process [35]. The branching criterion and dead-roots eliminating criterion are listed in
(16)
N
i
=
N
i
+
w
i
,
if
F
i
>
Branch
G
,
N
i
=
N
i
-
1
,
if
F
i
<
N
mority
,
where w
i is the branching number defined by (10), Branch
G is the branching threshold value defined by (9), and N
mority is the death threshold value parameter.


## 2.6. ARFO Algorithm

We implement an instantiation algorithm of the proposed root foraging model, namely, artificial root foraging optimizer (ARFO), for numerical function optimization. The pseudocode of the detailed procedure is given as follows.Step 1 (initialization).
Initialize the positions of root tips population and evaluate the nutrient (fitness) of the population. Set maximum number of cycles (LimitC).Step 2.
Iteration = 0.Step 3 (auxin concentration calculation).
(
1
) Calculate the nutrient concentration and auxin concentration values of the population by (2) and (3).
(
2
) Divide the population into two groups: mainroots group and lateral-roots group by (5).Step 4 (mainroots operation).
Loop over each mainroot tip(1)
Select a random portion of mainroots group to take regrowing operator by (6), (7), and (8).
(2)
Evaluate the fitness of the renewal mainroots and apply greedy selection.
(3)
Adjust the corresponding nutrient concentration value by (4).
(4)
If the condition of branching determined by (9) is met, continue; otherwise, go to Step 5.
(5)
Calculate the branching number by (10) and the branching new roots by (11).
(6)
The population size increases by the corresponding number calculated by (10).Step 5 (lateral-roots operation).
Loop over each root tip of lateral-roots(1)
Lateral-root takes regrowing operator by (14) and (15).
(2)
Evaluate the fitness of the renewal lateral-root and apply greedy selection.
(3)
Adjust the corresponding nutrient concentration value by (4).Step 6 (dead roots eliminating).
Remove the dead individuals from the population according to their auxin concentration values by (16).Step 7.
Iteration = Iteration + 1;Step 8.
If the iteration is greater than LimitC, stop the procedure; otherwise, go to Step3.Step 9.
Output the best solution achieved.By using the regrowing and branching operators, ARFO not only alleviates the premature convergence in the entire solution space, but also improves the accuracy of optimum in the local solution space, at least to some extent.(1)
Exploring strategy. The regrowing operator implemented by mainroots employs a large elongate-length unit to explore the previously unscanned regions in the search space.
(2)
Exploiting strategy. To perform fine-tuning exploitation of the global optimum, both the branching operator and the lateral-root’s regrowing operator employ a small elongate-length unit and a random angle.


## 3. Benchmark Test

### 3.1. Test Functions

In the benchmark testing studies, according to the no free lunch (NFL) theorem [36], a set of 10 benchmark functions (consisting of basic benchmark functions f
1
~
f
5 and CEC’05 functions f
6
~
f
10) was employed, which are listed as follows. The dimensions, initialization ranges, global optimum, and the criterion of each function are summarized in Table 1. f
1
~
f
5 were adopted widely in evolutionary computation domain to show solution quality and convergence rate. f
6
~
f
10 are shifted and rotated functions selected from CEC2005 functions; their global optimum values are different to each other. In order to compare the different algorithms fairly, we decide to use the number of function evaluations (FEs) as a time measure substitute the number of iterations.Table 1
Parameters of the test functions.
f
Dimensions
Initial range
x
*
f
(
x
*
)
f
1
20
[−100, 100]
D
[
0,0
,
…
,
0
]
0
f
2
20
[−30, 30]
D
[
1,1
,
…
,
1
]
0
f
3
20
[−5.12, 5.12]
D
[
0,0
,
…
,
0
]
0
f
4
20
[−500, 500]
D
[
420.9867
,
…
,
420.9867
]
0
f
5
20
[−600, 600]
D
[
0,0
,
…
,
0
]
0
f
6
20
[−100, 100]
D
[
0,0
,
…
,
0
]
-
450
f
7
20
[−100, 100]
D
[
0,0
,
…
,
0
]
390
f
8
20100
[−100, 100]
D
[
0,0
,
…
,
0
]
-
450
f
9
20
No bounds
[
0,0
,
…
,
0
]
-
180
f
1
0
20
[−5, 5]
D
[
0,0
,
…
,
0
]
-
330(
1) Sphere Function
(17)
f
1
(
x
)
=
∑
i
=
1
n
x
i
2
.(
2) Rosenbrock Function
(18)
f
2
(
x
)
=
∑
i
=
1
n
100
×
(
x
i
+
1
-
x
i
2
)
2
+
(
1
-
x
i
)
2
.(
3) Rastrigin Function
(19)
f
3
(
x
)
=
∑
i
=
1
n
100
×
(
x
i
+
1
-
x
i
2
)
2
+
(
1
-
x
i
)
2
.(
4) Schwefel Function
(20)
f
4
(
x
)
=
D
*
418.9829
+
∑
i
=
1
D
-
x
i
sin
⁡
(
|
x
i
|
)
.(
5) Griewank Function
(21)
f
5
(
x
)
=
1
4000
∑
i
=
1
n
x
i
2
-
∏
i
=
1
n
cos
⁡
(
x
i
i
)
+
1
.(
6) Shifted Sphere Function
(22)
f
6
(
x
)
=
∑
i
=
1
D
z
i
2
+
f
bias
1
,
z
=
x
-
o
.(
7) Shifted Rosenbrock’s Function
(23)
f
7
(
x
)
=
∑
i
=
1
D
-
1
(
100
(
z
i
2
-
z
i
+
1
)
2
+
(
z
i
2
-
1
)
2
)
+
f
bias
.(
8) Shifted Schwefel’s Problem 1.2
(24)
f
8
(
x
)
=
∑
i
=
1
D
(
∑
j
=
1
i
z
j
)
2
+
f
bias
,
z
=
x
-
o
.(
9) Shifted Rotated Griewank’s Function without Bounds
(25)
f
9
(
x
)
=
∑
i
=
1
D
z
i
2
4000
-
∏
i
=
1
D
cos
⁡
(
z
i
i
)
+
1
+
f
bias
.(
10) Shifted Rastrigin’s Function
(26)
f
10
(
x
)
=
∑
i
=
1
D
(
z
i
2
-
10
cos
⁡
(
2
π
z
i
)
+
10
)
+
f
bias
,
z
=
x
-
o
.


### 3.2. Parameters Settings for the Involved Algorithms

Experiments were conducted to compare the performance of ARFO with four successful evolutionary algorithms (EAs) on ten benchmark functions with 20 dimensions and 100 dimensions: standard particle swarm optimization (PSO) [8], real-coded genetic algorithm (GA) [5], different evolution algorithm (DE) [3], and artificial bee colony algorithm (ABC) [4].The experiments were run 30 times, respectively, for each algorithm on each benchmark function and the maximum evaluation number is set as 100000. The mean values and standard deviation of the results are presented. Population size in PSO, GA, DE, and ABC was 20. The PSO algorithm we used is the standard one (i.e., the global version with inertia weight). The parameters were given by the default setting of [8]; the acceleration factors c
1 and c
2 were both 2.0; a decaying inertia weight ω starting at 0.9 and ending at 0.4 was used. For DE, F was set to 0.5 and crossover rate was chosen to be 0.9 as recommended in [3]. The GA algorithm we executed is a real-coded genetic algorithm with intermediate crossover and Gaussian mutation. The control parameters, for example, mutation rate and crossover rate, were set to be the same as [5]. For ABC, the limit parameter is set to be S
N
×
D, where D is the dimension of the problem and S
N is the number of employed bees, which is half of ABC population size. The parameter setting of the ARFO algorithm uses our empirical values summarized in Table 2.Table 2
Parameters of ARFO algorithm for simulation and optimization progress.
Simulation setting
Values
Optimization setting
Values
The number of initial roots
1
The number of initial population
20
The maximum number of root tip
200
The maximum number of population
100
Cr
0.1
Cr
0.1
BranchG
100
BranchG
100
Nmority
−1
Nmority
−1
S
max
4
S
max
4
S
min
1
S
min
1
α
max
π/4
α
max
π/4
β
max
π
β
max
π


### 3.3. Results for Benchmark Functions

In Table3 reporting the means and variances of the five algorithms on the ten 30-dimensional test functions for 30 runs is illustrated where best results among involved algorithms are shown in bold. Figures 1 and 2 demonstrate the performance characteristics in terms of convergence rate regarding the best fitness value of the median run of each algorithm on basic benchmarks and CEC 2005 benchmarks, respectively.Table 3
Comparison of results with 30 dimensions obtained by ARFO, ABC, PSO, GA, and DE.
Function
ARFO
ABC
PSO
GA
DE
f
1
Mean
3.85
E
-
13
1.01
E
-
600
6.86
E
-
03
5.24
E
-
01
4.44
E
+
00
Std
5.29
E
-
13
4.33
E
-
600
5.37
E
-
03
1.86
E
-
01
3.89
E
+
00
Min
6.70
E
-
14
2.31
E
-
600
2.05
E
-
03
4.14
E
-
01
1.93
E
+
00
Max
1.18
E
-
12
1.21
E
-
500
1.54
E
-
02
6.73
E
-
01
1.09
E
+
01
f
2
Mean
6.54
E
-
01
3.22
E
+
00
5.63
E
+
01
2.16
E
+
02
4.83
E
+
04
Std
1.43
E
+
00
3.06
E
+
00
5.54
E
+
01
1.15
E
+
02
5.40
E
+
04
Min
0
1.56
E
-
01
1.68
E
+
01
7.17
E
+
01
1.34
E
+
04
Max
3.25
E
+
00
7.51
E
+
00
1.46
E
+
02
3.91
E
+
02
1.38
E
+
05
f
3
Mean
1.94
E
+
00
4.02
E
+
01
6.15
E
+
01
1.08
E
+
02
2.21
E
+
02
Std
1.10
E
+
00
1.06
E
+
01
2.72
E
+
01
3.09
E
+
02
9.41
E
+
01
Min
8.19
E
-
01
2.35
E
+
01
4.44
E
+
01
5.89
E
+
01
1.09
E
+
02
Max
3.25
E
+
00
5.04
E
+
01
1.10
E
+
02
1.34
E
+
02
3.63
E
+
02
f
4
Mean
6.88
E
-
04
2.32
E
+
03
3.62
E
+
03
2.43
E
+
02
2.43
E
+
03
Std
5.74
E
-
04
2.24
E
+
02
1.82
E
+
02
1.78
E
+
02
4.91
E
+
02
Min
2.15
E
-
04
2.22
E
+
03
3.40
E
+
03
1.22
E
+
02
1.74
E
+
03
Max
1.38
E
-
03
2.48
E
+
03
3.92
E
+
03
3.51
E
+
02
2.95
E
+
03
f
5
Mean
4.51
E
-
03
6.81
E
-
02
7.54
E
-
01
3.23
E
+
00
2.95
E
+
01
Std
4.42
E
-
03
5.96
E
-
02
2.33
E
-
01
6.27
E
-
02
2.19
E
+
01
Min
0
1.20
E
-
02
5.27
E
-
01
1.88
E
+
00
8.04
E
+
00
Max
1
.
04
E
-
02
1.72
E
-
01
9.86
E
-
01
3.76
E
+
00
4.89
E
+
01
f
6
Mean
6.80
E
-
14
3.29
E
-
14
5.12
E
+
02
7.72
E
+
01
5.25
E
+
03
Std
2.77
E
-
14
2.84
E
-
14
6.36
E
+
02
3.21
E
+
01
3.24
E
+
03
Min
4.85
E
-
14
1.49
E
-
14
1.95
E
+
02
6.40
E
+
01
8.39
E
+
02
Max
1.01
E
-
13
5.10
E
-
14
1.27
E
+
03
1.20
E
+
02
1.23
E
+
04
f
7
Mean
1.32
E
+
00
6.88
E
+
00
5.71
E
+
06
1.99
E
+
03
2.30
E
+
09
Std
1.83
E
+
00
3.47
E
+
00
1.88
E
+
07
5.89
E
+
03
2.50
E
+
09
Min
4.85
E
-
04
3.93
E
+
00
3.82
E
+
02
3.15
E
+
03
2.60
E
+
08
Max
3.38
E
+
00
1.21
E
+
01
3.13
E
+
07
3.05
E
+
03
3.60
E
+
09
f
8
Mean
1.76
E
+
02
7.36
E
+
02
1.58
E
+
02
1.71
E
+
04
1.20
E
+
04
Std
3.31
E
+
02
4.74
E
+
02
4.54
E
+
01
3.96
E
+
03
3.61
E
+
04
Min
1.76
E
+
01
5.29
E
+
02
1.39
E
+
02
1.11
E
+
04
8.71
E
+
03
Max
3.39
E
+
02
1.03
E
+
03
1.84
E
+
02
1.69
E
+
04
1.46
E
+
04
f
9
Mean
1.44
E
+
03
1.70
E
+
03
4.26
E
+
03
2.06
E
+
03
1.75
E
+
03
Std
6.29
E
-
01
6.07
E
-
13
4.25
E
+
02
3.14
E
+
02
4.97
E
+
01
Min
1.41
E
+
03
1.64
E
+
03
4.01
E
+
03
1.97
E
+
03
1.64
E
+
03
Max
1.48
E
+
03
1.77
E
+
03
4.52
E
+
03
2.16
E
+
03
1.81
E
+
03
f
1
0
Mean
6.12
E
+
00
4.79
E
+
01
2.83
E
+
02
5.41
E
+
01
1.64
E
+
02
Std
5.81
E
-
01
1.44
E
+
01
5.34
E
+
01
4.70
E
+
00
2.92
E
+
01
Min
3.80
E
+
00
3.15
E
+
01
7.08
E
+
01
5.09
E
+
01
1.48
E
+
02
Max
7.14
E
+
00
5.81
E
+
01
3.76
E
+
02
6.17
E
+
01
1.79
E
+
02Convergence results of ARFO, ABC, PSO, GA, and DE. (a)–(e) correspond to 30-dimensionalf
1
–
f
5, respectively.
(a)
(b)
(c)
(d)
(e)Convergence results of ARFO, ABC, PSO, GA and DE. (a)–(e) correspond to30-dimensionalf
6
–
f
10, respectively.
(a)
(b)
(c)
(d)
(e)On the two unimodal functions, namely, Sphere and Rosenbrock benchmarks, both ARFO and ABC demonstrated competent performance. It is worthy noted ABC algorithm can consistently obtain the global minimum every time on Sphere. However, for the complex unimodal problem-Rosenbrock, ARFO markedly performed more powerfully than other cases.In terms of complex multimodal benchmarks including variable-separable Rastrigin, variable-separable Schwefel, and nonseparable Griewank, as shown in Figures1(c)–1(e), the convergence profiles of ARFO and ABC on these functions were significantly better than the other three algorithms. IN particular, on Schwefel and Griewank, the search performance order is ARFO
>
ABC
>
PSO
>
GA
>
DE. But ABC showed better performance than ARFO on Rastrigin.Seven shifted and rotated functionsf
6
~
f
10 are considered more difficult to be optimized. From Table 3, ARFO obtained best results on f
7, f
9, and f
10 and the search performance order is ARFO
>
ABC
>
PSO
>
GA
>
DE. PSO also achieved better results on f
8. From the graph in terms of convergence rate in Figure 2, it is clearly visible that our proposed ARFO and ABC algorithms achieved significantly better performance on most CEC 2005 benchmarks than PSO, DE, and GA. Meantime ARFO performed a little better than ABC except f
6 when other algorithms miss the global optimum basin. The shifted rotated Griewank function is a good example that ARFO successfully avoids falling into local optima and continues to find better results even after the PSO, GA, and DE seem to have stagnated.In general, ARFO exhibits more significant advantages in terms of accuracy and convergence rate on most of the tested benchmark functions. This can be explained that ARFO can achieve an appropriate balance between exploration and exploitation by incorporating the regrowing operator and the branching operator because the regrowing operation employs a large elongate-length unit to explore the previously unscanned regions while the branching uses a small elongation-length unit to perform fine-tuning exploitation of the global optimum. Furthermore, the branching operator can generate new individuals to maintain its population diversity, which is a key contributing factor.


## 3.1. Test Functions

In the benchmark testing studies, according to the no free lunch (NFL) theorem [36], a set of 10 benchmark functions (consisting of basic benchmark functions f
1
~
f
5 and CEC’05 functions f
6
~
f
10) was employed, which are listed as follows. The dimensions, initialization ranges, global optimum, and the criterion of each function are summarized in Table 1. f
1
~
f
5 were adopted widely in evolutionary computation domain to show solution quality and convergence rate. f
6
~
f
10 are shifted and rotated functions selected from CEC2005 functions; their global optimum values are different to each other. In order to compare the different algorithms fairly, we decide to use the number of function evaluations (FEs) as a time measure substitute the number of iterations.Table 1
Parameters of the test functions.
f
Dimensions
Initial range
x
*
f
(
x
*
)
f
1
20
[−100, 100]
D
[
0,0
,
…
,
0
]
0
f
2
20
[−30, 30]
D
[
1,1
,
…
,
1
]
0
f
3
20
[−5.12, 5.12]
D
[
0,0
,
…
,
0
]
0
f
4
20
[−500, 500]
D
[
420.9867
,
…
,
420.9867
]
0
f
5
20
[−600, 600]
D
[
0,0
,
…
,
0
]
0
f
6
20
[−100, 100]
D
[
0,0
,
…
,
0
]
-
450
f
7
20
[−100, 100]
D
[
0,0
,
…
,
0
]
390
f
8
20100
[−100, 100]
D
[
0,0
,
…
,
0
]
-
450
f
9
20
No bounds
[
0,0
,
…
,
0
]
-
180
f
1
0
20
[−5, 5]
D
[
0,0
,
…
,
0
]
-
330(
1) Sphere Function
(17)
f
1
(
x
)
=
∑
i
=
1
n
x
i
2
.(
2) Rosenbrock Function
(18)
f
2
(
x
)
=
∑
i
=
1
n
100
×
(
x
i
+
1
-
x
i
2
)
2
+
(
1
-
x
i
)
2
.(
3) Rastrigin Function
(19)
f
3
(
x
)
=
∑
i
=
1
n
100
×
(
x
i
+
1
-
x
i
2
)
2
+
(
1
-
x
i
)
2
.(
4) Schwefel Function
(20)
f
4
(
x
)
=
D
*
418.9829
+
∑
i
=
1
D
-
x
i
sin
⁡
(
|
x
i
|
)
.(
5) Griewank Function
(21)
f
5
(
x
)
=
1
4000
∑
i
=
1
n
x
i
2
-
∏
i
=
1
n
cos
⁡
(
x
i
i
)
+
1
.(
6) Shifted Sphere Function
(22)
f
6
(
x
)
=
∑
i
=
1
D
z
i
2
+
f
bias
1
,
z
=
x
-
o
.(
7) Shifted Rosenbrock’s Function
(23)
f
7
(
x
)
=
∑
i
=
1
D
-
1
(
100
(
z
i
2
-
z
i
+
1
)
2
+
(
z
i
2
-
1
)
2
)
+
f
bias
.(
8) Shifted Schwefel’s Problem 1.2
(24)
f
8
(
x
)
=
∑
i
=
1
D
(
∑
j
=
1
i
z
j
)
2
+
f
bias
,
z
=
x
-
o
.(
9) Shifted Rotated Griewank’s Function without Bounds
(25)
f
9
(
x
)
=
∑
i
=
1
D
z
i
2
4000
-
∏
i
=
1
D
cos
⁡
(
z
i
i
)
+
1
+
f
bias
.(
10) Shifted Rastrigin’s Function
(26)
f
10
(
x
)
=
∑
i
=
1
D
(
z
i
2
-
10
cos
⁡
(
2
π
z
i
)
+
10
)
+
f
bias
,
z
=
x
-
o
.


## 3.2. Parameters Settings for the Involved Algorithms

Experiments were conducted to compare the performance of ARFO with four successful evolutionary algorithms (EAs) on ten benchmark functions with 20 dimensions and 100 dimensions: standard particle swarm optimization (PSO) [8], real-coded genetic algorithm (GA) [5], different evolution algorithm (DE) [3], and artificial bee colony algorithm (ABC) [4].The experiments were run 30 times, respectively, for each algorithm on each benchmark function and the maximum evaluation number is set as 100000. The mean values and standard deviation of the results are presented. Population size in PSO, GA, DE, and ABC was 20. The PSO algorithm we used is the standard one (i.e., the global version with inertia weight). The parameters were given by the default setting of [8]; the acceleration factors c
1 and c
2 were both 2.0; a decaying inertia weight ω starting at 0.9 and ending at 0.4 was used. For DE, F was set to 0.5 and crossover rate was chosen to be 0.9 as recommended in [3]. The GA algorithm we executed is a real-coded genetic algorithm with intermediate crossover and Gaussian mutation. The control parameters, for example, mutation rate and crossover rate, were set to be the same as [5]. For ABC, the limit parameter is set to be S
N
×
D, where D is the dimension of the problem and S
N is the number of employed bees, which is half of ABC population size. The parameter setting of the ARFO algorithm uses our empirical values summarized in Table 2.Table 2
Parameters of ARFO algorithm for simulation and optimization progress.
Simulation setting
Values
Optimization setting
Values
The number of initial roots
1
The number of initial population
20
The maximum number of root tip
200
The maximum number of population
100
Cr
0.1
Cr
0.1
BranchG
100
BranchG
100
Nmority
−1
Nmority
−1
S
max
4
S
max
4
S
min
1
S
min
1
α
max
π/4
α
max
π/4
β
max
π
β
max
π


## 3.3. Results for Benchmark Functions

In Table3 reporting the means and variances of the five algorithms on the ten 30-dimensional test functions for 30 runs is illustrated where best results among involved algorithms are shown in bold. Figures 1 and 2 demonstrate the performance characteristics in terms of convergence rate regarding the best fitness value of the median run of each algorithm on basic benchmarks and CEC 2005 benchmarks, respectively.Table 3
Comparison of results with 30 dimensions obtained by ARFO, ABC, PSO, GA, and DE.
Function
ARFO
ABC
PSO
GA
DE
f
1
Mean
3.85
E
-
13
1.01
E
-
600
6.86
E
-
03
5.24
E
-
01
4.44
E
+
00
Std
5.29
E
-
13
4.33
E
-
600
5.37
E
-
03
1.86
E
-
01
3.89
E
+
00
Min
6.70
E
-
14
2.31
E
-
600
2.05
E
-
03
4.14
E
-
01
1.93
E
+
00
Max
1.18
E
-
12
1.21
E
-
500
1.54
E
-
02
6.73
E
-
01
1.09
E
+
01
f
2
Mean
6.54
E
-
01
3.22
E
+
00
5.63
E
+
01
2.16
E
+
02
4.83
E
+
04
Std
1.43
E
+
00
3.06
E
+
00
5.54
E
+
01
1.15
E
+
02
5.40
E
+
04
Min
0
1.56
E
-
01
1.68
E
+
01
7.17
E
+
01
1.34
E
+
04
Max
3.25
E
+
00
7.51
E
+
00
1.46
E
+
02
3.91
E
+
02
1.38
E
+
05
f
3
Mean
1.94
E
+
00
4.02
E
+
01
6.15
E
+
01
1.08
E
+
02
2.21
E
+
02
Std
1.10
E
+
00
1.06
E
+
01
2.72
E
+
01
3.09
E
+
02
9.41
E
+
01
Min
8.19
E
-
01
2.35
E
+
01
4.44
E
+
01
5.89
E
+
01
1.09
E
+
02
Max
3.25
E
+
00
5.04
E
+
01
1.10
E
+
02
1.34
E
+
02
3.63
E
+
02
f
4
Mean
6.88
E
-
04
2.32
E
+
03
3.62
E
+
03
2.43
E
+
02
2.43
E
+
03
Std
5.74
E
-
04
2.24
E
+
02
1.82
E
+
02
1.78
E
+
02
4.91
E
+
02
Min
2.15
E
-
04
2.22
E
+
03
3.40
E
+
03
1.22
E
+
02
1.74
E
+
03
Max
1.38
E
-
03
2.48
E
+
03
3.92
E
+
03
3.51
E
+
02
2.95
E
+
03
f
5
Mean
4.51
E
-
03
6.81
E
-
02
7.54
E
-
01
3.23
E
+
00
2.95
E
+
01
Std
4.42
E
-
03
5.96
E
-
02
2.33
E
-
01
6.27
E
-
02
2.19
E
+
01
Min
0
1.20
E
-
02
5.27
E
-
01
1.88
E
+
00
8.04
E
+
00
Max
1
.
04
E
-
02
1.72
E
-
01
9.86
E
-
01
3.76
E
+
00
4.89
E
+
01
f
6
Mean
6.80
E
-
14
3.29
E
-
14
5.12
E
+
02
7.72
E
+
01
5.25
E
+
03
Std
2.77
E
-
14
2.84
E
-
14
6.36
E
+
02
3.21
E
+
01
3.24
E
+
03
Min
4.85
E
-
14
1.49
E
-
14
1.95
E
+
02
6.40
E
+
01
8.39
E
+
02
Max
1.01
E
-
13
5.10
E
-
14
1.27
E
+
03
1.20
E
+
02
1.23
E
+
04
f
7
Mean
1.32
E
+
00
6.88
E
+
00
5.71
E
+
06
1.99
E
+
03
2.30
E
+
09
Std
1.83
E
+
00
3.47
E
+
00
1.88
E
+
07
5.89
E
+
03
2.50
E
+
09
Min
4.85
E
-
04
3.93
E
+
00
3.82
E
+
02
3.15
E
+
03
2.60
E
+
08
Max
3.38
E
+
00
1.21
E
+
01
3.13
E
+
07
3.05
E
+
03
3.60
E
+
09
f
8
Mean
1.76
E
+
02
7.36
E
+
02
1.58
E
+
02
1.71
E
+
04
1.20
E
+
04
Std
3.31
E
+
02
4.74
E
+
02
4.54
E
+
01
3.96
E
+
03
3.61
E
+
04
Min
1.76
E
+
01
5.29
E
+
02
1.39
E
+
02
1.11
E
+
04
8.71
E
+
03
Max
3.39
E
+
02
1.03
E
+
03
1.84
E
+
02
1.69
E
+
04
1.46
E
+
04
f
9
Mean
1.44
E
+
03
1.70
E
+
03
4.26
E
+
03
2.06
E
+
03
1.75
E
+
03
Std
6.29
E
-
01
6.07
E
-
13
4.25
E
+
02
3.14
E
+
02
4.97
E
+
01
Min
1.41
E
+
03
1.64
E
+
03
4.01
E
+
03
1.97
E
+
03
1.64
E
+
03
Max
1.48
E
+
03
1.77
E
+
03
4.52
E
+
03
2.16
E
+
03
1.81
E
+
03
f
1
0
Mean
6.12
E
+
00
4.79
E
+
01
2.83
E
+
02
5.41
E
+
01
1.64
E
+
02
Std
5.81
E
-
01
1.44
E
+
01
5.34
E
+
01
4.70
E
+
00
2.92
E
+
01
Min
3.80
E
+
00
3.15
E
+
01
7.08
E
+
01
5.09
E
+
01
1.48
E
+
02
Max
7.14
E
+
00
5.81
E
+
01
3.76
E
+
02
6.17
E
+
01
1.79
E
+
02Convergence results of ARFO, ABC, PSO, GA, and DE. (a)–(e) correspond to 30-dimensionalf
1
–
f
5, respectively.
(a)
(b)
(c)
(d)
(e)Convergence results of ARFO, ABC, PSO, GA and DE. (a)–(e) correspond to30-dimensionalf
6
–
f
10, respectively.
(a)
(b)
(c)
(d)
(e)On the two unimodal functions, namely, Sphere and Rosenbrock benchmarks, both ARFO and ABC demonstrated competent performance. It is worthy noted ABC algorithm can consistently obtain the global minimum every time on Sphere. However, for the complex unimodal problem-Rosenbrock, ARFO markedly performed more powerfully than other cases.In terms of complex multimodal benchmarks including variable-separable Rastrigin, variable-separable Schwefel, and nonseparable Griewank, as shown in Figures1(c)–1(e), the convergence profiles of ARFO and ABC on these functions were significantly better than the other three algorithms. IN particular, on Schwefel and Griewank, the search performance order is ARFO
>
ABC
>
PSO
>
GA
>
DE. But ABC showed better performance than ARFO on Rastrigin.Seven shifted and rotated functionsf
6
~
f
10 are considered more difficult to be optimized. From Table 3, ARFO obtained best results on f
7, f
9, and f
10 and the search performance order is ARFO
>
ABC
>
PSO
>
GA
>
DE. PSO also achieved better results on f
8. From the graph in terms of convergence rate in Figure 2, it is clearly visible that our proposed ARFO and ABC algorithms achieved significantly better performance on most CEC 2005 benchmarks than PSO, DE, and GA. Meantime ARFO performed a little better than ABC except f
6 when other algorithms miss the global optimum basin. The shifted rotated Griewank function is a good example that ARFO successfully avoids falling into local optima and continues to find better results even after the PSO, GA, and DE seem to have stagnated.In general, ARFO exhibits more significant advantages in terms of accuracy and convergence rate on most of the tested benchmark functions. This can be explained that ARFO can achieve an appropriate balance between exploration and exploitation by incorporating the regrowing operator and the branching operator because the regrowing operation employs a large elongate-length unit to explore the previously unscanned regions while the branching uses a small elongation-length unit to perform fine-tuning exploitation of the global optimum. Furthermore, the branching operator can generate new individuals to maintain its population diversity, which is a key contributing factor.


## 4. Multilevel Threshold for Image Segmentation by ARFO

### 4.1. Entropy Criterion Based Fitness Measure

The Otsu multithreshold entropy measure [37] proposed by Otsu has been popularly employed in determining whether the optimal threshold method can provide image segmentation with satisfactory results. Here, it is used as the objective function for the involved algorithms and its process can be described as followsLet the gray levels of a given image range over [0,L
-
1] and h
(
i
) denote the occurrence of gray-level i.Let(27)
N
=
∑
i
=
0
L
-
1
h
(
i
)
,
P
(
i
)
=
h
(
i
)
N
,
for
0
≤
i
≤
L
-
1
Maximize
f
(
t
)
=
w
0
w
1
(
u
0
-
u
1
)
2
,
where
(28)
w
0
=
∑
i
=
0
t
-
1
P
i
;
u
0
=
∑
i
=
0
t
=
1
i
×
P
i
w
0
(29)
w
1
=
∑
i
=
t
L
-
1
P
i
;
u
1
=
∑
i
=
t
L
=
1
i
×
P
i
w
1
and the optimal threshold is the gray level that maximizes (27). Then (28) can also be written as
(30)
f
(
t
)
=
δ
2
-
w
0
×
δ
0
2
-
w
1
×
δ
1
2
,
where w
0, w
1, u
0, and u
1 are the same as given in (9) and (10), and
(31)
δ
0
=
∑
i
=
0
t
-
1
(
i
-
u
0
)
2
×
P
i
w
0
,
δ
1
=
∑
i
=
t
L
-
1
(
i
-
u
1
)
2
×
P
i
w
1
,
δ
0
=
∑
i
=
0
L
-
1
(
i
-
u
)
2
×
P
i
w
,
w
=
∑
i
=
0
L
-
1
P
i
,
u
=
∑
i
=
0
L
=
1
i
×
P
i
w
.
Expanding this logic to multilevel threshold,
(32)
f
12
(
t
)
=
w
0
×
δ
0
2
+
w
1
×
δ
1
2
+
w
2
×
δ
2
2
+
⋯
+
w
1
×
δ
N
2
,
where N is the number of thresholds.Equation (32) is considered as the objective function for the proposed ARFO based procedure which is to be optimized (minimized). A close look into this equation will show it is very similar to the expression for uniformity measure.


### 4.2. Experiment Setup

As illustrated above, the performance evaluation based on entropy criterion is one of key indicators in image segmentation scenarios, our experimental studies are conducted on a suit of widely used image datasets [38–43], which contain avion.ppm, house.ppm, lena.ppm, peppers.ppm, safari04.ppm, and hunter.pgm (available at http://decsai.ugr.es/cvg/dbimagenes/). The size of each involved image is 512 * 512, with a unique grey-level histogram. The intrinsic multimodality of these tested images’ histograms will result in a more difficult optimization process. The segmentation results of the proposed algorithm and other algorithms are evaluated by the traditional Otsu method. The parameters of these algorithms (namely ARFO, ABC, PSO, GA, and DE) can be referred in Section 3.2, which have been proven as the optimal setting.The numbers of thresholdsM
-
1 investigated in the experiments were 2, 3, 4, 5, 7, and 9. The population size is 20. The maximum number of FEs is 2000, which has sufficient time to search on high dimension problems. Since the population-based algorithms have the randomized characteristics, all the experiments were repeated 30 times for each image for each M
-
1 value. In all experiments, to ensure that the initial values of each random algorithm are equal, we used the MATLAB command rand(
X
,
Y
)
*
(
max
⁡
L
-
min
⁡
L
)
+
min
⁡
L. The max
⁡
L and min
⁡
L are the maximum and minimum grays of the tested images. Figure 3 presents four original images and their histograms.Figure 3
Test images and their histograms.


### 4.3. Experimental Results of Multilevel Threshold

It is worthy noted the involved population-based algorithms for image segmentation, including ARFO, may only reduce the computation consumption but fail to significantly improve the accuracy of the optimal result. Therefore, we strived to utilize the proposed algorithm to obtain multiple thresholds with larger fitness values and fast computation ability. Due to the exhaustive search characteristic of classical Otsu method, its result can be regarded as the basis for comparison with other EA-based algorithms. We utilized (32) as the fitness function to evaluate involved algorithms.Case 1 (multilevel threshold results and efficiency and different methods with M − 1 = 2, 3, 4).
In tables reporting the fitness, mean computation time, and corresponding optimal thresholds (withM
-
1 = 2, 3, 4) obtained by Otsu are presented in Table 4. Due to the long consumption of CPU time when M
-
1
>
4, the correlative values are not worthy being listed in our experiment. However, other results attained by the population-based algorithms were compared. It should be noted that in the real-time image proceeding scenarios an appropriate balance of running time and high accuracy is demanded seriously; the termed CPU time of the population-based methods should be analyzed in detail. As shown in Table 5, there are no apparent differences about CPU times among the involved population-based methods, which are significantly superior in terms of time complexity for high-dimensional image segmentation problems.
From Table6, it can be observed that the termed fitness values result achieved by the proposed ARFO is generally closed to that by Otsu when M
-
1 = 2, 3, 4, whereas the performance regarding computation time of ARFO is significantly better than its counterpart Otsu. Furthermore, the ARFO-based method also performed most powerful among the population-based methods in most involved datasets. This can be explained by the fact that the ARFO has an appropriate balance between exploration and exploitation. Compared with other population-based algorithms, the ARFO is endowed with branching operations to vary population according to the nutrient concentration change to facilitate better exploitation search in complex higher dimensional search space as well as the branching operation will be used to fine tune the best solutions.Table 4
Objective values and thresholds by the Otsu method.
Image
M
-
1
=
2
M
-
1
=
3
M
-
1
=
4
Objective values
Optimal thresholds
Objective values
Optimal thresholds
Objective values
Optimal thresholds
Avion
4.11
E
+
04
114, 175
4.35
E
+
04
95, 151, 192
3.91
E
+
04
85, 131, 172, 208
House
3.04
E
+
04
109, 175
3.15
E
+
04
87, 136, 183
2.69
E
+
04
76, 121, 155, 191
Lena
1.77
E
+
04
95, 156
2.11
E
+
04
84, 123, 175
1.96
E
+
04
76, 121, 155, 191
Hunter
1.12
E
+
04
53, 118
1.13
E
+
04
37, 88, 144
1.07
E
+
04
33, 76, 121, 156
Mean CPU time
3.99
E
+
00
8.69
E
+
01
2.42
E
+
03Table 5
The mean CPU time of the compared population-based methods on Otsu algorithm.
Dim
Algorithm
ARFO
ABC
PSO
DE
GA
2 Dim
2.54
E
-
01
3.01
E
-
01
1.04
E
+
00
9.55
E
-
01
3.01
E
+
00
3 Dim
2.92
E
-
01
3.44
E
-
01
1.10
E
+
00
9.96
E
-
01
3.76
E
+
00
4 Dim
2.86
E
-
01
3.64
E
-
01
1.11
E
+
00
1.05
E
+
00
5.26
E
+
00
5 Dim
3.26
E
-
01
3.76
E
-
01
1.00
E
+
00
1.10
E
+
00
1.41
E
+
00
7 Dim
3.58
E
-
01
5.26
E
-
01
1.21
E
+
00
1.19
E
+
00
2.04
E
+
00
9 Dim
5.18
E
-
01
6.46
E
-
01
1.23
E
+
00
1.11
E
+
00
2.81
E
+
00Table 6
Objective value and standard deviation by the compared population-based methods on Otsu algorithm.
Image
M
-
1
Objective values (standard deviation)
ARFO
ABC
PSO
DE
GA
Avion
2
3.84
E
+
04
3.32
E
-
01
3.70
E
+
04
8.49
E
-
12
3.83
E
+
04
1.58
E
-
02
3.82
E
+
04
1.58
E
-
01
3.79
E
+
04
3.92
E
-
01
3
3.84
E
+
04
1.82
E
-
01
3.75
E
+
04
2.04
E
-
01
3.73
E
+
04
1.38
E
-
01
3.84
E
+
04
1.58
E
+
00
3.79
E
+
04
2.82
E
-
01
4
3.84
E
+
04
5.23
E
-
01
3.79
E
+
04
5.39
E
-
01
3.81
E
+
04
8.18
E
-
01
3.73
E
+
04
3.92
E
-
01
3.77
E
+
04
1.63
E
+
00
House
2
3.16
E
+
04
3.36
E
-
02
3.04
E
+
04
4.80
E
-
01
3.11
E
+
04
4.74
E
-
02
3.11
E
+
04
3.78
E
-
03
3.11
E
+
04
1.05
E
+
00
3
3.18
E
+
04
2.99
E
-
02
3.11
E
+
04
3.72
E
-
02
3.07
E
+
04
2.04
E
-
01
3.15
E
+
04
4.72
E
-
01
3.11
E
+
04
7.19
E
-
01
4
3.23
E
+
04
9.14
E
-
01
3.13
E
+
04
5.42
E
-
01
3.15
E
+
04
7.04
E
-
01
3.15
E
+
04
3.61
E
-
01
3.12
E
+
04
5.08
E
+
05
Lena
2
1.96
E
+
04
9.90
E
-
03
1.92
E
+
04
4.17
E
-
12
1.92
E
+
04
2.40
E
-
02
1.93
E
+
04
3.84
E
-
12
1.91
E
+
04
5.18
E
-
01
3
1.98
E
+
04
8.12
E
-
03
1.94
E
+
04
5.94
E
-
02
1.96
E
+
04
3.60
E
-
01
1.96
E
+
04
4.82
E
-
02
1.94
E
+
04
1.28
E
+
00
4
1.98
E
+
04
2.34
E
-
01
1.94
E
+
04
1.09
E
+
00
1.93
E
+
04
4.04
E
+
00
1.95
E
+
04
3.92
E
-
01
1.93
E
+
04
2.15
E
+
00
Peppers
2
1.93
E
+
04
3.78
E
-
01
1.84
E
+
04
4.25
E
-
12
1.89
E
+
04
3.58
E
-
01
1.87
E
+
04
4.71
E
-
12
1.86
E
+
04
1.71
E
-
01
3
1.94
E
+
04
9.03
E
-
01
1.87
E
+
04
1.50
E
-
02
1.90
E
+
04
1.84
E
-
01
1.89
E
+
04
5.07
E
-
01
1.89
E
+
04
3.82
E
-
01
4
1.95
E
+
04
8.79
E
-
01
1.88
E
+
04
1.50
E
-
01
1.93
E
+
04
1.51
E
+
00
1.93
E
+
04
7.17
E
-
12
1.91
E
+
04
4.27
E
+
00
Safari04
2
2.26
E
+
04
4.81
E
-
01
2.23
E
+
04
9.61
E
-
01
2.23
E
+
04
2.50
E
-
01
2.14
E
+
04
3.66
E
-
01
2.23
E
+
04
3.93
E
-
01
3
2.28
E
+
04
8.34
E
-
01
2.24
E
+
04
3.79
E
-
02
2.23
E
+
04
1.62
E
-
01
2.15
E
+
04
6.49
E
-
12
2.23
E
+
04
5.03
E
-
01
4
2.27
E
+
04
5.70
E
-
01
2.23
E
+
04
1.11
E
+
00
2.25
E
+
04
3.86
E
-
01
2.24
E
+
04
1.46
E
+
00
2.23
E
+
04
5.03
E
-
01
Hunter
2
1.02
E
+
04
1.25
E
-
01
1.01
E
+
04
3.56
E
-
02
1.00
E
+
04
4.74
E
-
01
1.01
E
+
04
2.22
E
-
04
1.01
E
+
04
6.26
E
-
01
3
1.05
E
+
04
1.34
E
+
00
1.02
E
+
04
7.27
E
-
01
1.02
E
+
04
7.07
E
-
01
1.02
E
+
04
2.57
E
+
00
1.03
E
+
04
1.62
E
+
00
4
1.05
E
+
04
9.16
E
-
01
1.03
E
+
04
3.60
E
+
00
1.03
E
+
04
1.99
E
+
00
1.04
E
+
04
3.63
E
-
01
1.04
E
+
04
3.82
E
+
00Case 2 (multilevel threshold results and efficiency and different methods withM
-
1 = 5, 7, 9).
In this section, these algorithms are conducted on image segmentation withM
-
1 = 5, 7, 9, which can be considered as higher dimensional problems. In Table 7 reporting the termed experimental results of average fitness value and standard deviation obtained by each population-based algorithm is illustrated, where the larger values and smaller standard deviations indicate the better achievement. From Table 7, it is clearly visible that there are statistically significant differences among experimental results achieved by these population-based approaches, in terms of both accuracy (fitness values) and stability (standard deviation). Depending on the population-varying model and the fast convergence ability, ARFO exhibits more efficient performances than other classical population-based algorithms, especially on high dimensional problems. The results as illustrated in Figure 4 prove that the ARFO-based algorithm is one of the suitable resolutions for multilevel image segmentation problems.Table 7
Objective value and standard deviation by the compared population-based methods on Otsu algorithm.
Image
M
-
1
Objective values (standard deviation)
ARFO
ABC
PSO
DE
GA
Avion
5
3.78
E
+
04
3.73
E
+
04
3.76
E
+
04
3.77
E
+
04
3.73
E
+
04
3.53
E
-
01
4.63
E
+
00
1.47
E
+
00
2.23
E
+
00
2.41
E
+
00
7
3.77
E
+
04
3.74
E
+
04
3.76
E
+
04
3.76
E
+
04
3.75
E
+
04
3.78
E
-
01
4.75
E
+
00
1.47
E
+
00
3.73
E
+
00
2.83
E
+
00
9
3.78
E
+
04
3.76
E
+
04
3.75
E
+
04
3.77
E
+
04
3.74
E
+
04
1.45
E
+
00
2.35
E
+
00
2.68
E
+
00
3.36
E
+
00
1.59
E
-
01
House
5
3.14
E
+
04
3.08
E
+
04
3.11
E
+
04
3.12
E
+
04
3.09
E
+
04
1.46
E
+
00
1.70
E
+
00
1.58
E
+
00
2.22
E
+
00
3.65
E
+
00
7
3.15
E
+
04
3.12
E
+
04
3.10
E
+
04
3.13
E
+
04
3.10
E
+
04
1.45
E
-
01
5.97
E
+
00
1.46
E
+
00
2.99
E
+
00
3.76
E
+
00
9
3.16
E
+
04
3.11
E
+
04
3.09
E
+
04
3.13
E
+
04
3.10
E
+
04
1.12
E
+
00
5.85
E
+
00
1.25
E
+
00
1.03
E
+
00
1.29
E
+
00
Lena
5
1.97
E
+
04
1.94
E
+
04
1.91
E
+
04
1.95
E
+
04
1.92
E
+
04
1.58
E
-
01
2.13
E
+
00
1.39
E
+
00
2.78
E
+
00
2.60
E
+
00
7
1.97
E
+
04
1.94
E
+
04
1.95
E
+
04
1.95
E
+
04
1.93
E
+
04
1.15
E
+
00
5.74
E
+
00
2.59
E
+
00
1.79
E
-
01
1.48
E
+
00
9
1.96
E
+
04
1.93
E
+
04
1.94
E
+
04
1.95
E
+
04
1.93
E
+
04
1.57
E
+
00
6.83
E
+
00
3.76
E
+
00
1.23
E
+
00
3.66
E
+
00
Peppers
5
1.93
E
+
04
1.89
E
+
04
1.89
E
+
04
1.92
E
+
04
1.88
E
+
04
5.06
E
-
01
1.57
E
+
00
1.79
E
+
00
2.19
E
+
00
4.31
E
+
00
7
1.93
E
+
04
1.88
E
+
04
1.89
E
+
04
1.92
E
+
04
1.88
E
+
04
1.39
E
+
00
2.67
E
+
00
4.01
E
+
00
5.07
E
-
01
2.83
E
+
00
9
1.92
E
+
04
1.89
E
+
04
1.91
E
+
04
1.92
E
+
04
1.90
E
+
04
1.20
E
+
00
2.69
E
+
00
2.01
E
+
00
5.07
E
-
01
4.39
E
+
00
Safari04
5
2.35
E
+
04
2.20
E
+
04
2.21
E
+
04
2.04
45
E
4
2.20
E
+
04
1.03
E
-
01
6.73
E
+
00
1.46
E
+
00
2.23
E
+
04
2.19
E
+
00
7
2.34
E
+
04
2.20
E
+
04
2.22
E
+
04
3.54
E
+
00
2.20
E
+
04
1.48
E
+
00
4.75
E
+
00
1.11
E
+
00
2.23
E
+
04
3.74
E
+
00
9
2.25
E
+
04
2.20
E
+
04
2.21
E
+
04
2.30
E
+
00
2.20
E
+
04
9.68
E
-
01
1.43
E
+
00
5.18
E
-
01
3.65
E
+
00
2.34
E
+
00
Hunter
5
1.03
E
+
04
1.02
E
+
04
1.02
E
+
04
1.03
E
+
04
1.02
E
+
04
7.91
E
-
01
4.32
E
+
00
3.65
E
+
00
4.95
E
-
01
3.61
E
+
00
7
1.05
E
+
04
1.02
E
+
04
1.02
E
+
04
1.03
E
+
04
1.03
E
+
04
1.83
E
+
00
3.98
E
+
00
2.35
E
+
00
1.76
E
+
00
2.66
E
+
00
9
1.04
E
+
04
1.03
E
+
04
1.03
E
+
04
1.03
E
+
04
1.03
E
+
04
1.37
E
+
00
5.85
E
+
00
2.45
E
+
00
4.94
E
-
01
3.98
E
+
00Figure 4
The comparison of the thresholded images by ARFO-based and ABC-based algorithms (M
-
1
=
9).


## 4.1. Entropy Criterion Based Fitness Measure

The Otsu multithreshold entropy measure [37] proposed by Otsu has been popularly employed in determining whether the optimal threshold method can provide image segmentation with satisfactory results. Here, it is used as the objective function for the involved algorithms and its process can be described as followsLet the gray levels of a given image range over [0,L
-
1] and h
(
i
) denote the occurrence of gray-level i.Let(27)
N
=
∑
i
=
0
L
-
1
h
(
i
)
,
P
(
i
)
=
h
(
i
)
N
,
for
0
≤
i
≤
L
-
1
Maximize
f
(
t
)
=
w
0
w
1
(
u
0
-
u
1
)
2
,
where
(28)
w
0
=
∑
i
=
0
t
-
1
P
i
;
u
0
=
∑
i
=
0
t
=
1
i
×
P
i
w
0
(29)
w
1
=
∑
i
=
t
L
-
1
P
i
;
u
1
=
∑
i
=
t
L
=
1
i
×
P
i
w
1
and the optimal threshold is the gray level that maximizes (27). Then (28) can also be written as
(30)
f
(
t
)
=
δ
2
-
w
0
×
δ
0
2
-
w
1
×
δ
1
2
,
where w
0, w
1, u
0, and u
1 are the same as given in (9) and (10), and
(31)
δ
0
=
∑
i
=
0
t
-
1
(
i
-
u
0
)
2
×
P
i
w
0
,
δ
1
=
∑
i
=
t
L
-
1
(
i
-
u
1
)
2
×
P
i
w
1
,
δ
0
=
∑
i
=
0
L
-
1
(
i
-
u
)
2
×
P
i
w
,
w
=
∑
i
=
0
L
-
1
P
i
,
u
=
∑
i
=
0
L
=
1
i
×
P
i
w
.
Expanding this logic to multilevel threshold,
(32)
f
12
(
t
)
=
w
0
×
δ
0
2
+
w
1
×
δ
1
2
+
w
2
×
δ
2
2
+
⋯
+
w
1
×
δ
N
2
,
where N is the number of thresholds.Equation (32) is considered as the objective function for the proposed ARFO based procedure which is to be optimized (minimized). A close look into this equation will show it is very similar to the expression for uniformity measure.


## 4.2. Experiment Setup

As illustrated above, the performance evaluation based on entropy criterion is one of key indicators in image segmentation scenarios, our experimental studies are conducted on a suit of widely used image datasets [38–43], which contain avion.ppm, house.ppm, lena.ppm, peppers.ppm, safari04.ppm, and hunter.pgm (available at http://decsai.ugr.es/cvg/dbimagenes/). The size of each involved image is 512 * 512, with a unique grey-level histogram. The intrinsic multimodality of these tested images’ histograms will result in a more difficult optimization process. The segmentation results of the proposed algorithm and other algorithms are evaluated by the traditional Otsu method. The parameters of these algorithms (namely ARFO, ABC, PSO, GA, and DE) can be referred in Section 3.2, which have been proven as the optimal setting.The numbers of thresholdsM
-
1 investigated in the experiments were 2, 3, 4, 5, 7, and 9. The population size is 20. The maximum number of FEs is 2000, which has sufficient time to search on high dimension problems. Since the population-based algorithms have the randomized characteristics, all the experiments were repeated 30 times for each image for each M
-
1 value. In all experiments, to ensure that the initial values of each random algorithm are equal, we used the MATLAB command rand(
X
,
Y
)
*
(
max
⁡
L
-
min
⁡
L
)
+
min
⁡
L. The max
⁡
L and min
⁡
L are the maximum and minimum grays of the tested images. Figure 3 presents four original images and their histograms.Figure 3
Test images and their histograms.


## 4.3. Experimental Results of Multilevel Threshold

It is worthy noted the involved population-based algorithms for image segmentation, including ARFO, may only reduce the computation consumption but fail to significantly improve the accuracy of the optimal result. Therefore, we strived to utilize the proposed algorithm to obtain multiple thresholds with larger fitness values and fast computation ability. Due to the exhaustive search characteristic of classical Otsu method, its result can be regarded as the basis for comparison with other EA-based algorithms. We utilized (32) as the fitness function to evaluate involved algorithms.Case 1 (multilevel threshold results and efficiency and different methods with M − 1 = 2, 3, 4).
In tables reporting the fitness, mean computation time, and corresponding optimal thresholds (withM
-
1 = 2, 3, 4) obtained by Otsu are presented in Table 4. Due to the long consumption of CPU time when M
-
1
>
4, the correlative values are not worthy being listed in our experiment. However, other results attained by the population-based algorithms were compared. It should be noted that in the real-time image proceeding scenarios an appropriate balance of running time and high accuracy is demanded seriously; the termed CPU time of the population-based methods should be analyzed in detail. As shown in Table 5, there are no apparent differences about CPU times among the involved population-based methods, which are significantly superior in terms of time complexity for high-dimensional image segmentation problems.
From Table6, it can be observed that the termed fitness values result achieved by the proposed ARFO is generally closed to that by Otsu when M
-
1 = 2, 3, 4, whereas the performance regarding computation time of ARFO is significantly better than its counterpart Otsu. Furthermore, the ARFO-based method also performed most powerful among the population-based methods in most involved datasets. This can be explained by the fact that the ARFO has an appropriate balance between exploration and exploitation. Compared with other population-based algorithms, the ARFO is endowed with branching operations to vary population according to the nutrient concentration change to facilitate better exploitation search in complex higher dimensional search space as well as the branching operation will be used to fine tune the best solutions.Table 4
Objective values and thresholds by the Otsu method.
Image
M
-
1
=
2
M
-
1
=
3
M
-
1
=
4
Objective values
Optimal thresholds
Objective values
Optimal thresholds
Objective values
Optimal thresholds
Avion
4.11
E
+
04
114, 175
4.35
E
+
04
95, 151, 192
3.91
E
+
04
85, 131, 172, 208
House
3.04
E
+
04
109, 175
3.15
E
+
04
87, 136, 183
2.69
E
+
04
76, 121, 155, 191
Lena
1.77
E
+
04
95, 156
2.11
E
+
04
84, 123, 175
1.96
E
+
04
76, 121, 155, 191
Hunter
1.12
E
+
04
53, 118
1.13
E
+
04
37, 88, 144
1.07
E
+
04
33, 76, 121, 156
Mean CPU time
3.99
E
+
00
8.69
E
+
01
2.42
E
+
03Table 5
The mean CPU time of the compared population-based methods on Otsu algorithm.
Dim
Algorithm
ARFO
ABC
PSO
DE
GA
2 Dim
2.54
E
-
01
3.01
E
-
01
1.04
E
+
00
9.55
E
-
01
3.01
E
+
00
3 Dim
2.92
E
-
01
3.44
E
-
01
1.10
E
+
00
9.96
E
-
01
3.76
E
+
00
4 Dim
2.86
E
-
01
3.64
E
-
01
1.11
E
+
00
1.05
E
+
00
5.26
E
+
00
5 Dim
3.26
E
-
01
3.76
E
-
01
1.00
E
+
00
1.10
E
+
00
1.41
E
+
00
7 Dim
3.58
E
-
01
5.26
E
-
01
1.21
E
+
00
1.19
E
+
00
2.04
E
+
00
9 Dim
5.18
E
-
01
6.46
E
-
01
1.23
E
+
00
1.11
E
+
00
2.81
E
+
00Table 6
Objective value and standard deviation by the compared population-based methods on Otsu algorithm.
Image
M
-
1
Objective values (standard deviation)
ARFO
ABC
PSO
DE
GA
Avion
2
3.84
E
+
04
3.32
E
-
01
3.70
E
+
04
8.49
E
-
12
3.83
E
+
04
1.58
E
-
02
3.82
E
+
04
1.58
E
-
01
3.79
E
+
04
3.92
E
-
01
3
3.84
E
+
04
1.82
E
-
01
3.75
E
+
04
2.04
E
-
01
3.73
E
+
04
1.38
E
-
01
3.84
E
+
04
1.58
E
+
00
3.79
E
+
04
2.82
E
-
01
4
3.84
E
+
04
5.23
E
-
01
3.79
E
+
04
5.39
E
-
01
3.81
E
+
04
8.18
E
-
01
3.73
E
+
04
3.92
E
-
01
3.77
E
+
04
1.63
E
+
00
House
2
3.16
E
+
04
3.36
E
-
02
3.04
E
+
04
4.80
E
-
01
3.11
E
+
04
4.74
E
-
02
3.11
E
+
04
3.78
E
-
03
3.11
E
+
04
1.05
E
+
00
3
3.18
E
+
04
2.99
E
-
02
3.11
E
+
04
3.72
E
-
02
3.07
E
+
04
2.04
E
-
01
3.15
E
+
04
4.72
E
-
01
3.11
E
+
04
7.19
E
-
01
4
3.23
E
+
04
9.14
E
-
01
3.13
E
+
04
5.42
E
-
01
3.15
E
+
04
7.04
E
-
01
3.15
E
+
04
3.61
E
-
01
3.12
E
+
04
5.08
E
+
05
Lena
2
1.96
E
+
04
9.90
E
-
03
1.92
E
+
04
4.17
E
-
12
1.92
E
+
04
2.40
E
-
02
1.93
E
+
04
3.84
E
-
12
1.91
E
+
04
5.18
E
-
01
3
1.98
E
+
04
8.12
E
-
03
1.94
E
+
04
5.94
E
-
02
1.96
E
+
04
3.60
E
-
01
1.96
E
+
04
4.82
E
-
02
1.94
E
+
04
1.28
E
+
00
4
1.98
E
+
04
2.34
E
-
01
1.94
E
+
04
1.09
E
+
00
1.93
E
+
04
4.04
E
+
00
1.95
E
+
04
3.92
E
-
01
1.93
E
+
04
2.15
E
+
00
Peppers
2
1.93
E
+
04
3.78
E
-
01
1.84
E
+
04
4.25
E
-
12
1.89
E
+
04
3.58
E
-
01
1.87
E
+
04
4.71
E
-
12
1.86
E
+
04
1.71
E
-
01
3
1.94
E
+
04
9.03
E
-
01
1.87
E
+
04
1.50
E
-
02
1.90
E
+
04
1.84
E
-
01
1.89
E
+
04
5.07
E
-
01
1.89
E
+
04
3.82
E
-
01
4
1.95
E
+
04
8.79
E
-
01
1.88
E
+
04
1.50
E
-
01
1.93
E
+
04
1.51
E
+
00
1.93
E
+
04
7.17
E
-
12
1.91
E
+
04
4.27
E
+
00
Safari04
2
2.26
E
+
04
4.81
E
-
01
2.23
E
+
04
9.61
E
-
01
2.23
E
+
04
2.50
E
-
01
2.14
E
+
04
3.66
E
-
01
2.23
E
+
04
3.93
E
-
01
3
2.28
E
+
04
8.34
E
-
01
2.24
E
+
04
3.79
E
-
02
2.23
E
+
04
1.62
E
-
01
2.15
E
+
04
6.49
E
-
12
2.23
E
+
04
5.03
E
-
01
4
2.27
E
+
04
5.70
E
-
01
2.23
E
+
04
1.11
E
+
00
2.25
E
+
04
3.86
E
-
01
2.24
E
+
04
1.46
E
+
00
2.23
E
+
04
5.03
E
-
01
Hunter
2
1.02
E
+
04
1.25
E
-
01
1.01
E
+
04
3.56
E
-
02
1.00
E
+
04
4.74
E
-
01
1.01
E
+
04
2.22
E
-
04
1.01
E
+
04
6.26
E
-
01
3
1.05
E
+
04
1.34
E
+
00
1.02
E
+
04
7.27
E
-
01
1.02
E
+
04
7.07
E
-
01
1.02
E
+
04
2.57
E
+
00
1.03
E
+
04
1.62
E
+
00
4
1.05
E
+
04
9.16
E
-
01
1.03
E
+
04
3.60
E
+
00
1.03
E
+
04
1.99
E
+
00
1.04
E
+
04
3.63
E
-
01
1.04
E
+
04
3.82
E
+
00Case 2 (multilevel threshold results and efficiency and different methods withM
-
1 = 5, 7, 9).
In this section, these algorithms are conducted on image segmentation withM
-
1 = 5, 7, 9, which can be considered as higher dimensional problems. In Table 7 reporting the termed experimental results of average fitness value and standard deviation obtained by each population-based algorithm is illustrated, where the larger values and smaller standard deviations indicate the better achievement. From Table 7, it is clearly visible that there are statistically significant differences among experimental results achieved by these population-based approaches, in terms of both accuracy (fitness values) and stability (standard deviation). Depending on the population-varying model and the fast convergence ability, ARFO exhibits more efficient performances than other classical population-based algorithms, especially on high dimensional problems. The results as illustrated in Figure 4 prove that the ARFO-based algorithm is one of the suitable resolutions for multilevel image segmentation problems.Table 7
Objective value and standard deviation by the compared population-based methods on Otsu algorithm.
Image
M
-
1
Objective values (standard deviation)
ARFO
ABC
PSO
DE
GA
Avion
5
3.78
E
+
04
3.73
E
+
04
3.76
E
+
04
3.77
E
+
04
3.73
E
+
04
3.53
E
-
01
4.63
E
+
00
1.47
E
+
00
2.23
E
+
00
2.41
E
+
00
7
3.77
E
+
04
3.74
E
+
04
3.76
E
+
04
3.76
E
+
04
3.75
E
+
04
3.78
E
-
01
4.75
E
+
00
1.47
E
+
00
3.73
E
+
00
2.83
E
+
00
9
3.78
E
+
04
3.76
E
+
04
3.75
E
+
04
3.77
E
+
04
3.74
E
+
04
1.45
E
+
00
2.35
E
+
00
2.68
E
+
00
3.36
E
+
00
1.59
E
-
01
House
5
3.14
E
+
04
3.08
E
+
04
3.11
E
+
04
3.12
E
+
04
3.09
E
+
04
1.46
E
+
00
1.70
E
+
00
1.58
E
+
00
2.22
E
+
00
3.65
E
+
00
7
3.15
E
+
04
3.12
E
+
04
3.10
E
+
04
3.13
E
+
04
3.10
E
+
04
1.45
E
-
01
5.97
E
+
00
1.46
E
+
00
2.99
E
+
00
3.76
E
+
00
9
3.16
E
+
04
3.11
E
+
04
3.09
E
+
04
3.13
E
+
04
3.10
E
+
04
1.12
E
+
00
5.85
E
+
00
1.25
E
+
00
1.03
E
+
00
1.29
E
+
00
Lena
5
1.97
E
+
04
1.94
E
+
04
1.91
E
+
04
1.95
E
+
04
1.92
E
+
04
1.58
E
-
01
2.13
E
+
00
1.39
E
+
00
2.78
E
+
00
2.60
E
+
00
7
1.97
E
+
04
1.94
E
+
04
1.95
E
+
04
1.95
E
+
04
1.93
E
+
04
1.15
E
+
00
5.74
E
+
00
2.59
E
+
00
1.79
E
-
01
1.48
E
+
00
9
1.96
E
+
04
1.93
E
+
04
1.94
E
+
04
1.95
E
+
04
1.93
E
+
04
1.57
E
+
00
6.83
E
+
00
3.76
E
+
00
1.23
E
+
00
3.66
E
+
00
Peppers
5
1.93
E
+
04
1.89
E
+
04
1.89
E
+
04
1.92
E
+
04
1.88
E
+
04
5.06
E
-
01
1.57
E
+
00
1.79
E
+
00
2.19
E
+
00
4.31
E
+
00
7
1.93
E
+
04
1.88
E
+
04
1.89
E
+
04
1.92
E
+
04
1.88
E
+
04
1.39
E
+
00
2.67
E
+
00
4.01
E
+
00
5.07
E
-
01
2.83
E
+
00
9
1.92
E
+
04
1.89
E
+
04
1.91
E
+
04
1.92
E
+
04
1.90
E
+
04
1.20
E
+
00
2.69
E
+
00
2.01
E
+
00
5.07
E
-
01
4.39
E
+
00
Safari04
5
2.35
E
+
04
2.20
E
+
04
2.21
E
+
04
2.04
45
E
4
2.20
E
+
04
1.03
E
-
01
6.73
E
+
00
1.46
E
+
00
2.23
E
+
04
2.19
E
+
00
7
2.34
E
+
04
2.20
E
+
04
2.22
E
+
04
3.54
E
+
00
2.20
E
+
04
1.48
E
+
00
4.75
E
+
00
1.11
E
+
00
2.23
E
+
04
3.74
E
+
00
9
2.25
E
+
04
2.20
E
+
04
2.21
E
+
04
2.30
E
+
00
2.20
E
+
04
9.68
E
-
01
1.43
E
+
00
5.18
E
-
01
3.65
E
+
00
2.34
E
+
00
Hunter
5
1.03
E
+
04
1.02
E
+
04
1.02
E
+
04
1.03
E
+
04
1.02
E
+
04
7.91
E
-
01
4.32
E
+
00
3.65
E
+
00
4.95
E
-
01
3.61
E
+
00
7
1.05
E
+
04
1.02
E
+
04
1.02
E
+
04
1.03
E
+
04
1.03
E
+
04
1.83
E
+
00
3.98
E
+
00
2.35
E
+
00
1.76
E
+
00
2.66
E
+
00
9
1.04
E
+
04
1.03
E
+
04
1.03
E
+
04
1.03
E
+
04
1.03
E
+
04
1.37
E
+
00
5.85
E
+
00
2.45
E
+
00
4.94
E
-
01
3.98
E
+
00Figure 4
The comparison of the thresholded images by ARFO-based and ABC-based algorithms (M
-
1
=
9).


## 5. Conclusions

This paper proposes a new type of bionic algorithm, namely, artificial root foraging optimizer (ARFO) inspired by the plant root behaviors including the tropism and growth behaviors. Through these ecological behaviors, ARFO model can set up the dynamics mechanism of root growing rapidly towards the global optima, in which the auxin concentration is utilized to determine how to select new growing points and branching number of roots. Then, we introduced three root growing operators to perform the exploring and exploiting progress. Then, against a suit of benchmark functions, ARFO appeared to perform more powerful than other state-of-the-art algorithms in most tested cases. Finally, the ARFO algorithm is applied for the real-world image segmentation problems. The experimental results of ARFO on each dataset indicate its significant improvement compared to other population-based methods. Furthermore, as an effective image segmentation method, ARFO can also be incorporated into other popular threshold segmentation methods based on optimizing the fitness function.



---

*Source: 471209-2014-07-16.xml*