```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```

# Introducción.


En el post anterior, ya se ha presentado los elementos básicos de un análisis de regresión LOGIT o PROBIT. También se expuso un pequeño ejemplo y las conclusiones que se pueden sacar del mismo. En el presente post voy a seguir utilizando este modelo, pero utilizando ya un ejemplo algo más complicado.

Antes de entrar en materia, voy a explicar dos conceptos estadísticos, que van a servir para medir el nivel de dependencia de unas variables sobre la variable dependiente que se quiere predecir con el modelo. Estos dos conceptos son los denominados ** WOE (Weight of Evidence)** y **VI ( information value)**.

**El valor de WOE** nos va a facilitar información del poder predictivo de una variable independiente en relación a otra variable dependiente.

**El valor de Information value (IV)* es una técnica fácil y muy usada para seleccionar variables importantes en el modelo predictivo. A continuación se muestra una tabla que contiene los criterios más utilizados a la hora de aplicar esta técnica.

| Information value (IV) | Capacidad Predicción       |
|-------------------|----------------------------|
| <=0.02            | No sirve para predicción   |
| entre 0.02 y 0.1  | Bajo nivel de predicción   |
| entre 0.1 y 0.3   | Nivel predictivo medio     |
| entre 0.3 y 0.5   | Alto nivel predicción      |
| >= 0.5            | Excelente nivel predicción |


En concreto la fórmula empleada para calcular el valor de WOE es la siguiente:

\\[ WOE= \frac{\% No-events \, o \, fracasos}{\% Events \, o\,  Exitos} \\]

Es decir y a efectos prácticos, tendremos dos tipos de observaciones, una correspondiente a una variable que toma dos valores (1 y 0) y la otra correspondiente a una variable categórica. Entones a cada valor de la variable categórica le asignamos "Éxito", si el valor de la variable dicotómica toma el valor 1 o "Fracaso" si el valor de la dicotómica es cero. Con el código que sigue se explica mejor la forma de calcular este indicador estadístico.

El valor de Information Value (IV) se obtiene mediante la aplicación de la siguiente fórmula:

\\[ IV=\sum (\% de \, no-events\,o\,fracasos\,-\,\%events \,o\, éxitos)*WOE \\]


```{r}
woe_iv=function(Target,VariableCategorica){
  Tabla_WOE=table(VariableCategorica,Target)
  DF_WOE=data.frame(FRACASOS=Tabla_WOE[,1],EXITOS=Tabla_WOE[,2])
  DF_WOE$EXITOS_PORC=DF_WOE$EXITOS/sum(DF_WOE$EXITOS)
  DF_WOE$FRACASOS_PORC=DF_WOE$FRACASOS/sum(DF_WOE$FRACASOS)
  DF_WOE$WOE=log(DF_WOE$EXITOS_PORC/DF_WOE$FRACASOS_PORC)
  DF_WOE$IV=(DF_WOE$EXITOS_PORC-DF_WOE$FRACASOS_PORC)*DF_WOE$WOE
  DF_WOE
}
a<-c(1,1,1,0,0,1,0,1,1,1,0,0,1)
b<-c("a","b","c","a","a","b","b","c","a","b","c","c","c")
woe_iv(a,b)



```

Otra forma de calcula los valores anteriores es utilizando el paquete *InformationValue*.

```{r}
if (!require("InformationValue")) install.packages("InformationValue")
require("InformationValue")
print("Valores de WOE")
WOE(X=as.factor(b),Y=a)
print("Tabla de valores de WOE y de IV")
WOETable(X=as.factor(b),a)

```

En base a todos los conceptos expuestos anteriormente, vamos a elaborar un ejemplo con 15 variables. Unas son de tipo continuo y otras categóricas (tipo factor en el argot de R), y sobre ellas vamos aplicar los conceptos anteriormente expuestos ( en resumen cálculo de IV), para elegir las variables a priori más influyentes con las que se hará el análisis de regresión logit. Comenzamos por una presentación de los resultados. 

```{r }
#Leemos el fichero
datos2 <- read.csv("./data/adult.csv")
head(datos2)

```

La variable denominada *ABOVE50k* toma solo valores 1 ó 0 y es la que se va a utilizar como variables dependiente para hacer el análisis.


```{r}
colnames(datos2)
```
```{r}
str(datos2)
```

Observamos que básicamente hay dos tipos de variables, unas de tipo int y otras de tipo Factor.

En un problema de clasificación lo ideal es que la variable dependiente ( en este caso ABOVE50 ) esté equilibrada, es decir que aproximadamente tenga el mismo número de casos de uno y otro valor. Veamos en esta ocasión cómo se comporta este factor.

```{r}
table(datos2$ABOVE50K)
```

Como puede verse existe un sesgo claro hacia el valor cero donde se clasifican la mayor parte de los casos. Esta consideración hay que tenerla en cuenta a la hora de proceder a partir la muestra en elementos de los conjuntos train y test. Veamos a continuación cómo lo hacemos en este caso.

```{r}

# Creamos datos de training y test 
input_unos <- datos2[which(datos2$ABOVE50K == 1), ]  # Donde están los 1
input_ceros <- datos2[which(datos2$ABOVE50K == 0), ]  # Donde están los 0
set.seed(100)  # Generamos una semilla
input_unos_training <- sample(1:nrow(input_unos), 0.7*nrow(input_unos))  # sacamos el 70% de los unos para trainig

input_ceros_training <- sample(1:nrow(input_ceros), 0.7*nrow(input_ceros))  # Sacamos una muestra del 70% de ceros

training_unos <- input_unos[input_unos_training, ]  
training_ceros <- input_ceros[input_ceros_training, ]
trainingData <- rbind(training_unos, training_ceros)  # junto los unos y los ceros del training 

# Creamos los datos del test
test_unos <- input_unos[-input_unos_training, ]
test_ceros <- input_ceros[-input_ceros_training, ]
testData <- rbind(test_unos, test_ceros)  # Junto los unos y los ceros del subconjunto test

```

Veamos que de esta manera las proporciones se mantienen

```{r}
print("Proporción en la muestra total")
table(datos2$ABOVE50K)/nrow(datos2)
print("------------------");print("Proporcion en los datos train")
table(trainingData$ABOVE50K)/nrow(trainingData)
print("------------------");print("Proporcion en los datos test")
table(testData$ABOVE50K)/nrow(testData)
```

A continuación vamos a distinguir y separar las variables continuas de las variables de tipo factor.

```{r}
#Primero los nombre de las variables de tipo factor
Variables_Factor<-c ("WORKCLASS", "EDUCATION", "MARITALSTATUS", "OCCUPATION", "RELATIONSHIP", "RACE", "SEX", "NATIVECOUNTRY")
#A continuación los nombres de las variables de tipo continuo
Variables_Continua<-c("AGE", "FNLWGT","EDUCATIONNUM", "HOURSPERWEEK", "CAPITALGAIN", "CAPITALLOSS")

#Creamos a continuación un dataframe vacio que contiene el nombre de las variables independientes, #y además la columna IV que contendrá el valor de este indicador

iv_df <- data.frame(VARS=c(Variables_Factor, Variables_Continua), IV=numeric(14))  
```

A  continuación se muestra el código para tener el valor final de IV para cada una de las variables, información que nos servirá para elegir las variables que entrarán en el análisis.

```{r, warning=FALSE}
# El paquete smbbinning (https://cran.r-project.org/web/packages/smbinning/smbinning.pdf)
#Sirve entre otras cosas para clacular el IV.
#Además la función smbinning convierte una variable continua en categórica usando el método
#denominado "recursive partitioning"

if (!require("smbinning")) install.packages("smbinning"); require("smbinning")

# Calculamos IV para las variables categóricas
for(factor_var in Variables_Factor){
  smb <- smbinning.factor(trainingData, y="ABOVE50K", x=factor_var)  # WOE table
  if(class(smb) != "character"){ # heck if some error occured
    iv_df[iv_df$VARS == factor_var, "IV"] <- smb$iv
  }
}

# Calculamos IV para las variables continuas
for(continuous_var in Variables_Continua){
  smb <- smbinning(trainingData, y="ABOVE50K", x=continuous_var)  # WOE table
  if(class(smb) != "character"){  # any error while calculating scores.
    iv_df[iv_df$VARS == continuous_var, "IV"] <- smb$iv
  }
}

#Ordenamos las variables por el valor de IV obtenido

(iv_df <- iv_df[order(-iv_df$IV), ] )

```

Elegimos las variables con mayor valor de IV para hacer la regresión logit

```{r}
logitMod <- glm(ABOVE50K ~ RELATIONSHIP + MARITALSTATUS+AGE + EDUCATIONNUM, data=trainingData, family=binomial(link="logit"))
logitMod
```

```{r}
summary(logitMod)
```

A continuación calculamos las predicciones mediante dos fórmulas.

```{r}
predicted1 <- plogis(predict(logitMod, testData))  # Prediccion de los datos test
head(predicted1)
# o
predicted2 <- predict(logitMod, testData, type="response")  # Prediccion de los datos test
head(predicted1)
```

Como podemos observar el valor devuelto tanto por uno como por otro método es una probabilidad. Para clasificar las respuestas como 0 ó 1 (posibles valores de la variable dependiente), inicialmente se podría pensar que si p<0.5 se asigna cero y en caso contrario un uno. Pero esta es una decisión opinática y fácil de usar pero sin ningún criterio claro. Para optimizar la elección de este valor de corte, se puede usar la función **optimalCuttof** del paquete *InformationValue* que previamente ya se ha cargado en este ejemplo.

```{r}
if (!require("InformationValue")) install.packages("InformationValue")
require("InformationValue")
(PuntoCorte<-optimalCutoff(testData$ABOVE50K,predicted2)[1])
```

Si ahora se quiere obtener la predicción real ( es decir valores cero ó uno ), se haría de la siguiente manera.

```{r}
predic<-ifelse(predicted2>PuntoCorte,1,0)
head(predic)
```

Si quisiéramos sacar la curva de ROC para el modelo obtenido se haría de la siguiente manera. El gráfico también muestra el área bajo la curva, de tal manera que cuanto más área, mejor clasificación se obtiene. En este caso concreto, quizá se podría mirar si introduciendo alguna variable se puede conseguir mejorar la clasificación obtenida.

```{r}
plotROC(testData$ABOVE50K,predicted1)
```

