Selection of Third-Party Reverse Logistics Service Provider Based on Intuitionistic Fuzzy Multi-Criteria Decision Making

: The scientiﬁc selection of a third-party reverse logistics service provider (3PRLP) is helpful for enterprises to obtain the expected ecological and economic beneﬁts. Different enterprises have different requirements for 3PRLP selection and, thus, should adopt personalized and simpliﬁed evaluation criteria. However, there is a lack of research on criteria screening. Therefore, this paper proposes a criteria screening method based on a rough set for the ﬁrst time. The 3PRLP selection is a multi-criteria decision making (MCDM) problem in essence, and different criteria can be expressed in different forms. The existing research mostly uses one method for 3PRLP selection, lacking the comprehensive application of various methods. In this paper, various criteria values are transformed into intuitionistic fuzzy numbers (IFNs) for the comparison and combination of various intuitionistic fuzzy MCDM methods. In terms of criteria weighting, a subjective weighting method based on an analytical network process (ANP) is proposed due to the possible correlation between the criteria at the same level. Meanwhile, an objective weighting method based on intuitionistic fuzzy entropy is proposed. The subjective and objective weights are integrated to form the more scientiﬁc combination weights. Combining the modeling principles of different intuitionistic fuzzy MCDM methods, the representative methods under each principle are chosen to build a combination evaluation idea that integrates multiple single evaluation models, and the speciﬁc evaluation steps are given, including the single evaluation, Kendall compatibility test, combination evaluation, and Spearman consistency test. An illustrative example of 3PRLP selection is provided to verify the feasibility of the methods of criteria screening and weighting and the combination evaluation idea.


Introduction
With the continuous development of the global economy, more and more resources are consumed. However, resource waste, improper disposal of a large amount of waste, and the resulting environmental damage also occur frequently. In order to improve the utilization rate of resources, save costs, and reduce environmental damage, many enterprises have taken a series of measures, such as resource conservation and recycling, reverse logistics (RL), and green material supply [1]. Reverse logistics can be defined as "the process of planning, implementing, and controlling the efficient, cost-effective flow of raw materials, in-process inventory, finished goods, and related information from the point of consumption to the point of origin for the purpose of recapturing value or proper disposal" [2]. In the past few decades, RL has been highly valued by academics and entrepreneurs for three main reasons. First, more and more enterprises have realized that a win-win situation can be achieved in corporate image and interests through RL activities, such as recycling of the back-flow products and waste materials [3]. Second, the development of e-commerce has led to a sharp increase in online sales, meaning that the return flow of goods has continued to increase [4]. Third, environmental protection laws and policies, corporate social responsibility, and other factors induce the development of RL [5]. Due to the complexity of RL related activities, an enterprise must have a professional team responsible for this part of the business and must also provide the necessary infrastructure and expertise to maintain its normal operation. However, most enterprises tend to focus on their core competitiveness, and they do not handle RL business by themselves, but outsource it to a trusted third-party RL provider (3PRLP) [6]. In reality, there are many 3PRLPs to choose from. These 3RPLPs differ in professional level, service content, service cost, output benefit, etc. Therefore, an enterprise should select the best provider from these 3PRLPs according to its own requirements. This means that the selection process of the best 3PRLP is a complex multi-criteria decision making (MCDM) problem. However, different enterprises have different requirements for RL. For example, an enterprise with strong economic strength may not consider the cost factor but might instead pay more attention to cultural compatibility and try to highlight the exchange of corporate culture through cooperation with 3PRLPs. Another enterprise, limited by capital, may emphasize cost more than cultural compatibility. This shows that different enterprises may use different evaluation criteria when selecting the best 3PRLP. How to screen the evaluation criteria that meet the personalized requirements of one enterprise from many alternative criteria is an important part of the 3PRLP election. In addition, the expression forms of these criteria are different, which can be expressed as certain crisp numbers, uncertain interval numbers, linguistic variables, fuzzy numbers, hesitant fuzzy numbers, intuitionistic fuzzy numbers (IFNs), Fermatean fuzzy numbers, interval-valued IFNS, bipolar complex fuzzy numbers (BCFNs), linear Diophantine fuzzy numbers (LDFNs), etc. There are also many applicable MCDM models. How to scientifically choose MCDM models, compare and integrate the results of various models, and give a more reliable and consistent conclusion are also important aspects that should be paid attention to in the 3PRLP selection.
However, the existing research on screening evaluation criteria that meet the personalized requirements of enterprises has yet to be seen. The rough set method can be used to analyze whether there is redundancy among condition attributes relative to the decision attribute in a knowledge representation system, and an realize effective reduction in conditional attributes on the premise of meeting decision support [7]. We use the rough set method to screen the 3RPLPs evaluation criteria that meet the personalized requirements of enterprises. Due to the complexity of the 3PRLP selection problem, most research only considers one expression of the criteria and one MCDM model. In fact, there are many models applicable to a certain type of expression, and there are no strict pros and cons among them. Even some models that consider hybrid expressions only give the subjective weights of the criteria, but fail to give the objective weights by the data dispersion degree of the criteria. Thus, IFN, through the concepts of membership, non-membership, and hesitation [8], can better reflect the uncertainty of information than the special case forms, such as the crisp number and interval number. Compared with the generalized forms, such as interval-valued IFN, BCFN, neutrosophic numbers, and LDFN, IFN is easier for evaluators to understand and assign. Moreover, the intuitionistic fuzzy entropy reflecting information difference is relatively mature in theory, and there are many intuitionistic fuzzy MCDM models, which have been used for 3PRLP selection [9,10]. Therefore, we apply IFNs to unify various expression forms of all criteria, calculate the criteria's objective weights by means of a intuitionistic fuzzy entropy method, and determine the combination weights with an ANP subjective weighting method, which can reflect the association among the criteria at the same level. Based on various single intuitionistic fuzzy MCDM models, we propose a combination evaluation idea for their comparison and integration. The key novelties of this research can be presented as follows: (1) A rough set method is proposed for criteria screening, which can reflect the comprehensive judgment of experts on the condition attributes (i.e., evaluation criteria) that affect the 3PRLP evaluation of an enterprise; (2) A combined weighting method is proposed, which can combine the subjective weighting based on ANP (analysis network process) and the objective weighting based on Hesitant fuzzy linguistic terms AO Linguistic rating variables [35] IFNs TOPSIS Fuzzy entropy [9] Linguistic variables quantified with IFNs GRA ANP [10] Interval-valued IFNs projection Entropy [36] Fermatean fuzzy numbers EDAS CRITIC [37] Single-valued neutrosophic numbers CoCoSo CRITIC [38] Interval Pythagoras hesitant fuzzy numbers AO BWM [18] BCFNs AO CRITIC [39] LDFNs AO Proportion of expectation score of LDFN [40] Crisp numbers and linguistic variables Neighborhood rough set-TOPSIS-VIKOR No [41] Crisp numbers, IFN S and hesitant fuzzy numbers Optimization of the weighted distance measures with ideal solutions Optimization method [42] Crisp numbers, intervals, and linguistic terms CPT Principal component analysis and AHP [25] According to the literature and the meaning of specific criteria, the criteria for 3PRLP evaluation can be expressed in hybrid forms. For example, the cost criterion can be expressed as crisp numbers or as interval numbers, while the capacity criterion can be expressed as linguistic variables, fuzzy numbers, IFNs, or other forms. The concepts of the crisp number, interval number, fuzzy number and linguistic variable are special cases of IFNs, and they can be transformed into IFNs with little information distortion. The concepts, such as the federal fuzzy number, interval-valued IFN, BCFN, neutrosophic number, and LDFN are the extension or improvement of IFN, but compared with them, the membership and non-membership of an IFN make it easier to assign values to the 3PRLP evaluation criteria. In addition, the theoretical basis of intuitionistic fuzzy entropy based on the intuitionistic fuzzy set (IFS) is relatively mature. Therefore, intuitionistic fuzzy MCDM methods are used for evaluating 3PRLPs in this research.

Intuitionistic Fuzzy MCDM
The core idea of the fuzzy set theory proposed by Zadeh [43] in 1965 is to expand the characteristic function whose value is 0 or 1 to the membership function which can take any value in the closed interval [0, 1]. In 1986, Atanassov [44] generalized the fuzzy set theory and proposed the concept of IFS. In addition to membership, IFS also considered the information of non-membership and hesitation, which can more accurately describe uncertainty and fuzziness. Therefore, IFS has been widely used in pattern recognition, medical diagnosis, image processing, MCDM, and other fields [45,46]. For the MCDM problem under the intuitionistic fuzzy environment, scholars have proposed a variety of Systems 2022, 10, 188 6 of 27 models. According to the modeling principle, the models can be divided into the following four categories: (1) models based on AOs, including various AOs and WASPAS (weighted aggregates sum product assessment); (2) models based on criteria preferences, including ELECTRE (elimination et choice translating reality) and PROMETHEE (preference ranking organization method for enrichment evaluation); (3) models based on evidential reasoning (ER); (4) models based on distance, trend, or utility from reference points, including TOPSIS, EADS, MABAC (multi-attributive border approximation area comparison), GRA, MULTIMOORA (multiplicative MOORA), CPT, MARCOS (measurement of alternatives and ranking according to the compromise solution), and VIKOR. Typical intuitionistic fuzzy MCDM models are shown in Table 3.  [47]; neutral averaging AO [48]; geometric Heronian mean AOs [49]; Einstein weighted averaging AO [50]; weighted Heronian mean AO [51]; WASPAS [52] The aggregation of decision information is realized through the weighting operator, and the alternatives are ranked according to the aggregated scores.
Criteria preferences ELECTRE [53][54][55] According to the harmony and disharmony indexes of the criteria set, the preference relationship on the alternative set is constructed, and the alternatives are ranked accordingly.
PROMETHEE [56,57] According to the preference function of each criterion given by the decision-maker, the priority relationship between alternatives and the complete ranking of alternatives are determined.
Evidential reasoning ER [58][59][60] Each criterion is regarded as an evidence, and the ER algorithm is applied to aggregate the basic reliability allocation of each criterion to obtain the comprehensive evaluation value, and then the alternatives are ranked accordingly.
Reference points TOPSIS [61][62][63][64] The alternatives are ranked according to the relative distance between each alternative and the positive and negative ideal points.
EADS [52,65] The evaluation score of alternatives is calculated according to the positive and negative distance between each alternative and the average solution, and the alternatives are ranked accordingly.
MABAC [66] According to the distance between the criterion function of each alternative and the border approximation area, the comprehensive value of each alternative is calculated, and the alternatives are ranked accordingly.
GRA [67,68] The alternatives are ranked according to the trend correlation between each alternative and the ideal reference point.
MULTIMOORA [69] The alternatives are ranked by calculating the additive utility function value (the ratio system to the ideal point) of each alternative.
CPT [70][71][72] According to the value function result relative to the reference point, the cumulative prospect value of each alternative is calculated, and the alternatives are ranked accordingly.
MARCOS [73] According to the utility value relative to the positive and negative ideal points, the comprehensive utility value of each alternative is obtained, and the alternatives are ranked accordingly.
VIKOR [74][75][76][77][78] According to the group utility value and individual regret value relative to the reference points, the benefit ratio value of each alternative is obtained by compromise, and the alternatives are ranked accordingly.
and v A (x) are the membership degree and non-membership degree that x belongs to the IFS A, respectively.
is defined as the score function of x, here s A ∈ [−1,1] reflecting the net membership degree or non-ambiguity degree of x belonging to A. π A (x) = u A (x) + v A (x) is defined as the accuracy function of x, reflecting the degree of accuracy or non-hesitation of x belonging to A.
Definition 4 [79]. For an IFS A, real function E: A→R + is an intuitionistic fuzzy entropy, if it meets the following conditions:

Definition 5. For an IFS
n; x i ∈ X}, an intuitionistic fuzzy entropy can be defined as follows: , so E(A) = E(A c ).
Definition 6 [80]. Let x 1 = u 1 , v 1 and x 2 = u 2 , v 2 be two IFNs, the intuitionistic fuzzy distance of them is as follows:

Transformation of Hybrid Expressions into IFNs
In order to transform crisp numbers, interval numbers, and linguistic variables into IFNs, their values need to be standardized first. Subsequent intuitionistic fuzzy AO involves the multiplication of multiple IFNs. As long as the membership degree of one criterion of an evaluation object equals 1 or the non-membership degree equals 0, the membership degree or the non-membership degree in the comprehensive evaluation result is 0, which is not conducive to the comparison between evaluation objects. Therefore, normalization and vector normalization methods are preferred in the standardization. The value by normalization is usually smaller than that by vector normalization, and the latter is preferred in the case that its standardization result is in the range of [0, 1]. Let the values of n criteria C 1 , C 2 , . . . , C n in m alternatives form the evaluation matrix [x ij ] m×n , where x ij > 0, i = 1, 2, . . . , m; j = 1,2, . . . , n. The value of x ij after intuitionistic fuzzy transformation is [z ij ] m×n . 1 For a crisp number x ij , vector normalization and normalization formula are (3) and (4), respectively, as follows: where C Benefit and C Cost are benefit-type and cost-type criteria sets, respectively. Then transform y ij into the IFN z ij = u ij , v ij = y ij , 1 − y ij .
2 If x ij is an interval number [x ij L , x ij R ], vector normalization and normalization formula are (5) and (6), respectively, as follows.
Then transform the interval number [y ij can be first transformed into standard values d j L and d ij R by the bipolar ratio method, respectively [81]. Taking seven-level linguistic variables as an example, the standard values of the bipolar ratio method are shown in Table 4

Criteria Screening Based on Rough Set
It is difficult for different enterprises to establish a consistent index system because of different objects and requirements for 3PRLP evaluation. In order to fully reflect the personalized requirements of an enterprise's 3PRLP selection and to simplify the criteria, considering the advantages of rough set theory in conditional attribute reduction, we apply a rough set method to screen the criteria. The process depends on the individual experience and judgment of experts in the field. Supposing that an enterprise has established a preliminary evaluation criteria system with t criteria and the corresponding sub-criteria, the organizers consult h experts with high theoretical level and practical experience about 3PRLP, and each expert evaluates the importance of the t criteria and the corresponding sub-criteria according to Likert's five-point scale. The values 1, 2, 3, 4, and 5 represent unimportant, general, important, very important and especially important, respectively. With each criterion as the decision attribute and the corresponding sub-criteria as the condition attributes, t decision tables can be obtained. Taking one decision table as an example, as shown in Table 5, the condition attributes C 1 , C 2 , . . . , C s are the corresponding sub-criteria, the decision attribute D is the criterion, and x i (j) and d i are, respectively, Likert values given by the ith expert on the importance of C j and D relative to 3PRLP evaluation. Table 5. Decision table of one criterion.

Condition Attributes
Decision The steps of criteria screening are as follows [81]: Step 1: According to the decision attribute D, divide the domain U = {1, 2, . . . , h} into q equivalent classes about D: U/D = {H 1 , H 2 , . . . , H q }; Step 2: Calculate the lower approximation of the kth equivalence class H k with respect to the conditional attribute set C = {C 1 , C 2 , . . . , Cs}: Step 3: Remove the attribute C j from C, j = 1, 2, . . . s, and calculate pos(C − C j , D). If pos(C − C j , D) = pos(C, D), it indicates that C j is a redundant attribute. Delete C j from C, and get the reduced set of conditional attributes; Step 4: According to Step 2 and Step 3, test whether there is a redundant attribute in the reduced condition attribute set. Repeat the above steps until all the attributes are non-redundant. Then, we obtain the reduced secondary index set C r corresponding to the primary index.

Criteria Weighting
The weights reflect the importance of criteria and have an important impact on 3PRLP evaluation. There are many methods of criteria weighting, including subjective weighting, objective weighting, and combination weighting. Subjective weighting methods mainly include AHP [82], ANP [82], BWM [83], SWARA [84], PIPRECIA (pivot pairwise relative criteria importance assessment) [85], FUCOM (full consistency method) [86], LBWA (levelbased weight assessment) [87], etc. For the 3PRLP evaluation problem, there are usually some correlations among the criteria of the same level. For example, timeliness of response is not only associated with customer satisfaction, value recovery ratio, and environmental protection effect, which are all the sub-criteria of service quality, but also associated with value-added service capability and network coverage rate, which are both in the sub-criteria of service capability. In the subjective weighting methods, the ANP method can take the relations among the criteria at the same level into consideration, so we apply it for criteria weighting subjectively. Objective weighting methods mainly include CRITIC [8], deviation coefficient, entropy, etc. Considering that intuitionistic fuzzy entropy can better reflect the degree of dispersion among criterion data, and has a mature theoretical basis compared with other concepts, such as the correlation coefficient and standard deviation of IFNs, we apply the intuitionistic fuzzy entropy method for criteria weighting objectively.

Subjective Weighting Based on ANP
Step 1: Construct the typical structure of the ANP. A simple schematic diagram is shown in Figure 1. The first part is the control layer that represents the evaluation based on the network structure criteria system. The second part is the network layer, where each element group or element is not independent, one element may affect any element in the entire network system, and vice versa. In Figure 1, C i represents the ith element group, e ik represents the kth element of C i , and the directed line represents the relationship between elements.  Step 2: Construct a weightless hypermatrix. Taking the element ejl in a certain element group Cj as criterion, in terms of the degree that ejl affects each element eik (k = 1, 2, …, ni) in the element group Ci, we can construct a judgment matrix. After the consistency test, we output the eigenvector ⋯ that satisfies the consistency. Summarizing the eigenvectors of each judgment matrix into a matrix Wij, which reflects the influence relationship between the elements in the element group Ci and those in Cj, and we can obtain the weightless hypermatrix as follows: Step 2: Construct a weightless hypermatrix. Taking the element e jl in a certain element group C j as criterion, in terms of the degree that e jl affects each element e ik (k = 1, 2, . . . , n i ) in the element group C i , we can construct a judgment matrix. After the consistency test, we output the eigenvector W jl i1 W jl i2 · · · W jl in i T that satisfies the consistency. Summarizing the eigenvectors of each judgment matrix into a matrix W ij , which reflects the influence relationship between the elements in the element group C i and those in C j , and we can obtain the weightless hypermatrix as follows: Step 3: Construct weighted hypermatrix. With element group C j as criterion, we conduct pairwise comparison of element groups and construct the judgment matrix as follows:  After the consistency test, we output the eigenvector a 1j a 2j · · · a Nj T that satisfies the consistency. Therefore, the weight matrix A reflecting the relationship between elements is as follows: Multiplying the weightless hypermatrix W S by the weight matrix A, we obtain the weighted hypermatrix W: W = AW S .
Step 4: Calculate the limit hypermatrix. Due to the introduction of feedback and interdependence in ANP, it is relatively complicated to determinate element priority. The two elements can be compared directly or indirectly, and the stable element priority is ensured by solving the limit hypermatrix, as follows: In the limit hypermatrix W l , the value of element in each column is its weight, which reflects the limit relative priority of each element to the corresponding elements in this column. We denote the subjective weights of n sub-criteria by ANP method as η j , j = 1, 2, . . . , n.

Objective Weighting Based on Intuitionistic Fuzzy Entropy
By evaluating n sub-criteria of m 3PRLPs (A 1 , A 2 , . . . , A m ), we obtain the evaluation matrix [x ij ] m×n , and transform it into the intuitionistic fuzzy matrix [z ij ] m×n . The intuitionistic fuzzy entropy of the jth criterion is as follows: Systems 2022, 10, 188

of 27
The objective weight of the jth criterion is as follows: By synthesizing the subjective and objective weights with the weight α and 1 − α, respectively, we can obtain the combination weight of each criterion as follows: w j = αη j + (1 − α)τ j , j = 1, 2, . . . , n.

Evaluation Model Based on Intuitionistic Fuzzy MCDM
There is no strict distinction between the intuitionistic fuzzy MCDM methods in Table 3, and they can be theoretically used for 3PRLP evaluation. Since models based on criteria preferences require more subjective parameter values, we do not consider such methods. Among the models based on AOs, the hybrid weighted averaging operator (HWAO) considers the importance and position of each criterion simultaneously and has some good properties, such as idempotency and boundedness [8], so we apply HWAO for 3PRLP evaluation. Among the three models based on the distance from the reference points, namely TOPSIS, EADS, and MABAC, we choose TOPSIS because of its clear meaning, simple calculation, because it has the most extensive application. Among the models based on the utility from the reference points, namely MULTIMOORA, CPT, MARCOS, and VIKOR, we choose VIKOR because it considers group utility and individual regret simultaneously. In this way, we use five single MCDM methods to evaluate 3PRLPs at the same time, namely HWAO, ER, TOPSIS, GRA, and VIKOR. Because the results of single evaluation methods may not be the same, we need to test their compatibility using Kendall's concordance coefficient method and obtain the compatible models [81]. Since the purpose of 3PRLP evaluation is to determine the reasonable and consistent ranking of multiple 3PRLPs, we apply the ranking-based combination evaluation models, and include Borda count, comprehensive Borda, Copeland, and fuzzy Borda models [88], to combine the results of compatible models. Finally, we output the best combination evaluation results that pass the Spearman consistency test [81]. The evaluation process is shown in Figure 2, which includes the single evaluation, Kendall compatibility test, combination evaluation, and Spearman consistency test. Systems 2022, 10, x FOR PEER REVIEW 13 of 27 Figure 2. The process of 3PRLP evaluation.

HWAO Model
The steps of the HWAO model for 3PRLP evaluation are as follows: Step 1: Calculate the weighted matrix, as follows: = × , where = = 〈 , 〉.
Step 2: Reorder the values of n criteria of each 3PRLP from large to small. Let ( ) be the jth IFN, j = 1, 2, …, n, and the corresponding value before weighting is ziσ(j), and the criterion weight is wiσ(j).
Step 3: Apply a normal distribution method to calculate the position weight ωj, j = 1, 2, …, n, and then the comprehensive value of the ith 3PRLP is as follows: Step 4: Rank m 3PRLPs according to their comprehensive values.

TOPSIS
The steps of TOPSIS model are as follows: Step 1: Determine the positive and negative ideal points as follows: where = 〈max , min 〉 , = 〈min , max 〉.
Step 2: Calculate the distances of each 3PRLP from the positive and negative ideal points as follows: Step 3: Calculate the proximity of each 3PRLP as follows:

HWAO Model
The steps of the HWAO model for 3PRLP evaluation are as follows: Step 1: Calculate the weighted matrix, as follows: Z = z ij m×n , where z ij = w j z ij = w j u ij , v ij .
Step 2: Reorder the values of n criteria of each 3PRLP from large to small. Let z iσ(j) be the jth IFN, j = 1, 2, . . . , n, and the corresponding value before weighting is z iσ(j) , and the criterion weight is w iσ(j) .

TOPSIS
The steps of TOPSIS model are as follows: Step 1: Determine the positive and negative ideal points as follows: where Step 2: Calculate the distances of each 3PRLP from the positive and negative ideal points as follows: Step 3: Calculate the proximity of each 3PRLP as follows: Step 4: Rank m 3PRLPs according to their proximities from large to small.

VIKOR
The steps of the VIKOR model are as follows: Step 1: Calculate the group utility value P i and the individual regret value N i of the ith 3PRLP as follows: Step 2: Calculate the benefit ratio value Q i as follows: where γ is the compromise coefficient between the group utility and individual regret, 0 ≤ γ ≤ 1.
Step 3: Rank m 3PRLPs according to their benefit ratio values from small to large.

GRA
The GRA method evaluates each object based on its relation degree with the reference sequence. Its steps are as follows: Step 1: Take the positive ideal point as reference point, and calculate the relation coefficient between the jth criterion of each 3PRLP and reference point as follows: where ρ is the distinguishing coefficient, ρ ∈ [0, 1].
Step 2: Calculate the relation degree of each 3PRLP as follows: Step 3: Rank m 3PRLPs according to their relation degrees from large to small.

ER
As an uncertain reasoning method, ER takes the index values as the evidence source and the evaluation grades as the identification frame and uses the ER algorithm to fuse the index information to obtain the confidence distribution belonging to each grade and expected utility of each evaluation object [89]. We apply the ER software IDS, developed by Professor Yang of the University of Manchester, for the 3PRLP evaluation [90].
Step 1: Establish an indicator hierarchy structure composed of target and n sub-criteria. Set the best and worst grades for the target, and the utility value is 1 and 0, respectively. According to the rule method, the confidence levels of the best and worst grades of each sub-criterion belonging to the (best, worst) are set as 1, 0 and 0, 1 , respectively.
Step 2: Add m 3PRLPs, and the confidence level of the jth sub-criterion of the ith 3PRLP belonging to the best and worst grades is u ij , v ij .
Step 3: Input the weights of n sub-criteria and evaluate the confidence level u i , v i (i = 1, 2, · · · , so that the target of each 3PRLP belongs to the (best, worst).
Step 4: Rank m 3PRLPs according to their confidence levels.

Compatibility Test
We test the compatibility of single evaluation models using Kendall's concordance coefficient method. Let r ik be the order value of the ith 3PRLP in the kth single model, i = 1, 2, . . . , m; k = 1, 2, . . . , g, and we can calculate the value of statistical indicator as follows: Given the significance level α, if the value of χ 2 is not less than the critical value χ 2 α (m − 1), then the g models are compatible. In the case of incompatibility, we eliminate one single model and calculate the value of the statistical indicator of the remaining models. The set of compatible models with the largest statistical value can be obtained.

Borda Count
By converting the ranking value p ik of A i in the kth single evaluation model into a score, namely f ik = m − p ik + 1, we calculate the average score of each 3PRLP as follows: We then rank m 3PRLPs according to the average scores from large to small. If there are multiple 3PRLPs whose average sores are same, we calculate the variance of scores in different evaluation models and prefer the 3PRLP with a smaller variance.

Comprehensive Borda
Based on the results of each 3PRLP in g compatible single evaluation models, we count the number of models in which the ranking value of A i is less than that of A l and the number of models in which the ranking value of A l is less than that of A i . If the former is larger than the latter, b il = 1; otherwise, b il = 0. Calculate the Borda score of each 3PRLP as follows: and rank m 3PRLPs according to the Borda scores from large to small.

Copeland
The Copeland model is an extension of the Borda model. In addition to defining the value of superior order, it also defines the values of equivalent and inferior orders. If the frequency of A i ranking superior to A l in g compatible single evaluation models is greater than that of A l ranking superior to A i , c il = 1; if the former is equal to the latter, c il = 0; if the former is smaller than the latter, c il = −1. Calculate the Copeland score of each 3PRLP as follows: and rank m 3PRLPs according to the Copeland scores from large to small.

Fuzzy Borda
The fuzzy Borda model takes the difference in scores and in ranking order of various models into consideration, and the calculation steps are as follows: Step 1: Calculate the membership degree of A i that is excellent in the kth model as follows: where, if the kth model is HWAO or ER, f ik is the score of IFN; if it is TOPSIS or GRA, f ik is the value of proximity or relation degree, respectively. If the model is VIKOR, as the benefit ratio is a cost-type index, the following membership function is used: Step 2: Count the fuzzy frequency of A i ranking the hth in g compatible single evaluation models, as follows: where, Then, calculate the fuzzy frequency ratio as follows: Step 3: Calculate the fuzzy Borda score of each 3PRLP as follows: and rank m 3PRLPs according to the fuzzy Borda scores from large to small.

Spearman Consistency Test
The Spearman rank correlation coefficient between the lth combination evaluation model and the kth single evaluation model is as follows: where p il and p ik are the ranking values of A i in the lth combination model and the kth single model, respectively. We calculate the value of statistical indicator as follows: where ρ l is the average correlation coefficient as follows: Given the significance level α, if the value of t l is not less than the critical value t α (m − 2), it means that the lth combination evaluation model is consistent with all the compatible single evaluation models. We output the rank result of the combination evaluation model that satisfies the consistency test and has the maximum value of the statistical indicator.

An Illustrative Example
Company H is a manufacturer of cold chain equipment. In order to realize the recycling and reuse of waste products, the decision-makers of Company H plan to select the best 3PRLP from six 3PRLPs including A 1 , A 2 , . . . , and A 6 . Based on the criteria systems proposed in a variety of sources in the literature, a preliminary evaluation criteria system is constructed, as shown in Table 6. By consulting ten experts with high theoretical level and practical experience, four decision tables are established, which include service quality, service competency, service cost, and cooperative alliance decision tables. For example, the service capacity decision table is shown in Table 7. According to the reduction steps of the rough set method, we obtain 15 reduction sets of conditional attributes as shown in Table 8. By consulting experts, we select the elements in the last one with the longest length as the secondary indexes of service capacity. Similarly, we can obtain the sub-criteria of the other criteria. The final evaluation criteria system for 3PRLP selection of company H is shown in Table 9, where the content in parentheses represents the code for the criterion.

Service quality (B 4 )
Value recovery ratio (C 12 ) Timeliness of response (C 13 ) Customer satisfaction (C 14 ) Environmental protection effect (C 15 ) Through a questionnaire survey on the correlations of evaluation criteria, we obtain the correlation matrix as shown in Table 10, where the number 1 means that the influencing factor on the left affects the affected factor on the top, and the number 0 means that the factor on the left has no impact relationship with the factor at the top. Table 10. Correlation matrix of evaluation criteria.

Affected Factors
Then, we apply the software SuperDecisions (www.superdecisions.com) to draw the network hierarchy, as shown in Figure 3. Then, we apply the software SuperDecisions (www.superdecisions.com) to draw the network hierarchy, as shown in Figure 3. Inputting comparison judgement matrixes of criteria and sub-criteria, we conduct consistency tests and obtain the weightless hypermatrix WS that satisfies consistency as follows: The weighted hypermatrix W is as follows: 0000 0.2237 0.0597 0.1763 0.0773 0.2575 0.0202 0.0630 0.0000 0.0630 0.0000 0.0582 0.1742 0.6125 0.1506  0.1166 0.0000 0.1790 0.0625 0.2318 0.0515 0.0100 0.0223 0.0000 0.0223 0.0000 0.1644 0.0581 0.1225   Inputting comparison judgement matrixes of criteria and sub-criteria, we conduct consistency tests and obtain the weightless hypermatrix W S that satisfies consistency as follows: The weighted hypermatrix W is as follows: The limit hypermatrix W l is: According to the limit hypermatrix, we get the weights of the sub-criteria C 1 , C 2 , . . . ,  Table 11. The criteria C 3 and C 13 are expressed as interval-value fuzzy numbers; C 5 , C 6 , C 10 and C 12 are expressed as general interval numbers; C 14 is expressed as crisp number, and the other five criteria are expressed as linguistic variables with seven-level linguistic term set {highest, very high, high, average, low, very low, lowest}. By standardizing the values of various variables and transforming them into IFNs, we get the intuitionistic fuzzy matrix as follows According to Formulas (7) and (8), we calculate the intuitionistic fuzzy entropy weights of the criteria C 1 , C 2 , . . . , C 15 as 0.0803, 0.0460, 0.0752, 0.0517, 0.0781, 0.0677, 0.0580, 0.0595, 0.0670, 0.0594, 0.0648, 0.0402, 0.0843, 0.0930 and 0.0748, respectively. By synthesizing the subjective and objective weights with the same weight of 0.5, we obtain the combination weights of 15 criteria as 0.1352, 0.0681, 0.0450, 0.0410, 0.0740, 0.0561, 0.0337, 0.0687, 0.0408, 0.0614, 0.0390, 0.0617, 0.0775, 0.1399, and 0.0580, respectively. It can be seen from Figure 4 that, among the 15 criteria, customer satisfaction (C 14 ) and corporate reputation (C 1 ) are the two most important criteria for Company H to select 3PRLP, and the sum of their weights is close to 30%. 0.0595, 0.0670, 0.0594, 0.0648, 0.0402, 0.0843, 0.0930 and 0.0748, respectively. By synthesizing the subjective and objective weights with the same weight of 0.5, we obtain the combination weights of 15 criteria as 0.1352, 0.0681, 0.0450, 0.0410, 0.0740, 0.0561, 0.0337, 0.0687, 0.0408, 0.0614, 0.0390, 0.0617, 0.0775, 0.1399, and 0.0580, respectively. It can be seen from Figure 4 that, among the 15 criteria, customer satisfaction (C14) and corporate reputation (C1) are the two most important criteria for Company H to select 3PRLP, and the sum of their weights is close to 30%. We use the HWAO, TOPSIS, VIKOR, GRA, and ER models to evaluate the six 3PRLPs; their evaluation results are shown in Table 12, where the compromise coefficient in VIKOR and the distinguishing coefficient in GRA are both 0.5. The result of the compatibility test using Kendall's concordance coefficient method is shown in Table 13. The value of Asymp. sig. is 0.0024, which indicates that the evaluation results of five models have a relatively significant consistency, and the five single evaluation models are compatible.  We use the HWAO, TOPSIS, VIKOR, GRA, and ER models to evaluate the six 3PRLPs; their evaluation results are shown in Table 12, where the compromise coefficient in VIKOR and the distinguishing coefficient in GRA are both 0.5. The result of the compatibility test using Kendall's concordance coefficient method is shown in Table 13. The value of Asymp. sig. is 0.0024, which indicates that the evaluation results of five models have a relatively significant consistency, and the five single evaluation models are compatible. By inputting the results of the above single evaluation models into four combination evaluation models, we calculate the results of combination evaluation, which can be seen in Table 14. The ranking results of the Borda count and fuzzy Borda models are identical, and those of the comprehensive Borda and Copeland models are also identical. We calculate the values of statistical indicator in two cases, and the former and the latter are 3.2468 and 2.9598, respectively. They are greater than the critical value t 0.05 (4) = 2.132, and the former is larger, so we output the former result as the final ranking result. Therefore, A 5 ranks the first in the six 3PRLPs, and Company H can select it as the 3PRLP partner. According to the scores of IFNs of 6 3PRLPs' 15 criteria (see Figure 5), A 5 ranks second in the 2 most important criteria C 14 and C 1 , which is the main reason for it ranking first. In addition, A 5 ranks the sixth in cultural and strategic compatibility (C 3 ) and information level (C 7 ), the fifth in communication level (C 4 ), implicit cost (C 6 ), and value recovery ratio (C 12 ), and the fourth in benefit-risk sharing level (C 2 ), which indicate that A 5 has certain gaps in the above six criteria compared with the other 3PRLPs. Therefore, Company H can enhance communication with A 5 to improve the strategic compatibility, and require it to improve its information level, communication level, value recovery ratio, and benefit-risk sharing level, and lower implicit cost. By inputting the results of the above single evaluation models into four combination evaluation models, we calculate the results of combination evaluation, which can be seen in Table 14. The ranking results of the Borda count and fuzzy Borda models are identical, and those of the comprehensive Borda and Copeland models are also identical. We calculate the values of statistical indicator in two cases, and the former and the latter are 3.2468 and 2.9598, respectively. They are greater than the critical value t0.05(4) = 2.132, and the former is larger, so we output the former result as the final ranking result. Therefore, A5 ranks the first in the six 3PRLPs, and Company H can select it as the 3PRLP partner. According to the scores of IFNs of 6 3PRLPs' 15 criteria (see Figure 5), A5 ranks second in the 2 most important criteria C14 and C1, which is the main reason for it ranking first. In addition, A5 ranks the sixth in cultural and strategic compatibility (C3) and information level (C7), the fifth in communication level (C4), implicit cost (C6), and value recovery ratio (C12), and the fourth in benefit-risk sharing level (C2), which indicate that A5 has certain gaps in the above six criteria compared with the other 3PRLPs. Therefore, Company H can enhance communication with A5 to improve the strategic compatibility, and require it to improve its information level, communication level, value recovery ratio, and benefit-risk sharing level, and lower implicit cost.

Conclusions
The evaluation or selection of 3PRLPs plays an important role in helping enterprises to carry out RL effectively and achieve the goal of environmental protection and the associated economic benefits. In view of the personalized requirements of enterprises in the process of selecting 3PRLP, this paper puts forward a criteria screening method based on a rough set. In terms of the interactive relationship of the criteria in the same layer, an ANP-based subjective weighting method is put forward. Combined with the hybrid attribute characteristics of the criteria, this paper puts forward the unified IFN transformation methods, and constructs five single evaluation models based on intuitionistic fuzzy decision making, including HWAO, TOPSIS, VIKOR, GRA, and ER. Kendall's

Conclusions
The evaluation or selection of 3PRLPs plays an important role in helping enterprises to carry out RL effectively and achieve the goal of environmental protection and the associated economic benefits. In view of the personalized requirements of enterprises in the process of selecting 3PRLP, this paper puts forward a criteria screening method based on a rough set. In terms of the interactive relationship of the criteria in the same layer, an ANP-based subjective weighting method is put forward. Combined with the hybrid attribute characteristics of the criteria, this paper puts forward the unified IFN transformation methods, and constructs five single evaluation models based on intuitionistic fuzzy decision making, including HWAO, TOPSIS, VIKOR, GRA, and ER. Kendall's concordance coefficient method is applied for compatibility test on the single evaluation models. Based on the ranking methods, Borda count, comprehensive Borda, Copeland, and fuzzy Borda models are applied for combination evaluation on the results of compatible single evaluation models. Meanwhile, the Spearman rank correlation coefficient method is applied to conduct a consistency test and obtain the best 3PRLP. The results of an illustrative example verify that the criteria screening method and the evaluation models are feasible and effective.
Compared with the existing literature, this paper makes the following contributions: (1) The problem of criteria screening based on a rough set method is discussed, which can meet the personalized requirements of an enterprise and achieve the reduction in criteria; (2) considering the subjective guidance of the important criteria and the objective differences between the evaluation objects, a combined weighting method based on ANP subjective weighting and intuitionistic fuzzy entropy objective weighting is proposed; (3) by combining the principles of different intuitionistic fuzzy MCDM methods and based on the representative models, a combination evaluation idea integrating multiple intuitionistic fuzzy MCDM models is put forward, which enhances the consistency and credibility of decision-making results.
There are some shortcomings in this research, as follows: (1) Only one criteria screening method is provided. In fact, there are other methods, such as the direct scoring and screening of the importance of the criteria by experts. Are these methods compatible? Can they be used together? These questions must be considered; (2) there are many criteria weighting methods. For example, the CRITIC method considers the correlation and difference between criteria at the same time. When a scientific correlation coefficient of IFN S is defined and CRITIC objective weighting is adopted, is it more advantageous than the intuitionistic fuzzy entropy weighting? This must be determined; (3) the single MCDM models can replace each other. For example, VIKOR can be replaced by MARCOS, another utility-based model, while an intuitionistic fuzzy MCDM model can be replaced by its generalized form, such as an interval-valued intuitionistic fuzzy MCDM model. Under these circumstances, how will the evaluation results change? Combined with the above shortcomings, further research should be carried out on the following aspects: (1) In terms of criteria screening, other screening methods can be considered; (2) in terms of criteria weighting, other subjective, objective and combination weighting methods can be adopted; (3) in terms of evaluation models, other single evaluation models can be explored.