Archivo de la categoría: Modelos

KNN con SAS. Mejorando K-Means

Imagen de previsualización de YouTube

La clasificación por k vecinos más cercanos es EL MÉTODO supervisado no paramétrico. El KNN, si empleamos las siglas en inglés, clasifica las observaciones en función de su probabilidad de pertenecer a uno u otro grupo, en el video que encabeza la entrada queda muy bien explicado. El caso es que tenemos la posibilidad de realizar esta clasificación con SAS STAT y el PROC DISCRIM y me parece interesante dedicarle unas líneas. Hace años ya hablamos de segmentación con SAS  y vamos a emplear los mismos datos para ilustrar esta entrada. Primero generamos un conjunto de datos con datos simulados de 3 esferas que clasificamos en 3 grupos:

data pelota;
do i = 1 to 1000;
a=0; b=5; x = a+(b-a)*ranuni(34);
a=0; b=5; y = a+(b-a)*ranuni(14);
grupo=1;
distancia = sqrt(((x-2.5)**2)+((y-2.5)**2));
if distancia < 2.5 then output;
end;
run;

data pelota1;
set pelota;
grupo=1;
run;

data pelota2;
set pelota;
x = x+4.5;
grupo=2;
run;

data pelota3;
set pelota;
x = x+2.5;
y = y+3.5;
grupo=3;
run;

data datos;
set pelota1 pelota2 pelota3;
run;

proc gplot data=datos;
	plot y * x = grupo; 
run;quit;

KNN_SAS1

Si realizamos un análisis mediante k-means sin asignar centroides obtenemos esta clasificación Seguir leyendo KNN con SAS. Mejorando K-Means

Valor atípico o pocos registros. Animación con R

outlier¿Cómo influye un solo punto en una recta de regresión? Evidentemente cuanto menos observaciones tengo más puede “descolocar” la recta de regresión. Sin embargo, cuantos más puntos tengo más complicado es encontrar ese punto con una recta de regresión, sin analizar los residuos podríamos hasta pasarlo por alto, aunque puede ser que nos interese ese punto. El código de R que genera la animación es:

library(animation)
saveGIF(
  for (i in c(100,50,25,10,5,1)){
    x <- seq(-500,500, by = i )
    y=sin(x)+x/100
    y[10]=y[10]+10
    plot (y,x,main=paste(“Regresión lineal con “,1000/i,” observaciones”))
    reg <- lm(y~x)
    points( fitted.values(reg),x, type=”l”, col=”red”, lwd=2)},
  interval = .85, ,movie.name=”/Users/raulvaquerizo/Desktop/R/animaciones/outlier.gif”)

Interpretación de los parámetros de un modelo GLM

Muchos estudiantes  terminarán trabajando con GLM que siguen buscando relaciones lineales en multitud de organizaciones a lo largo del planeta. Y hoy quería ayudar a esos estudiantes  a interpretar los parámetros resultantes de un GLM, más concretamente los resultados de un PROC GENMOD de SAS aunque lo que vaya a contar ahora se puede extrapolar a otras salidas de SAS o R. En la línea de siempre no entro en aspectos teóricos y os remito a los apuntes del profesor Juan Miguel Marín. Con un GLM al final lo que buscamos (como siempre) es distinguir lo que es aleatorio de lo que es debido al azar a través de relaciones lineales de un modo similar a como lo hace una regresión lineal, sin embargo los GLM nos permiten que nuestra variable dependiente no sólo siga una distribución normal, puede seguir otras distribuciones como Gamma, Poisson o Binomial. Además un GLM puede trabajar indistintamente con variables categóricas y numéricas pero yo recomiendo trabajar siempre con variables categóricas y en la práctica cuando realizamos un modelo de esta tipo siempre realizaremos agrupaciones de variables numéricas. Si disponemos de variables agrupadas, de factores, los parámetros de los modelos nos servirán para saber como se comporta nuestra variable dependiente a lo largo de cada nivel del factor.

El modelo siempre fija un nivel base del factor, un nivel que promedia nuestros datos y el resto de niveles corrigen el promedio en base al coeficiente estimado. Imaginemos que modelizamos el número de abandonos en 3 carreras de coches, cuando la carrera se disputa en un circuito A el número de abandonos es 5, sin embargo en el circuito B son 10 y en el circuito C son 15. Si fijamos como nivel base de nuestro factor circuito el B tendríamos un modelo de este modo abandonos = 10 + 0.5*es circuito A + 1*es circuito B + 1.5*es circuito C + Error. El nivel base promedia nuestro modelo por lo que va multiplicado por 1 y el resto de niveles se corrigen por su multiplicador. Esta es la base de la modelización multivariante en el sector asegurador. Veamos en un ejemplo como se articulan los parámetros de estos modelos. Simulamos unos datos con la probabilidad de tener un siniestro por edad, zona y edad:

data datos_aleatorios;
do idcliente = 1 to 2000;
if ranuni(1) >= 0.75 then sexo = “F”;
else sexo=”M”;
edad = ranpoi(45,40);
if ranuni(8)>=0.9 then zona=1;
else if ranuni(8)>0.7 then zona=2;
else if ranuni(8)>0.4 then zona=3;
else zona=4;
output;end;
run;

data datos_aleatorios;
set datos_aleatorios;
if zona=1 then incremento_zona = 0.1+(0.5-0.1)*ranuni(8);
if zona=2 then incremento_zona = 0.1+(0.7-0.1)*ranuni(8);
if zona=3 then incremento_zona = 0.1+(0.2-0.1)*ranuni(8);
if zona=4 then incremento_zona = 0.1+(0.9-0.1)*ranuni(8);

incremento_edad=exp(1/edad*10)-1;

sini = (ranuni(9) – sum(incremento_zona, -incremento_edad))>0.8 ;
run;

Se da una probabilidad aleatoria de tener un siniestro que se ve incrementada o decrementada por la zona y la edad, el sexo, aunque aparece, no influye. El número de siniestros suponemos que sigue una distribución de poisson. Para entender mejor como funciona un GLM vamos a agregar los datos por los factores en estudio y contamos el número de clientes sumando el número de siniestros:

proc sql;
create table datos_agregados as select
sexo,
case
when edad<= 30 then “1 menos 30”
when edad<= 40 then “2 31-40”
when edad<= 50 then “3 41-50”
else “4 mas 50” end as edad,
zona,
log(count(idcliente)) as exposicion,
sum(sini) as sini
from datos_aleatorios
group by 1,2,3;
quit;

Nuestros datos tienen que ir ponderados por el logaritmo del número de clientes, será nuestro offset, ya que no es lo mismo tener un siniestro en un grupo de 2 clientes que un siniestro en un grupo de 20 clientes. Ponderados por el logaritmo porque siempre cuesta menos trabajar con números pequeños y además tienen unos superpoderes de los que no somos conscientes hasta que trabajamos con ellos. Ahora estos datos son los que emplearemos para el modelo:

proc genmod data=datos_agregados;
class sexo edad zona;
model sini = sexo edad zona / dist = poisson
link = log
offset = exposicion;
run;

GENMOD como todos los procedimientos de SAS necesita que le indiquemos las variables categóricas, en model especificamos el modelo y las opciones son los aspectos más interesantes. En dist especificamos la distribución de nuestra variable dependiente en link la función de enlace que vamos a emplear y como offset para ponderar la variable exposición que hemos creado a la hora de agregar los datos. Y los parámetros que estima este modelo son:

parametros_modelo_GENMOD1

El intervalo de confianza de algunos estimadores contiene el valor 0, esos factores no son significativos como es el caso del sexo o bien puede haber agrupación de niveles de factores como es el caso de la edad entre 31 y 50 años o la zona 2 que puede unirse con otra zona. Vemos que el estimador del nivel base siempre es el último nivel del factor (esto puede cambiarse) y toma valores 0 y no 1 como habíamos usado en el ejemplo, para transformarlo en 1 sólo hemos de realizar el exponencial:

Parameter   Estimate
Intercept 0.062
sexo F 0.932
sexo M 1.000
edad 1 menos 30 2.951
edad 2 31-40 2.062
edad 3 41-50 1.494
edad 4 mas 50 1.000
zona 1 1.854
zona 2 1.218
zona 3 3.091
zona 4 1.000
Scale 1.000

Vemos que la zona 3 tiene casi el triple de siniestralidad que la zona 4 y lo mismo sucede con las edades jóvenes frente a las mayores edades, en cuanto al sexo que no fue significativo tenemos que las mujeres tienen un 7% menos de siniestralidad. Algunos resultados, aunque no salgan estadísticamente significativos, es evidente que pueden interesarnos comercialmente ya que mi producto puede dar un descuento a las mujeres y aunque sea pequeño se puede mantener, igual razonamiento para algunas zonas o grupos de edad susceptibles de unirse entre ellos. Este ejemplo es muy burdo pero aquellos que empiecen a trabajar con GLM se van a encontrar situaciones de este tipo, es necesario interpretar los parámetros estimados para describir como funciona el modelo pero igual de  importante es la agrupación de factores y el posterior suavizado de los parámetros.

Truco para EMB Emblem. Cambiar el nivel base de un factor

Un buen truco que me han descubierto hoy para los usuarios de EMB Emblem, como cambiar el nivel base de un factor de datos sin necesidad de pasar por los datos (habitualmente SAS) o sin hacerlo a posteriori (habitualmente Excel y lo que hacía el ahora escribiente). Cuando se generan los datos se genera el fichero binario *.BID y el fichero que se emplea para leer ese fichero *.FAC; para alterar el nivel base debemos abrir este archivo *.FAC con un block de notas o cualquier editor de texto plano. Al abrirlo tendremos lo siguiente:

XXXXXXXXXX        –> ES EL NOMBRE DEL ARCHIVO
99                                 –>ES EL NÚMERO DE FACTORES
XXXXXXXXXX.Bid –> ES EL NOMBRE DEL ARCHIVO BINARIO CON LOS DATOS
0
99999999                   –> ES EL NÚMERO DE REGISTROS
No. Levels Base         –> ES LA DESCRIPCION DE LO QUE VIENE A CONTINUACIÓN
XXXXXXXXXXXX  –> NOMBRE DE UN FACTOR
999                 99         –> EL PRIMERO ES EL NÚMERO DE NIVELES Y EL SEGUNDO ES EL NIVEL BASE

ASÍ PARA TODOS LOS FACTORES. Si deseamos cambiar el nivel base de un factor sólo tenemos que cambiar ese segundo número que hay tras cada definición de los factores. Un truco muy sencillo y que será de utilidad para aquellos que usáis EMB Emblem. Saludos.

Juego de modelos de regresión con R

Rplot

Os propongo un juego con R. El juego parte de unos datos aleatorios que he generado con R (los que veis arriba) que dividimos en entrenamiento y test. Sobre el conjunto de datos de entrenamiento he realizado varios modelos y valoro las predicciones gráficamente sobre los datos de test. El juego consiste en asociar cada resultado gráfico de test a cada código de R correspondiente y justificar brevemente la respuesta.

Los gráficos de los datos de test son:

Figura A:
Rplot01

Figura B:
Rplot02

Figura C:
Rplot03

Figura D:
Rplot05

Figura E:
Rplot07

Figura F:
Rplot08

Figura G:
Rplot06

Los códigos R que tenéis que asociar a cada figura son:

Código 1: Red neuronal con una sólo capa y 2 nodos:
mejor.red {
mejor.rss for(i in 1:50){
modelo.rn linout=T, trace=F,decay=0.1)
if(modelo.rn$value < mejor.rss){
mejor.modelo mejor.rss

return(mejor.modelo)
}}
}

mejor.red(2)

Código 2: Regresión lineal
lm(dep ~ indep,entrenamiento)

Código 3: Máquina de vector de soporte con un margen muy alto
svm(dep ~ indep ,entrenamiento, method=”C-classification”,
kernel=”radial”,cost=100,gamma=100)

Código 4: Árbol de regresión
rpart(dep~indep,entrenamiento)

Código 5: Regresión LOESS
loess (dep ~ indep, data = entrenamiento)

Código 6: Máquina de vector de soporte con un margen bajo
svm(dep ~ indep ,entrenamiento, method=”C-classification”,
kernel=”radial”,cost=10,gamma=10)

Código 7: K vecinos más cercanos K-nn
train.kknn(dep ~ indep, data = entrenamiento,
k = 4, kernel = c(“rectangular”))

Por ejemplo la figura A irá con el código 2 porque se trata de una estimación lineal. Y ahora os toca a vosotros asociar figuras a modelos de R.

Medir la importancia de las variables en una red neuronal con R

Sigo a vueltas con esta gran web y hoy vamos a medir la importancia de las variables en una red neuronal. Al igual que sucede en un modelo de regresión los parámetros obtenidos pueden servirnos para determinar la importancia de una variable dentro del modelo. En el caso de una red neuronal los pesos de la red pueden ser utilizados para determinar cómo influye una variable en el modelo. Para ilustrar este tipo de tareas el gran @beckmw realizó esta entrada:

http://beckmw.wordpress.com/2013/08/12/variable-importance-in-neural-networks/

El método que emplea para determinar esta importancia fue propuesto por Garson en 1991 y parte de la idea de buscar todas las conexiones que tiene una variable dentro de una red neuronal y se ponderan para obtener un valor único que describe el grado de asociación de nuestra variable dependiente con cada una de sus regresoras. Garson realizó este método para poder obtener un valor de 0 a 1, pero @beckmw lo ha modificado para que el signo tenga su importancia. La función que ha creado este genio la tenéis en Github

source(url(‘https://gist.githubusercontent.com/fawda123/6206737/raw/d6f365c283a8cae23fb20892dc223bc5764d50c7/gar_fun.r’))

Por cierto, qué manía tienen algunos con usar devtools. Haciendo sólo:

#Dejamos los nombres de los coeficientes de una forma más adecuada
tr=nnet.fit$coefnames
tr=substr(tr,16,30)
nnet.fit$coefnames=tr

#Pintamos la importancia de las variables
gar.fun(‘medv’,nnet.fit)

Obtenemos el gráfico con el que comienza esta entrada al blog. Ni se os ocurra comenzar a tocar este gráfico con los temas de ggplot2, somos gente productiva. Vemos como a la izquierda del gráfico se sitúan las variables con mayor peso negativo y a la derecha las variables con mayor peso positivo. Podemos eliminar algunas variables y seguramente el comportamiento predictivo del modelo no empeoraría.

Hay que dejar de pensar en las redes neuronales como una caja negra sin sentido. Saludos.

Representación de redes neuronales con R

En la última entrada realizamos un modelo de regresión con redes neuronales. Hoy quería mostraros como representar gráficamente la red neuronal creada en esa entrada. A la modelización con redes neuronales siempre se le ha achacado un comportamiento de “caja negra”, nosotros pasamos unas variables de entrada por una capa oculta y obtenemos una salida. No hay parámetros ni inferencia sobre los mismos, no sabemos lo que hace la red por dentro. En el caso concreto de R y continuando con la entrada anterior  si hacemos summary(bestnn):

a 12-4-1 network with 57 weights
options were – linear output units decay=1e-05
b->h1 i1->h1 i2->h1 i3->h1 i4->h1 i5->h1 i6->h1 i7->h1 i8->h1 i9->h1
-5.61 -3.80 -1.03 0.74 -2.81 2.83 2.37 2.86 6.72 4.68
i10->h1 i11->h1 i12->h1
1.65 0.87 -8.16 Seguir leyendo Representación de redes neuronales con R

Regresión con redes neuronales en R

La última técnica que me estoy estudiando este verano es la regresión con redes neuronales. El ejemplo que os voy a poner es completamente análogo a este: http://heuristically.wordpress.com/2011/11/17/using-neural-network-for-regression/ Vamos a trabajar con el paquete nnet, si dispusiera de tiempo replicaría este ejemplo en otra entrada con neuranet. Para realizar el ejemplo tenemos el conjunto de datos housing que contiene el censo de 1970 de 506 hogares de Boston. Empecemos a trabajar con la consola de RStudio (¡!)

#install.packages(“mlbench”)
library(mlbench)
data(BostonHousing)
summary(BostonHousing$medv)

Como variable dependiente vamos a emplear el valor mediano de las vivendas ocupadas en dólares.  El primer paso será realizar un modelo de regresión lineal:

lm.fit <- lm(medv ~ ., data=BostonHousing)
lm.predict <- predict(lm.fit)
mean((lm.predict – BostonHousing$medv)^2)

plot(BostonHousing$medv, lm.predict,
main=”Regresión lineal”,
xlab=”Actual”)

El ajuste ofrece una suma de cuadrados de 21.9  Vamos a realizar este modelo con redes neuronales Seguir leyendo Regresión con redes neuronales en R

Modelos lineales dinámicos (DLM) con R

Otro de los modelos que está tocando estudiar este verano son los Dinamic Linear Models (DLM). Para estudiar este tipo de modelos es imprescindible leer este documento.  Estos métodos parten de una idea: “la vida no es fácil cuando tienes que hacer estimaciones sobre una serie temporal”. Una serie temporal es un vector de datos aleatorios, una sucesión de observaciones de la forma Yt con t=1,2,.. Si sobre esta  sucesión tenemos una característica que puede influir estamos ante un modelo de espacio estado. Estos modelos tienen  una cadena de Markov (http://es.wikipedia.org/wiki/Cadena_de_M%C3%A1rkov) porque esa característica que afecta a la serie es una cadena de Markov y eso nos permite que los Yt sean independientes ya que dependen sólo de esa característica. El más importante modelo de espacio estado es el modelo lineal dinámico.

Los modelos lineales dinámicos vienen representados por una ecuación de observación Yt=F’·Bt + Vt con Vt según N(0,vt) y una ecuación de sistema o ecuación de estado Bt = Gt·Bt-1 + wt  con wt según N(0,wt). Luego está caracterizado por F y G que son matrices conocidas, B es un vector de parámetros desconocidos , Vt es la varianza de la ecuación de observaciones y Wt la varianza de la ecuación del sistema, ambas conocidas. El modelo más sencillo de este tipo es el paseo aleatorio, un modelo que no presenta estacionalidad ninguna. Sin entrar mucho más en temas matemáticos si es importante comentar que esta técnica la emplearemos tanto para filtrar, suavizar y predecir Seguir leyendo Modelos lineales dinámicos (DLM) con R

Primeros pasos con regresión no lineal (nls) con R

La regresión no lineal se da cuando tenemos que estimar Y a partir de una función del tipo Y=f(X,Beta) + Error donde Beta son Beta1, Beta2,…, Beta n. Unos datos X e Y se relacionan mediante una función no lineal respecto a unos parámetros Beta desconocidos. Y cómo obtenemos estos Beta desconocidos, a través de mínimos cuadrados o bien con otros métodos como máxima verosilimilitud. Este cálculo llevará asociada su inferencia estadística habitual. La función que asocia los pares de datos (x1,y1), (x2, y2),…, (yn, xn) será una función conocida. Por eso esta técnica es muy utilizada en ciencias químicas, geodinámica,… donde ya se conoce la relación que hay entre las variables independientes y la variable dependiente pero es necesario realizar modelos con los pares de datos disponibles de cara a obtener estimaciones.

Para la realización de este monográfico con R vamos a emplear el conjunto de datos Thurber Son datos de un estudio NIST de movilidad de electrones en semiconductores la variable respuesta es la movilidad del electrón y la variable regresora es el logaritmo de la densidad. El modelo es:

Nuestra variable está relacionada con su regresora por un modelo racional con siete parámetros Beta1, Beta2,…, Beta7 y cúbicas. Comenzamos el trabajo con los datos en R Seguir leyendo Primeros pasos con regresión no lineal (nls) con R