martes, 8 de mayo de 2012

knitr: Un paquete avanzado para la edición de informes dinámicos

El paquete knitr de Yihui Xie se ha diseñado para superar algunos problemas de Sweave y avanzar en la edición de informes dinámicos con R.

Una fórmula del propio autor resume el conjunto de aplicaciones incluidas en el paquete:

knitr ≈ Sweave + cacheSweave + pgfSweave + weaver + R2HTML::RweaveHTML + highlight::HighlightWeaveLatex + 0.2 * brew + 0.1 * SweaveListingUtils + more

Las ventajas sobre Sweave son múltiples. En primer lugar y lo más evidente es el coloreado del código, pero también la habilidad para formatear el código automáticamente.
Este paquete incluye el soporte para varios dispositivos gráficos, entre ellos el tikzDevice que traduce el gráfico a código LaTeX y, por lo tanto, escribir fórmulas en LaTeX dentro del gráfico en la misma fuente que el texto.
Entre otras opciones, se permite controlar el tamaño de los gráficos o incluir dos gráficos de diversos tipos (normal y ggplot2, por ejemplo) en la misma figura.
El paquete knitr recoge lo mejor de otros paquetes que ya utilizaban el almacenamiento en caché para los objectos, de forma que se evitan cálculos repetidos innecesarios con el consiguiente ahorro de tiempo, incluso con los gráficos.

Además, los viejos archivos .Rnw de Sweave se pueden reciclar fácilmente.
Lo mejor es ver los ejemplos que incluye el propio manual y juzgar entonces sus muchas ventajas.

Para los que no saben LaTeX o no necesitan documentos que incluyan muchas fórmulas, el paquete knitr y la nueva versión 0.96 de RStudio (hoy todavía en desarrollo) incluyen la edición de documentos Markdown.
Un documento Markdown no es más que un documento de texto plano que incluye algunas marcas para darle formato y convertirlo fácilmente a otro formato como HTML.
En esta nueva versión de RStudio ya podremos crear archivos con formato Markdown y extensión .Rmd. Un botón MD en la barra de herramientas mostrará una guía rápida para aprender a escribir en este formato en 3 minutos.
El formato Markdown permite escribir títulos con # y listas con el guión -, de manera que es fácil de escribir y de recordar, al mismo tiempo que es legible.
Con el paquete knitr podemos intercalar código R y cuando pulsemos el botón Knit HTML obtendremos un documento HTML con el texto y los resultados. ¿No es fantástico? Una forma sencilla de redactar informes dinámicos.



Enlaces:

knitr: http://yihui.name/knitr/

Markdown + knitr + RStudio

Uso del formato Markdown para escribir páginas web

jueves, 8 de marzo de 2012

El mapa de Barcelona

Hace ya días que tenía ganas de publicar este mapa. En realidad desde que vi el mapa de Berlín o esta infografía del New York Times o los mapas de Oscar Perpiñán Lamigueiro.

Bueno, pues sin intentar ser tan espectacular, ya que no tengo los conocimientos suficientes, vamos a explicar de forma sencilla cómo se puede utilizar el mapa de Barcelona o de cualquier otra ciudad, región o estado.

En primer lugar debemos obtener el mapa propiamente dicho. En el caso de Barcelona nos vamos al Geoportalbcn del Ayuntamiento de Barcelona y concretamente a la pestaña de descargas. Allí nos podemos descargar varios archivos en diferentes formatos. A nosotros nos interesan los archivos de las Divisiones administrativas en formato SHP (ESRI shapefile). Estas divisiones administrativas son de diversos tipos: distritos, barrios, secciones censales,... Para trabajar con pocos datos tomaremos como unidades los 10 distritos. De modo que descomprimimos el archivo y nos preparamos para su lectura en R con la carga de algunos paquetes.

library(sp)
library(rgdal)
library(classInt)
library(RColorBrewer)

setwd("~/R/DIVADM_SHP")
ogrInfo(".", "BCN_Districte_SHP")

Para utilizar el paquete rgdal es preciso disponer de los programas adecuados que en linux se instalan con la siguiente instrucción:

$ sudo apt-get install libgdal1-dev libproj-dev

Las siguientes instrucciones leen el mapa y nos informan de su contenido:

ogrInfo(".", "BCN_Districte_SHP")
map.bcn <- readOGR(".",layer="BCN_Districte_SHP")
names(map.bcn)

En particular, los nombres de los distritos no son correctos seguramente por problemas de codificación, luego hay que arreglarlo:

map.bcn$NDistric <- c("Ciutat Vella","Eixample",
"Sants-Montjuïc",

"Les Corts", "Sarrià-Sant Gervasi", "Gràcia",
"Horta-Guinardó", "Nou Barris", "Sant Andreu",
"Sant Martí")

Ahora necesitamos una variable con los datos que queremos representar: demográficos, políticos,...
En este caso vamos a utilizar la proporción de mujeres en cada distrito ya que son datos incluidos en el propio archivo que hemos leído, pero que conste que podemos utilizar cualquier vector de datos, tantos como unidades territoriales.

plot_var <-
map.bcn@data[, "Dones"]/(map.bcn@data[,"Dones"]+map.bcn@data[,"Homes"])


A continuación preparamos los intervalos y los colores:

cuts <- 8
plot_pal <- brewer.pal(cuts, "Purples")
plot_intvl <- classIntervals(plot_var, cuts, style="pretty")
plot_colors <- findColours(plot_intvl, plot_pal)

y ya podemos dibujar el mapa:

plot(map.bcn, col=plot_colors, border="black")
title("Percentatge de dones a Barcelona")
text(map.bcn$Coord_X,map.bcn$Coord_Y,map.bcn$NDistric)

Finalmente podemos añadir la leyenda:

plot_lbls <- paste(round(100*plot_intvl$brks, 0), "%", sep="")[seq(2, length(plot_intvl$brks)-1)]
plot_lbls <- c(
paste("menys del", head(plot_lbls[1],1)),
paste(head(plot_lbls,-1), "-", tail(plot_lbls,-1)),
paste("més del", tail(plot_lbls,1))
)
legend("bottomright",
legend=plot_lbls, fill=attr(plot_colors, "palette"),
bty="n", border=attr(plot_colors, "palette"), cex=0.9
)


Referencias

http://www.nceas.ucsb.edu/scicomp/usecases/ReadWriteESRIShapeFiles


Maps with R de Oscar Perpiñán

domingo, 12 de febrero de 2012

Una animación con R



La idea es sencilla. Se trata de componer un vídeo mediante la unión de imágenes como en una película de dibujos animados.
En este ejemplo, vamos a utilizar el parámetro de amplitud h de una función coseno

f(t) = h*cos(omega(t-rho))

con omega=1 y rho=0 fijos.

En primer lugar fijamos los parámetros y hacemos un gráfico preliminar:

h <- 1
rho <- 0
omega <- 1

x <- seq(-1.5*pi+0.5, 2*pi+0.5, length=250)
y <- h*cos(omega*(x-rho))
plot(x, y, type="l", ylim=c(-1.5,1.5))
abline(h=0, v=0, col="grey75")



Ahora vamos a crear una secuencia de imágenes variando el parámetro h entre los valores 1 y 1.5 primero, después bajándolo hasta 0.5 y luego subiéndolo otra vez hasta 1.
Para un vídeo de 8 segundos necesitamos al menos 8*25=200 imágenes.

hnew <- c(seq(1, 1.5, by=0.01),
seq(1.48, 0.5, by=-0.01),
seq(0.52, 1, by=0.01))

También debemos crear una carpeta donde guardaremos todas las imágenes separadas del resto de documentos. La carpeta debe existir realmente para que la podamos utilizar.

# imagesDir <- "C:/Mis documentos/R/video"
imagesDir <- "/home/francesc/R/video"

A continuación abrimos el dispositivo para crear múltiples imágenes en formato PNG, dentro de la carpeta escogida y con una numeración de tres dígitos correlativa.

png(file.path(imagesDir, "Rplot%03d.png"))

Ahora ya podemos generar las imágenes con un for.

for (i in 1:length(hnew)){
plot(x, y, type="l", ylim=c(-1.5,1.5),
main=expression(f(t)==h*cos(omega*t-omega*rho)), lwd=2)
abline(h=0, v=0, col="grey75")
h <- hnew[i]
lines(x, h*cos(omega*(x-rho)), col=2)
lines(c(0,0),c(0,h), col=2, lwd=2)
lines(c(-0.05,0.05),c(h,h), col=2, lwd=2)
lines(c(-0.05,0.05),c(0,0), col=2, lwd=2)
lines(c(0,1.8), c(h/2,0.85), col=2, lwd=2)
text(1.9,0.87, paste("h=", round(h,2), sep=""), adj=0, col=2)
text(1.9,1.1, "Changing the amplitude:", font=2, adj=0, col=2)
}
dev.off()

La última instrucción cierra el dispositivo png.
El último paso es crear el vídeo con algún programa específico. Una opción es ImageMagick que se puede ejecutar desde la consola con la instrucción:

convert /home/francesc/R/video/*.png video.mpeg

En este caso no ha funcionado (otras veces sí) y he utilizado el programa OpenShot que permite crear vídeos con la importación de una secuencia de imágenes con el mismo nombre que nosotros le hemos dado: Rplot%03d.png.

¿Os animáis a crear vuestros vídeos?

miércoles, 8 de febrero de 2012

Sweave desde un editor de texto


En algún otro artículo ya he resaltado la importancia de utilizar Sweave para mezclar código de R con texto y producir informes automáticos con resultados estadísticos y gráficos dinámicos. Sweave se puede utilizar desde LibreOffice, no sin ciertas dificultades, con el paquete odfWeave, pero sin duda los mejores resultados se obtienen con el texto en formato LaTeX.

Supongamos que tenemos una buena instalación de R y de algún programa de LaTeX para nuestro sistema operativo. Entonces se puede ejecutar Sweave desde la consola de R o incluso desde una GUI como RStudio, pero es mucho mejor ejecutarlo desde el propio editor de LaTeX que seguramente es donde hemos escrito el documento .Rnw que queremos compilar. Este es el objetivo: ejecutar Sweave sobre el documento que estamos editando sin salir del editor.

Sweave en un sistema Windows

Supongamos que utilizamos Texmaker con la distribución MikTeX de LaTeX, aunque el procedimiento es similar si utilizamos otro editor y otra distribución de LaTeX.
Sweave viene con la distribución estándar de R y se puede ejecutar desde su consola, pero yo prefiero compilar el documento directamente desde el editor. Para ello necesitamos un archivo llamado Sweave.bat que se puede obtener en el CRAN dentro del archivo batchfiles_x.x-x.zip. Para instalar el archivo Sweave.bat debemos dejarlo en alguna carpeta que sea accesible en nuestro PATH o añadir dicha carpeta al PATH.
Ahora ya deberíamos ser capaces de ejecutar este .bat desde la consola, sólo necesitamos indicarle al editor donde encontrarlo.
Abrimos Texmaker y seleccionamos Opciones -> Configurar Texmaker e incluir la orden:

Sweave.bat --pdf %.Rnw

en alguna de las categorías, por ejemplo como Montaje rápido.

Un detalle importante: en Windows hay una dificultad añadida ya que necesitamos que LaTeX sea capaz de hallar el archivo Sweave.sty. Para ello debemos añadir el siguiente path

C:\Archivos de programa\R\R-x.x.x\share\texmf

en las MikTeX Options: Roots como se indica en esa página.

Sweave en un sistema Mac OS X

Los pasos para poner a punto un ordenador Mac OS X son incluso más sencillos que los explicados para un ordenador con Windows. Necesitamos el script Sweave.sh que hallaremos en el CRAN.
Para instalarlo, debemos copiarlo en

/usr/local/bin

y, dentro del programa Terminal, escribimos

sudo chmod +x /usr/local/bin/Sweave.sh

para hacerlo ejecutable.

Ahora ya podemos ejecutar este script desde la consola. Sin embargo, para que podamos ejecutarlo desde TEXshop, sólo hay que indicarle a TEXShop donde hallar dicho Sweave.sh.
Abrimos TEXShop y buscamos TEXShop→Preferences. Hacemos click en la pestaña Misc y escribimos

/usr/local/bin/Sweave.sh −ld

en el campo LaTeX Progam Personal Script.

Sweave en un sistema Linux

Para poner a punto un ordenador con Linux, el procedimiento es también muy sencillo, aunque los detalles pueden variar ligeramente dependiendo de la distribución y el editor utilizados. El ejemplo que explicaré se basa en Kubuntu Linux y el editor es Kile.
Para ejecutar Sweave desde el editor necesitamos el script Sweave.sh que podemos hallar en el CRAN. Luego instalaremos ese script en

/usr/local/bin

entonces abrimos la Konsole y escribimos

sudo chmod +x /usr/local/bin/Sweave.sh

para hacerlo ejecutable.
El siguiente paso es indicarle a Kile donde hallar este script. Abrimos Kile y hacemos click en Settings→Configure Kile. Nos vamos a Tools→Build y añadiremos una nueva herramienta que podemos llamar Sweave. En la ventana superior escribiremos la orden Sweave.sh y

−ld '%source'

en la de abajo.
Recordemos que el archivo Sweave.sty debe ser accesible.

Más información en el artículo:
http://www.scribd.com/doc/6451985/Learning-to-Sweave-in-APA-Style

sábado, 31 de diciembre de 2011

Sweave, FactoMineR y PDFs


El asunto es que estoy escribiendo un informe en LaTeX con Sweave. Lo hago así porque de esta forma me aseguro que añadir una variable o cambiar un dato no será un problema. El informe se volverá a generar en un minuto con todos los cambios.
El análisis que hago en concreto es un Análisis Factorial Múltiple o MFA (Multiple Factor Analysis) con el paquete FactoMineR. Se trata de un análisis de reducción de la dimensión cuando tenemos varios grupos de variables, tanto cuantitativas como cualitativas. Además este análisis permite estudiar los ejes principales mediante algunos gráficos y, como era de esperar, la representación de los individuos con las componentes principales.
En concreto utilizo PDFLaTeX directamente con el archivo .tex que resulta de procesar el archivo .Rnw con el Sweave. Así, los gráficos deben estar en formato PDF o PNG y no EPS. Ahora con Sweave esta es la opción por defecto.
Bueno, el caso es que los gráficos generados por el Sweave estaban vacíos. Eran archivos PDF con el nombre habitual por defecto, pero vacíos. ¿Sería un bug de Sweave? ¿Un error de mi instalación de R? ¿Un problema al cambiar de carpeta los archivos PDF?

Pues me costó un buen rato descubrir la razón del problema. Como siempre, lo primero fue aislar el error. No era un problema de Sweave, ni de mi instalación de R, ni de las opciones del Sweave (cambio de carpeta para las imágenes). Tuve que hacer unas cuantas pruebas, para descubrir que el problema estaba en los gráficos que generaba el paquete FactoMineR. En mi caso la función plot.MFA(). Tal vez, si hubiera leído mejor la ayuda de la función...

Bien, pues aquí va la explicación:
Resulta que para generar un gráfico en PDF, debemos desviar el dispositivo de salida gráfica hacia un archivo PDF. Eso se hace manualmente con la instrucción pdf() y seguramente lo debemos cerrar con la instrucción dev.off(). Eso es lo que hace automáticamente el Sweave.
Resulta que leyendo el código de la función plot.MFA() observé que, por defecto, hace una llamada a un nuevo dispositivo con la instrucción dev.new(). Ese era el problema.
Si abrimos el dispositivo pdf() y luego otro dispositivo con dev.new(), entonces el gráfico va al segundo dispositivo y al cerrarlo con dev.off(), cerramos el segundo y no el primero. La gente de FactoMineR tenía el problema previsto con un parámetro new.plot=TRUE que para mi gusto debería estar en FALSE por defecto. Ellos lo han puesto en TRUE para que podamos ver unos cuantos gráficos sin machacar el anterior, pero eso en mi caso era mortal.
En resumen, con la opción new.plot=FALSE, he podido generar el informe con PDFLaTeX sin problemas.
Moraleja: si pensáis en una función que genere varios gráficos seguidos, hay que pensar en un parámetro que regule la presentación uno a uno y que, mejor por defecto, no abra un nuevo dispositivo para cada uno.

sábado, 5 de noviembre de 2011

RKWard 0.5.7 para el nuevo R 2.14.0



Los desarrolladores de RKWard han lanzado una nueva versión (23 de octubre de 2011) de este magnífico software. Lo que me ha pasado es que al actualizar R a su versión 2.14.0, el RKWard dejó de funcionar. Pero como era de esperar, los amigos de RKWard ya estaban preparados.
Como yo utilizo Kubuntu (natty) he buscado el repositorio adecuado para actualizar RKWard que es

deb http://ppa.launchpad.net/rkward-devel/rkward-stable-cran/ubuntu natty main

y automáticamente se ha actualizado sin problema ninguno. Como veis todavía sigo con la versión natty o 11.04 de Kubuntu, ya que me gusta esperar algunos días antes de pasar a la nueva 11.10.

Además de nuevos detalles y corrección de errores, esta versión de RKWard viene con una importante mejora en la forma de añadir complementos, ya que ahora estos complementos se instalan como paquetes de R.

Aunque hay una versión para Windows, yo no la recomiendo ya que requiere un escritorio KDE para Windows y todavía no está demasiado desarrollada. Una opción alternativa y que comentaré próximamente es el programa multiplataforma RStudio que viene pegando fuerte.

martes, 12 de abril de 2011

El código de una función


Una de las ventajas de utilizar un programa como R de código libre es la posibilidad de aprender a través de la lectura del código de las funciones incorporadas, ya sea de los paquetes oficiales del CRAN o de los que proporcionan todas las personas que contribuyen de forma particular.

La forma más sencilla de acceder al código de una función es con la llamada a su nombre (sin paréntesis). Por ejemplo, para ver el código de la función matrix():

> matrix
function (data = NA, nrow = 1, ncol = 1, byrow = FALSE, dimnames = NULL)
{
data <- as.vector(data)
if (missing(nrow))
nrow <- ceiling(length(data)/ncol)
else if (missing(ncol))
ncol <- ceiling(length(data)/nrow)
.Internal(matrix(data, nrow, ncol, byrow, dimnames))
}
< environment: namespace:base >


Desgraciadamente los comentarios de estas funciones en los paquetes instalados se han eliminado para ahorrar memoria. Si queremos acceder al código fuente original, podemos descargar directamente el código de los paquetes del CRAN y descomprimirlos. Como los archivos están escritos con texto plano, podemos utilizar cualquier editor de texto para leer su contenido. Para los paquetes básicos de R, la carpeta es $R_HOME/src/library/NombrePaquete/R/.

En algunos casos, una función llama a otra función parcialmente oculta. Esa función oculta está en un namespace y se puede acceder a ella mediante la función getAnywhere("NombreFunción"). Esta función proporciona el namespace del que procede la función requerida y entonces podemos acceder a las fuentes de paquete correspondiente. Esto es así para los métodos S3 tales como plot.factor:

> plot.factor
Error: object 'plot.factor' not found
> getAnywhere("plot.factor")
A single object matching ‘plot.factor’ was found
It was found in the following places
registered S3 method for plot from namespace graphics
namespace:graphics
with value

function (x, y, legend.text = NULL, ...)
{
...
}
< environment: namespace:graphics >


El archivo que contiene el código de plot.factor es $R_HOME/src/graphics/R/plot.R.

Como ejemplo consideremos la función t.test(). Si queremos acceder a su código tenemos:

> t.test
function (x, ...)
UseMethod("t.test")
< environment: namespace:stats >

La llamada a UseMethod() oculta el código, pero nos indica que se trata de una función S3 genérica que llama a un método específico que será el apropiado para la clase del objeto. De modo que podemos preguntar cuales son los métodos apropiados así:

> methods(t.test)
[1] t.test.default* t.test.formula*

Non-visible functions are asterisked

De modo que ahora ya conocemos los métodos y podemos ver el código:

> getAnywhere(t.test.default)
A single object matching ‘t.test.default’ was found
It was found in the following places
registered S3 method for t.test from namespace stats
namespace:stats
with value

function (x, y = NULL, alternative = c("two.sided", "less", "greater"),
mu = 0, paired = FALSE, var.equal = FALSE, conf.level = 0.95,
...)
{
AQUÍ VEMOS EL CÓDIGO
}
< environment: namespace:stats >


En el caso de fuentes S4 es conveniente trabajar con los archivos fuente del paquete. Las funciones getClass(), getGeneric() y getMethod() nos pueden ayudar.

Cuando en el código de una función de R aparecen funciones del tipo .C(), .Call(), .Fortran(), .External() o .Internal() y .Primitive() es que llaman a código compilado. En estos casos deberemos mirar las fuentes originales (en C, C++ o Fortran) si queremos comprender todo el código.

Para més información consultar el artículo de Uwe Ligges:

R Help Desk en la Rnews_2006-4