“Drivers” del cambio de uso.

1. Cambio de coberturas y uso del suelo

La distribución actual de la cobertura terrestre en Colombia es un producto de la huella humana histórica explicada por múltiples factores y factores que desempeñaron diferentes roles en la transformación del ecosistema naturales (Etter et al., 2008). Desde el año 1700, la expansión agrícola seguida del aumento de la población después de 1900 redujo un tercio del bosque natural colombiano, donde el pastoreo de ganado es uno de los principales impulsores del cambio forestal durante este período; donde la industrialización, la urbanización y la globalización promovieron la expansión de la frontera agrícola en todo el país (Etter et al., 2008). El cultivo de coca comúnmente discutido como un impulsor de la deforestación, se propone solo como una actividad asociada a la agrícola (Dávalos et al. 2016), donde otros usos legales son más propensos a explicar las tasas de deforestación que el cultivo de coca (Dávalos et al. UNODC y MINAM 2011, Armenteras et al., 2013 y Chadid et al., 2015).

La deforestación en colombia es un problema vigente y muy frecuentemente se dicute a cerca de las dinamicas de cambio de estos ecosistemas en lo que se han enfocado multiples estudios (Armenteras et al. 2013, Sánchez-Cuervo et al. (2012), (???) 2011, (???) 2008). Sin embargo, se pocos estudios se han enfocado en las dinaámicas de cambio de cambio de ecosistemas diferentes a los bosques y las causas y agentes asociados a estos cambios.

Por lo general diferentes fuentes no oficiales varían en sus estimaciones sobre la misma (Armenteras et al. 2013, Sánchez-Cuervo et al. (2012)), la tendencia es preocupante, debido a los efectos directos que tiene la deforestación sobre la pérdida de biodiversidad (Brooks et al. 2002 Jha and Bawa (2006)).

Para poder actuar efectivamente frente al fenómeno de la deforestación y degradacioón/peérdida de otros ecosistemas, es necesario entender sus causas. (Geist and Lambin 2002) hacen una revisión de la literatura sobre causas y agentes de la deforestación en el trópico y proponen un marco conceptual para su análisis. El marco de análisis subdivide las causas en próximas y subyacentes. Las causas próximas son las actividades humanas que afectan directamente el bosque, como por ejemplo: la agricultura, la tala, la minería, y la expansión de infraestructura, mientras que las causas subyacentes se entienden como procesos fundamentales que afectan las causas próximas y que pueden ser de origen social, político, económico, y cultural (Figura 1).

Figura 1. Patrones de causalidad de las deforestación en el tropico propuesto por Geist & Lambin 2001. Tomado de Geist & Lambin 2001.

En cuanto a deforestación en Colombia las causas próximas varían por región y han sido descritas por (Armenteras et al. 2013) en las subregiones: piedemonte, norte, nororiente y sur. El piedemonte es la zona de más desarrollo económico de la región y concentra la mayor parte de la deforestación debida principalmente que es debida a la colonización y a actividades de ganadería, minería y cultivos ilícitos. La subregión norte carece de alternativas económicas viables y allí las zonas de reserva forestal han sido ocupadas con la expectativa de sustracción y titulación, fortalecidas por la perspectiva de yacimientos minero-energéticos. La subregión nororiental se caracteriza por su alto nivel de conservación de bosque que es amenazado por la reciente minería ilegal de oro y coltán, así como por el descubrimiento de yacimientos mineros dentro de los resguardos indígenas. La zona sur presenta bajos niveles de deforestación que son ocasionados por la minería ilegal y la extracción selectiva de maderas de alto valor que son traficadas ilegalmente.

2. Cambio de uso en Orinoquia

Estudios enfocados a entender las causas economicas historicas asociadas a los cambios de usos y coberturas sugieren que las actividades economicas han variado a lo largo de los diferentes periodos desde las actividades de caza/colecta (prehispanico) hasta la intesificación de ganadería en pasturas y agricultura a gran escala actual (Etter et al. 2010).

Romero et al. 2011 encontro que en la orinoquia colombiana entre el periodo de 1987 y 2007 el área experimento algún cambio de coberturas o uso de la tierra donde la mayoria de los cambios ocurrieron en la última decada. Una fracción importante de los cambios de coberturas/usos estaba asociada a la expansión de cultivos de palma.

Figura 2.Principales actividades economicas asociadas a la región de la Orinoquía desde la epoca pre hispanica hasta el 2000. Tomado de Etter 2010.

Para este análisis se siguió el marco conceptual propuesto por Geist and Lambin (2002), agrupando las variables como causas próximas y subyacentes a nivel municipal y por grupos como factores poblacionales, económicos, de infraestructura y políticos.

library (sp)  # classes for spatial data
library (raster)  # grids, rasters
library (rasterVis)
library (rasterVis)  # raster visualisation
library (maptools) # map simple operations
library (rgeos) # coord transform
library (rgdal) # coord transform
library (ggmap) # visualization
library (tidyverse) # date managemgemt
library (mapview) # Visualization widget in esri sat image and open street map
library (readxl)# read excel data
library (sf) # shp and vector map handling 
library (gfcanalysis) # Hansen et al. 2013 Global Forest Change dataset
library (rgdal)# used to read/write shapefiles and rasters
library (dplyr)
library (printr)
library (corrplot) # correlations
library (knitr) # to convert code in web pabe 
library (kableExtra) # nice tables
library (printr) # table do not work
library (DT) # show table by part 
library (xlsx)

Datos

Se recopilo una tabla de datos con 1200 variables (en la filas), para cada municipio de Meta (en las columnas). La base de datos se ensamblo de diversas fuentes como; DANE, DNP, Ministerio de Agricultura, IGAC, etc. Para ver una descripción detallada de las variables por favor consulte ESTE ARCHIVO y para ver la base de datos recopilada consulte ESTE ENLACE.

Estas variables fueron agrupadas en grandes grupos relacionados con el tipo de actividad agrícola, el tipo de uso y variables demográficas. Para agrupar las variables se realizó un análisis de correlación simple para seleccionar una o dos variables de cada grupo correlacionado.

##################################
### Read Data set from excel
##################################

# Specify sheet with a number or name
info_meta <- read_excel ("C:/Users/diego.lizcano/Box Sync/CodigoR/Cambio_Meta/data/NASCA000228-2017_ProductoB_municipal.xlsx", sheet = "BASE_DATOS")

meta_datos <- read_excel ("C:/Users/diego.lizcano/Box Sync/CodigoR/Cambio_Meta/data/NASCA000228-2017_ProductoB_municipal.xlsx", sheet = "VARIABLES")

x <- unique(meta_datos[4])
# knitr::kable(x, digits = 2, caption = "Tipos de variables y grupos")

#### show table by part
datatable(x, options = list(pageLength = 15))

Agrupando variables correlacionadas

Variables agropecuarias

# Correlation matrix from info_meta
# with mpg, cyl, and disp as rows 
# and hp, drat, and wt as columns 

########################################
# function to plot correlation #
########################################

f_plot_correla <- function(dataset=info_meta, inicio, final) {
  x1 <- dataset[inicio:final]
  y1 <- dataset[inicio:final]
  
  corMat1 <- cor(x1, y1)

  corrplot(corMat1, order = "FPC",type = "lower", tl.pos = "lower", tl.cex = 0.5 )#, title= titulo)
}

f_plot_correla(dataset = info_meta, 151, 200 )

Variables agropecuaria-tendencia

f_plot_correla(dataset = info_meta, 201, 230)

Variables agropecuaria-uso

f_plot_correla(dataset = info_meta, 231, 255)

Variables agropecuaria-vivienda

f_plot_correla(dataset = info_meta, 262, 285)

Variables agropecuaria-maquinaria

f_plot_correla(dataset = info_meta, 292, 333 )

Variables agropecuaria-asistencia tecnica

f_plot_correla(dataset = info_meta, 361, 396 )

Variables agropecuaria-agua

# f_plot_correla(dataset = info_meta, 397,465)
x8 <- info_meta[397:465]
y8 <- info_meta[397:465]
corMat8 <- cor(x8, y8)

# corrplot(corMat, method = "color", type="lower")#,
         # col = colorRampPalette(c("blue", "white", "red"))(100))

corrplot(corMat8, order = "original",type = "lower", tl.pos = "full", tl.cex = 0.5)

Variables agropecuaria-suelo

f_plot_correla(dataset = info_meta, 466,495)

Variables agropecuaria-ganado

f_plot_correla(dataset = info_meta, 541,594)

Análisis

Deforestación en Meta

Se usaron los datos de deforestación de Hansen et al. (2013), Global Forest Change Dataset, los cuales de descargaron del servidor de Google usando el paquete gfcanalysis (Zvoleff 2015).

Año deforestado desde 2000 a 2014

Col_full<-getData('GADM', country='COL', level=2)


###############################################################################
# Download data from Google server for a given AOI
###############################################################################
col_full_sf <- st_as_sf(Col_full) # convert to sf object

meta_municip <- filter(col_full_sf, NAME_1 == "Meta")
aoi.full <- as(meta_municip, "Spatial")

aoi <- as(meta_municip, "Spatial")

# plot(aoi, lty=2, col="#00ff0050") # full map

# plot(gfc_data[[4]], main="year of forest loss, 2000 to 2014")
# plot(aoi, add=TRUE, lty=2, col="#00ff0050")

path_gfc_data <- "C:/Users/diego.lizcano/Box Sync/CodigoR/Cambio_Meta/meta/meta_extract_thresholded.tif"

gfc_data <- raster(path_gfc_data, band=2) 

# color palete
library(colorspace)
myTheme <- rasterTheme(region=sequential_hcl(14, power=-1))


myTheme <- rasterTheme(region =c("#FFFBFC", "#C1545A" ))

# levelplot(Aug, par.settings = myTheme, contour = TRUE)
detach(package:ggmap, unload=TRUE)
detach(package:ggplot2, unload=TRUE)
# plot as rasterVis
EN_map <- levelplot(gfc_data, margin = list(FUN = 'median'), contour=FALSE, par.settings = myTheme)#BuRdTheme)
EN_map + layer(sp.lines(aoi, lwd=0.8, col='black'))

# plot(gfc_data, main="year deforested 2000 to 2014", col=rainbow)
# plot (aoi, add=T)
###############################################################################
# Performing thresholding and calculate basic statistics
###############################################################################
# 
# # Calculate and save a thresholded version of the GFC product
# gfc_thresholded <- threshold_gfc(gfc_data,
#                                  forest_threshold=forest_threshold, 
#                                  filename="C:/Users/diego.lizcano/Box Sync/CodigoR/Cambio_Meta/meta/meta_extract_thresholded.tif",
#                                  overwrite=TRUE)
# 
# # see http://azvoleff.com/articles/analyzing-forest-change-with-gfcanalysis/
# # for threshold details
# 
# # change gfc_thresholded to UTM to get meters
# 
# 
# # Calculate annual statistics on forest loss/gain
# 
# ###############################################################################
# # Make visualization of forest change
# ###############################################################################
# 
# # Calculate and save a thresholded annual layer stack from the GFC product 
# # (useful for simple visualizations, etc.)
# gfc_thresholded_annual <- annual_stack(gfc_thresholded)
# writeRaster(gfc_thresholded_annual, filename='C:/Users/diego.lizcano/Box Sync/CodigoR/Cambio_Meta/meta/meta_thresholded_annual.tif', overwrite=TRUE)
# 
# # Save a simple visualization of the thresholded annual layer stack (this is 
# # just an example, and is using the data in WGS84. The data should be projected 
# # for this).
# animate_annual(aoi, gfc_thresholded_annual)
# 
# # Stop the parallel processing cluster
#  if (require(snow)) endCluster()
# 
# # my.at <- seq(0, 6, 1)
# # 
# # myColorkey <- list(at=my.at, ## where the colors change
# #                    labels=list(
# #                      at=my.at ## where to print labels
# #                    ))
# # levelplot(r, at=my.at, colorkey=myColorkey)
# 
# size_scale = 1
# maxpixels = 50000
# 
# library(plyr)
# pyr<-list()
# aoi_tr <- spTransform(aoi, CRS=CRS("+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"))
# aoi_points <- fortify(aoi_tr, region = "id")
# aoi_df <- join(aoi_points, aoi@data, by = "id")
# 
# # paco<- ggplot(aoi_df) + 
# #   aes(long,lat,group=group) + 
# #   geom_polygon(alpha = 0.1) +
# #   geom_path(linetype = label, alpha = 0.7) +
# #   coord_equal() 
# 
# for (i in 1:13){
# p<-gplot(gfc_thresholded_annual[[7]], maxpixels = 50000) + 
#   geom_tile(aes(fill = factor(value, levels = c(1, 2, 3, 4, 5, 6, 0)))) + coord_fixed() + 
#   scale_fill_manual("Cover", breaks = c("1", "2", "3", 
#                                         "4", "5", "6", "0"), labels = c("Forest", "Non-forest", 
#                                                                         "Forest loss", "Forest gain", "Loss and gain", "Water", 
#                                                                         "No data"), values = c("#008000", "#ffffb3", "#ff0000", 
#                                                                                                "#0000ff", "#ff00ff", "#c0c0c0", "#101010"), drop = FALSE)  
# 
# 
# mp<-p + geom_path(data = aoi_df, aes(x = long, y = lat, group = group, alpha = 0.05) ) 
# 
# pyr[[i]] <- mp
#   }
# 
# 
# multiplot(multiplot(pyr[[1]], pyr[[2]], pyr[[3]], pyr[[4]],pyr[[5]], pyr[[6]],
#                     cols=3))

El área deforestada y el área de bosque se cuantifico para cada municipio

def_muni <- read.xlsx("C:/Users/diego.lizcano/Box Sync/CodigoR/Cambio_Meta/data/municipi_table.xlsx", sheetName = "Sheet1", header = TRUE)


ggplot(meta_municip) +
        geom_sf(aes(fill = NAME_2))

##### show table by tandas
datatable(def_muni, options = list(pageLength = 15))
# kable(beetle_sp_sitio3)
# kable(def_muni2, digits = 2, caption = "Deforestación por municipios") %>%
#   kable_styling(bootstrap_options = c("striped", "hover"))

Regresion: GLMM nested group effect by year

library (sjPlot)
library (lme4)
library (nlme)
library (arm)
             
### And finally, we can fit nested group effect terms through the following
# # varying-intercept group effect using the variable year
mod_lmer1<-lmer(cover ~ poblacion + coca + (1|year), data=def_muni)

# plot probability curve of fixed effects
sjp.lmer(mod_lmer1, type = "fe.pc")

# plot qq-plot of random effects
# sjp.lmer(mod_lmer1, type = "re.qq")

Regresion: GLMM nested group effect by municipio

### And finally, we can fit nested group effect terms through the following
# # varying-intercept group effect using the variable year
mod_lmer2<-lmer(cover ~ poblacion + coca + (1|municipio), data=def_muni)

# plot probability curve of fixed effects
sjp.lmer(mod_lmer2, type = "fe.pc")

# plot qq-plot of random effects
sjp.lmer(mod_lmer2, type = "re.qq")

# plot probability curves for each covariate
# grouped by random intercepts
sjp.lmer(mod_lmer2,
          type = "ri.pc",
          show.se = TRUE)

# plot probability curves for each covariate
# grouped by random intercepts
sjp.lmer(mod_lmer2,
          type = "ri.pc",
          facet.grid = FALSE)

# sort all predictors
sjp.lmer(mod_lmer2,
         facet.grid = FALSE,
         sort.est = "sort.all",
         y.offset = .4)

# plot effects
sjp.lmer(mod_lmer2, type = "eff")

Info R session

print(sessionInfo(), locale = FALSE)
## R version 3.4.0 (2017-04-21)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows 10 x64 (build 15063)
## 
## Matrix products: default
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] arm_1.9-3             MASS_7.3-47           nlme_3.1-131         
##  [4] lme4_1.1-15           Matrix_1.2-12         sjPlot_2.4.0         
##  [7] xlsx_0.5.7            xlsxjars_0.6.1        rJava_0.9-9          
## [10] DT_0.2                kableExtra_0.6.1.9000 knitr_1.18           
## [13] corrplot_0.84         printr_0.1            gfcanalysis_1.4      
## [16] sf_0.6-0              readxl_1.0.0          mapview_2.2.0        
## [19] leaflet_1.1.0         forcats_0.2.0         stringr_1.2.0        
## [22] dplyr_0.7.4           purrr_0.2.4           readr_1.1.1          
## [25] tidyr_0.7.2           tibble_1.4.1          tidyverse_1.2.1      
## [28] ggmap_2.6.1           ggplot2_2.2.1.9000    rgdal_1.2-16         
## [31] rgeos_0.3-26          maptools_0.9-2        rasterVis_0.43       
## [34] latticeExtra_0.6-28   RColorBrewer_1.1-2    lattice_0.20-35      
## [37] raster_2.6-7          sp_1.2-6             
## 
## loaded via a namespace (and not attached):
##   [1] backports_1.1.2    blme_1.0-4         plyr_1.8.4        
##   [4] lazyeval_0.2.1     TMB_1.7.12         splines_3.4.0     
##   [7] crosstalk_1.0.0    TH.data_1.0-8      digest_0.6.13     
##  [10] foreach_1.4.4      htmltools_0.3.6    magrittr_1.5      
##  [13] modelr_0.1.1       R.utils_2.6.0      sandwich_2.4-0    
##  [16] jpeg_0.1-8         colorspace_1.3-2   rvest_0.3.2       
##  [19] haven_1.1.0        crayon_1.3.4       RCurl_1.95-4.10   
##  [22] jsonlite_1.5       hexbin_1.27.1      bindr_0.1         
##  [25] survival_2.41-3    zoo_1.8-1          iterators_1.0.9   
##  [28] glue_1.2.0         gtable_0.2.0       emmeans_1.1       
##  [31] webshot_0.5.0      sjstats_0.14.0     sjmisc_2.6.3      
##  [34] maps_3.2.0         abind_1.4-5        scales_0.5.0.9000 
##  [37] mvtnorm_1.0-6      DBI_0.7            ggeffects_0.3.0   
##  [40] Rcpp_0.12.14       merTools_0.3.0     viridisLite_0.2.0 
##  [43] xtable_1.8-2       units_0.5-1        foreign_0.8-67    
##  [46] mapproj_1.2-5      stats4_3.4.0       prediction_0.2.0  
##  [49] survey_3.32-1      animation_2.5      htmlwidgets_0.9   
##  [52] httr_1.3.1         modeltools_0.2-21  geosphere_1.5-7   
##  [55] pkgconfig_2.0.1    R.methodsS3_1.7.1  nnet_7.3-12       
##  [58] labeling_0.3       tidyselect_0.2.3   rlang_0.1.6       
##  [61] reshape2_1.4.3     munsell_0.4.3      cellranger_1.1.0  
##  [64] tools_3.4.0        cli_1.0.0          sjlabelled_1.0.6  
##  [67] broom_0.4.3        evaluate_0.10.1    yaml_2.1.16       
##  [70] satellite_1.0.1    RgoogleMaps_1.4.1  bindrcpp_0.2      
##  [73] coin_1.2-2         mime_0.5           R.oo_1.21.0       
##  [76] xml2_1.1.1         compiler_3.4.0     bayesplot_1.4.0   
##  [79] rstudioapi_0.7     png_0.1-7          e1071_1.6-8       
##  [82] stringi_1.1.6      classInt_0.1-24    psych_1.7.8       
##  [85] nloptr_1.0.4       effects_4.0-0      stringdist_0.9.4.6
##  [88] pillar_1.1.0       pwr_1.2-1          lmtest_0.9-35     
##  [91] estimability_1.2   bitops_1.0-6       httpuv_1.3.5      
##  [94] R6_2.2.2           codetools_0.2-15   gdalUtils_2.0.1.7 
##  [97] assertthat_0.2.0   proto_1.0.0        rprojroot_1.3-2   
## [100] rjson_0.2.15       mnormt_1.5-5       multcomp_1.4-8    
## [103] parallel_3.4.0     hms_0.4.0          udunits2_0.13     
## [106] grid_3.4.0         glmmTMB_0.2.0      coda_0.19-1       
## [109] class_7.3-14       minqa_1.2.4        rmarkdown_1.8     
## [112] snakecase_0.7.1    carData_3.0-0      shiny_1.0.5       
## [115] lubridate_1.7.1    base64enc_0.1-3

Literatura Citada

Armenteras, D., N. Rodríguez, and J. Retana. 2013. Landscape dynamics in northwestern Amazonia: an assessment of pastures, fire and illicit crops as drivers of tropical deforestation. PloS one 8:e54310.

Brooks, T. M., R. A. Mittermeier, C. G. Mittermeier, G. A. B. da Fonseca, A. B. Rylands, W. R. Konstant, P. Flick, J. Pilgrim, S. Oldfield, G. Magin, and C. Hilton-Taylor. 2002. Habitat loss and extinction in the hotspots of biodiversity. Conservation Biology 16:909–923.

Geist, H. J., and E. F. Lambin. 2002. Proximate causes and underlying driving forces of tropical deforestation. BioScience 52:143–150.

Hansen, M. C., P. V. Potapov, R. Moore, M. Hancher, S. A. Turubanova, A. Tyukavina, D. Thau, S. V. Stehman, S. J. Goetz, T. R. Loveland, A. Kommareddy, A. Egorov, L. Chini, C. O. Justice, and J. R. G. Townshend. 2013. High-Resolution Global Maps of 21st-Century Forest Cover Change. Science 342:850–853.

Jha, S., and K. S. Bawa. 2006. Population growth, human development, and deforestation in biodiversity hotspots. Conservation Biology 20:906–912.

Sánchez-Cuervo, A. M., T. M. Aide, M. L. Clark, and A. Etter. 2012. Land cover change in Colombia: surprising forest recovery trends between 2001 and 2010. PloS one 7:e43943.

Zvoleff, A. 2015. Gfcanalysis: Tools for working with hansen et al. global forest change dataset.