An upper bound of the mutation probability in the genetic algorithm for general 0-1 knapsack problem

Abstract

As an important part of genetic algorithms (GAs), mutation operators is widely used in evolutionary algorithms to solve $\mathcal{NP}$-hard problems because it can increase the population diversity of individual. Due to limitations in mathematical tools, the mutation probability of the mutation operator is primarily empirically set in practical applications.

In this paper, we propose a novel reduction method for the 0-1 knapsack problem(0-1 KP) and an improved mutation operator (IMO) based on the assumption $\mathcal{NP}\neq\mathcal{P}$, along with the utilization of linear relaxation techniques and a recent result by Dey et al. (Math. Prog., pp 569-587, 2022). We employ this method to calculate an upper bound of the mutation probability in general instances of the 0-1 KP, and construct an instance where the mutation probability does not tend towards 0 as the problem size increases. Finally, we prove that the probability of the IMO hitting the optimal solution within only a single iteration in large-scale instances is superior to that of the traditional mutation operator.

Keywords: genetic algorithm; mutation operator; mutation probability; upper bound; reduction method

1 Introduction

Definition 1.

[1, 2, 3, 4] Given a set $N$ consisting of $n$ items, where the $j$-th item has a profit and a weight, and denoted by $p_{j}$ and $w_{j}$, respectively. The goal of the 0-1 knapsack problem(0-1 KP) is to find an optimal subset of the set $N$ such that its total profit is maximized without exceeding the knapsack capacity $C$. The problem can be formulated as follows.

$$\max f(\textbf{X})=\sum\limits_{j=1}^{n}x_{j}p_{j},$$ (1.1)

s.t.

$$g(\textbf{X})=\sum\limits_{j=1}^{n}x_{j}w_{j}\leq C,$$ (1.2)

$$x_{j}\in\{0,1\},$$ (1.3)

where $x_{j}=1$ represents the $j$-th item is packed in the knapsack while $x_{j}=0$ indicates not. The objective function (1.1) of the model aims to maximize the total profit of the selected items. Inequality (1.2) is a capacity constraint that ensures the total weight of selected items do not exceed the knapsack capacity $C$. The decision variables $x_{j}$ are binary, as indicated in (1.3).

A substantial number of commercial decision-making and industrial manufacturing problems can be formulated as the 0-1 KP, such as budget control and cargo loading. This has resulted in a significant focus on developing efficient algorithms to tackle these problems. However, the finding for a polynomial-time algorithm to exactly solve the 0-1 KP, which is widely recognized as an $\mathcal{NP}$-hard problem[5, 6, 7], poses a substantial challenge. Existing algorithms like dynamic programming [8, 9] can only solve the problem exactly within a pseudo-polynomial time complexity of $O(nC)$, and even the fastest exact algorithms require $O(2^{n/2})$ [10, 11]. Additionally, a significant body of research indicates that unless one intends to demonstrate $\mathcal{NP}=\mathcal{P}$, we can assume $\mathcal{NP}\neq\mathcal{P}$ [12, 13]. If $\mathcal{NP}\neq\mathcal{P}$, for a general instance of the 0-1 KP, there are items whose selection in the optimal solution cannot be determined without exhaustively enumerating all feasible solutions.

Coincidentally, in the 1960s, Holland proposed a genetic algorithm(GA) [14] that simulates biological evolution, consisting of selection, crossover, and mutation operators, primarily used to search the binary solution space. Genetic algorithms(GAs) do not require the solution space to be continuous or need to calculate the gradient of the objective function during iteration, and as a result, the algorithm has been widely applied to solving various $\mathcal{NP}$-hard problems, such as the task scheduling problem [15], the Traveling Salesman Problem [16], and the Vehicle Routing Problems [17].

It is worth noting that the convergence analysis of GAs is primarily limited to convex problems or $\mathcal{P}$ problems [18, 19, 20, 21, 22, 23]. However, there is a lack of theoretical results for $\mathcal{NP}$-hard problems, especially for the optimal parameter control of operators such as the mutation operator (MO)[24, 25, 26]. Despite the widespread and efficient performance of GAs in practical industrial and application processes[27], there is still widespread skepticism among researchers regarding GA research[28, 29, 30, Seretal2019]. Moreover, the No Free Lunch (NFL) theorem states that there is no superior parameterization for every problem[27, 31, 32, 33, 34], but it is still important to determine if there exists an optimal parameterization for a given 0-1 KP instance.

1.1 Our results

In this paper, we build upon the assumption of $\mathcal{NP}\neq\mathcal{P}$ [5, 6, 7], improve the recent theoretical results of the branch and bound algorithm (B&B) [35, 36], integrating generating function approaches[37, 38], introduce a new reduction method and an improved mutation operator (IMO). Simultaneously, the method is applicable in computing an upper bound of the mutation probability for the 0-1 KP, as well as in constructing a counterexample where the mutation probability does not tend toward 0 with an increasing number of decision variables [25]. The contributions of this paper can be summarized as follows.

1. We propose an improved reduction method for the 0-1 KP, and provide the upper bound of leaves in the search tree (Theorem 3.6).

2. We propose the IMO and provide a general mathematical formulation for the upper bound of the mutation probability in solving the 0-1 KP (Theorem 3.7).

3. We construct examples where the mutation probability of the IMO does not tend towards 0 with the increase in problem size (Theorem 3.9).

4. We demonstrate that for large-scale instances, the performance of the IMO is superior to that of the MO (Theorem 3.11).

1.2 Organization of remainder paper

The remainder of this paper is organized as follows. In Section 2, we introduce the preliminary and relevant concepts of GAs. Section 3 presents the theoretical results, including a novel reduction method, the IMO proposed from this reduction method, and the associated theoretical findings. Additionally, we provide a convergence analysis of both the IMO and the MO. Finally, in the last section, we summarize the work of this paper and discuss future research directions.

2 Preliminary and existing conclusions

2.1 Preliminary

From the objective function (1.1) and inequality (1.2), a primitive greedy algorithm can be easily identified: the likelihood of an item appearing in the optimal solution increases as its profit increases and its weight decreases. Consequently, we introduce the classic concept of profit density[1, 2] to characterize the return on profit per unit weight of the items.

Definition 2.

For $j\in N$, the profit density of item $j$ is defined by

$$e_{j}=p_{j}/w_{j}.$$ (2.4)

Without loss of generality, we always assume that the items are arranged in non-increasing order, denoted as $e_{1}\geq e_{2}\geq\dots\geq e_{n}$. The greedy algorithm selects items according to their profit density in descending order until the first item that cannot be packed in the knapsack. The item is known as the break item $b$ and the solution obtained employing this greedy algorithm is the break solution. The characteristics of the break item and the break solution are defined as follows.

Definition 3.

The break item $b$ is defined to be the first item that cannot be packed in the knapsack with profit density descending order, that is,

$$\sum\limits_{j=1}^{b-1}w_{j}\leq C\quad\text{and}\quad\sum\limits_{j=1}^{b}w_{j}>C$$ (2.5)

Correspondingly, the residual capacity $r$ is defined as

$$r=C-\sum\limits_{j=1}^{b-1}w_{j}$$ (2.6)

Definition 4.

The break solution is donated as $\textbf{X}^{*}=(x^{*}_{1},x^{*}_{2},\dots,x^{*}_{n})\in\{0,1\}^{n}$, where $x^{*}_{j}=1$ if $j\in\{1,2,\dots,b-1\}$ and $x^{*}_{j}=0$ otherwise. For the sake of comparison, let Y be the optimal solution.

To further verify that items with higher profit density are more likely to be selected in the optimal solution, Pisinger[39] conducted a comparison between the optimal solution Y and the break solution $\textbf{X}^{*}$ for 1000 instances. Each instance had a data size of 1000, with the break item being the 500th. Pisinger found that the primary concentrations of differences between Y and $\textbf{X}^{*}$ were around the break item.

Based on empirical evidence derived from the greedy algorithm and a significant number of instances, it can be empirically observed that items with higher profit density tend to be more frequently selected for inclusion in the optimal solution Y. However, it is regrettable that this observation may be challenged by a classic counterexample, which can be shown as follows.

Example 2.1.

[2] Let $n=2$ and $C=M$. The first item is given by $w_{1}=1$ and $p_{1}=2$ whereas the second item is defined by $w_{2}=p_{2}=M$, where $M$ is a large constant.

Clearly, the break solution is $\textbf{X}^{*}=(1,0)$, but the optimal solution for the counterexample is $\textbf{Y}=(0,1)$. From this counterexample, we can observe that the profit density of an item does not necessarily determine whether it will be selected in the optimal solution, and also implies that in the process of solving through GA, the optimal probability distribution for the MO cannot be determined for a given instance.

Although it is currently impossible to determine all decision variables in polynomial time complexity, with the continuous development of solving algorithms, a reduction algorithm has been discovered that can exactly solve a part of decision variables in polynomial time complexity. Given an instance $\mathscr{P}$, the reduction algorithm treats the solving process as a search tree, where each branch represents a decision variable in the problem. We employ the notation $(\mathscr{P}|x_{j}=\alpha)$ to denote the subproblem where the $j$-th decision variable is fixed at $\alpha\in\{0,1\}$. To facilitate subsequent discussions, $\overline{v}(\mathscr{P}|x_{j}=\alpha)$ is introduced as an upper bound for this subproblem. Furthermore, we introduce $\underline{v}(\mathscr{P})$ as a lower bound for the original problem $\mathscr{P}$, typically computed using the greedy algorithm or a meta-heuristic algorithm.

Despite the subproblem $(\mathscr{P}|x_{j}=\alpha)$ is also known to be $\mathcal{NP}$-hard, relaxing the solution space to the real field converts it into a convex problem[40]. Consequently, the optimal solution for the relaxed subproblem, $\overline{v}(\mathscr{P}|x_{j}=\alpha)$, can be obtained within polynomial time complexity. These findings lead to several significant conclusions.

Theorem 2.2.

For the 0-1 KP, if $\overline{v}(\mathscr{P}|x_{j}=\alpha)\leq\underline{v}(\mathscr{P})$, then $x_{j}$ is set to $1-\alpha$ in the optimal solution, where $\alpha\in\{0,1\}$.

The decision variable $x_{j}$ is fixed, and the branch $x_{j}=\alpha$ is pruned. To compute $\overline{v}(\mathscr{P}|x_{j}=\alpha)$, the most common approach is to relax the solution space from the integer field to the real field. By doing this, we can obtain the well-known Dantzig upper bound[41], which is denoted as $U$ and is given by

$$U=\sum\limits_{j=1}^{b-1}p_{j}+\frac{rp_{b}}{w},$$

where $r$ represents the residual capacity. It is clear that the decision variables dominated by the Dantzig upper bound $U$ can be expressed as follows:

Theorem 2.3.

[39] For the $j$-th item, if $e_{j}\geq e_{b}$(resp. $e_{j}<e_{b}$) and $p_{j}>(w_{j}+r)p_{b}/w_{b}$(resp. $p_{j}<(w_{j}-r)p_{b}/w_{b}$), then we have $x_{j}=1$(resp. $x_{j}=0$) in the optimal solution.

The time complexity of computing the Dantzig bound is $O(n)$. It is evident that we can solve a decision variable exactly within polynomial time complexity. A significant amount of research has been conducted to improve the upper bounds obtained within the time complexity of $O(n^{2})$ [39, 42, 43, 44, 45], thereby enabling the exact fixing of more decision variables within polynomial time complexity. However, in the worst-case scenario, the method cannot solve all decision variables. Consequently, for a preprocessed $\mathcal{NP}$-hard problem, to the best of our knowledge, there is currently no relevant research available regarding which decision variables are more likely to be selected in the optimal solution. Moreover, if an $\mathcal{NP}$-hard problem can be exactly solved within polynomial time complexity, heuristic algorithms, such as GAs, that can only provide locally optimal solutions within polynomial time complexity lose their significance. Based on the current state of research, we make the following assumptions.

Assuption 1.

$\mathcal{NP}\neq\mathcal{P}$.

If Assumption 1 holds, then for an unpruned feasible solution, we cannot determine its status as an exact solution within polynomial time complexity. Consequently, each unpruned feasible solution has an equal probability of being the optimal solution.

2.2 Genetic algorithm Framework

With the prosperity of commercial activities, a vast number of business operations can be described as $\mathcal{NP}$-hard problems. Considering Assumption 1, the development of exact algorithms for solving such problems has been slow. Meanwhile, heuristic algorithms have flourished [50, 51]. As a classic algorithm in the field of heuristics, GAs have been extensively applied to solve $\mathcal{NP}$-hard problems, particularly the knapsack problem and its variants [46], over the past decades. In GAs, each solution in the search space can be represented by an individual composed of several genes. By simulating natural selection, genetic mutation, and crossover operations, GAs evolve and select individuals within the current population, aiming to find the optimal solution. The flowchart of the GA can be illustrated in Figure 1. It is worth noting that Rudolph pointed out, through Markov chain analysis, that the GA lacks global convergence, whereas the Elitist Genetic Algorithm (EGA), which selects the current population and the optimal individual from the previous generation, possesses the property of converging to the global optimum as time approaches infinity [22]. The main distinction between EGA and GA lies in the selection operators, while the same mutation operator is employed. Since this paper primarily focuses on the theoretical results of mutation operators, we will still adopt the representation of the GA in this study.

Figure 1: Flowchart of a basic genetic algorithm

The GA begins by initializing the population and generating initial solutions, and the fitness value of each individual is calculated by the objective function. The algorithm then checks if the termination condition is satisfied based on the specified number of iterations. If the termination condition is met, the results are output; otherwise, the next generation population is obtained by applying the crossover, mutation, and selection operators sequentially. Common crossover operators include single-point crossover, k-point crossover, and uniform crossover. Frequent mutation operators include flip bit mutation, swap mutation, and inversion mutation. Selection operators mainly consist of roulette wheel selection, stochastic universal sampling, and rank-based selection. Various works on these operators are available in the references [47, 48]. In this study, we utilize the widely employed single-point crossover, flip bit mutation, and roulette wheel selection operators, along with presenting pseudo-code for the GA.

Let $O_{t}^{k}(j)$ represent the value of the $j$-th decision variable in the $k$-th individual of the $t$-th $(1\leq k\leq\text{pop})$ generation where pop represents the population size, and let $\min f(O_{t})$ (resp. $\max f(O_{t})$) denote the minimum (resp. maximum) fitness value of individuals in the $t$-th generation population. Additionally, LOS refers to the local optimal solution, while $I$ and $t$ represent the maximum number of iterations and the current iteration number, respectively. Moreover, $p_{m}$ and $p_{c}$ represent the mutation probability and the crossover probability, respectively.

Algorithm 1 Genetic algorithm

  Input: $n,W,P,C,\text{pop},I,p_{c},p_{m}$

  Output: LOS

  Generating initial population $O_{1}$ and let $t=1$.

  while $t\leq I$ do

     //Crossover Operator

     for $k=1:2:\text{pop}$ do

        if ${\rm rand}\leq p_{c}$ then

           $r_{t}=\text{rand}*(n-1)+1,T=O_{t}^{k}$,

           $O_{t}^{k}(r_{t}+1:n)=O_{t}^{k+1}(r_{t}+1:n),$

           $O_{t}^{k+1}(r_{t}+1:n)=T(r_{t}+1:n).$

        end if

     end for

     //Mutation Operator

     $O_{t}=|O_{t}-[{\rm rand}(\text{pop},n)\leq p_{m}]|.$

     //Selection Operator

     The selection probability of the $k$-th individual is determined as $\frac{f(O_{t}^{k})-\min f(O_{t})+1}{\sum_{k^{\prime}=1}^{\text{pop}}f(O_{t}^{k^{\prime}})-\text{pop}\min f(O_{t})+\text{pop}}$, and resulting in $O_{t+1}$.

     $t=t+1.$

  end while

  $\text{LOS}=\max f(O_{I})$

The GA starts by generating initial solutions and then enters a loop. Next, the algorithm performs the crossover and mutation operators on the individuals determined by random numbers. These operations involve randomly generating crossover points, denoted as $t_{r}$, crossing two individuals from position $t_{r}$ to $n$, and flipping their bits. Afterward, the GA calculates the probability of each individual being selected using the roulette wheel selection operator. Particularly, as problem size increases, the difference between the objective function values of sub-optimal and local optimal solutions decreases, increasing the likelihood of sub-optimal solutions being selected. To address this, we assign a value of 1 to the objective function value of the worst solution in the current generation, then subtract this value from the objective function values of the other solutions, and add 1. Finally, the algorithm terminates the loop after a limited number of iterations and returns the local optimal solution (LOS).

It is worth noting that GAs were initially perceived as a robust search algorithm, leading to a lack of extensive attention to parameter control in the early stages of its development [22]. As the research on algorithm theory progressed, however, numerous researchers acknowledged the impact of parameter control on search outcomes, and it has since become one of the most crucial issues in the field of evolutionary algorithms [49].

3 Mainly theoretical result

3.1 A novel reduction method

Recently, a new reduction method was proposed by Dey et al.[36] for analyzing the performance of the B&B in solving random integer linear programming problems. This method divides the decision variables into different regions and constrains the branching generated by decision variables within each region. In this section, we combine traditional reduction methods [39] with the approach by Dey et al. [36] to solve the 0-1 KP. We introduce a colored regions partitioning method for 0-1 KP, as demonstrated in Theorem 3.1. This method determines the minimum or maximum number of items that should be included in the optimal solution for each colored region.

Theorem 3.1.

Given a set of items $N^{\prime}$ where the $j$-th item satisfies the inequality

$$\frac{p_{j}}{w_{j}+r/i}>\frac{p_{b}}{w_{b}}$$ (3.7)

for $i\in\mathbb{N}$, we can conclude that $\sum\limits_{j\in N^{\prime}}|x_{j}-y_{j}|\leq i-1$, where Y represents the optimal solution.

It is evident that for the 0-1 KP, the number of leaves in the search tree for items in $N^{\prime}$ is reduced from $2^{|N^{\prime}|}$ to $\sum\limits_{k=1}^{i-1}\binom{|N^{\prime}|}{k}$. Finding ways to effectively reduce the algorithmic complexity and improve the solving speed, particularly by exploring more efficient prune strategies between different colored regions, remains an intriguing question. Regrettably, existing research has not investigated the relationships between items within different colored regions. Therefore, we present the following definition based on Theorem 3.1.

Definition 5.

The item set $N_{1}=\{1,2,\dots,b-1\}$, consisting of items with a higher profit density than the break item $b$, can be divided into $m$ disjoint subsets, denoted as $N_{1}=\sqcup_{i=1}^{m}N_{1i}$. Each $N_{1i}$ subset is defined as follows:

$$N_{1i}=\{j|\lfloor\frac{rp_{b}}{p_{j}w_{b}-p_{b}w_{j}}\rfloor+1=i,e_{j}>e_{b}\}.$$ (3.8)

If no item satisfies the formula (3.8) for any integer $i(1\leq i\leq m)$, then $N_{1i}=\emptyset$.

Theorem 3.2.

Let $D_{1i}$ and $s_{1i}$ denote the item set and the number of items that are not packed into the optimal solution in $N_{1i}$, respectively, where $1\leq i\leq m$. Furthermore, we have $s_{1i}=|D_{1i}|=\sum\limits_{j\in N_{1i}}|y_{j}-1|$. All $s_{1i}$ form the vector $S=(s_{1i})_{i=1}^{m}$. Therefore, we can conclude that:

$$\sum\limits_{i=1}^{m}\frac{s_{1i}}{i}\leq 1.$$

Proof.

Given a vector $\textbf{Y}^{\prime}$, we can determine the corresponding $D^{\prime}_{1i}$ and $S^{\prime}=(s^{\prime}_{1i})^{m}_{i=1}$ of $\textbf{Y}^{\prime}$ for $i=1,2,\dots,m$ using the definitions of $D_{1i}$ and $S$.

Assume, for the sake of contradiction, that $\textbf{Y}^{\prime}$ is the optimal solution and satisfies

$$\sum\limits_{i=1}^{m}s^{\prime}_{1i}/i>1.$$

It is clear that the optimal solution $\textbf{Y}^{\prime}$ must not exceed the Dantzig bound. Therefore, we have

$$\sum\limits_{i=1}^{m}\sum\limits_{j\in D^{\prime}_{1i}}p_{j}\leq(\sum\limits_{i=1}^{m}\sum\limits_{j\in D^{\prime}_{1i}}w_{j}+r)p_{b}/w_{b}.$$ (3.9)

Using the formulas (3.8) and Theorem 3.1, we can deduce the following inequality:

$$\frac{\sum\limits_{i=1}^{m}\sum\limits_{j\in D^{\prime}_{1i}}p_{j}}{\sum\limits_{i=1}^{m}\sum\limits_{j\in D^{\prime}_{1i}}(w_{j}+r/i)}=\frac{\sum\limits_{i=1}^{m}\sum\limits_{j\in D^{\prime}_{1i}}p_{j}}{\sum\limits_{i=1}^{m}\sum\limits_{j\in D^{\prime}_{1i}}w_{j}+r\sum\limits_{i=1}^{m}s^{\prime}_{1i}/i}>\frac{p_{b}}{w_{b}}$$

and

$$\sum\limits_{i=1}^{m}\sum\limits_{j\in D^{\prime}_{1i}}p_{j} >(\sum\limits_{i=1}^{m}\sum\limits_{j\in D^{\prime}_{1i}}w_{j}+r\sum\limits_{i=1}^{m}s^{\prime}_{1i}/i)p_{b}/w_{b}$$

$$>(\sum\limits_{i=1}^{m}\sum\limits_{j\in D^{\prime}_{1i}}w_{j}+r)p_{b}/w_{b}.$$ (3.10)

Inequality (3.9) and inequality (3.10) contradict each other, proving the proposition. Hence, Theorem 3.2 is established. ∎

It is worth noting that $s_{1i}$ represents the items in $N_{1i}$ that are not selected by the optimal solution, so $s_{1i}$ belongs to the set of natural numbers, i.e., $s_{1i}\in\mathbb{N}$. Even when the integer constraints of the original problem are relaxed, i.e., $x_{j}\in[0,1]$ and $s_{1i}\in\mathbb{R}$, the conclusion still holds. This means that

Corollary 3.3.

Given a solution $\textbf{X}\in[0,1]^{n}$ and $s_{1i}=\sum\limits_{j\in D_{1i}}|x_{j}-1|$. If

$$\sum\limits_{i=1}^{m}s_{1i}/i>1,$$ (3.11)

then X is not the optimal solution.

Furthermore, we establish the following equivalent corollary to facilitate algorithmic computation, considering the presence of numerous empty sets $N_{1i}$ in Theorem 3.2.

Corollary 3.4.

Consider a vector $H=\{h_{1},h_{2},\dots,h_{n}\}$, where each element is defined as:

$$\begin{split}h_{j}=\left\{\begin{array}[]{ll}\lfloor\frac{rp_{b}}{p_{j}w_{b}-p_{b}w_{j}}\rfloor+1,&\text{if}\ e_{j}>e_{b},\\ \infty,&\text{otherwise}.\end{array}\right.\end{split}$$

If X is the optimal solution, then the following inequality holds:

$$\sum\limits_{j=1}^{n}\frac{1-x_{j}}{h_{j}}\leq 1.$$

Based on Corollary 3.4, we can derive the following corollary for items with a profit density lower than the break item $b$.

Corollary 3.5.

Consider a vector $L=\{l_{1},l_{2},\dots,l_{n}\}$, where each element is defined as:

$$\begin{split}l_{j}=\left\{\begin{array}[]{ll}\lfloor\frac{rp_{b}}{p_{b}w_{j}-p_{j}w_{b}}\rfloor+1,&\text{if}\ e_{j}<e_{b},\\ \infty,&\text{otherwise}.\end{array}\right.\end{split}$$

If X is the optimal solution, then the following inequality holds:

$$\sum\limits_{j=1}^{n}\frac{x_{j}}{l_{j}}\leq 1.$$

The number of leaves pruned in the search tree of the 0-1 KP by Corollaries 3.4 and 3.5 can be expressed using generating functions[37, 38] as follows.

Theorem 3.6.

Consider the number of leaves in the search tree, denoted as $\omega$, and the number of items in $N_{1i}$, denoted as $n_{1i}$, respectively. Additionally, let $\mathcal{X}(\lambda)$ be the polynomial representing a potential optimal solution, given by:

$$\mathcal{X}(\lambda):=\prod\limits_{i=1}^{m}\left(\sum\limits_{j=0}^{\min\{n_{1i},i\}}{n_{1i}\choose j}\lambda^{j/i}\right).$$ (3.12)

Then, the number of leaves $\omega$ is equal to the sum of the coefficients of terms with an exponent no greater than 1.

It is evident from Theorem 3.6 that Corollary 3.4 and Corollary 3.5 significantly enhance the efficiency of B&B algorithms in solving problems such as the 0-1 KP.

3.2 An improvement mutation operator and the upper bound of the mutation probability

As one of the most commonly used tools to enhance population diversity in various heuristic algorithms, the flip bit mutation operator has been widely applied due to its simplicity. Specifically, for $O_{t}^{k}(j)$, the flip bit mutation compares the mutation probability $p_{m}$ (typically set as $p_{m}\in[0.001,0.01]$) with a randomly generated number. If the generated random number is less than $p_{m}$, the value of $O_{t}^{k}(j)$ is flipped, changing 1 to 0 and 0 to 1. This procedure effectively increases the diversity of the population, enabling better exploration of the search space. To facilitate subsequent discussions, we first provide the mathematical representation of the flip bit mutation.

$$O_{t}^{k}(j)=\begin{cases}1-O_{t}^{k}(j),&\text{if }\text{rand}<p_{m},\\ O_{t}^{k}(j),&\text{otherwise},\end{cases}$$ (3.13)

where rand is a randomly generated number within the interval $[0,1]$.

From Theorem 3.6, it becomes evident that items positioned farther from the line formed by the break item and the origin are more likely to be selected in the optimal solution. Building upon this observation and combining it with a greedy algorithm, we propose the IMO. Given a mutation probability $p_{m}$, the IMO initially employs the break item to partition the set of items into two parts. Mutation operations are then performed on items belonging to different parts based on the value of $O_{t}^{k}(j)$. The mathematical representation of the IMO can be described using formula (3.14).

$$O_{t}^{k}(j)=\begin{cases}1-O_{t}^{k}(j),&\text{if }e_{j}>e_{b}\text{ and }\text{rand}<O_{t}^{k}(j)p_{m}+(1-O_{t}^{k}(j))(1-p_{m}),\\ 1-O_{t}^{k}(j),&\text{if }e_{j}\leq e_{b}\text{ and }\text{rand}<(1-O_{t}^{k}(j))p_{m}+O_{t}^{k}(j)(1-p_{m}),\\ O_{t}^{k}(j),&\text{otherwise}.\end{cases}$$ (3.14)

For a given mutation probability $p_{m}$, the IMO yields an expectation of $p_{m}n_{1i}$ for the number of unselected items in each set $N_{1i}$. However, due to the constraint imposed by Corollary 3.3, there is a maximum limit on the number of unselected items among all sets. Solutions that do not meet this constraint cannot be considered optimal. Consequently, we can further reduce the ineffective search space of the mutation operator and determine an upper bound for the mutation probability $p_{m}$.

Theorem 3.7.

Let $\overline{p}_{m}$ denote the upper bound of the mutation probability, which can be expressed as:

$$\overline{p}_{m}=\min\{\frac{1}{\sum_{j=1}^{n}1/h_{j}},\frac{1}{\sum_{j=1}^{n}1/l_{j}}\}.$$ (3.15)

Proof.

Without loss of generality, we assume that $0<\frac{1}{\sum_{j=1}^{n}1/h_{j}}\leq\frac{1}{\sum_{j=1}^{n}1/l_{j}}$. It is evident that, according to Corollary 3.3 and Corollary 3.4, if $\overline{p}_{m}>\frac{1}{\sum_{j=1}^{n}1/h_{j}}$, the individual processed by the IMO cannot be an optimal solution in terms of mathematical expectation. The proof of Theorem 3.7 is thereby established. ∎

In the past decades, numerous researchers have conducted extensive research on the optimal mutation probability for convex functions or $\mathscr{P}$ problems such as OneMax and LeadingOnes. The commonly obtained value for this probability is $\frac{1}{n}$ [18, 19, 20]. These studies have led to the following conclusions.

Proposition 1.

[25] With the increase in problem size, the mutation probability gradually tends to 0, i.e.,

$$\lim\limits_{n\rightarrow\infty}p_{m}=0.$$

As the problem size increases, it becomes apparent that for a general instance of the 0-1 KP, the number of items with profit density greater than the break item also increases. This raises an intriguing question of whether the upper bound $\overline{p}_{m}$ for the mutation probability in the IMO gradually tends to 0.

Theorem 3.8.

Let $R$ be a constant such that $R\in\mathbb{N}^{+}$ and consider an instance $\mathscr{P}$ of the 0-1 KP with $p_{j},w_{j}\in\{1,2,\dots,R\}$. We can conclude that

$$\lim\limits_{n\rightarrow\infty}\overline{p}_{m}=0.$$

Proof.

For an instance of the 0-1 KP, without loss of generality, let the profit density of the break item be constant at $p_{b}/w_{b}$, and let the residual capacity $r$ also be constant. The knapsack capacity $C$ can then be expressed as:

$$C=\sum\limits_{j=1}^{b-1}w_{j}+r.$$

Consider the item set $N_{1}$ consisting of items with profit density greater than the break item $b$, and let $R$ be any constant. It is evident that there exists a constant $M\in\mathbb{N}^{+}$ such that $m<M$ and $N_{1}$ can be expressed as the disjoint union of sets $N_{1i}$, as defined in Definition 5. In other words, for any $j\in N_{1}$, we have:

$$\frac{p_{j}}{w_{j}+r/M}>\frac{p_{b}}{w_{b}}.$$

Furthermore, by leveraging Corollary 3.4, we are able to acquire a vector $H$ that is associated with the given instance. Consequently, the following conclusion can be derived:

$$\lim\limits_{n\rightarrow\infty}\overline{p}_{m}\leq\lim\limits_{n\rightarrow\infty}\frac{1}{\sum_{j=1}^{n}1/h_{j}}\leq\lim\limits_{n\rightarrow\infty}\frac{M}{n}=0.$$

The proof of Theorem 3.8 is complete. ∎

It is notable that when the profit and weight values of items in a 0-1 KP instance are limited, as exemplified in Theorem 3.8, the instance can be solved within $O(Rn^{2})$ time complexity by Dynamic Programming[8, 9]. Since $R$ is a constant, the instance belongs to the $\mathcal{P}$ complexity class. Consequently, it becomes an interesting question of whether the upper bound of the mutation operator in the IMO still tends towards 0 when the instance of the 0-1 KP no longer limits the profit and weight values of items.

Theorem 3.9.

For any given constant $\theta$ within the interval $(0,1)$, it is feasible to construct an instance of the 0-1 KP with $p_{j},w_{j}\in\mathbb{N}^{+}$, and we have

$$\lim\limits_{n\rightarrow\infty}\overline{p}_{m}=\theta.$$

Proof.

Following the same approach as the proof of Theorem 3.8, we assume that the profit density of the break item remains constant at $p_{b}/w_{b}$, while maintaining a constant residual capacity $r$. Furthermore, the knapsack capacity $C$ must satisfy

$$C=\sum\limits_{j=1}^{b-1}w_{j}+r.$$

Without loss of generality, we set $\theta=\frac{1}{2}$. Next, we construct the profit and weight of items in the instance, ensuring that their profit densities exceed that of the break item. Let $p_{j}$ and $w_{j}$ represent the profit and weight, respectively, of the $j$-th item, satisfying

$$\lfloor\frac{rp_{b}}{p_{j}w_{b}-p_{b}w_{j}}\rfloor+1=2^{j-1}.$$ (3.16)

In other words, we have $h_{j}=2^{j-1}$. Since $p_{j},w_{j}\in\mathbb{N}^{+}$, it is clear that there exist $p_{j}$ and $w_{j}$ that satisfy equality (3.16). Consequently, we can conclude that

$$\lim\limits_{n\rightarrow\infty}\overline{p}_{m}=\lim\limits_{n\rightarrow\infty}\frac{1}{\sum_{j=1}^{n}1/h_{j}}=\frac{1}{2}=\theta.$$

The proof of Theorem 3.9 is complete. ∎

3.3 Comparison analysis

In $\mathcal{NP}$-hard problems, given two solutions $\textbf{X}^{\prime}$ and $\textbf{X}^{\prime\prime}$, even if $\|\textbf{X}^{\prime}-\textbf{Y}\|<\|\textbf{X}^{\prime\prime}-\textbf{Y}\|$, we cannot conclude that $f(\textbf{X}^{\prime})\geq f(\textbf{X}^{\prime\prime})$, where Y is the optimal solution. In other words, solutions with a smaller Hamming distance to the optimal solution Y do not necessarily dominate solutions with a larger Hamming distance to Y. Therefore, the existing convergence analysis of GAs is limited to specific scenarios. For a general instance of $\mathcal{NP}$-hard problems, the expected runtimes of the GA often exceed $O(2^{n})$, which lacks strong theoretical guarantees for the performance.

To demonstrate the performance of the IMO, we evaluate the algorithm by considering the probability $\tau(\text{Algorithm})$ of hitting the optimal solution Y from the zero vector 0 within a single iteration for a given instance $\mathscr{P}$. This evaluation is conducted to assess the performance of the IMO and the MO on general large-scale instances.

Example 3.10.

Given an instance $\mathscr{P}$ of the 0-1 KP with $n$ items, let $\lambda_{1}$ and $\lambda_{2}$ denote the number of selected and unselected items within the first $b-1$ items, respectively. Similarly, $\lambda_{3}$ and $\lambda_{4}$ represent the number of selected and unselected items after the $(b-1)$-th item, respectively. We have $\lambda_{1}=\sum\limits_{j=1}^{b-1}|y_{j}-1|$, $\lambda_{2}=\sum\limits_{j=1}^{b-1}|y_{j}|$, $\lambda_{3}=\sum\limits_{j=b}^{n}|y_{j}|$, and $\lambda_{4}=\sum\limits_{j=b}^{n}|y_{j}-1|$.

Theorem 3.11.

If $\lambda_{1}<\lambda_{2}$ and $p_{m}<0.5$, then $\tau({\rm IMO})>\tau({\rm MO})$.

Proof.

$$\frac{\tau({\rm IMO})}{\tau({\rm MO}} =\frac{(1-p_{m})^{\lambda_{1}}p_{m}^{\lambda_{2}}p_{m}^{\lambda_{3}}(1-p_{m})^{\lambda_{4}}}{p_{m}^{\lambda_{1}+\lambda_{3}}(1-p_{m})^{\lambda_{2}+\lambda_{4}}}$$

$$=\frac{p_{m}^{\lambda_{2}+\lambda_{3}}(1-p_{m})^{\lambda_{1}+\lambda_{4}}}{p_{m}^{\lambda_{1}+\lambda_{3}}(1-p_{m})^{\lambda_{2}+\lambda_{4}}}$$

$$=(\frac{1-p_{m}}{p_{m}})^{\lambda_{2}-\lambda_{1}}$$

$$>1.$$

This completes the proof of Theorem 3.11. ∎

Theorem 3.11 demonstrates that when the number of items in the optimal solution exceeds the number of items not in the optimal solution before the break item, the probability of the IMO hitting the optimal solution in a single iteration is higher, which aligns with empirical experimental data. Moreover, since GAs are primarily employed for solving large-scale $\mathcal{NP}$-hard problems, which rarely exhibit the scenario of Example 3.10 with $\lambda_{1}\geq\lambda_{2}$, and the mutation probability $p_{m}$ does not exceed $0.5$, the IMO generally outperforms the MO in large-scale instances.

4 Conclusion and future work

Heuristic algorithms with parameters exhibit significant performance differences depending on the value of the parameters. Due to limitations in mathematical tools, the problem of finding the optimal parameters for heuristic algorithms has received extensive attention in recent decades. Based on the assumption that $\mathcal{NP}\neq\mathcal{P}$, we propose a novel reduction method and apply it to compute the upper bound of the mutation probability.

Furthermore, our method not only proves that in the 0-1 KP, when the weight and profit of items are limited, the mutation probability tends toward 0, but also demonstrates that the mutation probability can tend toward a constant within the open interval $(0,1)$ when the weight and profit values of items are unrestricted.

For future work, we can approach it from two perspectives. Firstly, the upper bound of the method deserves further improvement. Our research only utilized the linear relaxation technique, and it would be worthwhile to explore better upper bound computation methods, such as the Lagrangian relaxation technique. Secondly, Dey et al. have demonstrated that such reduction methods can apply to not only one-dimensional cases but also multidimensional problems. Therefore, the application scope of these reduction methods can be further expanded, such as multidimensional knapsack problems.