El próximo 1 de Febrero a las 19:oo h tendremos nueva reunión del “Grupo de Usuarios de R de Almería”, en la que contaremos con dos charlas muy interesantes.
El lugar elegido para las presentaciones son las instalaciones El Cable – Andalucía Open Future C/ Magistral Domínguez, 23 Almería.
¡Os esperamos!
Posted por @mirallesjm
a 30 de Junio, 2017
Podría resulta interesante tanto a educadores, padres y alumnos conocer en que media afectan características demográficos, sociales y relacionadas con la propia Escuela a la hora de explicar el éxito o fracaso escolar. Esto es precisamente lo que me ha motivado a realizar este estudio analítico, para ello nos centraremos en la asignatura de Matemáticas.
¿En que medida afecta al suspenso en matemáticas factores como?
A través de nuestro estudio daremos respuesta a estas preguntas. Para ello utilizaremos el algoritmo XGBoost, una técnica de Machine Learning muy popular en cualquier competición de Data Science. Nos apoyaremos en el conjunto de datos del consumo de alcohol en estudiantes de la UCI Machine Learning Repositorio
¿Qué es XGBoost?
Técnica supervisada que nos permite realizar tanto regresión como clasificación. Su capacidad de hacer computación en paralelo hace que sea 10 veces más rápido que el gradient boosting.
Su alto poder predictivo la convierte en una técnica ideal para competiciones. Cuenta con funciones adicionales para realizar validación cruzada y búsqueda de variables importante, además existen muchos parámetros que nos permiten optimizar el modelo.
Exploración analítica de los datos
Vamos a leer en los datos y exploraremos su estructura:
student <- read.table(“student-mat.txt”, header = TRUE, sep = “\t”, quote = “\””, dec = “.”, na.strings = “NA”, fill = TRUE, comment.char = “”,stringsAsFactors = T)
str(student)
‘data.frame’: 395 obs. of 33 variables:
$ school : Factor w/ 2 levels “GP”,”MS”: 1 1 1 1 1 1 1 1 1 1 …
$ sex : Factor w/ 2 levels “F”,”M”: 1 1 1 1 1 2 2 1 2 2 …
$ age : int 18 17 15 15 16 16 16 17 15 15 …
$ address : Factor w/ 2 levels “R”,”U”: 2 2 2 2 2 2 2 2 2 2 …
$ famsize : Factor w/ 2 levels “GT3″,”LE3”: 1 1 2 1 1 2 2 1 2 1 …
$ Pstatus : Factor w/ 2 levels “A”,”T”: 1 2 2 2 2 2 2 1 1 2 …
$ Medu : int 4 1 1 4 3 4 2 4 3 3 …
$ Fedu : int 4 1 1 2 3 3 2 4 2 4 …
$ Mjob : Factor w/ 5 levels “at_home”,”health”,..: 1 1 1 2 3 4 3 3 4 3 …
$ Fjob : Factor w/ 5 levels “at_home”,”health”,..: 5 3 3 4 3 3 3 5 3 3 …
$ reason : Factor w/ 4 levels “course”,”home”,..: 1 1 3 2 2 4 2 2 2 2 …
$ guardian : Factor w/ 3 levels “father”,”mother”,..: 2 1 2 2 1 2 2 2 2 2 …
$ traveltime: int 2 1 1 1 1 1 1 2 1 1 …
$ studytime : int 2 2 2 3 2 2 2 2 2 2 …
$ failures : int 0 0 3 0 0 0 0 0 0 0 …
$ schoolsup : Factor w/ 2 levels “no”,”yes”: 2 1 2 1 1 1 1 2 1 1 …
$ famsup : Factor w/ 2 levels “no”,”yes”: 1 2 1 2 2 2 1 2 2 2 …
$ paid : Factor w/ 2 levels “no”,”yes”: 1 1 2 2 2 2 1 1 2 2 …
$ activities: Factor w/ 2 levels “no”,”yes”: 1 1 1 2 1 2 1 1 1 2 …
$ nursery : Factor w/ 2 levels “no”,”yes”: 2 1 2 2 2 2 2 2 2 2 …
$ higher : Factor w/ 2 levels “no”,”yes”: 2 2 2 2 2 2 2 2 2 2 …
$ internet : Factor w/ 2 levels “no”,”yes”: 1 2 2 2 1 2 2 1 2 2 …
$ romantic : Factor w/ 2 levels “no”,”yes”: 1 1 1 2 1 1 1 1 1 1 …
$ famrel : int 4 5 4 3 4 5 4 4 4 5 …
$ freetime : int 3 3 3 2 3 4 4 1 2 5 …
$ goout : int 4 3 2 2 2 2 4 4 2 1 …
$ Dalc : int 1 1 2 1 1 1 1 1 1 1 …
$ Walc : int 1 1 3 1 2 2 1 1 1 1 …
$ health : int 3 3 3 5 5 5 3 1 1 5 …
$ absences : int 6 4 10 2 4 10 0 6 0 0 …
$ G1 : int 5 5 7 15 6 15 12 6 16 14 …
$ G2 : int 6 5 8 14 10 15 12 5 18 15 …
$ G3 : int 6 6 10 15 10 15 11 6 19 15 …dim(student)
[1] 395 33
En nuestro caso contamos con un dataset que contiene 395 alumnos con 33 características medibles. Nuestro objetivo es clasificar a los alumnos en función de la variable respuesta G3 – nota final (numérica: de 0 a 20, target).
Veamos la distribución de la nota final, en este caso utilizaremos la libraria ggplot2.
# Distribucion Target
library(ggplot2)
ggplot(student, aes(x=G3)) +
geom_histogram(binwidth = 1, color=”white”, fill=rgb(0.2,0.7,0.1,0.4))
# discretizaremos la variable objetivo
student$Result <- as.factor(ifelse(student$G3 <= 9.9, ‘0’,
ifelse(student$G3 > 9.9 & student$G3 <=13.9, ‘1’,
ifelse(student$G3 > 13.9 & student$G3 <=16.9, ‘2’,
ifelse(student$G3 > 16.9, ‘3’, ‘0’)))))table(student$Result)
0 1 2 3
130 165 76 24
El siguiente paso es realizar un estudio bivariante de nuestra target (Result) frente a algunos de las principales variables de nuestro dataset. Con ello podremos observar como se comportan las distribuciones y hacernos una idea de aquellos posibles predictores.
library(GGally)
library(ggplot2)levels(student$Result) <- c(“Suspenso”,”Aprobado”,”Notable”,”Sobresaliente”)
ggpairs(student[,c(“absences”,”failures”,”health”,”goout”,”Medu”,”Result”)],
aes(col = Result))
library(scales)
library(gridExtra)plot1 <- ggplot(student, aes(x = Result, y = age, fill= Result)) +
scale_x_discrete(“Target”) +
scale_y_continuous(“Edad”, breaks = seq(0,22, 2)) +
geom_boxplot()plot2 <- ggplot(student, aes(x = Result, y = studytime , fill= Result)) +
scale_x_discrete(“Target”) +
scale_y_continuous(“Tiempo Estudio”, breaks = seq(0,4, 1)) +
geom_boxplot()plot3 <- ggplot(student, aes(x= sex, group = Result )) +
geom_bar(aes(y = ..prop.., fill = factor(..x..)), stat = “count”) +
geom_text(aes(label = scales::percent(..prop..), y = ..prop..), stat = “count”, vjust = -.5) +
labs(y = “Percent”, fill = “Sexo”) +
facet_grid(~Result) +
scale_y_continuous(labels = percent)plot4 <- ggplot(student, aes(x= school, group = Result )) +
geom_bar(aes(y = ..prop.., fill = factor(..x..)), stat = “count”) +
geom_text(aes(label = scales::percent(..prop..), y = ..prop..), stat = “count”, vjust = -.5) +
labs(y = “Percent”, fill = “School”) +
facet_grid(~Result) +
scale_y_continuous(labels = percent)grid.arrange(plot1,plot2,plot3,plot4,ncol=2,nrow=2)
Este análisis descriptivo podemos observar como las mujeres presentan mejores notas que los hombre, como afecta el absentismo, salud, salir fuera con los amigos, nivel de estudio de la madre al éxito/fracaso en el examen de matemáticas.
Construcción del Modelo Xgboost
Ya estamos preparados para construir nuestro modelo de clasificación supervisado. Necesitaremos cargar la librería xgboost.
library(xgboost)
# creamos nuestro train y test
set.seed(1976)
trainIndex <- createDataPartition(student$Result, p= .7, list = F, times = 1)train <- student[trainIndex,]
test <- student[trainIndex,]# utilizaremos la data.table para ser mas eficientes en memoria
library(data.table)
train = data.table(train)
test = data.table(test)# preparando los datos usuando XGBoost
library(Matrix)# necesitamos convertir los datos en una matriz dumificada
sparse.matrix.train <- sparse.model.matrix(Result ~ .-1, train)
sparse.matrix.test <- sparse.model.matrix(Result ~ .-1, test)# xgboost necesita que la variable objetivo sea numerica
target_train = as.numeric(train[,Result])
target_test = as.numeric(test[,Result])# convertimos al formato xgb.DMatrix que nos permite utilizar funciones avanzadas
dtrain <- xgb.DMatrix(data = sparse.matrix.train, label=target_train)
dtest <- xgb.DMatrix(data = sparse.matrix.test, label=target_test)
watchlist <- list(train = dtrain, test =dtest)# construimos nuestro machine learning model
# primero definimos los parámetros xgboost
param <- list(“objective” = “multi:softprob”, # multiclass classification
“num_class” = 4, # numero de clases
“eval_metric” = “merror”, # evaluacion de la metrica
“nthread” = 8, # numero de hilos a utilizar
“max_depth” = 16, # profundidad maxima del arbol
“eta” = 0.3, # tamaño de cada paso
“gamma” = 0, # reduccion minima de perdidas
“subsample” = 1, # Parte de las instancias de datos para crecer el arbol
“colsample_bytree” = 1, # Sub-muestra de columnas al construir cada árbol
“min_child_weight” = 12 # Cantidad mínima de peso de la instancia necesaria en un child
)
Una forma de medir el progreso en el aprendizaje de un modelo es proporcionar a XGBoost un segundo conjunto de datos ya clasificada. Por lo tanto, se puede aprender en el primer conjunto de datos y poner a prueba su modelo en la segunda. Algunas métricas se miden después de cada ronda durante el aprendizaje, para evitar este sobreaprendizaje utilizaremos el parámetro watchlist que nos permitirá observar el error de calculo en cada ronda tanto para el train como para el test.
# entrenamiento del modelo real con datos completos
xgb <- xgb.train(param=param,
data=dtrain,
nrounds=10,
verbose=1,
watchlist = watchlist)[1] train-merror:0.464029 test-merror:0.464029
[2] train-merror:0.413669 test-merror:0.413669
[3] train-merror:0.402878 test-merror:0.402878
[4] train-merror:0.381295 test-merror:0.381295
[5] train-merror:0.377698 test-merror:0.377698
[6] train-merror:0.352518 test-merror:0.352518
[7] train-merror:0.323741 test-merror:0.323741
[8] train-merror:0.323741 test-merror:0.323741
[9] train-merror:0.320144 test-merror:0.320144
[10] train-merror:0.316547 test-merror:0.316547
Hemos calculado para cada ronda la misma métrica de error promedio (nrounds a 10). En este caso, train-merror esta asociado a la muestra entrenamiento y test-merror al conjunto de prueba.
Ambos indicadores relacionados con los errores de entrenamiento y prueba son muy similares. Por lo tanto, el aprendizaje realizado en el train coincide con las observaciones del test.
Por último, observamos el gráfico que nos muestra de manera ordenada aquellas variables mas importantes y que finalmente contribuyen a explicar las notas finales en matemáticas.
# obtenemos el modelo
modelo= xgb.dump(xgb, with_stats=TRUE)
# obtenemos los nombres de las variables
names = dimnames(dtrain)[[2]]
# generamos la matriz de variables importantes
importance_matrix = xgb.importance(names, model=xgb)
# generamos el grafico de importancia
xgb.plot.importance(importance_matrix)
Atributos como el absentismo, fracasos en cursos anteriores, salida con los amigos y nivel de estudio de la madre, entre otras, hacen bastante explicativo el modelo que nos permite determinar la calificación el la asignatura de matemáticas.
Hasta aquí el final del post, espero que os haya parecido interesante.
Ya tenemos disponibles en nuestro canal en Vimeo los dos vídeos de nuestra reunión del pasado 12 de enero. Recordad que además las presentaciones de estas dos charlas están disponibles en Slideshare en www.slideshare.net/AlmeriaRUsers
Posted por @mirallesjm
a 19 de Enero, 2016
Motivación
La mayoría de las veces que pasamos por una estación de servicio, no podemos evitar girar la cabeza para revisar el panel de precios de la gasolina e intentar realizar una comparativa mental del resto de gasolineras conocidas. Establecer esta comparativa de precios resulta difícil principalmente porque no tenemos un único valor a comparar (Precio gasolina 95/98, gasóleo A/B, etc.) y dichos precios se mueven con un margen muy pequeño.
Nuestro objetivo pasa por visualizar esta información mediante un mapa que nos permita ver una primera distribución de precios y posteriormente aplicar técnicas de Machine Learning que nos permitan segmentar las estaciones de servicio vía precios y poder determinar grupos homogéneos de clasificación.
Análisis
Para realizar el estudio, disponemos de una base de datos de todas las estaciones de servicio de España (http://www.geoportalgasolineras.es/#/Descargas).
Lo primero que hacemos es descargarnos la información y leerla:
# Lectura del Dataset
Fuel <- read.table(“Gasolineras.txt”, header = TRUE, sep = “\t”, quote = “\””, dec = “,”,
na.strings = “NA”, fill = TRUE, comment.char = “”,stringsAsFactors = T)# Dimension
dim(Fuel)
El dataset esta formado por 9.951 estaciones de servicio y 25 variables.
Observamos la estructura del dataset
# Recodificamos el nombre de los campos “Precio”
FuelPrecio <- select(Fuel, contains(“Precio”))
names(FuelPrecio) <- c(“PGasolina95″,”PGasoleoA”,”PGasoleoB”,”PBioEtanol”,”PNuevoGasoleoA”,
“PBioDiesel”,”PGasolina98″,”PGasNaturalComp”,”PGasNaturalLic”,
“PGasLicuPetro”)
# Estructura del dataset
str(FuelPrecio)
‘data.frame’: 9951 obs. of 10 variables:
$ PGasolina95 : num NA 1.17 1.23 1.16 1.22 …
$ PGasoleoA : num 1.06 1.05 1.12 1.03 1.12 …
$ PGasoleoB : num 0.726 NA NA NA NA NA 0.785 NA NA 0.8 …
$ PBioEtanol : num NA NA NA NA NA …
$ PNuevoGasoleoA : num NA 1.13 1.18 NA 1.18 …
$ PBioDiesel : num NA NA NA NA NA …
$ PGasolina98 : num NA 1.25 1.34 1.26 1.33 …
$ PGasNaturalComp: num NA NA NA NA NA NA NA NA NA NA …
$ PGasNaturalLic : num NA NA NA NA NA NA NA NA NA NA …
$ PGasLicuPetro : num NA NA NA NA NA NA NA NA NA NA …summary(FuelPrecio)
PGasolina95 PGasoleoA PGasoleoB PBioEtanol
Min. :0.779 Min. :0.749 Min. :0.544 Min. :1.149
1st Qu.:1.195 1st Qu.:1.079 1st Qu.:0.729 1st Qu.:1.277
Median :1.239 Median :1.129 Median :0.789 Median :1.345
Mean :1.214 Mean :1.106 Mean :0.772 Mean :1.353
3rd Qu.:1.259 3rd Qu.:1.155 3rd Qu.:0.823 3rd Qu.:1.410
Max. :1.360 Max. :1.229 Max. :1.039 Max. :1.599
NA’s :476 NA’s :126 NA’s :7515 NA’s :9939
PNuevoGasoleoA PBioDiesel PGasolina98 PGasNaturalComp
Min. :0.814 Min. :0.909 Min. :0.890 Min. :0.689
1st Qu.:1.170 1st Qu.:1.039 1st Qu.:1.329 1st Qu.:0.853
Median :1.204 Median :1.089 Median :1.360 Median :0.898
Mean :1.184 Mean :1.081 Mean :1.340 Mean :0.859
3rd Qu.:1.220 3rd Qu.:1.126 3rd Qu.:1.389 3rd Qu.:0.898
Max. :1.329 Max. :1.232 Max. :1.555 Max. :0.995
NA’s :2834 NA’s :9873 NA’s :3429 NA’s :9914
PGasNaturalLic PGasLicuPetro
Min. :0.652 Min. :0.469
1st Qu.:0.689 1st Qu.:0.582
Median :0.689 Median :0.612
Mean :0.771 Mean :0.603
3rd Qu.:0.898 3rd Qu.:0.638
Max. :0.898 Max. :0.727
NA’s :9932 NA’s :9462
Podemos observamos la fuerte presencia de valores perdidos, dichos valores pueden ser un problema cuando se analiza un conjunto de datos.
En algunos casos, si el numero de datos perdidos es relativamente pequeño la eliminación de ellos puede ser una estrategia a seguir. En nuestro caso debido a la presencia significativa de valores perdidos, optamos por la opción de estimarlos mediante técnicas de Machine Learning.
Tratamiento de valores Missing
Una primera visualización de los valores perdidos nos la proporciona el paquete VIM.
library(VIM)
aggr_plot <- aggr(Fuel, numbers=TRUE, sortVars=TRUE, labels=names(Fuel), cex.axis=.5, gap=3, ylab=c("Histograma de valores perdidos","Pattern"))
La visualización nos permite ver que casi un 75% de los registros no presentan missing, el 22% de los missing es debido a la variable “Precio Bio Etanol”.
Igualmente podemos realizar otra aproximación visual comparando dos variables:
marginplot(FuelPrecio[,c(“PGasolina95″,”PGasolina98”)])
Aquí podemos ver la distribución de valores perdidos del “Precio de la gasolina 98” con falta de valores del “Precio de la gasolina 95”. La parte azul nos muestra la distribución del resto de puntos.
Imputación de Misssing
A través de la librería MICE podemos asignar valores perdidos mediante técnicas predictivas.
library(mice)
Fuelmice<- mice(FuelPrecio,m=5,maxit=50,meth='rf',seed=500)completedFuel <- complete(Fuelmice,1) summary(completeFuel)PGasolina95 PGasoleoA PGasoleoB
Min. :0.779 Min. :0.749 Min. :0.544
1st Qu.:1.189 1st Qu.:1.079 1st Qu.:0.720
Median :1.235 Median :1.129 Median :0.789
Mean :1.213 Mean :1.106 Mean :0.768
3rd Qu.:1.259 3rd Qu.:1.155 3rd Qu.:0.820
Max. :1.360 Max. :1.229 Max. :1.039
PBioEtanol PNuevoGasoleoA PBioDiesel
Min. :1.149 Min. :0.814 Min. :0.909
1st Qu.:1.283 1st Qu.:1.155 1st Qu.:1.064
Median :1.345 Median :1.194 Median :1.118
Mean :1.346 Mean :1.174 Mean :1.108
3rd Qu.:1.399 3rd Qu.:1.219 3rd Qu.:1.155
Max. :1.599 Max. :1.329 Max. :1.232
PGasolina98
Min. :0.890
1st Qu.:1.319
Median :1.359
Mean :1.338
3rd Qu.:1.384
Max. :1.555
En este caso utilizamos como modelo de imputación un Random Forest.
Fuel <- select(Fuel, -contains(“Precio”), completeFuel)
Clustering
Determinaremos una clusterización que nos permita comparar estaciones de servicio similares y diferentes entre si, de cara a poder seleccionar y encontrar estaciones de servicio adaptadas a nuestras necesidades.
Para ello determinaremos con las variables numéricas, una vez escaladas, las componentes principales y determinaremos el número optimo de segmentos.
# Seleccionamos los atributos numéricos y escalamos
Fuel_num
Fuel_scale # Aplicamos análisis de componentes principales
pc plot(pc)
plot(pc, type = “l”)
summary(pc)
# Obtenemos componentes principales
pc summary(pc)# Primeras componentes principales
comp# Plot
plot(comp, pch = 16, col= rgb(0,0,0,0.5))
wssfor (i in 2:15)
wss[i]
plot(1:15, wss, type = “b”, xlab = “Number of Cluster”, ylab = “Within groups sum of squares”)
# K-Means Clustering con 4 clusters
Fuel_Cluster Fuel_Cluster
aggregate(Fuel[,16:22], by=list(Fuel_Cluster$cluster), FUN = mean)library(ggplot2)
Fuel$cluster
Fuel$PC1 Fuel$PC2
ggplot(Fuel, aes(PC1, PC2, color = Fuel$cluster)) +
geom_point()
El siguiente paso es ponerle cara a los clusters:
aggregate(Fuel[,16:22], by=list(Fuel_Cluster$cluster), FUN = mean)
Group.1 PGasolina95 PGasoleoA PGasoleoB PBioEtanol PNuevoGasoleoA
1 1 0.9402179 0.8586645 0.6861961 1.291664 0.9427451
2 2 1.2504438 1.1443182 0.8026884 1.472251 1.2108244
3 3 1.1654730 1.0503863 0.6936435 1.304441 1.1177545
4 4 1.2453041 1.1389383 0.7956955 1.301777 1.2053461
PBioDiesel PGasolina98
1 1.045405 1.049529
2 1.136626 1.372269
3 1.044700 1.299260
4 1.132222 1.368115
En este caso los cluster1 corresponden con estaciones de servicio económicas frente al cluster2 donde los precios son mucho mas altos.
Visualización
Una de las formas mas atractiva de visualizar los datos es a través de los mapas. En este caso utilizaremos todo el potencial de las librerías ggmap y ggplot2.
En nuestro Dataset, ya disponemos tanto de la longitud como de la latitud que nos permitirá georefenciar la estación de servicio y de esta forma visualizar su Cluster en mapa.
library(ggmap)
library(ggplot2)
Spain.map = get_map(location = “Spain”, zoom = 6, color=”bw”)
p p p p + theme(axis.title=element_blank(),
axis.text=element_blank(),
axis.ticks=element_blank()) + ggtitle(“Clusters Estaciones Servicio”)
Podemos visualizar los clusters de manera individual:
Podemos hacer zoom y ver el caso particular de la provincia de Almería.
El próximo 12 de Enero a las 20:oo h recibiremos el año con una nueva reunión del “Grupo de Usuarios de R de Almería”, en la que tendremos dos presentaciones muy interesantes.
El lugar elegido para las presentaciones son las instalaciones Workspace Coworking en Plaza Administración Vieja 1, 1º izq.
¡Os esperamos!
a 11 de Octubre, 2016
Esta claro que el mundo de fútbol levanta pasiones, el fútbol y todo lo que lo rodeada, entre ellas algunas web que nos permiten gestionar de manera virtual a nuestro equipo, además de recoger los puntos de la jornada sirve para confeccionar el once inicial.
En este caso, nuestro objetivo es poder llegar a segmentar a los jugadores de la liga profesional de fútbol en función a una serie de inputs que la web pone a nuestra disposición, y de esta forma poder ayudar en la toma decisiones a la hora de gestionar mi equipo. http://www.comuniazo.com/
En primer lugar necesitamos extraer la información necesaria para seleccionar los atributos que determinaran la segmentación. Para ello, realizamos scraping a la web “comuniazo” a través de la siguiente url:
url <- "http://www.comuniazo.com/comunio/jugadores"
La libreria XML nos permite extraer información de la página web mediante la función readHTMLTable
library(XML)
jugadores= readHTMLTable(jugadores,
stringsAsFactors = T,
colnames = c("Posicion","Equipo","Jugador","Puntos","Media","Puntos_Casa","Media_Casa","Puntos_Fuera","Media_fuera", "Valor"),
colClasses = c("character","character","character","FormattedNumber","FormattedNumber","FormattedNumber","FormattedNumber","FormattedNumber","FormattedNumber"))
El dataset esta formado por 504 registros y 14 variables candidatas.
Análisis de componentes principales
A continuación vamos a aplicar una de las técnicas de reducción de dimensionalidad , aplicamos análisis de componentes principales de los datos en R, para ello previamente escalaremos el dataset:
# Seleccionamos los atributos numéricos y escalamos jugadores_num <- jugadores[,4:10] jugadores_scale <- scale(jugadores_num) # Aplicamos análisis de componentes principales pc <- princomp(jugadores_scale) plot(pc) plot(pc, type = "l") summary(pc)
Seleccionando dos componentes principales conseguimos explicar el 86% de la varianza.
# Obtenemos componentes principales pc <- prcomp(jugadores_scale) summary(pc) Importance of components: Comp.1 Comp.2 Comp.3 Comp.4 Standard deviation 2.2477665 0.9923867 0.71057691 0.50588385 Proportion of Variance 0.7232142 0.1409699 0.07227477 0.03663246 Cumulative Proportion 0.7232142 0.8641840 0.93645881 0.97309127 Comp.5 Comp.6 Comp.7 Standard deviation 0.37326927 0.22058428 0 Proportion of Variance 0.01994385 0.00696488 0 Cumulative Proportion 0.99303512 1.00000000 1 # Primeras componentes principales comp <- data.frame(pc$x[,1:2])
Análisis Cluster
Vamos a segmentar con la técnica más simple, para ello utilizaremos k-medias en R.
En primer lugar, ¿cómo determinamos el número de grupos? El método más simple es mirar a la suma de cuadrados dentro de los grupos y escoger el “punto” en la gráfica, de forma similar a como ocurre con el gráfico de sedimentación que hicimos para el PCA previamente. A continuación mostramos el código utilizado:
# Plot plot(comp, pch = 16, col= rgb(0,0,0,0.5)) wss <- (nrow(jugadores_scale)-1)*sum(apply(jugadores_scale,2,var)) for (i in 2:15) wss[i] <- sum(kmeans(jugadores_scale, centers =i)$withinss) plot(1:15, wss, type = "b", xlab = "Number of Cluster", ylab = "Within groups sum of squares")
# K-Means Clustering con 6 clusters jugadores_Cluster <- kmeans(comp, 6, nstart = 25, iter.max = 1000) jugadores_Cluster aggregate(jugadores[,4:10], by=list(jugadores_Cluster$cluster), FUN = mean) K-means clustering with 6 clusters of sizes 72, 114, 71, 64, 144, 39 Cluster means: PC1 PC2 1 -1.51933814 -1.1588460 2 0.81871056 -0.4593748 3 -0.09203462 0.9141368 4 -2.37569832 0.9286046 5 2.47658186 0.1106443 6 -4.66639217 -0.1143984 Clustering vector: [1] 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 4 6 6 6 6 6 6 6 6 4 6 4 6 6 6 6 6 [34] 6 6 4 4 6 6 1 4 4 6 6 6 4 1 4 1 4 4 4 1 4 1 4 4 1 4 4 4 1 4 4 6 1 [67] 4 1 4 1 1 1 4 4 4 4 4 4 4 1 6 1 4 4 4 1 1 4 4 4 1 4 4 4 4 1 4 4 4 [100] 1 4 4 4 4 1 4 4 1 4 1 4 4 1 1 4 4 1 1 4 1 1 1 1 1 1 4 4 1 1 4 4 1 [133] 4 3 3 3 1 3 1 1 3 3 1 3 1 4 1 1 1 3 1 1 1 3 1 4 3 3 1 1 3 3 1 1 1 [166] 1 2 1 3 3 2 1 3 3 1 3 1 3 3 1 3 1 1 1 1 3 3 3 2 3 2 3 3 2 4 1 1 3 [199] 1 2 3 3 1 2 1 2 2 3 2 3 3 3 2 2 2 2 3 3 2 1 1 3 2 3 3 3 3 3 4 4 3 [232] 1 2 2 2 3 3 2 2 3 3 2 2 2 3 2 2 2 2 1 3 3 2 2 3 3 2 2 2 2 2 3 2 2 [265] 2 3 2 3 2 2 2 2 2 2 2 2 2 2 2 3 3 2 3 2 2 2 2 2 2 2 2 2 2 3 2 2 2 [298] 3 3 3 2 2 2 2 2 3 3 3 3 2 1 3 2 2 2 2 2 2 2 2 5 2 2 3 5 2 2 2 3 2 [331] 3 2 2 2 2 2 2 2 2 2 2 2 2 2 2 5 2 5 5 2 2 5 3 5 2 2 2 2 5 5 5 5 2 [364] 5 2 2 2 2 2 5 2 5 2 2 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 [397] 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 [430] 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 [463] 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 [496] 5 5 5 5 5 5 5 5 5 Within cluster sum of squares by cluster: [1] 76.21358 75.80518 55.22837 67.37217 75.01832 118.73351 (between_SS / total_SS = 84.6 %) Available components: [1] "cluster" "centers" "totss" "withinss" [5] "tot.withinss" "betweenss" "size" "iter" [9] "ifault"
Una vez calculados los clusters correspondientes los visualizaremos en función de las dos componentes principales utilizadas:
library(ggplot2) jugadores$cluster <- as.factor(jugadores_Cluster$cluster) jugadores$PC1 <- comp$PC1 jugadores$PC2 <- comp$PC2 ggplot(jugadores, aes(PC1, PC2, color = jugadores$cluster)) + geom_point()
El siguiente paso sería identificar los segmentos y establecer una etiqueta que los identifique y poder gestionar la estrategia.
Hasta aquí llega este post. ¡Espero que le haya gustado! Cualquier pregunta o comentario, no dude en dejar un comentario o acercarse a mí en Twitter: @mirallesjm
No puedes perderte el evento nacional de Usuarios de R 2016
Albacete 17 y 18 de noviembre de 2016
Agrupación Politécnica Superior – Universidad de Castilla la Mancha Visita la web
