domingo, 26 de abril de 2009

Curso de análisis de datos en R - Sesión 1

Una introducción a R

R es un lenguaje orientado a objetos. Aunque existen algunas interfaces gráficas para R como Rcommander, es muy recomendable aprender R como un lenguaje en vez de tratarlo como un programa estadístico convencional. Como ambiente de trabajo, R ofrece una serie de ventajas:
  • Sus posibilidades gráficas son excelentes
  • Es muy flexible. Los procedimientos estadísticos estándar se pueden aplicar con sólo utilizar el comando apropiado. Además, existen multitud de librerías (a los que llamaremos paquetes de ahora en adelante) programadas por los usuarios de todo el mundo para llevar a cabo procedimientos específicos.
  • Es libre. Libre no quiere decir gratuito (aunque R también lo es). Libre significa que podemos acceder al código escrito por otros usuarios y modificarlo libremente. A pesar de que R viene sin garantía alguna (al iniciar la sesión de R saldrá la siguiente advertencia “R is free software and comes with ABSOLUTELY NO WARRANTY”), la mayor parte del código de R, o por lo menos, el código más comúnmente utilizado por los usuarios, ha sido meticulosamente supervisado por estadísticos y académicos de mucho prestigio de todo el mundo (el llamado “R Core team”).
  • Podemos además programar nuestras propios procedimientos y aplicaciones
  • En la misma página desde la que se puede bajar el programa, existe abundante documentación sobre cómo utilizarlo
  • Es gratuito.

Esta documentación forma parte del Curso de Doctorado de Análisis de Datos Ecológicos en R, que se impartirá entre el 5 y el 13 de mayo de 2009 en el Departamento de Ecología de la Universidad de Alcalá. La primera sesión constituye una introducción somera a R, en dónde el alumno podrá explorar los conceptos más básicos del manejo de esta herramienta, desde la instalación del software a la manipulación de datos y la programación de funciones sencillas. La documentación para esta primera sesión se puede descargar aquí. El índice de contenidos se detalla a continuación:
  1. ¿ Qué es R?
  2. ¿Cómo instalar R?
  3. CRAN y paquetes
  4. Tipos de objetos en R y la función str()
  5. Funciones y argumentos
  6. El menú de ayuda: Aprendiendo a ser autosuficientes
  7. Como leer datos en R
  8. Recomendaciones para organizar una sesión de trabajo
  9. Funciones básicas para la manipulación de datos
  10. Como repetir un procedimiento con el comando for()
  11. Citando R en los trabajos científicos
  12. Referencias

lunes, 13 de abril de 2009

Estima de los atributos funcionales de ecosistemas forestales a partir del NDVI

En 2006, Domingo Alcaraz publicaba un trabajo en Global Ecology and Biogeography sobre la caracterización de los tipos funcionales de ecosistemas de la Península Ibérica a partir del uso de índices de vegetación obtenidos de imágenes de satélite. Este trabajo consiste, básicamente, en definir unos atributos que permitan distinguir desde un punto funcional distintos tipos de ecosistemas. En su trabajo, Alcaraz et al. (2006) definen con estos fines tres atributos: la integral del índice de verdor NDVI (NDVI-I), la diferencia entre el máximo y el mínimo NDVI (RREL) y el mes de máximo NDVI (MMAX) (Figura 1). El primero se relaciona funcionalmente con la productividad de los ecosistemas, mientras que el segundo y el tercero están más vinculados a la estacionalidad. Utilizando estos tres atributos, calculados para una serie temporal larga obtenida a partir de imágenes NOAA/AVHRR a una resolución de 1 x 1 km, Alcaraz et al. (2006) clasificaron todos los tipos de ecosistemas en la Península Ibérica a lo largo de gradientes de productividad y estacionalidad.

Figura 1. Los tres atributos del NDVI empleados en la caracterización de los tipos funcionales de ecosistemas por Alcaraz et al. (2006): la integral del NDVI (NDVI.I), la diferencia entre el máximo y el mínimo NDVI (RREL) y el mes de máximo NDVI (MMAX). Ejemplo tomado de un pinar de Pinus nigra en la Sierra de los Filabres, Almería.

Partiendo de esta idea, pensamos en utilizar estos atributos, no para distinguir distintos tipos funcionales de ecosistemas, sino para intentar detectar el decaimiento en plantaciones de Pinus sylvestris y P. nigra en la Sierra de los Filabres. Para ello seguimos los siguientes pasos:
  1. Se promediaron en cada fecha los valores de las series temporales de NDVI obtenidos a partir de las imágenes MODIS (250 x 250 m). Esto supone unos 8 a 9 datos por fecha (2000-2008). Cómo estamos trabajando con valores promedios no hace falta interpolar los datos faltantes como hicimos anteriormente.
  2. Se calculó la integral total anual (NDVI-I) como la suma del área bajo la curva en cada intervalo.
  3. Se calculó la diferencia entre el mínimo y el máximo NDVI (RREL) y el mes del máximo NDVI (MMAX). Estos cálculos son bastante intuitivos así que no hace falta explicarlos (ver Figura 1).
Para los datos de Filabres, usamos el siguiente código en R, que permite calcular estos atributos a partir de los datos de las series temporales de NDVI obtenidos previamente de las imágenes MODIS, y obtener gráficos de dichos atributos para cada una de las 76 parcelas de estudio. Como en ocasiones anteriores, el código permite acceder directamente a los datos y ejecutar las distintas funciones diseñadas para tal fin.

load(url("http://archivos-para-subir.googlegroups.com/web/R+step+1.Rob?gda=GQGm0T4AAACkOfwqHAVd4YqgfIB09GDRsFEbW00qADF89i9HmLDg8UMvt5QpYAf3GGSwB4Eu5X3jsKXVs-X7bdXZc5buSfmx"))
ndvi <- points_ndvi@data[,-c(1:4)]
julday <- as.numeric(substr(labels(ndvi)[[2]], 5, 7)) year <- as.numeric(substr(labels(ndvi)[[2]], 1, 4))

NDVI.I <- rep(0,dim(ndvi)[1])
RREL <- rep(0, dim(ndvi)[1])
MMAX <- rep(0, dim(ndvi)[1])
for (i in 1:dim(ndvi)[1]) {
ndvi.mean <- tapply(as.numeric(ndvi[i, ]), as.factor(julday), mean)/10000
NDVI.i <- abs(ndvi.mean[1] - ndvi.mean[length(ndvi.mean)])/2 + min(c(ndvi.mean[1], ndvi.mean[length(ndvi.mean)]))
for (j in 1:(length(ndvi.mean)-1)) {
NDVI.i <- NDVI.i + abs(ndvi.mean[j+1] - ndvi.mean[j])/2 + min(c(ndvi.mean[j+1], ndvi.mean[j]))
}
NDVI.I[i] <- NDVI.i
RREL[i] <- max(ndvi.mean) - min(ndvi.mean)
MMAX[i] <- dimnames(ndvi.mean)[[1]][ndvi.mean == max(ndvi.mean)]
png(as.character(paste("NDVI.I", i, ".png", sep = ""), width = 550))
dimnames(ndvi.mean)[[1]] <- c("Jan", "", "Feb", "", "Mar", "", "Apr", "", "May", "", "Jun", "", "Jul", "", "Aug", "", "Sep", "", "Oct", "", "Nov", "", "Dec")
plot(c(0, ndvi.mean), type = "n", xlab = "", ylab = "NDVI", axes = F, ylim = c(0.1,0.9))
axis(1, at = c(2:24), labels = dimnames(ndvi.mean)[[1]], las = 2)
axis(2)
polygon(x = c(2:24, 24, 2), y = c(ndvi.mean, 0.1, 0.1), col = "grey90")
segments(x0 = 1, y0 = min(ndvi.mean), x1 = 1, y1 = max(ndvi.mean), col = "red", lwd = 3)
pos <- pmatch(ndvi.mean[ndvi.mean == max(ndvi.mean)], ndvi.mean)
segments(x0 = pos +1, y0 = 0.1, x1= pos+1, y1 = as.numeric(ndvi.mean[ndvi.mean == max(ndvi.mean)]), col = "blue", lwd = 3, lty = 3)
box()
title(main = paste(points_ndvi@data$Species[i], points_ndvi@data$Id[i]), sub =paste("Level of damage =", points_ndvi@data$Damage[i]))
text(x = 12, y = mean(c(max(ndvi.mean), 0.1)), labels = "NDVI.I", cex = 1.5)
text(x = 1.5, y =max(ndvi.mean) + 0.05, labels = "RREL", cex = 1.2)
text(x = ifelse(pos < y ="ifelse(pos" labels = "MMAX" cex =" 1.2)
dev.off()
}

El resultado de aplicar este código son tres vectores con los valores de estos atributos para cada una de las 76 parcelas (NDVI.I, RREL, MMAX) y 76 gráficas con la representación de dichos valores.

Estos valores se utilizarán junto a otras variables (p.e. insolación, altitud, etc) para intentar discriminar masas forestales con distintos grados de decaimiento. Pronto esperamos poder mostrar avances en esta línea, en la que llevamos trabajando ya varios meses.

Referencias:

Alcaraz, D., Paruelo, J.M. & Cabello, J. 2006. Identification of current ecosystem functional types in the Iberian Peninsula. Global Ecology and Biogeography 15: 200-212.

Una nueva revista de conservación tropical: Tropical Conservation Science

Hace poco descubrí una nueva revista de investigación centrada en la conservación de los ecosistemas tropicales: Tropical Conservation Science. Esta revista, que se viene publicando desde el año 2008 (tan sólo cuenta con cinco volúmenes publicados hasta el momento), contiene artículos muy interesantes sobre distintos aspectos relacionados con la conservación de los ecosistemas tropicales, como la política territorial, la conservación de comunidades y especies, o la conservación genética de poblaciones, entre otros muchos.

Tropical Conservation Science es una revista de libre acceso (open source), lo que favorece la diseminación de los trabajos científicos frente a revistas de suscripción. El Dr. Alejandro Estrada, un eminente primatólogo de la Universidad Nacional Autónoma de México, es
junto con Rhett Butler, de mongabay.com, uno de los dos editores ejecutivos de la revista. La revista cuenta además con una larga lista de editores asociados de la talla y calidad del Dr. Rodolfo Dirzo (Stanford University, USA), la Dra. Sarah Boyle (Arizona State University, USA) o el Dr. William F. Laurance (The Smithsonian Tropical Research Institute, Panama), entre otros.

Tropical Conservation Science ofrece, por tanto, una alternativa sólida y con futuro de proyección a revistas científicas ya consolidadas en el ámbito de la ecología tropical, como Biotropica o Journal of Tropical Ecology. Por todo ello, espero que esta breve reseña contribuya a su difusión, para que sea conocida por un público amplio y llegue a ser en el futuro una revista de referencia en el campo de la conservación de los ecosistemas tropicales.

Buscar entradas