1. Introduction

This notebook illustrates land cover change in Montes de Maria from 2002 to 2012. Land cover categories are represented by two raster datasets extracted from the Mapas de Cobertura de la Tierra (IDEAM). The diffeR library is used to estimate metrics of difference for comparing the two land cover maps.

Preliminaries

## Loading required package: sp
## Loading required package: rgdal
## rgdal: version: 1.4-8, (SVN revision 845)
##  Geospatial Data Abstraction Library extensions to R successfully loaded
##  Loaded GDAL runtime: GDAL 2.4.2, released 2019/06/28
##  Path to GDAL shared files: /Library/Frameworks/R.framework/Versions/3.6/Resources/library/rgdal/gdal
##  GDAL binary built with GEOS: FALSE 
##  Loaded PROJ.4 runtime: Rel. 5.2.0, September 15th, 2018, [PJ_VERSION: 520]
##  Path to PROJ.4 shared files: /Library/Frameworks/R.framework/Versions/3.6/Resources/library/rgdal/proj
##  Linking to sp version: 1.3-2
## Loading required package: ggplot2
## Loading required package: reshape2

The IDEAM land cover datasets, in vector format, were obtained previously from its geoportal. Then, the datasets were clipped to match the extent of the study area. Finally, the cropped datasets were converted to raster format.

ruta <- "../IDEAM/"
files <- list.files(path = ruta, pattern = "tif", all.files = TRUE)
files
## [1] "LC_2002.tif" "LC_2012.tif"
(file1 <- paste(ruta,files[1],sep=""))
## [1] "../IDEAM/LC_2002.tif"
lc1 <- raster(paste(ruta,files[1],sep=""))
lc2 <- raster(paste(ruta,files[2],sep=""))
lc1
## class      : RasterLayer 
## dimensions : 1840, 1863, 3427920  (nrow, ncol, ncell)
## resolution : 5e-04, 5e-04  (x, y)
## extent     : -75.70386, -74.77236, 9.225483, 10.14548  (xmin, xmax, ymin, ymax)
## crs        : +proj=longlat +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +no_defs 
## source     : /Users/ivanlizarazo/Documents/ivan/PR/montes/IDEAM/LC_2002.tif 
## names      : LC_2002
lc2
## class      : RasterLayer 
## dimensions : 1840, 1862, 3426080  (nrow, ncol, ncell)
## resolution : 5e-04, 5e-04  (x, y)
## extent     : -75.70386, -74.77286, 9.225483, 10.14548  (xmin, xmax, ymin, ymax)
## crs        : +proj=longlat +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +no_defs 
## source     : /Users/ivanlizarazo/Documents/ivan/PR/montes/IDEAM/LC_2012.tif 
## names      : LC_2012
(lc1n <- crop(lc1,  extent(lc2)))
## class      : RasterLayer 
## dimensions : 1840, 1862, 3426080  (nrow, ncol, ncell)
## resolution : 5e-04, 5e-04  (x, y)
## extent     : -75.70386, -74.77286, 9.225483, 10.14548  (xmin, xmax, ymin, ymax)
## crs        : +proj=longlat +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +no_defs 
## source     : /private/var/folders/l_/1nv5l0dx78g_d3k3_t85lrjr0000gn/T/RtmpXKssDY/raster/r_tmp_2020-06-02_122809_20709_78110.grd 
## names      : LC_2002 
## values     : 11, 99  (min, max)
unique(values(lc1))
##  [1] NA 24 31 52 51 23 32 33 42 11 41 21 22 12 99 13 14
unique(values(lc2))
##  [1] NA 24 41 31 52 23 32 33 42 51 11 22 21 13 12
mlc <- stack(lc1n,lc2)

3. Visualize Land Cover Change

This section visualizes changes using the color palette provided by IDEAM.

LC 2002

library(knitr)
include_graphics("yLC2002.png")

LC 2012

library(knitr)
include_graphics("YLC2012.png")

4. Another visualization of Land Cover Change

This time, a customized color palette is used to emphasize changes.

library(leaflet)
library(RColorBrewer)
pal <- colorFactor(c('#152106', '#225129', '#369b47', '#30eb5b', '#387242', 
    '6a2325','#c3aa69', '#b76031' , '#d9903d', '#91af40', 
    '111149', '#3486eb', '#344feb', '#33280d', '#d7cdcc',
    '#f7e084' ), 
    c(11, 12, 13, 14, 21, 22, 23, 24, 31, 32, 33, 41,  42,  51, 52, 99), na.color = "transparent")

LC 2002

leaflet() %>% addTiles() %>%
  addRasterImage(mlc[[1]], colors = pal, opacity = 1.0, method="ngb") %>%
  addLegend(pal = pal, values = c(11, 12, 13, 14, 21, 22, 23, 24, 31, 32, 33, 41,  42,  51, 52, 99),
    title = "Montes de María - LandCover 2002")

LC 2012

leaflet() %>% addTiles() %>%
  addRasterImage(mlc[[2]], colors = pal, opacity = 1.0, method="ngb") %>%
  addLegend(pal = pal, values =  c(11, 12, 13, 14, 21, 22, 23, 24, 31, 32, 33, 41,  42,  51, 52, 99),
    title = "Montes de María - LandCover 2012")

5. Metrics of land cover change from 2002 to 2012

comp <- mlc[[2]] # 2012 Land Cover
ref <-  mlc[[1]] # 2002 Land Cover
ctmatCompRef <- crosstabm(comp, ref)
ctmatCompRef
##       11  12  13    21   22     23     24     31    32   33    41 42    51   52
## 11 15035 157   0    24    0    224   1025     65    94   81     8  0     0    0
## 12     0 172   0     0    0      0    179      0     0    0     0  0     0    0
## 13    39   0 385     0    0    549      0    101    19   46     0  0     0    0
## 21     0   0   0  3930    0    756    477     20     0    0    94  0     0    0
## 22     0   0   0  4126 2258   3636   7953      0     0    0     0  0    13    0
## 23   113  32   0  9423  133 434475 156627  39440 44956 3159 10946  0   805    4
## 24    56   0   0 11692  199 144277 240588  24722 48599  572  6627  0  2360    0
## 31     0   0  42   174   51  36103  66971 201788 26408  206  2234  0  1374  804
## 32     0   0   0   481    0  64813  92105  32977 98758  432  2374  0   315   68
## 33     0   0   0   186    0   2759   2572    992   225 1738     0  0  1740    0
## 41    30   0   0   569    0   1831   2510   3652  1618    0 15144  0  2255  782
## 42     0   0   0     0    0     81      0    248     0    0     0 46     0    0
## 51    35   0   0    33    0    228    914    885   185    0  7023  0 33272    0
## 52     0   0   0     0    0    649    255   1510   148  378    20  0     0 5433
## 14     0   0   0     0    0      0      0      0     0    0     0  0     0    0
## 99     0   0   0     0    0      0      0      0     0    0     0  0     0    0
##    14   99
## 11 88    0
## 12  1    0
## 13  0    0
## 21  0    0
## 22  0    0
## 23  0    0
## 24  0  247
## 31  0 1139
## 32  0  590
## 33  0    0
## 41  0    0
## 42  0    0
## 51  0    0
## 52  0    0
## 14  0    0
## 99  0    0
categoryComponentsPlot(ctmatrix = ctmatCompRef, units = "pixels")

differenceMR(comp, ref, eval='original', percent=TRUE)
##   Resolution Quantity Exchange    Shift  Overall
## 1      5e-04 6.991132 35.46072 3.317866 45.76971
diffTablej(ctmatCompRef)
##    Category Omission Agreement Comission Quantity Exchange Shift
## 1        11      273     15035      1766     1493      354   192
## 2        12      189       172       180        9        0   360
## 3        13       42       385       754      712       84     0
## 4        21    26708      3930      1347    25361     2694     0
## 5        22      383      2258     15728    15345      664   102
## 6        23   255906    434475    265638     9732   462320 49492
## 7        24   331588    240588    239351    92237   444652 34050
## 8        31   104612    201788    135506    30894   182848 26376
## 9        32   122252     98758    194155    71903   244118   386
## 10       33     4874      1738      8474     3600     7524  2224
## 11       41    29326     15144     13247    16079    21140  5354
## 12       42        0        46       329      329        0     0
## 13       51     8862     33272      9303      441     8934  8790
## 14       52     1658      5433      2960     1302     1792  1524
## 15       14       89         0         0       89        0     0
## 16       99     1976         0         0     1976        0     0
## 17  Overall   888738   1053022    888738   135751   688562 64425

It is evident what are the most important land cover changes. Isn’t it?