R

El jueves pasado y durante un breve receso de mi gripe, me entrevistaron en Canal Extremadura Radio. Durante una hora larga (que luego hubo que recortar a los 30 minutos que dura el programa de divulgación científica Principio de Incertidumbre) hablé de estadística, big data y R con Jorge Solís Bejarano.

A quien tengo que agradecer, primero, que contase conmigo; pero además y sobre todo, lo bien documentado que estuvo (lo cual me lleva a pensar que habrá que estar atentos a otras grabaciones de su programa).

Tenía en mente escribir estas líneas desde hace un tiempo. La reciente noticia de la adquisición de Revolution Analytics por parte de Microsoft la ha adelantado, como mucho, unos pocos días.

S, el lenguaje del que R es una implementación libre, vivió su ciclo propietario completo: nació en los laboratorios Bell, creció con Insightful, se reprodujo (R fue su vástago) y creo que ha muerto sin pena ni gloria en manos de Tibco. Como casi cualquier otro producto similar.

Preguntaba el otro día Emilio Torres esto en R-help-es. Resumo la pregunta. Se trata de una simulación de unos datos y su ajuste mediante una regresión logística para ver si los coeficientes obtenidos son o no los esperados (teóricamente y por construcción).

El código de Emilio (cuyos resultados no podemos reproducir porque no nos ha contado qué similla usa) es

logisticsimulation <- function(n){
  dat <- data.frame(x1=sample(0:1, n,replace=TRUE),
                    x2=sample(0:1, n,replace=TRUE))
  odds <- exp(-1 - 4 * dat$x1 + 7*dat$x2 - 1 *dat$x1* dat$x2 )
  pr <- odds/(1+odds)
  res <- replicate(100, {
    dat$y <- rbinom(n,1,pr)
    coef(glm(y ~ x1*x2, data = dat, family = binomial()))
  })
  t(res)
}

res <- logisticsimulation(100)
apply(res,2,median)
## (Intercept)          x1          x2       x1:x2
## -1.0986123 -18.4674562  20.4823593  -0.0512933

Efectivamente, los coeficientes están lejos de los esperados, i.e., -1, -4, 7 y 1.

En todos mis Linux, desde el principio de los tiempos, R guardaba los paquetes en

• /usr/lib/R/library
• /usr/lib/R/site-library (¡a veces y no sé por qué!)
• /usr/local/lib/R/site-library

Bajo /usr/lib deberían instalarse solo aquellos que vienen de serie con la instalación de R (o que se instalan usando el sistema de actualización de paquetes de la distribución de Linux) mientras que bajo /usr/local vivirían los instalados posteriormente por el usuario (véase esto).

Por supuesto, para escribir /usr/local/lib/R/site-library hacen falta permisos de superusuario y los paquetes ahí instalados están disponibles para todos los usuarios de la máquina. Pero de un tiempo a esta parte y por culpa, creo, de RStudio (tanto en versión de escritorio como de servidor), se me han comenzado a instalar paquetes en ~/R, bajo mi directorio personal. ¡Anatema!

Creo que todos sabéis la historia de las admisiones de la Universidad de Berkeley y la paradoja de Simpson. Con palabras, muchas palabras, está contado, por ejemplo, aquí. Y si buscáis ubc admissions simpson en Google la encontraréis también en modo --verbose en muchos más sitios.

En R puede resumirse en

library(reshape2)
library(plyr)

data(UCBAdmissions)

raw <- as.data.frame(UCBAdmissions)

dat <- dcast(raw, Gender + Dept ~ <a href="http://inside-r.org/packages/cran/AdMit">Admit)

mod.0 <- glm(cbind(Admitted, Rejected) ~ Gender, data = dat, family = binomial)
mod.1 <- glm(cbind(Admitted, Rejected) ~ Gender + Dept, data = dat, family = binomial)

Echad un vistazo a los coeficientes de Gender en ambos modelos y veréis.

library(plyr)

dat <- data.frame( a = sample(c("x", "y"),    100, replace = T),
                    b = sample(c(TRUE, FALSE), 100, replace = T))

ddply(dat, .(a), summarize, b = sum(b), no.b = sum(!b))
ddply(dat, .(a), summarize, no.b = sum(!b), b = sum(b))

Huele a bicho, ¿verdad?

Tengo por delante otro proyecto que tiene mucho de análisis exploratorio de datos. Sospecho que más de un árbol construiré. Los árboles son como la Wikipedia: prácticamente nunca el último pero casi siempre el primer recurso.

Esta vez, además, por entretenerme un poco, probaré el paquete [evtree](http://cran.r-project.org/web/packages/evtree/index.html). Aunque no porque espere sorprendentes mejoras con respecto a los tradicionales, ctree y rpart.

¿Qué tiene aquél que los diferencie de los otros dos? Que la optimización es global. Tanto ctree como rpart utilizan algoritmos recursivos: al definir un nuevo corte del espacio, el algoritmo solo tiene en cuenta la región definida por los cortes anteriores. La optimización es local. evtree utiliza un algoritmo global de la familia de los evolucionarios (¡qué tufillo a lentorro!). Los detalles están aquí.

Acabo de subir a Github una nueva versión de rPython-win, que soluciona uno de mis bugs históricos: ha pasado tanto tiempo en estado “pendiente” que casi le cojo cariño. Tiene (o tenía) que ver con particularidades no documentadas de las APIs para C de Python en distintas versiones de Windows y creo que no afecta al paquete en otras plataformas.

Y aprovechando que el Pisuerga pasa por Valladolid, un enlace: Calling Python from R with rPython.

Suelo trabajar un servidor con ocho CPUs. Cuando quiero paralelizar código en R, suelo utilizar [parallel::mclapply](https://stat.ethz.ch/R-manual/R-devel/library/parallel/html/mclapply.html) (como aquí). Pero no tengo una máquina. Tengo varias. Y antes, de hecho, muchas.

¿Cómo paralelizar en distintas máquinas?

Se puede usar Spark (y SparkR), por ejemplo. Pero una ruta que no había ensayado jamás es la de la vieja escuela, i.e., MPI, snow y demás.

Pero si

• tienes varios servidores corriendo un sistema operativo decente,
• instalas R y snow (y todo lo que necesites) en todos ellos y
• configuras los servidores para poder acceder a través de ssh sin contraseña desde uno central,

y, entonces, ejecutas

Me tengo que morder la lengua por no contarlo todo. Escribiré hasta donde pueda hacerlo. Que es casi nada. La cosa es que ha llegado a mis oídos que una muy importante empresa española con muchos, muchos empleados planea una migración muy seria de SAS a R.

Lo cual no deja de ser un cotilleo empresarial más. Que, como tal, no tendría cabida aquí. Salvo por el hecho de que me consta que me leen muchos estudiantes, muchos profesionales que se replantean sus carreras, muchos desempleados que se están formando de cara a su reincorporación.