Estadística

Los actuarios fueron tal vez los primeros científicos de datos. Aparentemente, la primera tabla de mortalidad fue creada por John Graunt en 1662.

Los actuarios hablan de la esperanza de vida. Pero no son tan pendejos como los científicos de datos de hoy en día en eso de pretender que la esperanza de vida es la vida. Los actuarios saben que en una cohorte con una esperanza de vida de 78 años habrá quien muera a los tres, a los quince, a los cincuenta y a los noventa. Nadie les exige que acierten cuándo va a morir cada cual. En lugar de eso, estudian la distribución de los fallecimientos a lo largo del tiempo y calculan indicadores útiles para sus clientes.

El primer enlace de la selección de esta semana es The evolution of the modern census. Todos sabemos que lo que llevó a José y María a Belén hace más de 2000 años fue dizque tenían que censarse. Hay noticias de censos anteriores. Desde entonces hasta ahora ha habido muchos, muchísimos censos, pero su mismo concepto y finalidad ha ido cambiando a lo largo de la historia: ya no se trata solamente de contar, medir la riqueza o el poderío militar. Ahora nos interesan otros aspectos relacionados ya no tanto con el cuántos sino con el cómo somos.

Después de hablar con tirios y troyanos sobre mi entrada sobre los efectos de binarizar una variable objetivo continua, he decidido tomarme la justicia por mi mano y llamar a la caballería. Es decir, utilizar random forests.

Aquí va el código:

library(randomForest)

set.seed(1234)

my.coefs <- -2:2
n <- 200
train.n <- floor(2*n/3)

test.error <- function(){
  X <- matrix(rnorm(n*5), n, 5)
  Y <- 0.2 + X %*% my.coefs + rnorm(n)
  Y.bin <- factor(Y>0)

  train <- sample(1:n, train.n)

  X <- as.data.frame(X)
  X$Y <- Y

  modelo <- randomForest(Y ~ .,
    data = X[train,])
  pred <- predict(modelo, X[-train,])
  error.cont <- length(pred) -
    sum(diag(table(pred >0, Y[-train]>0)))

  X$Y <- Y.bin
  modelo <- randomForest(Y ~ .,
    data = X[train,])
  pred <- predict(modelo, X[-train,])
  error.bin <- length(pred) -
    sum(diag(table(pred, Y.bin[-train])))

  data.frame(error.cont = error.cont,
    error.bin = error.bin)
}

errores <- do.call(rbind,
  replicate(1000, test.error(), simplify = F))

sapply(errores, fivenum)

El resultado, si te interesa, en tu pantalla.

Supongamos que queremos construir un modelo para predecir quién ganará un determinado partido de baloncesto basándonos en datos diversos. Y en un histórico, por supuesto.

Podemos utilizar una regresión logística así:

set.seed(1234)

my.coefs <- -2:2
n <- 200
train.n <- floor(2*n/3)

test.error.glm <- function(){
  X <- matrix(rnorm(n*5), n, 5)
  Y <- (0.2 + X %*% my.coefs + rnorm(n)) > 0

  train <- sample(1:n, train.n)

  X <- as.data.frame(X)
  X$Y <- Y

  mod.glm <- glm(Y ~ ., data = X[train,],
    family = binomial)

  glm.pred <- predict(mod.glm, X[-train,],
    type = "response")

  error <- length(glm.pred) -
    sum(diag(table(glm.pred > 0.5, Y[-train,])))
}

errores.glm <- replicate(1000, test.error.glm())

El código anterior hace lo siguiente:

Acá va otra selección de cuatro enlaces relevantes –que no necesariamente nuevos— de la semana. El primero, Using Metadata to find Paul Revere recoge a modo de historia, que algunos encontrarán amena, una aplicación de rudimentos del álgebra lineal al análisis de redes sociales. Dada una matriz de incidencia A (personas que pertenecen a clubes) es posible calcular índices de proximidad entre personas (o entre clubes) calculando no más AA'. El resto hasta ganar el premio de Netflix es pura heurística.

Voy a comentar y recomendar hoy un artículo, Statistics and data mining: intersecting disciplines (lo siento, he perdido el enlace para su libre descarga), del siempre recomendable David Hand. Trata de un asunto que para muchos de los que seáis estadísticos y trabajéis en el asunto rodeados de gente procedente de otras disciplinas —¡ay, esos ingenieros!—, seguro, os produce dolores de cabeza: esa brecha que separa los mundos de la estadística y de la llamada minería de datos (y de otras maneras más recientemente).

Primero, canuto es como llamo yo a los streams: lugares por donde nos llegan datos —de sujetos que sacan dinero del cajero, de sensores de lo que sea, de clientes que llaman por teléfono—, típicamente en forma asíncrona. Es decir, que los eventos suceden cuando se les antoja, sin una periodicidad preestablecida.

Y uno quiere medir cosas. Por ejemplo, la frecuencia. O la intensidad de uso —¿cuánto se utiliza una tarjeta de crédito?—. Son estos fines para los que la gente todavía utiliza técnicas mandadas a recoger. Y ni siquiera. Podría contar la anécdota de un muy desavisado —por no decir cosas peores— que trabajaba conmigo no hace tanto. Se ufanaba de que tirando de su mucho ingenio y utilizando tecnologías punterísimas era capaz de asociar a cada sujeto de nuestra base de datos la serie temporal de sus actividades (asíncronas) durante los últimos meses. Y la verdad, a fuerza de tesón, lo consiguió durante un fin de semana en el que, me da la sensación, durmió poco. Llegó la reunión del lunes y nos contó con prolijidad de detalles sus muchos méritos y logros. Pero entonces se encogió un poco, se deshinchó un mucho y preguntó con mucha humildad: _y ahora que tengo todo esto, ¿qué más puedo hacer? _(entiéndase: para realizar un análisis estadístico del comportamiento de los sujetos).

Nunca he entendido eso de los odds. Me refiero a eso que mencionan las películas: ocho contra uno a favor de tal, cinco contra tres a favor de cual. Y no creo que sea el único al que le son ajenos. De hecho, la página de la Wikipedia en español correspondiente a la inglesa para odds se refiere a ellas como cuotas, término que jamás hasta hoy había visto así usado. Tampoco lo han visto, se concoce, los lexicógrafos de la RAE.

Un cascarrabias escribió el otro día el siguiente alegato en contra de la homeopatía: Contra la estupidez… y la estafa. Menciona incluso a gentes aún más cascarrabias que han propuesto el manifiesto No sin evidencia.

Sus razones tendrán, pero es manifiesto que nunca han leído el articulito que G.C. Smith y J.P. Pell escribieron en diciembre de 2003 en BMJ y que lleva por título Parachute use to prevent death and major trauma related to gravitational challenge: systematic review of randomised controlled trials.

No conocía el paquete gbm. Pero como ahora ando rodeado de data scientists que no son estadísticos…

Bueno, la cuestión es que había que ajustar un modelo para el que yo habría hecho algo parecido a

dat <- read.csv("http://www.ats.ucla.edu/stat/data/poisson_sim.csv")
summary(m.glm <- glm(num_awards ~ prog + math, family = "poisson", data = dat))
# Call:
#   glm(formula = num_awards ~ prog + math, family = "poisson", data = dat)
#
# Deviance Residuals:
#   Min       1Q   Median       3Q      Max
# -2.1840  -0.9003  -0.5891   0.3948   2.9539
#
# Coefficients:
#   Estimate Std. Error z value Pr(>|z|)
# (Intercept) -5.578057   0.676823  -8.242   <2e-16 ***
#   prog         0.123273   0.163261   0.755     0.45
# math         0.086121   0.009586   8.984   <2e-16 ***
#   ---
#   Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#
# (Dispersion parameter for poisson family taken to be 1)
#
# Null deviance: 287.67  on 199  degrees of freedom
# Residual deviance: 203.45  on 197  degrees of freedom
# AIC: 385.51
#
# Number of Fisher Scoring iterations: 6

como en esta página.