Credibilistic Mean-Semi-Entropy Model for Multi-Period Portfolio Selection with Background Risk

In financial markets, investors will face not only portfolio risk but also background risk. This paper proposes a credibilistic multi-objective mean-semi-entropy model with background risk for multi-period portfolio selection. In addition, realistic constraints such as liquidity, cardinality constraints, transaction costs, and buy-in thresholds are considered. For solving the proposed multi-objective problem efficiently, a novel hybrid algorithm named Hybrid Dragonfly Algorithm-Genetic Algorithm (HDA-GA) is designed by combining the advantages of the dragonfly algorithm (DA) and non-dominated sorting genetic algorithm II (NSGA II). Moreover, in the hybrid algorithm, parameter optimization, constraints handling, and external archive approaches are used to improve the ability of finding accurate approximations of Pareto optimal solutions with high diversity and coverage. Finally, we provide several empirical studies to show the validity of the proposed approaches.


Introduction
As a research field, portfolio selection is used to accomplish the investments in financial markets by spreading investors' capital among several different assets considering return and risk. Since the pioneering work of Markowitz [1] in single-period investment problems, the mean-variance portfolio selection problem has attracted much attention and has become a research hotspot. By introducing different risk measures, a large variety of portfolio selection models have been presented, such as the mean-variance-skewness model [2], the mean-conditional value at risk (CVaR) model [3], the mean-value at risk (VaR) model [4], the mean-semi-variance model [5] and the minimax risk model [6]. In addition, entropy can also be used as a risk measure because it does not depend on symmetric membership functions and can be calculated from non-metric data. Philippatos and Wilson [7] first replaced variance with entropy as a risk measure. Later, Rödder et al. [8] provided a new and efficient method for determining the portfolio weights on the basis of a rule inference mechanism with both maximum entropy and minimum relative entropy. Nawrocki and Harding [9] provided two alternative weighted computations of entropy to measure portfolio risk. Usta and Kantar [10] presented a multi-objective model founded on mean, variance, skewness and entropy to adequately diversify the portfolio. Yu et al. [11] discussed the performance of the models with diverse entropy measures by comparing the mean-variance efficiency, portfolio values, and diversity.
Traditionally, researchers dealt with the uncertainty of portfolio selection problems by applying probability theory. For example, Beraldi et al. [12] proposed a mean-CVaR model considering a complex If the set function Cr satisfies the aforementioned four axioms, the {Θ, P, Cr} will be credibility space. Definition 1. Let ξ be a fuzzy variable defined on the credibility space {Θ, P, Cr} with membership function µ{x}. For any set A of real numbers, the credibility is defined as Credibility measure is an increasing function of set A. It is obvious that the credibility measure is self-dual. (2) Theorem 1. Let ξ be a fuzzy variable with a finite expected value; let µ and ν be any given two real numbers. Then, Theorem 2. Suppose that ξ and η are two independent fuzzy variables. The expected value of these variables are finite. Then, for any numbers µ and ν, Example 1. The expected value for the trapezoidal fuzzy variable ξ = (ξ a , ξ b , ξ c , ξ d ) is given by For the sake of determining the credibility of a fuzzy event, the trapezoidal fuzzy variable ξ has a membership function illustrated below: Then, the credibility of fuzzy event {ξ ≤ r} is given as below:

Notation
At the beginning of the investment, we assume that the investor's initial wealth is W 1 . The investor allocates W 1 among n risky assets and a risk-free asset at the start of T − 1 period and acquires the ultimate wealth at the final period T. As a matter of convenience, we list all the symbols used below: i: the exponents for the n risky assets, i = 1, 2, ...n. t: the exponents for the T investment period, t = 1, 2, ...T. W t : the wealth accumulated at the start of the t-th investment period.
x it : the proportion of the whole wealth that investor spreads to the i-th risky asset during the t-th investment period. r it : fuzzy variables that represent the return rate on the i-th risky asset during the t-th investment period, . r f : the variable that represents the return rate on the risk-free asset. r b : the fuzzy variable that represents the return rate on background asset, r b = (ξ a , ξ b , ξ c , ξ d ). ub it : the upper limit that can be assigned to the i-th risky asset during the t-th investment period. lb it : the lower limit that can be assigned to the i-th risky asset during the t-th investment period. f it : the cost on transaction of the i-th risky asset during the t-th investment period. L it : the fuzzy variables that represent the turnover rates on the i-th risky asset during the t-th investment period, . L t : the acceptable minimum expected liquidity during the t-th investment period. m it : the 0-1 variables expressing whether the i-th risky asset is chosen for the portfolio during the t-th investment period or not: if the i-th risky asset is chosen to the portfolio during the t-th period, 0, otherwise, Z t : the desired number of risky assets that can be chosen for each investment interval.

Maximize Ultimate Wealth
According to Equations (3) and (5), the expected value of the portfolio x t during the t-th investment period is Moreover, from Equation (5), the expected value of the background asset is Additionally, we apply a V-shaped function that expresses the differences between the two diverse portfolios of the two adjacent periods. Then, the expense on transaction of the i-th risky asset during the t-th investment period is f it |x it − x it−1 |. Furthermore, from Equations (8) and (9), the net return rate at period t can be denoted as Then, the expected value of the wealth at the beginning of the period t + 1 is expressed as W t+1 = W t (1 + R t ). Thus, after accomplishing the investment through the entirety of investment periods, from Equation (10), the ultimate wealth at the end of the period T is denoted as

Minimize Risk
Fuzzy entropy has been extensively applied to characterize uncertainty since Luca and Termini [58] first defined a non-probabilistic entropy in the framework of fuzzy set entropy. Since then, various definitions for fuzzy entropy have been proposed; see, for example, Li and Liu [59], Zhou et al. [60], Qin et al. [61], and Xu et al. [62]. Fuzzy entropy is more convenient than fuzzy variance because it does not depend on symmetric membership functions and can be calculated from non-metric data. It is used to express the uncertainty of both low and high extreme returns. However, what investors really dislike is the downside uncertainty. Therefore, fuzzy semi-entropy introduced by Zhou et al. [60] matches reality more exactly as the downside risk measure. In this section, we used the semi-entropy to quantify the portfolio downside risk.

Definition 3.
Assume that there is a continuous fuzzy variable δ whose expected value E[δ] is finite. The function o(x) is equal to Cr{δ = x}. Then, the semi-entropy of δ is defined as [60] where Because of S e (0) = 0, the semi-entropy of δ can be transformed into Theorem 3. Suppose there is a continuous fuzzy variable δ whose expected value E[δ] is finite. Then, for these two real numbers λ and ω with λ > 0, where Furthermore, according to Equation (16), we obtain the cumulative portfolio risk with background risk as follows: In Equation (17), according to the definition of semi-entropy and Equation (16), where Similarly, the semi-entropy of background asset where

• Liquidity
In the process of making a portfolio decision, one of the key elements that should be considered is liquidity for investors. It measures the degree of probability that investors will convert an asset into income. Investors prefer assets with higher liquidity because their returns tend to rise over time. Generally, liquidity is measured by the turnover rate of assets. Because turnover rates cannot be precisely predicted, we suppose that the turnover rates of risky assets are fuzzy variables characterized by trapezoidal numbers. On account of the former discussion, by Equation (5), the constraint of the portfolio liquidity is expressed as • The desired number of risky assets that are selected into the portfolio during the t-th investment period is expressed as • The risk-free asset constrained in each period is • The lower and upper limits that can be assigned to the i-th risky asset during the t-th investment period are given as • Whether the i-th risky asset is selected into the portfolio during the t-th investment period is shown as • No short selling of assets during any investment period

The Proposed Model
Over the entire investment horizons, investor intends to obtain the greatest final wealth and minimize the risk at the same time to find a first-rank invest strategy. Then, we supply the multi-objective model for multi-period portfolio selection problems in the following: s e (r it ) and s e (r b ) in the proposed model are defined by Equations (18) and (19), respectively.

Standard Dragonfly Algorithm (DA)
The static and dynamic swarming behaviors of dragonflies inspire the DA algorithm. These two behaviors represent the exploration phase and the exploitation phase, which are two major phases of the meta-heuristic algorithm. Five diverse operators determine the movement of swarm dragonflies:

• Separation
For the individual i, its separation is calculated as S i = − N ∑ k=1 (P − Pn k ). Pn k denotes the k-th adjacent individual's position. P denotes the current individual's position. N is the number of neighboring individuals.

• Alignment
For the individual i, its alignment is given as

• Cohesion
For the individual i, the cohesion is calculated as Pn k /N − P.

• Attraction and Distraction
For the individual i, the index for an individual being attracted by a food source is evaluated as F i = P + − P, where P + is the food source's position. In addition, the index for an individual fleeing an enemy is calculated as E i = P − + P, where P − is the enemy's position.
In order to find some new individuals in the search space, two vectors are employed. The step vector ∆P is used to update the locations of individuals, and the position vector P is introduced for simulating movements of the individuals. The movement directions of the individuals are given by the ∆P. If an individual has at least one neighbor, then ∆P is evaluated as In Equation (27), the separation weight is indicated by s, the alignment weight is shown by a, the cohesion weight is represented by c, and the food element and the enemy element are denoted as f and e, respectively. Furthermore, t is the iteration counter. According to ∆P in Equation (27), P is given as If an individual has no neighbors, the Lévy Flight equation will be applied to update P. This equation can improve the randomness, global search capacity and chaotic behavior of individuals. P is calculated as In Equation (29), the equation of Lévy flight is In Equation (30), η 1 and η 2 are two random numbers taking values in [0, 1], and ϑ is a constant, γ is calculated as where Γ(χ) = (χ + 1)!.

The Hybrid DA-GA for the Proposed Model
A good metaheuristic algorithm should better balance exploration and exploitation processes. The exploration process is used to investigate the new search space to find great global optima, while the exploitation process is used to focus on the search of local areas. Excessive exploitation results in premature convergence, while overmuch exploration leads to slow convergence. DA has advantages in exploring the global search space by using the food source and enemy source. However, the use of Lévy Flight results in a large movement that leads to local convergence and pushes the algorithm apart from the global optimum [56]. In addition, NSGA II, developed by Deb et al. [63], is a well-known meta-heuristic approach for solving multi-objective optimization problems. It has an improved mechanism that depends on the non-domination rank and the crowding distance and conducts constraints by using an adapted explanation of dominance instead of the penalty functions. Thus, NSGA II has a good ability to attain diverse and uniformly distributed Pareto solutions. In this paper, for the sake of solving the proposed model efficiently, a novel hybrid algorithm named HDA-GA is developed by combining the strengths of DA and NSGA II.

Parameter Optimization
In the static swarm of DA, the probability of alignments is low, while the probability of cohesion is high. In order to enhance the information exchange of the dragonflies from global exploration to local exploitation, dragonflies are assigned with higher alignment weights and lower cohesion weights when the global space is explored and designed on the contrary when the local area is exploited. Therefore, the exponential function is introduced to adjust the swarming elements a and c. The factors a and c are given as follows: where h is adaptively decreased as the iteration increases. Moreover, in order to enhance the randomness, the standard DA selected the positions of food source P + and enemy P − by using a roulette-wheel mechanism. However, in the global search space, it may lead to poor exploration ability. Inspired by the ideas in [51], we propose a new method for choosing food sources and enemies. P gbest and P gworst are defined as the best and the worst solutions in each iteration. The selections of P + and P − are given as follows:

Constraints Handling
Note that the standard DA only considered the non-constrained situation. However, there exist constraints in the proposed model that cannot be ignored. In this paper, to handle the constraints, we employ the constrained domination approach proposed by Deb et al. [63].
If any of the conditions below is true, a solution S k is constrained-dominated by another solution S j . (1) Both solutions are feasible, and solution S k is dominated by solution S j ; (2) The feasible solution is S j , but the infeasible one is S k ; (3) Both are infeasible, but comparing the constrained violations of these two solutions, the violation solution S j has is smaller.
For the t − th inequality constraint g t (s) ≤ 0 and equality constraint h t (s) = 0, the constrained violation is estimated as where ι is a tolerance coefficient that violates the equality constraints. After the normalization of cv t , the constrained violation of solution S j is given as For the purpose of drifting the solutions towards the Pareto front and making the Pareto-optimal set as diverse as possible, a joint strategy combining the constrained non-dominated sorting and crowding distance assignment is implemented. In the strategy, how close a solution is to its neighbors is measured by crowding distance distance k . Diversity improves with larger distance k . In the proposed algorithm, the crowding distance distance k measure introduced by Deb et al. [63] is employed and calculated as follows: In Equation (38), the maximum and minimum of the first objective function is shown as F 1max and F 1min , respectively. Similarly, the maximum and minimum of the second objective function are illustrated as F 2max and F 2min , respectively. The constrained non-dominated sorting pseudo-code is summarized as Algorithm 1.
1: Classify feasible and infeasible groups in the population by Equation (37) 2: For p = 1 to f easible_population do 3: Calculate S p , a set of solutions that the p − th individual dominates 4: Calculate n p , the number of individuals that dominate the p − th individual 5: End for 6: Create first front whose n p = 0 7: While (n p > 0) 8: Create subsequent fronts by traversing S p 9: Crowding distance assignment by Equation (38) 10: End While 11: For q = 1 to in f easible_population do 12: Sort infeasible individual by Equations (36) and (37) 13: End for 14: Combine the feasible and infeasible solutions

External Archive
An external archive is widely used to solve multi-objective problems and to maintain the Pareto optimal solutions during optimization. The standard MODA applies an archive to retain the best elite solutions and updates the archive with respect to the non-dominated sorting. However, the updating progress deletes the infeasible solutions directly. It did not consider the constrained situation either. Based on constrained dominate rules and crowding distance, an external archive is used to improve the speed of convergence and retain the diversity of the solution set. The archive is divided into two subsets, Archive 1 and Archive 2 . Archive 1 saves solutions obtained by DA, while Archive 2 saves solutions solved by NSGA II. Finally, Archive 1 and Archive 2 make up a new set New_Archive for the next generation. Initially, this archive is empty. As the iteration goes by, feasible and infeasible solutions enter the archive, and the size of the archive may be huge. If the archive is full, one or more than one solution may be deleted. The progress of this method is summarized as pseudo-code shown in Algorithm 2.
Through the above discussions, Algorithm 3 describes the proposed hybrid algorithm. In the hybrid algorithm, both DA and NSGA II start with the same initial population. The external archive is divided into two parts, where one retains feasible solutions and the other saves infeasible solutions during each iteration. Each of the two parts is evolved by a respective algorithm and then recombined in the updating archive process. Algorithm 2 Update archive.
1: Classify the population by Equation (37) 2: Divide the archive to Archive 1 whose CV t = 0 and Archive 2 Whose CV t = 0 3: While (N Archive 1 > 0) 4: estimate the rank of each solution according to the Equation (37) 5: Constrained non-dominated sorting by Algorithm 1 6: Calculate the crowding distance by Equation (38) 7: End While 8: While (N Archive 2 > 0) Sort by Equation (37) Set the distance to in f Select individuals from the particular front based on crowding distance by Equation (38) 12: end if 13: For i = 1 to Archive1 do 14: Find their neighbors with respect to the Euclidean distance 15: Calculate S, A, C, F and E 16: If an individual has one neighbor at least 17: Update ∆X t by Equation (27) and X t+1 by Equation (28) 18: end if 19: If an individual has no neighbor 20: Update X t+1 by Equation (29) Table 1, and Table 2 shows the fuzzy turnover rates of these 10 risky assets. The background asset returns are given by experts' estimations. In this empirical study, we hypothetically set the initial wealth as W 1 = 1, the lower and upper bounds are set as u it = 0.1 and l it = 0.5, respectively, the unit transaction cost is f t = 0.003, and the desired number of risky assets chosen for the portfolio during the t-th investment period is Z t = 5. In addition, we assume that n = 10 and T = 3. The fuzzy variable r b = (0.080, 0.090, 0.109, 0.121) is the return rate on a background asset, the return rate on risk-free assets is r f = 0.01, and the accepted minimum expected liquidities during each investment interval are designed as L 1 = 0.0045, L 2 = 0.0035, and L 3 = 0.0025.
The parameters in NSGA II and MODA are equal to those in HDA-GA. MOPSO: The modulus of personal learning c1 is 1, the modulus of global learning c2 is 2, and the initial weight w is 0.5.
MOSSA: The initial range r is 0.2, and the initial max velocity V max is 0.04.
MOABC: The f ood_Number is 200, and the limit is 50.
In addition, each algorithm independently runs 30 times, and the average results are obtained after running.

Performance Measure Metrics
Five performance metrics, GD, Spacing, Diversity, CM and MPFE, are selected to compare the performances of the algorithms. Generation Distance (GD): This convergence metric is employed to compute the distance between the approximated Pareto frontier and the true Pareto frontier. It is calculated as [63] where N is the number of the obtained solutions, and d m is the minimum Euclidean distance between each of the obtained solutions and the true Pareto frontier. A smaller value of GD means that the obtained Pareto frontier is closer to the true Pareto frontier. Spacing: This diversity metric is applied to measure the propagate of the obtained values. It is evaluated as [64] where d k is the minimum distance between the k − th solution and its adjacent solutions, where where d m is the Euclidean distance between the solution obtained with the algorithm and the nearest solution on the Pareto frontier. The smaller the value of this metric is, the better the convergence toward the true Pareto frontier. Maximum Pareto front error (MPFE): This convergence-diversity metric is employed to measure the quality of the obtained solutions in terms of diversity and convergence on a single scale. It is expressed as [64] where Q is the number of objective functions and P is the number of the Pareto solutions. MPFE aims to find the maximum minimum distance between each solution obtained with the algorithm and the corresponding nearest solution on the Pareto frontier. The convergence and the diversity of the algorithm improve with smaller values of this metric.

Experimental Results Based on the Zdt Functions
In this section, we select four ZDT functions as benchmarks and present a comparison of these functions to verify the validity of the proposed HDA-GA. The details of the four ZDT functions are in Appendix A. Tables 3 and 4 show the best (Best), mean (Mean) and standard deviation (SD) of the five performance metrics. The bold fonts indicate better results. It can be easily observed that the proposed HDA-GA is superior to the other five algorithms within the five performance metrics.
ZDT1 is a relatively easier problem than the other three ZDT problems. From Table 3, MOPSO and MOSSA have better SD than HDA-GA. However, HDA-GA has the smallest Mean of the five metrics among the six algorithms, which means that HDA-GA converges to the Pareto frontier with the best distribution, spread, and diversity. Five disjoint curves make up the Pareto front of ZDT3. With respect to GD and Diversity, although MOPSO can obtain the Best, HDA-GA performs better between Mean and SD. In addition, HDA-GA has the smallest Best and Mean of two metrics, CM and MPFE. Moreover, HDA-GA owns a better Spacing than others, which means solutions produced by HDA-GA have a better distribution than others.
ZDT6 is another difficult problem for many multi-objective optimization algorithms to achieve a set of solutions with good convergence and diversity. From Table 4, for GD and CM, although MOPSO has smaller SD, HDA-GA performs better in Best and Mean than the others. For diversity metrics Spacing and Diversity, solutions produced by HDA-GA spread out better over the Pareto frontier with a better distribution. The results of MPFE demonstrate a superior convergence and diversity ability of HDA-GA. For ZDT2, although MOPSO and MOSSA perform more stably than HDA-GA with respect to GD and Spacing, HDA-GA has the smallest Mean, Best and SD of Diversity, CM and MPFE, which indicates that HDA-GA finds a better distribution and spread with a smaller convergence metric than others.
Based on the above discussion, HDA-GA has a superior convergence and diversity ability with a better distribution and spread. It indicates that HDA-GA outperforms the other algorithms in most of the performance metrics.

Experimental Results Based on the Proposed Model
This section presents three cases with different cardinality constraints. For the proposed model, the minimum (Min), maximum (Max), mean (Mean), standard deviation (SD) and range (Range) of the results found by six different algorithms are revealed in Table 5. The bold fonts indicate better results. Given the comparisons among the six algorithms, HDA-GA can own the smallest mean value in all the cases. In addition, according to the comparison of min and max index, we can see that HDA-GA can acquire a set of non-dominated solutions with better distribution. Finally, the comparison of Range index illustrates that HDA-GA can search space reliably and extensively. Although the MODA is more stable than the HDA-GA in terms of SD index, it is easier for MODA to fall into local optimization. These results indicated that HDA-GA performs better than the other algorithms. Moreover, for a fair comparison of the performances among the algorithms, GD, Spacing, Diversity, CM, and MPFE are employed as the performance measurement metrics. Table 6 presents some results in terms of the five metrics above. For GD and CM, the index values indicate results obtained by the proposed HDA-GA are closer to the Pareto front than the other algorithms in the three cases. Meanwhile, for Spacing and Diversity, HDA-GA performs better than the other algorithms, which means that it finds a better spread and distribution metric than others. Moreover, for MPFE, HDA-GA has a superior convergence and diversity ability. Figures 1-3 display the Pareto front and the efficient frontiers of the six algorithms under the three cases above. It can be seen that the proposed HDA-GA can obtain a set of non-dominated solutions that approach the Pareto front properly. Moreover, we can see that the proposed HDA-GA performs better with accurate convergence, preferable coverage, and better diversity.

Experimental Results with and without Background Risk
We present four cases to analyze the impact of background risk in the proposed model. Case 1: Without background risk asset (BR); Case 2: With background risk asset BR1 whose fuzzy return is r b 1 = (0, 024, 0.027, 0.0327, 0.0363); Case 3: With background asset BR2 whose fuzzy return is r b 2 = (0, 040, 0.045, 0.0545, 0.0605); Case 4: With background asset BR3 whose fuzzy return is r b 3 = (0, 080, 0.090, 0.109, 0.121). The experimental results indicate that the background risk has a significant impact on the portfolio selection.
From Table 7, it can be observed that cases considering background risk have higher returns and risk than that without background risk. Ignoring background risk will cause the underestimation of risk and the reduction of return in the actual investment. In addition, Figure 4 shows the Pareto frontiers of the above four cases. The shapes of the Pareto frontiers are approximately the same, and the Pareto frontier moves right as the background risk is concerned. It can be observed that there is a positive correlation between the background asset return and portfolio return. When the risk is the same, a portfolio with background risk can obtain a higher return than that without background risk. It indicates that considering background risk avoids the reduction of return in the actual investment and the ignorance of the potential income in the actual investment. Moreover, the risk of background assets is positively correlated with portfolio risk. When the return is the same, a portfolio with background risk is riskier than one without background risk. Considering background risk can prevent investors from underestimating the investment risk and ignoring the potential risk.

Conclusions
In the real world, investors usually need to optimize the portfolio strategies from time to time. In this paper, we proposed a mean-semi-entropy model based on the credibility theory by taking buy-in thresholds, cardinality, liquidity, and transaction costs into account. In particular, background risk is also considered in the proposed model. To solve the proposed multi-objective model, a hybrid algorithm, HDA-GA, combining the advantages of dragonfly algorithm (DA) and non-dominated sorting genetic algorithm II (NSGA II), is developed. Finally, we conducted a series of experiments to demonstrate the effectiveness of the proposed model and the hybrid algorithm. The numerical results showed that (1) the proposed algorithm HDA-GA is superior to the other five algorithms, namely, NSGA II, MODA, MOPSO, MOSSA, and MOABC, with accurate convergence, preferable coverage, and better diversity; (2) the mean-semi-entropy model can lead to more distributive investments; and (3) considering background risk will prevent investors from the underestimation of risk in the actual investment.
Future research directions include but are not limited to the following: (1) considering a more general transaction cost structure as in Beraldi et al. [12]; (2) extending the proposed model by adding other constraints of real markets such as minimum transaction lots, skewness, and class constraints; and (3) applying other metaheuristic algorithms such as the estimation of distribution algorithm (EDA), the krill herd (KH) algorithm, and bacterial foraging optimization (BFO) for solving the proposed model. Author Contributions: Conceptualization, formal analysis, visualization, and supervision, J.Z.; validation, investigation, data curation, writing-original draft preparation, Q.L.; methodology, software, resources, writing-review and editing, and project administration, J.Z. and Q.L.
Funding: This research was supported by the Humanity and Social Science Foundation of Ministry of Education of China (No. 19YJAZH005).

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.