Modelización: modelos de clasificación
Joan Claverol _ Curso Juan Gabriel Gomila _ Udemy
23 de agosto de 2018
- Generando matrices de confusión
- Diagramas mosaicos
- Análisis de componentes principales
- Diagramas ROC
- Los árboles de decisión
- Bosques aleatorios o “random forest”
- Máquinas de soporte vectorial
- Naïve Bayes
- K Nearest Neighbors
- Redes neuronales para clasificar
- Análsisi del discriminante lineal
- La regresión logística
- Sentimental analysis
Generando matrices de confusión
Usaremos las notas del colegio para evaluar estos modelos
cp <- read.csv("../DataSets/college-perf.csv")
head(cp)Podemos en “Perf” como el valor real y “Pred” para ver la predicción que ha hecho un modelo.
Primero de todo vamos a dar valor a los tres niveles de la variable “Perf”:
cp$Perf <- ordered(cp$Perf,
levels = c("Low", "Medium", "High")) # convertira un valor de factores para darle un orden, con este orden damos a Low un valor mas bajo que High
cp$Pred <- ordered(cp$Pred,
levels = c("Low", "Medium", "High"))
str(cp)## 'data.frame': 4000 obs. of 7 variables:
## $ SAT : int 1380 1100 1110 1180 1240 1140 970 1100 1230 1280 ...
## $ GPA : num 2.53 3.18 2.73 2.49 2.89 2.85 2.37 2.67 3.01 3.08 ...
## $ Projects : int 1 1 2 3 3 2 1 2 2 0 ...
## $ Community: int 0 5 10 0 5 0 0 0 0 5 ...
## $ Income : int 41800 37600 34800 24100 56000 50800 47000 50900 37500 60200 ...
## $ Perf : Ord.factor w/ 3 levels "Low"<"Medium"<..: 1 1 2 1 2 1 2 2 2 1 ...
## $ Pred : Ord.factor w/ 3 levels "Low"<"Medium"<..: 1 1 2 3 2 1 2 2 2 1 ...Si nos fijamos ahora “cp” tiene las variables “Perf” y “Pred”, ordenadas por factores.
Vamos ahora a generar una matriz de confusión o tabla de doble entrada con los valores actuales y los predichos:
table <- table(cp$Perf, cp$Pred, dnn = c("Actual", "Predicho"))
table## Predicho
## Actual Low Medium High
## Low 1150 84 98
## Medium 166 1801 170
## High 35 38 458Con esta tabla de doble entrada podemos ver todos los cruzes posibles entre dos variables.
Podemos ver las proprciones o los porcentajes en vez de los valores absolutos en la tabla:
prop.table(table)## Predicho
## Actual Low Medium High
## Low 0.28750 0.02100 0.02450
## Medium 0.04150 0.45025 0.04250
## High 0.00875 0.00950 0.11450Pero para interpretar fila por fila o columna por columna, escalaremos los valores según lo que busquemos. Por ejemplo, normalizamos por fila:
round(prop.table(table, 1)*100, 2) # anadimos 1 para decir que escale por fila## Predicho
## Actual Low Medium High
## Low 86.34 6.31 7.36
## Medium 7.77 84.28 7.96
## High 6.59 7.16 86.25O lo podemos normalizar por columnas:
round(prop.table(table, 2)*100, 2) # lo hacemos por columna## Predicho
## Actual Low Medium High
## Low 85.12 4.37 13.50
## Medium 12.29 93.66 23.42
## High 2.59 1.98 63.09Diagramas mosaicos
Vamos a representar la performance del modelo de forma visual para detectar rápidamente donde estan los errores:
barplot(table, legend = TRUE, xlab = "Nota predecida por el modelo")
Este gráfico nos muestra las predicciones. Podemos ver claramente como el modelo predice de la mayoria de las notas “medias” correctamente, y tiene más errores en los extremos.
Ahora usaremos el mosaic plot para la eficiencia del modelo:
mosaicplot(table, main = "Eficiencia del modelo")
Aquí podemos ver la tabla en formato de gráfico. Mosaic plot toma el eje de las X de la tabla como los valores “Actuales”, y el eje de las Y como los valores del modelo predictivo.
Para comparar la performance del modelo para diferentes variables usamos la función summary sobre la tabla que hemos creado:
summary(table)## Number of cases in table: 4000
## Number of factors: 2
## Test for independence of all factors:
## Chisq = 4449, df = 4, p-value = 0El valor Chisq es muy alto, luego tenemos 4 grados de libertad, y un p-value. Un p-value bajo indica que las proporciones para las diferentes clases son significativamente diferentes y por lo tanto todas y cada una de ellas estan relacionadas, es decir, nos indica dependencia entre el valor actual y el predicho, por lo tanto el modelo es eficiente.
Análisis de componentes principales
Esta técnica de clasificación estadística para reducir la dimensionalidad, busca las variables dentro del dataset que son principales, y que definen la estructura el mismo. Busca reemplazar un numero grande de variables correlacionades entre si, para solo tomar un grupo de variables que no tengan relación alguna entre ellas.
usarrests <- read.csv("../DataSets/USArrests.csv", stringsAsFactors = F)
head(usarrests)Recomendación, en este caso, ya que cada estado está asociado a una observación, podemos convertir el numero de la columna en el nombre del estado:
rownames(usarrests) <- usarrests$X
usarrests$X <- NULL
head(usarrests)Ahora vamos a calcular para cada columna la variabilidad (la varianza):
apply(usarrests, 2, var)## Murder Assault UrbanPop Rape
## 18.97047 6945.16571 209.51878 87.72916Esta varianza nos muestra que los asaltos es la que varía mas de un estado a otro. La magnitud de las variables puede afectar mucho en el cálculo del acp, es por eso que debemos normalizar los resultados (escalar y centrar).
Usaremos el metodo prcomp, al que podemos pedir que nos centre (nos centre en la media) y que nos escale (que nos divida por la desviación típica) para perder así esta variabilidad desde una variable a la otra.
acp <- prcomp(usarrests,
center = TRUE,
scale = TRUE)
acp## Standard deviations (1, .., p=4):
## [1] 1.5748783 0.9948694 0.5971291 0.4164494
##
## Rotation (n x k) = (4 x 4):
## PC1 PC2 PC3 PC4
## Murder -0.5358995 0.4181809 -0.3412327 0.64922780
## Assault -0.5831836 0.1879856 -0.2681484 -0.74340748
## UrbanPop -0.2781909 -0.8728062 -0.3780158 0.13387773
## Rape -0.5434321 -0.1673186 0.8177779 0.08902432Podemos ver que nos ha devuelto la desviación standar, visualizemos los resultados:
plot(acp, type = "l")
Aquí estamos viendo las desviaciones típicas para cada uno de los componentes. Ahora tenemos que elegir hasta que “codo” queremos quedarnos. En este caso una buena opcion seria elegir la 1a y 2a analisis principal. Hagamos un resumen del acp para ver mas información:
summary(acp)## Importance of components:
## PC1 PC2 PC3 PC4
## Standard deviation 1.5749 0.9949 0.59713 0.41645
## Proportion of Variance 0.6201 0.2474 0.08914 0.04336
## Cumulative Proportion 0.6201 0.8675 0.95664 1.00000Con este resumen ver información muy relevante, como la desviacion estandar, la proporción de la varianza y la proporción acumulada. Ejemplo, para la primera componente principal, solo se explicará el 62% de los datos. De esta forma podemos elegir el numero de componentes que nos interesa.
Hay diferentes puntos que nos debemos fijar en la matriz de rotacion de los componentes principales:
Para la primera componente principal, que representa el 62% de los datos, las variables de “Muerder”, “Assault” y “Rape” tienen un orden de magnitud muy parecido, indicando que PC1 pondera por igual estas 3 variables, por lo tanto estan muy correlacionadas las unas con las otras.
En la PC2, da un peso muy grande en la variable “UrbanPop”, de modo que las otras variables están menos relacionada con la población urbana.
acp## Standard deviations (1, .., p=4):
## [1] 1.5748783 0.9948694 0.5971291 0.4164494
##
## Rotation (n x k) = (4 x 4):
## PC1 PC2 PC3 PC4
## Murder -0.5358995 0.4181809 -0.3412327 0.64922780
## Assault -0.5831836 0.1879856 -0.2681484 -0.74340748
## UrbanPop -0.2781909 -0.8728062 -0.3780158 0.13387773
## Rape -0.5434321 -0.1673186 0.8177779 0.08902432Podemos ver estas observaciones que hemos comentado antes, con un diagrama:
biplot(acp, scale = 0)
Este diagrama nos muestra de forma resumida como la primera y segunda componente principal nos exlica toda la información necesaria. Nos muestra donde se colocarían cada uno de los estados analizados dentro del diagrama.
Explicación del gráfico al final del video:
https://www.udemy.com/r-data-science/learn/v4/t/lecture/8961386?start=0
Ahora vamos a calcular para cada estado su componente principal numero 1:
pc1 <- apply(acp$rotation[,1]*usarrests, 1, sum)
pc1## Alabama Alaska Arizona Arkansas California
## -109.70674 -200.93002 -198.67595 -154.13078 -141.66518
## Colorado Connecticut Delaware Florida Georgia
## -181.98644 -87.23615 -190.34537 -153.19304 -173.82313
## Hawaii Idaho Illinois Indiana Iowa
## -60.04080 -109.02150 -125.99912 -114.92576 -52.53493
## Kansas Kentucky Louisiana Maine Maryland
## -116.69053 -67.92231 -194.84993 -63.43161 -233.71833
## Massachusetts Michigan Minnesota Mississippi Missouri
## -93.89436 -207.86176 -65.68138 -193.67727 -99.69902
## Montana Nebraska Nevada New Hampshire New Jersey
## -102.55422 -81.94814 -224.43659 -49.63545 -152.84113
## New Mexico New York North Carolina North Dakota Ohio
## -192.57846 -216.11699 -129.31137 -54.54813 -98.67129
## Oklahoma Oregon Pennsylvania Rhode Island South Carolina
## -140.26444 -90.91451 -111.36422 -122.84293 -209.74003
## South Dakota Tennessee Texas Utah Vermont
## -53.63715 -158.88941 -147.16953 -128.55049 -34.79666
## Virginia Washington West Virginia Wisconsin Wyoming
## -137.72208 -113.16669 -76.95290 -54.51128 -134.92677Ahora la segunda componente:
pc2 <- apply(acp$rotation[,2]*usarrests, 1, sum)
pc2## Alabama Alaska Arizona Arkansas California
## -194.7112752 -40.5473153 59.0145557 29.5446483 -234.5123542
## Colorado Connecticut Delaware Florida Georgia
## -24.4602694 -19.4446345 34.6766247 -280.3381377 -25.0709671
## Hawaii Idaho Illinois Indiana Iowa
## -49.3852136 15.7574333 -199.2179979 -8.8481365 -23.5265544
## Kansas Kentucky Louisiana Maine Maryland
## 12.9551436 -83.5608566 -30.0747857 -8.3751863 48.5206319
## Massachusetts Michigan Minnesota Mississippi Missouri
## -106.8894943 -32.3535756 -23.6218683 41.8469753 -146.9363486
## Montana Nebraska Nevada New Hampshire New Jersey
## -9.8906062 -8.3126180 40.4256132 -33.7453036 -11.6274411
## New Mexico New York North Carolina North Dakota Ohio
## 61.5587283 36.4081282 -283.9333766 -0.3290096 -12.7011753
## Oklahoma Oregon Pennsylvania Rhode Island South Carolina
## 19.6635681 -134.2822019 -5.5095502 -2.6673081 46.2369768
## South Dakota Tennessee Texas Utah Vermont
## -65.8260256 -22.3842141 13.8088319 12.9422324 -37.3683410
## Virginia Washington West Virginia Wisconsin Wyoming
## -16.1242322 4.3864873 9.4959571 -26.9978277 -16.9910326Añadimos las variables que hemos creado al dataset usarrest:
usarrests$pc1 <- pc1
usarrests$pc2 <- pc2
head(usarrests)Ahora podría eliminar las 4 primeras columnas y quedarme con las 2 componentes principales que explican el 86% de los datos:
usarrests[,1:4] <- NULL
head(usarrests)Diagramas ROC
Las curvas ROC nos muestra la fiabilidad en nuestro modelo. Nos dirá el porcentaje de éxito de una observación aplicando un determinado algoritmo, si será cierto o fracasará el resultado. Los éxitos se miden por la sensibilidad, o también llamada tasa de verdaderos positivos. El fracaso de falsos positivos, o también 1 - Especifidad. Los dos conceptos aportant la probabilidad de que la observación acierte (sensibilidad) o que falle(1-Especifidad).
ROC curve
Los diagramas ROC, nos permiten usarlos como técnica para saber las probabilidades que para cierto caso pertenezca a una o varias clases y ser nosotros los que determinemos la probabilidad de corte. Este diagrama nos permite visualizar el numero de verdaderos y falsos positivos de los diferentes niveles de corte en los que obtenemos.
library(ROCR)## Loading required package: gplots##
## Attaching package: 'gplots'## The following object is masked from 'package:stats':
##
## lowessdata1 <- read.csv("../DataSets/roc-example-1.csv")
data2 <- read.csv("../DataSets/roc-example-2.csv")
head(data1)head(data2)Tenemos dos datasets, uno con variables categóricas y el otro con variables numéricas.
Vamos el generar el objeto predicción a partir de las probabilidades y las clases, indicaremos que nos prediga la probabilidad en función de la clase:
pred1 <- prediction(data1$prob, data1$class)Ahora vamos a medir la performance de la predicción:
“tpr”: true positive rate
“fpr”: false positive rate
perf1 <- performance(pred1, "tpr", "fpr")
plot(perf1)
lines(par()$usr[1:2], par()$usr[3:4]) # anadimos la diagonal para ver el umbral de la curva ROC
La linea central indica un caso totalmente aleatorio.
La pregunta debería ser, a partir de que valores yo estoy seguro de que siempre me podría dar verdaderos positivos, por lo tanto debemos decidir donde cortamos. Para eso deberemos crear un objeto adicional, la probabilidad de los cortes.
prob.cuts.1 <- data.frame(cut = perf1@alpha.values[[1]], # @ debido a que es un objeto de la estructura interna de R
fpr = perf1@x.values[[1]],
tpr = perf1@y.values[[1]])
# analizaremos los primeros y los ultimos valores de la tabla
head(prob.cuts.1)tail(prob.cuts.1)Esta tabla nos permitirá encontrar cual es la mejor tasa de falsos positivos respecto la tasa de verdadores positivos:
head(prob.cuts.1[prob.cuts.1$tpr >= 0.8,])Podemos ver que el primer valor, en corte en “0.49”, acertaremos en el 81% de los casos y fallaremos en el 21% de lo casos. Podría ser un buen indicador que el modelo funciona correctamente. Es decir, que la probabilidad de exito en nuestro caso deberia ser por encima del 49%. Es decir que si cogemos una observación y el algoritmos nos dice que tiene un (ejemplo) 53% de exito, entonces es con un 80 de probabilidades un éxito.
Vamos a verlo con el ejemplo de etiquetas:
pred2 <- prediction(data2$prob, data2$class, label.ordering = c("non-buyer","buyer"))
perf2 <- performance(pred2, "tpr", "fpr")
plot(perf2)
lines(par()$usr[1:2], par()$usr[3:4])
Los árboles de decisión
Vamos a usar las siguientes librerias:
library(caret)## Loading required package: lattice## Loading required package: ggplot2library(rpart)
library(rpart.plot)
banknote <- read.csv("../DataSets/banknote-authentication.csv")
head(banknote)Para crear un modelo debemos dividir los datos entre los datos de entrenamiento y los datos de test.
set.seed(2018) # definimos una semillas para que no varie el resultado cada vez que ejecutemos este codigo
training.ids <- createDataPartition(banknote$class,
p = 0.7,
list = F)Ahora vamos a crear un árbol aleatorio:
mod <- rpart(class ~., data = banknote[training.ids,], # defino la variable dependiente y el dataset
method = "class", # especifico que es un metodo de classificacion
control = rpart.control(minsplit = 20, cp = 0.01)) # defino un control del modelo, solo debe conisderar nodos con 20 casos en su interior, y un parametro para ajustar si nivel de complejidad
mod## n= 961
##
## node), split, n, loss, yval, (yprob)
## * denotes terminal node
##
## 1) root 961 427 0 (0.555671176 0.444328824)
## 2) variance>=0.297405 500 53 0 (0.894000000 0.106000000)
## 4) curtosis>=-4.3856 474 33 0 (0.930379747 0.069620253)
## 8) variance>=1.5968 347 2 0 (0.994236311 0.005763689) *
## 9) variance< 1.5968 127 31 0 (0.755905512 0.244094488)
## 18) curtosis>=-1.7388 97 8 0 (0.917525773 0.082474227) *
## 19) curtosis< -1.7388 30 7 1 (0.233333333 0.766666667)
## 38) skew>=5.28905 7 0 0 (1.000000000 0.000000000) *
## 39) skew< 5.28905 23 0 1 (0.000000000 1.000000000) *
## 5) curtosis< -4.3856 26 6 1 (0.230769231 0.769230769)
## 10) variance>=2.293 7 1 0 (0.857142857 0.142857143) *
## 11) variance< 2.293 19 0 1 (0.000000000 1.000000000) *
## 3) variance< 0.297405 461 87 1 (0.188720174 0.811279826)
## 6) skew>=7.5653 69 12 0 (0.826086957 0.173913043)
## 12) variance>=-4.8358 57 0 0 (1.000000000 0.000000000) *
## 13) variance< -4.8358 12 0 1 (0.000000000 1.000000000) *
## 7) skew< 7.5653 392 30 1 (0.076530612 0.923469388)
## 14) variance>=-0.458565 60 19 1 (0.316666667 0.683333333)
## 28) curtosis>=0.297461 22 5 0 (0.772727273 0.227272727) *
## 29) curtosis< 0.297461 38 2 1 (0.052631579 0.947368421) *
## 15) variance< -0.458565 332 11 1 (0.033132530 0.966867470) *Ahora vamos a visualizar la información de forma grafica con un diagrama de árbol:
prp(mod, type = 2, extra = 104, nn = TRUE, # type = 2 nos permite ver el pocentaje de datos por nuivel
fallen.leaves = TRUE, faclen = 4, varlen = 8,
shadow.col = "gray", roundint=FALSE)
diferencia entre nodos hoja (finales) y nodos intermediorios.
La poda de árboles de decisión
Esta técnica nos permite reducir el tamaño del árbol. Vamos a ver como se dividen las componentes principales dentro de un árbol de decisión:
mod$cptable## CP nsplit rel error xerror xstd
## 1 0.67213115 0 1.00000000 1.0000000 0.03607406
## 2 0.10538642 1 0.32786885 0.3489461 0.02627735
## 3 0.03278689 2 0.22248244 0.2435597 0.02255369
## 4 0.02810304 3 0.18969555 0.1990632 0.02061446
## 5 0.01873536 4 0.16159251 0.1756440 0.01947412
## 6 0.01639344 6 0.12412178 0.1592506 0.01861618
## 7 0.01405152 7 0.10772834 0.1334895 0.01714873
## 8 0.01170960 9 0.07962529 0.1194379 0.01627482
## 9 0.01000000 10 0.06791569 0.1194379 0.01627482En esta tabla podemos ver:
Las componentes principales (CP)
El error relativo
Número de particiones
Error total
Desviación estándard
Podemos elegir el mejor componente principal, eligiendo el que tenga un error más bajo. Ejemplo; en el caso de la primera componente principal, podríamos elegirlo dado que es el mas alto, sin embargo, nos debemos asegurar que:
El error total no debe ser mayor que el error minimo permitido (xerror del último componente principal) más la desciación típica
Para poder podar este árbol también debemos evitar que el modelo se ajuste demasiado a la información con la que le hemos entrenado, evitando así el “overfitting”. Para evitarlo, usaremos la técnica de validación cruzada o “cross-validation”, y miraremos la complejidad generada a través de la validación cruzada, que justamene eslo que nos indica “xerror”.
- Supongamos que queremos elegor la primera componente:
Estaria esperando un error no suporior a 0.1194379 + 0.03607406, es decir que estaría esperando un error de validación cruzada no superior a 0.155512, sin embargo tiene un error de 1.
Si vamos a la octava componente, 0.1194379 + 0.01627482, es decir 0.1357127, sí que se adecua al que estamos buscando.
Miremos que pasa cuando recortamos hasta la 3 componente principal:
mod.pruned <- prune(mod, mod$cptable[3, "CP"])
prp(mod.pruned, type = 2, extra = 104, nn = TRUE,
fallen.leaves = TRUE, faclen = 4, varlen = 8,
shadow.col = "gray", roundint=FALSE)
Y ahora hasta la 8 componente principal:
mod.pruned <- prune(mod, mod$cptable[8, "CP"])
prp(mod.pruned, type = 2, extra = 104, nn = TRUE,
fallen.leaves = TRUE, faclen = 4, varlen = 8,
shadow.col = "gray", roundint=FALSE)
Ahora podemos usar la función predict() para predecir el resto de la información que no está introducido:
pred.pruned <- predict(mod, banknote[-training.ids,],
type = "class")
table(banknote[-training.ids,]$class, pred.pruned, dnn = c("Actual", "Predicho"))## Predicho
## Actual 0 1
## 0 223 5
## 1 11 172Ahora vamos a ver como funciona el modelo recortado (pruned):
pred.pruned <- predict(mod.pruned, banknote[-training.ids,],
type = "class")
table(banknote[-training.ids,]$class, pred.pruned, dnn = c("Actual", "Predicho"))## Predicho
## Actual 0 1
## 0 219 9
## 1 11 172Es muy similar al modelo inicial, pero siendo un árbol menos pesado con el que trabajar.
Podemos hacerlo mirando la probabilidad:
pred.pruned2 <- predict(mod.pruned, banknote[-training.ids,],
type = "prob")
head(pred.pruned2)## 0 1
## 1 0.9942363 0.005763689
## 11 1.0000000 0.000000000
## 14 1.0000000 0.000000000
## 16 0.9942363 0.005763689
## 19 0.9175258 0.082474227
## 27 0.9175258 0.082474227Ahora tenemos la probabilidad para cada una de las observaciones. Ahora podemos generar undiagrama ROC para saber si nuestro modelo es bueno o no:
# del paquete ROC usamos prediction
pred <- prediction(pred.pruned2[,2], banknote[-training.ids, "class"])
perf <- performance(pred, "tpr", "fpr")
plot(perf)
Podemos ver que el diagrama ROC es muy bueno, con una tasa de aciertos positivo smuy alta.
Bosques aleatorios o “random forest”
UNa muy buena técnica pero muy costosa computacionalmente. Usaremos las siguientes librerias:
library(randomForest)## randomForest 4.6-14## Type rfNews() to see new features/changes/bug fixes.##
## Attaching package: 'randomForest'## The following object is masked from 'package:ggplot2':
##
## marginlibrary(randomForest)Usaremos el mismo dataset que en el anterior tema “banknote”:
banknote$class <- factor(banknote$class)Para el random forest no necesitamos partir los datos, ya que el propio modelo lo hace de forma inherente. Aún así lo separaremos para hacer el tesing:
set.seed(2018)
training.ids <- createDataPartition(banknote$class,
p = 0.7,
list = F)Ahora usaremos la técnica del random forest, que se trata de construir varios árboles de clasificación (un bosque):
mod <- randomForest( x = banknote[training.ids, 1:4], # indicamos las variables independientes
y = banknote[training.ids, 5], # variables independientes
ntree = 500, # número de árboles a construir
keep.forest = TRUE)# si queremos retener todos los árboles de salida)Veamos como el modelo predice ahora los valores que hemos reservado para hacer el testing:
pred <- predict(mod, banknote[-training.ids,])
table(banknote[-training.ids, "class"], pred, dnn = c("Actual", "Predicho"))## Predicho
## Actual 0 1
## 0 223 5
## 1 3 180Este modelo es muy bueno, como podemos ver en la tabla.
El modelo de random forest en general no se queda con la información de los árboles que ha utilizado, y por lo tanto, no se puede usar de cara a futuro. Con el parámetro keep.forest = TRUE, nos aseguramos que guarde esta información para así usarla a futuro.
Vamos a generar un diagrama ROC:
probs <- predict(mod, banknote[-training.ids,], type = "prob")
head(probs)## 0 1
## 1 1.00 0.00
## 2 1.00 0.00
## 5 0.92 0.08
## 9 1.00 0.00
## 12 1.00 0.00
## 13 1.00 0.00pred <- prediction(probs[,2], banknote[-training.ids,"class"])
perf <- performance(pred, "tpr", "fpr")
plot(perf)
Nos sale un modelo mucho mejor que solo haciendo un árbol de decisión. Por lo tanto, aunque tenga un coste computacional muy alto, realmente sale a cuenta.
Máquinas de soporte vectorial
Este conjunto de algoritmos llamados “support vector machines” son útilies tanto en campos de clasificación como en regresión. SVM intenta elavorar un modelo que nos represente los puntos de muestra en el espacio, e intentar separar la muestras en los espcios de forma más amplia posible a través de un hiperplano de separación, definido como el vector entre los dos puntos de las dos clases que esten más cerca el uno del otro, también llamado el vector de soporte. El mejor vector es aquel que consiga un margen máximo entre los elementos a separar.
SVM
¿Cómo lo utlizamos directamente en R?
library(e1071)
library(caret)Usaremos el dataset “banknote” usado anteriormente para crear el modelo:
banknote <- read.csv("../DataSets/banknote-authentication.csv")
banknote$class <- factor(banknote$class)
set.seed(2018)
t.ids <- createDataPartition(banknote$class, p = 0.7, list = F)
mod <- svm(class ~., data = banknote[t.ids,])Comprobamos los resultados:
table(banknote[t.ids, "class"], fitted(mod), dnn = c("Actual", "Predicho")) # fitted: podemos elegir las varaibles que el modelo a predecido de forma automática, evitando que deba crear la predicción## Predicho
## Actual 0 1
## 0 534 0
## 1 0 427Nos da un resultado muy bueno respecto el conjunto de entrenamiento, ahora vamos a ver como es respecto el conjunto de testing:
pred <- predict(mod, banknote[-t.ids,])
table(banknote[-t.ids, "class"], pred, dnn = c("Actual", "Predicho"))## Predicho
## Actual 0 1
## 0 228 0
## 1 0 183Nos da también un resultado muy bueno para los resultados de validación.
Visualizamos un los gráficos del conjunto de entrenamiento a través de solo 2 de las 4 variables (no podemos ver un espacio de 4 dimensiones):
plot(mod, data = banknote[t.ids,], skew ~ variance)
Ahora hacemos la visualización respecto al conjunto de validación:
plot(mod, data = banknote[-t.ids,], skew ~ variance)
El modelo usa lo que llamamos un nucleo de convolución, también llamado kernel, para tranformar información que en principio no es lineal, en espacios de dimensión superior donde el espacio se puede separar más facilmente con una frontera o hiperplanos que maximicen la anchura, el margen entre los diferentes casos.
El SVM predice que tipo de modelo será, basado en la naturaleza de la variable a predecir. En el caso en que dicha variable es un factor, SVM constuye un modelo que clasfica.
Parámetros de un SVM
Por defect escala todas las variables para que tenga media 0 y desviación típca unitaria, antes de construir el modelo para que tengan mejor resultados.
El type() es el algoritmo que se utiliza para llevar a cabo la predicción, se puede poner algun tipu u otro en reación si es de regresión o de clasificación.
En vez de separar por un hiperplano, también se puede separar por un nucleo de convolución en el parámetro kernel().
Para balancear las clases, se puede añadir peso con el parámetro class.weights(). Ejemplo:
mod <- svm(class ~., data = banknote[t.ids,], class.weights = c("0" = 0.3, "1" = 0.7))Hay un parámetro llamado cost(), que nos permite que un coste bajo va a permitir menos clasificaciones por error, que un coste elevado:
mod <- svm(class ~., data = banknote[t.ids,], cost=1)Para ajustar los parametros la libreia svm trae una función para ajustar los parámetros lo mejor posible:
tuned <- tune.svm(class ~., data = banknote[t.ids,],
gamma = 10^(-6:-1), cost = 10^(1:2))
summary(tuned)##
## Parameter tuning of 'svm':
##
## - sampling method: 10-fold cross validation
##
## - best parameters:
## gamma cost
## 0.1 10
##
## - best performance: 0
##
## - Detailed performance results:
## gamma cost error dispersion
## 1 1e-06 10 0.44429768 0.04595582
## 2 1e-05 10 0.44429768 0.04595582
## 3 1e-04 10 0.06972723 0.01635636
## 4 1e-03 10 0.01977019 0.01509369
## 5 1e-02 10 0.01457259 0.01223107
## 6 1e-01 10 0.00000000 0.00000000
## 7 1e-06 100 0.44429768 0.04595582
## 8 1e-05 100 0.06972723 0.01635636
## 9 1e-04 100 0.01977019 0.01509369
## 10 1e-03 100 0.01457259 0.01223107
## 11 1e-02 100 0.00000000 0.00000000
## 12 1e-01 100 0.00000000 0.00000000Deberíamos quedarnos con un coste y una gamma menor posible que minimize el error, y obtener los mejores resultados posibles.
Naïve Bayes
Se trata de un clasificado probabilistico, asume que la presencia o la ausencia de una característica particular nunca está relacionada con la presencia o la auséncia de cualquier otra característica en la variable elegida:
Naïve Bayes
Requiereo que todas las variables sean categórica para poder trabajar.
Vamos a usarlo en R:
library(e1071)
library(caret)Trabajaremos con un nuevo dataset:
ep <- read.csv("../DataSets/electronics-purchase.csv")
head(ep)Vamos a crear el modelo:
set.seed(2018)
t.ids <- createDataPartition(ep$Purchase, p = 0.67, list = F)
mod <- naiveBayes(Purchase ~., data = ep[t.ids,])Podemos ver que el modelo ha echo una predicción a través de la discretización de las variables. Vamos a probar el modelo con los datos de validación:
pred <- predict(mod, ep[-t.ids,])
tab <- table(ep[-t.ids,]$Purchase, pred, dnn = c("Actual", "Predicha"))
confusionMatrix(tab)## Confusion Matrix and Statistics
##
## Predicha
## Actual No Yes
## No 1 1
## Yes 2 0
##
## Accuracy : 0.25
## 95% CI : (0.0063, 0.8059)
## No Information Rate : 0.75
## P-Value [Acc > NIR] : 0.9961
##
## Kappa : -0.5
## Mcnemar's Test P-Value : 1.0000
##
## Sensitivity : 0.3333
## Specificity : 0.0000
## Pos Pred Value : 0.5000
## Neg Pred Value : 0.0000
## Prevalence : 0.7500
## Detection Rate : 0.2500
## Detection Prevalence : 0.5000
## Balanced Accuracy : 0.1667
##
## 'Positive' Class : No
## K Nearest Neighbors
Este algoritmo no paramétrico que sirve tanto como regresión como clasificación. En ambos casos se basa en elegir dentro de un conjunto de K elementos los que se encuentran para catalogar o para clasificar una regresión. En el caso de clasificación, el objeto se va a clasificar por mayoría de sus vecinos. La k indica cuantos vecinos existen dentro de la clasificación. Para la regresión cada valor cercano se le asignará un valor, y el valor predicho será el promedio de todos los valores cercanos.
KNN
En esta imágen podemos ver que si tenemos k igual a 3, nos sale que el triangulo rojo sería la opción para la redonda verde, pero con k igual a 5 saldría un cuadrado azul. Por lo tanto es muy importante saber elegir la k.
library(class)
library(caret)Usaremos un dataset llamado vacation trip.
vac <- read.csv("../DataSets/vacation-trip-classification.csv")Importante el KNN requiere que las varaibles independientes o predictoras sean numéricas, y que las variables dependientes sea una categoría, por lo tanto deberé estandarizar las varaibles para que sean del tipo que queremos. En este caso no hace falta.
Para evitar que los grandes números dominen sobre los pequeños, vámos a estandarizar los datos:
vac$Income.z <- scale(vac$Income)
vac$Family_size.z <- scale(vac$Family_size)
head(vac)Ahora vamos a hacer la partición de los datos para hacer el modelo, però haremos 3 particiones:
set.seed(2018)
t.ids <- createDataPartition(vac$Result, p = 0.5, list = F)
train <- vac[t.ids,]
temp <- vac[-t.ids,] # conjunto temporal que ahora partiremos
v.ids <- createDataPartition(temp$Result, p = 0.5, list = F)
val <- temp[v.ids,]
test <- temp[-v.ids,]Vamos a generar predicciones para k igual a 1, solo el vecino más cercano:
pred1 <- knn(train[,4:5], # del conjunto de entrenamiento
val[,4:5], # quiero predecir los del conjunto de validación
train[,3], # tomo como categoria de la variables elegida
k = 1) # 1 solo vecino
errmat1 <- table(val$Result, pred1, dnn = c("Actual","Predichos"))
errmat1## Predichos
## Actual Buyer Non-buyer
## Buyer 3 2
## Non-buyer 2 3Vamos a incrementar el número de vecinos a ver si tenemos una mejora:
pred1 <- knn(train[,4:5], val[,4:5], train[,3], k = 5)
errmat1 <- table(val$Result, pred1, dnn = c("Actual","Predichos"))
errmat1## Predichos
## Actual Buyer Non-buyer
## Buyer 3 2
## Non-buyer 1 4También podemos usar el valor de k para generar predicciones y matriz de confusión utilizando la partición de testing junto con la de entrenamiento, en lugar de la de validacion.
pred2 <- knn(train[,4:5], test[,4:5], train[,3], k = 1)
errmat2 <- table(test$Result, pred2, dnn = c("Actual","Predichos"))
errmat2## Predichos
## Actual Buyer Non-buyer
## Buyer 3 2
## Non-buyer 1 3Eligiendo el mejor número de vecinos para la decisión
Crearemos una función para elegir la mejor k:
knn.automate <- function(tr_predictors, val_predictors, tr_target,
val_target, start_k, end_k){
for(k in start_k:end_k){
pred <- knn(tr_predictors, val_predictors, tr_target, k)
tab <- table(val_target, pred, dnn = c("Actual", "Predichos"))
cat(paste("Matriz de confusión para k = ",k,"\n"))
cat("====================================\n")
print(tab)
cat("------------------------------------\n")
}
}
knn.automate(train[,4:5], val[4:5], train[,3], val[,3], 1, 8)## Matriz de confusión para k = 1
## ====================================
## Predichos
## Actual Buyer Non-buyer
## Buyer 3 2
## Non-buyer 2 3
## ------------------------------------
## Matriz de confusión para k = 2
## ====================================
## Predichos
## Actual Buyer Non-buyer
## Buyer 3 2
## Non-buyer 2 3
## ------------------------------------
## Matriz de confusión para k = 3
## ====================================
## Predichos
## Actual Buyer Non-buyer
## Buyer 3 2
## Non-buyer 2 3
## ------------------------------------
## Matriz de confusión para k = 4
## ====================================
## Predichos
## Actual Buyer Non-buyer
## Buyer 3 2
## Non-buyer 0 5
## ------------------------------------
## Matriz de confusión para k = 5
## ====================================
## Predichos
## Actual Buyer Non-buyer
## Buyer 3 2
## Non-buyer 1 4
## ------------------------------------
## Matriz de confusión para k = 6
## ====================================
## Predichos
## Actual Buyer Non-buyer
## Buyer 3 2
## Non-buyer 1 4
## ------------------------------------
## Matriz de confusión para k = 7
## ====================================
## Predichos
## Actual Buyer Non-buyer
## Buyer 3 2
## Non-buyer 0 5
## ------------------------------------
## Matriz de confusión para k = 8
## ====================================
## Predichos
## Actual Buyer Non-buyer
## Buyer 2 3
## Non-buyer 0 5
## ------------------------------------Con 6 vecinos catalogamos correctamente el 80% de los casos, por lo tanto parece ser una buena opción para tomar.
Podemos también determinar el número de vecinos de forma automatica usando la siguiente función del paquete caret:
trcontrol <- trainControl(method = "repeatedcv",
number = 10, # numero de iteraciones de remuestreo
repeats = 3) # conjunto de repeticiones de la validaci´n cruzada
caret_knn_fit <- train(Result ~ Family_size + Income,
data = train,
method = "knn", # tipo de operación que queremos que haga
trControl = trcontrol, # tecnicas de muestreo para el entreno
preProcess = c("center", "scale"),
tuneLength = 10)
caret_knn_fit## k-Nearest Neighbors
##
## 21 samples
## 2 predictor
## 2 classes: 'Buyer', 'Non-buyer'
##
## Pre-processing: centered (2), scaled (2)
## Resampling: Cross-Validated (10 fold, repeated 3 times)
## Summary of sample sizes: 19, 19, 18, 19, 19, 19, ...
## Resampling results across tuning parameters:
##
## k Accuracy Kappa
## 5 0.7666667 0.53333333
## 7 0.6888889 0.38000000
## 9 0.8222222 0.64666667
## 11 0.7888889 0.58000000
## 13 0.8666667 0.72666667
## 15 0.8388889 0.68000000
## 17 0.8166667 0.64666667
## 19 0.4722222 -0.04333333
## 21 0.4944444 0.00000000
## 23 0.5166667 0.03000000
##
## Accuracy was used to select the optimal model using the largest value.
## The final value used for the model was k = 13.Nos muestra en consola el valor de las muestras y cual es la mejor k.
También podemos usar el knn para calcular las probabilidades de una clase a través de pesos: (ejemplo, a veces es más dificil determinar que es un comprador o no, por lo tanto ponderaremos un si lo es):
pred5 <- knn(train[,4:5], val[,4:5], train[,3], k = 5, prob = T)
pred5## [1] Buyer Non-buyer Buyer Non-buyer Non-buyer Buyer Buyer
## [8] Non-buyer Non-buyer Non-buyer
## attr(,"prob")
## [1] 0.8 0.6 0.8 1.0 1.0 0.8 0.6 1.0 1.0 0.6
## Levels: Buyer Non-buyerAquí tenemos las probabilidades respecto si es una categoría u otra.
Redes neuronales para clasificar
Usaremos el paquet nnet() para redes neuronales:
library(nnet)
library(caret)Una red neuronal requiere que las variables predictoras sean numéricas, ya las dependientes sean variables de salida 0 o 1. Sin embargo, la función neural network se encarga de hacer todo el trabajo de generar variables dummies para corregir los errores que tengan en el dataframe de salida, y dar la variable categorica como valor final.
Usaremos el csv de banknote:
bn <- read.csv("../DataSets/banknote-authentication.csv")
bn$class <- factor(bn$class) # en principio no haría falta però por si acaso lo podemos convertirCreamos una partición para hacer el modelo:
set.seed(2018)
t.id <- createDataPartition(bn$class, p = 0.7, list = F)Constriumos la red neuronal
mod <- nnet(class ~ ., data = bn[t.id,],
size = 3, # argumento para especificar el número de unidades de la capa más interna que la propia red neuronal ha de utilizar. Otra forma de calular este valor es poner un valor cercano al promedio al número de unidades en la entrada y en la salida. En este caso tenemos 4 variables de entrada, así que podríamos haber indicado un sizze 3 o 4. Un valor muy superior puede dar muy buen resultado, pero el coste computacional será mayor.
maxit = 10000, # número de iteraciones máximas que debe tratar de llevar a cabo para intentar converger. El agoritmo se parará si converge antes, y sin no ha encontrado convergencia después del número máximo indicado también se parará.
decay = .001, # controlar el overfitting, evitar que se ajuste demasiado a nuestros datos de entrenamiento.
rang = 0.05, # especifica el rango de pesos aleatorios iniciales que hay que asignar a la red neuronal. Si por defecto los parámetros de entrada tienen valores muy grandes, intentar que el rango se ajuste de modo que el rango por el máximo de las variables de dentro del dataset sera lo más cercano a uno possible.
na.action = na.omit, # los NAs puede hacer que el modelo no performe nada bien. Es por eso que con este parámetro poddemos omitir estos valores.
skip = T) # añade una capa adicional para separar los ndos de entrada de los nodos de salida.## # weights: 23
## initial value 606.896634
## iter 10 value 44.133919
## iter 20 value 3.286036
## iter 30 value 0.874201
## iter 40 value 0.766866
## iter 50 value 0.605232
## iter 60 value 0.494422
## iter 70 value 0.411223
## iter 80 value 0.382613
## iter 90 value 0.375601
## iter 100 value 0.373759
## iter 110 value 0.372608
## iter 120 value 0.371611
## iter 130 value 0.371244
## iter 140 value 0.370991
## iter 150 value 0.370940
## iter 160 value 0.370921
## iter 170 value 0.370905
## iter 180 value 0.370900
## iter 190 value 0.370898
## iter 200 value 0.370898
## iter 210 value 0.370898
## iter 220 value 0.370897
## final value 0.370897
## convergedEl modelo ha empezado a hacer iteraciones, llegando a la más óptima. Vamos a ver su predicción sobre e conjunto de entrenamiento:
pred <- predict(mod, newdata = bn[-t.id,], type = "class")Vamos a hacer una confusion matrix:
table(bn[-t.id,]$class, pred,dnn = c("Actual", "Predichos") )## Predichos
## Actual 0 1
## 0 228 0
## 1 0 183El algoritmo suele ir muy bien para casos de segmentación, reconocimiento de imágenes o de classificación.
Atención si hay NAs por medio de los datos, el algoritmo puede no performar nada bien.
Calculando cuál sería el máximo del parámetro rang dentro del modelo de la red neuronal:
#rang * max(|variables|) ~ 1
apply(bn, 2, max)## variance skew curtosis entropy class
## " 6.8248000" " 12.95160000" "17.92740000" " 2.4495000" "1"Viendo estos resultados vemos que la variables más grande (la curtosis) está alrededor de vente. Por lo tanto, deberíamos ve qué número multiplicado por 20 se acerca a 1. En este caso sería 0.05, que es el valor que hemos usado para el rang dentro del modelo de redes neuronales.
Vamos a ver la performance de la red neuronal usando ROC curve:
library(ROCR)
pred2 <- predict(mod, newdata = bn[-t.id,], type = "raw") # con el type "raw" no da el valor básico de probabiliades
perf <- performance(prediction(pred2, bn[-t.id,"class"]),
"tpr", "fpr")
plot(perf)
Un modelo ROC perfecto. En este caso, núnca sacara falsos positivos, siempre verdaderos positivos.
Análsisi del discriminante lineal
Está muy relacionad con el análisis de la varianza y de la regresión, en el sentido que ambos intentan expresar las variables dependientes como combinación lineal del resto de características o medidas.
Sin embargo, la ANOVA, utiliza variables independientes que son categóricas y la variable dependiente continua.
En el caso del análisis del discriminante lineal, las variables independientes son continuas, y la dependiente es categórica.
La idea es analizar similitudes y diferencias con los diferentes objetos que queremos observar. Ver si las medias son parecidas, o las varianzas son parecidas entre ellas. La función lo mirará por nosotros. Usaremos el paquete MASS:
library(MASS)
library(caret)
bn<- read.csv("../DataSets/banknote-authentication.csv")
bn$class <- factor(bn$class)
set.seed(2018)
t.id <- createDataPartition(bn$class, p=0.7, list = F)Constriumos la función en base a “lda” del paquete MASS:
mod <- lda(bn[t.id,1:4], # parámetros de entrada
bn[t.id,5]) # parámetros de salida, variable categórica
# también se puede usar esta metadología
#mod <- lda(class ~., data = bn[t.id,])En este caso vamos a crear la predicción pero vamos a hacerlo creando una nueva columna con la predicción dentro del dataset original:
bn[t.id, "Pred"] <- predict(mod, # modelo
bn[t.id, 1:4])$class # solo nos quedaremos con la categoría de classe
postResample(bn[t.id, "class"], bn[t.id, "Pred"])## Accuracy Kappa
## 0.9781478 0.9559639table(bn[t.id, "class"], bn[t.id, "Pred"], dnn = c("Actual", "Predichos"))## Predichos
## Actual 0 1
## 0 513 21
## 1 0 427Vamos a ver la performance sobre nuestros datos de testing:
bn[-t.id, "Pred"] <- predict(mod, bn[-t.id, 1:4])$class
table(bn[-t.id, "class"], bn[-t.id, "Pred"], dnn = c("Actual", "Predichos"))## Predichos
## Actual 0 1
## 0 217 11
## 1 0 183La regresión logística
Es un tipo de análisis de la regresión que se utiliza para predecir el resultado de una variable categórica en función de variables independientes o predictoras. Es muy útil para modelar la probabilidad de un evento que ocurre en función de otros factores, y forma parte del conjunto de modelos lineales generalizados (glm).
En general, la regresión logística parte de poder trabajar con probabilidades conocidas de diferentes eventos y modelarlas con una transformación logarítmica:
En la gráfica superior, cada vez que tengamos un resultado por encima de 0.1, la categoría será 1, mientras que si está por debajo, la categoría será 0. Por lo tanto para este modelo necesitaremos que sean variables dicotómicas (0 i 1).
Atención restricciones para las variables de regresión:
Las variables de salida deben ser 0 o 1 y categóricas
Requierque todas las variables predictoras sean numéricas
Usaremos funciones própias del paquete stats (paquete glm):
library(caret)
bh <- read.csv("../DataSets/boston-housing-logistic.csv")
bh$CLASS <- factor(bh$CLASS, levels = c(0,1)) # modificamos la variable para que sea factor, con los niveles 0 i 1
set.seed(2018)
t.id <- createDataPartition(bh$CLASS, p=0.7, list = F)
mod <- glm(CLASS ~ ., data = bh[t.id, ],
family = binomial) # indicamos que el resultado de salida es 0 i 1
summary(mod)##
## Call:
## glm(formula = CLASS ~ ., family = binomial, data = bh[t.id, ])
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.47008 -0.44341 0.05835 0.37747 2.76863
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 27.978520 4.180802 6.692 2.20e-11 ***
## NOX -21.550773 4.297748 -5.014 5.32e-07 ***
## DIS -0.350211 0.150864 -2.321 0.0203 *
## RAD 0.183897 0.076947 2.390 0.0169 *
## TAX -0.004861 0.004096 -1.187 0.2353
## PTRATIO -0.889614 0.157880 -5.635 1.75e-08 ***
## B 0.006623 0.002938 2.254 0.0242 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 353.03 on 254 degrees of freedom
## Residual deviance: 154.22 on 248 degrees of freedom
## AIC: 168.22
##
## Number of Fisher Scoring iterations: 6Lo que busca esta función, es aplicar una función lineal, pero con forma logarítmica. Encima podemos ver los parámetros que ha calculado nuestro modelo.
Vámos a calcular las probabilidades de éxito (1) o de fallo (0). Primero lo haremos para las variables de entrenamiento:
bh[-t.id, "PROB_SUCCESS"] <- predict(mod, newdata = bh[-t.id,], type="response") # miramos que probabilidad de éxito tienen la variables de entrada
bh[-t.id, "PRED_50"] <- ifelse(bh[-t.id, "PROB_SUCCESS"]>=0.5, 1, 0) # vamos a dividir en casos de éxito o fracaso a partir de que esten por encima o por debajo de 0.5
# Analizamos los resultados con una matriz de confusión:
table(bh[-t.id,"CLASS"], bh[-t.id,"PRED_50"], dnn=c("Actual","Predicho"))## Predicho
## Actual 0 1
## 0 39 12
## 1 3 54Sentimental analysis
Conectandonos a la API de Twitter
Deberemos instalar varias librerias que no vienen de serie. Vamos a ver como lo podemos hacer:
# install.packages(c("twitteR", "RColorBrewer", "plyr",
# "ggplot2", "devtools", "httr"))
# require(devtools)
# install_url("https://cran.r-project.org/src/contrib/Archive/Rstem/Rstem_0.4-1.tar.gz")
# install_url("https://cran.r-project.org/src/contrib/Archive/slam/slam_0.1-37.tar.gz")
# install_url("https://cran.r-project.org/src/contrib/Archive/sentiment/sentiment_0.2.tar.gz")
library(slam)
library(sentiment)## Loading required package: tm## Loading required package: NLP##
## Attaching package: 'NLP'## The following object is masked from 'package:ggplot2':
##
## annotate## Loading required package: Rstemlibrary(twitteR)Ahora vamos a introducir las claves para acceder a la API de twitter:
- api_key
- api_secret
- api_token
- api_access_secret
Vamos a ver si hemos podido congurarlo todo correctamente:
setup_twitter_oauth(api_key, api_secret,
access_token, access_token_secret)## [1] "Using direct authentication"1## [1] 1Limpiando los tweets descargados
Descarga de tweets:
tweets <- searchTwitter("mallorca", # la palabra que buscaremos
n = 1500, # nos descargamos 1500 tweets
lang = "es") # especificamos el lenguajeVamos a crear una función para limpiar espacios inútiles (dígitos, espacios extras, puntuaciones, …).
texts <- sapply(tweets, function(x) x$getText()) # extraemos los textos aplicando la función x$getText()
head(texts, 2)## [1] "@grinwall31 Desde Mallorca <U+0001F4AA>! #irph https://t.co/5HAbx8ZLOS"
## [2] "RT @rfuertespuch: GH | Bea deja los estudios y vive una experiencia terrorífica en Mallorca https://t.co/YsnzHKMniu"Limpiamos los tweets:
clean.data <- function(text){
#eliminar re-tweets y @ del texto original
text = gsub("(RT|VIA)((?:\\b\\W*@\\w+)+)", "", text)
text = gsub("@\\w+", "", text)
#eliminar signos de puntuación y dígitos del 0 al 9
text = gsub("[[:punct:]]", "", text)
text = gsub("[[:digit:]]", "", text)
#eliminar links html, tabulaciones y espacios adicionales
text = gsub("http\\w+","",text)
text = gsub("[ \t]{2,}", "", text)
text = gsub("^\\s+|\\s+$", "", text)
}
texts <- clean.data(texts)
head(texts, 2)## [1] "Desde Mallorca <U+0001F4AA> irph"
## [2] "GHBea deja los estudios y vive una experiencia terrorífica en Mallorca"Ahora vamos a hacer otra función para gestionar NAs dentro del texto:
handle.error <- function(x){
#crear el valor omitido
y = NA
#try_catch error
try_error <- tryCatch(tolower(x), error=function(e) e)
#si no hay error
if(!inherits(try_error, "error"))
y = tolower(x)
#devolvemos el resultado
return(y)
}
texts = sapply(texts, handle.error)
head(texts, 2)## Desde Mallorca <U+0001F4AA> irph
## "desde mallorca <U+0001F4AA> irph<U+7378><U+5F73><U+6572><U+6F70><U+7472><U+2E00><U+6F63><U+2F6D><U+2F69><U+7378><U+5F73><U+6572><U+6F70><U+7472><U+6300><U+3D6B>0<U+6400><U+313D><U+3630><U+3934><U+3336><U+3330><U+3733><U+3932><U+3630><U+3934><U+3935><U+3600><U+3934><U+3336><U+3330><U+3733><U+3932><U+3630><U+3934><U+3935><U+6900><U+7474><U+7265><U+632E><U+6D6F>"
## GHBea deja los estudios y vive una experiencia terrorífica en Mallorca
## "ghbea deja los estudios y vive una experiencia terrorífica en mallorca"Eliminar los NAs y NULL la columna de los nombres:
texts <- texts[!is.na(texts)]
names(texts) <- NULLAnálisis de sentimiento
La función classify_emotion no ayuda a analizar el texto y a classificarlo en diferentes categorias de emociones (enfado, angusti, miedo,…)
# analisis de sentimiento
class_emo <- classify_emotion(texts, algorithm = "bayes", prior = 1)
head(class_emo)## ANGER DISGUST FEAR
## [1,] "1.46871776464786" "3.09234031207392" "2.06783599555953"
## [2,] "1.46871776464786" "3.09234031207392" "2.06783599555953"
## [3,] "1.46871776464786" "3.09234031207392" "2.06783599555953"
## [4,] "1.46871776464786" "3.09234031207392" "7.34083555412328"
## [5,] "1.46871776464786" "3.09234031207392" "7.34083555412328"
## [6,] "1.46871776464786" "3.09234031207392" "2.06783599555953"
## JOY SADNESS SURPRISE BEST_FIT
## [1,] "1.02547755260094" "1.7277074477352" "2.78695866252273" NA
## [2,] "1.02547755260094" "1.7277074477352" "2.78695866252273" NA
## [3,] "1.02547755260094" "1.7277074477352" "2.78695866252273" NA
## [4,] "1.02547755260094" "7.34083555412328" "2.78695866252273" "fear"
## [5,] "1.02547755260094" "7.34083555412328" "2.78695866252273" "fear"
## [6,] "1.02547755260094" "1.7277074477352" "2.78695866252273" NAEstos tweets nos muestran como podemos usar para predecir el sentimiento, y que valor le da a cada uno de los sentimientos que considera.
La última columna indica cuál es el valor que ha ganado por peso, siendo NA los tweets que no ha pedido predecir.
emotion <- class_emo[,7] # guardamos los valores resultantes como un nuevo vector de emociones
emotion[is.na(emotion)]<-"unknown" # le damos a todos los valores de NA el valor de "unknown"
head(emotion)## [1] "unknown" "unknown" "unknown" "fear" "fear" "unknown"Ahora haremos el análisis de la polaridad para ver como es classifcado un tweet (positivo, negativo, neutro).
#analisis de la polaridad
class_pol <- classify_polarity(texts, algorithm = "bayes") # usamos el algoritmo de naivebayes para predecir cuál es la probabilidad de que cada palabra signifique alegría y tristeza, en este caso que sea positiva, negativa o neutra.
polarity <- class_pol[,4] #extraemos la cuarta columna
head(class_pol)## POS NEG POS/NEG BEST_FIT
## [1,] "1.03127774142571" "0.445453222112551" "2.31512017476245" "positive"
## [2,] "1.03127774142571" "0.445453222112551" "2.31512017476245" "positive"
## [3,] "1.03127774142571" "0.445453222112551" "2.31512017476245" "positive"
## [4,] "1.03127774142571" "0.445453222112551" "2.31512017476245" "positive"
## [5,] "1.03127774142571" "0.445453222112551" "2.31512017476245" "positive"
## [6,] "1.03127774142571" "0.445453222112551" "2.31512017476245" "positive"Este algoritmo, si tiene dudas entre positivo y neutro, se queda en positivo. y lo mismo con negativo y neutro. Delante de positivo y neutro no sabe decidir.
Creamos un data frame que resuma los sentimientos (tanto sentiment como polaridad):
sent_df <- data.frame(text = texts,
emotion = emotion, polarity = polarity, stringsAsFactors = F)
sent_df <- within(sent_df, emotion <- factor(emotion,
levels = names(sort(table(emotion),
decreasing= T))))
library(RColorBrewer)
library(ggplot2)
ggplot(sent_df, aes(x=emotion))+
geom_bar(aes(y = ..count.., fill=emotion))+
scale_fill_brewer(palette = "Set2")+
labs(x="Categorías de emocion", y = "Número de Tweets")+
labs(title = "Análisis de Sentimiento acerca de Machine Learning")
ggplot(sent_df, aes(x=polarity))+
geom_bar(aes(y = ..count.., fill = polarity))+
scale_fill_brewer(palette = "Set3")+
labs(x="Categorías de polaridad", y = "Número de Tweets")+
labs(title="Análisis de Sentimiento acerca de Machine Learning")