Clasificación binaria….Naive Bayes

Un poco sobre análisis de textos

Lo principal de esta entrada es dar un ejemplo de como clasificar correos, esto lo hago siguiendo el ejemplo de libro Machine Learning for Hackers [1]. La técnica es conocida como Naive Bayes, pero creo conveniente comentar un poco de la librería de R para text mining.

La librería tm   permite realizar text mining en R project. El text mining se aboca sobre la conversión textos a datos que puedan ser analizados, ya sea por herramientas estadísticas, o por técnicas de procesamiento natural del lenguaje (NLP, en ingles).

La librería es bastante conocida en la red y existe mucha información, en varios blog se muestran ejemplos interesantes. En particular se mostraban análisis de textos de Twitter, un buen ejemplo es el de Yanchang Zhaog, este se puede encontrar en sus notas data mining con R. Lo desafortunado es que la API de twitter fue actualizada y ya no se puede replicar el ejemplo con el código de Yanchang.

Para ejemplificar el funcionamiento de tm tomo los títulos de las primeras planas de un periódico desde el día 1 al 23 de Marzo.

Los pasos son los siguientes:

  1. Extraemos los HTML de cada portada desde la página del periódico.
  2. Tomo el texto de la página web o HTML.
  3. Le quito a los textos la mayoría de palabras que forman parte de los títulos y secciones del periódico.
  4. Proceso con tm los archivos.
  5. Construyo una nubes de palabras.

Podríamos hacer otro tipo de análisis  después del punto 4, pero implica hablar de otros conceptos, como cluster, PCA  o redes. Observación: No pongo todo el código, pero prácticamente solo falta el proceso del paso 2.

#Paso 0
#Cargamos las librerías 
library(tm)
library(wordcloud)
library(stringi)

#Paso 1
#Extraemos los datos
direc<-""http://www.eluniversal.com.mx/hemeroteca/edicion_impresa_201503"
comp_dir<-".html"
dias<-c("01", "02", "03", "04", "05", "06", "07", "08", "09", "10", "11", "12", "13", "14", "15", "16", "17", "18", "19", "20", "21", "22", "23")
#Paso 2
#Extraemos los datos
edi_impresa<-character(length(dias))

for (i in 1:length(dias)){
     edi_impresa[i]=htmlToText(stri_paste(direc,dias[i],comp_dir))
    }
#Paso 3 y 4
#Construimos el Corpus
Doc<-Corpus(VectorSource(edi_impresa))

#Visualizamos Doc
Doc[[1]]

#Procesamos Doc con la librería tm

Doc1<-tm_map(Doc, function(x) stri_replace_all_regex(as.character(x), "<.+?>", " "))
Doc2 <- tm_map(Doc1, function(x) stri_replace_all_fixed(x, "\t", " "))
Doc3 <- tm_map(Doc2, PlainTextDocument)
Doc4 <- tm_map(Doc3, stripWhitespace)
Doc5 <- tm_map(Doc4, removeWords, stopwords("spanish"))
Doc6 <- tm_map(Doc5, removePunctuation)
Doc7<-tm_map(Doc6,removeNumbers)

#Hacemos una lista con todas las palabras de títulos o secciones del periódico
L<-#LO OMITO
Doc8 <- tm_map(Doc7, removeWords, L)
Doc9 <- tm_map(Doc8,content_transformer(tolower))

#Contruimos nuestra matriz de TDM
DocsTDM <- TermDocumentMatrix(Docs14)

Lo que hice hasta esta parte del código es la construcción de la matriz de términos de los documentos. Con ello puedo revisar las frecuencias de las palabras. Como observación, debido a que no hice bien la lista de palabras a remover aparecerán algunas que se refieren a títulos y secciones del periódico.

#Revisamos las frecuencias de las palabras
Grupo1<-rowSums(as.matrix(DocsTDM))
Grupo1<-subset(Grupo1,Grupo1>=25)

#Graficamos nuestras palabras
library(ggplot2)
qplot(names(Grupo1), Grupo1, geom="bar", xlab="Terms") + coord_flip()
barplot(Grupo1, las=2)

 Palabras_más_frecuentes

Por último hago la nube de palabras.

#Nube de palabras
Grupo1<-subset(Grupo1,Grupo1>=10)
m<-as.matrix(Grupo1)
word_freqs = sort(rowSums(m), decreasing=TRUE)
dm = data.frame(word=names(word_freqs), freq=word_freqs)
#Nube de palabras
wordcloud(dm$word, dm$freq, random.order=FALSE, colors=brewer.pal(8, "Dark2"))

Imagen_datos_nube

Bueno, la idea sería procesar correctamente los datos y presentar la nube de palabras con mayor frecuencia de apariciones en el texto o con algún criterio que de deseé. Otro ejemplo sería revisar cuales palabras persistieron por el paso de los días en el mes de Marzo. También se puede hacer otro tipo de análisis, por ejemplo buscar cluster  o tópicos en los datos, lo cual resulta más interesante.

Cómo observación, la librería tm cuenta con aproximadamente 54 funciones y cabe mencionar que el objeto principal es la construcción y procesamiento de Corpus, que es como un concentrado de todos los textos  que se analizan. Para procesar el contenido del Corpus se usa la función tm_map, la cual funciona análogamente a las funciones apply, lapply, sapplycon la diferencia de que se aplica solo a los Corpus. Para mayor detalle se puede revisar el manual de la librería [].

Si se tiene interés en aprender sobre procesamiento de textos siguiendo las metodologías de NLP, se tienen muy buenos cursos de la plataforma Coursera. Recomiendo el curso impartido por la Universidad de Stanford.

¿Web Scraping con R?

El procesar textos ha sido una área de investigación desde hace tiempo, ahora se tienen librerías y metodologías estándar, pero donde abunda mucho información para analizar y era complicado extraerla es en la web. Ha este tipo de programas se les llama web scraping.

Mario Annau liberó una librería en R llamada tm.plugin.webminig y por otro lado Hadley Wickham liberó la librería rvest. Las dos librerías se toman como modelo librerías de Python, las cuales son más robustas y se han mejorado con el tiempo.

No hago algún ejemplo, pero se puede revisar la documentación y replicar el ejemplo que comparten cada uno sobre el uso de su librería. Lo interesante es que las dos tratan de simplificar pasos que tm no hace por si sola, principalmente por la interacción con textos HTML.

Como recomendación ( y eso en cada entrada de esta categoría) creo importante revisar el repositorio de Hahley Wickham, el cual a contribuido a la comunidad de R project con librerías muy populares y de mucho utilidad, basta mencionar ggplot2 y reshape. Se ha convertido en el RockStar de R para muchos.

Clasificación de mensajes-detectando Spam

 

En la entrad sobre Máquina de Soporte Vectorial (SVM) hice un ejemplo comparativo para clasificar correos. Replico el código en esta entrada solo enfocándome en Naive Bayes.

#Se instala la librería kernlab para usar sus datos 
library(kernlab)
#Se extraen los datos
data(spam)
#Se genera una muestra de 100 índices
index=sample(1:nrow(spam),100)
#Se carga la librería e1071 para usar la función naiveBayes
library(e1071)
model=naiveBayes(type~.,data=spam[-index,])
model=naiveBayes(type~.,data=spam[-index,])
pred=predict(mol,spam[index,])
table(pred,spam[index,58])
 
pred   nonspam   spam
nonspam  38      0
 spam     29     33

Los datos que se usan en el código anterior son datos ya procesados, los cuales ya cuentan con la etiqueta o variable categórica de ser Spam o no. Lo recomendable es  cargar los datos y explorar la información, en este caso solo hago uso de los algoritmos Naive Bayes que se encuentra en la librería e1071.

Hago otro ejemplo con los datos de la librería kernlab, estos corresponden a 8993 personas de las cuales se toman 14 datos y uso el algoritmo de Naive Bayes para evaluar la clasificación de 1000 de estas personas por sexo.

#Se carga las librerías
library(kernlab)
library(e1071)
#Se extraen los datos
data("income")
#Se toman 1000 personas de los datos
index=sample(1:nrow(spam),1000)
#Se genera el modelo Naive Bayes de la librería e1071
modelo=naiveBayes(SEX~.,data=income[-index,])
#Se estima la predicción
pred=predict(modelo,income[index,])
#Se presentan los datos en modo de una tabla
table(pred,income[index,2])
 
pred  M   F
 M   234 162
 F   220 384

En este último ejemplo, al igual que el primero se hace uso de datos previamente procesados. Lo cual facilita la aplicación del algoritmo.

Las tablas de los dos ejemplos muestra dos categoría, los errores del algoritmo se pueden considerar como la cantidad de Hombres de la fila M de Male que son considerados como mujeres en la columna F de Famale, es decir; 384. Equivalente las 220 mujeres que son clasificadas como hombres en la segunda fila.

Estos errores son conocidos como falsos negativos y falsos positivos, el objetivo de las técnicas de clasificación es tener el mínimo de este tipo de errores.

El siguiente ejemplo consiste  también en clasificar correos, donde un conjunto son identificados como Spam y otro conjunto se sabe que no es spam. Estos son extraídos desde SpamAssasin, la dirección es http://spamassassin.apache.org/publiccorpus/ . Pero recomiendo bajar todos los datos y código libro Machine Learning for Hackers  dese el repositorio de Jhon M. White en Github, para replicar de manera completa el ejemplo.

A diferencia de los primero dos ejemplos de clasificación, en este caso los datos o textos se van a procesar para identificar las palabras que aparecen con mayor frecuencia en los Spam y las palabras de los que no son Spam.

Se hará uso de 3 funciones, la primera para extraer el mensaje del correo, la segunda para convertir el mensaje en una matriz de términos del documento, con la cual se revisa la frecuencias de las palabras y se aplica una tercera función para clasificar cualquier correo con respecto a las palabras que se encuentran con mayor frecuencia en los correos identificados como Spam.

Naive Bayes

Antes doy una breve explicación  de la idea a tras de la clasificación. Suponiendo que se tienen dos clase, hombres y mujeres; y una lista de personas, por ejemplo: Guadalupe, Ernesto, Sofía, Andrea, Guadalupe, Sofía, Luis, Miguel, Guadalupe y Jaime

Nombre     Sexo
Guadalupe- Mujer
Ernesto  - Hombre
Sofía    - Mujer
Andrea   - Mujer
Guadalupe- Hombre
Sofía    - Mujer
Luis     - Hombre
Miguel   - Hombre
Guadalupe- Mujer
Jaime    - Hombre

De la lista anterior al considerar dos categorías, hombres y mujeres, se puede tomar una persona, ejemplo Guadalupe,  pero uno puede preguntarse ¿cuál es la probabilidad de que sea mujer dado que la persona se llama Guadalupe?

P(Mujer|Guadalupe)=[P(Guadalupe|Mujer)P(Mujer)]/P(Guadalupe)

La ecuación anterior expresa P( | ) como la probabilidad condicional y la regla de Bayes. Con los datos del ejemplo anterior  se puede calcular el valor de P(Mujer|Guadalupe)=2/3.

Análogamente se puede calcular P(Hombre|Guadalupe)=1/3

Entonces  lo que se tiene es que tomando a una persona que se llama Guadalupe es más probable que sea mujer. Pero se pueden considerar más atributos, como edad, peso, color de ojos, etc. Así  se puede uno volver a preguntar, ¿cuál es la probabilidad de que sea Mujer dado que se llama la persona Guadalupe?, pero claro ahora se pueden considerar los nuevos atributos en el calculo de la probabilidad condicional.

Lo que se modificará en el calculo es la cantidad de factores, ejemplo:

P(Guadalupe|Mujer)=P(Peso|Mujer)*P(Color de ojos|Mujer)*P(Edad|Mujer)

En el ejemplo de los correos se tiene una lista de palabras que se encuentran en aquellos que son Spam, entonces el total de atributos será dado por el total de palabras. La pregunta que se hace es: ¿el correo es Spam dado que tienen estas palabras?

Lo que  no menciono en el calculo anterior es que los atributos  se consideran independientes condicionalmente. Los aspectos teóricos de manera breve se pueden consultar en la parte 4 de las notas de Andrew Ng o en las notas de Andrew Moore, con más detalles se pueden encontrar en la referencias [2,3,4].

Espero que la explicación pese a ser corta ayude a entender lo que se va hacer con los correos que se procesan.

Correos..¿spam o no?

Como antes mencioné, tomo un conjunto de correos identificados como Spam y otro que se tienen identificado que no lo es. Se procesan para formar una matriz de términos para cada conjunto de correos. Con la matriz se asignarán densidades a las palabras  y por último se construye una función para clasificar.

Los primeros 500 correos se consideran como el conjunto de entrenamiento para el proceso de clasificación.

#Cargamos los datos
library(tm)
library(ggplo2)
spam.path <- file.path("data", "spam")
easyham.path <- file.path("data", "easy_ham")

#Funciones que usaremos
######Primera función

get.msg <- function(path)
{
 con <- file(path, open = "rt")
 text <- readLines(con)
 msg <- text[seq(which(text =="")[1]+1, length(text), 1)]
 close(con)
 return(paste(msg, collapse = "\n"))
}

#####Segunda función

get.tdm2 <- function(doc.vec)
{
 library(stringi)
 doc.corpus<-Corpus(VectorSource(doc.vec))
 cor.doc<-tm_map(doc.corpus, function(x) stri_replace_all_regex(as.character(x), "<.+?>", " "))
 cor.doc <- tm_map(cor.doc, PlainTextDocument)
 cor.doc <- tm_map(cor.doc,content_transformer(tolower))
 cor.doc <- tm_map(cor.doc,removeWords,stopwords("english"))
 cor.doc<- tm_map(cor.doc, removePunctuation)
 cor.doc<- tm_map(cor.doc, removeNumbers)
 cor.doc<-tm_map(cor.doc, stripWhitespace)
 doc.dtm<- TermDocumentMatrix(cor.doc)
 return(doc.dtm)
}

######Tercer función
classify.email <- function(path, training.df, prior = 0.5, c = 1e-6)
{
 #Extrae el mensaje
 #Convierte los datos en una matrz de terminos del documento

 msg <- get.msg(path)
 msg.tdm <- get.tdm(msg)
 msg.freq <- rowSums(as.matrix(msg.tdm))
 msg.match <- intersect(names(msg.freq), training.df$términos)
 if(length(msg.match) < 1)
 {
 return(prior * c ^ (length(msg.freq)))
 }
 else
 {
 match.probs <- training.df$ocurrencia[match(msg.match, training.df$términos)]
 return(prior * prod(match.probs) * c ^ (length(msg.freq) - length(msg.match)))
 }
}

######Procesar los los Spam

spam.docs <- dir(spam.path)
spam.docs <- spam.docs[which(spam.docs != "cmds")]
all.spam <- sapply(spam.docs,function(p) get.msg(file.path(spam.path, p)))
spam.tdm <- get.tdm2(all.spam)

######Procesamos 500 correos que no son spam

easyham.docs<-dir(easyham.path)
easyhamo.docs500=easyham.docs[1:500]
easyham500<-sapply(easyham.docs500,function(p)get.msg(file.path(easyham.path,p))
easyham500.tdm<-get.tdm2(easyham500)

En el código anterior se obtiene las matrices de términos y se puede ver los valores de esta en R, solo colocando el nombre:

#Revisamos los valores de los corpus

spam.tdm

<<TermDocumentMatrix (terms: 23109, documents: 500)>>
Non-/sparse entries: 70715/11483785
Sparsity : 99%
Maximal term length: 298
Weighting : term frequency (tf)

#Corpus para correos que no son Spam
easyham500.tdm

<<TermDocumentMatrix (terms: 13118, documents: 500)>>
Non-/sparse entries: 50247/6508753
Sparsity : 99%
Maximal term length: 245
Weighting : term frequency (tf)

Con las matrices de términos se construye una matriz con las densidades de las palabras.

#Matriz de densidad

spam.matrix <- as.matrix(spam.tdm)
spam.counts <- rowSums(spam.matrix)
spam.df <- data.frame(cbind(names(spam.counts),as.numeric(spam.counts)),stringsAsFactors = FALSE)
names(spam.df) <- c("términos", "frecuencia")
spam.df$frecuencia <- as.numeric(spam.df$frecuencia)
spam.occurrence <- sapply(1:nrow(spam.matrix),
 function(i)
 {
 length(which(spam.matrix[i, ] > 0)) / ncol(spam.matrix)
 })
spam.density <- spam.df$frecuencia / sum(spam.df$frecuencia)

#Aplicamos el mismo proceso para easyham500.tdm

#Términos frecuentes en spam.df

head(spam.df[with(spam.df, order(-ocurrencia)),],n=15)
      términos frecuencia density ocurrencia
6019    email    837    0.007295962 0.574
15199   please   459    0.004001011 0.520
3539    click    350    0.003050880 0.454
11722   list     419    0.003652339 0.442
21468   will     843    0.007348262 0.442
7317    free     651    0.005674637 0.420
9990  information 364   0.003172915 0.374
13837   now      329    0.002867827 0.352
2899    can      518    0.004515302 0.344
7687    get      425    0.003704640 0.334
14323   one      376    0.003277517 0.322
13511   new      336    0.002928845 0.300
16452   receive  327    0.002850394 0.298
19437   time     316    0.002754509 0.296
12527   message  245    0.002135616 0.286

#Términos frecuentes en easy.df

 head(easy.df[with(easy.df, order(-ocurrencia)),],n=15)
       términos frecuencia density ocurrencia
5052    group   232      0.003446125 0.388
12391   use     272      0.004040284 0.380
12979   wrote   238      0.003535249 0.380
1577    can     348      0.005169187 0.368
6982    list    249      0.003698642 0.368
8258    one     358      0.005317727 0.338
6527    just    273      0.004055138 0.326
8120    now     231      0.003431271 0.324
4812    get     230      0.003416417 0.282
6924    like    232      0.003446125 0.282
3653    email   188      0.002792549 0.276
11337   subject 162      0.002406346 0.270
11854   time    188      0.002792549 0.258
12834   will    318      0.004723567 0.254
6112   information 162   0.002406346 0.232

 De los datos anteriores se observa que comparten los dos tipos de correos algunas palabras, pero se aprecia que las ocurrencias en los Spam tienen valores más altos.

Lo que sigue es comparar las palabras contra una muestra de correos que no son Spam, pero por el texto que contiene son más complicados para detectarse como no-Spam.

Usamos la función classify.email se prueba contra el conjunto de entrenamiento que son las matrices spam.df y easy.df

#Clasificación
hardham.path <- file.path("data", "hard_ham")
hardham.docs <- dir(hardham.path)
hardham.docs <- hardham.docs[which(hardham.docs != "cmds")]

hardham.spamtest <- sapply(hardham.docs, function(p) classify.email(file.path(hardham.path, p), training.df = spam.df))
 
hardham.hamtest <- sapply(hardham.docs,function(p) classify.email(file.path(hardham.path, p), training.df = easyham.df))

hardham.res <- ifelse(hardham.spamtest > hardham.hamtest,TRUE,FALSE)

summary(hardham.res)
   Mode FALSE TRUE NA's 
logical  242   7   0

Lo que se observa es que se clasifican 7 email como Spam, cuando no lo son. Es posible que se deba algún error en el procesamiento del TDM de los Spam, ya que el resultado obtenido es demasiado bueno.

Para cerrar esta entrada, hago lo que sería un paso de prueba del clasificador. Se inicia con un conjunto de entrenamiento y se prueba como responde para clasificar  2500 correos que no son Spam y que son fáciles de identificar (easyham).

#Función clasificador de Spam

spam.classifier <- function(path)
{
 pr.spam <- classify.email(path, spam.df)
 pr.ham <- classify.email(path, easy.df)
 return(c(pr.spam, pr.ham, ifelse(pr.spam > pr.ham, 1, 0)))
}
#Cargamos los 2500 correos

easyham.docs <- dir(easyham.path)
easyham.docs <- easyham.docs[which(easyham.docs != "cmds")]

#Los clasificamos

easyham.class <- suppressWarnings(lapply(easyham.docs,function(p){spam.classifier(file.path(easyham.path, p))}))

#Contruimos una matriz con los datos
easyham.matrix <- do.call(rbind, easyham.class)
easyham.final <- cbind(easyham.matrix, "EASYHAM")

class.matrix <- rbind(easyham.final)
class.df <- data.frame(class.matrix, stringsAsFactors = FALSE)

names(class.df) <- c("Pr.SPAM" ,"Pr.HAM", "Class", "Type")
class.df$Pr.SPAM <- as.numeric(class.df$Pr.SPAM)
class.df$Pr.HAM <- as.numeric(class.df$Pr.HAM)
class.df$Class <- as.logical(as.numeric(class.df$Class))
class.df$Type <- as.factor(class.df$Type)

#Hacemos la gráfica con la información de los correos clasificados

ggplot(class.df, aes(x = log(Pr.HAM), log(Pr.SPAM))) +
 geom_point(aes(shape = Type, alpha = 0.5)) +
 stat_abline(yintercept = 0, slope = 1) +
 scale_shape_manual(values = c("EASYHAM" =1),name = "Email Type") +
 scale_alpha(guide = "none") +
 xlab("log[Pr(HAM)]") +
 ylab("log[Pr(SPAM)]") +
 theme_bw() +
 theme(axis.text.x = element_blank(), axis.text.y = element_blank())

#Hacemos una tabla con las densidades 

easyham.col <- get.results(subset(class.df, Type == "EASYHAM")$Class)
class.res <- rbind(easyham.col)
colnames(class.res) <- c("NOT SPAM", "SPAM")
print(class.res)

          NOT SPAM SPAM
easyham.col 0.974 0.026

Se observamos que el clasificador identifica aproximadamente el 97.5% de los correos que no son Spam, lo cual es correcto y solo tienen un error (falsos positivos) de 2.5%. Este ejercicio se puede repetir con el resto de conjuntos de correos que se encuentran en los datos del texto Machine Learning for Hackers.

Grafica_Ham_vs_Spam

Espero que los ejemplos den una idea aproximada de lo que se hace con este algoritmo y como se usa para construir un clasificador, en el texto se tienen más detalles, pero en general traté de poner en esta entrada las ideas básicas.

También se debe de observar que la clasificación no se reduce solo a tener dos clases (Spam y no Spam), como mencioné en otras entradas existen varias técnicas de clasificación multi-clase y en muchas situaciones es probable que se tengan más de dos clases para analizar.

Referencias:

1.-http://shop.oreilly.com/product/0636920018483.do

2.-http://shop.oreilly.com/product/0636920018483.do

3.-http://www.amazon.com/Artificial-Intelligence-Modern-Approach-Edition/dp/0136042597

4.-http://statweb.stanford.edu/~tibs/ElemStatLearn/

5.-http://www.autonlab.org/tutorials/naive02.pdf

6.-http://cs229.stanford.edu/

Análisis de Cluster, un ejemplo sencillo.

El análisis de agrupamiento (cluster), es sencillo de explicar y de implementar en R project. Lo delicado es procesar la información y la definición de la métrica que consideramos apropiada para nuestros datos o experimento.

Hago un ejemplo sencillo antes de explicar qué se hace con esta técnica.

#Ejemplo de Cluster
data(iris)
#Generamos 5 grupos en nuestros datos
modelo5<-kmeans(iris[1:2],5)
plot(iris[1:2],col=modelo5$cluster,main="Ejemplo de 5 cluster",xlab="longitud de la Sápal",ylab="Ancho de la Sápal")
#Generamos 3 grupos en nuestros datos
modelo3<-kmeans(iris[1:4],3)
plot(iris[1:4],col=modelo3$cluster,main="Ejemplo de 3 cluster")

 Ejemplo_5_cluster

Ejemplo_3_cluster

Los dos gráficos muestran de colores distintos un conjunto de puntos, por lo que se ve en el código, usando una función kmeans para definir modelo5,modelo3. Luego indico al gráfico que tomara los colores partiendo de un valor $cluster. Esto podríamos deducirlo viendo solo el  código, la pregunta es cómo se define qué puntos pintar de cierto color. Si uno hace el “scatter plot” de los mismos datos se observa una “regla”.

#Scatter Plot 
data(iris)
plot(iris[1:4],col=c(2,3),main="Scatter Plot de parejas de vectores de datos")

Ejemplo_Scatter_Plot_Datos Iris

Si uno compara los 3 gráficos, la regla para determinar el color de los puntos tienen que ver con un modo de agruparlos, de tener un modo de verlos cuan cercanos o similares son. Eso se nota principalmente entre los dos últimos gráficos.

La idea de “similar”, puede ser pensada como un modo de comparar cosas, si se piensa en tres personas  donde dos estudiaron en la misma universidad, tienen casi la misma edad y alturas y tienen diferencia de peso de solo un par de kilos, son “casi iguales” con respecto al tercer individuo, que no estudió en la misma universidad, tiene más años, tiene sobrepeso y mide más que los otros dos. Entonces podemos pensar que los dos primeros son “similares” o al compararlos “son más cercanos o parecidos”.

La noción de “cercanía” está asociada con la noción de distancia, podemos definir que dos cafeterías tienen una distancia pequeña si las separa una calle, pensando que su unidad de medida son el número de calles entre ellas. Entonces la noción a través de todo esto es la de “métrica“, cómo medir similitudes.

En geometría, la métrica es la distancia entre puntos la cual tienen ciertas propiedades; medir del punto 1 al punto 2 es igual que iniciar midiendo desde el punto 2 hasta el punto 1 (simetría); la distancia entre los puntos es por lo menos cero (positiva) y  medir la distancia entre dos puntos es menor o igual que la suma de las distancia de los puntos a un punto intermedio (desigualdad del tringualo) y la última, si la distancia entre dos puntos es cero los puntos no pueden ser distintos (identidad de los indiscernibles).

Entonces la idea de la técnica es medir la distancia de cada pareja de puntos de los datos. Esto no necesariamente implica que la distancia sea para números o vectores, puede ser definida por las características de los datos o sobre cadenas de texto o mensajes, características de una persona,etc.

Aclarado el código anterior, lo que se hizo con el comando kmeans fue aplicar un algoritmo que permite separar entre los grupos de datos, tomando como parte inicial el  eligir la cantidad de grupos que se desean, después se calcula el centroide con respecto a todos los datos y se hace un proceso de separación e iteracciones para determinar los cúmulos, la técnica se llama k-Means Clustering. Para entender en detalles el algoritmo se puede revisar las notas de Andrew Ng o de Andrew Moore, las notas de éste último son acompañadas de varias  gráficas lo cual ayuda a visualizar el algoritmo de un modo sencillo.

Un ejemplo puede ser pensando en una matriz de datos, donde cada fila representa un individuo y cada columnas es una característica la cual puede tomar el valor de -1,0 o 1. Entonces lo que se hace es realizar el producto de la matriz y su transpuesta para obtener la correlación entre cada uno de los individuos y lo que se usa como métrica es la distancia euclidiana. Con ellos se puede tener una matriz con información sobre las distancia entre nuestros individuos.

#Ejemplos 
x.matriz<-matrix(sample(c(-1,0,1),24,replace=TRUE),nrow=4,ncol=6)
#Etiquetamos las columnas y las filas
row.names(ex.matriz)<-c('A','B','C','D')
colnames(ex.matriz)<-c('P1','P2','P3','P4','P5','P6')
#Multiplicamos por la transpuesta
ex.mult<-ex.matriz%*%t(ex.matriz)

#Calculamos la distancia
ex.dis<-dist(ex.mult)

Ahora usamos el comando cmdscale para calcular el escalamiento multidimensional, obtendremos como resultado los puntos que son cercanos o próximos, según la métrica.

#Cluster
ex.mds<-cmdscale(ex.dis)
plot(ex.mds,type='n')
text(ex.mds,c('A','B','C','D'))

Ejemplo_Escalado_Multidimensional

Se concluye del ejemplo anterior que la estrategia es; procesar nuestra información, definir o usar una métrica, construir la matriz de similaridad y aplicar el análisis de escalamiento multidimensional o un algoritmo de agrupamiento.

Polarización de un Senado

Muestro un ejemplo del texto Machine Learning for Hackers [1]. El ejemplo usa los datos del senado de EEUU, se encuentran en la plataforma http://www.voteview.com y en caso de tener una falla en la descargar se puede bajar del repositorio de GitHub en https://github.com/johnmyleswhite/ML_for_Hackers.

La idea es la siguiente, revisar el comportamiento de las votaciones del Senado y revisar el comportamiento del periodo 110 al 111 del congreso.

Lo que se hará es procesar la información, asignar ciertas variables categóricas que solo tendrán valores -1,0 y 1 para cada una de las características de los senadores analizados. Esto permite definir una métrica sencilla, como la del último ejemplo y basta considerar el productos de matrices para formar la matriz de similaridad.

#Análisis de votos 
library(foreign)
library(ggplot2)
#Tomamos la ruta de los archivos
data.dir 6, 0,no.pres[,i])
         no.pres[,i] 0 & no.pres[,i] < 4, 1, no.pres[,i])
         no.pres[,i] 1, -1, no.pres[,i])
 }
 
 return(as.matrix(no.pres[,10:ncol(no.pres)]))
}

#Aplicamos la función a cada entrada de rollcall.data
rollcall.simple

Con lo anterior pueden construir las gráficas que para apreciar la evolución del senado, donde la hipótesis inicial es que el congresos se polarizó en el 111 congreso.

Primero se revisa el comportamiento solo de las votaciones del 111 congreso y después se compara con respecto los congresos previos, partiendo desde el 101.

MDS_Votacion_Senado111_US

MDS_Votacion_Senado del 101 al 111_congreso_USEn resumen la última gráfica muestra que el Senado en distintos periodos se puede dividir en dos grupos, se puede decir que se polariza entre Republicanos y Demócratas. Los pocos senadores considerados como “independientes” parece, al revisar la última gráfica; que siempre votaran como los Demócratas.

La moraleja es que bajo una métrica sencilla se pueden extraer ciertas propiedades del comportamiento o dinámica del los datos. Sería interesante hacer un ejercicio similar para  el Senado de otros países.

Resumiendo, la técnica o análisis de agrupamientos (cluster) requiere como primer paso explorar nuestra información y procesarla si es necesario para poder definir una métrica que tenga coherencia con el tipo de información y datos.

Después se debe de calcular la matriz de distancias y por último aplicar escalamiento multidimensional o algún algoritmos para detectar aglomeraciones.

Como paso intermedio, tienen sentido plantearse una pregunta respecto a los posibles grupos que se busca encontrar en los datos y como siempre para cerrar el análisis, hacer usos de gráficos que permitan visualizar el resultado.

Ha esta entrada le falta mucho muchos ejemplos, que espero más adelante compartir otros. Respecto al tema falta mostrar varios detalles; tipos de métricas para tipos de datos, algoritmos de aglomeración para tipos de análisis y la creación de modelos basados en esta técnica.

Los aspectos teóricos se encontrarán en la categoría “Sobre Machine Learning”. Una disculpa para la gente de matemáticas por no haber puesto la definición formal de espacio métrico, métrica y mencionar algunas de sus propiedades, espero posteriormente hacer algunas entradas con un enfoque formal.

Actualización 19-05-2015

En la primera versión de esta entrada faltó mencionar muchos detalles sobre la técnica de K-medias, pero además no mencioné nada sobre otro algoritmo Hirarchical Clustering o Cluster jerárquico.

Comparto un ejemplo de esta técnica y en resumen lo que se hace en este algoritmo es construir un árbol que ayude aglomerar los datos en conjuntos los cuales al terminar el algoritmo muestra las aglomeraciones obtenidas.

El costo computacional es mucho mayor que los métodos k-medias, pero en general para muestras de datos relativamente grandes no resulta tan costoso implementarlo y además se han desarrollado algoritmos para reducir el costo en tiempo y procesamiento.

Los datos los tomo de la página del libro Elements of Stadistical Learning [2], en el código indico como descargarlo directamente.

Respecto al origen de los datos, cabe resaltar que corresponde al análisis de tumores cancerígenos, se analiza su valor genético a 64 tumores y por medio de las técnicas de Cluster se trata de detectar cuantos cánceres corresponde a qué tipo de tumor.

Algunas precisiones. Se descargaron los nombres de los tumores en un archivo txt, a que los nombres de las columnas no vienen por default en los datos. La matriz de datos está formada por 64 columnas y 6830 filas.

#Ejemplo de Cluster con datos de tumores 
#Ubicado en el directorio donde estan los nombres de los genes se procede a cargarlos en R
#Se carga como un vector
Nombres=read.table("nombres_tumores.txt")
#Visualizamos los datos
head(Nombres)
# V1
# CNS
# CNS
# CNS
# RENAL
# BREAST
# CNS

#Cargamos los datos de los genes

datos=scan("http://statweb.stanford.edu/~tibs/ElemStatLearn/datasets/nci.data")

#Visualizamos los datos
head(M)
#Tranformamos los datos en una matriz de 64 renglones por 6830 columnas
#Después invertimos la matrix
M=matrix(datos,nrow=64)
M2=t(M)

#Agremagamos los nombres a las columnas de M2

colnames(M2)=as.vector(t(Nombres))
 
head(M2)
#Así se visualiza la tabla de datos
# CNS CNS CNS RENAL BREAST CNS CNS BREAST NSCLC NSCLC RENAL RENAL RENAL RENAL RENAL RENAL RENAL BREAST NSCLC RENAL UNKNOWN
#[1,] 0.300 0.679961 0.940 2.800000e-01 0.485 0.310 -0.830 -0.190 0.460 0.760 0.270 -0.450 -0.030 0.710 -0.360 -0.210 -0.500 -1.060 0.150 -0.290 -0.200
#[2,] 1.180 1.289961 -0.040 -3.100000e-01 -0.465 -0.030 0.000 -0.870 0.000 1.490 0.630 -0.060 -1.120 0.000 -1.420 -1.950 -0.520 -2.190 -0.450 0.000 0.740
#[3,] 0.550 0.169961 -0.170 6.800000e-01 0.395 -0.100 0.130 -0.450 1.150 0.280 -0.360 0.150 -0.050 0.160 -0.030 -0.700 -0.660 -0.130 -0.320 0.050 0.080
#[4,] 1.140 0.379961 -0.040 -8.100000e-01 0.905 -0.460 -1.630 0.080 -1.400 0.100 -1.040 -0.610 0.000 -0.770 -2.280 -1.650 -2.610 0.000 -1.610 0.730 0.760
#[5,] -0.265 0.464961 -0.605 6.250000e-01 0.200 -0.205 0.075 0.005 -0.005 -0.525 0.015 -0.395 -0.285 0.045 0.135 -0.075 0.225 -0.485 -0.095 0.385 -0.105
#[6,] -0.070 0.579961 0.000 -1.387779e-17 -0.005 -0.540 -0.360 0.350 -0.700 0.360 -0.040 0.150 -0.250 -0.160 -0.320 0.060 -0.050 -0.430 -0.080 0.390 -0.080

#Aplicamos k-Means con 3 cluster a
KM3=kmeans(M2,3)

#Generamos una gráfica para visualizar los componentes

library(cluster)
clusplot(M2,KM3$clusteres,shape=TRUE,color=TRUE,labels=3)

#Aplicamos Hierarchical Clustering
H1<-hclust(dist(M),"average")
#Dendrograma
plot(H1,col=2,frame.plot=TRUE)

 Las gráficas que obtengo del código  son las  siguientes:

Dendrograma

Cluster_Genes

La primera gráfica es el modo visual en que se comportan los clusters en el método jerárquico y la segunda es una herramienta visual para observar el comportamiento del método k-medias.

Una observación entre los dos métodos es, k-medias pide como parámetro conocer cuantos clusters se buscamos en el algoritmo y el jerárquico no, de cierto modo la aglomeración por el jerárquico es cercana a lo que sería un método no supervisado.

Para el  algoritmo de k-medias resulta importante el análisis de la elección de la cantidad de clusters en los datos.

Existen varias técnicas para tratar de elegir el mejor valor de k, en R project, existen varias librerías entre las cuales están: “fpc” y “mclust”. Para los datos de este ejemplo probé con mclust, pero le resultaba muy costoso así que cancelé el proceso.

Pero otro modo muy sencillo y que no requiere ninguna librería es mediante el análisis de los cuadrados de las distancias de cada  dato aglomerado con respecto a la cantidad de clusters que se calculan.

#Elección de K para k-medias
wss <- (nrow(M2)-1)*sum(apply(M2,2,var))
 for (i in 2:10) wss[i] <- sum(kmeans(M2,centers=i)$withinss)
plot(1:10, wss, type="b", xlab="Númbero de Clusters",ylab="Suma de cuadrados dentro de cada cluster")

 La gráfica que se obtiene es la siguiente:

K-Medias

Esta última técnicas es sencilla, pero el problema es que no siempre indica de modo claro como elegir la cantidad de clusters de un modo óptimo. Por lo cual es bueno recurrir a otro tipo de algoritmo y compara nuestra selección y al final validar con alguna muestra de datos.

En un trabajo de investigación de Robert Tibshirani-G. Walther y Trevor Hastie[5] plantean un método para obtener un óptimo del número de clusters en el algoritmo k-meadias, sobre este reporte de investigación se basan varias librerías en R project para estimar el valor óptimo de k.

Espero que este complemente de una idea en general del tipo de análisis y de las técnicas básicas de como usar este tipo algoritmo. Cabe mencionar que las técnicas de Cluster están relacionadas con la Mezcla de Modelos y el algoritmo EM, lo cual se puede usar para la construcción de modelos. En otra entrada explico el uso y relación de esto.

Referencias:

1.-http://shop.oreilly.com/product/0636920018483.do

2.-http://statweb.stanford.edu/~tibs/ElemStatLearn/

3.-http://www.amazon.com/Cluster-Analysis-Brian-S-Everitt/dp/0470749911

4.-http://www.amazon.com/Machine-Learning-Probabilistic-Perspective-Computation/dp/0262018020

5.-http://web.stanford.edu/~hastie/Papers/gap.pdf