Procesos de Poisson no homogéneos: la historia de un fracaso

Partamos el tiempo en, p.e., días y contemos una serie de eventos que suceden en ellos. Es posible que esos recuentos se distribuyan según un proceso de Poisson de parámetro $\lambda$, que es un valor que regula la intensidad. Si los días son homogéneos, i.e., no hay variaciones de intensidad diaria, estimar $\lambda$ (por máxima verosimilitud), es tan fácil como calcular la media de los sucesos por día. Pero puede suceder que la intensidad varíe en el tiempo (p.e., se reduzca los fines de semana). O que fluctúe de cualquier manera. O que haya periodos de gran intensidad y otros de calma. Es decir, que el proceso no sea homogéneo y que $\lambda$ varíe en el tiempo. ...

8 de agosto de 2014 · Carlos J. Gil Bellosta

Coclustering con blockcluster

Guardo desde hace un tiempo el enlace al paquete blockcluster de R que igual puede ser del interés de alguno de mis lectores. No lo he probado pero sospecho que cualquier día me puede sacar de un apuro. Implementa lo que dice, el coclústering, concepto que se explica mejor, como el efecto de las dietas milagrosas, con la foto del antes y el después: Esto es: la entrada es una matriz y la salida es una matriz reorganizada tanto en sus filas como en sus columnas en la que se han detectado bloques homogéneos. ...

1 de agosto de 2014 · Carlos J. Gil Bellosta

Incrementalidad via particionamiento recursivo basado en modelos

Planteas un modelo tal como resp ~ treat y no encuentras diferencia significativa. O incluso puede ser negativa. Globalmente. La pregunta es, con el permiso del Sr. Simpson (o tal vez inspirados por él), ¿existirá alguna región del espacio en la que el tratamiento tiene un efecto beneficioso? Puede que sí. Y de haberla, ¿cómo identificarla? De eso hablo hoy aquí. E incluyo una protorespuesta. Primero, genero datos: n <- 20000 v1 <- sample(0:1, n, replace = T) v2 <- sample(0:1, n, replace = T) v3 <- sample(0:1, n, replace = T) treat <- sample(0:1, n, replace = T) y <- v1 + treat * v1 * v2 y <- exp(y) / (1 + exp(y)) y <- sapply(y, function(x) rbinom(1,1,x)) dat <- data.frame( y = y, treat = factor(treat), v1 = v1, v2 = v2, v3 = v3) Como puede apreciarse, solo las variables v1 y v2 (y no v3) interaccionan con el tratamiento: solo en la región donde v1 = v1 = 1 el efecto del tratamiento es positivo. ...

30 de julio de 2014 · Carlos J. Gil Bellosta

Estrategias escalables (con R)

Hay quienes preguntan cómo cargar con R un csv de 8GB en un portátil de 4GB de RAM. La verdad, he leído respuestas la mar de extravagantes a este tipo de cuestiones: p.e., recomendar SQLite. Yo recomendaría Scalable Strategies for Computing with Massive Data. Entre otras cosas, porque para eso lo escribieron sus autores: para que se lea. Y porque está cargado de razón y buenos consejos. Una cosa con la que tropezará enseguida quien lo hojee es: ...

9 de julio de 2014 · Carlos J. Gil Bellosta

Vectorización en R: un contraejemplo

No hay regla sin excepción, dicen. Para la recomendación casi única para quienes se quejan de la lentitud de R, es decir, ¡vectoriza!, he encontrado hoy una. Sí, el artículo deja R por los suelos. En el fondo, no tanto, porque viene a decir que R es malo para lo que la documentación de R dice que es malo: véase cómo en Writing R Extensions nos advierten que la convolución is hard to do fast in interpreted R code, but easy in C code. Y el problema que tratan de resolver los autores contiene una convolución (a través de una cadena de Markov, para pasar de un nivel de capital al del siguiente periodo). Es decir, en cierta medida solo viene a confirmar que la documentación de R es buena. ...

4 de julio de 2014 · Carlos J. Gil Bellosta

Disponible una nueva versión de MicroDatosEs

Acabo de subir a CRAN una nueva versión de MicroDatosEs, un paquete para procesar automáticamente en R ficheros de microdatos públicos españoles. A los cambios y mejoras a los que me referí el otro día, esta nueva versión añade otras, obra de Carlos Neira, que es ahora contribuidor oficial del paquete. Carlos también contribuyó a detectar y corregir un error inducido por el INE, que cambió el formato del fichero introduciendo una nueva variable sin aviso previo. ...

27 de junio de 2014 · Carlos J. Gil Bellosta

Grupo de usuarios de R de Portugal

Nuestros vecinos portugueses acaban de abrir un foro para sus usuarios de R un poco al estilo de nuestro r-help-es (¿todavía no estás dado de alta en él?). Espero que les sirva de base para organizar una comunidad vibrante de usuarios. Y que algún día podamos organizar unas Jornadas Ibéricas de Usuarios de R.

23 de junio de 2014 · Carlos J. Gil Bellosta

Agrupación de grafos por topología

Anuncio algo que no he conseguido hacer: agrupar grafos por topología. Pero no me he quedado lejos. Y espero que si alguien tiene alguna idea al respecto, nos lo haga saber al resto en la coda. Contexto (disfrazado). Hay usuarios que tienen correos electrónicos. La relación esperada es de uno a uno. Pero la realidad es, como siempre, mucho más compleja: hay usuarios que tienen varios correos y correos compartidos por varios usuarios. ...

13 de junio de 2014 · Carlos J. Gil Bellosta

A vueltas con el t-test

Me gustaría no tener que hacer más t-tests en la vida, pero no va a ser el caso. El problema al que me refiero le surgió a alguien en una galaxia lejana y, de alguna manera, me salpicó y me involucró. Es, simplificándolo mucho, el siguiente. Tiene una muestra $X = x_1, \dots, x_n$ y quiere ver si la media es o no cero. ¿Solución de libro? El t-test. Pero le salen cosas raras e inesperadas. De ahí lo del salpicón. ...

10 de junio de 2014 · Carlos J. Gil Bellosta

Validación cruzada en paralelo

Estoy sin tiempo, así que os suelto el código y me largo a casa a no cenar. Es así: library(parallel) cl <- makeCluster(8) # solo si hay aleatorización # clusterSetRNGStream(cl, 123) clusterEvalQ(cl, { # las librerías necesarias tienen que cargarse # en cada esclavo library(rpart) # en la práctica, hay que cargar los datos # (¿desde fichero?) en cada esclavo my.data <- iris # lo mismo con las funciones necesarias foo <- function(x, dat){ train <- 1:nrow(dat) %% 10 != 1 mod <- rpart(Species ~ ., data = dat[train,]) res <- predict(mod, dat[!train,]) } }) res <- parSapply(cl, 0:9, function(x) foo(x, my.data), simplify = F)

6 de junio de 2014 · Carlos J. Gil Bellosta