Clifford W Fong
email: cwfong@internode.on.net

Covid-19 drug repurposing: evaluation of inhibitors in SARS-CoV-2 infected cell lines

Keywords: Caco-2, VeroE6, VeroCCL81, HuH, Calu-3, COVID-2019, SARS-CoV-2; ACE2 receptor binding, spike serine proteases, S-RBD, TMPRSS2, IC 50, linear free energy relationships, HOMO-LUMO, quantum mechanics;

Structure activity relationships SAR, ΔG desolv,CDS free energy of water desolvation, ΔG lipo,CDS lipophilicity free energy, cavity dispersion solvent structure of the first solvation shell CDS, Dipole moment DM, Molecular Volume Vol, HOMO highest occupied molecular orbital, LUMO lowest unoccupied molecular orbital, HOMO-LUMO energy gap, linear free energy relationships LFER, Receptor binding domain of S protein of

Introduction

The susceptibility of various cell lines to repurposed antiviral drugs offers a rapid and effective screening method for evaluating host cell responses to viral infection such as Covid-19. Examination of different human cell lines may have relevance to therapeutic treatment of various human tissues. Inhibition of SARS-CoV-2 infected Caco-2 cells can be studied using the LFER method to separate Caco-2 cell entry processes involving ACE2, TMPRSS2 or S-RBD from intracellular inhibitory processes. [START_REF] Fong | Predicting drug permeability and relationship with COVID-19 infected Caco-2 cells, hal archives[END_REF] The extensive study by Ellinger [START_REF] Ellinger | Identification of inhibitors of SARS-CoV-2 in-vitro cellular toxicity in human (Caco-2) cells using a large scale drug repurposing collection[END_REF] likely predominantly involves the inhibition of Caco-2 cell entry processes involving ACE2, TMPRSS2 or S-RBD. Ellinger's data for Caco-2 cells is a valuable source for evaluating the efficacy of SARS-CoV-2 therapeutics.

Chu [3] evaluated the replicative capability of SARS-CoV and SARS-CoV-2 in 25 different cell lines including nine of human origin, including pulmonary (Calu-3), intestinal (Caco-2), hepatic (Huh-7), and neuronal (U251) cells. It was found that both human Calu-3 cells and Caco-2 cells showed the greatest replication and were best suited for studying SARS-CoV-2 replicative processes. It is also known that Caco-2 cells were the only human cell type of 13 tested refractory cell lines that supported efficient SARS-CoV replication and expression of the SARS-CoV receptor, ACE2. [START_REF] Ec Mossel | Exogenous ACE2 Expression Allows Refractory Cell Lines To Support Severe Acute Respiratory Syndrome Coronavirus Replication[END_REF] It is well known that choice of cell lines for identifying SARS-CoV-2 antiviral drug efficacy is crucial: for example chloroquine when tested in Vero cells which are derived from the kidney cells of the African Green monkey, showed antiviral properties; but when tested in Calu-3 human lung cells, it showed no activity. It has been shown that chloroquine can inhibit the virusactivating enzyme cathepsin L in Vero cells, so chloroquine has an inhibitory effect on this enzyme and can thus prevent the coronavirus from infecting the Vero cell. Cathepsin L requires an acidic environment to function. Both hydroxychloroquine and chloroquine decrease the acidity, which then disables the cathepsin L enzyme, blocking the virus from infecting the Vero cells. But using the more clinically appropriate human lung cells Calu-3, the virus is activated by the protease TRMPSS2 to initiate cellular infection. Chloroquine and hydroxychloroquine target a pathway for viral activation that is not active in Calu-3 lung cells, which have very low levels of cathepsin L. TMPRSS2 does not require an acidic environment to function, so chloroquine or hydroxychloroquine are not inhibitory in human lung cells. [START_REF] Hoffmann | Chloroquine does not inhibit infection of human lung cells with SARS-CoV-2[END_REF][START_REF] Hoffmann | SARS-CoV-2 Cell Entry Depends on ACE2 and TMPRSS2 and Is Blocked by a Clinically Proven Protease Inhibitor[END_REF][START_REF] Hoffmann | Nafamostat mesylate blocks activation of SARS-CoV-2:New treatment option 1for COVID-19[END_REF][START_REF] Ou | Hydroxychloroquine-mediated inhibition of SARS-CoV-2 entry is attenuated by TMPRSS2[END_REF] So virus entry inhibition is cell type dependent, as both pH-dependent and pH-independent pathways are available for entry into cells. The spike (S) protein of SARS-CoV-2, which mediates viral entry, is activated by the endosomal pH dependent cysteine protease cathepsin L in some cell lines. However entry into airway epithelial cells, which express low levels of cathepsin L, depends on the pH-independent transmembrane TMPRSS2. [START_REF] Hoffmann | Chloroquine does not inhibit infection of human lung cells with SARS-CoV-2[END_REF][START_REF] Hoffmann | SARS-CoV-2 Cell Entry Depends on ACE2 and TMPRSS2 and Is Blocked by a Clinically Proven Protease Inhibitor[END_REF][START_REF] Hoffmann | Nafamostat mesylate blocks activation of SARS-CoV-2:New treatment option 1for COVID-19[END_REF][START_REF] Ou | Hydroxychloroquine-mediated inhibition of SARS-CoV-2 entry is attenuated by TMPRSS2[END_REF][START_REF] Khan | Role of Endolysosomes in Severe Acute Respiratory Syndrome Coronavirus-2 Infection and Coronavirus Disease 2019 Pathogenesis: Implications for Potential Treatments[END_REF] Mykytyn et al [START_REF] Az Mykytyn | SARS-CoV-2 entry into human airway organoids is serine protease-mediated and facilitated by the multibasic cleavage site[END_REF] used human airway organoids (hAOs) as a better test vehicle than Calu-3 cells and found that the multibasic cleavage site (MBCS) of the SARS-CoV-2 S spike increased infectivity of the hAOs as well as forming syncytia (which are implicated in Covid-19 pathogenesis). The hAOs express ACE-2 and TMPRSS2. It was found that the MBCS increased entry speed and plasma membrane serine protease usage relative to cathepsin-mediated endosomal entry. Blocking serine proteases, but not cathepsins, effectively inhibited SARS-CoV-2 entry and replication in hAOs. Also it was shown that the SARS-CoV-2 MBCS facilitates serine protease-mediated entry in Calu-3 cells. SARS-CoV-2 infectivity requires binding to the ACE2 receptor followed by proteolytic activation of the S protein by host proteases at the cleavage site located at the S1/S2 boundary (S1/S2 cleavage site) and within the S2 domain (S2′ cleavage site). Several host proteases, including endosomal cathepsins, cell surface TMPRSS2 proteases, furin, and trypsin, have been identified to be responsible for S protein cleavage during virus entry or viral protein biogenesis in coronaviruses. Inhibitors of TMPRSS2 and furin have been identified as having antiviral potential against the SARS-CoV-2 virus. [START_REF] Hoffmann | A multibasic cleavage site in the spike protein of SARS-CoV-2 is essential for infection of human lung cells[END_REF]12] Dittmar et al [START_REF] Dittmar | Drug repurposing screens reveal FDA approved drugs active against SARS-Cov-2[END_REF] have evaluated a wide range of antivirals for activity against SARS-CoV-2 in the human HuH7.5, Vero and human Calu-3 cells. Major differences in drug sensitivity and entry pathways used by SARS-CoV-2 were found in these cells. Entry in lung epithelial Calu-3 cells is pH-independent and requires TMPRSS2, while entry in Vero and Huh7.5 cells requires low pH and is triggered by acid-dependent endosomal proteases. Ko et al [14] evaluated 24 FDA approved drugs (out of 28) for antiviral efficacy against SARS-CoV-2 in Vero and Calu-3 cells and found that 16 drugs were less effective in Calu-3 cells than in Vero cells, 6 drugs showed the same efficacy, and 6 drugs were more effective in Calu-3 cells, particularly nafamostat mesylate and camostat mesylate. Remdesivir, hydroxyprogesterone caproate, digitoxin, and cyclosporine showed modest increases only. Nafamostat can prevent the fusion of the envelope of the virus with surface membranes of host cells, and can inhibit membrane fusion at a concentration of less than one-tenth of that of camostat mesylate. Both nafamostat and camostat inhibit TMPRSS2. [START_REF] Hoffmann | SARS-CoV-2 Cell Entry Depends on ACE2 and TMPRSS2 and Is Blocked by a Clinically Proven Protease Inhibitor[END_REF][START_REF] Hoffmann | Nafamostat mesylate blocks activation of SARS-CoV-2:New treatment option 1for COVID-19[END_REF] Study objective:

Evaluate the application of a previously described quantitative LFER method used to describe the antiviral inhibition of repurposed drugs against the SARS-CoV-2 virus, focussing on different infected cell lines to identify predictive indicators of antiviral efficacy that may be relevant to various human tissues.

Results

We have previously shown that the general equation 1 can be applicable to the inhibition of the various proteases involved in the SARS-CoV-2, SARS-CoV and MERS-CoV replication process. Also eq 1 can be useful screening tools to evaluate the potential efficacy of therapeutic drugs that may be active against the SARS-CoV-2 virus. [START_REF] Fong | Predicting drug permeability and relationship with COVID-19 infected Caco-2 cells, hal archives[END_REF][START_REF] Fong | Predicting inhibition of SARS-CoV-2 in Caco-2 cells: viral cellular entry[END_REF][START_REF] Fong | Predicting inhibition of repurposed drugs for SARS-CoV-2 viral activity and cellular entry[END_REF] Eq 1 Inhibition COVID-19 = ΔG desolv,CDS + ΔG lipo,CDS + Dipole Moment + Molecular Volume + HOMO-LUMO With the following independent variables, or molecular specifiers: ΔG desolv,CDS is the free energy of water desolvation, ΔG lipo,CDS is the lipophilicity free energy in octane, DM is the dipole moment in water, Mol Vol is the molecular volume in water, and HOMO-LUMO is the energy gap in water. The general method is to test inhibition against all of the five molecular specifiers for the drugs totally unconstrained, then to identify which of the molecular specifiers best fits the inhibitory experimental data.

Many in silico studies seeking to find potential repurposed antivirals that may be effective against SARS-CoV-2 have used molecular docking screening of the ACE2 receptor, TMPRSS2 protein, as well as inhibitors of the M pro involved in virus replication. The extension of these searches using antivirals with virus infected cells is more problematical because it is not clear which site(s) is actively inhibited when the various drugs are found to be effective in SARS-CoV-2 infected cell lines. Since the ACE2 and the TMPRSS2 are found on the cell surface of the many cells and tissues that are actively infected by the virus, these targets are likely inhibited first. The inhibition of M pro requires the virus to penetrate into the infected cells as well as the drug to also be transported across the cell membrane.

Ellinger et al [START_REF] Ellinger | Identification of inhibitors of SARS-CoV-2 in-vitro cellular toxicity in human (Caco-2) cells using a large scale drug repurposing collection[END_REF] investigated the inhibition of viral induced cytotoxicity using the human epithelial colorectal adenocarcinoma cell line Caco-2 and a SARS-CoV-2 isolate obtained from an individual originally exposed to the virus in the Wuhan region of China. A total of 64 compounds with IC 50 <20 μM were identified in the final screening.

Analysis of Ellinger's IC 50 data for SARS-CoV-2 infected Caco-2 cells using the general equation 1 for 56 widely structurally diverse drugs is shown below in eq 2. No correlation was found with the dipole moment. The molecular volumes have been scaled by 35 times to allow a direct comparison of the relative magnitudes of the four molecular specifiers. It is noted that the P-values of the coefficients of the ΔG desolvCDS, ΔG lipoCDS (octane) and Molecular Volume molecular specifiers in eq 8(b) are significant at 0.015, 0.018 and 0.003 levels. The low regression coefficient is partly due to the low sensitivity of the coefficients to IC 50 , since the regression coefficient is dependent on the slope of the regression line. Eq 2 also includes 9 positively charged drugs and their neutral counterparts, since these drugs can be charged at neutral pH conditions in vitro. Eq 2 shows dependency on ΔG lipoCDS (octane) but not with ΔG lipoCDS (octanol) and no dependency on the HOMO-LUMO gap. We have previously shown that octanol better represents a cell membrane environment where some water interaction can occur, whereas octane is a better representative of a more hydrophobic environment found in deep pockets where drug-protein environments exists with very little water interaction. [START_REF] Fong | Predicting drug permeability and relationship with COVID-19 infected Caco-2 cells, hal archives[END_REF] The data of Dittmar [START_REF] Dittmar | Drug repurposing screens reveal FDA approved drugs active against SARS-Cov-2[END_REF] (see Table 1) which includes the inhibition of the SARS-CoV-2 infected Calu-3, Vero CCL81 (Dittmar found that SARS-CoV-2 was cytopathic in Vero E6, but not in Vero CCL81, so Vero CCL81 was used to focus on infection rather than toxicity) and human hepatocyte Huh7.5 cells by 22 neutral and 16 potentially charged (at neutral pH experimental levels) antiviral drugs has also been analysed using the general equation 1 to find out which of the drug molecular specifiers can describe strong LFERs, and to compare these results with that found for infected Caco-2 cells. The 38 drugs were tested against ΔG desolv,CDS , ΔG lipo,CDS , octane or octanol, DM, Molec Vol, and the HOMO-LUMO energy gap for virus infected HuH7. It is noted that Ko's data base was mainly comprised of neutral species with only five charged species under the experimental conditions so testing of charged species alone was not possible. The results from eq 4(b) for infected Calu-3 cells showing predominant dependency on the HOMO-LUMO gap and to a lesser extent on ΔG lipo,CDS in octanol is very different from those in eqs 3(c)(i,ii) which show predominant dependency on the lipophilicity, ΔG lipo,CDS in octane.

Discussion

It has been shown that the LFER method applies to a wide and diverse range of drugs that can inhibit the SARS-CoV-2 infected cell lines: Caco-2, HUH, Vero and Calu3, as shown in eqs 2-4. Since virus entry inhibition is cell type dependent, with both pH-dependent and pH-independent pathways are available for entry into cells. Entry in lung epithelial Calu-3 cells is pHindependent and requires TMPRSS2, while entry in Vero and Huh7.5 cells requires low pH and is triggered by acid-dependent endosomal proteases. It can be seen that the charged status of the drugs has a major influence on inhibitory behaviour, as shown in eqs 3(a)-(c).

Eqs 3(c)(i) and 3(c)(ii) are consistent with previous findings that entry in lung epithelial Calu-3 cells is pH-independent and requires TMPRSS2, since the neutral inhibitors show similar behaviour to those for the charged drugs, except for a significant water desolvation effect ΔG desolv,CDS in eq 3(c)(ii) which is not present for the neutral species in eq 3(c)(i).

However eq 4(b) for a different set of antiviral drugs with Calu-3 cells gives a different outcome with a strong dependence on the HOMO-LUMO gap, and a lesser dependence on the ΔG lipo,CDS lipophilicity in octanol rather than octane. This lipophilicity differs from previous dependency in all eqs 3(a-c). We have previously found that octanol better represents a cell membrane environment where some water interaction can occur, whereas octane is a better representative of a more hydrophobic environment found in deep pockets where drug-protein environments exists with very little water interaction. [START_REF] Fong | Predicting drug permeability and relationship with COVID-19 infected Caco-2 cells, hal archives[END_REF] TMPRSS2 is a type II transmembrane serine protease (TTSPs), with an extracellular region composed of a low-density lipoprotein (LDL) receptor class A domain, a scavenger receptor cysteine-rich (SRCR) domain and a peptidase S1 domain containing the catalytic triad. There is also a transmembrane sector and a cytoplasmic domain. [START_REF] Gn Brooke | Structural and functional modelling of SARS-CoV-2 entry in animal models[END_REF] Singh et al [START_REF] Singh | Structure-based drug repositioning over the human TMPRSS2 protease domain: search for chemical probes able to repress SARS-CoV-2 Spike protein cleavages[END_REF] have used docking studies to characterize drugs which are more likely to bind with the peptidase S1 domain containing the catalytic triad or a different set of drugs which can bind with the exosites. The exosites are more solvent-exposed than the catalytic site pocket and drugs binding to exposed pockets have different physico-chemical properties than drugs that bind to the deeper catalytic cavity. [START_REF] Singh | Structure-based drug repositioning over the human TMPRSS2 protease domain: search for chemical probes able to repress SARS-CoV-2 Spike protein cleavages[END_REF] Nafamostat and camostat are known to bind to the catalytic cavity site in human lung cells. [START_REF] Hoffmann | SARS-CoV-2 Cell Entry Depends on ACE2 and TMPRSS2 and Is Blocked by a Clinically Proven Protease Inhibitor[END_REF][START_REF] Hoffmann | Nafamostat mesylate blocks activation of SARS-CoV-2:New treatment option 1for COVID-19[END_REF] It is not clear why eq 4(b) is quite different from eq 3(c)(i) and 3(c)(ii) but a different antiviral drug set and possibly differences in the notional Calu-3 cell lines used in the Dittmar [START_REF] Dittmar | Drug repurposing screens reveal FDA approved drugs active against SARS-Cov-2[END_REF] and Ko [14] studies may be responsible. However the dependence on ΔG lipo,CDS ,octanol in eq 4(b) and ΔG lipo,CDS ,octane in eq 3(c)(i) and 3(c)(ii) is important, especially if ΔG lipo,CDS ,octanol better represents an exosite which is solvent exposed, and ΔG lipo,CDS ,octane better represents a buried site, such as the catalytic site of TMPRSS2 as defined by Singh. [START_REF] Singh | Structure-based drug repositioning over the human TMPRSS2 protease domain: search for chemical probes able to repress SARS-CoV-2 Spike protein cleavages[END_REF]

Conclusions

We have shown that the LFER method can be used to quantify the antiviral inhibition of various drugs in infected Caco-2, Vero, HuH and Calu-3 cell lines. A quantitative dependency can be found with various molecular specifiers such as the ΔG desolv,CDS , the free energy of water desolvation, ΔG lipo,CDS the lipophilicity free energy in octane or octanol, DM the dipole moment in water, Mol Vol the molecular volume in water, and HOMO-LUMO the energy gap in water.

It appears that the inhibition of the various infected cell lines is the result of factors specific to the inhibitors and the different cell types, which have implications for exploratory studies of potential therapeutics of the antiviral efficacy in various human organs and tissues.

There is some evidence that it may be possible to identify various drugs that may target the more solvent exposed exosite of TMPRSS2 in infected Calu-3 cells.

Materials and methods

All calculations were carried out using the Gaussian 09 package. Energy optimizations were at the DFT/B3LYP/6-31G(d) (6d, 7f) level of theory for all atoms in water. Selected optimizations at the DFT/B3LYP/6-311+G(d,p) (6d, 7f) level of theory gave very similar results to those at the lower level. Energy calculations were conducted at the DFT/B3LYP/6-31G(d,p) (6d, 7f) for neutral and cationic compounds with optimized geometries in water, using the IEFPCM/SMD solvent model. With the 6-31G* basis set, the SMD model achieves mean unsigned errors of 0.6 -1.0 kcal/mol in the solvation free energies of tested neutrals and mean unsigned errors of 4 kcal/mol on average for ions. [START_REF] Av Marenich | Universal Solvation Model Based on Solute Electron Density and on a Continuum Model of the Solvent Defined by the Bulk Dielectric Constant and Atomic Surface Tensions[END_REF] The 6-31G** basis set has been used to calculate absolute free energies of solvation and compare these data with experimental results for more than 500 neutral and charged compounds. The calculated values were in good agreement with experimental results across a wide range of compounds. [21,2 2] Adding diffuse functions to the 6-31G* basis set (ie 6-31+G**) had no significant effect on the solvation energies with a difference of less than 1% observed in solvents, which is within the literature error range for the IEFPCM/SMD solvent model. It is noted that high computational accuracy for each species in different environments is not the focus of this study, but comparative differences between various species is the aim of the study. Experimental errors in inhibitory and docking binding studies are substantially higher than those for calculated molecular specifiers.

  Eq 2 IC 50 = 0.76ΔG desolv,CDS + 0.83ΔG lipo,CDS + 1.71Mol Vol + 4.34 Where R 2 = 0.174, SEE = 5.55, SE(ΔG desolvCDS ) = 0.30, SE(ΔG lipoCDS ) = 0.34, SE(Mol Vol) = 0.54, F=3.64, Significance= 0.018

5 ,

 5 Vero CCL81 and Calu-3 cells. The best fit unconstrained relationships with EC 50 are shown in eqs 3(a-c) below. Also shown in eqs 3(a-c) are the subsets of the primary 38 drug subdivided into their 22 neutral and 17 ionized drug groups. Eq 3(a) Dittmar 38 drugs with infected HuH7.5 cells EC 50 = -0.050ΔG desolv,CDS -0.008DM -0.038Mol Vol + 0.793 Where R 2 = 0.192, SEE = 0.406, SE(ΔG desolvCDS ) = 0.023, SE(DM)= 0.003, SE(Mol Vol) = 0.0016, F=2.69, Significance= 0.061. P-values: ΔG desolvCDS = 0.040, DM= 0.027, Mol Vol=0.040; Mol Vol scaled by 20 times. Eq 4(b) Ko 28 drugs with infected Calu-3 cells IC 50 = -3.743ΔG lipo,CDS -9.687HOMO-LUMO + 57.108 Where R 2 = 0.414, SEE = 14.15, SE(ΔG lipo,CDS ) = 1.227, SE(HOMO-LUMO)= 2.996, F=8.828, Significance= 0.00125. Pvalues: ΔG lipo,CDS (octanol) = 0.0053, HOMO-LUMO=0.00343. No dependence was found with any other molecular specifiers, including ΔG lipo,CDS = octane.

Table 1 . Inhibition of SARS-CoV-2 infected HuH, VeroCCL81, Calu-3 cells [13] and molecular specifiers of inhibitors EC 50 HuH μM EC 50 Vero μM EC 50 Calu3 μM ΔGdesolv,CDS Kcal/mol ΔGlipo,CDS octane kcal/mol ΔGlipo,CDS octanol kcal/mol

 1

	Dipole	Molec	HOMO-
	Moment	Vol	LUMO
	D	cm	3 /mol	eV

Table 2 . Inhibition of SARS-CoV-2 infected Vero and Calu-3 cells [14] and molecular specifiers of inhibitors IC 50 Vero μM IC 50 Calu3 μM ΔGdesolv,CDS Kcal/mol ΔGlipo,CDS octane kcal/mol ΔGlipo,CDS octanol kcal/mol

 2

	Dipole	Molec	HOMO-
	Moment	Vol	LUMO
	D	cm	3 /mol	eV


Gen-Wen Hu

Cun-Jing Zheng

Wei-Xiang Zhong

Dan-Ping Zhuang

Ben-Heng Xiao

Dr Yì Xiáng

J Wáng

Usefulness of diffusion derived vessel density computed from a simplified IVIM imaging protocol: an experimental study with rat biliary duct blockage induced liver fibrosis

Keywords: Magnetic resonance imaging (MRI), intravoxel incoherent motion (IVIM), diffusion, perfusion, liver, fibrosis, viral hepatitis

Aim: Liver vessel density can be evaluated by DDVD (diffusion derived vessel density): DDVD(b0b1) = Sb0/ROIarea0 -Sb1/ROIarea1, where Sb0 and Sb1 refer to the liver signal when b is 0 or 1 (s/mm 2 ); ROIarea0 and ROIarea1 refer to the region-of-interest on b= 0 or 1 images.

Sb1 and ROIarea1 may be replaced by Sb25 and ROIarea25 or Sb50 and ROIarea50 (i.e., by images of b=25 or 50 s/mm 2 ). With rat biliary duct ligation (BDL) model, this study assessed the usefulness of liver DDVD computed from a simplified IVIM imaging protocol.

Materials and Methods: Male Sprague-Dawley rats were used, and 1, 2, 3, 4 weeks post BDL had rat number of 5, 5, 5,3 respectively. Partial recanalization (n=12) was performed with rats had BDL for 7 days and MRI was done mean 14 days after recanalization. Liver diffusion MR images were acquired using with a 3.0T magnet with a b-value distribution of 0, 30, 50, 75, 100, 150, 300, 700, 1000 s/mm 2 . DDVDmean (control rats n=6) was the mean of DDVD(b0b25) and DDVD(b0b50). IVIM fitting started from b = 0 with segmented fitting and a threshold b of 50 s/mm 2 (n=5 for control rats). Three 3-D space was constructed using PF, D slow , and D fast ; PF, D slow , and DDVD mean ; or PF, D fast , and DDVD mean .

Results: BDL week 1, 2, 3, 4 rats and partial recanalization rats had mean liver fibrosis METAVIR score of 1.5, 2.4, 3.25, 3.7 and 2.42 respectively. The control rats and BDL rats (n=18) had a liver DDVD mean of 84.0±26.2 and 44.7±14.4 au/pixel (p<0.001). Each diffusion MRI parameter alone could not totally separate healthy livers and all fibrotic livers (n=30, BDL and BDL+ recanalization rats). All 3-D spaces totally separated healthy livers and all fibrotic livers. The mean relative distance between healthy liver cluster and fibrotic liver cluster was 0.331 for PF, D slow , and D fast ; Dfast, 0.381 for or PF, D fast , and DDVD mean ; and 0.384 for PF, D slow , and DDVD mean .

Conclusion: A combination of PF, D slow , and D fast allows a total separation of all healthy livers and fibrotic livers, with the integration of DDVD improved the separation.

On liver diffusion weighted imaging, blood vessels show high signal when there is no diffusion gradient (b=0 s/mm 2 ), while show low signal even when very low b-values (such as b=1 s/mm 2 ) is applied. Recently Wang [START_REF] Wang | Living tissue intravoxel incoherent motion (IVIM) diffusion MR analysis without b=0 image: an example for liver fibrosis evaluation[END_REF] proposed that liver vessel density can be measured by a diffusion weighted imaging derived surrogate biomarker (DDVD: diffusion derived vessel density): DDVD (b0b1) = Sb0/ ROIarea0 -Sb1/ROIarea1 [START_REF] Wang | Living tissue intravoxel incoherent motion (IVIM) diffusion MR analysis without b=0 image: an example for liver fibrosis evaluation[END_REF] where Sb0 refers to the measured liver signal intensity when b=0 s/mm 2 , and Sb1 refers to the measured liver signal intensity when b=1 s/mm 2 . ROIarea0 and ROIarea1 refer to the region-of-interest (ROI) on b=0 s/mm 2 and b=1 s/mm 2 images, respectively, and '/' means 'divided by'.

Sb1 and ROIarea1 can also be approximated by other low b-value diffusion image's data [START_REF] Xiao | Diffusion MRI Derived per Area Vessel Density as a Surrogate Biomarker for Detecting Viral Hepatitis B-Induced Liver Fibrosis: A Proof-of-Concept Study[END_REF][START_REF] Huang | Age and gender dependence of liver diffusion parameters and the possibility that intravoxel incoherent motion modeling of the perfusion component is constrained by the diffusion component[END_REF].

While we advocate sampling 16 b-values for liver IVIM imaging [START_REF] Wang | Topics on quantitative liver magnetic resonance imaging[END_REF][START_REF] Wang | Diffusion-weighted MRI of the liver: Challenges and some solutions for the quantification of apparent diffusion coefficient and intravoxel incoherent motion[END_REF], simplified IVIM protocols using less b-values are often adopted clinically, and very low b-value images are not always acquired. With a rat biliary duct ligation model, this study assessed the usefulness of liver DDVD measure computed from a simplified IVIM imaging protocol, alone and in combination with other IVIM parameters, to assess liver fibrosis. Instead of DDVD computed from the signal difference between on b=0 s/mm 2 and b =2~15 s/mm 2 images as we previously reported [START_REF] Wang | Living tissue intravoxel incoherent motion (IVIM) diffusion MR analysis without b=0 image: an example for liver fibrosis evaluation[END_REF][START_REF] Xiao | Diffusion MRI Derived per Area Vessel Density as a Surrogate Biomarker for Detecting Viral Hepatitis B-Induced Liver Fibrosis: A Proof-of-Concept Study[END_REF][START_REF] Huang | Age and gender dependence of liver diffusion parameters and the possibility that intravoxel incoherent motion modeling of the perfusion component is constrained by the diffusion component[END_REF],

in this study DDVD was computed from the signal difference between on b=0 s/mm 2 and b =30 or 50 s/mm 2 images.

Material and Methods

Animals:

Male Sprague-Dawley rats (weight: 300 ± 20g) were used. For bile duct blockage ligation (BDL) surgery, rats' lower segment of common bile duct ligated. Recanalization operation was performed with rats common bile duct isolated and ligated for 7 days. Via a laparotomy, the ligated distal end of the common bile duct was dissected, and an anastomosis was established between the common bile duct and the jejunum, which was sectioned 5 cm from the duodenojejunal angle [START_REF] Dos Santos | Influence of biliary anastomosis on recovery from secondary biliary cirrhosis[END_REF]. gives each parameter equal weighting. A 3-D coordinate system was initially constructed with PF, D slow and D fast , and then D fast or D slow was replaced with DDVD mean . With support vector machine (SVM) approach, the best separating plane between health livers and fibrotic livers was defined as A*X+B*Y+C*Z+D=0. In the 3-d coordinate system, the distance between point ( and plane A*X+B*Y+C*Z+D=0 was calculated according to following equation:

(

The mean distance of points representing healthy livers to the plane and mean distance of points representing fibrotic livers to the line were calculated, and then these two distances were added up.

In total, two control rats failed MRI due to movement during the examination, with motion artifacts and low signal-to-noise ratio. In addition, for one additional control rat, only the DDVD measurement was available, with IVIM measure unusable due to artifact on high b-value images. MRI of all surgery rats could be analyzed.

After the MRI, the rats had histological examination of the liver with Hematoxylin-eosin staining as well as picrosirius staining. Histology was assessed independently by a histopathologist, with diagnosis for liver fibrosis was based on the METAVIR criteria. A mean liver fibrosis score for each experimental rat group was assigned by the sum of fibrosis staging value divided by the number of the rat in the group.

Results

MRI morphologically, all control rats livers had normal appearance, while all rats with BDL had dilated biliary duct system. Of the rats with re-cannulation operation, biliary duct system remained dilated on MRI, with none showing apparent resolvement of biliary duct dilatation.

The histopathological grading for liver fibrosis is shown in table-1. For the re-cannulation group, diverse liver fibrosis stages were noted. However, the mean of liver fibrosis score of the recannulation group was higher than the week-1 group (2.42 vs 1.6), suggesting overall progress of liver fibrosis, rather than overall recovery, following the re-cannulation surgery. The control rats had a liver DDVD mean of 84.0±26.2 au/pixel, which was significantly higher than that of the rats with BDL (44.7±14.4 au/pixel, p<0.001) (Fig- 3).

Fig- 1 .

 1 Fig-1. The experimental procedures for rats. In total 8 control rats and 30 rats with surgery were studied. *: Two of the control rats failed the MRI examination.

Fig- 3 .

 3 Fig-3. DDVD measures of control rats and rats with different time points post-BDL (biliary duct ligation). A shows the results of DDVD (b0b25) . B shows the results of DDVD (b0b50) . C shows the correlation between DDVD (b0b25) and DDVD (b0b50) . D shows the results of DDVD mean . All: all rats with BDL (n=18).

Fig- 4 .

 4 Fig-4. Separation of healthy livers (n=5) and all fibrotic livers (n=30, BDL rats and BDL+ recanalization rats) by diffusion MRI metrics. A: PF alone cannot totally separate healthy livers and all fibrotic livers. B: a combination of PF, D slow , and D fast totally separates healthy livers and all fibrotic livers. C: compared with (B), a combination of PF, D slow , and DDVD mean even better separates healthy livers and all fibrotic livers.

Table - 1

 - . METAVIR score of the experimental rats livers with surgery intervention.

	BDL duration	rat number	stage-1 stage-2 stage-3	stage-4	mean fibrosis score
	7 days	5	2	3			1.6
	14 days	5		3	2		2.4
	21 days*	4		1	1	2	3.25
	28 days	3			1	2	3.7
	Re-cannulation	12	2	5	3	2	2.42
	* One rat in this group missed histology examination.

Funding: This work was supported by grants from the Shenzhen Basic Research Plan Fund Project (No. JCYJ20180228164641207) and the Cultivating Fund Project of Shenzhen People's Hospital (No. SYKYPY201933), and partially supported by a grant from the Research Grants Council of Hong Kong SAR (Project Nos. 2141061).

All quantitative data analysis was implemented in MATLAB (MathWorks, Natick, MA, USA). For DDVD/area measurement, ROI for liver parenchyma was segmented on the b=0 s/mm 2 image (resulting in ROI area of area0) and the b=25 (or 50) s/mm 2 image (resulting in ROI area of area25 or area50) respectively (Figure-2) [START_REF] Xiao | Diffusion MRI Derived per Area Vessel Density as a Surrogate Biomarker for Detecting Viral Hepatitis B-Induced Liver Fibrosis: A Proof-of-Concept Study[END_REF].

Two parameters were obtained as

(2)

The average of area weighted DDVDs of included liver axial slice's was the value of the liver. For each rat, measurement was conducted twice, and the mean of the two measurement was adopted as the final result. DDVD mean is then the mean of DDVD(b0b25) and DDVD(b0b50).

The bi-compartmental IVIM modeling is based on Equation-4:

where SI (b) and SI (0) denote the signal intensity of images acquired with the b-factor value of b and b=0 s/mm 2 , respectively. D slow (or D) is the diffusion coefficient representing the slow 'pure' molecular diffusion (unaffected by perfusion). IVIM analysis follow the principle of our previous studies for human subjects [START_REF] Wang | Living tissue intravoxel incoherent motion (IVIM) diffusion MR analysis without b=0 image: an example for liver fibrosis evaluation[END_REF][START_REF] Wang | Topics on quantitative liver magnetic resonance imaging[END_REF][START_REF] Wang | Diffusion-weighted MRI of the liver: Challenges and some solutions for the quantification of apparent diffusion coefficient and intravoxel incoherent motion[END_REF][START_REF] Chevallier | Removal of evidential motion-contaminated and poorly fitted image data improves IVIM diffusion MRI parameter scan-rescan reproducibility[END_REF]. The standard analysis was segmented fitting with a threshold b of 50 s/mm 2 [START_REF] Wang | Dependence of intravoxel incoherent motion diffusion MR threshold b-value selection for separating perfusion and diffusion compartments and liver fibrosis diagnostic performance[END_REF]. For each rat, measurement was conducted twice, and the mean of the two measurement was adopted as the final result.

Discussion

In this study, a BDL model and a BDL-recanalization model were applied, and liver fibrosis of varying severity was included as expected (table-1). A longer duration of BDL was associated with a higher mean liver fibrosis score. The BDL-recanalization model applied in this study has been described by dos Santos et al [START_REF] Dos Santos | Influence of biliary anastomosis on recovery from secondary biliary cirrhosis[END_REF]. This type of operation is associated with complications.

The propulsive power of intestine wall can add to the biliary duct pressure, leading to bowelbiliary reflux. Free entering of the enteric content to the biliary duct may induce episodes of mechanical obstruction, leading to cholestasis and cholangitis. These were confirmed in this study by the continued existence of dilated biliary duct anatomically. The mean of liver fibrosis score of the re-cannulation group was higher than the week-1 group (2.42 vs 1.5), suggesting overall progress of liver fibrosis following the re-recanalization surgery.

The liver fibrosis caused with biliary duct blockage confirms what we have observed in our human patient studies that a combination of PF, D slow , and D fast allows a total separation of all healthy livers and fibrotic livers (Fig- 4B). Thus, this study further supports our previous results that diffusion MRI can be used to detect even stage-1 fibrosis (four stage-1 fibrosis livers in this study) [START_REF] Wang | Topics on quantitative liver magnetic resonance imaging[END_REF][START_REF] Wang | Diffusion-weighted MRI of the liver: Challenges and some solutions for the quantification of apparent diffusion coefficient and intravoxel incoherent motion[END_REF]. This study also confirmed that the integration of DDVD into the IVIM dataset improves the total separation of healthy livers and livers.

Declaration of Conflicting Interests:

The authors declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.


Mirpouya Mirmozaffari
email: mirpouya.mirmozaffari@mavs.uta.edu

Seyed Elham Shadkam

Mohammad Khalili

Kamyar Kabirifar

Reza Yazdani

Tayyebeh Asgari Gashteroodkhani

Seyed Mohammad Khalili

A Novel Artificial Intelligent Approach: Comparison of Machine learning Tools and Algorithms Based on Optimization DEA Malmquist Productivity Index for Ecoefficiency Evaluation

Keywords: Artificial Intelligent, Machine learning, Optimization, Data Envelopment Analysis, Two-stage additive models, Eco-efficiency, Association rules, Classification

published or not. The documents may come

The construction industry is considered as a main contributor to environmental pollution and degradation [START_REF] Nematzadeh | Post-fire compressive strength of recycled PET aggregate concrete reinforced with steel fibers: Optimization and prediction via RSM and GEP[END_REF], Shahmansouri et al., 2020b). Climate change as one of the biggest challenges facing humankind [START_REF] Shahmansouri | Artificial neural network model to predict the compressive strength of eco-friendly geopolymer concrete incorporating silica fume and natural zeolite[END_REF] and as a direct consequence of construction activities, is happening faster than previously predicted (Yazdani et al., 2020b). According to the recent researches, world's climate is experiencing a significant increase in average annual temperatures (Mahmoudi et al., 2019a). Meanwhile, the atmospheric concentrations of some greenhouse gases (GHGs), particularly CO2, engendered by construction activities, exacerbates climate change [START_REF] Kabirifar | Construction and demolition waste management contributing factors coupled with reduce, reuse, and recycle strategies for effective waste management: A review[END_REF], Yazdani et al., 2020a).

Construction industry is the largest source of GHG emissions in the world, responsible for up to 50% of emissions [START_REF] Khasreen | Life-cycle assessment and the environmental impact of buildings: a review[END_REF]. Many countries in the world are experiencing the consequences of climate change, therefore they urgently need to reduce GHGs as part of a strong global effort [START_REF] Fathollahi-Fard | A green home health care supply chain: New modified simulated annealing algorithms[END_REF]. For instance, the Intergovernmental Panel on Climate Change (IPCC) has suggested that the world must cut off emissions to net zero by no later than 2050 to have a chance of limiting warming to 1.5 C. However, there is growing consensus that the adaptation efforts today are insufficient in minimizing the rate and extent of projected climatic changes.

By transitioning away from conventional construction practices to novel technologies, methods, and materials, the construction industry has considerable potential to reduce GHG emissions (Shahmansouri et al., 2020a). Although there are many researches in the construction industry have focused on assessing the industry's preparedness and capacity to adapt to climate change risks (e.g. [START_REF] Hurlimann | Is the Australian construction industry prepared for climate change?[END_REF]), there is relatively little research in using ecofriendly materials in construction industry. The materials used in construction sector such as cement release a considerable amount of carbon dioxide (CO2) [START_REF] Tekin | Recycling zeolitic tuff and marble waste in the production of eco-friendly geopolymer concretes[END_REF]. Therefore, replacing some materials with eco-friendly alternatives can help to reduce the CO2 emission.

Many countries are investing on developing infrastructures (Chen et al., 2019, Kabirifar andMojtahedi, 2019), and in these projects cement is one of the most used construction materials in the world (Akbarzadeh Bengar et al., 2020, Akbarzadeh Bengar andShahmansouri, 2020). Although PC is beneficial in many aspects, its production emits approximately 7% of the entire human-related CO2 released into the atmosphere [START_REF] Metz | Carbon dioxide capture and storage: special report of the intergovernmental panel on climate change[END_REF][START_REF] Taylor | Energy efficiency and CO2 emissions from the global cement industry[END_REF]. Of the GHG emissions associated with cement production, around 50% comes from the calcination (CaO formation via deriving CO2 from CaCO3), with the other 50% from the energy used within the procedure [START_REF] Yang | Assessment of CO2 reduction of alkali-activated concrete[END_REF][START_REF] Ávalos-Rendón | Synthesis of belite cements at low temperature from silica fume and natural commercial zeolite[END_REF]. Approximately four billion tons of PC is made each year, which makes it the primary binder in the building industry [START_REF] Pacheco-Torgal | Durability of alkali-activated binders: a clear advantage over Portland cement or an unproven issue?[END_REF]. Based on estimates, around 6 billion tons of PC will be used annually in the next four decades [START_REF] Samimi | Influence of pumice and zeolite on compressive strength, transport properties and resistance to chloride penetration of high strength self-compacting concretes[END_REF]. Therefore, developing new binders requiring less energy to produce and leading to less GHG emissions has become a necessity. Therefore, it is crucial to develop a method capable of analyzing and assessing the efficiency of ecological systems with regard to cement production in cement companies.

Searching for a novel machine learning and optimization approaches to find the most appropriate model is still an open problem [START_REF] Yazdani | Invasive weed optimization algorithm for minimizing total weighted earliness and tardiness penalties on a single machine under aging effect[END_REF], Yazdani and Jolai, 2016, Yazdani et al., 2017a). Thus, in responding to this research gap, in current study the combination of the above-mentioned machine learning and optimization approaches/methods will be a new experiment in literature perception. DEA as one of the main concepts of the optimization method has been widely used for calculating energy and environmental efficiency and eco-efficiency since it was primarily proposed by [START_REF] Charnes | Measuring the efficiency of decision making units[END_REF]. It is a non-parametric frontier technique where the efficiency of a specific entity is calculated by its distance from the highest performance practice frontier created by the most exceptional performance entities inside the group. DEA is a general method for assessing the efficiency of ecological systems [START_REF] Li | Carbon emission abatement quota allocation in Chinese manufacturing industries: An integrated cooperative game data envelopment analysis approach[END_REF], Mirmozaffari et al., 2020).

Data mining is world-class expertise with exceptional potential to aid the companies. The data mining methods such as classification forecast upcoming trends and performances. Data mining is a process of extracting valuable information from large amounts of data. The volume of data is growing dramatically as the data created by commonplace actions. Thus, taking out association rules from a large quantity of data in the database is interesting in numerous industries. Association rule mining has an excellent reputation in data mining. Apriori is the critical algorithm in association rule mining. Apriori algorithm is the most classical and essential algorithm for mining frequent item sets proposed by [START_REF] Agarwal | Fast algorithms for mining association rules[END_REF]. Apriori is used to find all routine item sets in each database. The Recently, several open-source data mining tools and software are accessible such as the Waikato Environment for Knowledge Analysis (WEKA), Rapid miner, and Tanagra. These tools and software offer a set of approaches and algorithms that support in better evaluation of data.

The analysis of Eco-efficiency is separated into two stages. The first stage uses its inputs to produce outputs, in which outputs are consequently considered as the inputs of the second stage. The first stage is reflected as the production stage and the second stage as the pollution control stage. In current study, all the available companies' datasets are applied to three exclusive models, and their DMU's efficiency is compared to find the unfamiliar trends in cement companies. A dataset for 22 companies with one input, two intermediate products, and three outputs after converting two-stage to proposed single-stage model is used. MPI in DEA with CCR, BCC, and additive Models are applied to test and justify the alterations between companies. The use of DEA as a decision analysis tool is limitless in literature because DEA does not focus on finding a universal relationship for all unit's undervaluation in the sample. DEA authorities every group in the data to have its production function, and then it evaluates the efficiency of that single unit by comparing it to the efficiency of the other units in the dataset. After running the three MPI models in DEA SOLVER with every group in the data, data mining classification algorithms in three strong data mining tools which are WEKA, Tanagra and Rapid Miner have been applied to find the superior tool and algorithm. Many classification algorithms should be applied to analyze and classify data. These algorithms differ in their nature and performance. Thus, in order to choose the best algorithm, a comparison studies required to present the most appropriate approach to be implemented in certain domain. The rest of the paper is organized as follows:

Section 2 reviews the literature linked to the calculation of Eco-efficiency, energy use efficiency, and specifies the prospective role of the current study. Section 3 presents a clarification of the data set description and data sources. Section 4 discusses the nine parts of the research methodology (part1: CCR, part2; BCC part3: additive models, part4: radial and non-radial models part5:Proposed model, part6: MPI, part7: Apriori algorithm in association rules, part8: Three strong classification algorithms in machine learning and the final calculation of accuracy, part9: Assessment process of combining DEA and data mining. Section 5 presents an evaluation and result of MPI and data mining algorithms respectively followed by a discussion and conclusion of the experimental consequences in Section 6 and 7 correspondingly.

Background and literature review

In this section, first, the existing studies relevant to this paper are reviewed, then research gaps and the main contributions of this study are discussed. [START_REF] Graham | Developing a social perspective to farm performance analysis[END_REF] is the first to apply an ecologically complex DEA Malmquist productivity index to find how ecological facilities effect measured productivity in the long run. The capability to include ecological influences lacking price data into MPI makes the index attractive for This research. [START_REF] Ødegaard | Measuring the contribution of workers' health and psychosocial work-environment on production efficiency[END_REF] combine MPI and bootstrap DEA to examine the influence of labor quality attributes in the direction of firm productivity. Finally, [START_REF] Arabi | Power industry restructuring and ecoefficiency changes: A new slacks-based model in Malmquist-Luenberger Index measurement[END_REF] present a novel SBM for Malmquist-Lueberger (ML) Index, and in the next year, [START_REF] Arabi | A new slacks-based measure of Malmquist-Luenberger index in the presence of undesirable outputs[END_REF] proposed the same model integrating undesirable outputs. Meanwhile, in many production procedures, large amounts of undesirable outputs are produced. The side effects of undesirable outputs on efficiency cannot be handled by standard MPI to measure productivity change over time. Consequently, A DEA model and an algorithm were introduced which can eliminate a common problem. This model includes the best way among all other guidelines to increase desirable output and decrease undesirable outputs simultaneously. Finally, a simple example used to demonstrate the novel algorithm, and the actual application of steam power plants is applied to display the applicability of the projected model. In the next year, (Emrouznejad andYang, 2016a, Emrouznejad andYang, 2016b) developed the previous model given its summing of slacks of good and bad outputs as the objective function of their models. [START_REF] Ghulam | Efficiency and productivity of the cement industry: Pakistani experience of deregulation and privatisation[END_REF] applied the MPI to study productivity in cement companies in Pakistan. Results indicated that privatization and lower degrees of organizational involvement in such companies had an advanced effect on the company's productivity. This [START_REF] Zhang | Total-factor carbon emission performance of the Chinese transportation industry: A bootstrapped non-radial Malmquist index analysis[END_REF] using the MPI to assess the performance of events containing CO2 in the transportation industry in China. The investigation took place in several periods, showing that the performance of the Chinese transportation industry had fallen by 32.8 percentage. This declining performance was attributed to a low level of technology in the field Applying the DEA Slack Based Model (SBM) under the reflection of constant returns to scale and variable return to scale technologies, a composite index is considered from the efficiency scores of each model. These models reflect both desirable and undesirable outputs.

MPI and non-radial model based on energy, environmental and industrial ecological efficiency

Eco-efficiency evaluation with desirable or undesirable inputs and outputs

Increasing concerns around the topic of energy safety and global warming, the problem of energy efficiency has gained significant attention from researchers. According to the International Energy Agency (IEA) energy efficiency is a way of handling and restricting the growth in energy consumption. Something is more energy-efficient than other methods if it distributes more services for the same energy input or the services for less energy input.

The number of DEA applications regarding pollutant services, as well as undesired outputs, is noteworthy. Separating ecological and technical efficiencies for power plants, [START_REF] Korhonen | Eco-efficiency analysis of power plants: An extension of data envelopment analysis[END_REF] suggested an approach. They consider pollutants as the inputs in order to increase desirable outputs and decrease pollutants and inputs. Gomes and Lins utilized Zero Sum Gains DEA to calculate carbon dioxide emissions. According to technical aspects, [START_REF] Yang | The necessity of distinguishing weak and strong disposability among undesirable outputs in DEA: Environmental performance of Chinese coal-fired power plants[END_REF] have considered several impossibilities features as undesirable outputs. In an experimental study, [START_REF] Zhang | Eco-efficiency analysis of industrial system in China: A data envelopment analysis approach[END_REF] used DEA to assess the eco-efficiency of gross domestic products in China. [START_REF] Liu | DEA models with undesirable inputs and outputs[END_REF] came across an approach to combine desirable and undesirable inputs and outputs. [START_REF] Chu | Analysis of China's regional eco-efficiency: A DEA two-stage network approach with equitable efficiency decomposition[END_REF] emphasized on the eco-efficiency study of Chinese provincial-level regions, concerning each region as a two-stage network structure. The first stage is reflected in the production system, and the second stage is considered as the pollution control system. Regarding the pollution emissions as intermediate products, a two-stage DEA model is suggested to attain the eco-efficiency of the entire two-stage structure. [START_REF] Khalili-Damghani | Uncertain network data envelopment analysis with undesirable outputs to evaluate the efficiency of electricity power production and distribution processes[END_REF] have considered emissions as an undesirable output in the efficiency assessment of electric power production and distribution procedure. [START_REF] Oggioni | Eco-efficiency of the world cement industry: a data envelopment analysis[END_REF] .used DEA to evaluate energy as an input producing both desirable outputs (goods) and undesirable outputs (CO2 emissions). The exclusion of undesirable output does not appear to deliver a broad scale of the production procedure. Consequently, [START_REF] Zhou | Linear programming models for measuring economy-wide energy efficiency performance[END_REF] assess energy use efficiency in a combined production context of both desirable and undesirable output.

To comprehend whether this eco-efficiency is attributable to a sensible consumption of inputs or an actual CO2 lessening as a result of the ecological rule, they evaluate the circumstances where CO2 emissions can either be reflected as an input or as an undesirable output. Practical effects display that countries, where cement industries spend in scientifically innovative kilns and adopt substitute fuels and raw materials in their manufacturing processes, are ecoefficient. [START_REF] Environment | Eco-efficient cements: Potential economically viable solutions for a low-CO2 cement-based materials industry[END_REF] proposed two essential approaches, which have a positive influence in substantial additional reductions in CO2 emissions and increased the use of the low-CO2 supplement, including the more efficient use of cement clinker. This efficient use contributes a relative advantage to developed countries, such as India and China, which are encouraged to renovate their production processes. In terms of CO2 emissions, DEA, and ecoefficiency, some outstanding recent studies can be addressed in the literature [START_REF] Grigoroudis | Evaluation of National Environmental Efficiency Under Uncertainty Using Data Envelopment Analysis[END_REF][START_REF] Wang | Study on the decomposition of factors affecting energy-related carbon emissions in Guangdong province[END_REF][START_REF] Almutairi | Modeling, analysis, and evaluation of a carbon tax policy based on the emission factor[END_REF][START_REF] Bi | A clustering method for evaluating the environmental performance based on slacks-based measure[END_REF], Zeng et al., 2016a[START_REF] Li | Impact of energy conservation policies on the green productivity in China's manufacturing sector: Evidence from a three-stage DEA model[END_REF][START_REF] Jiang | Impact of carbon quota allocation mechanism on emissions trading: An agentbased simulation[END_REF][START_REF] Zeng | Forecasting the allocation ratio of carbon emission allowance currency for 2020 and 2030 in China[END_REF][START_REF] Zeng | Forecasting the allocative efficiency of carbon emission allowance financial assets in China at the provincial level in 2020[END_REF], Mahmoudi et al., 2019b[START_REF] Jin | Determining the optimal carbon tax rate based on data envelopment analysis[END_REF][START_REF] Khoshroo | Improving energy efficiency considering reduction of CO2 emission of turnip production: A novel data envelopment analysis model with undesirable output approach[END_REF][START_REF] Arabi | Eco-efficiency measurement and material balance principle: an application in power plants Malmquist Luenberger Index[END_REF][START_REF] Boyd | Estimating the linkage between energy efficiency and productivity[END_REF][START_REF] Wei | An empirical analysis of energy efficiency in China's iron and steel sector[END_REF][START_REF] Azadeh | An integrated DEA PCA numerical taxonomy approach for energy efficiency assessment and consumption optimization in energy intensive manufacturing sectors[END_REF][START_REF] Angulo-Meza | A multiobjective DEA model to assess the eco-efficiency of agricultural practices within the CF+ DEA method[END_REF][START_REF] Balitskiy | Energy efficiency and natural gas consumption in the context of economic development in the European Union[END_REF][START_REF] Bian | Energy efficiency analysis of the economic system in China during 1986-2012: a parallel slacks-based measure approach[END_REF][START_REF] Moya | Analysis of the Ecuadorian energy audit practices: a review of energy efficiency promotion[END_REF][START_REF] Mirmozaffari | Eco-Efficiency Evaluation in Two-Stage Network Structure: Case Study: Cement Companies[END_REF][START_REF] Lee | Eco-efficiency evaluation considering environmental stringency[END_REF][START_REF] Yu | Provincial allocation of carbon emission reduction targets in China: an approach based on improved fuzzy cluster and Shapley value decomposition[END_REF].

Machine learning algorithms in energy consumption

Making effective decisions is crucial in any field [START_REF] Khalili | A Bi-objective Model for Relief Supply Location in Post-Disaster Management[END_REF][START_REF] Aghapour | Capacity planning and reconfiguration for disaster-resilient health infrastructure[END_REF]. Innovative computational methods particularly machine learning and soft computing techniques have the potential to tackle a wide range of challenging problems (Yazdani et al., 2016, Yazdani et al., 2017b[START_REF] Yazdani | A modified imperialist competitive algorithm for a two-agent single-machine scheduling under periodic maintenance consideration[END_REF][START_REF] Yazdani | A genetic algorithm with modified crossover operator for a two-agent scheduling problem[END_REF][START_REF] Shahmansouri | Predicting compressive strength and electrical resistivity of eco-friendly concrete containing natural zeolite via GEP algorithm[END_REF][START_REF] Golilarz | Satellite image de-noising with harris hawks meta heuristic optimization algorithm and improved adaptive generalized gaussian distribution threshold function[END_REF], Golilarz et al., 2020), therefore they have widely been applied in different fields in recent years [START_REF] Yazdani | Reliability estimation using an integrated support vector regression-variable neighborhood search model[END_REF][START_REF] Azadeh | An integrated support vector regressionimperialist competitive algorithm for reliability estimation of a shearing machine[END_REF]. [START_REF] Yu | Provincial allocation of carbon emission reduction targets in China: an approach based on improved fuzzy cluster and Shapley value decomposition[END_REF] according to the relevant carbon emission features and decided that more carbon emission reduction quantity would be assigned to regions with large total emissions and high emission concentration. [START_REF] Emrouznejad | A novel inverse DEA model with application to allocate the CO2 emissions quota to different regions in Chinese manufacturing industries[END_REF] considered the same problem of this article, despite the fact that they applied an inverse DEA method and ignored the competitive and supportive relations between various sub-level trades and regions. [START_REF] Qing | A real-time monitoring method of energy consumption based on data mining[END_REF] extracted building energy consumption data by DBSCAN clustering and decision tree-based classification methods. Despite its capability to intensely comprehend the outline of energy consumption in constructions, the algorithm is too intricate to be appropriate for the rapid processing of data in the energy consumption monitoring platform.

Three strong data mining tools (WEKA, Rapid miner, and TANAGRA)

In this part, three strong data mining tools are introduced (Naik and Samant,2016):  WEKA is an extensively implemented tool for machine learning that was initially established at the University of Waikato in New Zealand. It includes many groups of advanced machine learning algorithms written in Java. This tool covers tools for regression, classification, clustering, association rules, visualization, and data pre-processing.

 Rapid miner (RM) is a user collaborative tool for machine learning procedures. It is written in Java. It shows an integrated method to design even very intricate problems -an integrated operative perception which lets the design of compound nested operator chains for a large number of learning problems. RM applied XML to define the operative trees modeling knowledge discovery (KD) procedures. It has flexible workers for data input and output in diverse file formats. It covers more than 100 learning systems for classification, regression and clustering algorithms.

 Tanagra is a free suite of machine learning software to study and theoretical methods established by [START_REF] Rakotomalala | TANAGRA: un logiciel gratuit pour l'enseignement et la recherche[END_REF] at the Lumière University, France. It provides numerous standard machine learning approaches such as: Visualization, Descriptive statistics, Instance selection, feature selection, feature construction, regression, factor analysis, clustering, classification and association rule learning. Tanagra supports an appropriate cooperation among the statistical methods (e.g. parametric and nonparametric statistical tests), the multivariate analysis approaches (e.g. factor analysis, correspondence analysis, cluster analysis, regression) and the machine learning techniques (e.g. neural network, support vector machine, decision trees, random forest).

Research gap analysis and contributions

DEA measure of energy use efficiency has two main benefits as compared to the old-style meaning, "the percentage of energy services to energy input." At the initial step, DEA provides somewhere to stay several inputs (energy and noenergy inputs) and several outputs (desirable and undesirable) in the production progression. Next, DEA can also put up the purposes of DMUs in evaluating energy use efficiency. However, far too little attention has been paid to developing an efficient solution method to cope with DEA and data mining associated with finding the superior model, algorithm and DMU. In conclusion, main contributions of this study can be summarized as follows:

 This research aims to study a comprehensive comparison of several efficiencies delivers insight into the firm's efficiency based on machine learning combined with innovative optimization models for Eco-efficiency Evaluation. This comparison is of considerable significance to cement companies' practitioners who desire to assess productivity and efficiency at a proper step of its progression.

 An exclusive, effective and easy to implement converting two-stage to single stage models in MPI are applied with DEA SOLVER, which eventually results in comparing several efficient and inefficient DMUs. This model is proposed to fix the efficiency of a two-stage process and prevent the dependency on various weights. In fact, by converting two-stage to single proposed stage, desirable and undesirable inputs and output can be evaluated with simple CCR, BCC and additive suggested model.

 After applying the abovementioned optimization part, best machine learning tool and algorithm are introduced.

DEA inputs and outputs as potential attributes for data mining classification algorithms in WEKA, Tanagra and Rapid miner are considered. In addition, data play the role of instances and finally efficient DMUs are applied for class yes and inefficient DMUs for class no. model, company, algorithm, and tool will be introduced. Thus, it can be beneficial for managers to remove unrelated data and conduct more effective processes.

Dataset description

The standard data set, collected in this study covers five years (2015-2019), which is collected from 22 cement companies. Consequently, the single input of the first stage, two intermediate elements, and three outputs of the second stage for the first company, which is 2015 to 2019, are presented in Table 1.

Research methodology

The objective of this study after converting two-stage to a single-stage is to compare companies' efficiency effectively.

Using a comparative DEA with window analysis methodology and clustering algorithms in data mining were established to determine the features of cement companies in terms of some DMUs and algorithms. Finally, the entire progression can be divided into five steps, as follows:

CCR Model:

The CCR models show a constant return to scale (CRS), which means that relative progress in all inputs is equal to the increase in outputs. The efficiency of an assumed DMU is calculated based on the model as follows:

CCRIO 𝑀𝑖𝑛𝜃 𝑆𝑡. ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑥 𝑖𝑗 ≤ 𝜃 𝑝 ,𝑖 = 1,…….,𝑚 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑦 𝑟𝑗 ≥ 𝑦 𝑟𝑝 , 𝑟 = 1,…….,𝑠 𝜆 𝑗 ≥ 0 ,𝑗 = 1,…,𝑛
Eq. 1

Where 𝜃𝑝 specifies the technical efficiency score of units DMU, λj indicates the dual variables that categorize the benchmarks for inefficient parts. If 𝜃𝑝is equal to one, then the DMU shows a technically efficient unit. It is on the efficiency frontier, which is collected from the set of efficient units. DEA calculates the efficiency of each observation based on the frontier that covers all the views. Inefficient DMUs can be improved (moved to the efficient frontier) with strategic directions for precision, which are the points along the frontier. The distance to the efficiency frontier distributes an amount of efficiency.  Pure Technical Efficiency (PTE): PTE provides the effect of scale size by just connecting a DMU to a unit of comparable scale and procedures on how a DMU develops its bases under the outer region.

4.2.BCC

 Scale Efficiency (SE): SE, shows how the scale size has positive or negative influences on efficiency. If after using both CRS and VRS models on similar data, there is a change in the two technical efficiencies, and based on these tags, DMU has a scaling efficiency and can be designed by:

SE =TE/ PTE Eq. 2
The is represented as follows:

BCCIO

𝑀𝑖𝑛𝜃 𝑆𝑡. ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑥 𝑖𝑗 ≤ 𝜃 𝑝 ,𝑖 = 1,…….,𝑚 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑦 𝑟𝑗 ≥ 𝑦 𝑟𝑝 , 𝑟 = 1,…….,𝑠 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 = 1 𝜆 𝑗 ≥ 0 ,𝑗 = 1,…,𝑛
Eq. 3

4.3.Additive Model (AM)

This model reflects the input surplus and the output shortage concurrently with incoming at a point on the efficient frontier. In addition, as a non-radial model, when dealing with undesirable outputs, slacks in undesirable outputs are accounted for the efficiency measure. In mathematical terms the Constant Return to Scale (CRS) for AM model or AM-CRS is given as:

𝑀𝑖𝑛𝑍 = ∑ 𝑚 𝑖 = 1 𝑆 - 𝑖 + ∑ 𝑠 𝑟 = 1 𝑆 + 𝑟 𝑆𝑡. ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑥 𝑖𝑗 + 𝑆 - 𝑖 = 𝑥 𝑖 ,𝑖 = 1,…….,𝑚 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑦 𝑟𝑗 -𝑆 + 𝑟 = 𝑦 𝑟 , 𝑟 = 1,…….,𝑠 𝜆 𝑗 ,𝑆 + 𝑟 , 𝑆 - 𝑖 ≥ 0 ,𝑗 = 1,…,𝑛
Eq. 4

After adding the following constrain to the abovementioned AM-CRS or relation number

∑ 𝑛 𝑗 = 1 𝜆 𝑗 = 1
Eq. 5

We can convert CRS to Variable Return to Scale (VRS) model. This is exactly like the converting CCR to BCC model. In technical relations, the AM-VRS model is given as:

𝑀𝑖𝑛𝑍 = ∑ 𝑚 𝑖 = 1 𝑆 - 𝑖 + ∑ 𝑠 𝑟 = 1
𝜆 𝑗 = 1 𝜆 𝑗 ,𝑆 + 𝑟 , 𝑆 - 𝑖 ≥ 0 ,𝑗 = 1,…,𝑛

Radial and non-radial Models

The main AM objective, as a radial model, is to achieve the highest value of the input and output slacks in the given DMUs. DEA measures the efficiency of each reflection relative to the frontier that covers all the explanations. Inefficient

DMUs can be value-added (moved to the efficient frontier) with proposed advice for development, which are the points along the frontier. The distance to the efficiency frontier delivers a measure of efficiency. To clarify more, Figure 1 shows the main differences between radial and no-radial models with two inputs and one output in Farrell efficiency border.

Fig. 1. Radial and non-radial models with two inputs and one output

In Figure 1, the efficient frontier represented by the solid line passes through A, , C and no other points. Therefore,

F ′ F ″ ,
for ,in order to be efficient DMUs, it is needed to bring each DMU onto the efficient frontier. Thus, we want to F ′ , F ″ bring F DMU onto efficient frontier based on radial models' rule. The value of which is available in the following 𝜃 * equation should be applied to inputs and finally will bring it to (onto AC efficient frontier).

𝐹

F ′ 𝜃 * = 𝑂𝐷 𝑂D ′ = 𝑂𝐵 𝑂𝐵 ′
Eq. 7

Based on the above calculation, represents radial efficient DMU. In fact, in radial models, the reducing amounts of F ′ inputs and outputs should be the same. On the other hand, in non-radial models, it is not necessary to apply the same reducing amounts for inputs or outputs and this is one of the advantages of radial models. For non-radial or AM models, based on figure 1, the following two separate steps should be considered:

1. Reducing the first input by at the first step 𝑆 1 -= 𝐹𝐿 2. Reducing the second input by at the second step 𝑆 2 -= LF ″ Based on the above calculation, shows non-radial efficient DMU. The major concern of input (output)-oriented F ″ models focuses on the input (output)-side efficiency, and output (input)-side is a minor subject in measuring efficiency.

In addition, only the non-radial models can capture the whole features of efficiency. Thus, in order to get the whole advantages of efficiency, AM model has been used in this study.

4.5.1.A new approach in DEA two-stage model

Many papers have proposed a novel two-stage model in DEA (Chen et al., 2009a, Chen et al., 2009b[START_REF] Kao | Efficiency decomposition in two-stage data envelopment analysis: An application to non-life insurance companies in Taiwan[END_REF][START_REF] Lotfi | A new two-stage data envelopment analysis (DEA) model for evaluating the branch performance of banks[END_REF]. . Decreasing pollutant material investment, as well as energy consumption (inputs of the single proposed model) and increasing waste material removed as well as cement production (outputs of the unique proposed model) have a positive influence on profits, and efficiency of cement companies.

4.5.2.Primal model in

: 𝑪𝑪𝑹 𝑰𝑶 𝑀𝑎𝑥ɸ = ∑ 𝑠 𝑟 = 1 𝑌 𝑟 𝑢 𝑟 + ∑ 𝑑 ℎ = 1 𝑀 ℎ 𝑙 ℎ + ∑ 𝑣 𝑡 = 1 𝐸 𝑡 𝑏 𝑡 + ∑ 𝑞 𝑧 = 1 𝐹 𝑧 𝑎 𝑧 𝑆𝑡. ∑ 𝑚 𝑖 = 1 𝑋 𝑖 𝑣 𝑖 + ∑ 𝐾 𝑐 = 1 𝑁 𝑐 𝑔 𝑐 = 1 ∑ 𝑠 𝑟 = 1 𝑌 𝑟 𝑢 𝑟 + ∑ 𝑑 ℎ = 1 𝑀 ℎ 𝑙 ℎ + ∑ 𝑣 𝑡 = 1 𝐸 𝑡 𝑏 𝑡 + ∑ 𝑞 𝑧 = 1 𝐹 𝑧 𝑎 𝑧 -∑ 𝑚 𝑖 = 1 𝑋 𝑖 𝑣 𝑖 -∑ 𝐾 𝑐 = 1 𝑁 𝑐 𝑔 𝑐 ≤ 0 𝑢 𝑟 ,𝑙 ℎ , 𝑣 𝑖 ,𝑔 𝑐 ,𝑏 𝑡 ,𝑎 𝑧 ≥ 0 , 𝑗 = 1,…,𝑛
Eq. 8

4.5.3.Primal model in

: 𝑩𝑪𝑪 𝑰𝑶 𝑀𝑎𝑥ɸ = ∑ 𝑠 𝑟 = 1 𝑌 𝑟 𝑢 𝑟 + ∑ 𝑑 ℎ = 1 𝑀 ℎ 𝑙 ℎ + ∑ 𝑣 𝑡 = 1 𝐸 𝑡 𝑏 𝑡 + ∑ 𝑞 𝑧 = 1 𝐹 𝑧 𝑎 𝑧 + 𝑤 𝑆𝑡. ∑ 𝑚 𝑖 = 1 𝑋 𝑖 𝑣 𝑖 + ∑ 𝐾 𝑐 = 1 𝑁 𝑐 𝑔 𝑐 = 1 ∑ 𝑠 𝑟 = 1 𝑌 𝑟 𝑢 𝑟 + ∑ 𝑑 ℎ = 1 𝑀 ℎ 𝑙 ℎ + ∑ 𝑣 𝑡 = 1 𝐸 𝑡 𝑏 𝑡 + ∑ 𝑞 𝑧 = 1 𝐹 𝑧 𝑎 𝑧 -∑ 𝑚 𝑖 = 1 𝑋 𝑖 𝑣 𝑖 -∑ 𝐾 𝑐 = 1 𝑁 𝑐 𝑔 𝑐 + 𝑤 ≤ 0 𝑢 𝑟 ,𝑙 ℎ , 𝑣 𝑖 ,𝑔 𝑐 ,𝑏 𝑡 ,𝑎 𝑧 ≥ 0 , 𝑗 = 1,…,𝑛
Eq. 9

4.5.4.Primal proposed model in

= -∑ 𝑚 𝑖 = 1 𝑆 𝑖 --∑ 𝑘 𝑐 = 1 𝑆 - 𝑐 -∑ 𝑠 𝑟 = 1 𝑆 + 𝑟 -∑ 𝑑 ℎ = 1 𝑆 + ℎ -∑ 𝑣 𝑡 = 1 𝑆 + 𝑡 -∑ 𝑞 𝑧 = 1 𝑆 + 𝑧 𝑆𝑡. ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑋 𝑖𝑗 + 𝑆 - 𝑖 = 𝑋 𝑖 ,𝑖 = 1,…,𝑚 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑁 𝑐𝑗 + 𝑆 - 𝑐 = 𝑁 𝑐 ,𝑐 = 1,…,𝑘 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑌 𝑟𝑗 -𝑆 + 𝑟 = 𝑌 𝑟 ,𝑟 = 1,…,𝑠 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑀 ℎ𝑗 -𝑆 + ℎ = 𝑀 ℎ ,ℎ = 1,…,𝑑 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝐸 𝑡𝑗 -𝑆 + 𝑡 = 𝐸 𝑡 ,𝑡 = 1,…,𝑣 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝐹 𝑧𝑗 -𝑆 + 𝑧 = 𝐹 𝑧 ,𝑧 = 1,…,𝑞 ∑ 𝑚 𝑗 = 1 𝜆 𝑗 = 1 𝜆 𝑗 , 𝑆 - 𝑖 , 𝑆 - 𝑐 , 𝑆 + 𝑟 , 𝑆 + ℎ , 𝑆 + 𝑡 ,𝑆 + 𝑧 ≥ 0, 𝑗 = 1,…,𝑛
Eq. 10

Dual proposed model in

: 𝑪𝑪𝑹 𝑰𝑶 𝑀𝑖𝑛 𝜃 𝑆𝑡. ∑ 𝑠 𝑟 = 1 𝑌 𝑟 𝑢 𝑟 + ∑ 𝑑 ℎ = 1 𝑀 ℎ 𝑙 ℎ + ∑ 𝑣 𝑡 = 1 𝐸 𝑡 𝑏 𝑡 + ∑ 𝑞 𝑧 = 1 𝐹 𝑧 𝑎 𝑧 -∑ 𝑚 𝑖 = 1 𝑋 𝑖 𝑣 𝑖 -∑ 𝐾 𝑐 = 1 𝑁 𝑐 𝑔 𝑐 ≤ 0 , 𝑗 = 1,…, 𝑛 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑋 𝑖𝑗 ≤ 𝜃 𝑝 𝑋 𝑖𝑝 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑁 𝑐𝑗 ≤ 𝜃 𝑝 𝑁 𝑐𝑝 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑌 𝑟𝑗 ≥ 𝑌 𝑟𝑝 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑀 ℎ𝑗 ≥ 𝑀 ℎ𝑝 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝐸 𝑡𝑗 ≥ 𝐸 𝑡𝑝 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝐹 𝑧𝑗 ≥ 𝐹 𝑧𝑝 𝜆 𝑗 ≥ 0 𝜃 𝑝 𝑓𝑟𝑒𝑒
Eq. 11

4.5.6.Dual proposed model in

: 𝑩𝑪𝑪 𝑰𝑶 𝑀𝑖𝑛 𝜃 𝑆𝑡. ∑ 𝑠 𝑟 = 1 𝑌 𝑟 𝑢 𝑟 + ∑ 𝑑 ℎ = 1 𝑀 ℎ 𝑙 ℎ + ∑ 𝑣 𝑡 = 1 𝐸 𝑡 𝑏 𝑡 + ∑ 𝑞 𝑧 = 1 𝐹 𝑧 𝑎 𝑧 -∑ 𝑚 𝑖 = 1 𝑋 𝑖 𝑣 𝑖 -∑ 𝐾 𝑐 = 1 𝑁 𝑐 𝑔 𝑐 + 𝑤 ≤ 0 , 𝑗 = 1,…, 𝑛 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑋 𝑖𝑗 ≤ 𝜃 𝑝 𝑋 𝑖𝑝 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑁 𝑐𝑗 ≤ 𝜃 𝑝 𝑁 𝑐𝑝 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑌 𝑟𝑗 ≥ 𝑌 𝑟𝑝 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 𝑀 ℎ𝑗 ≥ 𝑀 ℎ𝑝
𝑌 𝑟 𝑢 𝑟 + ∑ 𝑑 ℎ = 1 𝑀 ℎ 𝑙 ℎ + ∑ 𝑣 𝑡 = 1 𝐸 𝑡 𝑏 𝑡 + ∑ 𝑞 𝑧 = 1 𝐹 𝑧 𝑎 𝑧 -∑ 𝑚 𝑖 = 1 𝑋 𝑖 𝑣 𝑖 -∑ 𝐾 𝑐 = 1 𝑁 𝑐 𝑔 𝑐 , 𝑗 = 1,…, 𝑛 𝑆𝑡. ∑ 𝑠 𝑟 = 1 𝑌 𝑟 𝑢 𝑟 + ∑ 𝑑 ℎ = 1 𝑀 ℎ 𝑙 ℎ + ∑ 𝑣 𝑡 = 1 𝐸 𝑡 𝑏 𝑡 + ∑ 𝑞 𝑧 = 1 𝐹 𝑧 𝑎 𝑧 -∑ 𝑚 𝑖 = 1 𝑋 𝑖 𝑣 𝑖 -∑ 𝐾 𝑐 = 1 𝑁 𝑐 𝑔 𝑐 + 𝑤 ≤ 0 , 𝑗 = 1,…, 𝑛 ∑ 𝑚 𝑖 = 1 𝑣 𝑖 ≥ 1 ∑ 𝑘 𝑐 = 1 𝑔 𝑐 ≥ 1 ∑ 𝑐 𝑟 = 1 𝑢 𝑟 ≥ 1 ∑ 𝑑 ℎ = 1 𝑙 ℎ ≥ 1 ∑ 𝑣 𝑡 = 1 𝑏 𝑡 ≥ 1 ∑ 𝑞 𝑧 = 1 𝑎 𝑧 ≥ 1 ∑ 𝑛 𝑗 = 1 𝜆 𝑗 = 1 𝜆 𝑗 ≥ 0 𝜃 𝑝 𝑓𝑟𝑒𝑒
Eq. 13

Evaluation in MPI

The MPI is measured to assess productivity growth based on the reference technology. The following two main topics are used in the calculation of MPI development:



The first issue is that the number of productivity changes over the period.



On the other hand, second is to decompose changes in productivity into what are generally denoted as a catching-up result or technical efficiency change (TEC) and a frontier shift result or technological change.

MPI measures the total factor productivity change of a DMU between two periods. The idea of productivity is typically denoted as labor productivity. This idea is related to TFP, defined as the product of efficiency change (catch-up) and technological change (frontier-shift). If TFP value is more than one, this specifies a progressive TFP growth from period (t) to period (t+1). However, a value of less than one yields a reduction in TFP development or performance relative to the preceding year. In fact:

Malmquist Productivity Index (MPI) = TEC ×TC Eq. 14

The MPI can be stated through distance function(E) based on the two following equations applying the reflections at time t and t+1 [START_REF] Lee | Malmquist Productivity Index using DEA frontier in Stata[END_REF]. Eq. 15

𝑀𝑃𝐼 𝑡 + 1 𝐼 = 𝐸 𝑡 + 1 𝐼 (𝑥 𝑡 + 1 + 𝑦 𝑡 + 1 ) 𝐸 𝑡 + 1 𝐼 (𝑥 𝑡 + 𝑦 𝑡 )
Eq. 16

The geometric mean of two MPI in equations ( 14) and ( 15) provides the equation ( 16):

𝑀𝑃𝐼 𝐺 𝐼 = (𝑀𝑃𝐼 𝑡 𝐼 𝑀𝑃𝐼 𝑡 + 1 𝐼 ) 1/2 = [ ( 𝐸 𝑡 𝐼 (𝑥 𝑡 + 1 + 𝑦 𝑡 + 1 ) 𝐸 𝑡 𝐼 (𝑥 𝑡 + 𝑦 𝑡 ) ) . ( 𝐸 𝑡 + 1 𝐼 (𝑥 𝑡 + 1 + 𝑦 𝑡 + 1 ) 𝐸 𝑡 + 1 𝐼 (𝑥 𝑡 + 𝑦 𝑡 ) ) ] 1/2
Eq. 17

Input oriented geometric mean of MPI can be disintegrated using the abovementioned input-oriented TC and input oriented TEC as follows:

𝑀𝑃𝐼 𝐺 𝐼 = (𝑇𝐶 𝐼 )(𝑇𝐸𝐶 𝐺 𝐼 ) 1/2 = 𝐸 𝑡 + 1 𝐼 (𝑥 𝑡 + 1 + 𝑦 𝑡 + 1 ) 𝐸 𝑡 𝐼 (𝑥 𝑡 + 𝑦 𝑡 ) [ ( 𝐸 𝑡 𝐼 (𝑥 𝑡 + 𝑦 𝑡 ) 𝐸 𝑡 + 1 𝐼 (𝑥 𝑡 + 𝑦 𝑡 ) ) . ( 𝐸 𝑡 𝐼 (𝑥 𝑡 + 1 + 𝑦 𝑡 + 1 ) 𝐸 𝑡 + 1 𝐼 (𝑥 𝑡 + 1 + 𝑦 𝑡 + 1 ) ) ] 1/2
Eq. 18 MPI specified by equations ( 16) and ( 17) can be well-defined applying DEA, such as distance function. These are the elements of MPI which can be driven from the assessment of distance functions defined on frontier technology. Thus, the abovementioned MPI is the most common method between the numerous techniques that have been established

to evaluate a production technology (Färe et al., 1994, Chen and[START_REF] Chen | DEA Malmquist productivity measure: New insights with an application to computer industry[END_REF].

Evaluation in Apriori algorithm in association rules

After inserting the attributes and class by the WEKA software, the Apriori algorithm in association rules is applied. To select strong rules from the set of all possible rules, limitations on several methods of standing, and interest are used.

The best-known restrictions are minimum commencements on support and confidence. Support is a sign of how recurrently the item set looks as if in the dataset. For example, the rule has the support of one in a data set, since it happens in 100% of all transactions (one out of one transaction). Confidence is a sign of how often the rule has been found to be precise. For example, the rule has the confidence of one in a data set, since, for 100% of the transactions, the rule is correct. If the rule had a lift of one, it would propose that the scene of being and the subsequent are independent of each other. When two procedures are independent, no rule can be drawn connecting those two events.

If the lift is more than 1, two procedures are independent of one another and make those rules hypothetically respected for forecasting the following in upcoming data sets. Under this condition that the lift is greater than one, the more the support and confidence of a rule in the Apriori algorithm, the higher it indicates a fixed outline in the dataset. If these procedures are relatively small, then any indiscretion would be less strong than it would be for rules with high confidence, great support, and high lift. After applying distinct data sets, association rules are assessed, and strong rules are extracted.

Evaluation in classification algorithms in machine learning

In this part three strong classification algorithms are introduced, and the basic calculation pf accuracy is presented:

K-Nearest Neighbor: K-nearest neighbors' algorithm (k-NN) is a non-parametric technique applied for classification and regression. In both classification and regression, the input contains the k closest training instances in the feature space. The output is a class member for classification and is an average of neighbor's value for regression.

Decision tree: Decision trees are generally implemented in operations research and decision analysis, to support classify an approach in order to reach a goal. In addition, it is a strong tool in machine learning. The flow chart structure of decision tree contains branch and leaf. Branch shows the result of the test and leaf displays the class label. The efficiency of the classification algorithms based on the following accuracy formula was compared:

𝑇𝑟𝑢𝑒𝑝𝑜𝑠𝑖𝑡𝑖𝑣𝑒 + 𝑇𝑟𝑢𝑒𝑛𝑒𝑔𝑎𝑡𝑖𝑣𝑒 𝑇𝑟𝑢𝑒𝑝𝑜𝑠𝑖𝑡𝑖𝑣𝑒 + 𝑇𝑟𝑢𝑒𝑛𝑒𝑔𝑎𝑡𝑖𝑣𝑒 + 𝐹𝑎𝑙𝑠𝑒𝑝𝑜𝑠𝑖𝑡𝑖𝑣𝑒 + 𝐹𝑎𝑙𝑠𝑒𝑛𝑒𝑔𝑎𝑡𝑖𝑣𝑒

Result in Apriori algorithm

For the last part, three strong rules based on pseudocode have been proposed:

1.Waste water removed = If "(-inf-134021.2]" 33 ==> Then, Cement production ="(-inf-1964348.9]" 32

Based on the first strong rule if wastewater removed within the quantified range happens 33 times (33 out of 110), cement production within the quantified range will occur 32 times (32 out of 110).

2.Energy consumption= If "(952946.9-1549674.8]" 26 ==>Then, Cement production ="(-inf-1964348.9]" 24

Based on the second strong rule if energy consumption within the quantified range happens 26 times (26 out of 110), cement production within the quantified range will occur 24 times (24 out of 110).

Result in classification algorithms

Conclusion

The performance of the proposed approach provides us with a chance to recognize pattern recognition of the whole, combining DEA and data mining techniques during the selected period (five years from 2015 to-2019). Meanwhile, the cement industry is one of the foremost manufacturers of naturally harmful material using an undesirable by-product; specific stress is given to that pollution control investment or undesirable output while evaluating energy use efficiency. The significant concentration of the study is to respond to five preliminary questions. First, whether the conversion of the two-stage model to a novel single stage has any positive impact on Eco-efficiency? Secondly, among three proposed models (CCR, BCC, and AM models), which model has the highest average efficiency over 5-year periods for 22 DMUs? Thirdly, whether AM proposed model has any corroboration impact on Eco-efficiency? Fourthly, whether combining DEA and data mining have any strengthening influence on Eco-efficiency? and finally, whether comparison of three strong classification algorithms in three strong machine learning tools can be useful for future

studies?

To answer the first question, this model is proposed to fix the efficiency of a two-stage process and prevent the dependency on various weights. Converting undesirable outputs, and desirable inputs to final desirable inputs in a single-stage model to minimize inputs as well as turning desirable outputs to final desirable outputs in the single stage model to maximize outputs to have a positive effect on the efficiency of the whole process. Decreasing pollutant material investment, as well as energy consumption (inputs of the single proposed model) and increasing waste material removed as well as cement production (outputs of the unique proposed model) have a positive influence on profits, and efficiency of cement companies. To answer the second question, the AM model has a more strengthening impact on Eco-efficiency compare with other models such as CCR and BCC input-oriented and output-oriented. To answer the third question, the AM model decreases the value of inputs and increases the value of outputs at the same time. To answer the fourth question, after applying previous steps in MPI, the data mining Apriori algorithm extracts the strongest rules and remove week rules. Finally, for the last and the fifth question, in order to implement the best tool and algorithm, a comparison of strong algorithms and tools will be beneficial for future studies. In conclusion, it is worth declaring the Limitations of this research. The most critical limitation is associated with the incomplete contemplation of undesirable output. The current study measured only pollution control investment as an undesirable output produced by the cement industry, and aforesaid unique methods have been proposed to decrease these undesirable outputs. However, other undesirable severe outputs also could have been combined into the analysis if reliable data were accessible. Containing these variables may affect the efficiency marks of this study. However, the general inference remains the same as the proposed combining method in this study, can put up several inputs and outputs. Along with the results obtained from efficiency analysis, the managers of the 17th cement company have the highest efficiency in cement production centers during the five years. They should try to have better efficiency in the future. Managers of other companies with less efficiency need to use specific data mining Apriori algorithm to enhance their efficiency. Consistent with the proposed approach, based on the geometric average, results, and predictions derived from the period in MPI can assist using a practical instrument for the general practitioner to compare the efficiency of uncertain cases and instruct effectively accordingly.

Fig. 2 .

 2 Figure2shows the structure of the proposed model:

Fig

  Figure3shows the flowchart of proposed model

Fig. 5 ,

 5 and Fig.6 and Table 3. In Fig 4. Fig. 5, and Fig.6 the horizontal axis characterizes the number of cement companies, and the vertical axis describes the average efficiency scores. As shown in Fig. 4, Fig. 5, Fig.6 and

3.

  Waste water removed in inefficient DMUs= If "(-inf-134021.2]" 20 ==> Then, Cement production ="(-inf-1964348.9]" 20 Based on the third strong rule (sup:0.18, Conf: 1, lift: 1:43) in step C, if energy consumption in inefficient DMUs within the quantified range happens 20 times (20 out of 110), cement production within the quantified range will occur 20 times (20 out of 110).

  main idea of the Apriori algorithm is to make numerous authorizations over the database. It employs an iterative method recognized as a breadth-first search (level-wise search) through the search space, where k-item sets are used to explore (k+1)-item sets.

	e r n a t i
	o
	n a
	l
	J
	o u
	r
	n
	a l
	o
	f
	E
	n
	e
	r g y
	S
	e c
	t o
	r
	M
	a n
	3	a g e m
		e n
		t

  addressed the regional distribution of carbon emission reduction goals in China based on the constituent part swarm optimization algorithm, fuzzy c-means clustering algorithm, and Shapley decomposition approach. Consequently, they clustered all regions into four classes

	e r n a t i
	o
	n a
	l
	J
	o u
	r
	n
	a l
	o
	f
	E
	n
	e
	r g y
	S
	e c
	t o
	r
	M
	a n
	5	a g e m
		e n
		t

Table 1 .

 1 The inputs, intermediate elements, and outputs for the 1ST DMUs

	Period	Input (Energy Consumption (1000TCE))	Intermediate Elements (Cement Production (1 ton))	Intermediate Elements (Pollution Control Investment 100000IRR))	Output (Waste Material Removed (1 ton))	Output (Waste gas Removed (1 ton))	Output (Solid waste Removed (1 ton))
	2015	1063116	730865	8000	107237	28.92	22421.24
	2016	1425490	1216569	143427	289883	19.24	543.25
	2017	1965650	1560395	153450	254899	25.93	2341.23
	2018	2539977	1861560	134566	264522	28.93	1451.22
	2019	2805942	2008740	144000	381952	6.54	8739.98

Energy consumptions in the companies are the only input in the first stage. Cement production (outputs of the first stage) and pollution control investments (input of the second stage) are two intermediate elements. Waste material removed in companies is the single output of the second stage.

  The BCC model changed the Constant Return to Scale (CRS) impression to Variable Return to Scale (VRS). The DMU controls under VRS, and it is observed that growth in inputs does not result in a relative change in the outputs. The BCC model divides Technical Efficiency (TE) based on the CCR model into two parts:

	e r n a t i
	o
	n a
	l
	J
	o u
	r
	n
	a l
	o
	f
	E
	n
	e
		r g y
		S
		e c
		t o
		r
		M
		a n
	Model:	7	a g e m
			e n
			t

Table 2 : Table2. Explanation of dimensionless parameters in nomenclature for primal and dual proposed model

 2

	Dimensionless parameters	Units
	𝐷𝑀𝑈 𝑗	Decision making units
	𝑆 -𝑖 (𝑆 -1 ,…,𝑆 -𝑚 )	Vector of input (Energy consumption) surplus (consumption slack)
	𝑆 + 𝑟 (𝑆 + 1 ,…,𝑆 + 𝑠 )	Vector of output (Wastewater removed) shortfall (production slack)
	𝑆 + ℎ (𝑆 + 1 ,…,𝑆 + 𝑑 )	Vector of output (Cement production) shortfall (production slack)
	𝑆 + 𝑡 (𝑆 + 1 ,…,𝑆 + 𝑣 )	Vector of output (Waste gas removed) shortfall (production slack)
	𝑆 + 𝑧 (𝑆 + 1 ,…,𝑆 + 𝑞 )	Vector of output (Solid waste removed) shortfall (production slack)
	𝑆 -𝑐 (𝑆 -1 ,…,𝑆 -𝑘 )	Vector of input (Pollution control investment) surplus (consumption slack)
	𝜆 𝑗	Non-negative scalar (dual variables that categorize the benchmarks for inefficient parts)
	𝑋 𝑖𝑗	input (Energy consumption) for 𝑚 𝑡ℎ	DMU 𝑛 𝑡ℎ
	𝑌 𝑟𝑗	output (Wastewater removed) for 𝑠 𝑡ℎ	DMU 𝑛 𝑡ℎ
	𝑀 ℎ𝑗	output (Cement production) for 𝑑 𝑡ℎ	DMU 𝑛 𝑡ℎ
	𝐸 𝑡𝑗	output (Waste gas removed) for 𝑣 𝑡ℎ	DMU 𝑛 𝑡ℎ
	𝐹 𝑧𝑗	output (Solid waste removed) for 𝑞 𝑡ℎ	DMU 𝑛 𝑡ℎ
	𝑁 𝑐𝑗	input (Pollution control investment) for 𝑘 𝑡ℎ	DMU 𝑛 𝑡ℎ
	𝑤	Free of sign for variable return to scale

Table 3 .

 3

	e r n a t i
	o
	n a
	l
	J
	o u
	r
	n
	a l
	o
	f
	E
	n
	e
	r g y
	S
	e c
	t o
	r
	M
	a n
	15	a g e m
		e n
		t

Table 3 .

 3 Productivity measurement results based on MPI-AM, MPI-BCC and MPI-CCR for 22 companies

		1		2		3
	Companies	MPI-AM	Rank	MPI-BCC	Rank	MPI-CCR	Rank
	1	1.098	5	1.094	5	1.090	5
	2	1.059	10	1.055	10	1.051	10
	3	1.133	4	1.129	4	1.125	4
	4	1.067	7	1.063	7	1.059	7
	5	0.898	21	0.894	21	0.890	21
	6	0.952	19	0.948	19	0.944	19
	7	1.061	8	1.057	8	1.053	8
	8	1.219	2	1.215	2	1.211	2
	9	1.051	11	1.046	11	1.043	11
	10	1.018	14	1.014	14	1.010	14
	11	0.874	22	0.870	22	0.866	22
	12	1.030	13	1.026	13	1.022	13
	13	1.058	9	1.054	9	1.050	9
	14	0.912	20	0.908	20	0.904	20
	15	1.015	16	1.011	16	1.007	16
	16	1.016	15	1.012	15	1.008	15
	17	1.431	1	1.427	1	1.423	1
	18	0.970	18	0.966	18	0.962	18
	19	1.033	12	1.029	12	1.025	12
	20	1.139	3	1.135	3	1.131	3
	21	1.071	6	1.067	6	1.063	6
	22	0.971	17	0.967	17	0.963	17
	The efficiency process for 22 companies are obtainable in Fig 4.

Table 3

 3

	, the 17th company has the highest efficiency score, and the 11th

Fig. 6. Average efficiency over 5-year periods for 22 DMUs (Result in MPI-CCR)

Table 6 .

 6 Accuracy of classification algorithms (all numbers are in percent) 17 miner tool showed the lowest accuracy whereas WEKA tool estimated the best accuracy when compared to two other tools. For K-NN, Tanagra represented the highest accuracy. Finally, it can be concluded from table 6, each tool has a highest accuracy for three suggested strong classification algorithms.

	Algorithms	WEKA	Tanagra	Rapid miner
	K-NN	89.5122	93.0918	97.1121
	Decision tree	91.2341	86.2467	80.1222
	Naïve Bayes	79.2311	96.1423	89.9087

From above table WEKA tool predicts a lowest accuracy for Naive Bayes, however for the same algorithm Tanagra tool predicts the best accuracy when compared to WEKA and Rapid miner. In case of Decision tree evaluation, Rapid

International Journal of Energy Sector Management