# "C:/Users/David-PC/Desktop/ESCRITORIO/Landsat_procesamiento"
getwd()
[1] "C:/Users/David-PC/Desktop/ESCRITORIO/Landsat_procesamiento"
library(raster)
# Costera
b1 = raster('./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B1.tif')
# Azul
b2 = raster('./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B2.tif')
# Verde
b3 = raster('./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B3.tif')
# Rojo
b4 = raster('./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B4.tif')
# Infrarojo cercano (NIR)
b5 = raster('./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B5.tif')
# Infrarrojo de Onda Corta 1 (SWIR 1)
b6 = raster('./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B6.tif')
# Infrarrojo de Onda Corta 2 (SWIR 2)
b7 = raster('./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B7.tif')
# Pancromatica
b8 = raster('./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B8.tif')
# Cirros (Cirrus)
b9 = raster('./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B9.tif')
# Sensor Térmico Infrarrojo 1 (TIRS 1)
b10 = raster('./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B10.tif')
# Sensor Térmico Infrarrojo 2 (TIRS 2)
b11 = raster('./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B11.tif')
b2
class : RasterLayer
dimensions : 7741, 7581, 58684521 (nrow, ncol, ncell)
resolution : 30, 30 (x, y)
extent : 446085, 673515, 362985, 595215 (xmin, xmax, ymin, ymax)
crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
source : C:/Users/David-PC/Desktop/ESCRITORIO/Landsat_procesamiento/data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B2.tif
names : LC08_L1TP_008057_20150104_20170415_01_T1_B2
values : 0, 65535 (min, max)
# Sistema de coordenadas de referencia (CRS)
crs(b2)
CRS arguments:
+proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
## CRS arguments:
## +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84
## +towgs84=0,0,0
# Numero de celdas, filas, columnas
ncell(b2)
[1] 58684521
## [1] 58684521
# Dimension
dim(b2)
[1] 7741 7581 1
## [1] 7741 7581 1
# Resolucion espacial
res(b2)
[1] 30 30
## [1] 30 30
# Numero de badas
nlayers(b2)
[1] 1
## [1] 1
# ¿Las bandas tienen la misma extension, numero de filas y columnas, proyeccion, resolucion y origen?
compareRaster(b2,b3)
[1] TRUE
## [1] TRUE
# Creamos el RasterStack
s2015 = stack(b5, b4, b3)
# Revisamos las propiedades del RasterStack
s2015
class : RasterStack
dimensions : 7741, 7581, 58684521, 3 (nrow, ncol, ncell, nlayers)
resolution : 30, 30 (x, y)
extent : 446085, 673515, 362985, 595215 (xmin, xmax, ymin, ymax)
crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
names : LC08_L1TP_008057_20150104_20170415_01_T1_B5, LC08_L1TP_008057_20150104_20170415_01_T1_B4, LC08_L1TP_008057_20150104_20170415_01_T1_B3
min values : 0, 0, 0
max values : 65535, 65535, 65535
## class : RasterStack
## dimensions : 7741, 7581, 58684521, 3 (nrow, ncol, ncell, nlayers)
## resolution : 30, 30 (x, y)
## extent : 446085, 673515, 362985, 595215 (xmin, xmax, ymin, ymax)
## crs : ++proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
## names : LC08_L1TP_008057_20180317_20180403_01_T1_B5, LC08_L1TP_008057_20180317_20180403_01_T1_B4, LC08_L1TP_008057_20180317_20180403_01_T1_B3
## min values : 0, 0, 0
## max values : 65535, 65535, 65535
# Primero creamos una lista de los raster layers usando:
filenames = paste0('./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B', 1:7, ".tif")
filenames
[1] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B1.tif"
[2] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B2.tif"
[3] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B3.tif"
[4] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B4.tif"
[5] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B5.tif"
[6] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B6.tif"
[7] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B7.tif"
## [1] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B1.tif"
## [2] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B2.tif"
## [3] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B3.tif"
## [4] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B4.tif"
## [5] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B5.tif"
## [6] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B6.tif"
## [7] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B7.tif"
###### no lee la banda 8 porque tiene diferente extension, la banda 9 es cirro, y las bandas 10 y 11 son termales
## [8] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B8.tif"
## [9] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B9.tif"
## [10] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B10.tif"
## [11] "./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B11.tif"
# Creo el RasterStack de las bandas seleccionadas, es decir de la 1 a la 7 (1:7)
landsat2015= stack(filenames)
landsat2015
class : RasterStack
dimensions : 7741, 7581, 58684521, 7 (nrow, ncol, ncell, nlayers)
resolution : 30, 30 (x, y)
extent : 446085, 673515, 362985, 595215 (xmin, xmax, ymin, ymax)
crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
names : LC08_L1TP//5_01_T1_B1, LC08_L1TP//5_01_T1_B2, LC08_L1TP//5_01_T1_B3, LC08_L1TP//5_01_T1_B4, LC08_L1TP//5_01_T1_B5, LC08_L1TP//5_01_T1_B6, LC08_L1TP//5_01_T1_B7
min values : 0, 0, 0, 0, 0, 0, 0
max values : 65535, 65535, 65535, 65535, 65535, 65535, 65535
## class : RasterStack
## dimensions : 7741, 7581, 58684521, 7 (nrow, ncol, ncell, nlayers)
## resolution : 30, 30 (x, y)
## extent : 446085, 673515, 362985, 595215 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=10 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
## names : LC08_L1TP//5_01_T1_B1, LC08_L1TP//5_01_T1_B2, LC08_L1TP//5_01_T1_B3, LC08_L1TP//5_01_T1_B4, LC08_L1TP//5_01_T1_B5, LC08_L1TP//5_01_T1_B6, LC08_L1TP//5_01_T1_B7
## min values : 0, 0, 0, 0, 0, 0, 0
## max values : 65535, 65535, 65535, 65535, 65535, 65535, 65535
# Trazando capas indivuduales (aZUL, VERDE, ROJO, INFRAROJO CERCANO) a partir del RasterStack
par(mfrow = c(2,2))
plot(b2, main = "Azul", col = gray(0:100 / 100), axes = FALSE)
plot(b3, main = "Verde", col = gray(0:100 / 100), axes = FALSE)
plot(b4, main = "Rojo", col = gray(0:100 / 100), axes = FALSE)
plot(b5, main = "NIR", col = gray(0:100 / 100), axes = FALSE)
# Ploteamos la combinacion de bandas 4, 3 y 2 que para L8 es color verdadero
landsatRGB2015 = stack(b4, b3, b2)
plotRGB(landsatRGB2015, axes = TRUE, stretch = "lin", main = "Composicion de color verdadero - Landsat 8")
# Ploteamos la combinacion de bandas 5, 4 y 3 que para L8 es color falso
par(mfrow = c(1,2))
plotRGB(landsatRGB2015, axes=TRUE, stretch="lin", main="Composicion de color verdadero")
landsatFCC2015 = stack(b5, b4, b3)
plotRGB(landsatFCC2015, axes=TRUE, stretch="lin", main="Composicion de color falso")
## solucion a la pregunta 1, ploteamos la combinacion de bandas 5, 6 y 4 que para L8 es color falso
landsatFCC2015_p1 = stack(b5, b6, b4)
plotRGB(landsatFCC2015_p1, axes=TRUE, stretch="lin", main="Composicion de color falsa - Agua / Tierra")
# Seleccionamos unicamente las 3 primeras bandas
landsat2015sub01 = subset(landsat2015, 1:3)
# Mismo
landsat2015sub02 = landsat2015[[1:3]]
# Numero de bandas en los datos originales y nuevos.
nlayers(landsat2015)
[1] 7
## [1] 7
nlayers(landsat2015sub01)
[1] 3
## [1] 3
nlayers(landsat2015sub02)
[1] 3
## [1] 3
## Eliminando las ultimas cuatro bandas, es decir, seleccionaremos para nuestro ejercicio las bandas de 1 a la 7 (1:7)
landsat2015 = subset(landsat2015, 1:7)
plot(landsat2015)
# Aca vemos que nombres tienen las bandas
names(landsat2015)
[1] "LC08_L1TP_008057_20150104_20170415_01_T1_B1" "LC08_L1TP_008057_20150104_20170415_01_T1_B2"
[3] "LC08_L1TP_008057_20150104_20170415_01_T1_B3" "LC08_L1TP_008057_20150104_20170415_01_T1_B4"
[5] "LC08_L1TP_008057_20150104_20170415_01_T1_B5" "LC08_L1TP_008057_20150104_20170415_01_T1_B6"
[7] "LC08_L1TP_008057_20150104_20170415_01_T1_B7"
## [1] "LC08_L1TP_008057_20180317_20180403_01_T1_B1" "LC08_L1TP_008057_20180317_20180403_01_T1_B2"
## [3] "LC08_L1TP_008057_20180317_20180403_01_T1_B3" "LC08_L1TP_008057_20180317_20180403_01_T1_B4"
## [5] "LC08_L1TP_008057_20180317_20180403_01_T1_B5" "LC08_L1TP_008057_20180317_20180403_01_T1_B6"
## [7] "LC08_L1TP_008057_20180317_20180403_01_T1_B7"
# Aca las renombramos y las ploteamos
names(landsat2015) = c('ultra-azul', 'azul', 'verde', 'rojo', 'NIR', 'SWIR1', 'SWIR2')
names(landsat2015)
[1] "ultra.azul" "azul" "verde" "rojo" "NIR" "SWIR1" "SWIR2"
## [1] "ultra.azuk" "azul" "verde" "rojo" "NIR" "SWIR1" "SWIR2"
## [1] "ultra.azul" "azul" "verde" "rojo" "NIR" "SWIR1" "SWIR2"
# Usando Extension
extent(landsat2015)
class : Extent
xmin : 446085
xmax : 673515
ymin : 362985
ymax : 595215
## class : Extent
## xmin : 446085
## xmax : 673515
## ymin : 362985
## ymax : 595215
# Cargamos un objeto espacial que se encuentra en un shapefile. Esta extension es el area de Sibate
sibate2015 = shapefile('C:/Users/David-PC/Desktop/ESCRITORIO/SIBATE_RURAL/rural.shp')
# Define que la extension es Sibate
e = extent(sibate2015)
# Corta la imagen landsat con el area de sibate
landsatcrop2015 = crop(landsat2015, e)
# Plotea a Sibate
plotRGB(landsatcrop2015, r= 6, g= 5, b=2, axes = TRUE, stretch = "lin", main = "Composicion de color falsa - Sibaté")
# Creacion de un compuesto de color verdadero, para landsat 8 la combinacion de color es 4,3,2
landsatcrop2015pregunta = crop(stack(b4, b3, b2), e)
# Ploteamos la combinacion
plotRGB(landsatcrop2015pregunta, axes = TRUE, stretch = "lin", main = "Sibaté 2015 - Color Verdadero")
# Aca guardanmos nuestro raster de Sibate y lo imprimimos para ver lo que se ha guardado
x = writeRaster(landsatcrop2015, filename="Sibaté_banda.tif", overwrite=TRUE)
plot(x)
# Lo guardamos tambien en formato .grd ya que este formato si guarda los nombres de las capas
writeRaster(landsatcrop2015, filename="Sibate_banda.grd", overwrite=TRUE)
# Comparacion del reflejo en la longitud de onda Ulta azul vs azul
pairs(landsatcrop2015[[1:2]], main = "Ultra-Azul versus Azul")
# Comparacion del reflejo en la longitud de onda Roja vs NIR
pairs(landsatcrop2015[[4:5]], main = "Roja versus NIR")
# Cargar los poligonos con informacion sobre el uso del suelo
sibatecobert = shapefile ('C:/Users/David-PC/Desktop/ESCRITORIO/PERCEPCION REMOTA/TRABAJO_COBERTURAS/cobertu.shp')
Z-dimension discarded
# Generar muestras de 300 puntos a partir de los poligonos
ptsibatecobert = spsample(sibatecobert,300, type='regular')
# Agregar la clase de cobertura del suelo a los puntos
ptsibatecobert$Name = over(ptsibatecobert, sibatecobert)$Name
# Extraer valores con puntos
df = extract(landsat2015, ptsibatecobert)
# Para ver algunos de los valores de reflectancia
head(df)
ultra.azul azul verde rojo NIR SWIR1 SWIR2
[1,] 7928 7237 6660 6081 12490 7668 6064
[2,] 7957 7278 6700 6101 13199 7933 6123
[3,] 7952 7288 6733 6158 13073 8111 6255
[4,] 7951 7297 6744 6105 14049 7740 6011
[5,] 7953 7273 6723 6076 13661 8025 6186
[6,] 7959 7262 6716 6121 12563 8113 6218
# ultra.azul azul verde rojo NIR SWIR1 SWIR2
#[1,] 7928 7238 6623 6033 12696 7520 5997
#[2,] 7957 7278 6700 6101 13199 7933 6123
#[3,] 7973 7285 6741 6148 13482 8072 6219
#[4,] 7962 7262 6710 6114 13413 8116 6199
#[5,] 7931 7246 6671 6060 13541 7655 6014
#[6,] 7943 7242 6623 6024 12437 7529 6039
numero= length(ptsibatecobert)
numero
[1] 298
df_samples = as (ptsibatecobert, "SpatialPointsDataFrame")
df_samples
class : SpatialPointsDataFrame
features : 298
extent : 577367, 585597.8, 486272.4, 501061.1 (xmin, xmax, ymin, ymax)
crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
variables : 1
names : Name
min values : Agua
max values : zona desnuda
df_samples@data=data.frame(ID=1:numero,size=1)
df_samples
class : SpatialPointsDataFrame
features : 298
extent : 577367, 585597.8, 486272.4, 501061.1 (xmin, xmax, ymin, ymax)
crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
variables : 2
names : ID, size
min values : 1, 1
max values : 298, 1
plot(landsatcrop2015)
plot(sibatecobert, add= TRUE)
plot(df_samples,pch=1, cex=(df_samples$size)/4, add=TRUE)
df_samples$Name =over(df_samples,sibatecobert)$Name
df_samples
class : SpatialPointsDataFrame
features : 298
extent : 577367, 585597.8, 486272.4, 501061.1 (xmin, xmax, ymin, ymax)
crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
variables : 3
names : ID, size, Name
min values : 1, 1, Agua
max values : 298, 1, zona desnuda
df1=raster::extract(landsatcrop2015, df_samples)
ms = aggregate(df1, list(ptsibatecobert$Name), mean)
# En lugar de la primera columna, usamos nombres de fila
rownames(ms) = ms[,1]
ms = ms[,-1]
ms
## ultra.blue blue green red NIR SWIR1 SWIR2
# Agua 8251.053 7501.298 6840.658 6066.404 5762.921 5081.289 5050.825
# Bosque 8026.761 7326.134 6773.410 6172.284 11207.410 7962.134 6412.746
# Cultivos 8417.864 7800.318 7603.682 7197.955 18156.818 14063.045 11198.864
# Fresa 8541.667 7870.667 7464.667 6902.000 13157.667 9832.333 7806.000
# Invernadero 11522.083 11035.833 11442.333 11179.083 18368.333 14785.000 11300.083
# pastos 8481.000 7926.727 8267.909 7331.636 21093.273 12449.000 8507.273
# zona desnuda 10975.500 11035.500 12884.750 14428.000 17930.750 20671.750 16907.500
# Cree un vector de color para las clases de cobertura del suelo para usar en el trazado
mycolor = c('darkred', 'yellow', 'burlywood', 'cyan', 'blue', 'green','magenta', "darkgreen")
# transformar ms de un data.frame a una matriz
ms = as.matrix(ms)
# Primero crea una parcela vacia
plot(0, ylim=c(0,25000), xlim = c(1,7), type='n', xlab="Bandas", ylab = "Reflectancia")
# Agrega las diferentes clases
for (i in 1:nrow(ms)){
lines(ms[i,], type = "l", lwd = 3, lty = 1, col = mycolor[i])
}
# Titulo
title(main="Perfil espectral de Landsat", font.main = 2)
# Leyenda
legend("topleft", rownames(ms),
cex=0.8, col=mycolor, lty = 1, lwd =3, bty = "n")
### solo puedo cargar las bandas de la 1 a la 7 porque la 8 no tiene la misma extension de las demas
raslistp3 = paste0('./data2015/LC08_L1TP_008057_20150104_20170415_01_T1_B', 1:7, ".tif")
##raslistp3= practica 3
# Stack de las bandas cargadas
landsatp3 = stack(raslistp3)
plot(landsatp3)
# Combinacion de las bandas 4,3,2
landsatRGBp3 = landsatp3[[c(4,3,2)]]
## landsatRGBp3 = practica 3
plot(landsatRGBp3)
# Combinacion de las bandas en falso color 5,4,3
landsatFCCp3 = landsatp3[[c(5,4,3)]]
## landsatFCCp3 = practica 3
plot(landsatFCCp3)
## vi2015= indice de vegetacion para la imagen 2015
vi2015 = function(img, k, i) {
bk = img [[k]]
bi = img [[i]]
vi2015=(bk-bi)/(bk+bi)
return(vi2015)
}
## para landsat NIR = 5, red = 4
## ndvi2015= indice de vegetacion de direfencia normalizada (NDVI) para la imagen landsat 2015
ndvi2015 = vi2015(landsatp3, 5, 4)
plot(ndvi2015, col=rev(terrain.colors(9)), main= "Landsat 2015 - NVDI")
# Seleccionando el shape del area de sibate, cargandolo y recortandolo con respecto del NDVI de la imagen Landsat
NDVI1_sibate_2015 = shapefile('C:/Users/David-PC/Desktop/ESCRITORIO/SIBATE_RURAL/rural.shp')
eNDVI1_sibate = extent(NDVI1_sibate_2015)
NDVI_recorte_1 = crop(ndvi2015, eNDVI1_sibate)
plot(NDVI_recorte_1, col=rev(terrain.colors(10)), main = "Landsat 2015 - NDVI")
vi2_2015 = function(x, y) {
(x-y)/(x+y)
}
ndvi2_2015 = overlay(landsatp3[[5]], landsatp3[[4]], fun = vi2_2015)
plot(ndvi2_2015, col=rev(terrain.colors(10)), main = "Landsat 2015 - NDVI v2")
# Seleccionando el shape del area de sibate, cargandolo y recortandolo con respecto del NDVI de la imagen Landsat
NDVI2_sibate_2015 = shapefile('C:/Users/David-PC/Desktop/ESCRITORIO/SIBATE_RURAL/rural.shp')
eNDVI2_sibate = extent(NDVI2_sibate_2015)
NDVI_recorte_2 = crop(ndvi2_2015, eNDVI2_sibate)
plot(NDVI_recorte_2, col=rev(terrain.colors(10)), main = "Landsat 2015 - NDVI")
## p1 es pregunta 1
ndvi2015_p1 = vi2015(landsatp3, 3, 6)
plot(ndvi2015_p1, col=rev(terrain.colors(10)), main= "Landsat 2015 - NDWI")
vi2_2015_p1 = function(x, y) {
(x-y)/(x+y)
}
ndvi2_2015_p1 = overlay(landsatp3[[3]], landsatp3[[6]], fun = vi2_2015_p1)
plot(ndvi2_2015_p1, col=rev(topo.colors(10)), main = "Landsat 2015 - NDWI")
### PARA INDICE DE AGUA DE DIFERENCIA NORMALIZADA - NDWI
# Seleccionando el shape del area de sibate, cargandolo y recortandolo con respecto del NDwI de la imagen Landsat
NDWI_sibate = shapefile('C:/Users/David-PC/Desktop/ESCRITORIO/SIBATE_RURAL/rural.shp')
eNDWI = extent(NDWI_sibate)
NDWI_recorte = crop(ndvi2_2015_p1, eNDWI)
plot(NDWI_recorte, col=rev(topo.colors(10)), main = "Landsat 2015 - NDWI")
### PARA INDICE INCORPORADO DE DIFERENCIA NORMALIZADA - NDBI
ndbi1_2015_p1 = overlay(landsatp3[[6]], landsatp3[[5]], fun = vi2_2015_p1)
plot(ndbi1_2015_p1, col=rev(heat.colors(10)), main = "Landsat 2015 - NDBI")
# Seleccionando el shape del area de sibate, cargandolo y recortandolo con respecto del NDBI de la imagen Landsat
NDBI_sibate = shapefile('C:/Users/David-PC/Desktop/ESCRITORIO/SIBATE_RURAL/rural.shp')
eNDBI = extent(NDBI_sibate)
NDBI_recorte = crop(ndbi1_2015_p1, eNDBI)
plot(NDBI_recorte, col=rev(heat.colors(10)), main = "Landsat 2015 - NDBI")
## ver histograma de datos
hist(NDVI_recorte_1,
main= "Dsitribucion de valores del NDVI",
xlab= "NDVI",
ylab= "Frecuencia",
col= terrain.colors(50),
xlim= c (-0.1,0.7),
breaks = 30,
xaxt= 'n')
axis(side = 1, at=seq(-0.5,1,0.05), labels = seq(-0.5,1,0.05))
## respuesta a la pregunta 2: ver histograma de datos
hist(NDWI_recorte,
main= "Dsitribucion de valores del NDWI",
xlab= "NDWI",
ylab= "Frecuencia",
col= topo.colors(50),
xlim= c (-0.5,0.3),
breaks = 30,
xaxt= 'n')
axis(side = 1, at=seq(-0.5,1,0.05), labels = seq(-0.5,1,0.05))
## respuesta a la pregunta 2: ver histograma de datos
hist(NDBI_recorte,
main= "Dsitribucion de valores del NDBI",
xlab= "NDBI",
ylab= "Frecuencia",
col= heat.colors(50),
xlim= c (-0.5,0.3),
breaks = 30,
xaxt= 'n')
axis(side = 1, at=seq(-0.5,1,0.05), labels = seq(-0.5,1,0.05))
veg = reclassify(NDVI_recorte_1, cbind(-Inf, 0.4, NA))
plot(veg, main = "Vegetacion")
land = reclassify(NDVI_recorte_1, c(-Inf, 0.25, NA, 0.25, 0.3, 1, 0.3, Inf, NA))
plot(land, main = " Que es esto?")
plotRGB(landsatcrop2015, r=7, g=6, b=4, axes = TRUE, stretch="lin", main = "Sibate composicion en falso color")
plot(land, add = TRUE, legend= FALSE)
vegc = reclassify(NDVI_recorte_1, c(-Inf,0.25,1, 0.25,0.3,2, 0.3,0.4,3, 0.4,0.5,4, 0.5,Inf, 5))
plot(vegc, col = rev(terrain.colors(4)), main = 'Umbral basado en NDVI')
# Umbral, en este caso para agua
agua = reclassify(NDWI_recorte, cbind(-Inf, 0.1, NA))
plot(agua, main = "Agua", col= topo.colors(10))
#Reclasifico el Umbral, en este caso para agua
aguac = reclassify(NDWI_recorte, c(-Inf,0.025,1, 0.025,0.03,2, 0.03,0.1,3, 0.1,0.2,4, 0.2,Inf, 5))
plot(aguac, col = rev(topo.colors(10)), main = 'Umbral basado en NDWI')
set.seed(1)
sr = sampleRandom(landsatp3, 100000)
plot(sr[,c(4,5)], main = "NIR - ROJO")
pca = prcomp(sr, scale = TRUE)
pca
Standard deviations (1, .., p=7):
[1] 2.59169224 0.42707865 0.29072035 0.10439979 0.06868262 0.02120325 0.01227403
Rotation (n x k) = (7 x 7):
PC1 PC2 PC3 PC4
LC08_L1TP_008057_20150104_20170415_01_T1_B1 -0.3828686 0.146466958 -0.285016623 0.64295434
LC08_L1TP_008057_20150104_20170415_01_T1_B2 -0.3824643 0.251163610 -0.246675390 0.25077406
LC08_L1TP_008057_20150104_20170415_01_T1_B3 -0.3831958 0.239695066 -0.165224729 -0.23270432
LC08_L1TP_008057_20150104_20170415_01_T1_B4 -0.3785953 0.424481057 0.002225169 -0.62514538
LC08_L1TP_008057_20150104_20170415_01_T1_B5 -0.3622638 -0.745210227 -0.440008806 -0.23096773
LC08_L1TP_008057_20150104_20170415_01_T1_B6 -0.3766513 -0.349896939 0.516775974 0.01974348
LC08_L1TP_008057_20150104_20170415_01_T1_B7 -0.3792861 -0.007760711 0.608233062 0.15820288
PC5 PC6 PC7
LC08_L1TP_008057_20150104_20170415_01_T1_B1 -0.10311614 0.074615517 -0.566686791
LC08_L1TP_008057_20150104_20170415_01_T1_B2 -0.04410742 0.257250131 0.773807743
LC08_L1TP_008057_20150104_20170415_01_T1_B3 -0.03046576 -0.843888776 0.034355651
LC08_L1TP_008057_20150104_20170415_01_T1_B4 -0.05067863 0.455034481 -0.275762605
LC08_L1TP_008057_20150104_20170415_01_T1_B5 0.23955731 0.093210905 -0.019919926
LC08_L1TP_008057_20150104_20170415_01_T1_B6 -0.68244282 -0.008319047 0.049611469
LC08_L1TP_008057_20150104_20170415_01_T1_B7 0.67882929 -0.017110515 0.002057162
screeplot(pca)
pci = predict(landsatp3, pca, index = 1:2)
plot(pci[[1]])
## Recortado a sibate
pci_sibate = shapefile('C:/Users/David-PC/Desktop/ESCRITORIO/SIBATE_RURAL/rural.shp')
epci = extent(pci_sibate)
pci_recorte = crop(pci, epci)
plot(pci_recorte, col=rev(terrain.colors(20)), main = "pci Sibate")
landsatFCC2015 = stack(b5, b4, b3)
## recortado a sibate
pc2 = shapefile('C:/Users/David-PC/Desktop/ESCRITORIO/SIBATE_RURAL/rural.shp')
epc2 = extent(pc2)
pc2_rec= crop(landsatFCC2015, epc2)
pc2 = reclassify(pci_recorte[[2]], c(-Inf,0,1,0,Inf,NA))
par(mfrow = c(1,2))
plotRGB(pc2_rec, r = 1, g = 2, b = 3, axes = TRUE, stretch = "lin", main = "FFC")
plotRGB(pc2_rec, r = 1, g = 2, b = 3, axes = TRUE, stretch = "lin", main = "FFC")
plot(pc2, legend = FALSE, add = TRUE)
library(raster)
landsat8 = stack(b2,b3,b4,b5,b6,b7)
names(landsat8) = c ('blue', 'green', 'red', 'NIR', 'SWIR1','SWIR2')
plot(landsat8)
landsat8FCC_imagen = stack(b5, b6, b4)
plotRGB(landsat8FCC_imagen, axes=TRUE, stretch="lin", main="Landsat composicion de color falsa - Agua / Tierra")
# recortando a sibate
sibatextend = shapefile('C:/Users/David-PC/Desktop/ESCRITORIO/SIBATE_RURAL/rural.shp')
ext = extent(sibatextend)
landsatfcc = crop(landsat8FCC_imagen, ext)
plotRGB(landsatfcc, axes = TRUE, stretch = "lin", main = "Landsat True Color Composite")
ndviLandSat8=(landsat8[['NIR']]-landsat8[['red']])/(landsat8[['NIR']]+landsat8[['red']])
plot(ndviLandSat8)
# recortando a sibate
extent(landsat8)
class : Extent
xmin : 446085
xmax : 673515
ymin : 362985
ymax : 595215
landsat8crop2015 = crop(ndviLandSat8, ext)
landsat8crop2015
class : RasterLayer
dimensions : 617, 351, 216567 (nrow, ncol, ncell)
resolution : 30, 30 (x, y)
extent : 576795, 587325, 484005, 502515 (xmin, xmax, ymin, ymax)
crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
source : memory
names : layer
values : -0.1625031, 0.6681765 (min, max)
plot(landsat8crop2015)
# convertir el raster a vector / matriz
nr=getValues(landsat8crop2015)
str(nr)
num [1:216567] 0.471 0.441 0.435 0.396 0.378 ...
#El algoritmo k-means resuelve un problema de optimizacion, siendo la funcion a optimizar (minimizar) la suma de las distancias cuadraticas de cada objeto al centroide de su cluster. # #
# Es importante configurar el generador de semillas porque `kmeans` inicia los centros en ubicaciones aleatorias, para eso utilizaremos set.seed
set.seed (99)
# Queremos crear 10 grupos, permitir 500 iteraciones, comenzar con 5 conjuntos aleatorios utilizando el método "Lloyd"
kmncluster = kmeans(na.omit(nr), centers = 10, iter.max = 500, nstart = 5, algorithm = "Lloyd")
# kmeans devuelve un objeto de la clase "kmeans"
str(kmncluster)
List of 9
$ cluster : int [1:216567] 6 10 10 3 3 3 3 10 6 3 ...
$ centers : num [1:10, 1] 0.285 0.099 0.376 0.548 0.331 ...
..- attr(*, "dimnames")=List of 2
.. ..$ : chr [1:10] "1" "2" "3" "4" ...
.. ..$ : NULL
$ totss : num 3096
$ withinss : num [1:10] 6 4.46 6.17 7.92 5.96 ...
$ tot.withinss: num 57.6
$ betweenss : num 3039
$ size : int [1:10] 32049 5743 33820 9132 34882 20727 25928 6570 19149 28567
$ iter : int 163
$ ifault : NULL
- attr(*, "class")= chr "kmeans"
# Use el objeto ndvi (en este caso se llama "landsat8crop2015") para establecer los valores del cluster en un nuevo raster
knr = setValues(landsat8crop2015, kmncluster$cluster)
# Tambien puedes hacerlo asi
knr= raster(landsat8crop2015)
values(knr) = kmncluster$cluster
knr
class : RasterLayer
dimensions : 617, 351, 216567 (nrow, ncol, ncell)
resolution : 30, 30 (x, y)
extent : 576795, 587325, 484005, 502515 (xmin, xmax, ymin, ymax)
crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
source : memory
names : layer
values : 1, 10 (min, max)
#Defina un vector de color para 10 grupos (mas informacion sobre como configurar el color mas adelante)
mycolor = c ("#fef65b","#ff0000", "#daa520","#0000ff","#0000ff","#00ff00","#cbbeb5",
"#c3ff5b", "#ff7373", "#00ff00", "#808080")
par(mfrow = c(1,2))
plot(landsat8crop2015, col = rev(terrain.colors(10)), main = 'Landsat-NDVI')
plot(knr, main = 'Unsupervised classification', col = mycolor )
landsat8p2 = stack(b6, b3, b2)
#plotRGB(landsat8p2, axes=TRUE, stretch="lin", main="Landsat composicion de color falsa ")
#crop landsat by the extent
landsatsibatep2 = crop(landsat8p2, ext)
plotRGB(landsatsibatep2, axes = TRUE, stretch = "lin", main = "Landsat FCC - SIBATE")
# convertir el raster a vector / matriz
nr_p2=getValues(landsatsibatep2)
str(nr_p2)
int [1:216567, 1:3] 13348 12473 12269 13387 14576 14992 14842 12425 11290 9894 ...
- attr(*, "dimnames")=List of 2
..$ : NULL
..$ : chr [1:3] "LC08_L1TP_008057_20150104_20170415_01_T1_B6" "LC08_L1TP_008057_20150104_20170415_01_T1_B3" "LC08_L1TP_008057_20150104_20170415_01_T1_B2"
# Es importante configurar el generador de semillas porque `kmeans` inicia los centros en ubicaciones aleatorias, para eso utilizaremos set.seed
set.seed (99)
# Queremos crear 10 grupos, permitir 500 iteraciones, comenzar con 5 conjuntos aleatorios utilizando el método "Lloyd"
kmncluster_p2 = kmeans(na.omit(nr), centers = 10, iter.max = 500, nstart = 5, algorithm = "Lloyd")
# kmeans devuelve un objeto de la clase "kmeans"
str(kmncluster_p2)
List of 9
$ cluster : int [1:216567] 6 10 10 3 3 3 3 10 6 3 ...
$ centers : num [1:10, 1] 0.285 0.099 0.376 0.548 0.331 ...
..- attr(*, "dimnames")=List of 2
.. ..$ : chr [1:10] "1" "2" "3" "4" ...
.. ..$ : NULL
$ totss : num 3096
$ withinss : num [1:10] 6 4.46 6.17 7.92 5.96 ...
$ tot.withinss: num 57.6
$ betweenss : num 3039
$ size : int [1:10] 32049 5743 33820 9132 34882 20727 25928 6570 19149 28567
$ iter : int 163
$ ifault : NULL
- attr(*, "class")= chr "kmeans"
# Use el objeto ndvi (en este caso se llama "landsatsibatep2") para establecer los valores del cluster en un nuevo raster
knr_p2 = setValues(landsatsibatep2, kmncluster_p2$cluster)
# Tambien puedes hacerlo asi
knr_p2= raster(landsatsibatep2)
values(knr_p2) = kmncluster_p2$cluster
knr_p2
class : RasterLayer
dimensions : 617, 351, 216567 (nrow, ncol, ncell)
resolution : 30, 30 (x, y)
extent : 576795, 587325, 484005, 502515 (xmin, xmax, ymin, ymax)
crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
source : memory
names : layer
values : 1, 10 (min, max)
#Defina un vector de color para 10 grupos (mas informacion sobre como configurar el color mas adelante)
#mycolor_p2 = c (topo.colors(10))
plot(knr_p2, main = 'Unsupervised classification', col = rev(topo.colors(10)))
# Directorio de trabajo
setwd("C:/Users/David-PC/Desktop/ESCRITORIO/Landsat_procesamiento/")
# Importar las imagenes
B2 = raster("LC08_L1TP_008057_20150104_20170415_01_T1/LC08_L1TP_008057_20150104_20170415_01_T1_B2.TIF")
B3 = raster("LC08_L1TP_008057_20150104_20170415_01_T1/LC08_L1TP_008057_20150104_20170415_01_T1_B3.TIF")
B4 = raster("LC08_L1TP_008057_20150104_20170415_01_T1/LC08_L1TP_008057_20150104_20170415_01_T1_B4.TIF")
B5 = raster("LC08_L1TP_008057_20150104_20170415_01_T1/LC08_L1TP_008057_20150104_20170415_01_T1_B5.TIF")
B6 = raster("LC08_L1TP_008057_20150104_20170415_01_T1/LC08_L1TP_008057_20150104_20170415_01_T1_B6.TIF")
B7 = raster("LC08_L1TP_008057_20150104_20170415_01_T1/LC08_L1TP_008057_20150104_20170415_01_T1_B7.TIF")
# Datos imagen MTL
dia <- 004
SUN_ELEVATION = 51.89134539
RADIANCE_MULT_BAND_2 = 1.3298E-02
RADIANCE_MULT_BAND_3 = 1.2254E-02
RADIANCE_MULT_BAND_4 = 1.0333E-02
RADIANCE_MULT_BAND_5 = 6.3236E-03
RADIANCE_MULT_BAND_6 = 1.5726E-03
RADIANCE_MULT_BAND_7 = 5.3006E-04
RADIANCE_ADD_BAND_2 = -66.49141
RADIANCE_ADD_BAND_3 = -61.27126
RADIANCE_ADD_BAND_4 = -51.66738
RADIANCE_ADD_BAND_5 = -31.61786
RADIANCE_ADD_BAND_6 = -7.86307
RADIANCE_ADD_BAND_7 = -2.65028
REFLECTANCE_MAXIMUM_BAND_2 = 1.210700
REFLECTANCE_MAXIMUM_BAND_3 = 1.210700
REFLECTANCE_MAXIMUM_BAND_4 = 1.210700
REFLECTANCE_MAXIMUM_BAND_5 = 1.210700
REFLECTANCE_MAXIMUM_BAND_6 = 1.210700
REFLECTANCE_MAXIMUM_BAND_7 = 1.210700
RADIANCE_MAXIMUM_BAND_2 = 805.01141
RADIANCE_MAXIMUM_BAND_3 = 741.81116
RADIANCE_MAXIMUM_BAND_4 = 625.53699
RADIANCE_MAXIMUM_BAND_5 = 382.79742
RADIANCE_MAXIMUM_BAND_6 = 95.19823
RADIANCE_MAXIMUM_BAND_7 = 32.08690
# Conversion
SUN_ELEVATION_R <- SUN_ELEVATION*pi/180 # radianes
d <- 1+0.0167*(sin((2*pi*(dia-93.5))/365)) # distancia del sol a la tierra
# Determinacion ESUN
ESUN2 <- pi*(d^2)*RADIANCE_MAXIMUM_BAND_2/REFLECTANCE_MAXIMUM_BAND_2
ESUN3 <- pi*(d^2)*RADIANCE_MAXIMUM_BAND_3/REFLECTANCE_MAXIMUM_BAND_3
ESUN4 <- pi*(d^2)*RADIANCE_MAXIMUM_BAND_4/REFLECTANCE_MAXIMUM_BAND_4
ESUN5 <- pi*(d^2)*RADIANCE_MAXIMUM_BAND_5/REFLECTANCE_MAXIMUM_BAND_5
ESUN6 <- pi*(d^2)*RADIANCE_MAXIMUM_BAND_6/REFLECTANCE_MAXIMUM_BAND_6
ESUN7 <- pi*(d^2)*RADIANCE_MAXIMUM_BAND_7/REFLECTANCE_MAXIMUM_BAND_7
# Elimar valores nulos
B2[B2==0] <- NA
B3[B3==0] <- NA
B4[B4==0] <- NA
B5[B5==0] <- NA
B6[B6==0] <- NA
B7[B7==0] <- NA
# Determinacion de radiancia sensor
L2 <- RADIANCE_MULT_BAND_2*B2 + RADIANCE_ADD_BAND_2
L3 <- RADIANCE_MULT_BAND_3*B3 + RADIANCE_ADD_BAND_3
L4 <- RADIANCE_MULT_BAND_4*B4 + RADIANCE_ADD_BAND_4
L5 <- RADIANCE_MULT_BAND_5*B5 + RADIANCE_ADD_BAND_5
L6 <- RADIANCE_MULT_BAND_6*B6 + RADIANCE_ADD_BAND_6
L7 <- RADIANCE_MULT_BAND_7*B7 + RADIANCE_ADD_BAND_7
# Determinacion de la radiancia minima
Lmin2 <- RADIANCE_MULT_BAND_2*min(getValues(B2), na.rm = TRUE)+ RADIANCE_ADD_BAND_2
Lmin3 <- RADIANCE_MULT_BAND_3*min(getValues(B3), na.rm = TRUE)+ RADIANCE_ADD_BAND_3
Lmin4 <- RADIANCE_MULT_BAND_4*min(getValues(B4), na.rm = TRUE)+ RADIANCE_ADD_BAND_4
Lmin5 <- RADIANCE_MULT_BAND_5*min(getValues(B5), na.rm = TRUE)+ RADIANCE_ADD_BAND_5
Lmin6 <- RADIANCE_MULT_BAND_6*min(getValues(B6), na.rm = TRUE)+ RADIANCE_ADD_BAND_6
Lmin7 <- RADIANCE_MULT_BAND_7*min(getValues(B7), na.rm = TRUE)+ RADIANCE_ADD_BAND_7
# Determinacion radiancia del objeto oscuro
LDOS2 <- 0.01*ESUN2*sin(SUN_ELEVATION_R)/pi*(d^2)
LDOS3 <- 0.01*ESUN3*sin(SUN_ELEVATION_R)/pi*(d^2)
LDOS4 <- 0.01*ESUN4*sin(SUN_ELEVATION_R)/pi*(d^2)
LDOS5 <- 0.01*ESUN5*sin(SUN_ELEVATION_R)/pi*(d^2)
LDOS6 <- 0.01*ESUN6*sin(SUN_ELEVATION_R)/pi*(d^2)
LDOS7 <- 0.01*ESUN7*sin(SUN_ELEVATION_R)/pi*(d^2)
# Determinar el efecto bruma
LP2 <- Lmin2-LDOS2
LP3 <- Lmin3-LDOS3
LP4 <- Lmin4-LDOS4
LP5 <- Lmin5-LDOS5
LP6 <- Lmin6-LDOS6
LP7 <- Lmin7-LDOS7
# Banda Reflectancia Superficie Dark Object Substrction (DOS)
RS_dos_B2 <- (pi*(L2-LP2)*(d)^2)/(ESUN2*sin(SUN_ELEVATION_R))
RS_dos_B3 <- (pi*(L3-LP3)*(d)^2)/(ESUN3*sin(SUN_ELEVATION_R))
RS_dos_B4 <- (pi*(L4-LP4)*(d)^2)/(ESUN4*sin(SUN_ELEVATION_R))
RS_dos_B5 <- (pi*(L5-LP5)*(d)^2)/(ESUN5*sin(SUN_ELEVATION_R))
RS_dos_B6 <- (pi*(L6-LP6)*(d)^2)/(ESUN6*sin(SUN_ELEVATION_R))
RS_dos_B7 <- (pi*(L7-LP7)*(d)^2)/(ESUN7*sin(SUN_ELEVATION_R))
# Combinacion de bandas
L8_B234567_RS_DOS <- stack(RS_dos_B2,RS_dos_B3,RS_dos_B4,RS_dos_B5,RS_dos_B6,RS_dos_B7)
# Exportacion imagen reflectancia Superfie DOS
#setwd("C:/Users/David-PC/Desktop/ESCRITORIO/Landsat_procesamiento/SHAPE/")
writeRaster(L8_B234567_RS_DOS,"L8_B234567_RS_DOS1.tif", drivername="Gtiff", overwrite=TRUE)
## Tema: Temperatura de brillo de Landsat 8 TIRS
# Importar las imagenes
B10 = raster("LC08_L1TP_008057_20150104_20170415_01_T1/LC08_L1TP_008057_20150104_20170415_01_T1_B10.TIF")
B11 = raster("LC08_L1TP_008057_20150104_20170415_01_T1/LC08_L1TP_008057_20150104_20170415_01_T1_B11.TIF")
# K imagen Landsat 8 TIRS
K1_CONSTANT_BAND_10 = 774.8853
K2_CONSTANT_BAND_10 = 1321.0789
K1_CONSTANT_BAND_11 = 480.8883
K2_CONSTANT_BAND_11 = 1201.1442
# Datos de radiancia del sensor imagen Landsat
RADIANCE_MULT_BAND_10 = 3.3420E-04
RADIANCE_MULT_BAND_11 = 3.3420E-04
RADIANCE_ADD_BAND_10 = 0.10000
RADIANCE_ADD_BAND_11 = 0.10000
# Elimar valores nulos
B10[B10==0] <- NA
B11[B11==0] <- NA
# Determinancion de la Radianca
L10 <- RADIANCE_MULT_BAND_10*(na.omit(B10)) + RADIANCE_ADD_BAND_10
L11 <- RADIANCE_MULT_BAND_11*(na.omit(B11)) + RADIANCE_ADD_BAND_11
# Determinacion Temperatura de brillo Celsius
TB_B10 <- K2_CONSTANT_BAND_10/log(K1_CONSTANT_BAND_10/L10+1) - 273.15
TB_B11 <- K2_CONSTANT_BAND_11/log(K1_CONSTANT_BAND_11/L11+1) - 273.15
# Exportacion imagen Temperatura de brillo Celsius
writeRaster(TB_B10,"L8_TB_B10.tif", drivername="Gtiff",overwrite=TRUE)
writeRaster(TB_B11,"L8_TB_B11.tif", drivername="Gtiff",overwrite=TRUE)
##Reproyectar y recorte de images con la zona de estudio Landsat 8 OLI
# Importar las imagen
L8_B234567_RS_DOS = stack("L8_B234567_RS_DOS1.tif")
# Verificar sistemas de coordenada
st_crs(L8_B234567_RS_DOS)
Coordinate Reference System:
User input: +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
wkt:
PROJCS["UTM Zone 18, Northern Hemisphere",
GEOGCS["WGS 84",
DATUM["WGS_1984",
SPHEROID["WGS 84",6378137,298.257223563,
AUTHORITY["EPSG","7030"]],
AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich",0,
AUTHORITY["EPSG","8901"]],
UNIT["degree",0.0174532925199433,
AUTHORITY["EPSG","9122"]],
AUTHORITY["EPSG","4326"]],
PROJECTION["Transverse_Mercator"],
PARAMETER["latitude_of_origin",0],
PARAMETER["central_meridian",-75],
PARAMETER["scale_factor",0.9996],
PARAMETER["false_easting",500000],
PARAMETER["false_northing",0],
UNIT["Meter",1]]
# Reproyectar imagen
L8_B234567_RS_DOS_py = projectRaster(L8_B234567_RS_DOS, crs = "+init=epsg:32618")
## importar shapefile
Zona_UTM_18n <- shapefile("rural.shp")
st_crs(Zona_UTM_18n)
Coordinate Reference System:
User input: +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
wkt:
PROJCS["UTM Zone 18, Northern Hemisphere",
GEOGCS["WGS 84",
DATUM["WGS_1984",
SPHEROID["WGS 84",6378137,298.257223563,
AUTHORITY["EPSG","7030"]],
AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich",0,
AUTHORITY["EPSG","8901"]],
UNIT["degree",0.0174532925199433,
AUTHORITY["EPSG","9122"]],
AUTHORITY["EPSG","4326"]],
PROJECTION["Transverse_Mercator"],
PARAMETER["latitude_of_origin",0],
PARAMETER["central_meridian",-75],
PARAMETER["scale_factor",0.9996],
PARAMETER["false_easting",500000],
PARAMETER["false_northing",0],
UNIT["Meter",1]]
# Recortar imagen segun shapefile Landsat 8 extension
clip_L8_B234567_RS_DOS_py <- crop(L8_B234567_RS_DOS_py,Zona_UTM_18n)
# Visualizacion del recorte de la zona de estudio
plotRGB(clip_L8_B234567_RS_DOS_py, 5,4,3, stretch = "lin")
# Guardar archivo por banda recortada y proyectada
L8_SR_B2 <- clip_L8_B234567_RS_DOS_py$L8_B234567_RS_DOS1.1
L8_SR_B3 <- clip_L8_B234567_RS_DOS_py$L8_B234567_RS_DOS1.2
L8_SR_B4 <- clip_L8_B234567_RS_DOS_py$L8_B234567_RS_DOS1.3
L8_SR_B5 <- clip_L8_B234567_RS_DOS_py$L8_B234567_RS_DOS1.4
L8_SR_B6 <- clip_L8_B234567_RS_DOS_py$L8_B234567_RS_DOS1.5
L8_SR_B7 <- clip_L8_B234567_RS_DOS_py$L8_B234567_RS_DOS1.6
# Exportacion imagen reflectancia Superfie DOS1 Proyectado y recortado
writeRaster(clip_L8_B234567_RS_DOS_py,"clip_L8_B234567_RS_DOS_py.tif", drivername="Gtiff",overwrite=TRUE)
writeRaster(L8_SR_B2,"L8_SR_B2.tif", drivername="Gtiff",overwrite=TRUE)
writeRaster(L8_SR_B3,"L8_SR_B3.tif", drivername="Gtiff",overwrite=TRUE)
writeRaster(L8_SR_B4,"L8_SR_B4.tif", drivername="Gtiff",overwrite=TRUE)
writeRaster(L8_SR_B5,"L8_SR_B5.tif", drivername="Gtiff",overwrite=TRUE)
writeRaster(L8_SR_B6,"L8_SR_B6.tif", drivername="Gtiff",overwrite=TRUE)
writeRaster(L8_SR_B7,"L8_SR_B7.tif", drivername="Gtiff",overwrite=TRUE)
##Clasificacion supervisada con Random Forest - Landsat 8 OLI
# Agregar raster multibandas de landsat 8
L8_2015 <- stack("clip_L8_B234567_RS_DOS_py.tif")
# Dimension de la banda
dim(L8_2015)
[1] 617 351 6
# Nombre de las bandas
names(L8_2015)
[1] "clip_L8_B234567_RS_DOS_py.1" "clip_L8_B234567_RS_DOS_py.2" "clip_L8_B234567_RS_DOS_py.3"
[4] "clip_L8_B234567_RS_DOS_py.4" "clip_L8_B234567_RS_DOS_py.5" "clip_L8_B234567_RS_DOS_py.6"
# cambiar nombre de las bandas
names(L8_2015) <- c("BLUE", "GREEN", "RED", "NIR", "SWIR1", "SWIR2")
# Agregar los shapefile del ROI 2019
roi <- shapefile("clasificacion.shp")
# Informacion de la tabla de atribu del shp
print(roi)
class : SpatialPointsDataFrame
features : 2236
extent : 576861.3, 587313.9, 484041, 502444.5 (xmin, xmax, ymin, ymax)
crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
variables : 1
names : Class
min values : Agua
max values : represa
head(roi@data,10)
# Ploteamos la imagen y las puntos de muestreo
plotRGB(L8_2015,5,4,3, stretch="lin")
plot(roi,add=TRUE,col="red",pch = 15, lwd=3)
# Extraer los valores del raster a tabla
Extraer_R_Shp <- extract(L8_2015,roi) # tenemos los valores Reflectancia Superficie
print(Extraer_R_Shp)
BLUE GREEN RED NIR SWIR1 SWIR2
[1,] 0.11948301 0.10713014 0.15300514 0.20682272 0.21726671 0.23338568
[2,] 0.14195211 0.14032544 0.18904631 0.28627995 0.25465634 0.25938857
[3,] 0.13981703 0.12992968 0.17885414 0.20590763 0.22130814 0.25041595
[4,] 0.13353890 0.12339737 0.17298281 0.20509428 0.23986316 0.25702471
[5,] 0.13272555 0.11833928 0.18025205 0.22367497 0.23935480 0.26587027
[6,] 0.11739877 0.10430878 0.15051429 0.18308212 0.20117721 0.22385383
[7,] 0.12504944 0.12171981 0.17222030 0.21627825 0.26657730 0.29423705
[8,] 0.12408359 0.11864429 0.16342606 0.20778859 0.23960897 0.24990755
[9,] 0.16721714 0.18783084 0.25812945 0.30389476 0.33190113 0.31569013
[10,] 0.13750404 0.12835379 0.19438384 0.25173664 0.23488127 0.24190079
[11,] 0.14421427 0.13440315 0.18312417 0.21787962 0.22941643 0.24884000
[12,] 0.11633123 0.11823761 0.16057937 0.19596913 0.21848677 0.23968942
[13,] 0.10212284 0.08895659 0.15107346 0.20140861 0.23887186 0.27967238
[14,] 0.12123682 0.11289991 0.16035064 0.19642667 0.21386071 0.23534288
[15,] 0.12009303 0.11099361 0.15963896 0.19932434 0.20892966 0.23008130
[16,] 0.13386934 0.12621871 0.17130531 0.21518528 0.22204526 0.24863666
[17,] 0.13432685 0.13193767 0.18236166 0.20326416 0.22977228 0.25041595
[18,] 0.14545973 0.13656366 0.20048392 0.23823957 0.26634854 0.28318012
[19,] 0.12446485 0.12393113 0.18406458 0.20867822 0.25341088 0.24474765
[20,] 0.24280888 0.30968258 0.34589398 0.34636855 0.37988999 0.35981625
[21,] 0.11048520 0.10702846 0.14987886 0.21734583 0.20770961 0.22479427
[22,] 0.11920342 0.10662179 0.16085896 0.18979251 0.20781128 0.22842909
[23,] 0.12624407 0.10896020 0.15999480 0.20537385 0.21879178 0.23651212
[24,] 0.11381490 0.10512215 0.15562309 0.19121596 0.20168558 0.22362503
[25,] 0.11955926 0.12797251 0.16243480 0.32898250 0.24881025 0.22329462
[26,] 0.20991859 0.20224260 0.23990552 0.38136935 0.30973679 0.31513089
[27,] 0.15661803 0.14357890 0.18261583 0.40424570 0.26855990 0.20860283
[28,] 0.15061949 0.12476992 0.15282722 0.32120451 0.19456859 0.15461437
[29,] 0.21157074 0.22529630 0.26933828 0.44804126 0.39295477 0.32021454
[30,] 0.13958828 0.11432330 0.13846667 0.37722617 0.18361349 0.14180356
[31,] 0.12118598 0.09121875 0.11948025 0.32567811 0.17522562 0.14208318
[32,] 0.13819031 0.12301610 0.15526725 0.34535182 0.21452157 0.17106001
[33,] 0.13412350 0.10763849 0.13905126 0.30547065 0.20331231 0.16874696
[34,] 0.15516923 0.12014393 0.14441423 0.37150708 0.19339935 0.14965782
[35,] 0.11727168 0.10865519 0.14388047 0.35803548 0.23122109 0.17975307
[36,] 0.16142194 0.13297978 0.15691935 0.38724095 0.23289867 0.21978688
[37,] 0.13142926 0.11592461 0.14108463 0.35373980 0.18661280 0.14383703
[38,] 0.12626949 0.10169078 0.13249370 0.29517630 0.17784365 0.14531130
[39,] 0.10507127 0.08481354 0.11782815 0.27532470 0.23195821 0.19607161
[40,] 0.15125492 0.17141111 0.23273796 0.32562730 0.33741680 0.31307203
[41,] 0.16205738 0.18755126 0.26628825 0.35061336 0.37147668 0.33706689
[42,] 0.21843347 0.24718082 0.30619279 0.34125948 0.33703554 0.32046872
[43,] 0.16254032 0.16856435 0.24547182 0.34059864 0.44353625 0.57050848
[44,] 0.16401453 0.15989697 0.21642032 0.23973925 0.23589797 0.26045614
[45,] 0.14253671 0.13010760 0.17717662 0.19957852 0.18780744 0.22893746
[46,] 0.14304507 0.13447942 0.18332751 0.19581662 0.20254979 0.23991817
[47,] 0.28368029 0.23523457 0.25718901 0.33658251 0.44587469 0.45945597
[48,] 0.18373853 0.20249678 0.26745746 0.27270666 0.27399930 0.28401890
[49,] 0.12586281 0.13539445 0.17737995 0.30445394 0.22743383 0.22611605
[50,] 0.10311412 0.09889485 0.14967553 0.18333629 0.16104247 0.17914303
[51,] 0.19512559 0.27142915 0.38285017 0.46720660 0.52271277 0.50396347
[52,] 0.16360785 0.19011843 0.26440740 0.29624388 0.30434820 0.31718978
[53,] 0.16198112 0.13984253 0.17209323 0.42310601 0.23724513 0.17873633
[54,] 0.09597179 0.07205392 0.11299895 0.11862161 0.11427364 0.15237758
[55,] 0.13897826 0.14002044 0.17961663 0.18514098 0.21759714 0.23178431
[56,] 0.17329192 0.15468636 0.19118132 0.39555272 0.23205988 0.17899053
[57,] 0.11007852 0.08056880 0.10674638 0.29260907 0.15628932 0.12777266
[58,] 0.11727168 0.10865519 0.14388047 0.35803548 0.23122109 0.17975307
[59,] 0.16213363 0.13577572 0.16606943 0.39677274 0.25541887 0.20252787
[60,] 0.17878212 0.13859706 0.16291772 0.41131201 0.21248816 0.16239238
[61,] 0.14891651 0.14230803 0.19268093 0.22835191 0.24748851 0.27542755
[62,] 0.14218086 0.13468276 0.18815672 0.21653244 0.23127194 0.24878915
[63,] 0.13430142 0.13125138 0.17290658 0.26149720 0.25003028 0.25738055
[64,] 0.13539438 0.13033636 0.19712889 0.29044852 0.29466400 0.31660518
[65,] 0.13023463 0.13445400 0.17966747 0.32888082 0.29550281 0.28295133
[66,] 0.12441401 0.11714465 0.16932280 0.24756806 0.25613058 0.25849897
[67,] 0.15984605 0.15608433 0.22114785 0.32916039 0.27588022 0.26648030
[68,] 0.16154903 0.11717007 0.13012993 0.35602742 0.17387846 0.13636403
[69,] 0.11228985 0.10006406 0.12494488 0.33459991 0.22740841 0.18694644
[70,] 0.12626949 0.11198489 0.14967553 0.26808050 0.22977228 0.19706292
[71,] 0.12184684 0.10024198 0.12573281 0.37549773 0.20875172 0.16086729
[72,] 0.12583739 0.10191954 0.12807116 0.36266160 0.21335237 0.16795899
[73,] 0.18533984 0.17822301 0.19900973 0.47473034 0.29303727 0.21536411
[74,] 0.15201746 0.12977718 0.15608059 0.30468270 0.21035306 0.18384542
[75,] 0.58754689 0.50305927 0.51372182 0.88973302 0.84374005 1.00053525
[76,] 0.65017569 0.66074973 0.83089924 1.22665107 1.43452680 1.43663681
[77,] 0.33850589 0.34951183 0.38991606 0.62424010 0.61609793 0.76833904
[78,] 0.12365148 0.10906187 0.13790751 0.32397512 0.20483738 0.16945866
[79,] 0.13516563 0.10715555 0.13691625 0.37171045 0.20387150 0.15921509
[80,] 0.11770378 0.10463922 0.13711958 0.34166616 0.24682765 0.20499343
[81,] 0.10763844 0.08438143 0.10753432 0.37076998 0.19253515 0.14681096
[82,] 0.14706104 0.11887304 0.13783126 0.37132919 0.20455779 0.15763915
[83,] 0.11439950 0.09576850 0.12413154 0.35806090 0.21520786 0.16376497
[84,] 0.11930509 0.08979537 0.11399020 0.35859463 0.19652575 0.14897153
[85,] 0.12253311 0.10298707 0.13107036 0.33724341 0.20994636 0.16905196
[86,] 0.11183233 0.08519480 0.11269394 0.27545178 0.17748781 0.14879359
[87,] 0.12367690 0.09429427 0.12392822 0.21366020 0.11897594 0.12881482
[88,] 0.10367330 0.08478811 0.11500689 0.28091669 0.15034156 0.14160022
[89,] 0.10446125 0.08059422 0.10832224 0.22202279 0.12093312 0.12604423
[90,] 0.12728620 0.10204662 0.13124827 0.26452199 0.17527644 0.15850338
[91,] 0.07670530 0.05700673 0.10184093 0.13509259 0.11572246 0.12792517
[92,] 0.07843369 0.06353904 0.10639055 0.18237041 0.13476042 0.13827042
[93,] 0.08206840 0.06803795 0.10377261 0.18089615 0.12395784 0.13219544
[94,] 0.11376406 0.10057241 0.15244599 0.22982617 0.21947806 0.21180555
[95,] 0.09899648 0.08178885 0.12120862 0.19009753 0.15867861 0.17299181
[96,] 0.09770019 0.08438143 0.13401872 0.22560674 0.17779282 0.16729811
[97,] 0.07579027 0.05110987 0.09256377 0.19792633 0.12459329 0.12754390
[98,] 0.10591005 0.09274381 0.12809658 0.21790501 0.14233494 0.14289655
[99,] 0.09988609 0.08936327 0.13734832 0.31914565 0.19418731 0.15814753
[100,] 0.08852445 0.07309604 0.10102759 0.17904063 0.13216780 0.15441102
[101,] 0.09004951 0.06292903 0.10046842 0.16663656 0.11747629 0.13130581
[102,] 0.09841187 0.08837198 0.12911326 0.20430629 0.17677610 0.18148151
[103,] 0.11157816 0.12591371 0.17407575 0.32816908 0.26741609 0.23272480
[104,] 0.09581929 0.10395294 0.15038720 0.34865615 0.28792828 0.23160641
[105,] 0.11203567 0.11978807 0.16558652 0.32183999 0.25549513 0.22728528
[106,] 0.15649094 0.17102985 0.25866321 0.37837005 0.39193806 0.34931850
[107,] 0.08610979 0.06236984 0.10435720 0.10301483 0.08346723 0.12640007
[108,] 0.11068854 0.10850269 0.15361516 0.15918903 0.14279246 0.16948406
[109,] 0.12548155 0.10049617 0.12524989 0.30961385 0.16571933 0.13257672
[110,] 0.07924705 0.06137855 0.11251602 0.22751309 0.32732594 0.27125892
[111,] 0.10522377 0.08549980 0.12443656 0.19474904 0.18208842 0.17759253
[112,] 0.10242785 0.09071040 0.13211246 0.18280253 0.15814483 0.16701850
[113,] 0.10044528 0.08288179 0.12166610 0.16096829 0.14012359 0.15784252
[114,] 0.05324491 0.03197045 0.05827636 0.22949573 0.09927713 0.09706736
[115,] 0.05192320 0.02889493 0.05560759 0.19032629 0.08644112 0.09239040
[116,] 0.05360076 0.03382593 0.05921679 0.23295259 0.10990179 0.10115971
[117,] 0.05403286 0.03158919 0.05908971 0.24525499 0.10736000 0.10235438
[118,] 0.05482080 0.03402928 0.05835262 0.24411118 0.10349648 0.09894831
[119,] 0.05403286 0.03186878 0.05860678 0.24395867 0.10059886 0.09401717
[120,] 0.04800890 0.01674534 0.04627959 0.06496382 0.02904765 0.06676877
[121,] 0.04841559 0.01613532 0.04615251 0.05639789 0.02452328 0.06547242
[122,] 0.05047441 0.02238804 0.05034629 0.12154470 0.05250832 0.07777487
[123,] 0.05301616 0.02930161 0.05581092 0.20316248 0.09465108 0.09475429
[124,] 0.05459205 0.03265673 0.05924221 0.20517051 0.09828583 0.09849078
[125,] 0.05245697 0.02782739 0.05649718 0.18511556 0.08659363 0.09302586
[126,] 0.05367702 0.03146210 0.05847970 0.23008035 0.10326773 0.10027008
[127,] 0.05227905 0.02838658 0.05555676 0.19858722 0.08890665 0.09498306
[128,] 0.05421078 0.03240255 0.05954721 0.23803625 0.10880881 0.10080384
[129,] 0.05601543 0.03659645 0.06178389 0.27855280 0.12248360 0.10708217
[130,] 0.04994064 0.02015129 0.04912628 0.08880607 0.04168032 0.07299624
[131,] 0.04844100 0.01641491 0.04661001 0.05840593 0.02594668 0.06618414
[132,] 0.04889852 0.01687243 0.04630500 0.06219324 0.02759884 0.06628581
[133,] 0.05332117 0.03105542 0.05731052 0.24139144 0.08173882 0.08715423
[134,] 0.05451579 0.03374968 0.05896262 0.23102081 0.09981090 0.09744864
[135,] 0.05271115 0.03075041 0.05743761 0.23623154 0.10049719 0.09640648
[136,] 0.05375327 0.03151293 0.05832720 0.23651116 0.10390317 0.10334568
[137,] 0.05474455 0.03496972 0.06188557 0.27163908 0.12667754 0.10797182
[138,] 0.05476997 0.03636769 0.06094513 0.24047638 0.11142685 0.10263397
[139,] 0.05395661 0.03245339 0.05901346 0.22985159 0.10919008 0.10382863
[140,] 0.05253322 0.02851367 0.05509925 0.18831825 0.08326389 0.08989941
[141,] 0.05016939 0.02780198 0.05710719 0.23313053 0.08687323 0.08934020
[142,] 0.04968646 0.02701403 0.05603968 0.21104212 0.08417894 0.08893352
[143,] 0.05192320 0.02874242 0.05761553 0.18145534 0.07622316 0.08756092
[144,] 0.04902561 0.02439602 0.05393007 0.17235565 0.08644112 0.09165327
[145,] 0.04983897 0.02355725 0.05260840 0.19022465 0.08600902 0.09813493
[146,] 0.05606626 0.03710480 0.06328350 0.29479501 0.12334782 0.10720927
[147,] 0.05149111 0.03019123 0.05886095 0.19881597 0.08369599 0.08771344
[148,] 0.04940687 0.02632776 0.05535342 0.19042797 0.08611069 0.09211080
[149,] 0.05205029 0.03280923 0.05769177 0.25239751 0.11834050 0.10344736
[150,] 0.04851726 0.02365891 0.05326923 0.16671281 0.06613228 0.08362109
[151,] 0.04818683 0.02335390 0.05212548 0.15466458 0.06191291 0.08153679
[152,] 0.05354993 0.03136043 0.05903887 0.20707689 0.09576946 0.09683859
[153,] 0.04811057 0.02129508 0.05082922 0.13732938 0.06142997 0.08146054
[154,] 0.05116068 0.02731904 0.05481967 0.18211621 0.08550067 0.09277167
[155,] 0.04887310 0.02391309 0.05202381 0.17624463 0.08730532 0.09396634
[156,] 0.05301616 0.03240255 0.05601426 0.26137012 0.10891049 0.10436241
[157,] 0.05446496 0.03324133 0.05959805 0.22019267 0.09188054 0.09404259
[158,] 0.05316867 0.03021664 0.05794595 0.20651768 0.09099092 0.09361047
[159,] 0.05380410 0.03268214 0.05881012 0.21477859 0.09060965 0.09305128
[160,] 0.05245697 0.02675986 0.05438758 0.15311408 0.06625936 0.08629001
[161,] 0.04864434 0.01631324 0.04663542 0.06257451 0.02838679 0.06681960
[162,] 0.05217738 0.02658193 0.05311674 0.14525986 0.06674231 0.08369735
[163,] 0.05421078 0.03087750 0.05542967 0.20438254 0.08450937 0.09066195
[164,] 0.04894935 0.01966836 0.04955837 0.09747367 0.04345956 0.07241162
[165,] 0.05334659 0.02942870 0.05624301 0.18928415 0.08684781 0.09068737
[166,] 0.04971188 0.02297264 0.05169339 0.13097484 0.06081995 0.08174013
[ reached getOption("max.print") -- omitted 2070 rows ]
# Seleccionar la columna de las clases
Class <- roi@data$Class
print(Class)
[1] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[8] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[15] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[22] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[29] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[36] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[43] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[50] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[57] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[64] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[71] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[78] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[85] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[92] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[99] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[106] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[113] "Artificial" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[120] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[127] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[134] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[141] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[148] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[155] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[162] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[169] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[176] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[183] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[190] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[197] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[204] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[211] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[218] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[225] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[232] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[239] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[246] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[253] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[260] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[267] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[274] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[281] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[288] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[295] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[302] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[309] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[316] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[323] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[330] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[337] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[344] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[351] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[358] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[365] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[372] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[379] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[386] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[393] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[400] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[407] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[414] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[421] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[428] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[435] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[442] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[449] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[456] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[463] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[470] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[477] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[484] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[491] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[498] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[505] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[512] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[519] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[526] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[533] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[540] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[547] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[554] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[561] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[568] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[575] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[582] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[589] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[596] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[603] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[610] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[617] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[624] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[631] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[638] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[645] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[652] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[659] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[666] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[673] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[680] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[687] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[694] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[701] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[708] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[715] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[722] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[729] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[736] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[743] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[750] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[757] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[764] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[771] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[778] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[785] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[792] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[799] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[806] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[813] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[820] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[827] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[834] "bosque" "bosque" "bosque" "bosque" "represa" "represa" "represa"
[841] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[848] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[855] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[862] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[869] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[876] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[883] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[890] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[897] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[904] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[911] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[918] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[925] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[932] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[939] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[946] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[953] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[960] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[967] "represa" "represa" "represa" "Agua" "Agua" "Agua" "Agua"
[974] "Agua" "Agua" "Agua" "Agua" "pastos" "pastos" "pastos"
[981] "pastos" "pastos" "pastos" "pastos" "pastos" "pastos" "pastos"
[988] "pastos" "pastos" "pastos" "pastos" "pastos" "pastos" "pastos"
[995] "pastos" "pastos" "pastos" "pastos" "pastos" "pastos"
[ reached getOption("max.print") -- omitted 1236 entries ]
# Crear un data frame con los datos reflectancia superficie y las clases
tabla <- data.frame(Extraer_R_Shp, Class)
summary(tabla$Class)
Agua Artificial bosque cultivos Desnudo pastos represa
8 113 724 747 23 489 132
# Seperando los datos en Entrenamiento y prueba
Muestra <- sample(1:2235,1000) # sample separa aleatoriamente los datos
Pruebas <- tabla[Muestra,]
Entrenamiento <- tabla[-Muestra,]
# Aplicando el algoritmo Random Forest
modelo <- randomForest(Class ~.,data=Entrenamiento, importance=TRUE)
prediccion<- predict(modelo, Pruebas[,-7])
# Definimos la paleta de colores segun las 7 clases
mycolor <- c("#00FF00","#0000FF","#228B22","#FF1493", "#1C1CE2", "#1C1CE2" ,"#5454AA")
## Visualizar prediccion a obtener
## Cantidad de punto que se encuentra en la clase
plot(prediccion, main = "Numero de pixeles identificados",
axes = TRUE, xlab = "Clases", ylab = "Numero de pixel", col= mycolor)
# Matriz de confusion
MC<-table(Pruebas[,7],prediccion)
MC
prediccion
Agua Artificial bosque cultivos Desnudo pastos represa
Agua 0 0 3 0 0 0 1
Artificial 0 53 0 2 1 1 0
bosque 0 0 326 1 0 1 0
cultivos 0 0 3 322 0 0 0
Desnudo 0 1 0 7 6 0 0
pastos 0 0 0 1 0 213 0
represa 0 0 0 0 0 0 58
# Indice de Kappa
kappa<-sum(diag(MC))/sum(MC)
kappa
[1] 0.978
# Ejecutamos el clasificador Random Forest
beginCluster()
8 cores detected, using 7
## 8 cores detected, using 7
rf_class<- clusterR(L8_2015, raster::predict,args = list(model =modelo))
endCluster()
# Ploteamos la clasificación con Support vector machine
plot(rf_class,main=paste("kappa = ",kappa,sep =""),col =mycolor,cex.lab=0.8,
cex.axis=0.8,cex.main=0.9)
# Exportacion de bandas
writeRaster(rf_class,"Class_L8_Random_Forest.tif", drivername="Gtiff",overwrite=TRUE)
## Clasificacion supervisada con Decision_Tree - Landsat 8 OLI
# Agregar raster multibandas de landsat 8
L8_2015 <- stack("clip_L8_B234567_RS_DOS_py.tif")
# Dimension de la banda
dim(L8_2015)
[1] 617 351 6
# Nombre de las bandas
names(L8_2015)
[1] "clip_L8_B234567_RS_DOS_py.1" "clip_L8_B234567_RS_DOS_py.2" "clip_L8_B234567_RS_DOS_py.3"
[4] "clip_L8_B234567_RS_DOS_py.4" "clip_L8_B234567_RS_DOS_py.5" "clip_L8_B234567_RS_DOS_py.6"
# cambiar nombre de las bandas
names(L8_2015) <- c("BLUE", "GREEN", "RED", "NIR", "SWIR1", "SWIR2")
# Agregar los shapefile del ROI 2019
roi <- shapefile("clasificacion.shp")
Warning in .Internal(gc(verbose, reset, full)) :
closing unused connection 9 (<-DESKTOP-QM41LE8:11036)
Warning in .Internal(gc(verbose, reset, full)) :
closing unused connection 8 (<-DESKTOP-QM41LE8:11036)
Warning in .Internal(gc(verbose, reset, full)) :
closing unused connection 7 (<-DESKTOP-QM41LE8:11036)
Warning in .Internal(gc(verbose, reset, full)) :
closing unused connection 6 (<-DESKTOP-QM41LE8:11036)
Warning in .Internal(gc(verbose, reset, full)) :
closing unused connection 5 (<-DESKTOP-QM41LE8:11036)
Warning in .Internal(gc(verbose, reset, full)) :
closing unused connection 4 (<-DESKTOP-QM41LE8:11036)
Warning in .Internal(gc(verbose, reset, full)) :
closing unused connection 3 (<-DESKTOP-QM41LE8:11036)
# Informacion de la tabla de atribu del shp
print(roi)
class : SpatialPointsDataFrame
features : 2236
extent : 576861.3, 587313.9, 484041, 502444.5 (xmin, xmax, ymin, ymax)
crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
variables : 1
names : Class
min values : Agua
max values : represa
head(roi@data,10)
# Ploteamos la imagen y las puntos de muestreo
plotRGB(L8_2015,5,4,3, stretch="lin")
plot(roi,add=TRUE,col="red",pch = 15, lwd=3)
# Extraer los valores del raster a tabla
Extraer_R_Shp <- extract(L8_2015,roi) # tenemos los valores Reflectancia Superficie
print(Extraer_R_Shp)
BLUE GREEN RED NIR SWIR1 SWIR2
[1,] 0.11948301 0.10713014 0.15300514 0.20682272 0.21726671 0.23338568
[2,] 0.14195211 0.14032544 0.18904631 0.28627995 0.25465634 0.25938857
[3,] 0.13981703 0.12992968 0.17885414 0.20590763 0.22130814 0.25041595
[4,] 0.13353890 0.12339737 0.17298281 0.20509428 0.23986316 0.25702471
[5,] 0.13272555 0.11833928 0.18025205 0.22367497 0.23935480 0.26587027
[6,] 0.11739877 0.10430878 0.15051429 0.18308212 0.20117721 0.22385383
[7,] 0.12504944 0.12171981 0.17222030 0.21627825 0.26657730 0.29423705
[8,] 0.12408359 0.11864429 0.16342606 0.20778859 0.23960897 0.24990755
[9,] 0.16721714 0.18783084 0.25812945 0.30389476 0.33190113 0.31569013
[10,] 0.13750404 0.12835379 0.19438384 0.25173664 0.23488127 0.24190079
[11,] 0.14421427 0.13440315 0.18312417 0.21787962 0.22941643 0.24884000
[12,] 0.11633123 0.11823761 0.16057937 0.19596913 0.21848677 0.23968942
[13,] 0.10212284 0.08895659 0.15107346 0.20140861 0.23887186 0.27967238
[14,] 0.12123682 0.11289991 0.16035064 0.19642667 0.21386071 0.23534288
[15,] 0.12009303 0.11099361 0.15963896 0.19932434 0.20892966 0.23008130
[16,] 0.13386934 0.12621871 0.17130531 0.21518528 0.22204526 0.24863666
[17,] 0.13432685 0.13193767 0.18236166 0.20326416 0.22977228 0.25041595
[18,] 0.14545973 0.13656366 0.20048392 0.23823957 0.26634854 0.28318012
[19,] 0.12446485 0.12393113 0.18406458 0.20867822 0.25341088 0.24474765
[20,] 0.24280888 0.30968258 0.34589398 0.34636855 0.37988999 0.35981625
[21,] 0.11048520 0.10702846 0.14987886 0.21734583 0.20770961 0.22479427
[22,] 0.11920342 0.10662179 0.16085896 0.18979251 0.20781128 0.22842909
[23,] 0.12624407 0.10896020 0.15999480 0.20537385 0.21879178 0.23651212
[24,] 0.11381490 0.10512215 0.15562309 0.19121596 0.20168558 0.22362503
[25,] 0.11955926 0.12797251 0.16243480 0.32898250 0.24881025 0.22329462
[26,] 0.20991859 0.20224260 0.23990552 0.38136935 0.30973679 0.31513089
[27,] 0.15661803 0.14357890 0.18261583 0.40424570 0.26855990 0.20860283
[28,] 0.15061949 0.12476992 0.15282722 0.32120451 0.19456859 0.15461437
[29,] 0.21157074 0.22529630 0.26933828 0.44804126 0.39295477 0.32021454
[30,] 0.13958828 0.11432330 0.13846667 0.37722617 0.18361349 0.14180356
[31,] 0.12118598 0.09121875 0.11948025 0.32567811 0.17522562 0.14208318
[32,] 0.13819031 0.12301610 0.15526725 0.34535182 0.21452157 0.17106001
[33,] 0.13412350 0.10763849 0.13905126 0.30547065 0.20331231 0.16874696
[34,] 0.15516923 0.12014393 0.14441423 0.37150708 0.19339935 0.14965782
[35,] 0.11727168 0.10865519 0.14388047 0.35803548 0.23122109 0.17975307
[36,] 0.16142194 0.13297978 0.15691935 0.38724095 0.23289867 0.21978688
[37,] 0.13142926 0.11592461 0.14108463 0.35373980 0.18661280 0.14383703
[38,] 0.12626949 0.10169078 0.13249370 0.29517630 0.17784365 0.14531130
[39,] 0.10507127 0.08481354 0.11782815 0.27532470 0.23195821 0.19607161
[40,] 0.15125492 0.17141111 0.23273796 0.32562730 0.33741680 0.31307203
[41,] 0.16205738 0.18755126 0.26628825 0.35061336 0.37147668 0.33706689
[42,] 0.21843347 0.24718082 0.30619279 0.34125948 0.33703554 0.32046872
[43,] 0.16254032 0.16856435 0.24547182 0.34059864 0.44353625 0.57050848
[44,] 0.16401453 0.15989697 0.21642032 0.23973925 0.23589797 0.26045614
[45,] 0.14253671 0.13010760 0.17717662 0.19957852 0.18780744 0.22893746
[46,] 0.14304507 0.13447942 0.18332751 0.19581662 0.20254979 0.23991817
[47,] 0.28368029 0.23523457 0.25718901 0.33658251 0.44587469 0.45945597
[48,] 0.18373853 0.20249678 0.26745746 0.27270666 0.27399930 0.28401890
[49,] 0.12586281 0.13539445 0.17737995 0.30445394 0.22743383 0.22611605
[50,] 0.10311412 0.09889485 0.14967553 0.18333629 0.16104247 0.17914303
[51,] 0.19512559 0.27142915 0.38285017 0.46720660 0.52271277 0.50396347
[52,] 0.16360785 0.19011843 0.26440740 0.29624388 0.30434820 0.31718978
[53,] 0.16198112 0.13984253 0.17209323 0.42310601 0.23724513 0.17873633
[54,] 0.09597179 0.07205392 0.11299895 0.11862161 0.11427364 0.15237758
[55,] 0.13897826 0.14002044 0.17961663 0.18514098 0.21759714 0.23178431
[56,] 0.17329192 0.15468636 0.19118132 0.39555272 0.23205988 0.17899053
[57,] 0.11007852 0.08056880 0.10674638 0.29260907 0.15628932 0.12777266
[58,] 0.11727168 0.10865519 0.14388047 0.35803548 0.23122109 0.17975307
[59,] 0.16213363 0.13577572 0.16606943 0.39677274 0.25541887 0.20252787
[60,] 0.17878212 0.13859706 0.16291772 0.41131201 0.21248816 0.16239238
[61,] 0.14891651 0.14230803 0.19268093 0.22835191 0.24748851 0.27542755
[62,] 0.14218086 0.13468276 0.18815672 0.21653244 0.23127194 0.24878915
[63,] 0.13430142 0.13125138 0.17290658 0.26149720 0.25003028 0.25738055
[64,] 0.13539438 0.13033636 0.19712889 0.29044852 0.29466400 0.31660518
[65,] 0.13023463 0.13445400 0.17966747 0.32888082 0.29550281 0.28295133
[66,] 0.12441401 0.11714465 0.16932280 0.24756806 0.25613058 0.25849897
[67,] 0.15984605 0.15608433 0.22114785 0.32916039 0.27588022 0.26648030
[68,] 0.16154903 0.11717007 0.13012993 0.35602742 0.17387846 0.13636403
[69,] 0.11228985 0.10006406 0.12494488 0.33459991 0.22740841 0.18694644
[70,] 0.12626949 0.11198489 0.14967553 0.26808050 0.22977228 0.19706292
[71,] 0.12184684 0.10024198 0.12573281 0.37549773 0.20875172 0.16086729
[72,] 0.12583739 0.10191954 0.12807116 0.36266160 0.21335237 0.16795899
[73,] 0.18533984 0.17822301 0.19900973 0.47473034 0.29303727 0.21536411
[74,] 0.15201746 0.12977718 0.15608059 0.30468270 0.21035306 0.18384542
[75,] 0.58754689 0.50305927 0.51372182 0.88973302 0.84374005 1.00053525
[76,] 0.65017569 0.66074973 0.83089924 1.22665107 1.43452680 1.43663681
[77,] 0.33850589 0.34951183 0.38991606 0.62424010 0.61609793 0.76833904
[78,] 0.12365148 0.10906187 0.13790751 0.32397512 0.20483738 0.16945866
[79,] 0.13516563 0.10715555 0.13691625 0.37171045 0.20387150 0.15921509
[80,] 0.11770378 0.10463922 0.13711958 0.34166616 0.24682765 0.20499343
[81,] 0.10763844 0.08438143 0.10753432 0.37076998 0.19253515 0.14681096
[82,] 0.14706104 0.11887304 0.13783126 0.37132919 0.20455779 0.15763915
[83,] 0.11439950 0.09576850 0.12413154 0.35806090 0.21520786 0.16376497
[84,] 0.11930509 0.08979537 0.11399020 0.35859463 0.19652575 0.14897153
[85,] 0.12253311 0.10298707 0.13107036 0.33724341 0.20994636 0.16905196
[86,] 0.11183233 0.08519480 0.11269394 0.27545178 0.17748781 0.14879359
[87,] 0.12367690 0.09429427 0.12392822 0.21366020 0.11897594 0.12881482
[88,] 0.10367330 0.08478811 0.11500689 0.28091669 0.15034156 0.14160022
[89,] 0.10446125 0.08059422 0.10832224 0.22202279 0.12093312 0.12604423
[90,] 0.12728620 0.10204662 0.13124827 0.26452199 0.17527644 0.15850338
[91,] 0.07670530 0.05700673 0.10184093 0.13509259 0.11572246 0.12792517
[92,] 0.07843369 0.06353904 0.10639055 0.18237041 0.13476042 0.13827042
[93,] 0.08206840 0.06803795 0.10377261 0.18089615 0.12395784 0.13219544
[94,] 0.11376406 0.10057241 0.15244599 0.22982617 0.21947806 0.21180555
[95,] 0.09899648 0.08178885 0.12120862 0.19009753 0.15867861 0.17299181
[96,] 0.09770019 0.08438143 0.13401872 0.22560674 0.17779282 0.16729811
[97,] 0.07579027 0.05110987 0.09256377 0.19792633 0.12459329 0.12754390
[98,] 0.10591005 0.09274381 0.12809658 0.21790501 0.14233494 0.14289655
[99,] 0.09988609 0.08936327 0.13734832 0.31914565 0.19418731 0.15814753
[100,] 0.08852445 0.07309604 0.10102759 0.17904063 0.13216780 0.15441102
[101,] 0.09004951 0.06292903 0.10046842 0.16663656 0.11747629 0.13130581
[102,] 0.09841187 0.08837198 0.12911326 0.20430629 0.17677610 0.18148151
[103,] 0.11157816 0.12591371 0.17407575 0.32816908 0.26741609 0.23272480
[104,] 0.09581929 0.10395294 0.15038720 0.34865615 0.28792828 0.23160641
[105,] 0.11203567 0.11978807 0.16558652 0.32183999 0.25549513 0.22728528
[106,] 0.15649094 0.17102985 0.25866321 0.37837005 0.39193806 0.34931850
[107,] 0.08610979 0.06236984 0.10435720 0.10301483 0.08346723 0.12640007
[108,] 0.11068854 0.10850269 0.15361516 0.15918903 0.14279246 0.16948406
[109,] 0.12548155 0.10049617 0.12524989 0.30961385 0.16571933 0.13257672
[110,] 0.07924705 0.06137855 0.11251602 0.22751309 0.32732594 0.27125892
[111,] 0.10522377 0.08549980 0.12443656 0.19474904 0.18208842 0.17759253
[112,] 0.10242785 0.09071040 0.13211246 0.18280253 0.15814483 0.16701850
[113,] 0.10044528 0.08288179 0.12166610 0.16096829 0.14012359 0.15784252
[114,] 0.05324491 0.03197045 0.05827636 0.22949573 0.09927713 0.09706736
[115,] 0.05192320 0.02889493 0.05560759 0.19032629 0.08644112 0.09239040
[116,] 0.05360076 0.03382593 0.05921679 0.23295259 0.10990179 0.10115971
[117,] 0.05403286 0.03158919 0.05908971 0.24525499 0.10736000 0.10235438
[118,] 0.05482080 0.03402928 0.05835262 0.24411118 0.10349648 0.09894831
[119,] 0.05403286 0.03186878 0.05860678 0.24395867 0.10059886 0.09401717
[120,] 0.04800890 0.01674534 0.04627959 0.06496382 0.02904765 0.06676877
[121,] 0.04841559 0.01613532 0.04615251 0.05639789 0.02452328 0.06547242
[122,] 0.05047441 0.02238804 0.05034629 0.12154470 0.05250832 0.07777487
[123,] 0.05301616 0.02930161 0.05581092 0.20316248 0.09465108 0.09475429
[124,] 0.05459205 0.03265673 0.05924221 0.20517051 0.09828583 0.09849078
[125,] 0.05245697 0.02782739 0.05649718 0.18511556 0.08659363 0.09302586
[126,] 0.05367702 0.03146210 0.05847970 0.23008035 0.10326773 0.10027008
[127,] 0.05227905 0.02838658 0.05555676 0.19858722 0.08890665 0.09498306
[128,] 0.05421078 0.03240255 0.05954721 0.23803625 0.10880881 0.10080384
[129,] 0.05601543 0.03659645 0.06178389 0.27855280 0.12248360 0.10708217
[130,] 0.04994064 0.02015129 0.04912628 0.08880607 0.04168032 0.07299624
[131,] 0.04844100 0.01641491 0.04661001 0.05840593 0.02594668 0.06618414
[132,] 0.04889852 0.01687243 0.04630500 0.06219324 0.02759884 0.06628581
[133,] 0.05332117 0.03105542 0.05731052 0.24139144 0.08173882 0.08715423
[134,] 0.05451579 0.03374968 0.05896262 0.23102081 0.09981090 0.09744864
[135,] 0.05271115 0.03075041 0.05743761 0.23623154 0.10049719 0.09640648
[136,] 0.05375327 0.03151293 0.05832720 0.23651116 0.10390317 0.10334568
[137,] 0.05474455 0.03496972 0.06188557 0.27163908 0.12667754 0.10797182
[138,] 0.05476997 0.03636769 0.06094513 0.24047638 0.11142685 0.10263397
[139,] 0.05395661 0.03245339 0.05901346 0.22985159 0.10919008 0.10382863
[140,] 0.05253322 0.02851367 0.05509925 0.18831825 0.08326389 0.08989941
[141,] 0.05016939 0.02780198 0.05710719 0.23313053 0.08687323 0.08934020
[142,] 0.04968646 0.02701403 0.05603968 0.21104212 0.08417894 0.08893352
[143,] 0.05192320 0.02874242 0.05761553 0.18145534 0.07622316 0.08756092
[144,] 0.04902561 0.02439602 0.05393007 0.17235565 0.08644112 0.09165327
[145,] 0.04983897 0.02355725 0.05260840 0.19022465 0.08600902 0.09813493
[146,] 0.05606626 0.03710480 0.06328350 0.29479501 0.12334782 0.10720927
[147,] 0.05149111 0.03019123 0.05886095 0.19881597 0.08369599 0.08771344
[148,] 0.04940687 0.02632776 0.05535342 0.19042797 0.08611069 0.09211080
[149,] 0.05205029 0.03280923 0.05769177 0.25239751 0.11834050 0.10344736
[150,] 0.04851726 0.02365891 0.05326923 0.16671281 0.06613228 0.08362109
[151,] 0.04818683 0.02335390 0.05212548 0.15466458 0.06191291 0.08153679
[152,] 0.05354993 0.03136043 0.05903887 0.20707689 0.09576946 0.09683859
[153,] 0.04811057 0.02129508 0.05082922 0.13732938 0.06142997 0.08146054
[154,] 0.05116068 0.02731904 0.05481967 0.18211621 0.08550067 0.09277167
[155,] 0.04887310 0.02391309 0.05202381 0.17624463 0.08730532 0.09396634
[156,] 0.05301616 0.03240255 0.05601426 0.26137012 0.10891049 0.10436241
[157,] 0.05446496 0.03324133 0.05959805 0.22019267 0.09188054 0.09404259
[158,] 0.05316867 0.03021664 0.05794595 0.20651768 0.09099092 0.09361047
[159,] 0.05380410 0.03268214 0.05881012 0.21477859 0.09060965 0.09305128
[160,] 0.05245697 0.02675986 0.05438758 0.15311408 0.06625936 0.08629001
[161,] 0.04864434 0.01631324 0.04663542 0.06257451 0.02838679 0.06681960
[162,] 0.05217738 0.02658193 0.05311674 0.14525986 0.06674231 0.08369735
[163,] 0.05421078 0.03087750 0.05542967 0.20438254 0.08450937 0.09066195
[164,] 0.04894935 0.01966836 0.04955837 0.09747367 0.04345956 0.07241162
[165,] 0.05334659 0.02942870 0.05624301 0.18928415 0.08684781 0.09068737
[166,] 0.04971188 0.02297264 0.05169339 0.13097484 0.06081995 0.08174013
[ reached getOption("max.print") -- omitted 2070 rows ]
# Seleccionar la columna de las clases
Class <- roi@data$Class
print(Class)
[1] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[8] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[15] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[22] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[29] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[36] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[43] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[50] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[57] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[64] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[71] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[78] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[85] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[92] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[99] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[106] "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial" "Artificial"
[113] "Artificial" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[120] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[127] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[134] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[141] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[148] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[155] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[162] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[169] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[176] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[183] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[190] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[197] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[204] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[211] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[218] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[225] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[232] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[239] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[246] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[253] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[260] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[267] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[274] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[281] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[288] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[295] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[302] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[309] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[316] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[323] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[330] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[337] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[344] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[351] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[358] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[365] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[372] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[379] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[386] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[393] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[400] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[407] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[414] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[421] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[428] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[435] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[442] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[449] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[456] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[463] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[470] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[477] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[484] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[491] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[498] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[505] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[512] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[519] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[526] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[533] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[540] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[547] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[554] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[561] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[568] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[575] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[582] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[589] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[596] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[603] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[610] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[617] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[624] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[631] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[638] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[645] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[652] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[659] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[666] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[673] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[680] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[687] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[694] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[701] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[708] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[715] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[722] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[729] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[736] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[743] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[750] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[757] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[764] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[771] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[778] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[785] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[792] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[799] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[806] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[813] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[820] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[827] "bosque" "bosque" "bosque" "bosque" "bosque" "bosque" "bosque"
[834] "bosque" "bosque" "bosque" "bosque" "represa" "represa" "represa"
[841] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[848] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[855] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[862] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[869] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[876] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[883] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[890] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[897] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[904] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[911] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[918] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[925] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[932] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[939] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[946] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[953] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[960] "represa" "represa" "represa" "represa" "represa" "represa" "represa"
[967] "represa" "represa" "represa" "Agua" "Agua" "Agua" "Agua"
[974] "Agua" "Agua" "Agua" "Agua" "pastos" "pastos" "pastos"
[981] "pastos" "pastos" "pastos" "pastos" "pastos" "pastos" "pastos"
[988] "pastos" "pastos" "pastos" "pastos" "pastos" "pastos" "pastos"
[995] "pastos" "pastos" "pastos" "pastos" "pastos" "pastos"
[ reached getOption("max.print") -- omitted 1236 entries ]
# Crear un data frame con los datos reflectancia superficie y las clases
tabla <- data.frame(Extraer_R_Shp, Class)
summary(tabla$Class)
Agua Artificial bosque cultivos Desnudo pastos represa
8 113 724 747 23 489 132
# Seperando los datos en Entrenamiento y prueba
Muestra <- sample(1:2235,1000) # sample separa aleatoriamente los datos
Pruebas <- tabla[Muestra,]
Entrenamiento <- tabla[-Muestra,]
# Aplicando el algoritmo DT
modelo <- rpart(Class ~.,data=Entrenamiento)
# Prediccion del modelo a obtener
prediccion <- predict(modelo, Pruebas[,-7],type="class")
# Definimos la paleta de colores segun las 7 clases
mycolor <- c("#00FF00","#0000FF","#228B22","#FF1493", "#1C1CE2", "#1C1CE2" ,"#5454AA")
## Visualizar prediccion a obtener
## Cantidad de punto que se encuentra en la clase
plot(prediccion, main = "Numero de pixeles identificados",
axes = TRUE, xlab = "Clases", ylab = "Numero de pixel", col= mycolor)
# Determinacion del indice de kappa de la clasificacion
MC <- table(Pruebas[,7],prediccion)
kappa<-sum(diag(MC))/sum(MC)
print(kappa)
[1] 0.952
# Graficando el arbol de decision
prp(modelo, extra=104, branch.type =1, box.col = c("pink","palegreen3")[modelo$frame$yval])
# Ejecutamos el clasificador DT
beginCluster()
8 cores detected, using 7
## 8 cores detected, using 7
tc_class<- clusterR(L8_2015, raster::predict,args = list(model =modelo,type="class"))
endCluster()
# Ploteamos la clasificación con Support vector machine
plot(tc_class,main=paste("kappa = ",kappa,sep =""),col =mycolor,cex.lab=0.8,
cex.axis=0.8,cex.main=0.9)
# Exportacion de bandas
writeRaster(tc_class,"Class_L8_Decision_Tree.tif", drivername="Gtiff",overwrite=TRUE)