1 Introduction

Supply chain management (SCM) has been considered as the most important and critical part of every business worldwide (Koberg & Longoni, 2019; Yazdani et al., 2021). Generally speaking, SCM is placed among the top priorities of every industry due to its high significance as well as the high complexity in processes, functions, and operations (Biswas, 2020). Recently, companies and organizations are not only aiming to address their SCM issues but also aim to consider environmental and social priorities along with economical features in order to transform the supply chain networks (SCNs) into a sustainable-based structure (Yazdani et al., 2021). Although aiming to achieve different targets in terms of modifying traditional SCNs into sustainable supply chains (SSCs) to increase the satisfaction of all stakeholders, companies and organizations soon understood greater challenges that they have to be resilient against all disruptions (Goli et al., 2020; Sangaiah et al., 2020). In this regard, resiliency also got high attention from companies and organizations, which could empower SCNs resilient against disruptions such as flood, storm, transportation risks, and pandemics such as COVID-19 (Elluru et al., 2019; Pavlov et al., 2019).

Since the rise of COVID-19 in China during late December 2019, healthcare supply chains (HSCs) have faced severe issues, risks, and challenges due to sudden increase in demand rate of medical products as a cause of the increase in COVID-19 infected cases (Queiroz et al., 2020). Thus, hospitals and medical centers faced significant shortages in usual products such as special waste bags, medical face masks, face shields, and gloves. According to a sudden increase in the demand rate of such products, the global market faced an unprecedented and unforeseen shock considering the low supply rate and unexpected restrictions in the transportation part of HSCs. Therefore, healthcare centers such as hospitals and medical clinics dealing with COVID-19 patients have to comprehensively consider all resiliency factors since any shortage caused by pandemic-related disruptions in HSCs would lead to irrecoverable consequences. Although resiliency against disruptions caused by the pandemic and its related restrictions play a significant role in HSC, the sustainability of HSCs is also an important matter. Since most of the daily and single-use products in hospitals are made of plastics, an increasing plastic waste generation has notorious impacts on the environment and society. Therefore, one of the critical problems healthcare centers are dealing with is developing a robust structure for a sustainable healthcare supply chain network to facilitate decision-making problems with higher reliability and accuracy.

In this regard, managerial sectors of healthcare centers are concerned about the availability of medical face masks and face shields due to their high demand in healthcare centers dealing with infected patients to the COVID-19 virus. As medical face masks and face shields play critical roles in protecting healthcare staff working in dangerous and infectious areas, they are continuously disposed of due to their single-use feature and short life cycle. Therefore, one of the important and critical decision-making problems of healthcare centers during the COVID-19 pandemic is to address the supply of face masks and face shields from the most decent suppliers. Thus, the supplier selection problem for face masks and face shields is of high significance for healthcare centers dealing with a high number of COVID-19 patients. According to the literature, there exist very limited number of studies which aimed to address supply chain management of primary healthcare products such as face masks (Orji & Ojadi, 2021; Setiawan et al., 2021; Zhu et al., 2021). However, these studies do not consider the uncertainty and subjectivity in information under a reliable environment while addressing the SC problem. On the other hand, none of these studies have not tried to address SCM of face masks and face shields which are accounted as primary and crucial materials during the pandemic.

Generally, sustainable supplier selection is one of the critical aspects of inefficient SCM. There are numerous studies (Milosevic et al., 2021; Petrovic et al., 2019) that argued that supplier performance fundamentally affects the efficiency or inefficiency of the supply chain. That is why the choice of suppliers is recognized as a strategic guideline (Samanta & Giri, 2021). Accordingly, the starting point for achieving sustainability in the SC is the choice of suppliers based on the principle of sustainability (Berisa, 2020; Milenkov et al., 2020). The choice of inadequate suppliers can damage the company's image in the eyes of the public. Therefore, every company or organization should select, identify and collaborate with the most appropriate suppliers to guarantee success in supply chain design (Zavadskas et al., 2020). Along with sustainability, healthcare centers should comprehensively consider the resiliency of suppliers as well. As the healthcare sector is directly dealing with human beings, any disruptions in the supply of such primary products could end up with irrecoverable consequences with high costs.

According to the aforementioned reasons, supplier selection for HSC during COVID-19 is a complex and difficult decision-making problem that is affected by many decision criteria. Based on the uncertain and critical circumstances during the COVID-19 pandemic in healthcare centers, addressing the supplier selection problem has become one of the top priorities of every hospital or medical clinic. In this regard, multi-criteria decision-making (MCDM) methods are reliable and straightforward soft computing-based tools that can empower healthcare managers to address the supplier selection problem considering sustainability and resiliency features (Alosta et al., 2021; Muhammad et al., 2021; Torkayesh et al., 2021; Yazdani et al., 2020; Yazdani et al., 2020). One of the primary conditions that multi-criteria models should fulfill is the possibility of perceiving subjective assessments, unreasonable arguments, and uncertainty in information (Liu et al., 2019; Pelissari et al., 2021; Torkayesh & Deveci, 2021; Torkayesh et al., 2021).

Therefore, in this paper, a new decision-making approach is presented that enables the processing of subjectivity and uncertainty in information during the selection of suppliers under sustainability and resiliency. In order to address the supplier selection problem considering the uncertainty in the information, a multi-criteria framework based on applying fuzzy rough numbers (FRNs) for information processing in the initial decision matrix is presented. The fuzzy rough extension of the CODAS (Combinative Distance-based Assessment) method (Keshavarz Ghorabaee et al., 2016) is used for supplier selection, while the fuzzy rough MACBETH (Measuring Attractiveness by a Categorical Based Evaluation Technique) technique (Bana e Costa and Vansnick, 1994) is used to determine the criteria weights. Although there are many MCDM methods for weighting and ranking purposes, this study adopts MACBETH and CODAS methods. Since MACBETH is one of the very first MCDM methods which has shown high reliability and accuracy to tackle complex MCDM problems; therefore, it is used to determine the weight coefficients of the decision criteria. On the other hand, due to solid structure of the CODAS based on two different distance measures, CODAS is adopted for evaluation of alternative suppliers with higher reliability and preciseness. Utilization of two distance measures increase the reliability of overall score determination which has high significance in prioritization problems. Based on the authors' best knowledge, applying FRN for the extension of CODAS and MACBETH multi-criteria techniques has not been considered in the literature so far. Therefore, FRN-based MACBETH-CODAS is the first integrated multi-criteria methodology based on defining a rough boundary interval using Bonferroni functions, thus enabling additional adaptability of the improved CODAS method. The presented FRN concept enables the definition of the degree of agreement in expert assessments using a rough boundary interval. In the case of complete consensus in expert estimates, the rough boundary interval has zero value, and the FRN is transformed into a classical fuzzy number. However, with the increasing discrepancy in expert estimates, the rough boundary interval increases, leading to an increase in the uncertainty footprint in the fuzzy rough number. The credibility of the presented multi-criteria methodology was presented through a case study in which the evaluation of five suppliers was performed concerning eighteen defined criteria.

Finally, overall contributions of this study can be summarized in development of a nvoel decision-making appraoch under FRNs to address supplier selection problem in HSC during COVID-19 using MACBETH and CODAS methods. This study is the first in its kind to develop FRN-MACBETH-CODAS decision-making approach. This study is also the first utilize the developed decision-making approach to address SS for face masks and face sheilds in the healthcare centers such as the hospitals. In general, this study aims to answer the following research questions:

  1. 1-

    Why HSC management is important during the COVID-19 pandemic?

  2. 2-

    What is the role of supplier selection in HSC management for medical products that are highly consumed during the COVID-19 pandemic?

  3. 3-

    What are the important decision criteria affecting supplier selection problem? How sustainability and resiliency goals can be considered in supplier selection problem?

  4. 4-

    What kind of decision-making approach is reliable and accurate to address supplier selection problem in HSC?

  5. 5-

    How can uncertainty and subjectivity in the information for importance of decision criteria and evaluation of supplier be considered in the decision-making process?

The rest of the paper is organized as follows. Section 2 presents a comprehensive background of the supplier selection problem, HSC management, and uncertain decision-making approaches. FRN MACBETH-CODAS approach and its preliminary definitions are explained in Sect. 3. A real-life case study of supplier selection problem for a hospital in Istanbul, Turkey is presented in Sect. 4. Section 4 also presents information about the results of the developed decision-making approach, extensive sensitivity analyses, and comparative analysis. Finally, we conclude in Sect. 5.

2 Literature review

This section presents a comprehensive background on supplier selection problem, sustainable supplier selection problem, resilient supplier selection problem, and applications of FRN-based decision-making approaches for different strategic decision-making problems. In this regard, Sect. 2.1 provides background of supplier selection problem in SCM with most recent studies. Most important criteria under different categories are identified from the literature for the supplier selection problem considering sustainability criteria. In Sect. 2.2, we discuss studies related to the FRN-based MCDM methods that are used in different decision-making problems.

2.1 Supplier selection in SCM

Supplier selection is one of the critical decision-making problems in SCMs where managers make decisions on choosing the best and most suitable supplier to satisfy the initial raw material, semi-products, or finished products. However, supplier selection is a very difficult, complicated, and multi-dimensional decision-making problem where supplier selection is performed by considering many decision criteria. Apart from the complexity of the supplier selection problem, managers mostly face significant challenges in terms of data availability while comparing the performance of suppliers against qualitative criteria. In order to address such issues, decision-makers and managers adopt linguistic opinion expression sets which empower them to express their opinions on the performance of suppliers against decision criteria with higher reliability and preciseness. In order to address the supplier selection problem in different applications and various environments where there exists uncertainty in data and decision-makers' opinions, many studies are conducted in the literature. In order to review these studies, this subsection is divided into two sections for applying MCDM methods for general supplier selection problem and supplier selection problem in HSC management.

Supplier selection and SCM studies are very rich in the literature recently. Due to the high significance of the decision-making problems in SCM, such as the supplier selection problem, many studies have aimed to tackle such vital decision-making problems in several applications and industries using MCDM methods. MCDM methods empower real-life decision-makers and managers to develop complex decision-making structures to determine the best supplier or sort a series of suppliers based on their performance over several decision criteria. Table 11 in Appendix summarizes the most important and recent studies that used different forms of MCDM methods to address the supplier selection problem.

MCDM methods have also become an important and straightforward decision-making tool to address the supplier selection problem in HSC. Since decision-making problems in HSCs are very critical due to their direct or indirect relationship with human beings and mortality, managers in such organizations show a high willingness to use reliable and straightforward tools to address their daily, weekly, or monthly SC problems such as supplier selection problem. In this regard, MCDM, as one of the potential decision-making tools, has shown high usability for such applications. Ishtiaq et al. (2018) developed AHP method to select the best recycling and supplier provider for healthcare waste in Pakistan. Yazdani et al. (2020) integrated DEMATEL, BWM, and modified EDAS methods to address the supplier selection problem in procurement operations of HSC of a Spanish hospital. Stevic et al. (2020) suggested a new MCDM method called MARCOS to address the sustainable Supplier selection problem for an eye surgical intervention in Bosnia and Herzegovina. Goh and Zhong (2020) used the fuzzy TOPSIS method for logistics service providers for several hospitals to optimize logistics operations of perishable healthcare products. Wei et al. (2019) used Probabilistic Linguistic MABAC to address healthcare supplier selection problems for China's medical consumption products. Bahadori et al. (2020) developed a new decision model using the artificial neural network and VIKOR under a fuzzy environment to address the supplier selection problem in military hospitals considering multi-aspect decision criteria. Gao et al. (2020) developed an improved version of traditional VIKOR method based on the Q-rung Interval-valued orthopair fuzzy information to address the healthcare supplier selection problem for medical consumption products.

2.2 FRN-based MCDM

FRNs are based on the traditional fuzzy logic which Zadeh (1969) introduced. Later, due to the high applicability of fuzzy logic in different industries and applications, new forms of fuzzy logic were introduced. FRNs are relatively considered one of the recent extensions of fuzzy logic that tries to put boundaries on the fuzzy numbers to increase the reliability in the decision-making process using linguistic terms (Dubois & Prade, 1990; Sun et al., 2008). Due to the high efficiency and preciseness of FRNs in decision-making problems, several studies have adopted FRNs and used them in different forms to address various problems. Pamučar et al. (2018b) developed a new modification of BWM and MABAC methods under interval-valued FRNs to deal with uncertainty and vagueness of the input data for the optimal selection of firefighting helicopters. Chen et al. (2020) introduced a new form of the DEMATEL method under FRNs empowered by graphic extensions to assess innovative value propositions for intelligent product-service systems. Mei and Chen (2021) used FRN-based BWM and DEA (data envelopment analysis) to assess occupational risks from human health and environmental aspects. Chen and Ming (2020) used a similar methodology to select the smart product service module. Liu et al. (2021a) utilized BWM and MABAC under a new extension of FRNs called q-rung orthopair. The applicability of the proposed methodology was assessed using a case study in the building sector in China. Zhu et al. (2020) developed a novel decision model based on AHP and TOPSIS under FRNs for design concept evaluation under uncertain and vague information. Zhang et al. (2020) applied a combination of TOPSIS and weighted arithmetic average (WAA) model with FRNs for evaluating a group of customers of a leading bank in China based on their demographic, income and transactional information. Deveci et al. (2020) utilized an interval-valued FRN based Delphi framework for solving the facility location selection problem for wind farms. Jia et al. (2019) extended the classical MABAC model using FRNs for enabling the process of supplier selection for medical devices. An integrated factor relationship (FARE) and multi-attributive border approximation area comparison (MABAC) models with FRNs was applied by Roy et al. (2020) for evaluation and selection of third-party logistics service providers on sustainability dimensions. Sun et al. (2020) proposed a VIKOR focused granular group decision analysis using FRN for evaluation of emergency preparedness. Sun et al. (2020) used multi granulation soft FRN in conjuction with TODIM algorithm in formulating policy decisions. Hou et al. (2021) put forth a FRN based model for identifying suitable indicators for measuring green productivity. In the context of smart city design, Wang et al. (2021) applied FRNs to measure energy efficiency.

2.3 Identification of criteria

With lapse of the time, complexity of decision-making problems has greatly increased. Unlike the past, supplier selection problem is not only affected by technical or economic criteria. With the rise of concepts such as sustainability, resiliency, and circular economy, environmental and social criteria attracted special attention. Since HSC is one of the critical global SCs, the decision-making process has to be built on a comprehensive framework in order to consider impact of environmental, and social aspects along with technical and economic aspects.

According to the literature, the most important and critical decision criteria required for supplier selection problem is presented in Table 1. Table 1 shows a hierarchical framework of the identified criteria under three main criteria, namely as technical, environmental, and social. Along with technical criteria, technical aspect includes economic costs, and resiliency criteria as well. Environmental dimension covers up crucial green and environmental decision criteria related to the manufacturing, packaging and disposal operations. Finally, social criteria aim to take into account various aspects of social concerns. Third column of Table 1 denotes type of each criterion whether it is cost or beneficial.

Table 1 Hierarchy of evaluation criteria

3 FRN MACBETH-CODAS approach

The following section presents an integrated FRN based MACBETH-CODAS model. The FRN MACBETH technique is used to determine the weight coefficients of the criteria, while the FRN CODAS method is used to evaluate the alternatives. FRNs are used to address uncertainties and inaccuracies in expert assessments of alternatives (Fig. 1).

Fig. 1
figure 1

FRN based MACBETH-CODAS model

Let us assume that the multi-criteria model has n criteria Cj = {C1, C2, …, Cn}, (j = 1, 2, …, n) which were used to evaluate l alternatives Ai = {A1, A2, …, Al}, (i = 1,2, …, l). Also, suppose that this is a group decision making model in which experts Eh = {E1, E2, …, En}, (h = 1, 2, …, e) evaluate alternatives Ai (i = 1, 2, …, l) in relation to the criteria Cj (j = 1, 2,., n) using a predefined fuzzy scale. Then, to evaluate alternatives using the FRN CODAS model, it is necessary to determine the weighting coefficients of the criteria. The following section of the paper presents the mathematical formulation of the FRN MACBETH-CODAS approach.

3.1 Preliminaries of FRNs

In this study, a new concept of hybrid fuzzy rough numbers is presented based on the adaptive definition of the rough boundary interval of a fuzzy number. The adaptability of the rough boundary interval is enabled by applying the Bonferroni function (Bonferroni, 1950). The application of the Bonferroni function enables flexible decision-making due to decision makers' risk attitudes. In addition, the Bonferroni function was applied as it allows the interrelationships between fuzzy sequences to be appreciated. One of the significant advantages of the fuzzy rough number concept is that the fuzzy rough number (FRN) limit values are defined based on the internal inaccuracies that exist in the original data set. The FRN concept is based on the settings of conventional rough numbers. As with classical rough numbers (RN), the lower limit, upper limit, and rough boundary interval are defined for FRN. The difference in relation to classical RNs is that FRNs enable the membership degree definition for each element in the fuzzy rough set. On the other hand, in classical RNs, each component of the rough boundary interval entirely belongs to the considered collection.

Assuming that \({\mathbb{N}}\) denotes a universe that contains all the decision-maker (DM) preferences that participate in decision making and that triangular fuzzy numbers represent these preferences \(\widetilde{\tau }_{i} = \left( {\tau_{i}^{(l)} ,\tau_{i}^{(m)} ,\tau_{i}^{(u)} } \right)\), with the mode, left endpoint, and right endpoint denoted by \(\tau_{i}^{(l)}\), \(\tau_{i}^{(m)}\), and \(\tau_{i}^{(u)}\), respectively. Then DMs preferences can be divided into x classes that satisfy the condition that \(\widetilde{\tau }_{1} \le \widetilde{\tau }_{2} \le ,..., \le \widetilde{\tau }_{x}\). If we assume that \(\widetilde{\Lambda }\) is a collection of \(\left( {\widetilde{\tau }_{1} ,\widetilde{\tau }_{2} ,...,\widetilde{\tau }_{x} } \right)\) and \(\vartheta\) is an arbitrary element of \({\mathbb{N}}\), then the lower approximation of class \(\widetilde{\tau }_{i}\) can be defined as follows:

$$\underline{Apr} \left( {\tau_{i}^{(l)} } \right) = \mathop \cup \limits_{1 \le i \le x} \left\{ {\vartheta \in {\mathbb{N}}/\widetilde{\Lambda }(\vartheta ) \le \tau_{i}^{(l)} } \right\}$$
(1)
$$\underline{Apr} \left( {\tau_{i}^{(m)} } \right) = \mathop \cup \limits_{1 \le i \le x} \left\{ {\vartheta \in {\mathbb{N}}/\widetilde{\Lambda }(\vartheta ) \le \tau_{i}^{(m)} } \right\}$$
(2)
$$\underline{Apr} \left( {\tau_{i}^{(u)} } \right) = \mathop \cup \limits_{1 \le i \le x} \left\{ {\vartheta \in {\mathbb{N}}/\widetilde{\Lambda }(\vartheta ) \le \tau_{i}^{(u)} } \right\}$$
(3)

Similarly, we can define the upper approximation as follows:

$$\overline{Apr} \left( {\tau_{i}^{(l)} } \right) = \mathop \cup \limits_{1 \le i \le x} \left\{ {\vartheta \in {\mathbb{N}}/\widetilde{\Lambda }(\vartheta ) \ge \tau_{i}^{(l)} } \right\}$$
(4)
$$\overline{Apr} \left( {\tau_{i}^{(m)} } \right) = \mathop \cup \limits_{1 \le i \le x} \left\{ {\vartheta \in {\mathbb{N}}/\widetilde{\Lambda }(\vartheta ) \ge \tau_{i}^{(m)} } \right\}$$
(5)
$$\overline{Apr} \left( {\tau_{i}^{(u)} } \right) = \mathop \cup \limits_{1 \le i \le x} \left\{ {\vartheta \in {\mathbb{N}}/\widetilde{\Lambda }(\vartheta ) \ge \tau_{i}^{(u)} } \right\}$$
(6)

Based on the defined lower and upper approximations of the FRN, we can represent the lower limit of \(\widetilde{\tau }_{i}\) as follows:

$$\underline{{Lim}} \left( {\tau _{i}^{{(l)}} } \right) = \left( {\frac{1}{{N_{{Ll}} }}\sum\limits_{{i,j = 1}}^{{N_{{Ll}} }} {\tau _{i}^{{(l)\delta _{1} }} } \left( {\prod\limits_{{j = 1}}^{{N_{{Ll}} }} {\tau _{j}^{{(l)\delta _{2} }} } } \right)^{{\frac{1}{{N_{{Ll}} - 1}}}} } \right)^{{\frac{1}{{\delta _{1} + \delta _{2} }}}} \left| {\tau _{i}^{{(l)\delta _{1} }} ,\tau _{j}^{{(l)\delta _{2} }} } \right. \in \underline{{Apr}} \left( {\tau _{i}^{{(l)}} } \right)$$
(7)
$$\underline{{Lim}} \left( {\tau _{i}^{{(m)}} } \right) = \left( {\frac{1}{{N_{{Lm}} }}\sum\limits_{{i,j = 1}}^{{N_{{Lm}} }} {\tau _{i}^{{(m)\delta _{1} }} } \left( {\prod\limits_{{j = 1}}^{{N_{{Lm}} }} {\tau _{j}^{{(m)\delta _{2} }} } } \right)^{{\frac{1}{{N_{{Lm}} - 1}}}} } \right)^{{\frac{1}{{\delta _{1} + \delta _{2} }}}} \left| {\tau _{i}^{{(m)\delta _{1} }} ,\tau _{i}^{{(m)\delta _{2} }} } \right. \in \underline{{Apr}} \left( {\tau _{i}^{{(m)}} } \right)$$
(8)
$$\underline{{Lim}} \left( {\tau _{i}^{{(u)}} } \right) = \left( {\frac{1}{{N_{{Lu}} }}\sum\limits_{{i,j = 1}}^{{N_{{Lu}} }} {\tau _{i}^{{(u)\delta _{1} }} } \left( {\prod\limits_{{j = 1}}^{{N_{{Lu}} }} {\tau _{j}^{{(u)\delta _{2} }} } } \right)^{{\frac{1}{{N_{{Lu}} - 1}}}} } \right)^{{\frac{1}{{\delta _{1} + \delta _{2} }}}} \left| {\tau _{i}^{{(u)\delta _{1} }} ,\tau _{i}^{{(u)\delta _{2} }} } \right. \in \underline{{Apr}} \left( {\tau _{i}^{{(u)}} } \right)$$
(9)

where \(N_{Ll}\), \(N_{Lm}\) and \(N_{Lu}\) represents a number of elements in \(\underline{Apr} (\tau_{i}^{(l)} )\), \(\underline{Apr} (\tau_{i}^{(m)} )\) and \(\underline{Apr} (\tau_{i}^{(u)} )\) respectively; \(\delta_{1} ,\delta_{2} \ge 0\) and \(\delta_{1} ,\delta_{2} \in \Re\), where \(\Re\) represents a set of real numbers.

Similarly, we can define an upper limit of \(\widetilde{\tau }_{i}\) as follows:

$$\overline{{Lim}} \left( {\tau _{i}^{{(l)}} } \right) = \left( {\frac{1}{{N_{{Ul}} }}\sum\limits_{{i,j = 1}}^{{N_{{Ul}} }} {\tau _{i}^{{(l)\delta _{1} }} } \left( {\prod\limits_{{j = 1}}^{{N_{{Ul}} }} {\tau _{j}^{{(l)\delta _{2} }} } } \right)^{{\frac{1}{{N_{{Ul}} - 1}}}} } \right)^{{\frac{1}{{\delta _{1} + \delta _{2} }}}} \left| {\tau _{i}^{{(l)\delta _{1} }} ,\tau _{j}^{{(l)\delta _{2} }} } \right. \in \overline{{Apr}} \left( {\tau _{i}^{{(l)}} } \right)$$
(10)
$$\overline{{Lim}} \left( {\tau _{i}^{{(m)}} } \right) = \left( {\frac{1}{{N_{{Um}} }}\sum\limits_{{i,j = 1}}^{{N_{{Um}} }} {\tau _{i}^{{(m)\delta _{1} }} } \left( {\prod\limits_{{j = 1}}^{{N_{{Um}} }} {\tau _{j}^{{(m)\delta _{2} }} } } \right)^{{\frac{1}{{N_{{Um}} - 1}}}} } \right)^{{\frac{1}{{\delta _{1} + \delta _{2} }}}} \left| {\tau _{i}^{{(m)\delta _{1} }} ,\tau _{j}^{{(m)\delta _{2} }} } \right. \in \overline{{Apr}} \left( {\tau _{i}^{{(m)}} } \right)$$
(11)
$$\overline{{Lim}} \left( {\tau _{i}^{{(u)}} } \right) = \left( {\frac{1}{{N_{{Uu}} }}\sum\limits_{{i,j = 1}}^{{N_{{Uu}} }} {\tau _{i}^{{(u)\delta _{1} }} } \left( {\prod\limits_{{j = 1}}^{{N_{{Uu}} }} {\tau _{j}^{{(u)\delta _{2} }} } } \right)^{{\frac{1}{{N_{{Uu}} - 1}}}} } \right)^{{\frac{1}{{\delta _{1} + \delta _{2} }}}} \left| {\tau _{i}^{{(u)\delta _{1} }} ,\tau _{j}^{{(u)\delta _{2} }} } \right. \in \overline{{Apr}} \left( {\tau _{i}^{{(u)}} } \right)$$
(12)

where \(N_{Ul}\), \(N_{Um}\) and \(N_{Uu}\) represents a number of elements in \(\overline{Apr} (\tau_{i}^{(l)} )\), \(\overline{Apr} (\tau_{i}^{(m)} )\) and \(\overline{Apr} (\tau_{i}^{(u)} )\) respectively. Then, based on the previously defined (expression (1)–(12) we can represent the FRN \(\widetilde{\tau }_{i}\) as follows:

$$\begin{aligned}FRN\left( {{{\widetilde \tau }_i}} \right) = & \left( \left[ {\underline {Lim} \left( {\tau _i^{(l)}} \right),\overline {Lim} \left( {\tau _i^{(l)}} \right)} \right],\left[ \underline {Lim} \left( {\tau _i^{(m)}} \right),\overline {Lim} \left( {\tau _i^{(m)}} \right) \right],\right.\\ &\quad\left.\left[ {\underline {Lim} \left( {\tau _i^{(u)}} \right),\overline {Lim} \left( {\tau _i^{(u)}} \right)} \right] \right)\\ = & \left( {\left[ {\tau _i^{(l) - },\tau _i^{(l) + }} \right],\left[ {\tau _i^{(m) - },\tau _i^{(m) + }} \right],\left[ {\tau _i^{(u) - },\tau _i^{(u) + }} \right]} \right)\end{aligned}$$
(13)

Rough boundary interval (RBI) of \(\tau_{i}^{(l)}\), \(\tau_{i}^{(m)}\), and \(\tau_{i}^{(u)}\) we can define it as follows:

$$\begin{gathered} RNBnd\left( {\tau_{i}^{(l)} } \right) = \underline{Lim} \left( {\tau_{i}^{(l)} } \right) - \overline{Lim} \left( {\tau_{i}^{(l)} } \right) \hfill \\ RNBnd\left( {\tau_{i}^{(m)} } \right) = \underline{Lim} \left( {\tau_{i}^{(m)} } \right) - \overline{Lim} \left( {\tau_{i}^{(m)} } \right) \hfill \\ RNBnd\left( {\tau_{i}^{(u)} } \right) = \underline{Lim} \left( {\tau_{i}^{(u)} } \right) - \overline{Lim} \left( {\tau_{i}^{(u)} } \right) \hfill \\ \end{gathered}$$
(14)

The RBI represents the degree of uncertainty that exists in the FRN. A higher RBI value represents greater uncertainty in the information, while a smaller interval represents less uncertainty (Pamucar et al., 2018; Pamucar et al., 2018). The value of \(RNBnd = 0\) represents the absence of uncertainty in the information, so then the FRN is transformed into a triangular fuzzy number.

In the following section, operations between two FRNs are defined \(FRN\left( {\widetilde{\tau }_{1} } \right) = \left( {\left[ {\tau_{1}^{(l) - } ,\tau_{1}^{(l) + } } \right],\left[ {\tau_{1}^{(m) - } ,\tau_{1}^{(m) + } } \right],\left[ {\tau_{1}^{(u) - } ,\tau_{1}^{(u) + } } \right]} \right)\) and \(FRN\left( {\widetilde{\tau }_{2} } \right) = \left( {\left[ {\tau_{2}^{(l) - } ,\tau_{2}^{(l) + } } \right],\left[ {\tau_{2}^{(m) - } ,\tau_{2}^{(m) + } } \right],\left[ {\tau_{2}^{(u) - } ,\tau_{2}^{(u) + } } \right]} \right)\), \(\theta > 0\):

$$\begin{aligned}FRN\left( {\widetilde{\tau }_{1} } \right) + FRN\left( {\widetilde{\tau }_{2} } \right) = &\left[ \left( {\tau_{1}^{(l) - } + \tau_{2}^{(l) - } ,\tau_{1}^{(l) + } + \tau_{2}^{(l) + } } \right),\left( {\tau_{1}^{(m) - } + \tau_{2}^{(m) - } ,\tau_{1}^{(m) + } + \tau_{2}^{(m) + } } \right),\right.\\&\left.\left( {\tau_{1}^{(u) - } + \tau_{2}^{(u) - } ,\tau_{1}^{(u) + } + \tau_{2}^{(u) + } } \right) \right]\end{aligned}$$
(15)
$$\begin{aligned}FRN\left( {\widetilde{\tau }_{1} } \right) \times FRN\left( {\widetilde{\tau }_{2} } \right) = &\left[ \left( {\tau_{1}^{(l) - } \times \tau_{2}^{(l) - } ,\tau_{1}^{(l) + } \times \tau_{2}^{(l) + } } \right),\left( {\tau_{1}^{(m) - } \times \tau_{2}^{(m) - } ,\tau_{1}^{(m) + } \times \tau_{2}^{(m) + } } \right),\right.\\&\left.\left( {\tau_{1}^{(u) - } \times \tau_{2}^{(u) - } ,\tau_{1}^{(u) + } \times \tau_{2}^{(u) + } } \right) \right]\end{aligned}$$
(16)
$$\begin{aligned}\theta \times FRN\left( {\widetilde{\tau }_{2} } \right) =& \left[ \left( {\theta \times \tau_{2}^{(l) - } ,\theta \times \tau_{2}^{(l) + } } \right),\left( {\theta \times \tau_{2}^{(m) - } ,\theta \times \tau_{2}^{(m) + } } \right),\right.\\&\left.\left( {\theta \times \tau_{2}^{(u) - } ,\theta \times \tau_{2}^{(u) + } } \right) \right] \end{aligned}$$
(17)

3.2 FRN-MACBETH

MACBETH is a reliable tool for determining the criteria weights. It is based on the application of a linear mathematical model and the methodology of pairwise comparisons of the criteria (Bana e Costa et al., 2001, 2004). Methodologically, the MACBETH method algorithm is reminiscent of the Analytic Hierarchy Process (AHP) method algorithm (Saaty, 1980), however, there are also significant differences between the two methodologies which are highlighted in Table 2.

Table 2 Characteristics MACBETH and AHP

The Crisp MACBETH methodology has been successfully used to determine the weighting coefficients of criteria and to evaluate alternatives in numerous MCDM studies (Bana e Costa and Vansnick, 1994), career management (Bana e Costa et al., 2004), energy storage (Montignac et al., 2009), healthcare management (Rodrigues, 2014), textile industry (Kundakcı & Tuş Işk, 2016), supply chain management (Ghosh & Biswas, 2016), off-road vehicle selection (Pamucar & Savin, 2020), and steam boiler systems (Kundakcı, 2018). However, according to the authors' best knowledge, the extension of the MACBETH method using RFNs has not been shown in the literature so far. This was one of the primary motives of the authors for the application of the FRN MACBETH method.

The weighting coefficients of the criteria are determined using the FRN MACBETH linear model, which is presented in the following section:

Step 1: Criteria from the set Cj = {C1, C2, …, Cn} are ranked according to their significance, so that the most influential criterion is in the first position, while the least influential criterion is in the last position.

Step 2: Forming a comparison matrix \({\mathbb{Z}} = \left[ {\widetilde{\varsigma }_{ij} } \right]_{n \times n}\), where \(\widetilde{\varsigma }_{ij} = \left( {\varsigma_{ij}^{(l)} ,\varsigma_{ij}^{(s)} ,\varsigma_{ij}^{(u)} } \right)\) represents the fuzzy value from the matrix \({\mathbb{Z}} = \left[ {\widetilde{\varsigma }_{ij} } \right]_{n \times n}\), which is defined based on the value from the fuzzy semantic scale (Table 3).

Table 3 Fuzzy semantic scale (Ecer and Pamucar, 2020)

In the matrix \({\mathbb{Z}}\), the comparison of criteria is performed only above the diagonal of the matrix, while the values of zero are entered on the diagonal of the matrix. Since e experts participate in the research, for each expert, we get a comparison matrix \({\mathbb{Z}}^{b} = \left[ {\widetilde{\varsigma }_{ij}^{b} } \right]_{n \times n}\), \(1 \le b \le e\).

Step 3: Formation of a fuzzy MACBETH linear model. A fuzzy linear model is formed based on values obtained in expert matrices \({\mathbb{Z}}^{b} = \left[ {\widetilde{\varsigma }_{ij}^{b} } \right]_{n \times n}\), \(1 \le b \le e\) considering the following rules:

  1. (1)

    Defining objective function fuzzy linear model, expression (18):

    $$\Im_{\min } = \wp (\widetilde{\chi }_{1} )$$
    (18)

    where \(\widetilde{\chi }_{1}\) represents the fuzzy value of the most influential criterion.

  2. (2)

    The first group of constraints is represented by the ordinal constraints, expression (19):

    $$\forall \widetilde{\chi }_{i} ,\widetilde{\chi }_{j} ,i,j \in \left\{ {1,2,..,n} \right\}:\widetilde{\chi }_{i} > \widetilde{\chi }_{j} \Rightarrow \wp (\widetilde{\chi }_{i} ) \ge \wp (\widetilde{\chi }_{j} ) + \widetilde{\varsigma }(i,j)$$
    (19)

    where \(\widetilde{\xi }(i,j)\) represents preference level difference between \(\widetilde{\delta }_{i} = \left( {\delta_{i}^{l} ,\delta_{i}^{s} ,\delta_{i}^{u} } \right)\) and \(\widetilde{\delta }_{j} = \left( {\delta_{j}^{l} ,\delta_{j}^{s} ,\delta_{j}^{u} } \right)\).

  3. (3)

    The second group of constraints is represented by semantic constraints, expression (20):

    $$\forall \widetilde{\chi }_{i} ,\widetilde{\chi }_{j} ,i,j,k,l \in \left\{ {1,2,..,n} \right\}:\wp (\widetilde{\chi }_{i} ) - \wp (\widetilde{\chi }_{j} ) \ge \wp (\widetilde{\chi }_{k} ) - \wp (\widetilde{\chi }_{l} ) + \widetilde{\varsigma }(i,j,k,l)$$
    (20)

    where the value of \(\widetilde{\chi }(i,j,k,l)\) is defined as the difference between \(\widetilde{\varsigma }(i,j)\) and \(\widetilde{\varsigma }(k,l)\).

  4. (4)

    The last constraint of the linear programming model is shown by the expression (21):

    $$\wp (\widetilde{\chi }_{n} ) = (1,1,2)$$
    (21)

    where \(\widetilde{\chi }_{n}\) represents the value assigned to the least influential criterion.

Step 4: By solving the fuzzy linear model, we obtain the fuzzy values of the weight coefficients of the criteria for each expert, i.e., for each criterion, we obtain the weight coefficient \(\widetilde{w}_{j}^{h} = \left( {w_{j}^{(l)h} ,w_{j}^{(m)h} ,w_{j}^{(u)h} } \right)\), (h = 1, 2, …, e). Using expression (1)–(2) fuzzy sequences \(w_{j}^{h(l)} = \left\{ {w_{j}^{1(l)} ,w_{j}^{2(l)} ,...,w_{j}^{e(l)} } \right\}\), \(w_{j}^{h(m)} = \left\{ {w_{j}^{1(m)} ,w_{j}^{2(m)} ,...,w_{j}^{e(m)} } \right\}\) and \(w_{j}^{h(m)} = \left\{ {w_{j}^{1(u)} ,w_{j}^{2(u)} ,...,w_{j}^{e(u)} } \right\}\) (h = 1,2, …, e) we transform into rough sequences and form a fuzzy rough weight coefficient \(\overline{w}_{j}^{h} = \left( {\left[ {\overline{w}_{j}^{h(l) - } ,\overline{w}_{j}^{h(l) + } } \right],\left[ {\overline{w}_{j}^{h(m) - } ,\overline{w}_{j}^{h(m) + } } \right],\left[ {\overline{w}_{j}^{h(u) - } ,\overline{w}_{j}^{h(u) + } } \right]} \right)\) (h = 1, 2, …, e). For FRN fusion of values \(\overline{w}_{j}^{h}\) (h = 1,2,, e) the FRN weighted geometric Bonferroni function was used, expression (22).

$$\begin{aligned} \overline{w}_{j} = &\left( {\frac{1}{e(e - 1)}\sum\limits_{\begin{subarray}{l} x,y = 1 \\ x \ne y \end{subarray} }^{e} {\overline{w}_{j(x)}^{{\delta_{1} }} \otimes \overline{w}_{j(y)}^{{\delta_{2} }} } } \right)^{{\frac{1}{{\delta_{1} + \delta_{2} }}}} \hfill \\ = &\left( \begin{gathered} \left[ {\left( {\frac{1}{e(e - 1)}\sum\limits_{\begin{subarray}{l} x,x = 1 \\ x \ne y \end{subarray} }^{e} {w_{j(x)}^{{(l) - }{\delta_{1} }} \otimes w_{j(y)}^{{(l) - }{\delta_{2} }} } } \right)^{{\frac{1}{{\delta_{1} + \delta_{2} }}}} ,\left( {\frac{1}{e(e - 1)}\sum\limits_{\begin{subarray}{l} x,x = 1 \\ x \ne y \end{subarray} }^{e} {w_{j(x)}^{{(l) + }{\delta_{1} }} \otimes w_{j(y)}^{{(l) + }{\delta_{2} }} } } \right)^{{\frac{1}{{\delta_{1} + \delta_{2} }}}} } \right], \hfill \\ \left[ {\left( {\frac{1}{e(e - 1)}\sum\limits_{\begin{subarray}{l} x,x = 1 \\ x \ne y \end{subarray} }^{e} {w_{j(x)}^{{(m) - }{\delta_{1} }} \otimes w_{j(y)}^{{(m) - }{\delta_{2} }} } } \right)^{{\frac{1}{{\delta_{1} + \delta_{2} }}}} ,\left( {\frac{1}{e(e - 1)}\sum\limits_{\begin{subarray}{l} x,x = 1 \\ x \ne y \end{subarray} }^{e} {w_{j(x)}^{{(m) + }{\delta_{1} }} \otimes w_{j(y)}^{{(m) + }{\delta_{2} }} } } \right)^{{\frac{1}{{\delta_{1} + \delta_{2} }}}} } \right], \hfill \\ \left[ {\left( {\frac{1}{e(e - 1)}\sum\limits_{\begin{subarray}{l} x,x = 1 \\ x \ne y \end{subarray} }^{e} {w_{j(x)}^{{(u) - }{\delta_{1} }} \otimes w_{j(y)}^{{(u) - }{\delta_{2} }} } } \right)^{{\frac{1}{{\delta_{1} + \delta_{2} }}}} ,\left( {\frac{1}{e(e - 1)}\sum\limits_{\begin{subarray}{l} x,x = 1 \\ x \ne y \end{subarray} }^{e} {w_{j(x)}^{{(u) + }{\delta_{1} }} \otimes w_{j(y)}^{{(u) + }{\delta_{2} }} } } \right)^{{\frac{1}{{\delta_{1} + \delta_{2} }}}} } \right] \hfill \\ \end{gathered} \right) \hfill \\ \end{aligned}$$
(22)

where e represents the number of experts participating in the research, while δ1, δ2 ≥ 0 are set of non-negative numbers.

Fuzzy rough sequences of weighting coefficients are defined for each expert. Expression (22) was used to aggregate experts' fuzzy rough sequences. The Bonferroni function was chosen to aggregate expert sequences since it enables the flexible fusion of weight coefficients and efficient representation of mutual relations between criteria. Also, the variation of the parameters δ1 and δ2 provides the possibility to check the robustness of the results and to check the influence of different levels of risk on the final decision.

3.3 FRN-CODAS

This section proposes an FRN extension of the CODAS methodology to deal with uncertain domain-based MCDM problems. The CODAS method is an efficient and updated decision-making methodology introduced by Keshavarz Ghorabaee et al. (2016). The CODAS technique applies Euclidean distance as the primary measure of estimating alternatives. If the Euclidean distances of the two alternatives are very close to each other, then the Taxicab distance is used to compare them. The degree of proximity of Euclidean distances is defined based on a threshold parameter. The use of two types of distances in the evaluation process contributes to increasing the accuracy of ranking results. Euclidean and Taxicab distances are measures for indifference spaces for l1-norm and l2-norm (Hassanpour & Pamucar, 2019). Due to its high reliability, different forms of CODAS have been used in various applications such as waste management (Karagoz et al., 2020), electric vehicles (Deveci et al., 2021; Simic et al., 2021), and transportation management (Seker & Aydin, 2020).

In the CODAS method, we first evaluate the alternatives in the indifference space for the l2-norm. If the alternatives are not comparable in this space, then the space of indifference to the l1-norm is considered. To perform this procedure, each pair of alternatives should be compared. The CODAS technique has significant advantages, which include: (1) a simple mathematical apparatus that is programmable and easy to use; and (2) Increasing the number of criteria/alternatives in the CODAS model does not complicate the evaluation algorithm.

For example, suppose we have l alternatives and n criteria, then we can present the algorithmic steps of the modified FRN CODAS method as follows:

Step 1. Construct the basic FRN decision-matrix (\(\overline{\overline{\aleph }}\)). Evaluation of the alternatives per each criterion by e expert is denoted as \(\widetilde{\theta }_{ij}^{b}\): \(1 \le b \le e\); i = 1,…,n; j = 1,…,l. The value of each pair \(\widetilde{\theta }_{ij}^{b}\) takes a fuzzy value from the predefined fuzzy scale. The judgment of e expert is presented as a matrix \(\overline{\overline{\aleph }}^{b} = \left[ {\widetilde{\theta }_{ij}^{b} } \right]_{l \times n}\), where \(1 \le b \le e\); i = 1,…,n; j = 1,…,l; and \(\widetilde{\theta }_{ij}^{b} = \left( {\theta_{ij}^{b(l)} ,\theta_{ij}^{b(m)} ,\theta_{ij}^{b(u)} } \right)\) represent linguistic variable taken from the preliminary defined linguistic scale used by expert e.

Based on response matrices \(\overline{\overline{\aleph }}^{b} = \left[ {\widetilde{\theta }_{ij}^{b} } \right]_{l \times n}\)(\(1 \le b \le e\)) obtained from each e expert, for each element \(\theta_{ij}^{e(l)}\), \(\theta_{ij}^{e(m)}\) and \(\theta_{ij}^{e(u)}\) of fuzzy number \(\widetilde{\theta }_{ij}^{b} = \left( {\theta_{ij}^{b(l)} ,\theta_{ij}^{b(m)} ,\theta_{ij}^{b(u)} } \right)\) we form matrices of the aggregated sequences of experts \(\overline{\overline{\aleph }}^{b(l)}\), \(\overline{\overline{\aleph }}^{b(m)}\) and \(\overline{\overline{\aleph }}^{b(u)}\)

$$\overline{\overline{\aleph }}^{b(l)} = \left[ {\begin{array}{*{20}c} {\theta_{11}^{1(l)} ,\theta_{11}^{2(l)} , \ldots ,\theta_{11}^{e(l)} ,\,\,\,\,\,\,\,\,\theta_{12}^{e(l)} ,\theta_{12}^{2(l)} , \ldots ,\theta_{12}^{e(l)} ,\,\,\,\,\,\, \ldots ,\,\,\,\,\,\,\,\theta_{1n}^{e(l)} ,\theta_{1n}^{2(l)} , \ldots ,\theta_{1n}^{e(l)} } \\ {\theta_{21}^{1(l)} ,\theta_{21}^{2(l)} , \ldots ,\theta_{21}^{e(l)} ,\,\,\,\,\,\,\,\,\theta_{22}^{e(l)} ,\theta_{22}^{2(l)} , \ldots ,\theta_{22}^{e(l)} ,\,\,\,\,\,\, \ldots ,\,\,\,\,\,\,\,\theta_{2n}^{e(l)} ,\theta_{2n}^{2(l)} , \ldots ,\theta_{2n}^{e(l)} } \\ { \ldots \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, \ldots \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, \ldots ,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, \ldots } \\ {\theta_{l1}^{1(l)} ,\theta_{l1}^{2(l)} , \ldots ,\theta_{l1}^{e(l)} ,\,\,\,\,\,\,\,\,\theta_{l2}^{e(l)} ,\theta_{l2}^{2(l)} , \ldots ,\theta_{l2}^{e(l)} ,\,\,\,\,\,\, \ldots ,\,\,\,\,\,\,\,\theta_{\ln }^{e(l)} ,\theta_{\ln }^{2(l)} , \ldots ,\theta_{\ln }^{e(l)} } \\ \end{array} } \right]$$
(23)
$$\overline{\overline{\aleph }}^{b(m)} = \left[ {\begin{array}{*{20}c} {\theta_{11}^{1(m)} ,\theta_{11}^{2(m)} , \ldots ,\theta_{11}^{e(m)} ,\,\,\,\,\,\,\,\,\theta_{12}^{e(m)} ,\theta_{12}^{2(m)} , \ldots ,\theta_{12}^{e(m)} ,\,\,\,\,\,\, \ldots ,\,\,\,\,\,\,\,\theta_{1n}^{e(m)} ,\theta_{1n}^{2(m)} , \ldots ,\theta_{1n}^{e(m)} } \\ {\theta_{21}^{1(m)} ,\theta_{21}^{2(m)} , \ldots ,\theta_{21}^{e(m)} ,\,\,\,\,\,\,\,\,\theta_{22}^{e(m)} ,\theta_{22}^{2(m)} , \ldots ,\theta_{22}^{e(m)} ,\,\,\,\,\,\, \ldots ,\,\,\,\,\,\,\,\theta_{2n}^{e(m)} ,\theta_{2n}^{2(m)} , \ldots ,\theta_{2n}^{e(m)} } \\ { \ldots \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, \ldots \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, \ldots ,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, \ldots } \\ {\theta_{l1}^{1(m)} ,\theta_{l1}^{2(m)} , \ldots ,\theta_{l1}^{e(m)} ,\,\,\,\,\,\,\,\,\theta_{l2}^{e(m)} ,\theta_{l2}^{2(m)} , \ldots ,\theta_{l2}^{e(m)} ,\,\,\,\,\,\, \ldots ,\,\,\,\,\,\,\,\theta_{\ln }^{e(m)} ,\theta_{ln}^{2(m)} , \ldots ,\theta_{ln}^{e(m)} } \\ \end{array} } \right]$$
(24)
$$\overline{\overline{\aleph }}^{b(u)} = \left[ {\begin{array}{*{20}c} {\theta_{11}^{1(u)} ,\theta_{11}^{2(u)} , \ldots ,\theta_{11}^{e(u)} ,\,\,\,\,\,\,\,\,\theta_{12}^{e(u)} ,\theta_{12}^{2(u)} , \ldots ,\theta_{12}^{e(u)} ,\,\,\,\,\,\, \ldots ,\,\,\,\,\,\,\,\theta_{1n}^{e(u)} ,\theta_{1n}^{2(u)} , \ldots ,\theta_{1n}^{e(u)} } \\ {\theta_{21}^{1(u)} ,\theta_{21}^{2(u)} , \ldots ,\theta_{21}^{e(u)} ,\,\,\,\,\,\,\,\,\theta_{22}^{e(u)} ,\theta_{22}^{2(u)} , \ldots ,\theta_{22}^{e(u)} ,\,\,\,\,\,\, \ldots ,\,\,\,\,\,\,\,\theta_{2n}^{e(u)} ,\theta_{2n}^{2(u)} , \ldots ,\theta_{2n}^{e(u)} } \\ { \ldots \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, \ldots \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, \ldots ,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, \ldots } \\ {\theta_{l1}^{1(u)} ,\theta_{l1}^{2(u)} , \ldots ,\theta_{l1}^{e(u)} ,\,\,\,\,\,\,\,\,\theta_{l2}^{e(u)} ,\theta_{l2}^{2(u)} , \ldots ,\theta_{l2}^{e(u)} ,\,\,\,\,\,\, \ldots ,\,\,\,\,\,\,\,\theta_{\ln }^{e(u)} ,\theta_{ln}^{2(u)} , \ldots ,\theta_{ln}^{e(u)} } \\ \end{array} } \right]$$
(25)

where \(\theta_{ij}^{e(l)}\), \(\theta_{ij}^{e(m)}\) and \(\theta_{ij}^{e(u)}\) represent sequences of fuzzy number \(\widetilde{\theta }_{ij}^{b} = \left( {\theta_{ij}^{b(l)} ,\theta_{ij}^{b(m)} ,\theta_{ij}^{b(u)} } \right)\).

Using expressions (1)–(2) each sequence \(\theta_{ij}^{e(l)}\), \(\theta_{ij}^{e(m)}\) and \(\theta_{ij}^{e(u)}\) is transformed into a rough sequence that together makes up FRN \(\overline{\theta }_{ij}^{b} = \left( {\left[ {\overline{\theta }_{ij}^{b(l) - } ,\overline{\theta }_{ij}^{b(l) + } } \right],\left[ {\overline{\theta }_{ij}^{b(m) - } ,\overline{\theta }_{ij}^{b(m) + } } \right],\left[ {\overline{\theta }_{ij}^{b(u) - } ,\overline{\theta }_{ij}^{b(u) + } } \right]} \right)\); \(1 \le b \le e\). For FRN fusion of values \(\overline{\theta }_{ij}^{b}\) (\(1 \le b \le e\)) the FRN weighted geometric Bonferroni function was used, expression (26):

$$\begin{aligned} \overline{\theta }_{ij} = & \left( {\frac{1}{e(e - 1)}\sum\limits_{\begin{subarray}{l} x,y = 1 \\ x \ne y \end{subarray} }^{e} {\overline{\theta }_{ij(x)}^{{\alpha_{1} }} \otimes \overline{\theta }_{ij(y)}^{{\alpha_{2} }} } } \right)^{{\frac{1}{{\alpha_{1} + \alpha_{2} }}}} \\ = & \left( \begin{gathered} \left[ {\left( {\frac{1}{e(e - 1)}\sum\limits_{\begin{subarray}{l} x,y = 1 \\ x \ne y \end{subarray} }^{e} {\theta_{ij(x)}^{{(l) - \alpha_{1} }} \otimes \theta_{ij(y)}^{{(l) - \alpha_{2} }} } } \right)^{{\frac{1}{{\alpha_{1} + \alpha_{2} }}}} ,\left( {\frac{1}{e(e - 1)}\sum\limits_{\begin{subarray}{l} x,y = 1 \\ x \ne y \end{subarray} }^{e} {\theta_{ij(x)}^{{(l) + \alpha_{1} }} \otimes \theta_{ij(y)}^{{(l) + \alpha_{2} }} } } \right)^{{\frac{1}{{\alpha_{1} + \alpha_{2} }}}} } \right], \hfill \\ \left[ {\left( {\frac{1}{e(e - 1)}\sum\limits_{\begin{subarray}{l} x,y = 1 \\ x \ne y \end{subarray} }^{e} {\theta_{ij(x)}^{{(m) - \alpha_{1} }} \otimes \theta_{ij(y)}^{{(m) - \alpha_{2} }} } } \right)^{{\frac{1}{{\alpha_{1} + \alpha_{2} }}}} ,\left( {\frac{1}{e(e - 1)}\sum\limits_{\begin{subarray}{l} x,y = 1 \\ x \ne y \end{subarray} }^{e} {\theta_{ij(x)}^{{(m) + \alpha_{1} }} \otimes \theta_{ij(y)}^{{(m) + \alpha_{2} }} } } \right)^{{\frac{1}{{\alpha_{1} + \alpha_{2} }}}} } \right], \hfill \\ \left[ {\left( {\frac{1}{e(e - 1)}\sum\limits_{\begin{subarray}{l} x,y = 1 \\ x \ne y \end{subarray} }^{e} {\theta_{ij(x)}^{{(u) - \alpha_{1} }} \otimes \theta_{ij(y)}^{{(u) - \alpha_{2} }} } } \right)^{{\frac{1}{{\alpha_{1} + \alpha_{2} }}}} ,\left( {\frac{1}{e(e - 1)}\sum\limits_{\begin{subarray}{l} x,y = 1 \\ x \ne y \end{subarray} }^{e} {\theta_{ij(x)}^{{(u) + \alpha_{1} }} \otimes \theta_{ij(y)}^{{(u) + \alpha_{2} }} } } \right)^{{\frac{1}{{\alpha_{1} + \alpha_{2} }}}} } \right] \hfill \\ \end{gathered} \right) \\ \end{aligned}$$
(26)

We thus obtain the basic FRN matrix of decision \(\overline{\overline{\aleph }} = \left[ {\overline{\theta }_{ij} } \right]_{l \times n}\); \(\overline{\theta }_{ij} = \left( {\left[ {\overline{\theta }_{ij}^{(l) - } ,\overline{\theta }_{ij}^{(l) + } } \right],\left[ {\overline{\theta }_{ij}^{(m) - } ,\overline{\theta }_{ij}^{(m) + } } \right],\left[ {\overline{\theta }_{ij}^{(u) - } ,\overline{\theta }_{ij}^{(u) + } } \right]} \right)\); i = 1,…,l; j = 1,…,n.

Step 2. Normalization of the basic matrix element (\(\overline{\overline{\aleph }}\)). Determine FRN normalized decision matrix \(\overline{\overline{\aleph }}^{N} = \left[ {\overline{\varphi }_{ij} } \right]_{l \times n}\) according to the type of each criterion using the following equations

$$\overline{\varphi }_{ij} = \left\{ \begin{gathered} \frac{{\overline{\theta }_{ij} }}{{\widehat{\theta }}}{ = }\left[ {\left( {\frac{{\overline{\theta }_{ij}^{(l) - } }}{{\widehat{\theta }}},\frac{{\overline{\theta }_{ij}^{(l) + } }}{{\widehat{\theta }}}} \right),\left( {\frac{{\overline{\theta }_{ij}^{(m) - } }}{{\widehat{\theta }}},\frac{{\overline{\theta }_{ij}^{(m) + } }}{{\widehat{\theta }}}} \right),\left( {\frac{{\overline{\theta }_{ij}^{(u) - } }}{{\widehat{\theta }}},\frac{{\overline{\theta }_{ij}^{(u) + } }}{{\widehat{\theta }}}} \right)} \right]{; }if \, j \in B, \hfill \\ \frac{{\widehat{{\widehat{\theta }}}}}{{\overline{\theta }_{ij} }}{ = }\left[ {\left( {\frac{{\widehat{{\widehat{\theta }}}}}{{\overline{\theta }_{ij}^{(u) + } }},\frac{{\widehat{{\widehat{\theta }}}}}{{\overline{\theta }_{ij}^{(u) - } }}} \right),\left( {\frac{{\widehat{{\widehat{\theta }}}}}{{\overline{\theta }_{ij}^{(m) + } }},\frac{{\widehat{{\widehat{\theta }}}}}{{\overline{\theta }_{ij}^{(m) - } }}} \right),\left( {\frac{{\widehat{{\widehat{\theta }}}}}{{\overline{\theta }_{ij}^{(l) + } }},\frac{{\widehat{{\widehat{\theta }}}}}{{\overline{\theta }_{ij}^{(l) - } }}} \right)} \right]{; }if \, j \in C \hfill \\ \end{gathered} \right.$$
(27)

where \(\widehat{\theta } = \mathop {\max }\limits_{1 \le i \le l} \left( {\overline{\theta }_{ij} } \right)\), \(\widehat{{\widehat{\theta }}} = \mathop {\min }\limits_{1 \le i \le l} \left( {\overline{\theta }_{ij} } \right)\), B and C represent the sets of benefit and cost criteria.

Step 3. Calculate FRN weighted normalized decision matrix (\(\overline{\overline{\Omega }}\)). The FRN weighted normalized matrix \(\overline{\overline{\Omega }} = [\overline{\tau }_{ij} ]_{l \times n}\) is calculated as follows

$$\begin{aligned}\overline{\tau }_{ij} = &\overline{w}_{j} \cdot \overline{\varphi }_{ij} = \left( {\left[ {w_{j}^{(l) - } ,w_{j}^{(l) + } } \right],\left[ {w_{j}^{(m) - } ,w_{j}^{(m) + } } \right],\left[ {w_{j}^{(u) - } ,w_{j}^{(u) + } } \right]} \right)\\ &\cdot \left( \left[ {\overline{\varphi }_{ij}^{(l) - } ,\overline{\varphi }_{ij}^{(l) + } } \right],\left[ {\overline{\varphi }_{ij}^{(m) - } ,\overline{\varphi }_{ij}^{(m) + } } \right],\left[ {\overline{\varphi }_{ij}^{(u) - } ,\overline{\varphi }_{ij}^{(u) + } } \right] \right)\end{aligned}$$
(28)

where \(\overline{w}_{j} = \left( {\left[ {w_{j}^{(l) - } ,w_{j}^{(l) + } } \right],\left[ {w_{j}^{(m) - } ,w_{j}^{(m) + } } \right],\left[ {w_{j}^{(u) - } ,w_{j}^{(u) + } } \right]} \right)\) represents the final FRN criteria weight of jth criterion.

Step 4. Determine FR negative-ideal solution. We obtain the FRN negative-ideal solution matrix \(\overline{\overline{NS}} = [\overline{\eta }_{j} ]_{1 \times n}\) as follows

$$\begin{aligned}\overline{\eta }_{j} =& \mathop {\min }\limits_{1 \le i \le l} \left\{ {\overline{\tau }_{ij} } \right\} = \left[ \left( {\mathop {\min }\limits_{1 \le i \le l} \left\{ {\overline{\tau }_{ij}^{(l) - } } \right\},\mathop {\min }\limits_{1 \le i \le l} \left\{ {\overline{\tau }_{ij}^{(l) + } } \right\}} \right),\right.\\&\left.\left( {\mathop {\min }\limits_{1 \le i \le l} \left\{ {\overline{\tau }_{ij}^{(m) - } } \right\},\mathop {\min }\limits_{1 \le i \le l} \left\{ {\overline{\tau }_{ij}^{(m) + } } \right\}} \right),\left( {\mathop {\min }\limits_{1 \le i \le l} \left\{ {\overline{\tau }_{ij}^{(u) - } } \right\},\mathop {\min }\limits_{1 \le i \le l} \left\{ {\overline{\tau }_{ij}^{(u) + } } \right\}} \right) \right] \end{aligned}$$
(29)

Step 5. Calculate the FRN weighted Euclidean (\(\overline{\overline{ED}}_{i}\)), and FRN weighted Hamming (\(\overline{\overline{HD}}_{i}\)) distances. We obtain \(\overline{\overline{ED}}_{i}\) and \(\overline{\overline{HD}}_{i}\) as follows:

  1. (1)

    FRN weighted Euclidean (\(\overline{\overline{ED}}_{i}\)) distances:

    $$\overline{\overline{ED}}_{i} = \sum\limits_{j = 1}^{n} {d_{E} \left( {\overline{\tau }_{ij} ;\overline{\eta }_{j} } \right)}$$
    (30)

    where \(d_{E} \left( {\overline{\tau }_{ij} ;\overline{\eta }_{j} } \right)\) we obtain as follows

    $$\begin{aligned} &d_{E} \left( {\overline{\tau }_{ij} ;\overline{\eta }_{j} } \right)\\ &\quad = \sqrt {\frac{{\left\{ {\overline{\tau }_{ij}^{(l) - } - \overline{\eta }_{j}^{(l) - } } \right\}^{2} + \left\{ {\overline{\tau }_{ij}^{(l) + } - \overline{\eta }_{j}^{(l) + } } \right\}^{2} + 2\left\{ {\overline{\tau }_{ij}^{(m) - } - \overline{\eta }_{j}^{(m) - } } \right\}^{2} + 2\left\{ {\overline{\tau }_{ij}^{(m) + } - \overline{\eta }_{j}^{(m) + } } \right\}^{2} + \left\{ {\overline{\tau }_{ij}^{(u) - } - \overline{\eta }_{j}^{(u) - } } \right\}^{2} + \left\{ {\overline{\tau }_{ij}^{(u) + } - \overline{\eta }_{j}^{(u) + } } \right\}^{2} }}{8}}\end{aligned}$$
    (31)
  2. (2)

    FRN weighted Hamming (\(\overline{\overline{HD}}_{i}\)) distances

    $$\overline{\overline{HD}}_{i} = \sum\limits_{j = 1}^{n} {d_{H} \left( {\overline{\tau }_{ij} ;\overline{\eta }_{j} } \right)}$$
    (32)

    where \(d_{H} \left( {\overline{\tau }_{ij} ;\overline{\eta }_{j} } \right)\) we obtain as follows

    $$\begin{aligned} &d_{H} \left( {\overline{\tau }_{ij} ;\overline{\eta }_{j} } \right) \\ &\quad= \frac{{\left| {\overline{\tau }_{ij}^{(l) - } - \overline{\eta }_{j}^{(l) - } } \right| + \left| {\overline{\tau }_{ij}^{(l) + } - \overline{\eta }_{j}^{(l) + } } \right| + 2\left| {\overline{\tau }_{ij}^{(m) - } - \overline{\eta }_{j}^{(m) - } } \right| + 2\left| {\overline{\tau }_{ij}^{(m) + } - \overline{\eta }_{j}^{(m) + } } \right| + \left| {\overline{\tau }_{ij}^{(u) - } - \overline{\eta }_{j}^{(u) - } } \right| + \left| {\overline{\tau }_{ij}^{(u) + } - \overline{\eta }_{j}^{(u) + } } \right|}}{8}\end{aligned}$$
    (33)

Step 6. Determine relative assessment matrix. By applying expression (30) we obtain elements of the relative assessment matrix \({\mathbb{Q}} = [\phi_{ik} ]_{l \times l}\)

$$\phi_{ik} = \left( {\overline{\overline{ED}}_{i} - \overline{\overline{ED}}_{k} } \right) + \left( {\beta \left( {\overline{\overline{ED}}_{i} - \overline{\overline{ED}}_{k} } \right) \times \left( {\overline{\overline{HD}}_{i} - \overline{\overline{HD}}_{k} } \right)} \right)$$
(34)

where \(k \in \{ 1,2, \ldots ,l\}\) and β is a threshold function that is defined as follows:

$$\beta (\varepsilon ) = \left\{ \begin{gathered} 1 \, if \, \left| \varepsilon \right| \ge \zeta \hfill \\ 0 \, if \, \left| \varepsilon \right| < \zeta \hfill \\ \end{gathered} \right.$$
(35)

The decision-maker can set the threshold parameter (\(\zeta\)) of this function. In this study, we use \(\zeta = 0.02\) for the calculations.

Step 7. Calculate the assessment score of each alternative. By applying expression (32) we obtain an assessment score

$$\Psi_{i} = \sum\limits_{k = 1}^{l} {\phi_{ik} }$$
(36)

The alternative with the highest assessment score is the most desirable alternative.

4 Case study

This section presents information on background of the case study, potential suppliers, and related numerical results on weight determination, supplier evaluation, sensitivity analyses, and comparative analysis.

In order to show the applicability of the developed novel decision model for real-life practices in HSC, a case study is investigated for a private hospital in European side of Istanbul, Turkey. This hospital is called as ABC hereafter due to confidentiality reasons. Due to its large infrastructure and many medical departments, ABC has been always one of the busiest hospitals in Istanbul. With the rise of COVID-19, the demand for medical services increased significantly as the ABC has had to take care of COVID-19 infected patients and patients diagnosed with other diseases. In this regard, ABC faced an unprecedented challenge in terms of single use medical products such as medical face mask and face shields for its working staff as well as its clients. Due to high demand of these products and the necessity of using them continuously, the problem brought up many challenges and risks to supply them in the right time, right amount, and of course in right price. Thus, one of the crucial decision-making problems for the ABC is to select the best supplier for face masks and face shields. In this regard, five experts from the managerial departments of the ABC are invited to participate in the following study.

Based on the interview with managerial roles in the ABC hospitals, 5 experts were invited to participate in the identification of potential suppliers and data collection. The expert team includes 3 male mangers and 2 female mangers of the hospital who has more than 10-year experience in HSC operations for different healthcare organizations and centers. Based on the background of SC operations for face mask and shields, the expert team decided to consider five potential companies as the suppliers. Based on the data privacy agreement and information disclosure of the hospital and the companies, these suppliers are called as (A1, A2, A3, A4, A5). All five suppliers are domestic and local Turkish companies that are producing different medical products including face masks, face shields, medical pads, and other primary products for healthcare centers. Three suppliers A2, A3, A5 have supplied different types of material to the hospital over recent decades. However, suppliers A1, and A4 are relatively new suppliers that almost emerged in era of the COVID-19 pandemic. All suppliers are located in within territory of Istanbul province except for A4 which is based in the Bursa province.

4.1 Results of FRN-MACBETH

This section presents information on the results of weighting calculations to determine weight coefficient of main criteria and sub-criteria using the FRN-MACBETH. In order to make the calculation straightforward, all processes and steps are itemized below.

Step 1: In the first step, the criteria were ranked according to expert preferences. For supplier evaluation, 18 criteria were defined, which were grouped into three clusters: Technical criterion (EC1), Environmental criterion (EN1), and social criterion (SO1) (Table 1).

The study involved five experts who ranked the criteria according to their preferences:

Expert 1: EC1 > SO1 > EN1; EC16 > EC13 > EC15 > EC17 > EC14 > EC18 > EC11 > EC12; EN15 > EN12 > EN11 > EN13 > EN14; SO14 > SO15 > SO11 > SO12 > SO13.

Expert 2: EC2 > SO1 > EN1; EC16 > EC13 > EC15 > EC12 > EC11 > EC17 > EC18 > EC14; EN12 > EN13 > EN15 > EN14 > EN11; SO13 > SO12 > SO14 > SO15 > SO11.

Expert 3: EC3 > SO1 > EN1; EC15 > EC12 > EC13 > EC16 > EC11 > EC14 > EC17 > EC18; EN11 > EN12 > EN14 > EN13 > EN15; SO15 > SO13 > SO14 > SO12 > SO11.

Expert 4: EC4 > SO1 > EN1; EC15 > EC17 > EC13 > EC16 > EC12 > EC14 > EC18 > EC11; EN14 > EN11 > EN12 > EN15 > EN13; SO14 > SO12 > SO15 > SO13 > SO11.

Expert 5: EC5 > SO1 > EN1; EC15 > EC17 > EC13 > EC16 > EC12 > EC14 > EC18 > EC11; EN14 > EN11 > EN12 > EN15 > EN13; SO14 > SO12 > SO15 > SO13 > SO11.

Step 2: Based on the defined criteria priorities, the experts made pairwise comparisons within each group/cluster. The comparison in pairs was performed using the fuzzy semantic scale shown in Table 3. The following section presents the expert matrices of cluster comparison, Table 4.

Table 4 Pairwise comparison matrices for clusters

Meanwhile, matrices for the remaining criteria were formed in the same way.

Step 3: Based on the pairwise comparison matrices of criteria and expressions (18)–(21), constraints for MACBETH linear models were formed. Four linear models were developed for each expert, i.e., one linear model for each criterion level. The next part shows a linear model created based on the preferences of Expert 1 for clusters level.

$$\begin{aligned} Expert{ 1} \\ { - - - - - - - - - - - } \\ \Im_{\min } = & \wp (\widetilde{\chi }_{SO1} ) \\ & \wp \left( {\widetilde{\chi }_{SO1} } \right) \ge \wp \left( {\widetilde{\chi }_{EN1} } \right) + VW;\wp \left( {\widetilde{\chi }_{SO1} } \right) \ge \wp \left( {\widetilde{\chi }_{EC1} } \right) + M; \\ & \wp \left( {\widetilde{\chi }_{EN1} } \right) \ge \wp \left( {\widetilde{\chi }_{EC1} } \right) + M; \\ & \wp \left( {\widetilde{\chi }_{SO1} } \right) - \wp \left( {\widetilde{\chi }_{EN1} } \right) \ge \wp \left( {\widetilde{\chi }_{EN1} } \right) - \wp \left( {\widetilde{\chi }_{EC1} } \right) + W - VW; \, \\ & \wp \left( {\widetilde{\chi }_{SO1} } \right) - \wp \left( {\widetilde{\chi }_{EC1} } \right) \ge \wp \left( {\widetilde{\chi }_{EN1} } \right) - \wp \left( {\widetilde{\chi }_{EC1} } \right) + M - VW; \\ & \wp \left( {\widetilde{\chi }_{EN1} } \right) - \wp \left( {\widetilde{\chi }_{EC1} } \right) \ge \wp \left( {\widetilde{\chi }_{EN1} } \right) - \wp \left( {\widetilde{\chi }_{EC1} } \right) + M - VW; \\ & \wp \left( {\widetilde{\chi }_{EC1} } \right) = (1,1,2); \, \wp \left( {\widetilde{\chi }_{j} } \right) \ge 0,\forall j = 1,2,3 \\ \end{aligned}$$

Similarly, MACBETH linear models were created for the remaining criteria from the considered set of criteria. Lingo 17.0 software on a personal computer was used to solve linear models.

Step 4: By solving MACBETH linear models, fuzzy values of weight coefficients are determined as in the Table 5.

Table 5 Fuzzy criteria weights

Step 5: By applying expressions (1)–(12), the fuzzy criteria weights were transformed into FRN weighting coefficients. The fusion of weights was performed using the FRN weighted geometric Bonferroni function, expression (22). The final values of the FRN weighting coefficients of the criteria are shown in Table 6.

Table 6 FRN weighting coefficients of the criteria

Local values of FRN weighting coefficients were obtained using the FRN weighted geometric Bonferroni function, expression (22). Figure 2 shows the local FRN values of the weighting coefficients of the cluster level. Global values of weighting coefficients were obtained by multiplying the weighting coefficients of the cluster by the values of the weighting coefficients of the criteria within the corresponding cluster. Finally, the global values of the FRN weighting coefficients of the criteria were used to evaluate the alternatives in the FRN CODAS model. Criteria are ranked based on their global weights in Table 6. Training about green practices for stakeholders (SO14), ensuring rights of the stakeholders (SO15), and job creation (SO11) are top three criteria with highest weight coefficients for sustainable supplier selection problem. Robustness (EC18), quality (EC14), and flexibility (SO17) are three criteria with the least weight coefficients for the sustainable supplier selection problem.

Fig. 2
figure 2

Local FRN values of cluster weighting coefficients

Based on the FRN weighting values (Table 6 and Fig. 2), we notice that the FRN weight coefficients have higher or lower values of the rough boundary interval. Based on that, we can conclude that there are different expert preferences when comparing the criteria. For weighting coefficients where there is a high consensus of experts when comparing the criteria, the length of the rough boundary interval is small. As the consensus of experts decreases, so does the length of the rough boundary interval. In the case of complete consensus in expert preferences, the FRN weighting coefficient is translated into a fuzzy value.

4.2 Results of FRN-CODAS

Similar to the previous sub-section, results of the ranking of suppliers considering all sub-criteria are presented here. For appropriate ranking of suppliers based on the experts' opinions, FRN-CODAS method is used. In other words, after calculating the weight coefficients of the criteria, expert evaluation of the alternatives \(A_{i} \left( {i = 1,2, \ldots ,5} \right)\) was carried out with the eighteen predefined evaluation criteria. For the ease of follow up, the calculations are itemized as different steps below.

Step 1: Construct the basic FRN decision-matrix (\(\overline{\overline{\aleph }}\)). To apply the proposed FRN CODAS multi-criteria model for evaluating the alternatives \(A_{i} \left( {i = 1,2, \ldots ,5} \right)\), the experts' Eh (h = 1,2,…,5) evaluated alternatives using the fuzzy scale shown in Table 7.

Table 7 Fuzzy scale for evaluating the alternatives (Tesfamariam & Sadiq, 2006)

The evaluation of alternatives was performed with the eighteen criteria shown in Table 1. Expert assessments of alternatives are shown in Table 8.

Table 8 Expert evaluation of the alternatives

Using expressions (1)–(12), the fuzzy sequences from Table 8 were transformed into FRN sequences. The FRN weighted geometric Bonferroni function, expression (26), was used for FRN sequence fusion. An example of the transformation of expert fuzzy estimates from Table 8 at position A1-EC12 is shown in the following section.

At position A1- EC12 in Table 8, the fuzzy values \(\widetilde{\theta }_{12}^{(1)} = \widetilde{\theta }_{12}^{(4)} = \widetilde{\theta }_{12}^{(5)} = \left( {7,9,10} \right)\) and \(\widetilde{\theta }_{12}^{(2)} = \widetilde{\theta }_{12}^{(3)} = \left( {5,7,9} \right)\) were obtained. Based on the presented fuzzy values, three sequences \(\theta_{12}^{(l)} = \left\{ {5,5,7,7,7} \right\}\), \(\theta_{12}^{(m)} = \left\{ {7,7,9,9,9} \right\}\) and \(\theta_{12}^{(u)} = \left\{ {9,9,10,10,10} \right\}\) were formed. Using the expression (1)–(12 and provided that δ1 = δ2 = 1, we can define the lower and upper limit sequences according to the following:

  1. (a)

    Lower limits:

$$\underline{Lim} (\theta_{12}^{(1)(l)} ) = \underline{Lim} (\theta_{12}^{(4)(l)} ) = \underline{Lim} (\theta_{12}^{(5)(l)} ) = \left( {\frac{1}{5}\left\{ {5 \cdot \left( {5 \cdot 7 \cdot 7 \cdot 5} \right)^{1/4} + 5 \cdot \left( {5 \cdot 7 \cdot 7 \cdot 5} \right)^{1/4} + \cdots + 7 \cdot \left( {5 \cdot 5 \cdot 7 \cdot 5} \right)^{1/4} } \right\}} \right)^{1/2} = 6.14;$$
$$\underline{Lim} (\theta_{12}^{(2)(l)} ) = \underline{Lim} (\theta_{12}^{(3)(l)} ) = 5;$$
$$\underline{Lim} (\theta_{12}^{(1)(m)} ) = \underline{Lim} (\theta_{12}^{(4)(m)} ) = \underline{Lim} (\theta_{12}^{(5)(m)} ) = \left( {\frac{1}{5}\left\{ {7 \cdot \left( {7 \cdot 9 \cdot 9 \cdot 7} \right)^{1/4} + 7 \cdot \left( {7 \cdot 9 \cdot 9 \cdot 7} \right)^{1/4} + \cdots + 9 \cdot \left( {7 \cdot 7 \cdot 9 \cdot 7} \right)^{1/4} } \right\}} \right)^{1/2} = 8.16;$$
$$\underline{Lim} (\theta_{12}^{(2)(m)} ) = \underline{Lim} (\theta_{12}^{(3)(m)} ) = 7;$$
$$\underline{Lim} (\theta_{12}^{(1)(u)} ) = \underline{Lim} (\theta_{12}^{(4)(u)} ) = \underline{Lim} (\theta_{12}^{(5)(u)} ) = \left( {\frac{1}{5}\left\{ {9 \cdot \left( {9 \cdot 10 \cdot 10 \cdot 9} \right)^{1/4} + 9 \cdot \left( {9 \cdot 10 \cdot 10 \cdot 9} \right)^{1/4} + \cdots + 10 \cdot \left( {9 \cdot 9 \cdot 10 \cdot 9} \right)^{1/4} } \right\}} \right)^{1/2} = 9.59;$$
$$\underline{Lim} (\theta_{12}^{(2)(u)} ) = \underline{Lim} (\theta_{12}^{(3)(u)} ) = 9.$$
  1. (b)

    Upper limits:

$$\overline{Lim} (\theta_{12}^{(1)(l)} ) = \overline{Lim} (\theta_{12}^{(4)(l)} ) = \overline{Lim} (\theta_{12}^{(5)(l)} ) = 7;$$
$$\overline{Lim} (\theta_{12}^{(2)(l)} ) = \overline{Lim} (\theta_{12}^{(3)(l)} ) = \left( {\frac{1}{5}\left\{ {5 \cdot \left( {5 \cdot 7 \cdot 7 \cdot 5} \right)^{1/4} + 5 \cdot \left( {5 \cdot 7 \cdot 7 \cdot 5} \right)^{1/4} + \cdots + 7 \cdot \left( {5 \cdot 5 \cdot 7 \cdot 5} \right)^{1/4} } \right\}} \right)^{1/2} = 6.14;$$
$$\overline{Lim} (\theta_{12}^{(1)(m)} ) = \overline{Lim} (\theta_{12}^{(4)(m)} ) = \overline{Lim} (\theta_{12}^{(5)(m)} ) = 9;$$
$$\overline{Lim} (\theta_{12}^{(2)(m)} ) = \overline{Lim} (\theta_{12}^{(3)(m)} ) = \left( {\frac{1}{5}\left\{ {7 \cdot \left( {7 \cdot 9 \cdot 9 \cdot 7} \right)^{1/4} + 7 \cdot \left( {7 \cdot 9 \cdot 9 \cdot 7} \right)^{1/4} + \cdots + 9 \cdot \left( {7 \cdot 7 \cdot 9 \cdot 7} \right)^{1/4} } \right\}} \right)^{1/2} = 8.16;$$
$$\overline{Lim} (\theta_{12}^{(1)(u)} ) = \overline{Lim} (\theta_{12}^{(4)(u)} ) = \overline{Lim} (\theta_{12}^{(5)(u)} ) = 10$$
$$\overline{Lim} (\theta_{12}^{(2)(u)} ) = \overline{Lim} (\theta_{12}^{(3)(u)} ) = \left( {\frac{1}{5}\left\{ {9 \cdot \left( {9 \cdot 10 \cdot 10 \cdot 9} \right)^{1/4} + 9 \cdot \left( {9 \cdot 10 \cdot 10 \cdot 9} \right)^{1/4} + \cdots + 10 \cdot \left( {9 \cdot 9 \cdot 10 \cdot 9} \right)^{1/4} } \right\}} \right)^{1/2} = 9.59.$$

Based on the defined lower and upper limits, the following FRNs are defined:

$$FRN\left( {\widetilde{\theta }_{12}^{(1)} } \right) = FRN\left( {\widetilde{\theta }_{12}^{(4)} } \right) = FRN\left( {\widetilde{\theta }_{12}^{(5)} } \right) = \left( {\left[ {5,6.14} \right],\left[ {7,8.16} \right],\left[ {9,9.59} \right]} \right);$$
$$FRN\left( {\widetilde{\theta }_{12}^{(2)} } \right) = FRN\left( {\widetilde{\theta }_{12}^{(3)} } \right) = \left( {\left[ {6.14,7} \right],\left[ {8.16,9} \right],\left[ {9.59,10} \right]} \right)$$

Using the FRN Bonferroni function (26), the fusion of FRNs was performed, and an aggregate value was obtained \(\overline{\theta }_{12} = \left( {\left[ {\overline{\theta }_{12}^{(l) - } ,\overline{\theta }_{12}^{(l) + } } \right],\left[ {\overline{\theta }_{12}^{(m) - } ,\overline{\theta }_{12}^{(m) + } } \right],\left[ {\overline{\theta }_{12}^{(u) - } ,\overline{\theta }_{12}^{(u) + } } \right]} \right) = \left[ {\left( {5.68,6.65} \right),\left( {7.69,8.66} \right),\left( {9.35,9.84} \right)} \right]\). The following section shows the application of the Bonferroni function (26) for the fusion of the FRN sequence \(\overline{\theta }_{12}^{(l) - }\):

$$\overline{\theta }_{12}^{(l) - } = \left( {\frac{1}{5(5 - 1)}\left( {5 \cdot 5 + 5 \cdot 5 + ... + 6.14 \cdot 6.14} \right)} \right)^{{\frac{1}{1 + 1}}} = 5.68$$

The remaining FRN sequences were obtained similarly. The aggregated FRN initial decision matrix is shown in Table 9.

Table 9 Aggregated initial decision matrix

Step 2. Normalization of the basic matrix element. After obtaining the aggregated decision matrix \(\overline{\overline{\aleph }}^{N} = \left[ {\overline{\varphi }_{ij} } \right]_{5 \times 18}\) (Table 9), we apply expression (23) to obtain the normalized matrix elements \(\overline{\overline{\aleph }}^{N}\) in the FRN CODAS methodology.

$$\overline{\overline{\aleph }}^{N} = \begin{array}{*{20}c} {A_{1} } \\ {A_{2} } \\ {A_{3} } \\ \begin{gathered} A_{4} \hfill \\ A_{5} \hfill \\ \end{gathered} \\ \end{array} \begin{array}{*{20}c} {EC11 \, EC12 \, {. }{\text{. }}{. }SO15} \\ {\left[ {\begin{array}{*{20}l} {\left( {\left[ {0.24,0.30} \right],\left[ {0.38,0.46} \right],\left[ {1.00,1.00} \right]} \right)} \hfill \\ {\left( {\left[ {0.13,0.18} \right],\left[ {0.19,0.28} \right],\left[ {0.30,0.65} \right]} \right)} \hfill \\ {\left( {\left[ {0.10,0.10} \right],\left[ {0.10,0.12} \right],\left[ {0.12,0.16} \right]} \right)} \hfill \\ {\left( {\left[ {0.12,0.16} \right],\left[ {0.16,0.23} \right],\left[ {0.24,0.44} \right]} \right)} \hfill \\ {\left( {\left[ {0.10,0.11} \right],\left[ {0.12,0.13} \right],\left[ {0.15,0.18} \right]} \right)} \hfill \\ \end{array} \begin{array}{*{20}l} {\left( {\left[ {0.10,0.11} \right],\left[ {0.12,0.13} \right],\left[ {0.15,0.18} \right]} \right)} \hfill \\ {\left( {\left[ {0.35,0.46} \right],\left[ {0.51,0.63} \right],\left[ {1.00,1.00} \right]} \right)} \hfill \\ {\left( {\left[ {0.10,0.11} \right],\left[ {0.12,0.14} \right],\left[ {0.16,0.19} \right]} \right)} \hfill \\ {\left( {\left[ {0.11,0.13} \right],\left[ {0.13,0.18} \right],\left[ {0.19,0.28} \right]} \right)} \hfill \\ {\left( {\left[ {0.10,0.10} \right],\left[ {0.10,0.12} \right],\left[ {0.12,0.16} \right]} \right)} \hfill \\ \end{array} {. }{\text{. }}{.}\begin{array}{*{20}l} {\left( {\left[ {0.55,0.73} \right],\left[ {0.76,0.90} \right],\left[ {0.95,1.00} \right]} \right)} \hfill \\ {\left( {\left[ {0.65,0.79} \right],\left[ {0.85,0.96} \right],\left[ {1.00,1.03} \right]} \right)} \hfill \\ {\left( {\left[ {0.59,0.77} \right],\left[ {0.79,0.94} \right],\left[ {0.97,1.02} \right]} \right)} \hfill \\ {\left( {\left[ {0.44,0.58} \right],\left[ {0.65,0.79} \right],\left[ {0.85,0.96} \right]} \right)} \hfill \\ {\left( {\left[ {0.59,0.69} \right],\left[ {0.80,0.90} \right],\left[ {0.97,1.02} \right]} \right)} \hfill \\ \end{array} } \right]} \\ \end{array}$$

Step 3. Calculate FRN weighted normalized decision matrix (\(\overline{\overline{\Omega }}\)). Applying the expression (24), the FRN weighted normalized matrix \(\overline{\overline{\Omega }} = [\overline{\tau }_{ij} ]_{5 \times 18}\) is calculated. Where the element \(\overline{\tau }_{11}\) we obtain as follows

$$\begin{aligned} \overline{\tau }_{11} = &\overline{w}_{1} \cdot \overline{\varphi }_{11} = \left( \left[ {0.002,0.008} \right],\left[ {0.028,0.075} \right],\left[ {0.077,0.080} \right] \right) \cdot \\ &\left( {\left[ {0.24,0.30} \right],\left[ {0.38,0.46} \right],\left[ {1.00,1.00} \right]} \right) \hfill \\ \, = &\left( {\left[ {0.0005,0.0024} \right],\left[ {0.0105,0.0349} \right],\left[ {0.0767,0.0803} \right]} \right). \hfill \\ \end{aligned}$$

The remaining elements of the matrix \(\overline{\overline{\Omega }} = [\overline{\tau }_{ij} ]_{5 \times 18}\) are calculated similarly.

Step 4. Determine FR negative-ideal solution. We obtain the FRN negative-ideal solution matrix \(\overline{\overline{NS}} = [\overline{\eta }_{j} ]_{1 \times 18}\) as follows:

$$\begin{aligned} \overline{\eta }_{EC11} = & \min \left\{ {\overline{\tau }_{1 - EC11} ;\overline{\tau }_{2 - EC11} ;\overline{\tau }_{2 - EC11} ;\overline{\tau }_{3 - EC11} ;\overline{\tau }_{4 - EC11} ;\overline{\tau }_{5 - EC11} } \right\} \\ = & \min \left\{ \begin{gathered} \left( {\left[ {0.0005,0.0024} \right],\left[ {0.0105,0.0349} \right],\left[ {0.0767,0.0803} \right]} \right); \hfill \\ \left( {\left[ {0.0003,0.0014} \right],\left[ {0.0051,0.0207} \right],\left[ {0.0233,0.0519} \right]} \right); \hfill \\ \left( {\left[ {0.0002,0.0008} \right],\left[ {0.0029,0.0091} \right],\left[ {0.0093,0.0125} \right]} \right); \hfill \\ \left( {\left[ {0.0003,0.0012} \right],\left[ {0.0044,0.0171} \right],\left[ {0.0184,0.0355} \right]} \right); \hfill \\ \left( {\left[ {0.0002,0.0008} \right],\left[ {0.0032,0.0098} \right],\left[ {0.0115,0.0141} \right]} \right) \hfill \\ \end{gathered} \right\} \\ = & \left( {\left[ {0.0002,0.0008} \right],\left[ {0.0029,0.0091} \right],\left[ {0.0093,0.0125} \right]} \right) \\ \end{aligned}$$

Applying the same calculations, we obtain other elements of the negative-ideal solution matrix \(\overline{\overline{NS}} = [\overline{\eta }_{j} ]_{1 \times 18}\)

$$\begin{gathered} \overline{\eta }_{EC12} = \left( {\left[ {0.0004,0.0011} \right],\left[ {0.0039,0.0096} \right],\left[ {0.0102,0.0132} \right]} \right);\\ \overline{\eta }_{EC13} = \left( {\left[ {0.0010,0.0013} \right],\left[ {0.0145,0.0236} \right],\left[ {0.0321,0.0393} \right]} \right); \hfill \\ \overline{\eta }_{EC14} = \left( {\left[ {0.0002,0.0007} \right],\left[ {0.0054,0.0151} \right],\left[ {0.0197,0.0266} \right]} \right);\\ \overline{\eta }_{EC15} = \left( {\left[ {0.0008,0.0014} \right],\left[ {0.0071,0.0120} \right],\left[ {0.0148,0.0170} \right]} \right); \hfill \\ \overline{\eta }_{EC16} = \left( {\left[ {0.0006,0.0017} \right],\left[ {0.0112,0.0354} \right],\left[ {0.0367,0.0564} \right]} \right);\\ \overline{\eta }_{EC17} = \left( {\left[ {0.0004,0.0018} \right],\left[ {0.0093,0.0347} \right],\left[ {0.0360,0.0561} \right]} \right); \hfill \\ \overline{\eta }_{EC18} = \left( {\left[ {0.0002,0.0009} \right],\left[ {0.0043,0.0154} \right],\left[ {0.0163,0.0240} \right]} \right);\\ \overline{\eta }_{EN11} = \left( {\left[ {0.0011,0.0035} \right],\left[ {0.0333,0.0737} \right],\left[ {0.0976,0.1256} \right]} \right); \hfill \\ \overline{\eta }_{EN12} = \left( {\left[ {0.0045,0.0092} \right],\left[ {0.0823,0.1748} \right],\left[ {0.2083,0.3195} \right]} \right);\\ \overline{\eta }_{EN13} = \left( {\left[ {0.0008,0.0025} \right],\left[ {0.0228,0.0468} \right],\left[ {0.0602,0.0800} \right]} \right); \hfill \\ \overline{\eta }_{EN14} = \left( {\left[ {0.0008,0.0053} \right],\left[ {0.0274,0.0822} \right],\left[ {0.0874,0.1250} \right]} \right);\\ \overline{\eta }_{EN15} = \left( {\left[ {0.0016,0.0097} \right],\left[ {0.0402,0.1051} \right],\left[ {0.1186,0.1502} \right]} \right); \hfill \\ \overline{\eta }_{SO11} = \left( {\left[ {0.0072,0.0100} \right],\left[ {0.1419,0.2342} \right],\left[ {0.3126,0.3778} \right]} \right);\\ \overline{\eta }_{SO112} = \left( {\left[ {0.0052,0.0068} \right],\left[ {0.077,0.1266} \right],\left[ {0.1440,0.1944} \right]} \right); \hfill \\ \overline{\eta }_{SO13} = \left( {\left[ {0.0024,0.0039} \right],\left[ {0.0696,0.1083} \right],\left[ {0.1535,0.1969} \right]} \right);\\ \overline{\eta }_{SO14} = \left( {\left[ {0.0014,0.0029} \right],\left[ {0.0496,0.0758} \right],\left[ {0.1073,0.1345} \right]} \right); \hfill \\ \overline{\eta }_{SO15} = \left( {\left[ {0.0032,0.0069} \right],\left[ {0.0618,0.1141} \right],\left[ {0.1280,0.1526} \right]} \right); \hfill \\ \end{gathered}$$

Step 5. Calculate the FRN weighted Euclidean (\(\overline{\overline{ED}}_{i}\)), and FRN weighted Hamming (\(\overline{\overline{HD}}_{i}\)) distances. We obtain the \(\overline{\overline{ED}}_{i}\) and \(\overline{\overline{HD}}_{i}\) distances of alternatives from the FRN negative-ideal solution, using Eqs. (26)–(29):

  1. (a)

    By applying Eqs. (26) and (27), we calculate \(\overline{\overline{ED}}_{i}\) for the first alternative:

    $$\overline{\overline{ED}}_{1} = \sum\limits_{j = 1}^{18} {\left( {0.0368 + 0.0013 + 0.0055 + \cdots + 0.00} \right) = 0.942}$$

    where each sequence is calculated using Eq. (27)

    $$\begin{aligned} d_{E} \left( {\overline{\tau }_{EC11} ;\overline{\eta }_{1} } \right) = & \sqrt {\frac{\begin{gathered} \left\{ {0.0005 - 0.0002} \right\}^{2} + \left\{ {0.00{24} - 0.0008} \right\}^{2} + 2\left\{ {0.0105 - 0.00{29}} \right\}^{2} + \hfill \\ 2\left\{ {0.0349 - 0.0091} \right\}^{2} + \left\{ {0.0767 - 0.0093} \right\}^{2} + \left\{ {0.0803 - 0.0125} \right\}^{2} \hfill \\ \end{gathered} }{8}} \\ = & 0.00368 \\ \end{aligned}$$

    In the same way, we calculate the Euclidean distances \(d_{E} \left( {\overline{\tau }_{ij} ;\overline{\eta }_{j} } \right)\) of the other alternatives:

    $$\overline{\overline{ED}}_{i} = \begin{array}{*{20}c} {A_{1} } \\ {A_{2} } \\ {A_{3} } \\ \begin{gathered} A_{4} \hfill \\ A_{5} \hfill \\ \end{gathered} \\ \end{array} \left[ {\begin{array}{*{20}r} \hfill {0.9417} \\ \hfill {0.9135} \\ \hfill {0.7858} \\ \hfill {0.6879} \\ \hfill {0.8050} \\ \end{array} } \right]$$
  2. (b)

    By applying Eq. (28) and (29) we calculate \(\overline{\overline{HD}}_{i}\) for the first alternative:

    $$\overline{\overline{HD}}_{1} = \sum\limits_{j = 1}^{18} {\left( {0.0872 + 0.0031 + 0.0128 + \cdots + 0.00} \right) = 2.419}$$

where each sequence is calculated using Eq. (29).

$$\begin{aligned} d_{H} \left( {\overline{\tau }_{EC11} ;\overline{\eta }_{1} } \right) = & \frac{\begin{gathered} \left| {0.0005 - 0.0002} \right| + \left| {0.00{24} - 0.0008} \right| + 2\left| {0.0105 - 0.00{29}} \right| + \hfill \\ 2\left| {0.0349 - 0.0091} \right| + \left| {0.0767 - 0.0093} \right| + \left| {0.0803 - 0.0125} \right| \hfill \\ \end{gathered} }{8} \\ = & 2.419 \\ \end{aligned}$$

In the same way, we calculate the Hamming distances \(d_{H} \left( {\overline{\tau }_{ij} ;\overline{\eta }_{j} } \right)\) of the other alternatives:

$$\overline{\overline{HD}}_{i} = \begin{array}{*{20}c} {A_{1} } \\ {A_{2} } \\ {A_{3} } \\ \begin{gathered} A_{4} \hfill \\ A_{5} \hfill \\ \end{gathered} \\ \end{array} \left[ {\begin{array}{*{20}r} \hfill {2.4190} \\ \hfill {2.2615} \\ \hfill {1.6830} \\ \hfill {1.2617} \\ \hfill {1.7367} \\ \end{array} } \right]$$

Step 6. Determine relative assessment matrix. By applying expression (30) we obtain the elements of the relative assessment matrix:

$${\mathbb{Q}} = \begin{array}{*{20}c} {A_{1} } \\ {A_{2} } \\ {A_{3} } \\ \begin{gathered} A_{4} \hfill \\ A_{5} \hfill \\ \end{gathered} \\ \end{array} \begin{array}{*{20}c} {A_{1} \, A_{2} \, A_{3} \, A_{4} \, A_{5} } \\ {\left[ {\begin{array}{*{20}l} {0.000} \hfill & {0.186} \hfill & {0.892} \hfill & {1.411} \hfill & {0.819} \hfill \\ { - 0.186} \hfill & {0.000} \hfill & {0.706} \hfill & {1.225} \hfill & {0.633} \hfill \\ { - 0.892} \hfill & { - 0.706} \hfill & {0.000} \hfill & {0.519} \hfill & { - 0.019} \hfill \\ { - 1.411} \hfill & { - 1.225} \hfill & { - 0.519} \hfill & {0.000} \hfill & { - 0.592} \hfill \\ { - 0.819} \hfill & { - 0.633} \hfill & {0.019} \hfill & {0.592} \hfill & {0.000} \hfill \\ \end{array} } \right]} \\ \end{array}$$

Element of the relative assessment matrix on the position A2-A1 we obtain by applying Eqs. (30) and (31)

$$\begin{aligned} \phi_{21} = & \left( {\overline{\overline{ED}}_{2} - \overline{\overline{ED}}_{1} } \right) + \left( {0.02 \cdot \left( {\overline{\overline{ED}}_{2} - \overline{\overline{ED}}_{1} } \right) \cdot \left( {\overline{\overline{HD}}_{2} - \overline{\overline{HD}}_{1} } \right)} \right) \\ = & \left( {0.9135 - 0.9417} \right) + \left( {0.02 \cdot \left( {0.9135 - 0.9417} \right) \cdot \left( {2.2615 - 2.4190} \right)} \right)\, = - 0.186 \\ \end{aligned}$$

Step 7. Calculate the assessment score (\(\Psi_{i}\)) of each alternative. Finally, we obtain the FRN CODAS criteria function (32) for the final ranking of the alternatives:

$$\begin{array}{*{20}c} { \, \begin{array}{*{20}c} { \, \Psi_{i} } & { \, Rank} \\ \end{array} } \\ {\begin{array}{*{20}c} \begin{gathered} A1 \hfill \\ A2 \hfill \\ A3 \hfill \\ A4 \hfill \\ A5 \hfill \\ \end{gathered} & {\left[ {\begin{array}{*{20}l} {3.308} \hfill \\ {2.379} \hfill \\ { - 1.098} \hfill \\ { - 3.748} \hfill \\ { - 0.841} \hfill \\ \end{array} } \right] \, \left[ {\begin{array}{*{20}c} 1 \\ 2 \\ 4 \\ 5 \\ 3 \\ \end{array} } \right]} \\ \end{array} } \\ \end{array}$$

The alternative should have the highest possible value \(\Psi_{i}\), so the following rank A1 > A2 > A5 > A3 > A4 is obtained. It is noticed from the initial decision matrix that for A1, the respondents are given higher ratings for most of the favorable attributes than non-beneficial criteria. A contrast is noticed in case of the worst performer, i.e., A4. Therefore, the final ranking order is in sync with the common notion of the decision makers.

4.3 Sensitivity analyses

In the following section, the validation of the results obtained in Sect. 3 was performed. The validation of the results was performed through four phases.

Phase 1 The parameters δ1 and δ2 are defined based on subjective assessments of the decision-maker, and the values δ1 = δ2 = 1 were adopted to obtain the initial results. The parameters δ1 and δ2 have a significant influence on the change of rough boundary interval of fuzzy rough numbers, so in some situations, it is expected that the change of the values of the parameters δ1 and δ2 affects the change of the initial results. Therefore, in the first phase, an experiment was conducted in which the change of parameters δ1 and δ2 in the interval 1 ≤ δ1, δ2 ≤ 70 was simulated.

Phase 2 The elements of the relative assessment matrix \({\mathbb{Q}} = [\phi_{ik} ]_{5 \times 5}\) are calculated using expression (30). Following the recommendations of Keshavarz Ghorabaee et al. (2016) for the calculation of the initial solution, the value β = 0.02 for the threshold function was adopted. By analyzing the mathematical form of expression (30), we notice that the threshold parameter β affects the final value of the elements relative assessment matrix. Therefore, in this phase, an experiment was conducted in which the change of the parameter β in the interval 0 ≤ β ≤ 1 was simulated.

Phase 3 In the third phase, an experiment was conducted in which 50 new vectors of weight coefficients of the criteria were formed. In each scenario, the decrease in the value of the most influential criterion by 2% was simulated, while at the same time the values of the remaining criteria were proportionally corrected.

Phase 4 In the fourth phase, a comparison was made with other multi-criteria techniques from the literature. The extant literature shows the utility of comparative analysis of results obtained by using various MCDM methods for validation purpose (Biswas et al., 2021; Pamučar et al., 2021). In this study three multi-criteria techniques were used for comparison: rough COPRAS (COmplex PRoportional ASsessment) method (Pamucar et al., 2018), rough MAIRCA (MultiAtributive Ideal-Real Comparative Analysis) method (Bozanic et al., 2020), and fuzzy MARCOS (Measurement of Alternatives and Ranking according to COmpromise Solution) methodology (Stankovic et al., 2020). There are no examples in the literature of applying the presented fuzzy rough methodology in multi-criteria decision-making. Therefore, extensions of COPRAS, MAIRCA and MARCOS methods in fuzzy and rough environment were used for comparison.

4.3.1 Change of values of parameters δ 1 and δ 2

The FRN approach presented in this study is based on the representation of inaccuracies and uncertainties using a rough boundary interval (RBI). A rough boundary interval is defined based on objective and subjective values. Objective values include a set of data based on which inaccuracy is defined, while subjective parameters include the parameters of the Bonferroni function used to calculate the RBI limit values. The values δ1 = δ2 = 1 were adopted for the calculation of the initial solution. Changes in the values of the parameters δ1 and δ2 in expressions (7)–(12) affect the change in RBI, which can also affect the final results of the decision-making model. Figure 3 shows the dependence of the assessment score on the change of the parameters δ1 and δ2 in the interval 1 ≤ δ1, δ2 ≤ 70.

Fig. 3
figure 3

Dependence of assessment score on parameter change in the interval 1 ≤ δ1, δ2 ≤ 70

Figure 3 shows the dependence of the assessment score of alternatives A1, A2, and A3 on the change of parameters 1 ≤ δ1, δ2 ≤ 70. Similarly, changes in the parameters δ1 and δ2 affect the change in the assessment score of the remaining alternatives. Increasing the values of the parameters δ1 and δ2 affects the increase in the RBI fuzzy rough numbers. RBI is a measure of uncertainty and inaccuracy in data. Adopting the value δ1 = δ2 = 1 in the initial scenario represents an optimistic scenario since the RBI is the smallest, and therefore the inaccuracy and risk in the information are the smallest. Figure 4 shows the change in assessment score alternatives through 70 scenarios. In each scenario, the value of the parameters δ1 and δ2 was increased by one.

Fig. 4
figure 4

Influence of parameter change 1 ≤ δ1, δ2 ≤ 70 on change of assessment scores

Decision-makers can set different values of the parameters δ1 and δ2 based on their preferences. However, the results from Fig. 4 indicate that the value of δ1 and δ2 should not be minimal, so it is recommended that the value be between four and eight. Also, from Fig. 4, we notice that the initial rank A1 > A2 > A5 > A3 > A4 is confirmed regardless of the change of parameters δ1 and δ2, so we can conclude that the initial solution is stable regardless of the level of risk when making decisions.

4.3.2 Change the value of the parameter β

As previously stated, when defining the initial solution, the value of the parameter β = 0.02 was adopted. It is evident that the value of the parameter β influences the definition of the elements of the relative assessment matrix elements. Therefore, in the next part, an experiment was conducted in which 100 scenarios were created. In the first scenario, the value β = 0.00 was adopted, while in each subsequent scenario, the value of β was defined using the expression βi = βi−1 + 0.01 (i = 1, 2,…,100). Simultaneously with the change in the value of the parameter β, the impact of these changes on the value of the assessment score alternative was analyzed, Fig. 5.

Fig. 5
figure 5

Influence of change of parameter 0 ≤ β ≤ 1 on change of assessment score alternative

The results in Fig. 5 show that increasing the value of the parameter β leads to a decrease in the value of the assessment score alternatives. This scenario is expected if we consider the influence of the parameter β in the linear function (30). The results show that the decline in the assessment score is constant in all alternatives and that there are no drastic changes that can affect the change in the ranking of alternatives. Also, this indicates to us the fact that the initial rank A1 > A2 > A5 > A3 > A4 was confirmed by both FRN weighted Euclidean and FRN weighted Hamming distances. Based on the above, we can conclude that the initial rank A1 > A2 > A5 > A3 > A4 is credible and that alternative A1 is the dominant alternative.

4.3.3 Analysis of changes in weight coefficients of criteria

In the following part, the influence of the change of the weight coefficient of the most important criterion (wSO14) on the ranking results is analyzed. Based on the recommendations of Kahraman (2002), new weight vectors were generated, and their influence on changes in the rankings of alternatives was analyzed. Since FRN weights were used in this study, for each of the six segments of the fuzzy rough number, the limits of change of the weighting factor of the wSO14 criterion were defined, according to the following: (1) lower bound − 0.0528 ≤ \(\Delta w_{s}^{(l) - }\) ≤ 0.163 and − 0.0731 ≤ \(\Delta w_{s}^{(l) + }\) ≤ 0.327; (2) modal value − 0.508 ≤ \(\Delta w_{s}^{(m) - }\) ≤ 1.789 and − 0.732 ≤ \(\Delta w_{s}^{(m) + }\) ≤ 2.827; and (3) upper bound − 0.7426 ≤ \(\Delta w_{s}^{(u) - }\) ≤ 2.965 and − 0.752 ≤ \(\Delta w_{s}^{(u) + }\) ≤ 3.047. Defined intervals are divided by 50 sequences within which new vectors of weight coefficients were generated. Thus, 50 new vectors of weight coefficients were formed. Figure 6 shows the new weighting vectors for the modal value \(w_{j}^{(m) + }\).

Fig. 6
figure 6

New vectors of criteria weights for modal value \(w_{j}^{(m) + }\)

After the calculation of new FRN vectors, their influence on the change of assessment score alternatives was analyzed. The influence of FRN vectors on the change of assessment score alternatives is shown in Fig. 7.

Fig. 7
figure 7

Influence of FRN vector on change of assessment score alternative

The results shown in Fig. 7 show that new vectors of criterion weight coefficients lead to a change in assessment score alternatives. These changes show that the FRN MACBETH CODAS model is sensitive to changes in the weight coefficients of the criteria, which is one of the significant features of the MCDM model. Also, the results show that a decrease in the significance of the most influential criterion leads to an increase in the assessment score of the dominant alternative (A1), and that the initial rank A1 > A2 > A5 > A3 > A4 was maintained through the scenarios.

4.4 Comparative analysis

In order to show the efficiency of the novel developed FRN approach, this section performs a comparative analysis test to compare the results of the developed approach with other FRN-based MCDM methods. The rough COPRAS method (Pamucar, Božanić, et al., 2018), the rough MAIRCA method (Bozanic et al., 2020) and the fuzzy MARCOS method (Stankovic et al., 2020) were chosen to compare the results. The presented MCMD techniques were used to compare the results, which, like the CODAS method, use a linear approach for normalizing the data. This criterion for selecting MCMD techniques was used because many studies (Pamucar & Cirovic, 2015; Yazdani, Torkayesh, et al., 2020) show that the application of different normalization techniques in MCMD models in some situations can lead to distortion of the results. These multi-criteria techniques were applied under the same conditions, and the results shown in Fig. 8 were obtained.

Fig. 8
figure 8

Ranks of alternatives based on different MCDM techniques

The results from Fig. 8 show that the application of FRN CODAS, rough COPRAS, and Rough MAIRCA techniques gave the same rank. Minor changes in the range of alternatives A2 and A5 occurred due to the fuzzy MARCOS method. However, by applying all multi-criteria techniques, alternative A1 has been proposed as the best solution, so we can conclude that it represents the dominant alternative from the considered set. In the next part (Table 10), a comparison of the applied multi-criteria techniques is performed and their specifics are pointed out.

Table 10 Comparisons of different methods

Based on the comparison of the results, it can be concluded that A1 is the best alternative according to all methodologies. However, unlike classical extensions of MCMD techniques in a fuzzy or rough environment, FRN CODAS methodologies allow flexible decision-making due to decision makers' risk attitudes and interrelationships between decision-making matrix parameters. These characteristics of the FRN CODAS model represent a significant advantage that affects the flexibility of the presented methodology when applied in real applications.

The Fuzzy rough methodology presented in this study is based on defining a rough boundary interval using Bonferroni functions. The application of the Bonferroni function has improved the adaptability of the multi-criteria framework as the rough boundary interval adjusts to the level of risk presented in the information. Applying the FRN approach in the CODAS methodology enables the presentation of the degree of agreement in expert assessments. A rough boundary interval defines the degree of agreement in the presented methodology. More significant discrepancies in expert estimates lead to an increase in the rough boundary interval of fuzzy rough numbers, while in the case of complete agreement, the fuzzy rough number is transformed into a fuzzy number. On the other hand, rough COPRAS, rough MAIRCA, and fuzzy MARCOS techniques do not have such adaptability, limiting their ability to process uncertainty and inaccuracy in real data.

5 Conclusions

This study presented a novel decision-making approach to address the supplier selection problem in HSC management during the COVID-19 pandemic. This study contributes to the literature in different directions as follows. First, this is the first study to develop FRN extension of MACBETH and CODAS methods. Also, this study for the first time developed an integrated decision-making approach by combining FRN MACBETH-CODAS methods. On the other hand, this study is first in its kind to use an FRN-based approach to address supplier selection problem of medical masks and face shields during COVID-19 pandemic. Both of the face masks and shields are highly consumed in hospitals and different places. The most important way to prevent the virus spread to use face masks and shields, specifically in medical departments that are dealing with infected people. Therefore, addressing supplier selection problem is one of the most critical operations for all hospitals in order to prevent any consequences caused from shortage of such products. The developed decision-making approach is applied for a real-life case study to address supplier selection in HSC management of a Turkish hospital in Istanbul during the COVID-19 pandemic. The developed methodology empowers real-life decision-makers in healthcare managerial sector to express their opinion against characteristics of the suppliers in flexible way through FRN environment. Thus, results of the developed methodology for the investigated case provide strong implications for suitability of suppliers based on their technical and sustainability performances.

Although this study develops a novel decision-making approach to address an important SC problem during COVID-19, there exist some limitations that should be considered in future studies. First, HSC operations like other SC operations are very depended the daily changes that may happen in the market. Therefore, due to high impact of future events on the decision-making problems, considering impacts of future events would help managers to make more reliable decisQ5ions. Stratification concept (Torkayesh et al., 2021) is one of the reliable concepts that can be integrated with the developed FRN-based approach to include such impacts in the calculations. Such impacts can be included in both of weighting and ranking operations of MCDM approaches.

This study can be extended in various future directions as follows. The most common direction is to utilize the developed methodology to investigate supplier selection problem in different SC applications such as transportation, waste management, perishable products, cosmetics, and dairy products. One may use other weighting method instead of MACBETH such as BWM which are more popular to address complex decision-making approaches. Another way to enhance reliability of the weight determination process is to integrate objective and subjective weight methods to determine optimal decision criteria with lower biasedness. The same can happen for the ranking method where CODAS can be integrated with other ranking MCDM methods such as CoCoSo, COPRAS, or MARCOS where finally a voting method can show the consensus final ranking.