En este ejercicio vamos a trabajar con las estructura de datos ‘array’ y ‘list’.
A partir de un fichero NetCDF ‘tas_obs_mm.nc’, que contiene las temperaturas medias mensuales observadas en una rejilla de 0.5ºx0.5º (base de datos observacionales E-OBS), vamos a extraer datos del fichero para manipularlos.
Borramos todos los objetos almacenados en nuestro entorno de trabajo
Para ver los objetos que tenemos cargados en memoria usamos ls().
Vamos a borrar todos los objetos antes de hacer nuestro ejercicio:
rm(list = ls())
ls()
Leemos el fichero netdcf con los datos de temperatura
Para manipular ficheros NetCDF primero tenemos que cargar el paquete ‘RNetCDF’, que ya está instalado:
library(RNetCDF)
Ahora ya podemos usar las funciones del paquete para trabajar con los datos del fichero:
fichero_obs <- paste0(getwd(),"/data/tas_obs_mm.nc")
ncfile_obs <- open.nc(fichero_obs) # abrimos la conexión con el fichero
class(ncfile_obs)
lobs <- read.nc(ncfile_obs) # leemos el fichero y devuelve un objeto de la clase list
class(lobs) # clase a la que pertenece el objeto lobs
str(lobs) # estructura del objeto lobs
print.nc(ncfile_obs)
file.inq.nc(ncfile_obs)
var.inq.nc(ncfile_obs, 'tas')
Extraemos los datos de temperatura:
tas_obs <- lobs$tas # devuelve un objeto de la clase array
class(tas_obs) # clase del objeto
dim(tas_obs) # dimensiones del objeto
typeof(tas_obs) # tipo de almacenamento interno
str(tas_obs) # estructura del objeto
Extraemos el array de longitudes:
lon_obs <- lobs$longitude #
class(lon_obs)
typeof(lon_obs)
dim(lon_obs)
Extraemos el array de latitudes:
lat_obs <- lobs$latitude
class(lat_obs)
typeof(lat_obs)
dim(lat_obs)
Extraemos el array de tiempo:
time_obs <- lobs$time #
class(time_obs)
typeof(time_obs)
dim(time_obs)
Leyendo sólo las variables que nos interesen del fichero:
No hace falta cargar todos los datos del fichero. Usando la función ‘var.get.nc’ extraemos el array del tiempo, por ejemplo:
var.get.nc(ncfile_obs, "time", start=NA, count=NA, na.mode=0, collapse=TRUE)
Añadimos el nombre a las dimensiones del array de temperaturas:
De esta forma tenemos bien identificado lo que representa cada dimensión cuando trabajemos con el array.
Se puede necesitar el nombre de las dimensiones para comprobar que existe dicha dimensión
names(dim(tas_obs)) <- c("lon", "lat", "time")
dim(tas_obs)
time_dim_name <- "time"
if(!(time_dim_name %in% names(dim(tas_obs)))) {
stop("Parameter 'tas_obs' must not have temporal dimension.")
}else {
print("Existe dimensión tenporal. ok.")
}
Se pueden utilizar los nombres de las dimensiones para reorganizar las dimensiones de los objetos array cuyas dimensiones no estén ordenadas igual.
Vamos a cambiar el orden de las dimensiones del array tas_obs2,
para que coincida con el orden de las dimensiones del array tas_obs, y
poder operar con ellos:
dim(tas_obs)
tas_obs2 <- rnorm(835*21*31, 15, 20)
dim(tas_obs2) <- c(time = 835, lon = 31, lat = 21)
dim(tas_obs2)
tas_obs2[5,15,10]
pos <- match(names(dim(tas_obs)), names(dim(tas_obs2)))
pos
tas_obs2 <- aperm(tas_obs2,pos)
dim(tas_obs2)
tas_obs2[15,10,5]
tas_sum <- tas_obs + tas_obs2
Clases ‘POSIXt’ y ‘Date’
Vamos a convertir el array ‘time_obs’ de tipo ‘character’ en un objeto de tipo ‘POSIXt’ para representar fechas y horas:
time_obs[1:24]
unitTime <- att.get.nc(ncfile_obs, "time", "units")
timeut_obs <- utcal.nc(unitTime, time_obs, type="c")
class(timeut_obs)
timeut_obs[1:24]
Paquete ‘lubridate’
Para facilitar el manejo de fechas y horas, existe el paquete ‘lubridate’. Facilita la manipulación de este tipo de datos.
library(lubridate)
paste0(year(timeut_obs), "-", month(timeut_obs, label = TRUE, abbr = TRUE))[1:24]
En el enlace 'https://estadistica-dma.ulpgc.es/cursoR4ULPGC/6h-Fechas.html' tenéis más información para practicar con fechas y horas.
Extrayendo datos de un punto de grid determinado
Extraemos la serie temporal en un punto de rejilla (por ej: [lon=1ºE, lat=42ºN]):
tas_ptoLonLat <- round(tas_obs[which(abs(lon_obs - 1.0) < 1e-10), which(abs(lat_obs - 42.0) < 1e-10), ],2)
tas_ptoLonLat[1:24]
Series temporales
R posee un tipo de clase para representar series temporales, ts dentro del paquete base. Podemos crear un objeto tipo ‘ts’ de R para representar nuestra serie:
ts_tas_ptoLonLat <- ts(tas_ptoLonLat, timeut_obs[1], timeut_obs[length(timeut_obs)], 1/12)
class(ts_tas_ptoLonLat)
attributes(ts_tas_ptoLonLat)
Para eliminar bucles 'for' en R existen varios funciones que hacen que los cálculos sean más rápidos.
tasmean_lonlat <- apply(tas_obs, c(1,2), mean)
tasmean_lonlat
tasmean_date <- apply(tas_obs,3, mean, na.rm=TRUE)
v_tas_obs <- as.vector(tas_obs)
a_tas_obs <- array(v_tas_obs, dim=c(lon=31, lat=21, month=12, year=70))
a_tas_mean_mm <- apply(a_tas_obs, c(1,2,3), mean)
Función ‘lapply’
Para trabajar con listas y devolver una lista
lapply(lobs, max, na.rm = TRUE)
Función ‘sapply’
Si queremos que devuelva un array en vez de una lista:
sapply(lobs, min, na.rm = TRUE)
Función ‘multiApply::Apply’
La función Apply del paquete multiApply es muy útil para tratar con varios arrays.
library(multiApply)
A <- array(1:30, c(5, 2, 2))
B <- array(1:30, c(5, 2, 2))
D <- Apply(data = list(A, B),
target_dims = c(2, 3),
fun = "%*%")$output1
D