06-Cross.Rmd

---
output:
  pdf_document: default
  html_document: default
---

```{r, echo=FALSE, eval=FALSE}
library(bookdown); library(rmarkdown); rmarkdown::render("06-Cross.Rmd", "pdf_book")
```

# Cross-Efficiency

## The Ratio Model

Put yourself as a competitor trying to argue that you are the best in converting inputs into outputs among a set of other units. You have data on what the competitors' inputs and outputs. You can

## The Linear Programs for DEA

On the other hand, what if it allowed for \index{Blending} blending of units. There are a few assumptions that we could make. Let's start by saying that we can compare any particular products by rescaling (up or down) any other product as well as any combination of units.

We'll start by creating a mathematical framework. Can you find a combination of units that produces at least as much output using less input? Let's define the proportion of input needed as $\theta$. A value of $\theta=1$ then means no input reduction can be found in order to produce that unit's level of output. The blend of other units is described by a vector $\lambda$. Another way to denote this is $\lambda_j$ is the specific amount of a unit *j* used in setting the target for for performance for unit *k*. Similarly, $x_j$ is the amount of input used by unit *j* and $y_j$ is the amount of output produced by unit *j*.

This can be easily expanded to the multiple input and multiple output case by defining $x_i,j$ to be the amount of the *i*'th input used by unit *j* and $y_{r,j}$ to be the amount of the *r*'th output produced by unit *j*. For simplicity, this example will focus on the one input and one output case rather than the *m* input and *s* output case but the R code explicitly allows for $m,s>1$. To make the code more readable and help avoid name space conflicts, I will use `NX` instead of *m* to refer to the number of inputs (x's) and `NY` to be the number of outputs (y's) instead of *s*. Also, *n* is used to denote the number of Decision Making Units (DMUs) and therefore I'll use `ND` to indicate that in the R code.

We have two important sets of variables now. The first is $u_r$ which is the weight on the *r*'th output. The second is $v_i$ which is the weight on the *i*'th input.

The multiplier model can be thought of as finding a weighting scheme for outputs over inputs that give you the best possible score while giving no one better than *1.0.*

$$
\begin{split}
 \begin{aligned}
    \text {max } & \frac{\sum_{r=1}^{N^Y} u_r y_{r,k}} {\sum_{i=1}^{N^X} v_i x_{i,k} } \\
    \text{s.t.:  } & \frac{\sum_{r=1}^{N^Y} u_r y_{r,j}} {\sum_{i=1}^{N^X} v_i x_{i,j} }
                          \leq 1 \forall \; j\\
                  & u_r, v_i\geq 0  \; \forall \; r,i
  \end{aligned}
\end{split}
(\#eq:Ch6LPCCRIOM-Ratio-1a)
$$

This isn't a linear program because we are dividing functions of variables by functions of variables. We need to make a few transformations. First, we clear the denominator of each of the constraints resulting in the following formulation.

$$
\begin{split}
 \begin{aligned}
    \text {max } & \frac{\sum_{r=1}^{N^Y} u_r y_{r,k}} {\sum_{i=1}^{N^X} v_i x_{i,k} } \\
    \text{s.t.:  } & \sum_{r=1}^{N^Y} u_r y_{r,j} - \sum_{i=1}^{N^X} v_i x_{i,j} 
                          \leq 0 \; \forall \; j\\
                   & u_r, v_i\geq 0  \; \forall \; r,i
  \end{aligned}
\end{split}
(\#eq:Ch6LPCCRIOM-Ratio-2a)
$$

Now we will convert the problem \index{Inputs} input and \index{Outputs} output constraints from inequalities into equalities by explicitly defining slack variables.

There are an infinite number of possible combinations of numerators and denominators that can give the same ratio. The next step is to select a normalizing value for the objective function. Let's set the denominator equal to one. In this case, we simply add a \index{Constraints} constraint, $\sum_{i=1}^{N^X} v_i x_{i,k}$, to the linear program.

$$
\begin{split}
 \begin{aligned}
    \text {max } & \sum_{r=1}^{N^Y} u_r y_{r,k} \\
    \text{s.t.:  } & \sum_{i=1}^{N^X} v_i x_{i,k} \\
    & \sum_{r=1}^{N^Y} u_r y_{r,j} - \sum_{i=1}^{N^X} v_i x_{i,j} 
                          \leq 0 \; \forall \; j\\
                  &  u_r, v_i\geq 0  \; \forall \; r,i
  \end{aligned}
\end{split}
(\#eq:Ch6LPCCRIOM-NoSlack)  
$$

## Creating the LP - The Algebraic Approach

We will implement this using the `ompr` package again.

We're going to use our data from earlier. For this example, we will use the dataset from Kenneth Baker's third edition of *Optimization Modeling with Spreadsheets*, pages 175-178, Example 5.3 titled "Hope Valley Health Care Association." In this case, a health care organization wants to benchmark six nursing homes against each other.

```{r, eval=FALSE}
library(kableExtra)           
library(dplyr)            
library(ROI)         
library(ROI.plugin.glpk)  
library(ompr)             
library(ompr.roi) 
```

```{r, include=FALSE, warning=FALSE, message=FALSE}
# This code chunk loads the above listed packages without
#     disruptive warnings and messages.
suppressPackageStartupMessages(library(kableExtra))           
suppressPackageStartupMessages(library(dplyr)) 
suppressPackageStartupMessages(library(ROI))
suppressPackageStartupMessages(library(ROI.plugin.glpk))
suppressPackageStartupMessages(library(ompr))
suppressPackageStartupMessages(library(ompr.roi))
```

```{r multiplier-model }
XBaker1 <- matrix(c(150, 400, 320, 520, 350, 320, 
                    0.2, 0.7, 1.2, 2.0, 1.2, 0.7),
                  ncol=2,dimnames=list(LETTERS[1:6],
                                       c("x1", "x2")))

YBaker1 <- matrix(c(14000, 14000, 42000, 28000, 19000, 14000, 
                    3500, 21000, 10500, 42000, 25000, 15000),
                  ncol=2,dimnames=list(LETTERS[1:6],
                                       c("y1", "y2")))

ND <- nrow(XBaker1); NX <- ncol(XBaker1); NY <- ncol(YBaker1); 
    # Define data size
xdata      <-XBaker1 [1:ND,]  # Call it xdata
dim(xdata) <-c(ND,NX)         # structure data correctly
ydata      <-YBaker1[1:ND,]
dim(ydata) <-c(ND,NY)
```

Remember the inputs are hard coded as "x1" and"x2" to represent the *staff hours per day* and the *supplies per day* respectively. The two outputs of *reimbursed patient-days* and *privately paid patient-days* are named "y1" and "y2".

```{r}
YBaker1 <- matrix(c(14000, 14000, 42000, 28000, 19000, 14000, 
                    3500, 21000, 10500, 42000, 25000, 15000),
                  ncol=2,dimnames=list(LETTERS[1:6],
                                       c("y1", "y2")))

ND <- nrow(XBaker1); NX <- ncol(XBaker1); NY <- ncol(YBaker1); 
# Define data size
```

Note that I'm naming the data sets based on their origin and then loading them into xdata and ydata for actual operation.

```{r Structure-Results}
# Need to remember to restructure the results matrices.

results.efficiency <- matrix(rep(-1.0, ND), nrow=ND, ncol=1)
results.lambda     <- matrix(rep(-1.0, ND^2), nrow=ND,ncol=ND)
results.vweight    <- matrix(rep(-1.0, ND*NX), nrow=ND,ncol=NX) 
results.uweight    <- matrix(rep(-1.0, ND*NY), nrow=ND,ncol=NY) 
results.xslack     <- matrix(rep(-1.0, ND*NX), nrow=ND,ncol=NX) 
results.yslack     <- matrix(rep(-1.0, ND*NY), nrow=ND,ncol=NY) 
```

We are now ready to do the analysis. In Baker, the analysis was done using the \index{Multiplier model} multiplier model. In chapter 2 we used the \index{Envelopment model} envelopment model to examine this case. Now we will use the multiplier model using the \index{DEAMultiplierModel} `DEAMultiplerModel` package by Aurobindh Kalathil Puthanpura.

```{r, include=FALSE, message=FALSE, warning=FALSE}
library (MultiplierDEA)
```

```{r, warning=FALSE, message=FALSE}
library (MultiplierDEA)
# Example from Kenneth R. Baker, Optimization Modeling with  
#    Spreadsheets, Third Edition, p. 176, John Wiley & Sons, Inc.
dmu <- c("A", "B", "C", "D", "E", "F")
x <- data.frame(c(150,400,320,520,350,320),
                c(0.2,0.7,1.2,2.0,1.2,0.7))
rownames(x) <- dmu
colnames(x)[1] <- c("StartHours")
colnames(x)[2] <- c("Supplies")
y <- data.frame(c(14,14,42,28,19,14),
                c(3.5,21,10.5,42,25,15))
rownames(y) <- dmu
colnames(y)[1] <- c("Reimbursed")
colnames(y)[2] <- c("Private")
# Calculate the efficiency score
result <- DeaMultiplierModel(x,y,"crs", "input")
# Examine the efficiency score for DMUs
```

```{r, echo=FALSE}
kbl (result$Efficiency, caption="CRS Efficiency Scores", 
     booktabs = T, digits = 3, align = 'c') |>
    kable_styling(latex_options = c("HOLD_position"))
```

The efficiency scores match the earlier results and that of Baker.

The package also includes a cross-efficiency calculation function so as to find the cross-efficiencies given the multiplier model. Cross-efficiency calculated in this way tends to suffer from frequent cases of \index{Multiple optima} multiple optima. Later in this chapter, we will explore how to handle these issues. For now, let's apply cross-efficiency to the same data set.

```{r}
# Example from Kenneth R. Baker: Optimization Modeling with Spreadsheets, 
#    Third Edition, p. 176, John Wiley & Sons, Inc.
dmu <- c("A", "B", "C", "D", "E", "F")
x <- data.frame(c(150,400,320,520,350,320),
                c(0.2,0.7,1.2,2.0,1.2,0.7))
rownames(x) <- dmu
colnames(x)[1] <- c("StartHours")
colnames(x)[2] <- c("Supplies")
y <- data.frame(c(14,14,42,28,19,14),
                c(3.5,21,10.5,42,25,15))
rownames(y) <- dmu
colnames(y)[1] <- c("Reimbursed")
colnames(y)[2] <- c("Private")
# Calculate the efficiency score
result <- CrossEfficiency(x,y,"crs", "input")
# Examine the cross-efficiency score for DMUs
kbl (result$ce_ave, 
     align = "c", booktabs = T, digits = 4,
     caption="Cross-Efficiency Scores")    |>
    kable_styling(latex_options = c("HOLD_position"))
```

As might be expected, the \index{Cross-efficiency} cross-efficiency values are all significantly lower than the original efficiency scores. In cross-efficiency, no nursing home is able to pick an extreme and unique weighting scheme to their advantage. It also tends to break the ties at 1.0 that frequently occur in DEA so it usually gives a unique ranking. These features come at a significant cost though, which will be discussed in more detail in the future.

## Implementing Cross-Efficiency

Let's revisit the Baker model again.

```{r Baker-cross, eval=TRUE}

multiplierIO <- function (x,y) 
  { 
  ND <- nrow(x); NX <- ncol(x); NY <- ncol(y); # Define data size
  
  results.efficiency <- matrix(rep(-1.0, ND), nrow=ND, ncol=1)
  results.lambda     <- matrix(rep(-1.0, ND^2), nrow=ND,ncol=ND)
  results.vweight    <- matrix(rep(-1.0, ND*NX), nrow=ND,ncol=NX) 
  results.uweight    <- matrix(rep(-1.0, ND*NY), nrow=ND,ncol=NY) 

    for (k in 1:ND) {

    result <- MIPModel()                                              |>
    add_variable(vweight[i], i = 1:NX, type = "continuous", lb = 0)   |>
      add_variable(uweight[r], r = 1:NY, type = "continuous", lb = 0) |>
    set_objective(sum_expr(uweight[r] * y[k,r], r = 1:NY), "max")     |>
    add_constraint(sum_expr(vweight[i] * x[k,i], i = 1:NX) == 1)      |>
    add_constraint((sum_expr(uweight[r] * y[j,r], r = 1:NY) -
                    sum_expr(vweight[i] * x[j,i], i = 1:NX)) 
                  <= 0, j = 1:ND)
    result

    result <- solve_model(result, with_ROI(solver = "glpk", 
                                           verbose = FALSE))
    results.efficiency[k] <- objective_value (result) 

    # Get the weights - Output weights
    tempvweight <- get_solution(result, vweight[i])
    results.vweight[k,] <- tempvweight[,3]

   # Get the weights- Input weights
    tempuweight <- get_solution(result, uweight[i])
    results.uweight[k,] <- tempuweight[,3]

  } # End of for k loop
  
  resultlist <- list(efficiency=results.efficiency, 
                     vweight=results.vweight, 
                     uweight=results.uweight)
 
  return(resultlist)  

} # End of function  

resfunc <- multiplierIO(x,y)
DMUnames <- list(c(LETTERS[1:ND]))
Vnames<- lapply(list(rep("v",NX)),paste0,1:NX)   
# Input weight names: v1, ...
Unames<- lapply(list(rep("u",NY)),paste0,1:NY)   
# Output weight names: u1, ...

dimnames(resfunc$efficiency) <- c(DMUnames,"CCR-IO")
dimnames(resfunc$vweight)    <- c(DMUnames,Vnames)
dimnames(resfunc$uweight)    <- c(DMUnames,Unames)
```

```{r, echo=FALSE}
kbl (cbind(resfunc$efficiency,resfunc$vweight,resfunc$uweight), 
        caption="Efficiency Scores and Weights", booktabs = T) |>
    kable_styling(latex_options = c("HOLD_position"))
```

You can think of the cross-efficiency by going back to revisit the original ratio model of the multiplier model. Let's extend the definition of the input and output weights to reflect for which unit the weights were determined. For example, $u_{r,k}$ is then the weight on output *r* when analyzed from the perspective of DMU *k* and a similar interpretation applies to input weight.

$$
\begin{split}
\begin{aligned}
    \text {max } & \theta_{k,k}=\frac{\sum_{r=1}^{N^Y} u_{r,k} y_{r,k}} {\sum_{i=1}^{N^X} v_i,k x_{i,k} } \\
    \text{s.t.:  } & \theta_{j,k}=\frac{\sum_{r=1}^{N^Y} u_{r,k} y_{r,j}} {\sum_{i=1}^{N^X} v_{i,k} x_{i,j} }
                          \leq 1 \; \forall \; j\\
            & u_{r,k}, v_{i,k}\geq 0  \; \forall \; r,i
  \end{aligned}
\end{split}
(\#eq:Ch6LPCCRIO-Cross1)
$$

Now the cross-efficiency for a unit j can be calculated as the following. I refer to the value $\theta_{k,j}$ as the cross-evaluation score for unit *j* from the perspective of the evaluation of unit *k*.

The cross-efficiency score for unit *j* is simply the average of all of the cross-evaluation scores given to *j* by the other units when they do their analysis.

$$
\begin{split}
 \begin{aligned}
    CE_j= \frac {\sum_{k=1}^{N^D}\theta_{k,j}} {N^D}
        = \frac {\sum_{k=1}^{N^D} {
                        \frac {\sum_{r=1}^{N^Y} u_{r,k} y_{r,j} } 
                       {\sum_{i=1}^{N^X} v_{i,k} x_{i,j} } }    }
                {N^D}
  \end{aligned}
\end{split}
(\#eq:Ch6Calc-Cross)
$$

```{r}
results.crosseval  <- matrix(rep(-1.0, ND*ND), nrow=ND,ncol=ND) 
results.crosseff   <- matrix(rep(-1.0, ND)) 

for (k in 1:ND) {
   for (j in 1:ND) {results.crosseval[k,j] <- 
     sum(resfunc$uweight[k,]*y[j,])/sum(resfunc$vweight[k,]*x[j,])
                   }
}
for (k in 1:ND) {
  results.crosseff[k]<-mean(results.crosseval[,k])
}
dimnames(results.crosseval) <- c(DMUnames,DMUnames)
dimnames(results.crosseff) <- c(DMUnames)
```

```{r, echo=FALSE}
colnames(results.crosseff)<-c("Cross-Efficiency")
kbl (rbind(results.crosseval,t(results.crosseff)), 
     booktabs = T, digits=5, escape=F,
     caption="Cross-Evaluation Matrix and Cross-Efficiency Scores from Direct Calculation") |>
    row_spec (7, bold=TRUE) |>
  add_header_above(c("Weights\nFrom"=1,"Peer Evaluation of:"=6)) |>
    kable_styling(latex_options = c("HOLD_position", "scale_down"))
```

Notice that the values on the diagonal, $\theta_{k,k}$ are the same as the original efficiency scores.

We are going to be doing the cross-efficiency calculations often so let's convert this to a function.

```{r}
calccrosseff <- function (x,y,vmatrix,umatrix) {

  ND <- nrow(x); NX <- ncol(y); NY <- ncol(y); # Define data size
  
  DMUnames <- rownames(x)
  
  if (is.null(DMUnames)) {DMUnames<-list(c(LETTERS[1:ND]))}  
  
  crosseval  <- matrix(rep(-1.0, ND*ND), nrow=ND, ncol=ND,
                       dimnames=list(DMUnames,DMUnames))
  crosseff   <- matrix(rep(-1.0, ND), nrow=ND, 
                       dimnames=list(DMUnames, "Cross-Eff") )

  for (k in 1:ND) {
     for (j in 1:ND) {crosseval[k,j]<-
       sum(umatrix[k,]*y[j,]) / sum(vmatrix[k,]*x[j,])
                   }
   }
  
   for (j in 1:ND) {crosseff[j]<-sum(crosseval[,j])/ND}

  resultlist <- list(crosseval=crosseval, crosseff=crosseff)
 
  return(resultlist)
}
```

```{r, echo=FALSE}
rescross <- calccrosseff (x, y, resfunc$vweight, resfunc$uweight)

kbl (rbind(rescross$crosseval,t(rescross$crosseff)), 
     booktabs = T, digits=5,
     caption="Cross-Evaluation Matrix and Cross-Efficiency Scores from Our Function") |>
    row_spec (7, bold=TRUE) |>
  add_header_above(c("Weights\nFrom"=1,"Peer Evaluation of:"=6)) |>
    kable_styling(latex_options = c("HOLD_position", "scale_down"))
```

## Dealing with Cross-Efficiency's Multiple Optima

\index{Multiple optima|(} As discussed in the chapter on the Multiplier model, DEA often has multiple optima, particularly for efficient DMUs. While this doesn't affect the efficiency scores in either the envelopment or the multiplier model, it can have a major impact on the cross-evaluation scores and therefore the \index{Cross-efficiency} cross-efficiency scores.

To examine this issue, we need to invoke a secondary objective function again in the same manner as we did for the envelopment model.

A variety of mechanisms exist for dealing with the issue of multiple optima.

```{r Baker-Cross-Eval-Sec, eval=TRUE}
# Implements secondary objective function to resolve multiple optima

multiplierIOSec <- function (x, y) 
  { 
  ND <- nrow(x); NX <- ncol(y); NY <- ncol(y); # Define data size
  
  results.efficiency <- matrix(rep(-1.0, ND), nrow=ND, ncol=1)
  results.lambda     <- matrix(rep(-1.0, ND^2), nrow=ND,ncol=ND)
  results.vweight    <- matrix(rep(-1.0, ND*NX), nrow=ND,ncol=NX) 
  results.uweight    <- matrix(rep(-1.0, ND*NY), nrow=ND,ncol=NY) 

    for (k in 1:ND) {
      result <- MIPModel()                                     |>
      add_variable(vweight[i], i = 1:NX, 
                   type = "continuous", lb = 0)        
      add_variable(uweight[r], r = 1:NY, 
                   type = "continuous", lb = 0)                |>
      set_objective(sum_expr(uweight[r] * y[k,r], 
                             r = 1:NY), "max")                 |>
      add_constraint(sum_expr(vweight[i] * x[k,i], 
                              i = 1:NX) == 1)                  |>
      add_constraint((sum_expr(uweight[r] * y[j,r], r = 1:NY)-
                    sum_expr(vweight[i] * x[j,i], i = 1:NX)) 
                  <= 0, j = 1:ND)
      result

      result2 <- MIPModel()                                    |>
      add_variable(vweight[i], i = 1:NX, 
                   type = "continuous", lb = 0)                |>
      add_variable(uweight[r], r = 1:NY, 
                   type = "continuous", lb = 0)                |>
      set_objective(sum_expr(uweight[r] * sum(y[,r]), 
                             r = 1:NY), "max")                 |>
          # Modified objective function for 
      add_constraint(sum_expr(vweight[i] * sum(x[,i]), 
                              i = 1:NX) == 1)                  |>
      add_constraint(result, sum_expr(uweight[r] * sum(y[,r]), 
                                      r = 1:NY)-
                    sum_expr(vweight[i] * sum(x[,i]), i = 1:NX)  
                  <= 0)
      add_constraint((sum_expr(uweight[r] * y[j,r], r = 1:NY)-
                    sum_expr(vweight[i] * x[j,i], i = 1:NX)) 
                  <= 0, j = 1:ND)

    result <- solve_model(result, with_ROI(solver = "glpk", 
                                           verbose = FALSE))
    results.efficiency[k] <- objective_value (result2) 

    # Get the weights - Output weights
    tempvweight <- get_solution(result2, vweight[i])
    results.vweight[k,] <- tempvweight[,3]

   # Get the weights- Input weights
    tempuweight <- get_solution(result2, uweight[i])
    results.uweight[k,] <- tempuweight[,3]

  } # End of for k loop
  
  resultlist <- list(efficiency=results.efficiency, 
                     vweight=results.vweight, 
                     uweight=results.uweight)
 
  return(resultlist)  

} # End of function  

```

Let's explore the impact of the secondary objective function. The \index{MultiplierDEA} `MultiplierDEA` package has a cross-efficiency function but the `Mal_Ben` is a specialized version of the cross-efficiency function that implements the secondary goal of either malevolently trying to minimize other people's scores or benevolently trying to maximize other people's scores based on the work of Doyle and Green.

This can be framed as a two-phase model in the same way as the envelopment model slack maximization model.

The first phase is still the standard DEA multiplier model as discussed earlier. In the second phase, we try to either minimize or maximize the collective population's score. Note that in this case, we can define the population input (or output) either including the unit being studied by summing over all of j or all of j excluding the unit studied, k.

$$
\begin{split}
\begin{aligned}
    \text {min/max } & \frac{\sum_{r=1}^{N^Y} u_{r,k} \sum_{j=1}^{N^D} y_{r,j}} {\sum_{i=1}^{N^X} v_{i,k} \sum_{j=1}^{N^D}x_{i,j} } \\
    \text{s.t.:  } & \theta_{j,k}=\frac{\sum_{r=1}^{N^Y} u_{r,k} y_{r,j}} {\sum_{i=1}^{N^X} v_{i,k} x_{i,j} }
                          \leq 1 \; \forall \; j\\
& \theta_{k,k}=\frac{\sum_{r=1}^{N^Y} u_{r,k} y_{r,k}} {\sum_{i=1}^{N^X} v_{i,k} x_{i,k} }
= \theta^*_{k,k}\\            
& u_{r,k}, v_{i,k}\geq 0  \; \forall \; r,i
  \end{aligned}
\end{split}
(\#eq:Ch6LPCCRIO-2Phase)
$$

This ratio model can be linearized in the same way as the original ratio DEA model.

$$
\begin{split}
\begin{aligned}
    \text {min/max } & \sum_{r=1}^{N^Y} u_{r,k} \sum_{j=1}^{N^D} y_{r,j}\\ 
       \text{s.t.:  } & \sum_{i=1}^{N^X} v_{i,k} \sum_{j=1}^{N^D}x_{i,j} = 1 \\
& \sum_{r=1}^{N^Y} u_{r,k} y_{r,j} - \sum_{i=1}^{N^X} v_{i,k} x_{i,j} 
                          \leq 0 \; \forall \; j\\
& \sum_{r=1}^{N^Y} u_{r,k} y_{r,k}-  \theta^*_{k,k} \cdot\sum_{i=1}^{N^X} v_{i,k} x_{i,k} \leq 0
\\            
& u_{r,k}, v_{i,k}\geq 0  \; \forall \; r,i
  \end{aligned}
\end{split}
(\#eq:Ch6LPCCRIO-2Phase)
$$

Let's compare these to the arbitrary cross-efficiency function which does neither. Note that this does not guarantee the existence of a unique optimal solution in terms of input and output weights, it does greatly decrease the potential.

```{r}
resfunc    <- CrossEfficiency (x, y, rts="crs", 
                               orientation="input")
resfuncmal <- Mal_Ben (x, y, rts="crs", 
                       orientation="input", phase="mal")
resfuncben <- Mal_Ben (x, y, rts="crs", 
                       orientation="input", phase="ben")
```

```{r, echo=FALSE}
rownames(resfunc$ce_ave)<-"Single-Phase"
rownames(resfuncmal$ce_ave)<-"Malevolent"
rownames(resfuncben$ce_ave)<-"Benevolent"

kbl (rbind(t(cbind(x,y)), 
           resfunc$ce_ave, resfuncmal$ce_ave, resfuncben$ce_ave),
     booktabs = T, digits=4, 
     caption = "Cross-Efficiency Data and Results for Baker" )  |>
     row_spec (5:7, bold=TRUE) |>
     pack_rows ("Inputs", 1, 2) |>
     pack_rows ("Outputs", 3, 4) |>
     pack_rows ("Cross-Efficiency Results", 5, 7) |>
    kable_styling(latex_options = c("HOLD_position", "scale_down"))
```

The cross-efficiency results are very different numerically but the relative ranking is similar. In this case, nursing home A stays the top scoring unit in all three cross-efficiency models but is joined by D in the Benevolent model. Also, G is the lowest ranked in all three models.

Next, let's review the implementation of multiple optima in Auro's multiplier model. To test it, we'll load data from Takamura and Tone's 2003 paper on sites for relocating the capital of Japan.

```{r}

XTT03 <- matrix(c(1,1,1,1,1,1,1),ncol=1,
                dimnames=list(LETTERS[1:7],"x"))

YTT03 <- matrix(c(5,7,8,4,9,10,4,
                  10,10,7,8,4,2,7,
                  3,3,5,3,4,10,7),
                ncol=3, dimnames=list(LETTERS[1:7], 
                                      c("C1","C2","C3")))

ceTT03    <- CrossEfficiency (XTT03, YTT03, rts="crs", 
                               orientation="input")
ceTT03mal <- Mal_Ben (XTT03, YTT03, rts="crs", 
                       orientation="input", phase="mal")
ceTT03ben <- Mal_Ben (XTT03, YTT03, rts="crs", 
                       orientation="input", phase="ben")
```

```{r, echo=FALSE}
rownames(ceTT03$ce_ave)<-"Single-Phase"
rownames(ceTT03mal$ce_ave)<-"Malevolent"
rownames(ceTT03ben$ce_ave)<-"Benevolent"

kbl (rbind(as.data.frame(t(cbind(XTT03, YTT03)),digits=1),
           as.data.frame(rbind(ceTT03$ce_ave, ceTT03mal$ce_ave, 
                               ceTT03ben$ce_ave), digits=5)), 
           booktabs = T, digits=4,
        caption = "Data and Cross-Efficiency Results for Takamura and Tone, 2003") |>
     row_spec (5:7, bold=TRUE) |>
     pack_rows ("Inputs", 1, 1) |>
     pack_rows ("Outputs", 2, 4) |>
     pack_rows ("Cross-Efficiency Results", 5, 7) |>
    kable_styling(latex_options = c("HOLD_position", "scale_down"))
```

In this case, site location B, is highest ranked in all three cross-efficiency models.

Comments about cross-efficiency functions from `DEAmultiplier` package.

-   Mal_and_Ben is not defined as a function.
-   `CrossEfficiency` function should allow secondary objective function choice of "Malevolent" or "Benevolent"
-   Results from `CrossEfficiency` do not match Takamura and Tone.
-   Naming convention from results of analysis. Would be good if column names follow convention like my "helper file"-For example lambdas could use L\_ as a prefix.
-   data validation. Ex. allow "Ben", "ben" "benevolent", "Benevolent". Also, RTS could do "CRS" or "crs"

Let's try the Spring case

```{r}
XSpringL <- matrix(c(1,1,1,1,1,1,1,1), ncol=1,
                   dimnames=list(LETTERS[1:8],"x"))

YSpringL <- matrix(c(77,10,88,13,6,35,75,3,
                     75,4,34,33,10,55,61,2,
                     80,60,2,44,27,33,50,90,
                     39,44,81,69,88,77,33,8),
                   ncol=4, dimnames=list(LETTERS[1:8],
                                         c("HW","MQ","Proj", "Exam")))
SpringLres <- CrossEfficiency(XSpringL, YSpringL)
```

```{r, echo=FALSE}
res <- rbind (SpringLres$ceva_matrix, SpringLres$ceva_max, 
              SpringLres$ce_ave)
resrows <-cbind (c("A", "B", "C", "D", "E", "F", "G", "H", 
         "Efficiency", "Cross-Efficiency"))
kbl (cbind(resrows, as.data.frame(res,digits=4)), 
     row.names=F, booktabs = T, digits=4, escape=F,
        caption = "Results from Class Grading Example")    |>
     row_spec (9:10, bold=TRUE) |>
     pack_rows ("Cross-Evaluation Matrix", 1, 8)           |>
     pack_rows ("Results", 9, 10) |>
         kable_styling(latex_options = c("HOLD_position", "scale_down"))
```

\index{Multiple optima|)}

## Future Issues for Cross-Efficiency Chapter

Issues to consider in the future for this chapter or other chapters:

-   Fixed Weighting nature of certain model dimensions (no-wait to talk about it in baseball chapter?)
-   Check to see if Auro is using population with or without k in mal_ben code.
-   Add reference to Doyle and Green.