Solving forest harvesting problem with artificial intelligence approaches: An example in Taiwan

Abstract

In this paper, we study the forest harvesting problem (FHP). A forest is assumed to be divided into several identical square units, and each unit has its harvesting value based on its type. Harvesting a unit will affect the growth and values of its neighboring units. In this FHP, the best harvesting plan of a unit must be identified to maximize three various objectives simultaneously. The FHP is a multiobjective mathematical and an NP-hard problem. We apply three artificial intelligence algorithms, namely, immune algorithm, genetic algorithm, and particle swarm optimization, for maximizing the weighted objective to solve the FHP. We also solve the following two sets of test problems: (i) a set of randomly generated FHP problems and (ii) a practical problem in Taiwan. Numerical results show the performance of the three algorithms for the test problems. Finally, we compare and discuss the effects of various weights for the three objectives.

Keywords

Forest harvesting problem optimization immune algorithm genetic algorithm particle swarm optimization

1 Introduction

The environmental organization Greenpeace International reported that the Indonesian rainforest, which is the most damaged among the world’s three major rainforests, is gradually disappearing at a speed of nearly 300 football fields per hour. In the 10 years since 2000, the Amazon rainforest has been destroyed by an area of 240,000 square kilometers, which is equivalent to the size of British territory, and the global forest area has continued to disappear at a rate of 7.3 million hectares per year. According to the Forest Vitality Series report published by the World Wide Fund for Nature in 2015, if the deforestation of forest continues, then a total of 170 million hectares of forest will disappear within 15 years. Particularly, if no action is taken, then more than 230 million hectares of forest will vanish before 2050.

Forests have multiple functions. For example, they can purify the air, reduce noise, store water and soil to regulate climate, prevent wind, and protect farmland. In addition, the spatial structure of forests affects the abundance of species, ecological processes [13], and sedimentation status [8].

Several aspects, such as biodiversity, landscape aesthetics, wood quality, water pollution, habitats of species, and erosion, should be considered in forest management and planning [14]. As mentioned in Reference [6], the adjacency of units in a forest is a key problem in forest management. The adequate adjacency of units in a forest can improve the condition of soil erosion and sedimentation and provide a good habitat for most forest species. In 2012, Fotakis et al. [6] considered a bi-objective forest planning problem to maximize timber harvest volume and minimize sedimentation. They proposed a spatial NSGA method to solve the problem and demonstrate its performance. In 2010, Billionnet [2] studied a three-objective forest harvesting problem (FHP), in which a forest is divided into several identical square units. The FHP aims to determine harvested and non-harvested units in forests to maximize three objectives related to the expected population of three species. Particularly, the first objective corresponds to the expected population of the first species in the harvested units; the second objective relates to the expected population of the second species to its neighboring non-harvested units; and the third objective is the expected population of the third species, which is proportional to the total edges between harvested and non-harvested units. They used a linear programming approach of bipartite graph to solve the three-objective FHP. In 2015, Liu and Lin [14] applied a cultural algorithm (CA) to study an FHP by maximizing the total amount of harvested wood subject to some related harvest schedules. They showed that CA outperforms simulated annealing (SA) for the test problems. In addition, several researchers have investigated other variants of the FHP and proposed several approaches to optimize different objectives. We refer to [3–5 , 16–18] for additional information regarding the variants and their approaches for maintaining the brevity of this paper.

As previously discussed, forests have numerous functions, and several aspects are involved in the FHP. On the one hand, forests should be maintained to purify the air, store water and soil to regulate climate, prevent wind and reduce disasters, and preserve habitats for species. On the other hand, matured trees should be harvested for wood, which is one of the most important resources. Therefore, harvesting matured trees and replanting trees in the forest such that the impact of harvesting can be reduced and the forest can be sustainable is a tradeoff problem.

In this study, we investigate a three-objective FHP, which is similar to that of [2]. In the FHP, a forest is assumed to be divided into several identical square units, and each unit has its harvesting value (volume) depending on its type. Additionally, harvesting a unit will affect the growth and values of its neighboring units. In the considered FHP, the three various objectives corresponding to three expected populations of three species are weighted into an objective. Thus, the considered FHP should determine the harvested and non-harvested units in the forest to maximize the weighted objective.

Considering that FHP is an NP-hard problem, the solving time for the optimal solution will drastically increase if the problem size is large, that is, many identical square units exist in the problem. We apply three artificial intelligence algorithms, namely, immune algorithm (IA), particle swarm optimization (PSO), and genetic algorithm (GA), to solve the FHP. We also solve the following two sets of test problems: (i) a set of randomly generated FHP problems and (ii) a practical problem in Taiwan. The numerical results show the performance of the three algorithms for the test problems. Finally, we compare and discuss the effects of various weights for the three objectives in the weighted objective. This study has multiple purposes:

Application of IA to solve the FHP;

Comparison of the numerical results of IA with those of GA and PSO;

Analysis of the effects of weight for the three various objectives in the FHP.

2 FHP and objectives

The notations and assumptions of the FHP are as follows.

Notations.

the number of rows for the FHP, 1≤i≤n

the number of columns for the FHP, 1≤j≤m

the unit of row i and column j in the FHP, 1≤i≤n, 1≤j≤m

the expected number of sides (each with length 3 km between the harvested and non-harvested units x_ij = 0 if unit C_ij is harvested; otherwise x_ij = 1, 1≤i≤n, 1≤j≤m

the expected population of species one harvested for unit C_ij, 1≤i≤n, 1≤j≤m

the expected population of species two for unit C_ij, 1≤i≤n, 1≤j≤m

the probability that unit C_ij is connected to unit C_kl, 1≤i≤n, 1≤j≤m, 1≤k≤n, 1≤l≤m

the connectivity of unit C_ij with the other units left in old growth, 1≤i≤n, 1≤j≤m, that is, the probability that species two can visit other neighboring non-harvested units

the minimal number of non-harvested units required in the forest

Assumptions.

The forest is divided into m×n identical square units C_ij, 1≤i≤m, 1≤j≤n (Fig. 1).

Three species are considered in the units of the forest.

Each forest unit C_ij is represented by a Boolean variable x_ij. If x_ij = 0, then this unit is harvested with the expected population of species one t_ij; otherwise, x_ij = 1 (Fig. 1).

Three objectives, namely, S₁, S₂, and S₃, are considered in the FHP and are defined as follows.

S₁ is the amount of the expected population of species one for the harvested units and can be expressed as: $S_{1} = \sum_{i = 1}^{m} \sum_{j = 1}^{n} t_{ij} (1 - x_{ij})$ (1)

S₂ is the sum of the expected population of species two to its neighboring non-harvested units. We assume that species two can only visit its neighboring non-harvested units within two layers of units (Fig. 2(b)). S₂ can be expressed as:

$S_{2} = \sum_{i = 1}^{m} \sum_{j = 1}^{n} p_{ij} {PR}_{ij} x_{ij},$ (2) where ${PR}_{ij} = 1 - \prod_{k = 1}^{m} \prod_{l = 1}^{n} (1 -$ pr_ijklx_kl).

S₃ is the expected population of the third species, which is proportional to the total edges between non-harvested and harvested units. It can be expressed as follows:

S_{3} = gL,

(3)

where L is the expected length of sides (km) between harvested and harvested units, and g is a constant.

The three objectives of the FHP are evaluated by a weighted objective, that is, we aim to maximize S₁w₁+ S₂w₂+ S₃w₃, where w₁, w₂, and w₃ are given weights for the three objectives in (1)–(3).

Fig. 1

Example of the FHP with m×n units.

Fig. 2

Example of FHP with 5×5 units.

Example: We consider an FHP with 5×5 units and a feasible solution in Fig. 2(a) (shaded units are non-harvested), which is similar to that in [2]. Suppose that pr_ijkl= 0.5 for all i, j, k, and l, and g = 1.26157. Then, the three objectives in (1)–(3) are computed as follows:

S₁ = 82 + 86 + (90×2) + 92 + 75 + (100×3)=815.

Fig. 2(b) shows the computation of pr_44kl, and Fig. 2(c) shows that PR₄₄ = 1-(1-0.5)⁵(1-0.15)⁴ = 0.98. Thus, S₂= (5×0.92) + (5.6×0.98) + (4.5×0.99) +(3.5×0.97) + (3.5×0.94) + (4.5×0.99) + (5×1)+ (4.5×0.99) + (3.5×0.94) + (4.5×1) + (5×0.99) + (2×0.62) + (4.5×0.99) + (5×0.98) +(4.5×0.92) + (5×0.91) = 67.25.

S₃ = gL = 1.26157 × 72 = 90.83, where L = 3×24 (24 = the expected number of sides between the harvested and non-harvested units).

Thus, the weighted objective of this FHP is to maximize S₁w₁+ S₂w₂+ S₃w₃, where w₁, w₂, and w₃ are given weights of the three objectives.

3 Approaches

The considered FHP is an NP-hard problem. Examples of approaches to solving this problem are as follows: (i) exact methods, such as linear programming approach; and (ii) approximate methods, such as heuristics methods. The exact methods can obtain global optimal solutions for the FHP but are time consuming for a large-size FHP, that is, the FHP has many units. Thus, we attempt to solve the FHP by using artificial evolutionary algorithms.

3.1 Three algorithms

In the past years, artificial intelligence evolutionary algorithms have been widely used to solve different complicated optimization problems. For example, [11] surveyed more than 40 novel optimization algorithms based on nature or animal behavior. Among various artificial intelligence algorithms, GA and PSO may be the most popular methods used to solve different types of complicated optimization problems due to their numerous successes in the literature, as presented in [1, 15]. In addition, IA, an evolutionary algorithm based on the immune mechanism, has attracted considerable attention, such as that presented in [12]. The IA is similar to GA in crossover and mutation, except for its unique memory mechanism. In the memory mechanism, although the objectives of IA are excellent, this algorithm can delete similar solutions. Moreover, IA can keep additional antibody diversity during the evolutionary process and has a higher probability to obtain optimal solutions than that of GA. Therefore, we adopt IA in this study to solve the FHP. In addition, the numerical results of IA are compared with those of GA and PSO. The main steps of IA are shown in Fig. 3, and the main steps of GA and PSO are referred to in [1, 15].

Fig. 3

Main steps of IA.

3.2 Encoding scheme and example

Our encoding of solution for an FHP with size m×n is a binary sequence of length m×n, where 0 indicates a unit to be harvested and 1 otherwise. Considering the presence of a constraint for the minimal number of non-harvested units for the habitat requirement (that is, NH) in the FHP, a random binary solution may be infeasible. However, we may use the following simple process to overcome the infeasibility of solutions. Suppose that A₁ = (a₁, a₂, ⋯ , a_m×n) is a binary sequence, where a_i ∈ {0, 1} , 1 ≤ i ≤ m × n.

If $\sum_{i = 1}^{m \times n} a_{i} ⩾ NH$ , then A₁ is a feasible solution for the FHP and conversion is not needed.

If $\sum_{i = 1}^{m \times n} a_{i} < NH$ and $\sum_{i = 1}^{m \times n} a_{i} ⩽ mn - NH$ , then A₁ is an infeasible solution for the FHP. We may replace the infeasible solution A₁ with B₁ = (b₁, b₂, …, b_m×n), where a_k+ b_k= 1, 1≤k≤m×n. Then, $\sum_{i = 1}^{m \times n} b_{i} = mn - \sum_{i = 1}^{m \times n} a_{i} ⩾ mn - (mn - NH) = NH$ , which implies that B₁ is a feasible solution for the FHP.

If $\sum_{i = 1}^{m \times n} a_{i} < NH$ and $\sum_{i = 1}^{m \times n} a_{i} > mn - NH$ , then A₁ is an infeasible solution for the FHP. We may add a penalty to the weighted objective. Particularly, the new weighted objective for this case is equal to the weighted objective $- R \times (NH - \sum_{i = 1}^{m \times n} a_{i})$ , where R > 0 is a given penalty weight.

Next, we use the following example to illustrate the encoding scheme.

Example: Consider an FHP example with 5×5 units in Fig. 1. Suppose that NH = 18.

If A₁= (1,0, ... ,0,1,0, ... ,0,1,0) is a binary sequence of length 25, then (a₁+ a₂+ ... + a₂₅) = 3 < 18 (= NH) and (a₁+ a₂+ ... + a₂₅)≤25 –18. Thus, we use B₁ =(0,1, ... ,1,0,1,1, ... ,1,0,1) to replace A₁.

If A₁= (1,..,1,0) is a binary sequence of length 25, then (a₁+ a₂+ ... + a₂₅) = 16. The new weighted objective = the weighted objective –R×(18 –16) because (a1 + a2+ ... + a25) <18 (= NH) and (a1 + a2 + ... + a25) >25 –18(= mn –NH).

Following the above encoding scheme, we may compute the weighted objective for a binary solution by using S₁, S₂, and S₃ in (1)–(3).

4 Numerical results and discussions

4.1 Parameter setting and equipment environment

In this study, the three algorithms, namely, IA, GA, and PSO, were used to solve the test problems of the FHP. Given that the parameter setting of algorithms affects their performance, we adopted the Taguchi method to find the optimal combination of parameters from a candidate set of suggested values in the literature. Finally, the adopted parameters for IA and GA in the experiment are as follows: population = 300, maximal iterations = 10000, crossover rate = 0.6, and mutation rate = 0.01. Meanwhile, the adopted parameters for PSO are as follows: population = 500, maximal iterations = 10000, c₁ = 2, and c₂ = 1.5. In addition, we set the three algorithms to stop in advance if no improvement is observed for 1000 consecutive iterations.

We used MATLAB 2009 to solve the problem with Windows 7 Intel(R) Core(TM) i5-4570 CPU 3.2 GHz and 4 GB RAM computer. Each test problem was solved 30 times using the three artificial intelligence algorithms.

4.2 Test problems

There are two main test problems as follows.

(1) Test Problem 1:

This problem comprises three test problems with various sizes, namely, 10×10, 10×20, and 20×20, where the corresponding p_ij was randomly selected from {1, 1.5, 2, 2.5, 3, 3.5, 4, 4.5, 5} and t_ij was randomly generated from the interval of [10, 50]. In addition, g = 1.26157.

(2) Test Problem 2:

One practical problem in Taiwan exists with a size of 30×30 (Fig. 4(a)). The corresponding p_ij, t_ij, and g of this problem were generated and set as those in Test Problem 1.

Fig. 4

Numerical results of Test problem 2 with NH = 0. (black unit = harvested unit).

We set N H = 0 and the weights for the three objectives in (1)–(3) as follows to analyze the effects of the weights for the three objectives: (1,1,1), (10,1,1), (1,10,1), (1,1,10), (10,10,1), (10,1,10), and (1,10,10). Therefore, we solved seven instances for each test problem by using the three adopted algorithms.

4.3 Numerical results and discussions

We used IA, GA, and PSO to solve two sets of test problems. IA, GA, and PSO were executed 30 times for each test problem. Numerical results, including the number of iterations for convergence, the weighted objective, and CPU time, are summarized in Tables 1 4. The best forest harvesting plans for Test Problem 2 are depicted in Fig. 4.

Table 1
Numerical results of Test Problem 1 (10×10)

w (1,1,1) (10,1,1) (1,10,1) (1,1,10) (10,10,1) (10,1,10) (1,10,10)

Avg 3026.8 29544 4092.66 8687.58 29451.2 30019.37 10028.9

OBJ Max 3043.80 29750 4110.10 9206.40 29750 30256 10570

GA Std 17.06 301.55 16.38 216.60 245.83 224.70 243.49

CI-H 3060.24 30135.04 4124.76 9112.12 29933.03 30459.78 10506.14

CI-L 2993.36 28952.96 4060.56 8263.04 28969.37 29578.96 9551.66

CPU Avg 21.86 17.74 26.30 32.08 16.29 21.28 25.60

Std 4.32 0.18 2.54 8.17 0.20 4.17 5.41

ITER Avg 165.60 28.37 106.63 425.90 29.73 266.77 230.03

Std 234.55 2.95 107.50 365.15 3.40 247.39 249.71

Avg 3041.96 29750 4106.49 8797.58 29750 30233.63 10095.41

OBJ Max 3043.80 29750 4110.10 9206.40 29750 30256 10570

Std 1.43 0 4.58 189.85 0 18.36 167.53

CI-H 3044.76 29750.00 4115.47 9169.69 29750.00 30269.62 10423.77

CI-L 3039.16 29750.00 4097.51 8425.47 29750.00 30197.64 9767.05

IA CPU Avg 78.93 72.08 81.88 110.78 69.81 86.33 86.62

Std 17.91 2.24 11.81 37.64 3.51 23.73 20.44

ITER Avg 206.63 121.57 244.43 776.80 132.93 372.30 391.50

Std 279.14 27.40 180.65 570.69 45.83 365.88 329.33

Avg 2991 27735.70 4029.15 8094.33 27905.70 29378.73 9457.08

OBJ Max 3035.70 28734 4088.70 8531.80 28863 30201 9800.60

Std 33.43 593.61 40 251.83 497.44 425.84 191.32

CI-H 3056.52 28899.18 4107.55 8587.92 28880.68 30213.38 9832.07

CI-L 2925.48 26572.22 3950.75 7600.74 26930.72 28544.08 9082.09

PSO CPU Avg 50.61 44.79 64.69 59.88 43.20 46.41 60.76

Std 9.49 3.23 13.25 9.90 1.41 5 10.82

ITER Avg 146.40 47.20 191.13 143.50 42.07 98.03 156.03

Std 204.41 70.12 250.57 182.38 26.50 123.28 195.44

		w	(1,1,1)	(10,1,1)	(1,10,1)	(1,1,10)	(10,10,1)	(10,1,10)	(1,10,10)
		Avg	3026.8	29544	4092.66	8687.58	29451.2	30019.37	10028.9
	OBJ	Max	3043.80	29750	4110.10	9206.40	29750	30256	10570
GA		Std	17.06	301.55	16.38	216.60	245.83	224.70	243.49
		CI-H	3060.24	30135.04	4124.76	9112.12	29933.03	30459.78	10506.14
		CI-L	2993.36	28952.96	4060.56	8263.04	28969.37	29578.96	9551.66
	CPU	Avg	21.86	17.74	26.30	32.08	16.29	21.28	25.60
		Std	4.32	0.18	2.54	8.17	0.20	4.17	5.41
	ITER	Avg	165.60	28.37	106.63	425.90	29.73	266.77	230.03
		Std	234.55	2.95	107.50	365.15	3.40	247.39	249.71
		Avg	3041.96	29750	4106.49	8797.58	29750	30233.63	10095.41
	OBJ	Max	3043.80	29750	4110.10	9206.40	29750	30256	10570
		Std	1.43	0	4.58	189.85	0	18.36	167.53
		CI-H	3044.76	29750.00	4115.47	9169.69	29750.00	30269.62	10423.77
		CI-L	3039.16	29750.00	4097.51	8425.47	29750.00	30197.64	9767.05
IA	CPU	Avg	78.93	72.08	81.88	110.78	69.81	86.33	86.62
		Std	17.91	2.24	11.81	37.64	3.51	23.73	20.44
	ITER	Avg	206.63	121.57	244.43	776.80	132.93	372.30	391.50
		Std	279.14	27.40	180.65	570.69	45.83	365.88	329.33
		Avg	2991	27735.70	4029.15	8094.33	27905.70	29378.73	9457.08
	OBJ	Max	3035.70	28734	4088.70	8531.80	28863	30201	9800.60
		Std	33.43	593.61	40	251.83	497.44	425.84	191.32
		CI-H	3056.52	28899.18	4107.55	8587.92	28880.68	30213.38	9832.07
		CI-L	2925.48	26572.22	3950.75	7600.74	26930.72	28544.08	9082.09
PSO	CPU	Avg	50.61	44.79	64.69	59.88	43.20	46.41	60.76
		Std	9.49	3.23	13.25	9.90	1.41	5	10.82
	ITER	Avg	146.40	47.20	191.13	143.50	42.07	98.03	156.03
		Std	204.41	70.12	250.57	182.38	26.50	123.28	195.44

95 % Confidnace interval = (CI-L, CI-H).

Table 2

Numerical results of Test Problem 1 (20×20)

		w	(1,1,1)	(10,1,1)	(1,10,1)	(1,1,10)	(10,10,1)	(10,1,10)	(1,10,10)
		Avg	12531.5	120497	16583.2	33604.03	120999.67	124386	39180.73
	OBJ	Max	12617	121820	16679	34350	122150	125190	39709
GA		Std	52.51	855.29	51.68	500.69	970.86	488.20	333.52
		CI-H	12634.42	122173.37	16684.49	34585.38	122902.56	125342.87	39834.43
		CI-L	12428.58	118820.63	16481.91	32622.68	119096.78	123429.13	38527.03
	CPU	Avg	60.23	57.05	99.96	124.89	60.27	54.54	144.92
		Std	16.47	19.03	32.69	47.08	18.38	15.88	48.83
	ITER	Avg	892.07	965	1050	1749.20	1063.40	779.93	2125.93
		Std	504.50	647.78	659.64	1015.22	636.84	523.65	1027.45
		Avg	12783.03	126267.67	16867.17	34829.77	126284.67	127343.67	40118.07
	OBJ	Max	12789	126320	16877	35924	126320	127420	40939
		Std	4.02	71.13	7.64	437.09	59.93	28.71	434.77
		CI-H	12790.91	126407.08	16882.14	35686.47	126402.13	127399.94	40970.22
		CI-L	12775.15	126128.26	16852.20	33973.07	126167.21	127287.40	39265.92
IA	CPU	Avg	576.56	482.99	566.49	768.65	477.58	522.93	767.88
		Std	93.63	71.02	111.63	212.13	71.85	91.18	227.69
	ITER	Avg	3291.07	2633.33	3086.97	4727.90	2602.97	2970.97	4688.67
		Std	699.67	538.69	803.53	1691.46	533.43	696.21	1657.34
		Avg	11974.57	107156.97	15832.97	28964.57	109462.73	116443.67	33995.43
	OBJ	Max	12431	116050	16251	30125	116360	121640	35584
		Std	286.31	4854.95	251.71	650.83	4097.81	2480.05	756.84
		CI-H	12535.74	116672.67	16326.32	30240.20	117494.44	121304.57	35478.84
		CI-L	11413.40	97641.27	15339.62	27688.94	101431.02	111582.77	32512.02
PSO	CPU	Avg	120.99	111.19	174.14	150.41	107.21	107.01	157.43
		Std	36.69	10.36	34.09	12.44	15.60	15.80	32.60
	ITER	Avg	162.13	69.20	238.27	131.63	81	110.83	183.07
		Std	332.28	98.07	216.04	76.52	148.95	179.58	228.13

95 % Confidnace interval = (CI-L, CI-H).

Table 3

Numerical results of Test Problem 1 (10×20)

		w	(1,1,1)	(10,1,1)	(1,10,1)	(1,1,10)	(10,10,1)	(10,1,10)	(1,10,10)
		Avg	6114.94	59751.37	8188.53	17525.30	59852.97	60837.57	20143.67
	OBJ	Max	6150.60	60370	8228.80	18084	60370	61197	20897
GA		Std	40.93	449.20	35.74	350.27	390.35	292.65	290.37
		CI-H	6195.16	60631.80	8258.58	18211.83	60618.06	61411.16	20712.80
		CI-L	6034.72	58870.94	8118.48	16838.77	59087.88	60263.98	19574.54
	CPU	Avg	30.92	20.50	40.73	75.48	20.58	39.16	62.65
		Std	7.93	0.25	9.27	25.07	0.25	11.74	23.30
	ITER	Avg	481.30	80.20	373.90	1724.77	80.93	991.73	1229.33
		Std	366.30	11.96	315.59	896.97	8.72	600.56	830.07
		Avg	6146.57	60370	8223.82	17362.10	60370	61156.40	20017.50
	OBJ	Max	6150.60	60370	8230.60	18154	60370	61196	20782
		Std	2.67	0	5.27	288.01	0	27.20	392.70
		CI-H	6151.80	60370.00	8234.15	17926.60	60370.00	61209.71	20787.19
		CI-L	6141.34	60370.00	8213.49	16797.60	60370.00	61103.09	19247.81
IA	CPU	Avg	143.73	123.79	145.79	171.14	119.95	132.44	168.52
		Std	25.28	11.02	19.34	58.06	11.60	21.32	42.14
	ITER	Avg	689.20	453.60	649.73	1006.47	439.30	587.70	1013
		Std	300.43	129.22	216.71	678.57	138.61	253.68	505.27
		Avg	5892.37	53642.67	7900.46	15161.57	54177.93	57849.77	17870.80
	OBJ	Max	6072.90	57138	8101.70	15904	56838	60341	18729
		Std	118.06	1594.83	110.59	422.28	1403.53	1114.59	447.86
		CI-H	6123.77	56768.54	8117.22	15989.24	56928.85	60034.37	18748.61
		CI-L	5660.97	50516.80	7683.70	14333.90	51427.01	55665.17	16992.99
PSO	CPU	Avg	56.51	53.46	86.34	79.20	53.01	64.47	80.20
		Std	5.62	1.39	14.48	10.65	1.49	15.68	9.91
	ITER	Avg	112.70	39.40	182.97	141.30	42.03	205.67	134.30
		Std	107.83	11.85	206.99	135.67	15.96	292.55	143.02

95 % Confidnace interval = (CI-L, CI-H).

Table 4

Numerical results of Test Problem 2 (30×30) with NH = 0

		w	(1,1,1)	(10,1,1)	(1,10,1)	(1,1,10)	(10,10,1)	(10,1,10)	(1,10,10)
		Avg	21187.73	187581	33723.80	56468.03	192926	202082	70698.23
	OBJ	Max	21221	189920	33843	57175	193490	203380	71738
GA		Std	93.46	1310.07	67.75	407.60	1342.11	779.87	427.19
		CI-H	21370.91	190148.74	33856.59	57266.93	195556.54	203610.55	71535.52
		CI-L	21004.55	185013.26	33591.01	55669.13	190295.46	200553.45	69860.94
	CPU	Avg	222.75	232.31	280.12	311.67	219.14	218.34	363.42
		Std	49.61	74.68	63.31	97.34	70.04	68.80	115.35
	ITER	Avg	843.03	992.30	898.97	1231.53	876.33	822.10	1531.67
		Std	414.57	646.63	432.64	683.03	603	588.18	787.87
		Avg	22365.87	206809.33	34913.17	63594.97	212191.33	217021.67	78048.83
	OBJ	Max	22399	207650	34941	64751	213050	217230	79615
		Std	14.54	425.08	17.68	459.36	479.55	109.99	752.48
		CI-H	22394.37	207642.49	34947.82	64495.32	213131.25	217237.25	79523.69
		CI-L	22337.37	205976.17	34878.52	62694.62	211251.41	216806.09	76573.97
IA	CPU	Avg	3718.54	3485.13	3729.15	3757.55	3532.15	3685.68	3631.95
		Std	219.10	442.73	268.56	133.93	379.93	200.25	279.67
	ITER	Avg	9465.93	8588.87	9172.23	9400.50	8816.30	9454.93	9595.43
		Std	848.71	1452.92	1018.28	666.93	1332.73	818.99	928.41
		Avg	21045.37	182618.33	32629.07	52388.07	190232.33	204477	66457.90
	OBJ	Max	22064	204880	33379	54127	209510	214010	68466
		Std	625.40	11726.64	413.99	960.84	9909.98	6483.34	1110.54
		CI-H	22271.15	205602.54	33440.49	54271.32	209655.89	217184.35	68634.56
		CI-L	19819.59	159634.12	31817.65	50504.82	170808.77	191769.65	64281.24
PSO	CPU	Avg	412.85	386.35	479.46	405.29	412.50	392.66	476.11
		Std	137.04	104.27	103.93	46.46	156.55	90.74	112.99
	ITER	Avg	467.27	333.20	277.17	138.40	467.57	370.57	305.80
		Std	505.41	388.04	258.39	123.03	578.94	329.51	306.06

95 % Confidnace interval = (CI-L, CI-H).

The following findings are obtained from Tables 1 4 and Fig. 4.

As the problem size (units) increased, the CPU time, the convergence iteration, and the standard deviation of the weighted objective also increased. For example, in Table 1, Test Problem 1 with weight (10,10,1) indicates that IA reported the following: CPU time = 69.81, convergence iteration = 132.93, and standard deviation of weighted objective = 0 for the instance of 10×10. Meanwhile, in Table 2, IA reported the following: CPU time = 477.58, convergence iteration = 2602.97, and standard deviation of weighted objective = 59.93 for the instance of 20×20. Similar results for GA and PSO are also obtained.

Among 28 test instances, IA, GA, and PSO achieved the best solutions of the three algorithms for 26, 12, and 0 times, respectively. For example, in Table 3, the best weighted objectives by IA, GA, and PSO are 18154, 18084, and 15904, respectively, for the 10×20 test problem with w = (1,1,10). Similar results for the other test problems are also obtained.

Among the three algorithms, the standard deviations of the objective function value for IA, GA, and PSO in Table 1 with w = (10,1,1) are 0, 301.55, and 593.61, respectively. Similar results for the other test problems are also obtained. Therefore, IA is more reliable than GA, and GA is more reliable than PSO.

Among the three algorithms, GA is faster than PSO, and PSO is faster than IA. For example, in Table 2, the average CPU times for IA, GA, and PSO are 482.99, 57.05, and 111.19, respectively, for the 20×20 test problem with w = (10,1,1). Similar results for the other test problems are also obtained. The main reason is that IA had to update the memory list for the diversity of solutions in each iteration of the algorithm.

Among the three algorithms, Table 3 shows that the average numbers of iterations for IA, GA, and PSO are 649.73, 373.9, and 182.97, respectively, for the 10×20 test problem with w = (1,10,1). Similar results for the other test problems are also obtained. Thus, IA requires more convergence iterations than GA and PSO. However, IA could continuously improve.

If the first objective function S₁ is emphasized, then additional units should be harvested to increase the expected population of species one. For example, Fig. 4(c) with w = (10,1,1)] shows more units to be harvested than Fig. 4(b) with w = (1,1,1).

If the objective function S₂ is emphasized, then fewer units are harvested. For example, Fig. 4(d) with w = (1,10,1) shows less units to be harvested than Fig. 4(b) with w = (1,1,1).

If the objective function S₃ is emphasized, then additional isolated units are harvested to increase the total length of sides between the harvested and non-harvested units. For example, Fig. 4(e) with w = (1,1,10) shows more isolated units to be harvested than Fig. 4(b) with w = (1, 1, 1).

To further evaluate the effect of NH, i.e., the minimal number of non-harvested units required in the forest, we also set NH = 0.7×mn for Test Problem 2 and experiments 30 trials for weights (1,1,1), (10,1,1), (1,10,1), (1,1,10). Compared with Fig. 4 (NH = 0), Fig. 5 illustrates that more units are non-harvested and reserved. The effects of various weights in Fig. 5 are similar to those in Fig. 4.

Fig. 5

Numerical results of Test problem 2 with NH = 0.7×mn. (black unit = harvested unit).

We used the following function to evaluate the corresponding values among the three algorithms to simplify the comparison of the three algorithms for the best solution, the CPU time, and the number of convergence iterations. $g (A) = V (A) / \max {V (A), V (B), V (C)}$ (4)

where V(i) denotes the corresponding value by using Algorithm i, i∈{IA, GA, PSO}. The corresponding g(A) values of the three algorithms for the two sets of test problems with various weights are depicted in Figs. 6 and 7, and the ratios of best solution for various weights in g(A) are shown in Table 5. The following results are observed in Figs. 6 and 7.

Fig. 6

Numerical results of Test Problem 1 ( = PSO, = GA, = IA).

Fig. 7

Numerical results of Test Problem 2 with NH = 0. ( = PSO, = GA, = IA).

Table 5

Ratios of best solution for various weights

Prob.	Method	(1,1,1)	(10,1,1)	(1,10,1)	(1,1,10)	(10,10,1)	(10,1,10)	(1,10,10)
10×10	PSO	0.9973	0.9658	0.9948	0.9267	0.9702	0.9982	0.9272
10×10	GA	1.0000	1.0000	1.0000	1.0000	1.0000	1.0000	1.0000
10×10	IA	1.0000	1.0000	1.0000	1.0000	1.0000	1.0000	1.0000
10×20	PSO	0.9720	0.9187	0.9629	0.8386	0.9212	0.9546	0.8692
10×20	GA	0.9866	0.9644	0.9883	0.9562	0.9670	0.9825	0.9700
10×20	IA	1.0000	1.0000	1.0000	1.0000	1.0000	1.0000	1.0000
30×30	PSO	0.9874	0.9465	0.9843	0.8761	0.9415	0.9860	0.8963
30×30	GA	1.0000	1.0000	0.9998	0.9961	1.0000	1.0000	1.0000
30×30	IA	1.0000	1.0000	1.0000	1.0000	1.0000	1.0000	0.9945

For Test Problem 1, Fig. 6(a) shows that the best solutions by IA are superior to those by GA and PSO, and their values of the best solution in (4) are g(IA) = 0.99974, g(GA) = 0.99099, and g(PSO) = 0.944550. In addition, for Test Problem 2, Fig. 7(a) also shows that the best solutions by IA are superior to those by GA and PSO, and their values of the best solution in (4) are as follows: g(IA) = 1.00000, g(GA) = 0.92273, and g(PSO) = 0.94164. Hence, IA outperformed GA and PSO.

For Test Problem 1, Fig. 6(b) shows that the average solutions of 30 trials by IA are better than those by GA and PSO, and their values of average solution in (4) are as follows: g(IA) = 0.99926, g(GA) = 0.98615, and g(PSO) = 0.91672. In addition, for Test Problem 2, Fig. 7(b) indicates that the average solutions by IA are better than those by GA and PSO, and their values of average solution in (4) are as follows: g(IA) = 1.00000, g(GA) = 0.92206, and g(PSO) = 0.89608. Thus, IA outperformed GA and PSO.

For Test Problem 1, Fig. 6(c) shows that the average CPU time of 30 trials by IA is larger than that by GA and PSO, and their values of average CPU time in (4) are as follows: g(IA) = 1.00000, g(GA) = 0.23001, and g(PSO) = 0.44348. In addition, for Test Problem 2, Fig. 7(c) indicates that IA is slower than GA and PSO, and their values of average CPU time in (4) are g(IA) = 1.00000, g(GA) = 0.07228, and g(PSO) = 0.11610. Therefore, IA consumed a considerable amount of time to solve the FHP because its memory mechanism increased the CPU time.

For Test Problem 1, Fig. 6(d) shows that the average number of convergence iterations of IA is larger than that of GA and PSO, and their values of average iteration in (4) are as follows: g(IA) = 0.95239, g(GA) = 0.50733, and g(PSO) = 0.20743. In addition, for Test Problem 2, Fig. 6(d) indicates that IA requires more convergence iterations than GA and PSO, and their values of average iteration in (4) are g(IA) = 1.00000, g(GA) = 0.11137, and g(PSO) = 0.03674. Hence, GA and PSO converged with less iteration than IA, but their solutions are worse than those by IA. Meanwhile, IA requires more convergence iterations than GA and PSO, but IA converged to a better solution compared with that of GA and PSO due to its unique memory mechanism of solution.

We used the following statistical hypothesis based on 30 trials to analyze the performance of IA, GA, and PSO algorithms for each test instance further:

H₀: Average(A) = Average(B),

H₁: Average(A) ≠ Average(B),

where Average(A) and Average(B) denote the average of objectives by using algorithms A and B, respectively, and A, B∈{IA, GA, PSO}, A≠B. The corresponding p-values of the statistical hypothesis are summarized in Table 6, which also show the following:

Table 6

P-values for statistics hypothesis test

w	IA v.s. GA				IA v.s. PSO				GA v.s. PSO
		P 1		P2		P 1		P2		P 1		P2
	10×10	20×20	10×20	30×30	10×10	20×20	10×20	30×30	10×10	20×20	10×20	30×30
(1,1,1)	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	(0.227)
(10,1,1)	0.001	0.000	0.000	0.001	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.028
(1,10,1)	0.001	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
(1,1,10)	0.041	0.000	(0.053)	0.041	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
(10,10,1)	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	(0.151)
(10,1,10)	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	(0.054)
(1,10,10)	0.223	0.000	(0.162)	0.223	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000

P-values in brackets indicate that the average of the weighted objective of the latter method is better than that of the former.

For Test Problem 1, (i) IA and GA outperform PSO for all test instances; (ii) except for Instance 10×10 with weight (1,10,10), and Instance 10×20 with weights (1,1,10) and (1,10,10), IA outperforms GA for other test instances.

For Test Problem 2, (i) IA outperforms PSO for all test instances; (ii) except for instance with weight (1,10,10), IA outperforms GA; (iii) except for instances with weights (1,1,1), (10,10,1), and (10,1,10), GA out performs PSO for other test instances.

By testing the statistical hypothesis, for all 28 test instances, IA outperformed GA for 24 instances, IA outperformed PSO for all 28 instances, and GA outperformed PSO for 25 instances.

5 Conclusions

In this study, we used IA, GA, and PSO to investigate the FHP with three various objectives. We solved the two main test problems to evaluate the performance of the three algorithms and the effects of weight on each objective. The main results and findings of this research are as follows.

With the increase in problem size (additional units in the forest), the CPU time, the convergence iteration, and the standard deviation of weighted objectives also increase.

For most FHP test problems, IA outperforms GA, and GA outperforms PSO.

IA is more robust than GA, and GA is more robust than PSO.

GA is faster than PSO, and PSO is faster than IA.

Different weights of objectives lead to different harvest planning. The results of this research can provide useful information for harvest planning.

In the future,

More practical constraints can be considered in the FHP. For example, there are several types of spe-cies or woods in each unit of forest. Different types of trees would have various mature periods for harvesting.

More algorithms can be used to test the FHP and compare their effectiveness and efficiency.

Footnotes

Acknowledgments

This research was supported in part by grants (No. MOST 106-2221-E-150-039 and MOST 107-2221-E-150-02) from Ministry of Science and Technology, Taiwan.

References

Aguilar-Rivera

, Valenzuela-Rendón

and Rodríguez-Ortiz

J.J.

, Genetic algorithms and Darwinian approaches in financial applications: a survey, Expert Systems with Applications 42 (2015), 7684–7697.

Billionnet

, Solving a cut problem in bipartite graphs by linear programming: application to a forest management problem, Applied Mathematical Modelling 34 (2010), 1042–1050.

Borges

J.G.

and Hoganson

H.M.

, Assessing the impact of management unit design and adjacency constraints on forest wide spatial conditions and timber revenues, Canadian Journal of Forest Research 29 (1999), 1764–1774.

Ducheyne

E.I.

, De Wulf

R.R.

and De Baets

, Single versus multiple objective genetic algorithms for solving the even-flow forest management problem, Forest Ecology and Management 201 (2004), 259–273.

Ducheyne

E.I.

, De Wulf

R.R.

and De Baets

, A spatial approach to forest management optimization: linking GIS and multiple objective genetic algorithms, International Journal of Geographical Information Science 20 (2006), 917–928.

Fotakis

D.G.

, Sidiropoulos

, Myronidis

and Ioannou

, Spatial genetic algorithm for multi-objective forest planning, Forest Policy and Economics 21 (2012), 12–19.

Heinonen

and Pukkala

, The use of cellular automaton approach in forest planning, Canadian Journal of Forest Resources 37 (2007), 2188–2200.

Hof

and Bevers

, Optimal timber harvest scheduling with spatially defined sediment objectives, Canadian Journal of Forest Research 30 (2000), 1494–1500.

Hof

J.G.

and Joyce

L.A.

, A mixed integer linear programming approach for spatially optimizing wildlife and timber in managed forest ecosystems, Forest Science 39 (1993), 816–834.

10.

Hoganson

H.M.

and Rose

D.W.

, A simulation approach for optimal timber management scheduling, Forest Science 44 (1984), 526–538.

11.

Jain

, Singh

and Rani

, A novel nature-inspired algorithm for optimization: Squirrel search algorithm, Swarm and Evolutionary Computation 44 (2019), 148–175.

12.

Jiang

, Hao

, Cai

and Ding

, An improved reinforcement-immune algorithm for agricultural resource allocation optimization, Journal of Computational Science 27 (2018), 320–328.

13.

Kurttil

, The spatial structure of forests in the optimization calculations of forest planning-a landscape ecological perspective, Forest Ecology and Management 142 (2001), 129–142.

14.

Liu

W.Y.

and Lin

C.C.

, Spatial forest resource planning using a cultural algorithm with problem-specific information, Environmental Modelling & Software 71 (2015), 126–137.

15.

Nouaouria

, Boukadoum

and Proulx

, Particle swarm classification: a survey and positioning, Pattern Recognition 46 (2013), 2028–2044.

16.

Pukkala

and Kangas

, A heuristic optimization method for forest planning and decision making, Scandinavian Journal of Forest Research 8 (1993), 560–570.

17.

Strange

, Meilby

and Thorsen

J.T.

, Optimizing land use in afforestation areas using evolutionary self-organization, Forest Science 48 (2002), 543–555.

18.

Zhu

and Bettinger

, Estimating the effects of adjacency and green-up constraints on landowners of different sizes and spatial arrangements located in the southeastern U.S, Forest Policy and Economics 10 (2008), 295–302.