U2A1 Clasificación lineal, cuadrática y de vecinos más próximos

Datos

En el ejemplo, se usan los datos del archivo de las muestras de cáncer mama.

Los datos del archivo Wisconsin18.RData proceden de un estudio sobre diagnóstico del cáncer de mama por imagen. Mediante una punción con aguja fina se extrae una muestra del tejido sospechoso de la paciente. La muestra se tiñe para resaltar los núcleos de las células y se determinan los límites exactos de los núcleos. Las variables consideradas corresponden a los valores medios de distintos aspectos de la forma de los núcleos de cada muestra.

Células

Tinción de los núcelos celulares

El fichero contiene un data.frame, llamado Wisconsin cuyas variables son:

Las 10 variables explicativas medidas en pacientes cuyos tumores fueron diagnosticados posteriormente.
La variable tipo que contiene el tipo de tumor (benigno o maligno).

Los datos son de acceso libre y fueron obtenidos de la universidad UCI (http://archive.ics.uci.edu/ml/datasets/Breast+Cancer+Wisconsin+%28Diagnostic%29 ), y son llamados “Wisconsin”, los datos consisten en:

ID number
Diagnosis (M = malignant, B = benign) 3-32

Se calculan diez características de valor real para cada núcleo celular

radius (mean of distances from center to points on the perimeter)
texture (standard deviation of gray-scale values)
perimeter
area
smoothness (local variation in radius lengths)
compactness (perimeter^2 / area - 1.0)
concavity (severity of concave portions of the contour)
concave points (number of concave portions of the contour)
symmetry
fractal dimension (“coastline approximation” - 1)

Importar datos

setwd("~/Estadistica Aplicada")
library(pacman)
p_load("MASS", "class")
load(url("http://verso.mat.uam.es/~joser.berrendero/datos/Wisconsin18.RData"))
head(Wisconsin)

##   radius texture perimeter  area smoothness compactness concavity concavepoints
## 1 13.540   14.36     87.46 566.3    0.09779     0.08129   0.06664      0.047810
## 2 13.080   15.71     85.63 520.0    0.10750     0.12700   0.04568      0.031100
## 3  9.504   12.44     60.34 273.9    0.10240     0.06492   0.02956      0.020760
## 4 13.030   18.42     82.61 523.8    0.08983     0.03766   0.02562      0.029230
## 5  8.196   16.84     51.71 201.9    0.08600     0.05943   0.01588      0.005917
## 6 12.050   14.63     78.04 449.3    0.10310     0.09092   0.06592      0.027490
##   symmetry fractal    tipo
## 1   0.1885 0.05766 benigno
## 2   0.1967 0.06811 benigno
## 3   0.1815 0.06905 benigno
## 4   0.1467 0.05863 benigno
## 5   0.1769 0.06503 benigno
## 6   0.1675 0.06043 benigno

¿Cuántos diagnósticos benignos y malignos tenemos?

table(Wisconsin$tipo)

## 
## benigno maligno 
##     357     212

Ahora, la forma en que se relacionan las variables por medio de una matriz de diagramas de dispersión

## COn el -11, podemos quitar la columna 11, y se colorea la columna entre benignos y malignos.
## Los benignos y malignos no se agrupan, por lo que son diferentes.
pairs(Wisconsin[,-11], col=Wisconsin$tipo)

Diagrama de estrellas

## Agrupa los datos que tenemos
stars(Wisconsin[-11], col.stars = Wisconsin$tipo, labels = NULL)

Se dispone en total de 569 imágenes, de las cuáles 357 corresponden a tumores benignos y 212 a malignos. Vemos que algunas variables están estrechamente relacionadas (por ejemplo, el radio, el perímetro y el área). También se observa que en términos generales, los tumores malignos corresponden a valores más altos de las variables.

Regla de clasificación lineal de Fisher

El comando básico es lda del paquete MASS. Tiene tres argumentos principales:

La expresión que determina la variable que queremos predecir (el tipo de tumor) y qué variables vamos a utilizar,
El data frame que contiene los datos,
El vector de probabilidades a priori de cada clase (por defecto, son las frecuencias relativas del vector que contiene las clases). Nosotros vamos a fijar probabilidades a priori iguales.

## Contiene los valores de las medias de cada una de las variables. 
resultado.lda <-lda(tipo ~ ., data=Wisconsin, prior = c(0.5,0.5) )

El objeto resultado.lda es una lista con todos los resultados del análisis. Veamos sus elementos más importantes. El vector resultado.lda$means contiene los vectores de medias x¯0 y x¯1 correspondientes a cada grupo:

## Con esto se pueden predecir los valoresy muestra que existe un valor medio para cada rango en benigno y maligno.
resultado.lda$means

##           radius  texture perimeter     area smoothness compactness  concavity
## benigno 12.14652 17.91476  78.07541 462.7902 0.09247765  0.08008462 0.04605762
## maligno 17.46283 21.60491 115.36538 978.3764 0.10289849  0.14518778 0.16077472
##         concavepoints symmetry    fractal
## benigno    0.02571741 0.174186 0.06286739
## maligno    0.08799000 0.192909 0.06268009

Como habíamos observado en las representaciones gráficas, los valores medios del grupo de tumores malignos son superiores a los del grupo de benignos. El vector resultado.lda$scaling contiene los coeficientes de la función discriminante lineal de Fisher:

## Sirven para predecir cuál es benigno o maligno.
resultado.lda$scaling

##                        LD1
## radius         2.173832578
## texture        0.097479319
## perimeter     -0.243883158
## area          -0.004235635
## smoothness     8.610211091
## compactness    0.431476344
## concavity      3.592356858
## concavepoints 28.529778564
## symmetry       4.489073661
## fractal       -0.529214778

Análisis Discriminante Lineal (ADL, o LDA por sus siglas en inglés) es una generalización del discriminante lineal de Fisher, un método utilizado en estadística, reconocimiento de patrones y aprendizaje de máquinas para encontrar una combinación lineal de rasgos que caracterizan o separan dos o más clases de objetos o eventos. La combinación resultante puede ser utilizada como un clasificador lineal, o, más comúnmente, para la reducción de dimensiones antes de la posterior clasificación.

Representando los datos de escalamiento de forma gráfica

barplot(as.vector(resultado.lda$scaling), names.arg = names(Wisconsin[-11]))

Si llamamos w al vector anterior, clasificamos una observación x en el grupo 1 siempre que

\[ w^\top \left(x-\frac{\bar{x}_0 + \bar{x}_1}{2}\right) > 0, \]

y en el grupo 0 en caso contrario. El valor del término de la izquierda de la desigualdad es la puntuación discriminante de x. Al aplicar plot a resultado.lda obtenemos histogramas de las puntuaciones discriminantes estandarizadas de las observaciones de cada grupo:

## Con esto, podemos ver como es que tienden a ser los tipos. Sin embargo, existe un punto en el que no se sabe diferenciar entre ambos, que es en el rango de -2 y 1.
plot(resultado.lda)

Vemos que, efectivamente, las puntuaciones en el grupo 0 tienden a ser negativas y en el grupo 1 tienden a ser positivas. Sin embargo, hay una zona de solapamiento de las puntuaciones discriminantes que llevará a un cierto porcentaje de errores de clasificación. Con el fin de calcular la tasa de error aparente, usamos el comando predict para clasificar los datos de la muestra y luego contamos la proporción de veces que nos hemos equivocado:

## Calcula la tasa de error
n <- resultado.lda$N
tipo.prediccion <- predict(resultado.lda)$class
1 - sum(tipo.prediccion == Wisconsin$tipo)/n

## [1] 0.05975395

Esto significa que nos hemos equivocado para aproximadamente un 5.97 % de las observaciones. Para calcular la tasa de error por validación cruzada, seleccionamos CV=TRUE en el comando lda, de la siguiente forma:

tipo.prediccion.cv <- lda (tipo ~ ., data=Wisconsin, prior = c(0.5,0.5),CV = TRUE)$class
1 - sum(tipo.prediccion.cv == Wisconsin$tipo)/n

## [1] 0.06326889

Si queremos clasificar nuevos vectores de observaciones, tenemos que crear previamente un data frame que los contenga y luego usar de nuevo el comando predict. El siguiente código genera aleatoriamente dos vectores de observaciones y los clasifica. También permite obtener las dos puntuaciones discriminantes correspondientes:

## Toma 10 valores al aleatorios con distribución normal 
x1 <- rnorm(10)
x2 <- rnorm(10)
nuevas.obs <- data.frame(rbind(x1,x2))
names(nuevas.obs)<- names(Wisconsin[1:10])
nuevas.obs

##       radius    texture   perimeter       area smoothness   compactness
## x1 0.3828126 -0.2998112 -0.02500274 -1.6613506   1.129039  0.0003731669
## x2 0.3798711  0.8058251 -0.17276578  0.7097629  -1.805693 -1.6135927767
##     concavity concavepoints   symmetry     fractal
## x1 -0.2453729     0.7131542  0.6293596  0.02074401
## x2  1.7958452    -2.2784414 -1.5159248 -1.36051351

Resultado de la clasificación

## Hace una predicción a partir de los datos
predict(resultado.lda, nuevas.obs)$class

## [1] maligno benigno
## Levels: benigno maligno

Puntuaciones discriminantes

## Hace el resultado a partir de una discriminación lineal
predict(resultado.lda, nuevas.obs)$x

##          LD1
## x1  21.67399
## x2 -91.07982

Ejercicios

Repite los cálculos utilizando únicamente las variables smoothness y concavepoints.

library(dplyr)

## 
## Attaching package: 'dplyr'

## The following object is masked from 'package:MASS':
## 
##     select

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

## Toma los datos a utilizar, que son smoothness y concavepoints
SCP <- select(Wisconsin,smoothness,concavepoints, tipo)

## Forma en que se relacionan smoothness y concavepoints por medio de una matriz de diagramas de dispersión.
## Se quita la 3ra columna o variable a utilizar que es tipo
pairs(SCP[-3], col=Wisconsin$tipo)

## Diagrama de estrellas
stars(SCP[-3])

## Contiene los valores de las medias de cada una de las variables. 
resultado.lda <-lda(tipo ~ ., data=SCP, prior = c(0.5,0.5) )
resultado.lda$means

##         smoothness concavepoints
## benigno 0.09247765    0.02571741
## maligno 0.10289849    0.08799000

resultado.lda$scaling

##                     LD1
## smoothness    -15.01049
## concavepoints  44.01965

## Se muestra el escalamiento en su forma gráfica
barplot(as.vector(resultado.lda$scaling), names.arg = names(SCP[-3]))

## Con esto, podemos ver como hacia donde tienden a ir los tipo, aunque, existe un punto en el ambos están en rango muy parecido, esto cuando casi tocan a 0. 
plot(resultado.lda)

A simple vista, se puede decir que el tipo benigno se agrupa hacia la parte izquierda y el maligno un poco más a la derecha, pero como se dijo, cuando casi tocan 0 ambos tienen un parecido, por ello se analizará el nivel de la tasa de error.

n <- resultado.lda$N
tipo.prediccion <- predict(resultado.lda)$class
1 - sum(tipo.prediccion == SCP$tipo)/n

## [1] 0.08787346

Se observa una tasa de error del 8.87%. A continuación se calcula la tasa de error por validación cruzada:

tipo.prediccion.cv <-lda(tipo ~ ., data=SCP, prior = c(0.5,0.5), CV =TRUE )$class
1 - sum(tipo.prediccion.cv == SCP$tipo)/n

## [1] 0.08787346

##  Genera aleatoriamente dos vectores de observaciones y los clasifica. También permite obtener las dos puntuaciones discriminantes correspondientes.
x1 <- rnorm(2)
x2 <- rnorm(2)
nuevas.obs <- data.frame(rbind(x1,x2))
names(nuevas.obs)<- names(SCP[1:2])
nuevas.obs

##    smoothness concavepoints
## x1  0.1615186    -0.4782110
## x2  0.6142505    -0.7416231