Introducción

El aprendizaje automatizado —machine learning (ML)— es una rama de la inteligencia artificial cuyo objetivo es conseguir que una “máquina” aprenda a partir de la experiencia. Básicamente, los algoritmos toman un conjunto de datos, los analizan para buscar en ellos ciertas pautas o patrones y una vez identificadas, las emplean para realizar predicciones. En otras palabras se trata de predecir el comportamiento futuro a partir de comportamientos pasados observados.

Dependiendo de cómo se aborde el problema del ML, los diferentes algoritmos se pueden agrupar en:

  • Aprendizaje supervisado

    Cuando las entradas al sistema (la base de conocimiento) están formadas por un conjunto de datos etiquetados a priori. Es decir, sabemos la clasificación correcta de un conjunto de datos, y a partir de ellos se va a generar una función predictora.

    • Algunos algoritmos supervisados: KNN, Random Forest, SVM, ANN, Naïve Bayes, …
  • Aprendizaje no supervisado

    Las entradas al sistema están formadas por un conjunto de datos de los que se desconoce su clasificación correcta. En este caso, se espera que el algoritmo sea capaz de reconocer patrones para poder etiquetar y clasificar nuevos datos de entrada.

    • Algunos algoritmos no supervisados: K-means, clustering jerárquico, ANN no supervisado, …

Pasos para aplicar ML

  • Recoger y almacenar datos

  • Explorar y preparar los datos

  • Entrenar al clasificador con nuestros datos

  • Evaluar el desempeño del modelo

  • Mejorarlo si fuese necesario

Una vez completados estos pasos, nuestro modelo puede ser implementado en la clasificación de nuevos datos.

KNN

K-Nearest Neighbours (Los \(k\) vecinos más próximos) es un algoritmo de clasificación supervisada,basado en distancias. Aunque no necesariamente la distancia euclídea generalmente se recurre a ella.

Siendo \(p\) y \(q\) dos observaciones, \(p_1\), \(p_2\), … , \(p_n\), y \(q_1\), \(q_2\), … , \(q_n\) los valores que toman esas observaciones en las \(n\) variables, la distancia euclídea entre estas dos observaciones se define como:

\[ dist(p,q) = \sqrt{(p_1-q_1)^2+(p_2-q_2)^2+...+(p_n-q_n)^2}\]

Su funcionamiento es muy simple: dado un conjunto de observaciones etiquetadas —clasificadas—, el algoritmo va a asignar una nueva observación a la clase —al grupo— más cercana.

Dado que esa distancia puede definir se de diferentes maneras y, en función de cómo se haga, los resultados podrían variar.

En el ejemplo de las siguientes imágenes se trata de clasificar un tomate a alguno de los grupos previamente definidos.

El principal problema que puede plantear este método es decidir sobre cuantos vecinos se calcula la distancia. Si nos fijamos en la imagen siguiente, al elegir \(k = 2\) —dos vecinos— nuestra observación —circulo verde— se clasificaría como triángulo rojo. Sin embargo si fijamos \(k = 3\) —tres vecinos— entonces nuestra observación será clasificada como cuadrado azul.

No hay una única regla, una formula general dice que el número de vecinos a tener en cuenta, está en función del tamaño de la muestra \(n\):

\[K = n^{\frac{1}{2}}\]

Otros métodos incluirían ensayos con diferentes números de vecinos y posteriormente evaluar el desempeño de la clasificación —curvas ROC, Cross–validation, …— para quedarnos con aquel número con el que se obtenga un mejor desempeño.

Implementación de KNN en R

Utilizaremos como ejemplo los datos de iris.

data( iris )

Preprocesar datos

El preprocesado de datos incluye normalizar si es necesario, o reducir el número de variables si fuese necesario mediante otras técnicas multivariantes.

Por ejemplo, estandarizamos, \(\dfrac{x-\overline{x}}{\sigma^2}\), mediante la función scale()

NormIris <- as.data.frame( scale( iris[, 1:4 ],
                                  scale = TRUE, 
                                  center = TRUE ) )
summary( NormIris )
  Sepal.Length       Sepal.Width       Petal.Length    
 Min.   :-1.86378   Min.   :-2.4258   Min.   :-1.5623  
 1st Qu.:-0.89767   1st Qu.:-0.5904   1st Qu.:-1.2225  
 Median :-0.05233   Median :-0.1315   Median : 0.3354  
 Mean   : 0.00000   Mean   : 0.0000   Mean   : 0.0000  
 3rd Qu.: 0.67225   3rd Qu.: 0.5567   3rd Qu.: 0.7602  
 Max.   : 2.48370   Max.   : 3.0805   Max.   : 1.7799  
  Petal.Width     
 Min.   :-1.4422  
 1st Qu.:-1.1799  
 Median : 0.1321  
 Mean   : 0.0000  
 3rd Qu.: 0.7880  
 Max.   : 1.7064  

Previsualizar datos

plot( iris[ ,1:4 ], col = iris$Species )

Generar las muestras

En una clasificación supervisada, normalmente se selecciona del total de datos, un subconjunto de aproximadamente de entre 2/3 y 3/4 del total para entrenar al clasificador y el resto de datos se reserva para testar el desempeño del clasificador, o dicho de otra forma, para comprobar cómo de bien o mal clasifica nuestra máquina·

#generar subconjunto de 100 datos para entrenar
seleccion <- sample( c( 1:150 ), 100 )
iris.train <- iris[ seleccion, 1:5 ]

#generar subconjunto de 50 datos para testar
iris.test <- iris[ -seleccion, 1:5 ]

Entrenar/testar clasificador

La función knn del paquete class entrena y clasifica en un único paso. Como parámetros le pasamos en conjunto de datos de entrenamiento, el conjunto de datos para el test, las clases a las que pertenece cada elemento, y el número de vecinos.

# Genera modelo y hace predicción a la vez
iris.pred <- knn( train = iris.train[ , 1:4 ],
                  test  = iris.test[  , 1:4 ], 
                  cl    = iris.train[ , 5 ],
                  k     = 2 )

# predicción
print( iris.pred )
 [1] setosa     setosa     setosa     setosa     setosa    
 [6] setosa     setosa     setosa     setosa     setosa    
[11] setosa     versicolor versicolor versicolor versicolor
[16] versicolor versicolor versicolor versicolor versicolor
[21] versicolor virginica  versicolor virginica  versicolor
[26] versicolor virginica  versicolor versicolor versicolor
[31] versicolor versicolor versicolor virginica  virginica 
[36] virginica  virginica  virginica  virginica  virginica 
[41] virginica  virginica  virginica  virginica  virginica 
[46] virginica  versicolor virginica  virginica  virginica 
Levels: setosa versicolor virginica
# realidad 
print( iris.test[ ,5 ] )
 [1] setosa     setosa     setosa     setosa     setosa    
 [6] setosa     setosa     setosa     setosa     setosa    
[11] setosa     versicolor versicolor versicolor versicolor
[16] versicolor versicolor versicolor versicolor versicolor
[21] versicolor versicolor versicolor versicolor versicolor
[26] versicolor versicolor versicolor versicolor versicolor
[31] versicolor versicolor versicolor virginica  virginica 
[36] virginica  virginica  virginica  virginica  virginica 
[41] virginica  virginica  virginica  virginica  virginica 
[46] virginica  virginica  virginica  virginica  virginica 
Levels: setosa versicolor virginica

Validar resultados

Podemos ver la diferencia entre lo clasificado y la realidad con table

table( Predic = iris.pred, 
       Realidad = iris.test[ , 5 ] )
            Realidad
Predic       setosa versicolor virginica
  setosa         11          0         0
  versicolor      0         19         1
  virginica       0          3        16

Una vez creado un clasificador y comprobado que clasifica bien, podemos probarlo clasificando un dato nuevo. Si encuentro un iris en el campo y mido pétalo y sépalo, puedo utilizar mi clasificador para determinar a que variedad pertenece.

miIris <- c( 5.0, 3.5, 1.3, 0.1 )

mi_predict <- knn( train = iris.train[ , 1:4 ], 
                   test = miIris, 
                   cl = iris.train[ , 5 ], 
                   k = 2, 
                   prob = TRUE ) 
mi_predict
[1] setosa
attr(,"prob")
[1] 1
Levels: setosa versicolor virginica

Máquina de soporte de vectores (SVM)

Algoritmo de clasificación supervisada, en el que dado un conjunto de datos etiquetados, el algoritmo construye un límite o frontera óptimo, llamado hiperplano, que separa los datos según su categoría para posteriormente asignar nuevos datos al grupo con mayor probabilidad de pertenencia.

El fundamento básico es el mismo que en el caso de KNN:

  • Conjunto de datos etiquetados
  • Construcción de un modelo
  • Dado un nuevo dato de etiqueta desconocida
  • Modelo es capaz de predecir a que clase pertenece

A diferencia de KNN, el algoritmo busca los límites o fronteras entra las clases. Esa frontera en un espacio de dos dimensiones —dos variables— sería una línea, en 3 dimensiones sería un plano, y genéricamente para un espacio multidimensional —multivariante— sería un hiperplano. Este hiperplano separa el hiperespacio en diferentes regiones, donde las observaciones que caigan dentro de la misma región pertenecerán a la misma clase.

La dificultad radica en definir cuál es el hiperplano óptimo. El algoritmo busca una frontera que maximize la distancia entre ella y los elementos marginales de cada clase.

En ocasiones definir esa frontera es intuitivo, pero es normal que cuando se trabaja con más de dos variables la definición de ese hiperplano se vuelva muy compleja y computacionalmente muy costosa.



Para simplificar el cálculo del hiperplano se pueden realizar una serie de transformaciones sobre los datos de manera que sea menos compleja la obtención de las diferentes regiones y por tanto hacer menos costoso —computacionalmente hablando— el aprendizaje.

Para poder llevar a cabo esto se utilizan distintas funciones kernel que pueden llevar a cabo diferentes tipos de transformaciones. Estas son algunas de ellas:

  • Polinomial: \(K(x_i, x_j) = (x_ix_j)^n\)

  • Perceptrón: \(K(x_i, x_j) = \left \|x_ix_j \right \|\)

  • Base radial Gaussiana: \(K(xi, xj)=e^{-\frac{(xi-xj)^2}{2\sigma^2}}\)

  • Sigmoidal: \(K(xi, xj)=tanh(xi\cdot xj-\Theta)\)

Gráficamente la función kernel aplica una transformación sobre los datos, simplificando enormemente la definición de la frontera.

Implementación de SVM en R

Seguiremos utilizando como ejemplo los datos de iris-

Preprocesar datos

Incluye normalizar si es necesario, reducir el número de variables si fuese necesario…

Para generar las muestras utilizamos el procedimiento anterior

#generar subconjunto de 100 datos para entrenar
seleccion <- sample( c( 1:150 ), 100 )
iris.train <- iris[ seleccion, 1:5 ]

#generar subconjunto de 50 datos para testar
iris.test <- iris[ -seleccion, 1:5 ]

Aplicar SVM

Con la función svm del paquete e1071 se puede generar un modelo mediante algoritmo de soporte de vectores. Como parámetros hay que proporcionar los datos en forma de fórmula de R, donde la variable respuesta son las clases frente al resto de variables predictoras Species ~ ., hay que definir sobre que datos se está refiriendo la fórmula, si queremos normalizar los datos y el tipo de función kernel que queremos utilizar.

# Esta función sólo genera el modelo
svm_model <- svm( Species ~ . , 
                  data   = iris.train, 
                  scale  = TRUE,
                  kernel = "linear", 
                  probability = TRUE)

summary( svm_model )

Call:
svm(formula = Species ~ ., data = iris.train, kernel = "linear", 
    probability = TRUE, scale = TRUE)


Parameters:
   SVM-Type:  C-classification 
 SVM-Kernel:  linear 
       cost:  1 

Number of Support Vectors:  22

 ( 9 1 12 )


Number of Classes:  3 

Levels: 
 setosa versicolor virginica

Una vez entrenado el modelo, podemos hacer predicción con el conjunto de datos de iris.test.

SVM_predict <- predict( svm_model, 
                        iris.test[,1:4], 
                        probability = TRUE, 
                        decision.values = TRUE )
summary(SVM_predict)
    setosa versicolor  virginica 
        18         15         17 

Contrastar la predicción con la realidad

# Versión simple
table( SVM_predict, 
       Realidad = as.factor( iris.test[ ,5 ] ) )
            Realidad
SVM_predict  setosa versicolor virginica
  setosa         18          0         0
  versicolor      0         15         0
  virginica       0          2        15

Igual que en el caso de KNN, podemos hacer una predicción sobre nuestro modelo con una muestra recogida en el campo, y obtendremos a qué clase ha sido asignada y con que probabilidad.

miIris <- matrix( c( 5.0, 3.5, 1.3, 0.1 ), 
                  nrow = 1, ncol = 4 )

colnames( miIris ) <- colnames( iris[ , 1:4 ] )

# predict requiere un matrix o data.frame no un vector
predict( svm_model, miIris, probability = TRUE )
     1 
setosa 
attr(,"probabilities")
    virginica    setosa versicolor
1 0.009314366 0.9763348 0.01435079
Levels: setosa versicolor virginica

Naïve Bayes

Es un clasificador probabilístico “ingenuo” basado en el teorema de Bayes, que asume que la probabilidad de cada variable es independiente de las demás. Básicamente, su lógica se basa en que dado un conjunto de datos de entrenamiento etiquetados, el clasificador calcula la probabilidad observada para cada clase, en función de los valores de sus variables. Cuando es usado posteriormente para predecir datos sin etiquetar, asigna estos datos a la clase con mayor probabilidad de pertenencia.

Refrescando la memoria bayesiana

Probabilidad de un evento

\[ P(X) = \frac{Casos\; favorables}{Casos\; posibles} \]

Probabilidad de que al lanzar un dado salga un 3:

\[ P(X=3) = \frac{1}{6} \]

Probabilidad condicionada

Probabilidad de que suceda un evento cuando ya ha sucedido otro:

\[ P(Y \mid X) \]

Probabilidad de que \(Y\) tome un valor determinado cuando \(X\) ya ha tomado uno dado.

Ejemplo: probabilidad de obtener un 10 entre dos dados cuando el primero ha salido 4

\[ P(Y=10 \mid X=4) \]

Teorema de Bayes

Basado en el Teorema de Bayes:

\[ P(Y \mid X) = \frac{ P(X \mid Y ) \, P(Y) }{ P(X) } \]

En castellano

\[p(a posteriori) = \frac{ P(probabilidad\ condicionada) \, P(a\ priori) }{ P(total) }\]

Veamos como funciona Naïve Bayes como predictor con un ejemplo:

Tenemos datos históricos sobre si se ha jugado o no al tenis dadas una serie de situaciones meteorológicas muy simples.
En base a esto, y empleando la lógica bayesiana determinaremos qué probabilidad tenemos de que haya partido un día nublado, frío, con alta humedad y sin viento.

La siguiente tabla, son nuestros datos históricos.

Cielo Temperatura Humedad Viento Jugar
Soleado Calor Alta Si No
Soleado Calor Alta No Si
Soleado Templado Baja Si Si
Nublado Frio Baja No Si
Nublado Calor Normal No Si
Nublado Frio Alta Si No
LLuvia Templado Alta No No

Probabilidades

¿Jugaremos al tenis (\(Y\)) un día nublado, frío, con humedad alta y sin viento (\(X\))?

\[ P(Y=Si \mid X_{v1,v2,v3,v4}) = \frac{P(X \mid Y_{viento}) \cdot P(Y)}{P(total)} \]

Probabilidades

Las probabilidades “observadas” se obtienen de las frecuencias de la tabla anterior.

\[\begin{array}{rlr} P(SI_{jugar}) &=& 4/7 & P(NO_{jugar}) &=& 3/7 \\ P(nublado \mid SI_{jugar}) &=& 2/4 & P(nublado \mid NO_{jugar}) &=& 1/3 \\ P(frio \mid SI_{jugar}) &=& 1/4 & P(frio \mid NO_{jugar}) &=& 1/3 \\ P(Alta_{humedad} \mid SI_{jugar}) &=& 1/4 & P(Alta_{humedad} \mid NO_{jugar}) &=& 3/3 \\ P(Si_{viento} \mid SI_{jugar}) &=& 1/4 & P(Si_{Viento} \mid NO_{jugar}) &=& 2/3 \\ \end{array}\]

La hora de la verdad ¿jugaremos?

\[ P(Y=Si \mid X_{v1,v2,v3,v4}) = \frac{2/4 \cdot 1/4 \cdot 1/4 \cdot 1/4 \cdot 4/7 =0.024}{0.024+0.21}=0.103 \] \[ P(Y=No \mid X_{v1,v2,v3,v4}) = \frac{1/3 \cdot 1/3 \cdot 3/3 \cdot 2/3 \cdot 3/7 = 0.21}{0.024+0.21}=0.897 \]

Implementación de Naïve Bayes en R

Recurriremos a la función naiveBayes de la librería del mismo nombre. Seleccionaremos, como en los casos anteriores un conjunto de datos de entrenamiento (iris.train) y de testeo (iris.test).

iris.NB <- naiveBayes( Species ~ ., 
                       data = iris.train )
iris.NB

Naive Bayes Classifier for Discrete Predictors

Call:
naiveBayes.default(x = X, y = Y, laplace = laplace)

A-priori probabilities:
Y
    setosa versicolor  virginica 
      0.32       0.33       0.35 

Conditional probabilities:
            Sepal.Length
Y                [,1]      [,2]
  setosa     5.034375 0.3755507
  versicolor 5.906061 0.5454009
  virginica  6.560000 0.6227076

            Sepal.Width
Y                [,1]      [,2]
  setosa     3.481250 0.3693390
  versicolor 2.818182 0.2909624
  virginica  2.960000 0.3050554

            Petal.Length
Y                [,1]      [,2]
  setosa     1.462500 0.1621429
  versicolor 4.290909 0.4888530
  virginica  5.531429 0.5454472

            Petal.Width
Y                [,1]      [,2]
  setosa     0.237500 0.1184578
  versicolor 1.339394 0.1999053
  virginica  2.008571 0.2559543
predNB <- predict( iris.NB, 
                   iris.test[ , 1:4] )
predNB
 [1] setosa     setosa     setosa     setosa     setosa    
 [6] setosa     setosa     setosa     setosa     setosa    
[11] setosa     setosa     setosa     setosa     setosa    
[16] setosa     setosa     setosa     versicolor versicolor
[21] versicolor versicolor versicolor versicolor versicolor
[26] versicolor versicolor versicolor versicolor versicolor
[31] versicolor versicolor versicolor versicolor versicolor
[36] virginica  virginica  virginica  virginica  virginica 
[41] virginica  virginica  virginica  virginica  versicolor
[46] virginica  virginica  virginica  virginica  virginica 
Levels: setosa versicolor virginica

Contrastar la predicción

table( Predicción = predNB, Realidad = iris.test[ , 5 ] )
            Realidad
Predicción  setosa versicolor virginica
  setosa         18          0         0
  versicolor      0         17         1
  virginica       0          0        14

Y con mi muestra encontrada en el campo … ¿que?

miIris <- matrix( c( 5.0, 3.5, 1.3, 0.1 ), 
                  nrow = 1,
                  ncol = 4 )

colnames( miIris ) <- colnames( iris[ , 1:4 ] )

# En este caso predict requiere un matrix o data.frame no un vector
pdr <- predict( iris.NB, miIris, probability = TRUE )
summary( pdr )
    setosa versicolor  virginica 
         1          0          0 

Lo que nos devuelve la predicción es la probabilidad de pertenencia a cada clase. En este caso es una probabilidad del 100\(\,\)% de pertenecer a la clase setosa, y por tanto 0\(\,\)% al resto de clases.

Redes Neuronales Artificiales (ANN)

Paradigma de aprendizaje automatizado inspirado en sistemas nerviosos biológicos, constituidos por una serie de nodos (neuronas) y una serie de conexiones entre los nodos (sinapsis).

La topología básica consta de:

  • Unos nodos de entrada al sistema, generalmente tantos nodos como variables.
  • Una serie de capas ocultas intermedias con un número variable de nodos.
  • Unos nodos de salida, uno por cada respuesta.

Cada nodo de entrada tiene asociado un peso \(W_i\) y en cada nodo: se aplica una función de activación (opcional) y una de propagación con la suma de las entradas ponderadas por sus pesos.

A pesar del nombre, su parecido con los sistemas biológicos se queda en la estructura física, es decir una serie de nodos interconectados formando una red con entradas y salidas.

Pero aunque parezca algo complejo, una especie de caja negra en la que le damos información y esperamos a que “eso” obre el milagro y nos de un resultado, el fundamento es bastante simple:

Dados una serie de parámetros de entrada, las redes neuronales buscan la forma de combinarlos para obtener el resultado esperado.

Aunque hay muchas topologías, la más común es la conocida como perceptrón

  • Componentes

    • \(N\) entradas, \(x_1\), …, \(x_n\)

    • Cada entrada con un peso \(x\), \(x_1\), …, \(x_n\)

    • Un nodo de entrada extra llamada bias (intercepto)

    • Suma de las entradas ponderada por sus pesos:

      \[y = \sum{x_0w_0,...,x_nw_n}\]

  • Función de activación p.e.:

    \(f_a(x) = 1\) si \(y>0\), \(f_a(x) = -1\) si \(x\leq0\)

Veamos su fundamento con un ejemplo:

Supongamos que tenemos la nota de dos exámenes parciales y la nota final. No sabemos como pondera el profesor las notas de los exámenes parciales para obtener la final. La red neuronal irá probando distintas combinaciones hasta que de con la correcta.

Nuestros parámetros de entrada son las dos notas y la salida será la nota final esperada.

Esto podría complicarse. Imaginemos ahora que tenemos las notas de dos exámenes más la nota de una práctica. Hemos suspendido un examen y sin embargo hemos aprobado la asignatura. En este caso ya no es tan trivial ponderar el peso de cada examen.

Pero podría ser peor. La retorcida mente de este profesor, podría tener en cuenta, la asistencia, la participación, la actitud con los compañeros, además de las notas de teoría y prácticas. Además a cada uno de estos parámetros le podría dar un peso diferente, incluso dar pesos diferentes a combinaciones de parámetros. Es más, podría ser que no todas las combinaciones fuesen lineales, algunas podrían ser cuadráticas, logarítmicas, …

Pero afortunadamente para nosotros, no tenemos que preocuparnos de todo esto, lo único que tendríamos que hacer básicamente, es decir cuales son nuestras notas de entrada y nuestra nota final. La red se encarga de generar el modelo y de esta forma la próxima vez que sepamos las notas iniciales, no tendremos que esperar a saber la nota final, nuestro modelo nos lo va a predecir.

Implementación de redes neuronales en R

Preprocesar datos

Incluye normalizar si es necesario, reducir el número de variables si fuese necesario, …

Generar las muestras

Este paso es igual lo que hemos visto en los ejemplos anteriores

data( iris )
seleccion  <- sample( c( 1:150 ), 100 )
iris.train <- iris[  seleccion, 1:5 ]
iris.test  <- iris[ -seleccion, 1:5 ]

Binarizar las categorias

En este caso hay que generar nuevas columnas de VERDADERO/FALSO para pertenencia de cada observación a cada clase —especie—.

iris.train <- cbind( iris.train, 
                     iris.train$Species == "setosa" )

iris.train <- cbind( iris.train, 
                     iris.train$Species == "versicolor" )

iris.train <- cbind( iris.train, 
                     iris.train$Species == "virginica" )

names( iris.train )[ 6 ] <- "setosa"
names( iris.train )[ 7 ] <- "versicolor"
names( iris.train )[ 8 ] <- "virginica"

head( iris.train )
    Sepal.Length Sepal.Width Petal.Length Petal.Width
140          6.9         3.1          5.4         2.1
131          7.4         2.8          6.1         1.9
82           5.5         2.4          3.7         1.0
112          6.4         2.7          5.3         1.9
49           5.3         3.7          1.5         0.2
70           5.6         2.5          3.9         1.1
       Species setosa versicolor virginica
140  virginica  FALSE      FALSE      TRUE
131  virginica  FALSE      FALSE      TRUE
82  versicolor  FALSE       TRUE     FALSE
112  virginica  FALSE      FALSE      TRUE
49      setosa   TRUE      FALSE     FALSE
70  versicolor  FALSE       TRUE     FALSE

Entrenar la red neuronal

La función neuralnet del paquete del mismo nombre, requiere la introducción de parámetros de forma similar a como se hizo en SVM, mediante una foŕmula, donde las variables dependientes —setosa, versicolor, virgínica— lo son como función de las predictoras. Podemos definir también la topología de la red, es decir, cuántas capas y cuántos nodos por capa.

Una vez entrenada nuestra red, podemos visualizar su topología, los pesos ponderados asignados a cada nodo. Con la función print podemos obtener bastante información acerca del modelo.

iris.nnt <- neuralnet( setosa + versicolor + virginica ~ Sepal.Length +
                         Sepal.Width +
                         Petal.Length +
                         Petal.Width, 
                       data = iris.train, hidden = c( 3, 2 ) )

# print(iris.nnt) # Ejecutar para ver los resultados
plot( iris.nnt, col.intercept = "blue" , rep = "best" )

Predicción

Con los datos que hemos reservamos para validar, vamos a hacer la predicción, en este caso en lugar de la función predict usada anteriormente, empleamos compute.

pred.nn <- compute( iris.nnt, iris.test[ 1:4 ] )
# print(pred.nn) # ejecutar para ver los resultados

El resultado de la predicción no es tan “amigable” como en los otros algoritmos vistos, pero con un poco de código lo podemos extraer y dejar más claro.

resultado <- 0
for ( i in 1:dim( pred.nn$net.result )[1] ){
  resultado[i] <- which.max( pred.nn$net.result[ i, ] )
}
resultado[ resultado == 1 ] <- "setosa"
resultado[ resultado == 2 ] <- "versicolor"
resultado[ resultado == 3 ] <- "virginica"

table( Predicción = resultado, Realidad = iris.test[ , 5 ] )
            Realidad
Predicción  setosa versicolor virginica
  setosa         18          0         0
  versicolor      0         14         1
  virginica       0          2        15

Evaluar la eficacia de un modelo

Vamos a estudiar algunos de los mecanismos de validación de la clasificación realizada por nuestro modelo predictivo, para ello generaremos un modelo nuevo con datos de pacientes que tienen un tumor, cuyo diagnóstico está etiquetado como benigno o maligno.

La primera columna contiene un identificador del sujeto, la segunda corresponde al diagnóstico y el resto son variables medidas al tumor.

wdbc <- read.csv( "recursos/datos/wdbc.csv" )

#eliminar columna de ID.
wdbc <- wdbc[ -1 ]

# Normalizar los datos
wdbc_N <- as.data.frame( scale( wdbc[ 2:31 ], 
                                scale  = TRUE, 
                                center = TRUE ) )

wdbc_N$diagnosis <- wdbc$diagnosis

Crear datos de entrenamiento y datos de test

# Datos
seleccion <- sample( nrow( wdbc_N ), 425 )
w_train <- wdbc_N[  seleccion, ]
w_test  <- wdbc_N[ -seleccion, ]

# Etiquetas
w_train_label <- w_train$diagnosis
w_test_label  <- w_test$diagnosis

Entrenar el modelo

w_predict<- knn( w_train[ , c( 1:( ncol( w_train ) - 1 ) ) ], 
                 w_test[  , c( 1:( ncol( w_test  ) - 1 ) ) ], 
                 cl = w_train_label, 
                 k = 21, 
                 prob = TRUE )
w_predict
  [1] M M B M M M B B M M B B M M M B M M B M B M B B B B B
 [28] B B B B B M B B B B B B B B B B M M M M M B M M B B B
 [55] B B M B B B M B B B B B M M M M M B B B B B M M B B B
 [82] B B M M B B B M B B B B B B M B M B B B B B M M B B M
[109] B B B B B M B B M B B M
 [ reached getOption("max.print") -- omitted 24 entries ]
attr(,"prob")
  [1] 1.0000000 1.0000000 0.9523810 1.0000000 0.9523810
  [6] 0.9047619 0.9047619 0.6666667 1.0000000 0.7619048
 [11] 1.0000000 1.0000000 1.0000000 0.9047619 0.9523810
 [16] 0.6190476 1.0000000 0.9523810 1.0000000 1.0000000
 [21] 1.0000000 0.6190476 1.0000000 0.8571429 0.9523810
 [26] 0.9523810 1.0000000 0.5238095 0.9523810 1.0000000
 [31] 0.9523810 0.9047619 0.7142857 1.0000000 0.5238095
 [36] 0.5714286 1.0000000 0.9047619 1.0000000 0.6190476
 [41] 1.0000000 1.0000000 1.0000000 0.7619048 0.9523810
 [46] 1.0000000 0.9523810 0.7619048 1.0000000 0.7619048
 [51] 0.8095238 0.9523810 0.7619048 0.9047619 0.9047619
 [56] 1.0000000 0.9047619 1.0000000 0.9047619 0.8095238
 [61] 1.0000000 1.0000000 1.0000000 1.0000000 0.9523810
 [66] 1.0000000 1.0000000 1.0000000 0.9523810 0.9523810
 [71] 0.7619048 0.9047619 0.9523810 1.0000000 1.0000000
 [76] 0.6818182 1.0000000 1.0000000 1.0000000 1.0000000
 [81] 0.8095238 1.0000000 1.0000000 0.7619048 1.0000000
 [86] 1.0000000 1.0000000 1.0000000 0.9047619 1.0000000
 [91] 0.9047619 1.0000000 1.0000000 1.0000000 1.0000000
 [96] 1.0000000 1.0000000 1.0000000 0.7619048 0.9523810
[101] 1.0000000 0.8571429 1.0000000 1.0000000 1.0000000
[106] 0.8095238 1.0000000 1.0000000 0.9523810 1.0000000
[111] 0.6190476 1.0000000 1.0000000 0.9523810 1.0000000
[116] 1.0000000 1.0000000 0.9047619 1.0000000 0.9523810
 [ reached getOption("max.print") -- omitted 24 entries ]
Levels: B M

Evaluar el modelo con tablas

Lo más sencillo, una tabla de contingencia con table:

table( Predicción = w_predict, Observación = w_test_label )
          Observación
Predicción  B  M
         B 90  5
         M  1 48

Una tabla de contingencia un poco más sofisticada con CrossTable

Calcula una tabla similar a la anterior pero más completa. Aporta proporciones de aciertos por filas, columnas y totales.

CrossTable( x = w_predict, 
            y = w_test_label, 
            dnn = c( "Predicción", "Observación" ),
            prop.chisq = FALSE )

 
   Cell Contents
|-------------------------|
|                       N |
|           N / Row Total |
|           N / Col Total |
|         N / Table Total |
|-------------------------|

 
Total Observations in Table:  144 

 
             | Observación 
  Predicción |         B |         M | Row Total | 
-------------|-----------|-----------|-----------|
           B |        90 |         5 |        95 | 
             |     0.947 |     0.053 |     0.660 | 
             |     0.989 |     0.094 |           | 
             |     0.625 |     0.035 |           | 
-------------|-----------|-----------|-----------|
           M |         1 |        48 |        49 | 
             |     0.020 |     0.980 |     0.340 | 
             |     0.011 |     0.906 |           | 
             |     0.007 |     0.333 |           | 
-------------|-----------|-----------|-----------|
Column Total |        91 |        53 |       144 | 
             |     0.632 |     0.368 |           | 
-------------|-----------|-----------|-----------|

 

Curvas ROC (Receiver Operating Characteristic) como evaluadores

Como vimos en el módulo anterior, es una representación gráfica de la \(sensibilidad\) frente a \(1 - especificidad\) para un clasificador binario, es decir sólo hay dos respuestas, positivo y negativo.

prob <- attr( w_predict, "prob" )
labels <- factor( w_test_label, 
                  labels = c( 1, 0 ) )

rocobj <- pROC::roc( labels, 
                     prob, 
                     auc = TRUE, 
                     ci = TRUE  )
Setting levels: control = 1, case = 0
Setting direction: controls > cases
print( rocobj )

Call:
roc.default(response = labels, predictor = prob, auc = TRUE,     ci = TRUE)

Data: prob in 91 controls (labels 1) > 53 cases (labels 0).
Area under the curve: 0.5522
95% CI: 0.4606-0.6439 (DeLong)
pROC::plot.roc( rocobj, 
                legacy.axes = TRUE, 
                print.thres = "best",
                print.auc   = TRUE, 
                auc.polygon = FALSE, 
                max.auc.polygon = FALSE,
                auc.polygon.col = "gainsboro", 
                col  = 2, 
                grid = TRUE )

Como se puede observar el la curva, el clasificador no parece ser muy bueno, ya que su AUC no es muy alta (0.5522496)

Vamos a repetir aumentando el número de vecinos a k=100.

w_predict<- knn( w_train[ , c( 1:( ncol(w_train) - 1 ) ) ], 
                 w_test[  , c( 1:( ncol(w_test)  - 1 ) ) ], 
                 cl = w_train_label, 
                 k = 100, 
                 prob = TRUE )

prob   <- attr( w_predict, "prob" )
labels <- factor( w_test_label, labels = c( 0, 1 ) )

rocobj <- pROC::roc( labels,
                     prob,
                     auc = TRUE, 
                     ci = TRUE  )
Setting levels: control = 0, case = 1
Setting direction: controls > cases
print( rocobj )

Call:
roc.default(response = labels, predictor = prob, auc = TRUE,     ci = TRUE)

Data: prob in 91 controls (labels 0) > 53 cases (labels 1).
Area under the curve: 0.7258
95% CI: 0.6362-0.8154 (DeLong)
pROC::plot.roc( rocobj,
                legacy.axes = TRUE,
                print.thres = "best",
                print.auc   = TRUE,
                auc.polygon = FALSE,
                max.auc.polygon = FALSE,
                auc.polygon.col = "gainsboro",
                col = 2, 
                grid = TRUE )

Ahora el AUC a aumentado a 0.7257931 por lo que al aumentar el número de vecinos, ha mejorado —en este caso— la eficiencia del predictor.

Validación cruzada, cross–validation

Como se ha comentado anteriormente, una vez que hemos generado el modelo, interesa conocer la precisión de este modelo en la predicción de los resultados. Dicho de otra manera, interesa estimar el error cometido en la predicción.

La validación cruzada aglutina una serie de técnicas, utilizadas para validar la precisión de un modelo predictivo. Para ello, comprueba la bondad predictiva del modelo sobre un conjunto de datos de prueba independientes de los datos empleados para generar el modelo.

Básicamente la secuencia sería:

  1. Construir el modelo con los datos de entrenamiento

  2. Aplicar el modelo con los datos de prueba, que serán independientes de los de entrenamiento

  3. Estimar el error en la predicción.

Hay distintas técnicas de validación cruzada para la estimación del error, aqui vamos a ver alguna de ellas.

Para todas las técnicas que vamos a ver, utilizaremos el paquete caret, uno de los más completos paquetes de machine learning que hay en R.

En primer lugar cargaremos y prepararemos los datos, esto nos servirá para todos los métodos que vamos a ver. A continuación iremos haciendo las validaciones con las diferentes técnicas.

La función CreateDataPartition nos permite crear un grupo de entrenamiento de forma similar a como lo hacíamos con sample.

set.seed(111)
train <- createDataPartition( iris$Species, 
                              p = 0.7, 
                              list = FALSE )
train.iris <- iris[ train, ]
test.iris  <- iris[ -train, ]

El método del conjunto de validación

En este caso la cuantificación del error se hace mediante la diferencia cuadrática media entre los valores observados y los predichos.

#Crear el modelo
model <- naiveBayes( Species ~ . , 
                     data = train.iris )

predictions <- predict( model, test.iris )
predictions <- as.numeric( predictions )

test <- as.numeric( as.factor( test.iris$Species ) )

data.frame( R2   = R2(   predictions, test ),
            RMSE = RMSE( predictions, test ),
            MAE  = MAE(  predictions, test ) )
         R2      RMSE        MAE
1 0.9363057 0.2108185 0.04444444
\(R^2\)
Es la correlación al cuadrado entre observaciones y predicciones. Cuanto más alto mejor será el modelo
RMSE, Root mean squared error
Error cuadrático medio, mide el error medio cometido por el modelo al predecir la clasificación de una observación. Cuanto más bajo sea, mejor será el modelo
MAE, Mean absolute Error
Error absoluto medio, es similar al RMSE, pero menos sensible a valores atípicos. Es la diferencia absoluta media entre observaciones y predicciones. Cuanto más bajo, mejor será el modelo

LOOCV, Leave one out cross validation

La validación cruzada dejando una observación fuera, consiste en dejar fuera del modelo una observación y contruir el modelo con el resto, a continuación comprobar el modelo con esa observación y anotar el error. Este proceso se repetirá con todas las observaciones. Finalmente se obtiene el error total cometido calculado como la media de los errores obtenidos para cada observación que se ha testado.

set.seed(111)
train <- trainControl( method = "LOOCV" )

model <- train( Species ~ ., 
                data = iris, 
                method = "naive_bayes", 
                trControl = train )
print( model )
Naive Bayes 

150 samples
  4 predictor
  3 classes: 'setosa', 'versicolor', 'virginica' 

No pre-processing
Resampling: Leave-One-Out Cross-Validation 
Summary of sample sizes: 149, 149, 149, 149, 149, 149, ... 
Resampling results across tuning parameters:

  usekernel  Accuracy   Kappa
  FALSE      0.9533333  0.93 
   TRUE      0.9600000  0.94 

Tuning parameter 'laplace' was held constant at a value
 of 0
Tuning parameter 'adjust' was held constant at
 a value of 1
Accuracy was used to select the optimal model using
 the largest value.
The final values used for the model were laplace =
 0, usekernel = TRUE and adjust = 1.

La ventaja del método es que utiliza todos los datos para entrenar y para testar el modelo.

La desventaja es que precisamente por utilizar todos los datos, tiene un alto coste computacional, especialmente cuando tenemos un conjunto de datos muy grande.

\(K\)-fold cross–validation y repeat \(K\)-fold cross-validation

Es la validación cruzada de \(k-\)iteraciones o de \(k-\)iteraciones repetidas.

En el primer caso, el conjunto de datos se divide en \(K\) subjuntos, uno de ellos se utiliza como datos de prueba y el resto \(K-1\) los datos de entrenamiento. La validación cruzada se realizará durante \(K\) iteraciones. En cada iteración un subconjunto será el conjunto de prueba y el resto de entrenamiento.

La idea es similar al método LOOCV, pero en lugar de hacerlo a nivel de observación, se hace a nivel de conjunto de observaciones.

El método de \(K\)-fold repetido es idéntico al método \(K\)-fold, sólo que en este caso, una vez terminadas las \(K\)-iteraciones, se vuelve a particionar en subconjuntos y se repite el proceso anterior. Estos ciclos se repetirán el número de veces que se indique.

\(K\)-fold

set.seed( 111 )
train <- trainControl( method = "cv", 
                       number = 10 )
model <- train( Species ~ ., 
                data=iris, 
                method = "naive_bayes",
                trControl = train )
print( model )
Naive Bayes 

150 samples
  4 predictor
  3 classes: 'setosa', 'versicolor', 'virginica' 

No pre-processing
Resampling: Cross-Validated (10 fold) 
Summary of sample sizes: 135, 135, 135, 135, 135, 135, ... 
Resampling results across tuning parameters:

  usekernel  Accuracy   Kappa
  FALSE      0.9533333  0.93 
   TRUE      0.9533333  0.93 

Tuning parameter 'laplace' was held constant at a value
 of 0
Tuning parameter 'adjust' was held constant at
 a value of 1
Accuracy was used to select the optimal model using
 the largest value.
The final values used for the model were laplace =
 0, usekernel = FALSE and adjust = 1.

Repeat \(K-\)fold

set.seed( 111 )
train <- trainControl( method = "repeatedcv", 
                       number = 10, 
                       repeats = 3 )
model <- train( Species ~ ., 
                data = iris, 
                method = "naive_bayes", 
                trControl = train )
print( model )
Naive Bayes 

150 samples
  4 predictor
  3 classes: 'setosa', 'versicolor', 'virginica' 

No pre-processing
Resampling: Cross-Validated (10 fold, repeated 3 times) 
Summary of sample sizes: 135, 135, 135, 135, 135, 135, ... 
Resampling results across tuning parameters:

  usekernel  Accuracy   Kappa    
  FALSE      0.9533333  0.9300000
   TRUE      0.9555556  0.9333333

Tuning parameter 'laplace' was held constant at a value
 of 0
Tuning parameter 'adjust' was held constant at
 a value of 1
Accuracy was used to select the optimal model using
 the largest value.
The final values used for the model were laplace =
 0, usekernel = TRUE and adjust = 1.

Manos a la obra

Ha llegado la hora de poner en práctica parte de lo aprendido en este módulo, para ello vamos ha crear un clasificador con un conjunto de datos biométricos de estudiantes de la Universidad de Murcia, y lo utilizaremos para intentar predecir el sexo de un individuo basado en sus datos biométricos, así que … manos a la obra.

Enlaces de interés

  • Un tutorial online en DataCamp (visitar).

  • Machine Learning Essentials de Alboukadel Kassambara es un buen libro práctico para iniciarse (visitar).

  • Una completa guía sobre el paquete Caret.

    • Caret Package – A Practical Guide to Machine Learning in R (visitar).
  • Un par de buenos libros al respecto:

    • Machine Learning with R. Brett Lantz. [PACKT]Publishing. 2013.
    • An Introduction to Statistical Learning: with Applications in R. Gareth James. (Springer Texts in Statistics). 2017.
  • Finalmente un enlace a un sitio donde aconsejan varios libros de ML en R (visitar)