Introducción:

Siendo la estructura del dosel del bosque una capa que determina las variaciones en los regímenes lumínicos del sotobosque, y que estas dinámicas cambiantes pueden ser capturadas a través de fotografías, podríamos pensar en que es factible que mediante el análisis de imágenes se pueda monitorear dinámicas temporales y espaciales de los ecosistemas forestales.

En este caso la utilización de índices de MIG y Anisotropía se hacen de manera experimental para la exploración de posibles patrones asociados a la estructura del dosel y su heterogeneidad.

Por otro lado, nuevos acercamientos de mediciones indirectas del dosel son metodologías que pueden ser probadas en diferentes ecosistemas y así lograr observar qué tan sensibles son estos métodos para detectar cambios en la estructura dentro de un mismo ecosistema.

¿Qué nos dicen los índices?

El MIG es una medida de la heterogeneidad y representa la aleatoriedad espacial de patrones en la imagen.

El MIG va de 0 a 1 donde valores cercanos a cero están asociados a patrones estructurados y valores cercanos a 1 están asociados con patrones aleatorios.

La Anisotropía es una medida que determina si los patrones de MIG son más heterogeneos a lo largo de una dirección (Horizontal o Vertical)

Valores bajos indican heterogeneidad hacia dirección vertical y valores altos indican heterogeneidad hacia dirección horizontal.

Paquetes a utilizar:

Se necesitan de varios paquetes para realizar el análisis. En caso de no tenerlos, instalar los siguientes:

  install.packages("raster")
  install.packages("rgdal")
  install.packages("devtools")
  library("devtools") #para el siguiente paso es necesario cargarlo 
  devtools::install_github("cmartin/LAI")
  devtools::install_github("cmartin/EXIFr")
  devtools::install_bitbucket("persican/imagemetrics", ref="36d34cc")

Recordar siempre cargar las librerías:

library(raster)
library(imagemetrics)
library(LAI)
library(EXIFr)

Indicar directorio donde tenemos los archivos

Es necesario indicarle a R dónde tenemos los archivos, para que pueda buscar y analizar las imágenes que hemos dentro de la carpeta. Así mismo, es bueno indicar el directorio de trabajo que utilizaremos durante la sesión:

#Directorio de trabajo:
setwd("C:camino/aldirectorio/detrabajo")
#Dirección a la carpeta con imágenes
images.folder = "C:/camino/folder/con/imagenes" 
#Archivo de salida:
results.file = "nombre.csv"

En ocasiones, el análisis toma mucho tiempo (depende de la cantidad de imágenes y el tamaño de cada una) por lo que es bueno tener una función ifelse que nos permite iniciar donde se cerró y se acabó una sesión. Esto para no tener que realizar el análisis desde cero cada vez que comenzamos una sesión.

if (file.exists(results.file)) {
        existing = read.csv(results.file)
        start = max(existing$i)+1
} else {
        start = 1
}

Índices MIG-Anisotropía-GAP Fraction

Lo siguiente es una función que a partir de los pixeles de cada una de las imágenes obtiene diferentes índices de acuerdo al arreglo que estos tengan.

files = dir(images.folder)
nbins = 15 # see MIG algorithm
nb.photos = length(files)

if (start <= nb.photos) {
        
        for (i in start:nb.photos) {
                file.to.analyze = files[[i]]
                path = paste(images.folder,file.to.analyze,sep="/")
                
                #índice GAP
                B <- raster(path, band = 3)
                binary_image <- unimodal_threshold(B)
                a = gap_fraction(binary_image)
                
                # Load 3 rasters from the target image, one for red, one for green and one for the blue channel
                R = raster(path, band=1)
                G = raster(path, band=2)
                B = raster(path, band=3)
                
                # Combine the RGB channels to create a grayscale image
                RGB = brick(R,G,B) 
                r.grey = mean(RGB)
                    
                # On the four channels, find either right (1), diagonal (2) or below (3) neighbors for histogram calculations
                #     
                v.grey.1 = getImagePixels(r.grey, side = 1)
                v.grey.2 = getImagePixels(r.grey, side = 2)
                v.grey.3 = getImagePixels(r.grey, side = 3)
                    
                # Calculate histograms from neighbor vectors
                prob.grey.1 = calculateHisto(reference_vector = v.grey.1$reference_vector,
                                             neighbour_vector = v.grey.1$neighbour_vector, nbins = nbins)
                   
                prob.grey.2 = calculateHisto(reference_vector = v.grey.2$reference_vector,
                                             neighbour_vector = v.grey.2$neighbour_vector, nbins = nbins)
                
                prob.grey.3 = calculateHisto(reference_vector = v.grey.3$reference_vector,
                                             neighbour_vector = v.grey.3$neighbour_vector, nbins = nbins)
                
               # Escribir archivo csv con los índices:
                write.table(
                        data.frame(
                                i = i,
                                ID = file.to.analyze,
                                MIG.grey = meanInformationGain(prob.grey.2),
                                Aniso.grey = meanInformationGain(prob.grey.1) / meanInformationGain(prob.grey.3),
                                GAP = a,
                                DateAndTime = read_exif_tags(path)[["DateTime"]]
                                
                        ),
                        file=results.file,
                        append=i!=1,
                        col.names = i==1,
                        row.names = FALSE,
                        sep=","
                )   
        }
}

Una vez que el proceso de detiene, debemos de encontrar un archivo .csv en el directorio de trabajo, que contiene los valores de los índices por imagen digital.

Preparación para el Análisis de los datos:

Ya con el archivo .csv con todos los índices de cada una de als imágenes que necesitamos analizar, podemos proceder a realizar un análisis exploratorio de los datos. Antes, debemos de arreglar un poco el excel:

Traer los datos a la sesión:

data <- read.csv("indices_pv.csv",h=T,stringsAsFactors = F)
#Revisemos el encabezado y las primeros 4 observaciones:
head(data,4)

##   i           ID  MIG.grey Aniso.grey       GAP         DateAndTime
## 1 1 IMG_6545.JPG 0.4948146  0.9986244 0.3288093 2016:09:30 09:05:58
## 2 2 IMG_6651.JPG 0.1635735  0.9644914 0.1687209 2016:09:30 16:30:45
## 3 3 IMG_6768.JPG 0.3716458  0.9764172 0.1893754 2016:10:01 10:11:44
## 4 4      E03.jpg 0.2157913  1.0359106 0.1305531 2017:05:13 09:41:04

Lo que observamos es que la columna DateAndTime tiene dos tipos de datos en una sola columna, por lo que sería recomendable separarla en una que contenga Date y otra que contenga Time

library(tidyr)
data<-separate(data,col=DateAndTime,into=c("Date","Time"),sep=" ")
head(data)

##   i           ID  MIG.grey Aniso.grey       GAP       Date     Time
## 1 1 IMG_6545.JPG 0.4948146  0.9986244 0.3288093 2016:09:30 09:05:58
## 2 2 IMG_6651.JPG 0.1635735  0.9644914 0.1687209 2016:09:30 16:30:45
## 3 3 IMG_6768.JPG 0.3716458  0.9764172 0.1893754 2016:10:01 10:11:44
## 4 4      E03.jpg 0.2157913  1.0359106 0.1305531 2017:05:13 09:41:04
## 5 5      E04.jpg 0.3363079  0.9896186 0.1852430 2017:05:13 09:51:27
## 6 6      E05.jpg 0.3648285  0.9854371 0.1987676 2017:05:13 10:01:11

Una vez separadas las columnas, tenemos de indicarle a R cómo debe de interpretar las fechas dándole un formato y almacenándolo en una nueva columna:

date <- data$Date
RDate <- strptime(as.character(date),"%Y:%m:%d")
data <- data.frame(data,RDate) #sobreescribimos objeto
data <- data[,-6] #Eliminamos columna innecesaria
head(data)

##   i           ID  MIG.grey Aniso.grey       GAP     Time      RDate
## 1 1 IMG_6545.JPG 0.4948146  0.9986244 0.3288093 09:05:58 2016-09-30
## 2 2 IMG_6651.JPG 0.1635735  0.9644914 0.1687209 16:30:45 2016-09-30
## 3 3 IMG_6768.JPG 0.3716458  0.9764172 0.1893754 10:11:44 2016-10-01
## 4 4      E03.jpg 0.2157913  1.0359106 0.1305531 09:41:04 2017-05-13
## 5 5      E04.jpg 0.3363079  0.9896186 0.1852430 09:51:27 2017-05-13
## 6 6      E05.jpg 0.3648285  0.9854371 0.1987676 10:01:11 2017-05-13

Análisis Imágenes Digitales Planas

Laboratorio de Ecología Funcional y Ecosistemas Tropicales. Por Ronny Hernández Mora

18 de mayo de 2017

¿Qué es esto?