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):

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 (¡!)

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.

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.

El tema que estoy estudiando estos días es la regresión por mínimos cuadrados parciales, partial least squares (PLS). Para documentarme teóricamente y conocer las principales posibilidades de R estoy empleando este documento. Para argumentar el uso de esta técnica de nuevo partimos del modelo lineal general Y = X • Beta + Error donde Beta = inv(X’X) * X’Y y ya analizamos los trastornos que nos provoca la inv(X’X) cuando hay columnas de X que son linealmente dependientes, cuando hay multicolinealidad. En ese caso empleábamos la regresión ridge. Bueno, imaginemos esta situación, tenemos más variables que observaciones. Entonces si que no somos capaces de tener una solución para la inv(X’X). Para este problema contamos con los mínimos cuadrados parciales.

Por lo visto no he estudiado lo suficiente. Tengo que redimirme y estudiar este verano determinadas técnicas avanzadas de predicción. Fundamentalmente tengo que trabajar con R y tener determinados conocimientos teóricos sobre estas técnicas. Así que he pensado que, a la vez que estudio yo, estudian todos mis lectores. Además es probable que genere debate.

En esta primera entrega vamos a tratar la regresión contraída o regresión ridge. En el blog ya hablamos del problema que suponía la multicolinealidad cuando tenemos este problema una de las posibles soluciones es la regresión contraída o regresión ridge. Como ya dijimos el modelo lineal se expresa como Y = X • Beta + Error la estimación de nuestros parámetros Beta por mínimos cuadrados ordinarios es Beta = inv(X’X) * X’Y cuando X’X no es invertible tenemos un problema. La regresión ridge plantea una solución a este problema con unos parámetros Beta_contraidos = inv(X’X + lambda*I) * X’Y si lambda es 0 estamos ante mínimos cuadrados ordinarios, en otro caso estamos ante un estimador sesgado de Beta. Este estimador sesgado es solución al problema de mínimos cuadrados penalizados y lo que hace es contraer los Betas en torno a 0. En resumen, metemos sesgo pero reducimos varianza.

Para determinar la distribución que sigue un vector de datos en R contamos con el paquete rriskDistributions. Este paquete de R nos permite realizar un test para las distribuciones siguientes:

• Normal
• Logística
• Uniforme
• Gamma
• Lognormal
• Weibull
• Cauchy
• Exponencial
• Chi-cuadrado
• F
• T-Student

Todos aquellos que estén trabajando con los modelos de supervisión de riesgos seguramente conocerán este paquete y si no lo conocen espero que lean estas líneas porque pueden ser de mucha ayuda para ellos, aunque se trate de software libre, no pasa nada, no receléis de R. La sintaxis es tan sencilla que se puede resumir en:

Buena noticia. Ya están en marcha las VI Jornadas de usuarios de R. Este año serán en Santiago de Compostela los días 23 y 24 de octubre y serán retrasmitidas por video streaming. Miguel Ángel Rodríguez Muiños debe tener mucha culpa de que este año se vayan a tierras gallegas donde hace mucho tiempo son referencia en formación sobre el uso de software libre.

A qué voy…

El modelo líneal se puede escribir de forma matricial como Y = X • Beta + Error. Donde Y es el vector con nuestra variable dependiente, X la matriz con las variables regresoras, Beta el vector de parámetros y el error esa parte aleatoria que tiene que tener todo modelo. La matriz con nuestras variables regresoras X ha de tener rango completo es decir, todas sus columnas tienen que ser linealmente independientes. Eso nos garantiza que a la hora de estimar por mínimos cuadrados ordinarios X’ X es invertible. Si no es invertible la estimación por mínimos cuadrados ordinarios “se vuelve inestable” ya que X’X =0 y 1/ X’X será muy complicado de calcular ya que los Beta son inversa(X’ X) •X’Y; por ello los resultados que arroja el modelo tienen una alta variabilidad. Cuando esto nos pasa tenemos colinealidad.

Los días 8, 9 y 10 de Abril de 2014.

Este curso pretenden introducir a los investigadores al análisis de datos que se están generando con las nuevas tecnologías de secuenciación, así como facilitar a aquellos que ya disponen de datos, su análisis e interpretación de resultados. El curso es un curso práctico en el que se utilizará una metodología similar a la que se ha usado en otros cursos anteriores y que tan buenos resultados ha obtenido.

La pasada semana, en un examen, me preguntaron cuál era el mejor punto para una prueba diagnóstica; era necesario razonar mi respuesta. Seguramente mi respuesta fue correcta pero mi razonamiento no lo fue y por eso quería redimirme. Para evaluar las pruebas diagnósticas con una respuesta binaria si/no contamos con la sensibilidad y la especificidad. La sensibilidad es la capacidad que tiene la prueba para acertar sobre los que de verdad tiene que acertar, la probabilidad de etiquetar como enfermos aquellos que verdaderamente están enfermos. La especificidad es una medida que nos indica cuanto nos hemos equivocado con los “unos”, la probabilidad de etiquetar enfermos a pacientes sanos. Una forma de medir cuanto acertamos y cuanto nos equivocamos con nuestra prueba. Para analizar el comportamiento de nuestra prueba diagnóstica debemos determinar un punto de corte. Para ilustrar como seleccionar el mejor punto de corte vamos a emplear unos datos sacados de la web de bioestadística del Hospital ramón y Cajal y vamos a elaborar una curva ROC con R y ggplot2.

Buena iniciativa en la que seguro que algún lector está interesado:

http://www.galpon.org/content/2a-xornada-dos-sabados-libres-na-altamar-estatistica-con-r-i

En este enlace podéis ver con más detalle el contenido del curso. Saludos.

Me llena de orgullo y satisfacción mostraros un ejemplo de uso de la librería wordcloud para la realización de nubes de palabras con R. Esta entrada no es muy innovadora porque ya tenemos alguna similar en el blog. Lo primero que tenéis que hacer es descargaros el discurso del Rey y ejecutad este código:

#Lectura del archivo
ubicacion="C:\\temp\\juancar.txt"
texto = read.table (ubicacion,sep="\r")

#Dejamos todas las palabras en mayúsculas
texto = toupper(textoV1)
#El texto lo transformamos en una lista separada por espacios
texto_split = strsplit(texto, split=" ")

#Deshacemos esa lista y tenemos el data.frame
texto_col = as.character(unlist(texto_split))
texto_col = data.frame(texto_col)
names(texto_col) = c("V1")

#Eliminamos algunos caracteres regulares
texto_colV1 = sub("([[:space:]])","",texto_colV1)
texto_colV1 = sub("([[:digit:]])","",texto_colV1)
texto_colV1 = sub("([[:punct:]])","",texto_colV1)
#Creamos una variable longitud de la palabra
texto_collargo = nchar(texto_colV1)

#Quitamos palabras cortas
texto_col = subset(texto_col,largo>4)

#Nube de palabras
#install.packages('wordcloud')
library(wordcloud)
library(RColorBrewer)
pesos = data.frame(table(texto_colV1))

#Paleta de colores
pal = brewer.pal(6,"RdYlGn")

#Realizamos el gráfico
png('C:\\temp\\Discurso del rey españa 2013.png', width=500, height=500)
wordcloud(pesosVar1,pesosFreq,scale=c(4,.2),min.freq=2,
max.words=Inf, random.order=FALSE,colors=pal,rot.per=.15)

dev.off()

Interesante el uso de la librería RColorBrewer. Particularmente me gusta mucho el resultado que nos da wordcloud para la realización de las nubes de palabras con una sintaxis sencilla. Considero imprescindible el uso de ramdom.order=FALSE. Espero que os sea de utilidad.