Interpolar datos en series temporales o espaciales

by Nacho

Aquí esta el código y alguna referencia para saber que hacer con aquellas series de datos en las que encontramos gaps, y que necesitamos rellenar para poder trabajar con ellas. El paquete [zoo] nos va a permitir obterner interpolaciones tanto con series temporales como con datos espaciales.

 
Temp <- read.table("temperatura 2010 2011.csv", header=TRUE, sep=";", na.strings="")
Temp$Date<-as.Date(Temp$Date) #caracterizamos el vector como datos temporales
x$Date<-seq(min(Temp$Date),max(Temp$Date),1) #creamos un nuevo vector con todas las fechas donde nos faltan datos
Mix<-merge(Temp,x, by="Date",all= T) # unimos el vector y la matriz de datos por la columna común "Date", nos apareceran 
                                     # fechas con NAs 
Data<-Mix[,c("Date","Temp")]         # Seleccionamos las columnas de interes y creamos una nueva matriz
library(zoo)                         # Cargamos el paquete zoo que nos va a permitir realizar las interpolaciones
z<-zoo(Data$Temp,Data$Date)          # zoo(objeto, index)
t1<-na.approx(z, na.rm = FALSE)      # Interpolación (na.rm= F no elimina aquellos gaps donde no puede realizar interpolación)
plain <- as.data.frame(coredata(t1)) # Convertimos en un dataframe
colnames(plain)<-"Temp1"             # denominamos la columna 
Data1<-c(Data,plain)                 # Creamos una nueva dataframe con la base con Temp (con Na's) y Temp1 (sin Na's)
Data1<-as.data.frame(Data1)

Métodos de interpolación:

na.approx(object, along = index(object), na.rm = TRUE, ...) #realiza una interpolación lineal para rellenar los NAs 
na.spline # realiza una interpolación cúbica por splines para rellenar los NAs 
na.locf # reemplaza NAs con el valor previo que no era NA
na.omit # ("zoo" object) elimina las observaciones incompletas
na.trim # recorta series de Nas lejos del inicio y al final pero no en el interior de la serie 
na.stinterp # (package [stinepack]) que realiza la interpolación de Stineman

Fuente:

Zoo Quick reference
stackoverflow.com

ANOVAs en R

by InesNaya

A veces uno tiene que empezar por lo básico para entender lo más complejo (aunque normalmente muchas veces esto nos lo saltamos y luego hay que volver a los inicios). Aquí os dejo una guía general de tipos de ANOVAs y como hacerlas en R que están bien explicadas y que por lo menos a mi me sirvieron para aclarar ciertas cosas de lo que hace la función aov() o no en R (por ejemplo: ANOVAs con efectos aleatorios):

R and Analysis of Variance

que encontré traducida al español al menos en parte aquí.

A muchos les parecerá muy básico, pero repito, a veces es mejor recordar la base antes de ponerse a añadir complejidad a las cosas.

GLS intervalos de confianza

by Nacho

Para aquellos que han leido Zuur (2009) “Mixed Effects Models and Extensions in Ecology with R” y en el capitulo 4 nos hace un plot de los predict del modelo sin intervalos de confianza al 95% y nos comenta que todavía no se ha desarrollado un código para calcularlos (predict.gls), pero que con un poco de código engorroso similar al de las ecuaciones de modelo lineales se podría obtener fácilmente. Ahí es donde por lo menos yo me quedé mirando la numeración de las páginas por si le faltaba alguna hoja, pero no, no se molesto en calcularlos ;)

Aquí esta la librería “rms”, que permite con la función Gls calcular los predict y graficarlos. El código cambia un poco pero esta basado en la librería “nlme”. Hay que fijarse en como se formulan los weights.

library(rms)
 
b3$BroodC <- NULL ##Eliminar una columna (datadist da problemas si en alguna columna hay NAs
dd <- datadist(b3); options(datadist='dd')  #opción de longitud necesaria para el Gls(  Gls {rms})
Form <- formula(log(Media) ~ logLc+fBrood2)
Ident<-varIdent( ~1 | fBrood2) 
M2<- Gls(Form,data=b3,weights = Ident)
 
###Ploteamos los residuales estandarizados
 
E <- resid(M, type = "normalized")
Fit <- fitted(M)
op <- par(mfrow = c(1, 3))
plot(x = Fit, y = E,
xlab = "Fitted values", ylab = "Residuals",
main = "Residuals versus fitted values")
text(Fit,E,b2$Cod,cex=0.6)
hist(E, nclass = 15)
plot(E ~fBrood2, data = b3, ylab = "Normalised residuals")
text(as.numeric(b3$fBrood2),E,b3$Cod,cex=0.6) ##Para saber identificar las observaciones
par(op) 
 
####Plot de los predict values and 95% confidence intervals in Gls
 
plot(Predict(M2))

Leer archivos de excell y spss

by noelas

Si queréis leer directamente archivos de excel podéis utilizar una de las dos siguientes librerías:

A)

library(RODBC)
 
table.lt<-odbcConnectExcel("datos_excel.xls")

B)

library(xlsReadWrite)
 
table.lt<-read.xls("datos_excel.xls")

También es posible leer directamente datos de archivos de SPSS

library (foreign)
 
table.lt<-read.spss("datos_spss.sav", to.data.frame=T)

R studio, una nueva interfaz de trabajo en R.

by Nacho

Hay varias herramientas para trabajar con R, editores como TinnR, StatET, Emacs, Vim, Notepad++, NppToR, etc., RStudio ™ es ​​un nuevo entorno de desarrollo integrado de R. RStudio combina una interfaz de usuario intuitiva R con un editor de R. A pesar de estar en una version 0.94, permite editar y llevar un seguimiento de todos los datos (bases, resultados) y los graficos realizados. Además ayuda a completar el codigo. También hay un libro que describe su uso y configuración. Por otro lado la posibilidad de uso en multiplataforma (Windows, linux, MacOS y soporte WEB) va a hacer que RStudio se convierta pronto en una de las herramientas más usadas de R.

Hay herramientas interesantes no incluidas en la instalación RStudio, como un paquete de R “manipulate”, que van a permitir, metiendo un poco de código, la modificación de los ejes, colores objetos, etc.

Fuente RStudio

Boxplot with outliers label & sample size

by Nacho

Una herramienta útil para hacer un boxplot e identificar los outliers. La función del boxplot es la siguiente:

boxplot.with.outlier.label <- function(y, label_name, ..., spread_text = T, data, plot = T, range = 1.5, label.col = "blue", 
										jitter_if_duplicate = T, jitter_only_positive_duplicates = F)
{	
	# jitter_if_duplicate - will jitter (Actually just add a bit of numbers) so to be able to decide on which location to plot the label when having identical variables...
	require(plyr) # for is.formula and ddply
 
	# a function to jitter data in case of ties in Y's
	jitter.duplicate <- function(x, only_positive = F)
	{
		if(only_positive) {
			ss <- x > 0
		} else {
			ss <- T
		}	
		ss_dup <- duplicated(x[ss])
		# ss <- ss & ss_dup
		temp_length <- length(x[ss][ss_dup])	
		x[ss][ss_dup] <- x[ss][ss_dup] + seq(from = 0.00001, to = 0.00002, length.out = temp_length)
		x
	}
	# jitter.duplicate(c(1:5))
	# jitter.duplicate(c(1:5,5,2))
	# duplicated(jitter.duplicate(c(1:5,5,2)))
	# jitter.duplicate(c(0,0,1:5,5,2))
	# duplicated(jitter.duplicate(c(0,0,1:5,5,2)))
 
 
 
	# handle cases where 
	if(jitter_if_duplicate) {
		# warning("duplicate jutter of values in y is ON")
		if(!missing(data)) {	#e.g: we DO have data
			# if(exists("y") && is.formula(y)) {		# F && NULL # F & NULL
			y_name <- as.character(substitute(y))	# I could have also used as.list(match.call())
												# credit to Uwe Ligges and Marc Schwartz for the help
												# https://mail.google.com/mail/?shva=1#inbox/12dd7ca2f9bfbc39
			if(length(y_name) > 1) {	# then it is a formula (for example: "~", "y", "x"
				model_frame_y <- model.frame(y, data = data)
				temp_y <- model_frame_y[,1]
				temp_y  <- jitter.duplicate(temp_y, jitter_only_positive_duplicates)	# notice that the default of the function is to work only with positive values...
				# the_txt <- paste(names(model_frame_y)[1], "temp_y", sep = "<<-") # wrong...
				the_txt <- paste("data['",names(model_frame_y)[1],"'] <- temp_y", sep = "")				
				eval(parse(text = the_txt))	# jutter out y var so to be able to handle identical values.
			} else {	# this isn't a formula
				data[,y_name] <- jitter.duplicate(data[,y_name], jitter_only_positive_duplicates)
				y <- data[,y_name]	# this will make it possible for boxplot(y, data) to work later (since it is not supposed to work with data when it's not a formula, but now it does :))
			}		
		} else {	# there is no "data"		 
			if(is.formula(y)) { # if(exists("y") && is.formula(y)) {		# F && NULL # F & NULL
				temp_y <- model.frame(y)[,1]
				temp_y  <- jitter.duplicate(temp_y, jitter_only_positive_duplicates)	# notice that the default of the function is to work only with positive values...
				environment(y) <- new.env()
				assign(names(model.frame(y))[1], temp_y, environment(y))
					# Credit and thanks for doing this goes to Niels Richard Hansen (2 Jan 30, 2011)
					# http://r.789695.n4.nabble.com/environment-question-changing-variables-from-a-formula-through-model-frame-td3246608.html
				# warning("Your original variable (in the global environemnt) was just jittered.")	# maybe I should add a user input before doing this....
				# the_txt <- paste(names(model_frame_y)[1], "temp_y", sep = "<<-")
				# eval(parse(text = the_txt))	# jutter out y var so to be able to handle identical values.
			} else {
				y <- jitter.duplicate(y, jitter_only_positive_duplicates)
			}		
		}
	}
	# the_txt <- paste("print(",names(model_frame_y)[1], ")")
	# eval(parse(text = the_txt))	# jutter out y var so to be able to handle identical values.
	# print(ls())
 
 
	# y should be a formula of the type: y~x, y~a*b
	# or it could be simply y
	if(missing(data)) {
			boxdata <- boxplot(y, plot = plot,range = range ,...)
		} else {
			boxdata <- boxplot(y, plot = plot,data = data, range = range ,...)
		}
	if(length(boxdata$names) == 1 && boxdata$names =="") boxdata$names <- 1	# this is for cases of type: boxplot(y) (when there is no dependent group)
	if(length(boxdata$out) == 0 ) stop("No outliers detected for this boxplot")
 
	if(!missing(data)) attach(data)	# this might lead to problams I should check out for alternatives for using attach here...
 
 
	# creating a data.frame with information from the boxplot output about the outliers (location and group)
	boxdata_group_name <- factor(boxdata$group)
	levels(boxdata_group_name) <- boxdata$names[as.numeric(levels(boxdata_group_name))]	# the subseting is for cases where we have some sub groups with no outliers
	boxdata_outlier_df <- data.frame(group = boxdata_group_name, y = boxdata$out, x = boxdata$group)
 
 
	# Let's extract the x,y variables from the formula:
	if(is.formula(y))
	{
		model_frame_y <- model.frame(y)
		y <- model_frame_y[,1]
		x <- model_frame_y[,-1]
		if(!is.null(dim(x))) {	# then x is a matrix/data.frame of the type x1*x2*..and so on - and we should merge all the variations...
			x <- apply(x,1, paste, collapse = ".")
		}
	} else {
		# if(missing(x)) x <- rep(1, length(y))
		x <- rep(1, length(y))	# we do this in case y comes as a vector and without x
	}	
 
	# and put all the variables (x, y, and outlier label name) into one data.frame
	DATA <- data.frame(label_name, x ,y)
	if(!missing(data)) detach(data)	# we don't need to have "data" attached anymore.
 
	# let's only keep the rows with our outliers 
	boxplot.outlier.data <- function(xx, y_name = "y")
	{
		y <- xx[,y_name]
		boxplot_range <- range(boxplot.stats(y, coef = range )$stats)
		ss <- (y < boxplot_range[1]) | (y > boxplot_range[2])
		return(xx[ss,])	
	}
	outlier_df <-ddply(DATA, .(x), boxplot.outlier.data)
 
 
	# create propor x/y locations to handle over-laping dots...
	if(spread_text) {
		# credit: Greg Snow
		require(TeachingDemos)		
		temp_x <- boxdata_outlier_df[,"x"]
		temp_y1 <- boxdata_outlier_df[,"y"]
		temp_y2 <- temp_y1
		for(i in unique(temp_x))
		{
			tmp <- temp_x == i
			temp_y2[ tmp ] <- spread.labs( temp_y2[ tmp ], 1.3*strheight('A'), maxiter=6000, stepsize = 0.05) #, min=0 )
		}
 
	}
 
 
	# max(strwidth(c("asa", "a"))
	# move_text_right <- max(strwidth(outlier_df[,"label_name"]))	
 
	# plotting the outlier labels :)  (I wish there was a non-loop wise way for doing this)
	for(i in seq_len(dim(boxdata_outlier_df)[1]))
	{
		# ss <- (outlier_df[,"x"]  %in% boxdata_outlier_df[i,]$group) & (outlier_df[,"y"] %in% boxdata_outlier_df[i,]$y)
 
		# if(jitter_if_duplicate) {
			# ss <- (outlier_df[,"x"]  %in% boxdata_outlier_df[i,]$group) & closest.number(outlier_df[,"y"]  boxdata_outlier_df[i,]$y)
		# } else {
		ss <- (outlier_df[,"x"]  %in% boxdata_outlier_df[i,]$group) & (outlier_df[,"y"] %in% boxdata_outlier_df[i,]$y)
		# }
 
		current_label <- outlier_df[ss,"label_name"]
		temp_x <- boxdata_outlier_df[i,"x"]
		temp_y <- boxdata_outlier_df[i,"y"]		
		# cbind(boxdata_outlier_df,		temp_y2)
		# outlier_df
 
 
		if(spread_text) {
			temp_y_new <- temp_y2[i] # not ss			
			move_text_right <- strwidth(current_label)
			text( temp_x+move_text_right, temp_y_new, current_label, col = label.col)			
			# strwidth
			segments( temp_x+(move_text_right/6), temp_y, temp_x+(move_text_right*.47), temp_y_new )
		} else {
			text(temp_x, temp_y, current_label, pos = 4, col = label.col)
		}		
	}
 
	# outputing some of the information we collected
	list(boxdata = boxdata, boxdata_outlier_df = boxdata_outlier_df, outlier_df=outlier_df)
}

Un ejemplo:

boxplot.with.outlier.label (Lip~Egon,Cod, varwidth=T, ylab="Lipidos(ug/mg)", xlab="Estado gonada", data=P1_98gon, ylim=c(50,350))
sample.size <- tapply(P1_98gon$Lip, P1_98gon$Egon, length) ####Codigo  para determinar el tamaño muestral
ss.ch <- paste("N=", sample.size, sep="")            ####Codigo  pegar N= a los datos de tamaño muestral
mtext(ss.ch, at=1:length(unique(P1_98gon$Lip)), line=2, side=3)   ####Codigo  poner una linea de texto en el eje x en cada boxplot


Si queremos que no se superpongan los labels instalar el paquete {teachingDemos} y {plyr}

boxplot.with.outlier.label (Lip~Egon,Cod, varwidth=T, ylab="Lipidos(ug/mg)", xlab="Estado gonada", data=A99gon, ylim=c(50,350), spread_text = F) #añadimos spread_text=F

Fuente: R-Staistics blog

Ideas for cooking ggplot2

by Nacho

Ideas para ggplot2

#####
library(ggplot2)
 
data = structure(list(Trt = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L), .Label = c("A",
"B"), class = "factor"), Clone = structure(c(1L, 2L,
3L, 4L, 5L, 6L, 7L, 8L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L), .Label =
c("1",
"2", "3", "4", "5", "6", "7", "8",
"Mean"), class = "factor"), x = c(-27.5996, -27.5333, -27.0267,
-27.6467, -26.7667, -28.07, -27.4608, -28.1867, -29.3833, -28.2196,
-29.6567, -28.9608, -29.6167, -30.1892, -28.5633, -30.2208),
    x.SE = c(0.3603, 0.3085, 0.3085, 0.3085, 0.3085, 0.3085,
    0.4831, 0.3085, 0.3085, 0.3601, 0.3085, 0.3603, 0.3085, 0.2402,
    0.3085, 0.3686), y = c(98.12, 69.84, 78.47, 68.03, 58.2,
    33.39, 46.57, 65.75, 40.01, 38.23, 34.09, 44.37, 31.92, 39.85,
    34.37, 41.27), y.SE = c(15.32, 12.51, 12.51, 12.51, 12.51,
    12.51, 15.32, 12.51, 12.51, 15.32, 12.51, 15.32, 12.51, 12.51,
    12.51, 15.32)), .Names = c("Trt", "Clone", "x", "x.SE",
"y", "y.SE"), class = "data.frame", row.names = 3:18)
 
# x and ylim for geom_errorbar with standard error
fig.xlim = aes(xmin = x - x.SE, xmax = x + x.SE, height=0)
fig.ylim = aes(ymin=y-y.SE, ymax=y+y.SE, height=0)
 
fig1 = ggplot(data, aes(x, y, shape=factor(Trt), fill=factor(Clone),
size=factor(Clone)) ) +
scale_shape_manual(values=c(23,21), "Treatment") +
scale_fill_brewer(palette="Set1", "Clone") +
theme_bw() +
coord_cartesian(xlim= c(-31,-26), ylim=c(0,120) ) +
scale_y_continuous(limits=c(0,120)) +
scale_x_continuous(limits=c(-31,-26)) +
geom_errorbar(fig.ylim, width=0, size=0.3,colour="black") +
geom_errorbarh(fig.xlim, size=0.3, colour="black") +
geom_point(aes(size=factor(Clone))) +
scale_size_manual(values=c(2,2,4,2,4,4,2,2), "Clone") +
opts(panel.grid.minor=theme_line(colour = NA, size = 0.0),
panel.grid.major=theme_line(colour = NA, size = 0.0),
strip.background = theme_rect(col="black",fill=NA))
fig1

############################ In Grey scale
 
fig2 = ggplot(data, aes(x, y, shape=factor(Trt))) +
scale_shape_manual(values=c(23,21), "Treatment") +
theme_bw() +
coord_cartesian(xlim= c(-31,-26), ylim=c(0,120) ) +
scale_y_continuous(limits=c(0,120)) +
scale_x_continuous(limits=c(-31,-26)) +
geom_errorbar(fig.ylim, width=0, size=0.3,colour="black") +
geom_errorbarh(fig.xlim, size=0.3, colour="black") +
geom_point(size=6, fill="white") +
geom_text(aes(label=Clone),size=4) +
opts(panel.grid.minor=theme_line(colour = NA, size = 0.0),
panel.grid.major=theme_line(colour = NA, size = 0.0),
strip.background = theme_rect(col="black",fill=NA))
 
fig2

#################Grey scale2 y conservar  grid lines
 
fig3 = ggplot(data, aes(x, y, shape=factor(Trt), colour = factor(Trt))) +
  geom_errorbar(fig.ylim, width=0, size=0.3) +
  geom_errorbarh(fig.xlim, size=0.3) +
  geom_point(size=6, fill="white") +
  geom_text(aes(label=Clone),size=3, colour = "grey50") +
  scale_shape_manual(values=c(23,21), "Treatment") +
  scale_colour_grey("Treatment", end = 0.6) +
  theme_bw() 
fig3

Fuente: Ideas for converting a colour to a gray graph

Trabajar con nuestros datos

by Nacho

Reshape es un paquete de R para facilitar la restructuración y la agregación de los datos. Un herramienta básica para empezar a trabajar con nuestros datos.

Algunas cosas que podremos hacer: transformar la estructura de nuestra tabla de datos de una forma sencilla, crear tablas con medias, SD, SE, sumas… de las variables según queramos agruparlas

Para empezar a ver lo que podemos hacer podemos comenzar con el práctico paper “Reshaping data with the reshape package”, publicado en el Journal of statistical software

     time 	treatment 	subject 	rep 	potato    buttery 	grassy 	rancid 	painty 	 
61 	1 	    1 	            3 	       1.00      2.90      0.00 	 0.00 	0.00   	5.50 	 
25 	1 	    1 	            3 	       2.00 	14.00      0.00 	 0.00 	1.10 	0.00 	 
62 	1 	    1 	           10 	       1.00 	11.00      6.40 	 0.00 	0.00 	0.00 	 
 
 
 
 
R> ffm <- melt(french_fries, id = 1:4, na.rm = TRUE)
R> head(ffm)
          time       treatment subject rep     variable     value
1          1             1       3       1     potato        2.9
2          1             1       3       2     potato       14.0
3          1             1      10      1      potato       11.0
4          1             1      10      2      potato        9.9
5          1             1      15      1      potato        1.2
6          1             1      15      2      potato        8.8

Espero que sea útil, más información en la web Reshape

Ggplot2 graficas de barras de error

by Nacho

Empezar con ggplot2….

 
#introducir tabla de datos 
 
 
   z<- structure(list(Intermoult = structure (c(1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 13L, 14L, 15L, 16L, 5L, 6L, 7L,
8L, 9L, 10L, 11L, 12L, 13L, 14L, 15L),.Label = c("Instar1","Instar2","Instar3 ","Instar4","Instar5 ","Instar6","Instar7","Instar8",
"Instar9","Instar10","Instar11","Instar12","Instar13","Instar14" ,"Instar15 ","Instar16"),class = "factor"), 
Way = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L , 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L), 
.Label = c(" Initial 200 larvae", " 80 larvae transferred from the mass culture"), class = "factor"),
Mean = c(8.4,7.6,9.4,11.3,16.7,10.3,17.3,14.2,13.4,14.4,14.2,16.4,19.0,22.3,27.0,NA,16.7,17.0,18.3,15.0,15.6,15.1,14.7,17.6,21.1,18.8,27.0),
SD = c(0.5,1.3,1.0,1.6,2.4,1.7,4.0,3.5,2.1,4.1,2.5,2.2,3.5,3.9,0.0,NA,5.0,6.7,5.8,4.7,4.0,2.6,2.1,3.2,4.5,1.3,0.0),
upp = c(11,19,13,15,20,15,22,22,17,25,19,20,25,28,27,NA,24,33,30,28,26,20,19,23,30,20,27),
low = c(8,6,8,9,13,9,12,10,11,12,10,13,14,19,27,NA,12,10,11,8,11,11,11,13,16,17,27)) ,
class = "data.frame",  ,row.names = c(NA, -27L) )
colnames(z)<- c("Intermoult", "Way", "Mean", "SD", "upper", "lower")
 
###Crear plot
 
library(ggplot2)
 
pdf(file = "intermuda2.pdf", width = 6, height = 6, dpi= 600)    #Guardarlo como pdf
 
 
 
p<-ggplot(z,aes(x=Intermoult,y=Mean,ymin=Mean-SD,ymax=Mean+SD,
             groups=Way,colour=Way))+  
  geom_point(position=position_dodge(width=0.5))+ 
  geom_pointrange(width=0.5,position=position_dodge(width=0.5))+
   labs(x="Intermoult",y="Intermoult Period (Days)")+ theme_bw()+ opts(legend.position = c(0.18,0.6))+ #colocar la leyenda 
 scale_colour_manual("", c(" Initial 200 larvae" = "#5CB0E2", " 80 larvae transferred from the mass culture" = "darkblue"))
 
   p + ylim(c(-1, 30))
 
 
 
 
 #Características del plot#########################
 
last_plot() + opts(panel.grid.minor = theme_line(colour = NA),
panel.grid.major = theme_line(colour = NA),
plot.background =  theme_rect(colour = NA, fill = NA),
axis.title.x = theme_text( face = "bold"),
axis.title.y = theme_text( face = "bold", angle = 90),
legend.key = theme_rect ( colour = NA, fill = NA))
 
 
dev.off
 
#colocar graficos en una página (creamos un nuevo grid y seleccionamos la forma)
library(gridExtra)
grid.newpage()
print(grid.arrange(plot1, plot2, plot3,plot4, plot5, nrow=1, ncol=5))   #ponerlos en columnas y filas

Link plots sobre áreas geográficas

by Nacho

Aquí un link con código para hacer diversos plots de datos en áreas