Glaciares en Chile

Tarapacá: Automatización del proceso de extracción de data por glaciar (análisis espacial) y en el tiempo.

VE-CC-CF

DataIntelligence

28-05-2021

Abstract

Éste trabajo contiene elas líneas para extraer el área considerada glaciar en cada límite administrativo entregado en un shp. Se propone una lógica para calcular el porcentaje glaciar. Se añade el código para automatizar la extracción en el tiempo. También avanzamos en torno a la determinación de las matrices de confusión para cada una de las técnicas de ML utilizadas, para mejorar las imágenes: RF, CART y SVM.

library(rgee)
# ee_Initialize("tarredwall@gmail.com", drive = TRUE)
ee_Initialize()

## -- rgee 1.0.9 --------------------------------------- earthengine-api 0.1.263 -- 
##  v email: not_defined
##  v Initializing Google Earth Engine:
 v Initializing Google Earth Engine:  DONE!
## 
 v Earth Engine user: users/rgee2 
## --------------------------------------------------------------------------------

R

A partir del 15 de febrero de 2018, se debe utilizar la siguiente denominación:

I Región de Arica y Parinacota. Lunes 24 II Región de Tarapacá. Martes 25 III Región de Antofagasta. Miercoles 26 IV Región de Atacama. Jueves 27 V Región de Coquimbo. Viernes 28 VI Región de Valparaíso. Lunes 31

Victor

VII Región Metropolitana de Santiago. Lunes 24 VIII Región del Libertador General Bernardo O’Higgins. Martes 25 IX Región del Maule. Miercoles 26 X Región del Ñuble. Jueves 27 XI Región del Biobío. Viernes 28 XII Región de La Araucanía. Lunes 31

Abner

XIII Región de Los Ríos. Miercoles 26 XIV Región de Los Lagos. Jueves 27 XV Región de Aysén del General Carlos Ibáñez del Campo. Viernes 28 XVI Región de Magallanes y la Antártica Chilena. Lunes 31

1. Cargamos el shp de la primera region y lo filterboundiamos con la imagen satelital sentinel. Cargamos el shp de glaciares, y obtenemos las muestras interiores de los limites administrativos de color rojo que significaran glaciares que tendren numero 0. Obtenemos muestras para los colores verdes fuera de los limites administrativos. Unimos ambas muestras y categorizamos nuevamente la imagen satelital sentinel.

I Región de Tarapacá.

1 Lo primero es dibujar el shp de la Región de Tarapacá

region_tara <- st_read("rrss/Region_01.shp", 
        quiet = TRUE) %>% 
        sf_as_ee()

## Registered S3 method overwritten by 'geojsonsf':
##   method        from   
##   print.geojson geojson

tara <- region_tara$geometry()
Map$setCenter(lon = -69.7522, lat = -20.18, zoom = 7)
Map$addLayer(tara)

2 lo segundo es asociar una imagen satelital Sentinel filtrada por fecha, nubosidad y vinculada a la region de estudio. Y establecemos que el parametro sea la media.

start <- ee$Date("2019-01-01")
finish <- ee$Date("2019-04-01")
cc <- 20

sentinel1 = ee$ImageCollection('COPERNICUS/S2')$filterDate(start, finish)$filterBounds(tara)$
  filter(ee$Filter$lt("CLOUDY_PIXEL_PERCENTAGE", cc))

first <- sentinel1$median()

vizParams <- list(
  bands = c("B8","B5" , "B3"),
  #  bands = c("B2", "B3"),
  min = 100,
  max = 1000,
  gamma = 2
)
Map$setCenter(lon = -69.7522, lat = -20.18, zoom = 7)
Map$addLayer(first, vizParams, "Landsat 8 image")

3 Lo tercero es construir un raster con NDGI mayor a 0 siempre vinculada a la primera imagen satelital que hemos extraído.

Construimos una imagen raster con una columna llamada nd con dos categorías para cada pixel: 0 si es rojo, 1 si es verde

getNDGI <- function(image) 
{
     image$normalizedDifference(c("B4", "B3")) > 0
}

ndgi <- getNDGI(first)

## Warning: Ops.ee.image.Image will be deprecated in rgee v.1.1.0. Please install
## rgeeExtra (https://github.com/r-earthengine/rgeeExtra). Deeply sorry for the
## inconveniences.

ndgiParams <- list(palette = c(
  "#d73027", "#f46d43", "#fdae61",
  "#fee08b", "#d9ef8b", "#a6d96a",
  "#66bd63", "#1a9850"
))

Map$setCenter(lon = -69.7522, lat = -20.18, zoom = 7)
Map$addLayer(ndgi, ndgiParams, "NDGI")

4 obtenemos los shp que nos indican la presencia de glaciares en la Region de Tarapaca.

glaciar_tara <- st_read("GlaciaresxRegion/Simple_TARAPACA.shp",
        quiet = TRUE) %>%
        sf_as_ee()
glaciar_tara <- glaciar_tara$geometry()
Map$setCenter(lon = -69.7522, lat = -20.18, zoom = 7)
Map$addLayer(glaciar_tara, ndgiParams, "NDGI")

5 Intersectamos las regiones glaciares de la region de Tarapaca con nuestro ndvi, y de esta forma obtenemos un universo muestral del cual poder extraer nuestro primer set de muestras (pixeles rojos interiores).

ESTA ES LA IMAGEN SOBRE LA QUE EXTRAEREMOS EL PRIMER SET DE MUESTRAS

glaciares_adm_con_ngvi <- ndgi$clip(glaciar_tara)
Map$setCenter(lon = -69.7522, lat = -20.18, zoom = 7)
Map$addLayer(glaciares_adm_con_ngvi, ndgiParams, "NDGI")

Con esto concluimos la primera parte de nuestra generacion de raster para, como siguiente paso, extraer las muestras.

6 Hay que hacer una diferencia (se le resta al shp de la región de Tarapaca los shps de los glaciares) para intersectarla con nuestra imagen ngdi y obtener el segundo universo muestral.

6.1 Primero obtenemos nuestros shapes sin zonas administrativas glaciares:

tara_sin_glacia <- tara$difference(glaciar_tara, ee$ErrorMargin(1))
Map$setCenter(lon = -69.7522, lat = -20.18, zoom = 7)
Map$addLayer(tara_sin_glacia, ndgiParams, "NDGI")

6.2 Necesitamos intersectar nuestro ndvi con la zona de Tarapaca que excluye los glaciares.

ndvi_tarapaca_sin_glaciares <- ndgi$clip(tara_sin_glacia)
Map$setCenter(lon = -69.7522, lat = -20.18, zoom = 7)
Map$addLayer(ndvi_tarapaca_sin_glaciares, ndgiParams, "NDGI")

7 Ahora tenemos dos regiones preparadas sobre las cuales podemos extraer nuestras muestras. A efectos de orden nombrémoslas:

7.1. pixeles rojos interiores

7.2. pixeles verdes exteriores

7.1. pixeles rojos interiores

Necesitamos intersectar las regiones administrativas glaciares con nuestro raster ndvi y de ahi extraer las muestras de color rojo que simbolizarán los glaciares.

limites_glaciares_con_rasters <- ndgi$clip(glaciar_tara)
Map$setCenter(lon = -68.69644, lat = -19.7581, zoom = 15)
Map$addLayer(limites_glaciares_con_rasters , ndgiParams, "NDGI")

1 Muestreo

# region_arica <- st_read("Regiones_separadas/Region_01.shp", quiet = TRUE)
# region_arica

Los rojos interiores son glaciares, 0 = rojo y son llamados muestras_in_gee_rojas.

muestras_in <- ndgi$sampleRegions(
    collection = ee$Feature(glaciar_tara),
    scale = 100,
    tileScale = 16,
    geometries = TRUE
  )

muestras_in_gee_rojas = muestras_in$filter(ee$Filter$eq('nd', 0))

Map$setCenter(lon = -68.69644, lat = -19.7581, zoom = 15)
Map$addLayer(
  eeObject =  muestras_in_gee_rojas,
  visParams = {},
  name = "puntos glaciares"
)

7.2. Pixeles verdes exteriores

ndvi_tara_sin_glaciares <- ndgi$clip(tara_sin_glacia)
Map$setCenter(lon = -68.69644, lat = -19.7581, zoom = 8)
Map$addLayer(ndvi_tara_sin_glaciares, ndgiParams, "NDGI")

Los verdes exteriores no son glaciares 1 = verde

muestras_in <- ndgi$sampleRegions(
    collection = ee$Feature(ndvi_tara_sin_glaciares),
    scale = 10000,
    tileScale = 16,
    geometries = TRUE
  )

muestras_out_gee_verdes = muestras_in$filter(ee$Filter$eq('nd', 1))

Map$setCenter(lon = -69.80584, lat = -17.67837, zoom = 7)
Map$addLayer(
  eeObject =  muestras_out_gee_verdes,
  visParams = {},
  name = "puntos no glaciares"
)

III debemos unir muestras y generar el random forest

union_total <- muestras_out_gee_verdes$merge(muestras_in_gee_rojas)

3 Aplicamos Random Forest

(volver al índice)

region_arica <- st_read("rrss/Region_01.shp", 
        quiet = TRUE) %>% 
        sf_as_ee()
arica <- region_arica$geometry()

start <- ee$Date("2019-01-01")
finish <- ee$Date("2019-04-01")
cc <- 20

sentinel1 = ee$ImageCollection('COPERNICUS/S2')$filterDate(start, finish)$filterBounds(arica)$
  filter(ee$Filter$lt("CLOUDY_PIXEL_PERCENTAGE", cc))

first <- sentinel1$median()

# Use these bands for prediction.
bands <- c("B2", "B3", "B4", "B8")

# This property of the table stores the land cover labels.
label <- "nd"

# Overlay the points on the imagery to get training.
training <- first$select(bands)$sampleRegions(
  collection = union_total,
  properties = list(label),
  scale = 10
)

# Train a CART classifier with default parameters.
trained <- ee$Classifier$smileRandomForest(10)$train(training, label, bands)

# Classify the image with the same bands used for training.
classified_rf <- first$select(bands)$classify(trained)

# Viz parameters.
# viz_img <- list(bands = c("B4", "B3", "B2"), max = 0.4)
vizParams <- list(
  bands = c("B8","B5" , "B3"),
  #  bands = c("B2", "B3"),
  min = 100,
  max = 1000,
  gamma = 2
)

viz_class <- list(palette = c("red", "green", "blue"), min = 0, max = 2)

Map$addLayer(first, vizParams, name = "image") +
Map$addLayer(classified_rf, viz_class, name = "classification")

18:29 18:30

Le aplicamos la matriz de confusion:

trainAccuracy = trained$confusionMatrix()
trainAccuracy$getInfo()

## [[1]]
## [1] 12  0
## 
## [[2]]
## [1]   0 451

18:40 18:43

4 Aplicamos CART

(volver al índice)

region_tara <- st_read("rrss/Region_01.shp", 
        quiet = TRUE) %>% 
        sf_as_ee()
tara <- region_tara$geometry()

start <- ee$Date("2019-01-01")
finish <- ee$Date("2019-04-01")
cc <- 20

sentinel1 = ee$ImageCollection('COPERNICUS/S2')$filterDate(start, finish)$filterBounds(tara)$
  filter(ee$Filter$lt("CLOUDY_PIXEL_PERCENTAGE", cc))

first <- sentinel1$median()

# Use these bands for prediction.
bands <- c("B2", "B3", "B4", "B8")

# This property of the table stores the land cover labels.
label <- "nd"

# Overlay the points on the imagery to get training.
training <- first$select(bands)$sampleRegions(
  collection = union_total,
  properties = list(label),
  scale = 10
)

# Train a CART classifier with default parameters.
trained <- ee$Classifier$smileCart(10)$train(training, label, bands)

# Classify the image with the same bands used for training.
classified <- first$select(bands)$classify(trained)

# Viz parameters.
# viz_img <- list(bands = c("B4", "B3", "B2"), max = 0.4)
vizParams <- list(
  bands = c("B8","B5" , "B3"),
  #  bands = c("B2", "B3"),
  min = 100,
  max = 1000,
  gamma = 2
)

viz_class <- list(palette = c("red", "green", "blue"), min = 0, max = 2)

Map$addLayer(first, vizParams, name = "image") +
Map$addLayer(classified, viz_class, name = "classification")

18:47 18:50

trainAccuracy = trained$confusionMatrix()
trainAccuracy$getInfo()

## [[1]]
## [1] 12  0
## 
## [[2]]
## [1]   0 451

18:51 18 53

4 Aplicamos SVM

(volver al índice)

region_arica <- st_read("rrss/Region_01.shp", 
        quiet = TRUE) %>% 
        sf_as_ee()
arica <- region_arica$geometry()

start <- ee$Date("2019-01-01")
finish <- ee$Date("2019-04-01")
cc <- 20

sentinel1 = ee$ImageCollection('COPERNICUS/S2')$filterDate(start, finish)$filterBounds(arica)$
  filter(ee$Filter$lt("CLOUDY_PIXEL_PERCENTAGE", cc))

first <- sentinel1$median()

# Use these bands for prediction.
bands <- c("B2", "B3", "B4", "B8")

# This property of the table stores the land cover labels.
label <- "nd"

# Overlay the points on the imagery to get training.
training <- first$select(bands)$sampleRegions(
  collection = union_total,
  properties = list(label),
  scale = 10
)
classifier = ee$Classifier$libsvm(
  kernelType = "RBF",
  gamma = 0.5,
  cost = 10
)
# Train a CART classifier with default parameters.
trained = classifier$train(training, label, bands)

# Classify the image with the same bands used for training.
classified <- first$select(bands)$classify(trained)

# Viz parameters.
# viz_img <- list(bands = c("B4", "B3", "B2"), max = 0.4)
vizParams <- list(
  bands = c("B8","B5" , "B3"),
  #  bands = c("B2", "B3"),
  min = 100,
  max = 1000,
  gamma = 2
)

viz_class <- list(palette = c("red", "green", "blue"), min = 0, max = 2)

Map$addLayer(first, vizParams, name = "image") +
Map$addLayer(classified, viz_class, name = "classification")

No se puede obtener la matriz de confusion para libsvm pues excede limites de memoria

trainAccuracy = trained$confusionMatrix()
trainAccuracy$getInfo()

## [[1]]
## [1] 12  0
## 
## [[2]]
## [1]   0 451

Analisis espacial

2. Una vez obtenido nuestro raster random forest interceptamos por limites glaciares, Y EN CADA UNO DE ELLOS determinamos la cantidad de pixeles rojos

mask_0 <- read_sf("GlaciaresxRegion/Simple_TARAPACA.shp")
#mask_0

reg_15_2 <- mask_0 %>% dplyr::select(COD_GLA,AREA_Km2)
varia <- reg_15_2 %>%
sf_as_ee()

region_0 <- varia$geometry()
sale_4 <- classified_rf$clip(region_0)

puntos_r <- sale_4$select('classification')$eq(0)


# Reducer$sum() returns a Reducer that computes the (weighted) sum of its inputs.

areas_rojas = puntos_r$reduceRegions(
reducer= ee$Reducer$sum(),
collection= varia,
scale= 30
)

rojas <- areas_rojas$getInfo()

rojas1 <- as.data.frame(rojas)
# rojas1


puntos_v <- sale_4$select('classification')$eq(1)

areas_verdes = puntos_v$reduceRegions(
reducer= ee$Reducer$sum(),
collection= varia,
scale= 30
)

verdes <- areas_verdes$getInfo()

verdes1 <- as.data.frame(verdes)
# verdes1



r_v <- rbind(colnames(rojas1),rojas1,verdes1)
r_v <- r_v[1:4,]
r_v <- r_v[-3,]
r_v <- t(r_v)
r_v <- as.data.frame(r_v)
names(r_v)[1] <- "col1"
names(r_v)[2] <- "col2"
names(r_v)[3] <- "col3"
r_v$col4 <- r_v$col3
r_v$col5 <- r_v$col3
# r_v


data1 <- filter(r_v, grepl("features.properties.COD_GLA",col1))
data1 <- as.data.frame(data1[,2])
# data1

data2 <- filter(r_v, grepl("features.properties.sum",col1))
data2 <- as.data.frame(data2[,c(2,3)])
# data2

data3 <- filter(r_v, grepl("features.properties.AREA_Km2",col1))
data3 <- as.data.frame(data3[,c(2,3)])
# data3

data4 <- cbind(data1,data2,data3)
colnames(data4) <- c("codigo_gla", "rojo", "verde", "area_Km2(shp)")
data4

##                             codigo_gla              rojo              verde
## features.properties.sum    CL101044025  49.2274509803922                  0
## features.properties.sum.1  CL101044018  14.4078431372549   98.0901960784314
## features.properties.sum.2  CL101044017  4.89411764705882    121.96862745098
## features.properties.sum.3  CL101044019                 0   97.7882352941177
## features.properties.sum.4  CL101044016  1.05098039215686   74.0392156862745
## features.properties.sum.5  CL101044015                 0   34.8666666666667
## features.properties.sum.6  CL101044014  8.13333333333333   21.7921568627451
## features.properties.sum.7  CL101044013                 0   40.1607843137255
## features.properties.sum.8  CL101044012                 0   13.8980392156863
## features.properties.sum.9  CL101044011  1.25882352941176   264.654901960784
## features.properties.sum.10 CL101044010                 0   5.82352941176471
## features.properties.sum.11 CL101044007                 0   57.3803921568627
## features.properties.sum.12 CL101044008 0.588235294117647   201.054901960784
## features.properties.sum.13 CL101044006                 0   221.235294117647
## features.properties.sum.14 CL101042002                 0   81.0078431372549
## features.properties.sum.15 CL101042001  2.75294117647059   77.1843137254902
## features.properties.sum.16 CL101042003                 0   63.2117647058824
## features.properties.sum.17 CL101042005  3.89411764705882   30.9882352941176
## features.properties.sum.18 CL101042004                 0    20.121568627451
## features.properties.sum.19 CL101610012                 0   154.098039215686
## features.properties.sum.20 CL101610013                 0   110.176470588235
## features.properties.sum.21 CL101610011                 0   6.61960784313725
## features.properties.sum.22 CL101610010                 0   94.2666666666667
## features.properties.sum.23 CL101610009                 0   72.7764705882353
## features.properties.sum.24 CL101610008                 0   41.6274509803922
## features.properties.sum.25 CL101610005                 0   78.9098039215686
## features.properties.sum.26 CL101610004                 0   63.0549019607843
## features.properties.sum.27 CL101610002                 0   113.752941176471
## features.properties.sum.28 CL101610003                 0   97.6431372549019
## features.properties.sum.29 CL101610001                 0   115.494117647059
## features.properties.sum.30 CL101611002                 0   159.976470588235
## features.properties.sum.31 CL101720004                 0   30.7843137254902
## features.properties.sum.32 CL101720003                 0                129
## features.properties.sum.33 CL101720001                 0   9.56470588235294
## features.properties.sum.34 CL101720002                 0   73.6862745098039
## features.properties.sum.35 CL101720008                 0   188.541176470588
## features.properties.sum.36 CL101720005                 0   22.3333333333333
## features.properties.sum.37 CL101720006                 0   15.2666666666667
## features.properties.sum.38 CL101720009                 0   13.2156862745098
## features.properties.sum.39 CL101720007                 0   9.65490196078431
## features.properties.sum.40 CL101730001                 0   129.050980392157
## features.properties.sum.41 CL101730002                 0   35.5254901960784
## features.properties.sum.42 CL101730009                 0    21.278431372549
## features.properties.sum.43 CL101730008                 0   19.0705882352941
## features.properties.sum.44 CL101730011                 0   53.4156862745098
## features.properties.sum.45 CL101730012  0.83921568627451   72.2509803921568
## features.properties.sum.46 CL101730013                 0   63.9176470588235
## features.properties.sum.47 CL101730010                 0   7.74117647058823
## features.properties.sum.48 CL101730014                 0   144.137254901961
## features.properties.sum.49 CL101730016                 0   25.1607843137255
## features.properties.sum.50 CL101730017                 0   9.94901960784314
## features.properties.sum.51 CL101730015                 0   9.06666666666666
## features.properties.sum.52 CL101730018                 0   11.5607843137255
## features.properties.sum.53 CL101730004                 0   33.2745098039216
## features.properties.sum.54 CL101730003                 0   79.1529411764706
## features.properties.sum.55 CL101730006                 0   27.6901960784314
## features.properties.sum.56 CL101730007                 0   24.5843137254902
## features.properties.sum.57 CL101730005                 0   13.2509803921569
## features.properties.sum.58 CL101044001  3.98823529411765   45.4196078431373
## features.properties.sum.59 CL101044002  1.10588235294118   25.9176470588235
## features.properties.sum.60 CL101044003 0.831372549019608   96.3333333333333
## features.properties.sum.61 CL101044004                 0   31.4117647058824
## features.properties.sum.62 CL101044009                 0   13.6901960784314
## features.properties.sum.63 CL101044005  11.5725490196078   10.5450980392157
## features.properties.sum.64 CL101040008                 0   51.1137254901961
## features.properties.sum.65 CL101611001                 0   128.905882352941
## features.properties.sum.66 CL101610024                 0   53.6862745098039
## features.properties.sum.67 CL101610023                 0    256.78431372549
## features.properties.sum.68 CL101610022                 0   229.074509803922
## features.properties.sum.69 CL101610021                 0   180.870588235294
## features.properties.sum.70 CL101610020                 0   106.007843137255
## features.properties.sum.71 CL101610018                 0    145.83137254902
## features.properties.sum.72 CL101610017                 0   119.317647058824
## features.properties.sum.73 CL101610016                 0   63.3176470588235
## features.properties.sum.74 CL101040006                 0   53.2901960784314
## features.properties.sum.75 CL101040005                 0   55.3372549019608
## features.properties.sum.76 CL101040004  4.65490196078431   8.56862745098039
## features.properties.sum.77 CL101040003  1.01960784313725    6.4156862745098
## features.properties.sum.78 CL101610015                 0   9.92156862745098
## features.properties.sum.79 CL101610014                 0   59.9921568627451
## features.properties.sum.80 CL101610019                 0   253.152941176471
## features.properties.sum.81 CL101040002                 0   50.8196078431372
## features.properties.sum.82 CL101040007  7.82352941176471    41.756862745098
## features.properties.sum.83 CL101040001  34.6627450980392   14.8117647058824
## features.properties.sum.84 CL101040009  55.3882352941176   46.2588235294118
## features.properties.sum.85 CL101044024  410.019607843137 0.0549019607843137
## features.properties.sum.86 CL101070001  22.3176470588235   7.57647058823529
## features.properties.sum.87 CL101044020  30.8039215686274                  0
## features.properties.sum.88 CL101044022  24.9058823529412                  0
## features.properties.sum.89 CL101044021  15.3411764705882                  0
## features.properties.sum.90 CL101044023   11.956862745098                  0
##                              area_Km2(shp)              NA
## features.properties.sum    0.0397270010711 0.0397270010711
## features.properties.sum.1   0.100876482862  0.100876482862
## features.properties.sum.2   0.109852809052  0.109852809052
## features.properties.sum.3  0.0873709651722 0.0873709651722
## features.properties.sum.4  0.0639482581884 0.0639482581884
## features.properties.sum.5  0.0390825631323 0.0390825631323
## features.properties.sum.6  0.0327605622869 0.0327605622869
## features.properties.sum.7  0.0382199791074 0.0382199791074
## features.properties.sum.8  0.0167628398049 0.0167628398049
## features.properties.sum.9   0.234736140147  0.234736140147
## features.properties.sum.10 0.0123072515574 0.0123072515574
## features.properties.sum.11 0.0490761861993 0.0490761861993
## features.properties.sum.12  0.181770894446  0.181770894446
## features.properties.sum.13  0.187241197168  0.187241197168
## features.properties.sum.14 0.0822324117868 0.0822324117868
## features.properties.sum.15 0.0811335051897 0.0811335051897
## features.properties.sum.16 0.0517317014976 0.0517317014976
## features.properties.sum.17 0.0342342838388 0.0342342838388
## features.properties.sum.18  0.025537434093  0.025537434093
## features.properties.sum.19   0.14556645363   0.14556645363
## features.properties.sum.20 0.0969989514597 0.0969989514597
## features.properties.sum.21 0.0133418051918 0.0133418051918
## features.properties.sum.22 0.0876500187874 0.0876500187874
## features.properties.sum.23 0.0722716222512 0.0722716222512
## features.properties.sum.24 0.0460606415757 0.0460606415757
## features.properties.sum.25 0.0863192881262 0.0863192881262
## features.properties.sum.26 0.0786716671606 0.0786716671606
## features.properties.sum.27 0.0996942476728 0.0996942476728
## features.properties.sum.28 0.0926249593749 0.0926249593749
## features.properties.sum.29  0.107310516435  0.107310516435
## features.properties.sum.30  0.125452215699  0.125452215699
## features.properties.sum.31 0.0370788979097 0.0370788979097
## features.properties.sum.32  0.121432954093  0.121432954093
## features.properties.sum.33 0.0156484103605 0.0156484103605
## features.properties.sum.34 0.0637414139919 0.0637414139919
## features.properties.sum.35  0.163304691787  0.163304691787
## features.properties.sum.36 0.0228456773695 0.0228456773695
## features.properties.sum.37 0.0161851633664 0.0161851633664
## features.properties.sum.38 0.0165129804368 0.0165129804368
## features.properties.sum.39 0.0198858464782 0.0198858464782
## features.properties.sum.40  0.121994691845  0.121994691845
## features.properties.sum.41 0.0391150993622 0.0391150993622
## features.properties.sum.42 0.0305984387632 0.0305984387632
## features.properties.sum.43 0.0225271588162 0.0225271588162
## features.properties.sum.44 0.0661863412709 0.0661863412709
## features.properties.sum.45 0.0641183674004 0.0641183674004
## features.properties.sum.46 0.0711972768362 0.0711972768362
## features.properties.sum.47 0.0127425943588 0.0127425943588
## features.properties.sum.48  0.130168840973  0.130168840973
## features.properties.sum.49 0.0425520003997 0.0425520003997
## features.properties.sum.50 0.0157886226892 0.0157886226892
## features.properties.sum.51 0.0155196563794 0.0155196563794
## features.properties.sum.52 0.0181153020415 0.0181153020415
## features.properties.sum.53 0.0346999506794 0.0346999506794
## features.properties.sum.54  0.071272988267  0.071272988267
## features.properties.sum.55 0.0311922006569 0.0311922006569
## features.properties.sum.56 0.0269613307305 0.0269613307305
## features.properties.sum.57 0.0182961256272 0.0182961256272
## features.properties.sum.58 0.0474476041426 0.0474476041426
## features.properties.sum.59  0.030925134966  0.030925134966
## features.properties.sum.60 0.0990839862833 0.0990839862833
## features.properties.sum.61 0.0316396276002 0.0316396276002
## features.properties.sum.62 0.0151640809662 0.0151640809662
## features.properties.sum.63 0.0261707721535 0.0261707721535
## features.properties.sum.64 0.0531454414338 0.0531454414338
## features.properties.sum.65  0.108143702435  0.108143702435
## features.properties.sum.66 0.0585930577914 0.0585930577914
## features.properties.sum.67  0.232784007153  0.232784007153
## features.properties.sum.68  0.189817275455  0.189817275455
## features.properties.sum.69  0.157189625326  0.157189625326
## features.properties.sum.70 0.0937245066846 0.0937245066846
## features.properties.sum.71  0.128313247654  0.128313247654
## features.properties.sum.72   0.11302848218   0.11302848218
## features.properties.sum.73 0.0630994631392 0.0630994631392
## features.properties.sum.74 0.0559590444049 0.0559590444049
## features.properties.sum.75 0.0578071886107 0.0578071886107
## features.properties.sum.76 0.0229483888793 0.0229483888793
## features.properties.sum.77 0.0126627103916 0.0126627103916
## features.properties.sum.78 0.0139952832236 0.0139952832236
## features.properties.sum.79 0.0548194450787 0.0548194450787
## features.properties.sum.80  0.224397707218  0.224397707218
## features.properties.sum.81 0.0476141336392 0.0476141336392
## features.properties.sum.82 0.0443321397721 0.0443321397721
## features.properties.sum.83 0.0490302735314 0.0490302735314
## features.properties.sum.84 0.0920502698629 0.0920502698629
## features.properties.sum.85  0.343230877798  0.343230877798
## features.properties.sum.86 0.0303463589963 0.0303463589963
## features.properties.sum.87 0.0304940656141 0.0304940656141
## features.properties.sum.88 0.0301132769341 0.0301132769341
## features.properties.sum.89 0.0145008622673 0.0145008622673
## features.properties.sum.90 0.0173113872952 0.0173113872952

Podemos obtener todos los glaciares de Chile y podemos obtener sus areas. A travaés de nuestros proceso de clasificación hemos podido determinar cantidades que debiesen ser glaciares en color rojo y cantidades que no son glaciares con valor de 1. Nuestro proceso de analisis nos llevo a la obtencion de estas cantidades cuyas magnitudes desconocemos. Pero lo anterior no es problema. Debemos seguir el siguiente procedimiento:

Calcular el porcentaje de rojo que hay dentro de la relación rojo-verde
Llamemos el porcentaje del punto 1 β.
debemos calcular el porcentaje que β representa en el area_Km2(shp). de cada glaciar.

\[ {\%rojo \over x} = {(rojo +verde) \over 100 } \]

Lo que nos interesa es obtener la superficie glaciar dentro de cada glaciar.

\[ { \%superficie\_roja \over 100} = {\%rojo \over region\_glaciar }\]

data4$porcentaje_rojo <- (100*as.numeric(data4$rojo)) / (as.numeric(data4$rojo) + as.numeric(data4$verde))
data4$porcentaje_area_glaciar <- (data4$porcentaje_rojo * as.numeric(data4$`area_Km2(shp)`)) / 100
data4

##                             codigo_gla              rojo              verde
## features.properties.sum    CL101044025  49.2274509803922                  0
## features.properties.sum.1  CL101044018  14.4078431372549   98.0901960784314
## features.properties.sum.2  CL101044017  4.89411764705882    121.96862745098
## features.properties.sum.3  CL101044019                 0   97.7882352941177
## features.properties.sum.4  CL101044016  1.05098039215686   74.0392156862745
## features.properties.sum.5  CL101044015                 0   34.8666666666667
## features.properties.sum.6  CL101044014  8.13333333333333   21.7921568627451
## features.properties.sum.7  CL101044013                 0   40.1607843137255
## features.properties.sum.8  CL101044012                 0   13.8980392156863
## features.properties.sum.9  CL101044011  1.25882352941176   264.654901960784
## features.properties.sum.10 CL101044010                 0   5.82352941176471
## features.properties.sum.11 CL101044007                 0   57.3803921568627
## features.properties.sum.12 CL101044008 0.588235294117647   201.054901960784
## features.properties.sum.13 CL101044006                 0   221.235294117647
## features.properties.sum.14 CL101042002                 0   81.0078431372549
## features.properties.sum.15 CL101042001  2.75294117647059   77.1843137254902
## features.properties.sum.16 CL101042003                 0   63.2117647058824
## features.properties.sum.17 CL101042005  3.89411764705882   30.9882352941176
## features.properties.sum.18 CL101042004                 0    20.121568627451
## features.properties.sum.19 CL101610012                 0   154.098039215686
## features.properties.sum.20 CL101610013                 0   110.176470588235
## features.properties.sum.21 CL101610011                 0   6.61960784313725
## features.properties.sum.22 CL101610010                 0   94.2666666666667
## features.properties.sum.23 CL101610009                 0   72.7764705882353
## features.properties.sum.24 CL101610008                 0   41.6274509803922
## features.properties.sum.25 CL101610005                 0   78.9098039215686
## features.properties.sum.26 CL101610004                 0   63.0549019607843
## features.properties.sum.27 CL101610002                 0   113.752941176471
## features.properties.sum.28 CL101610003                 0   97.6431372549019
## features.properties.sum.29 CL101610001                 0   115.494117647059
## features.properties.sum.30 CL101611002                 0   159.976470588235
## features.properties.sum.31 CL101720004                 0   30.7843137254902
## features.properties.sum.32 CL101720003                 0                129
## features.properties.sum.33 CL101720001                 0   9.56470588235294
## features.properties.sum.34 CL101720002                 0   73.6862745098039
## features.properties.sum.35 CL101720008                 0   188.541176470588
## features.properties.sum.36 CL101720005                 0   22.3333333333333
## features.properties.sum.37 CL101720006                 0   15.2666666666667
## features.properties.sum.38 CL101720009                 0   13.2156862745098
## features.properties.sum.39 CL101720007                 0   9.65490196078431
## features.properties.sum.40 CL101730001                 0   129.050980392157
## features.properties.sum.41 CL101730002                 0   35.5254901960784
## features.properties.sum.42 CL101730009                 0    21.278431372549
## features.properties.sum.43 CL101730008                 0   19.0705882352941
## features.properties.sum.44 CL101730011                 0   53.4156862745098
## features.properties.sum.45 CL101730012  0.83921568627451   72.2509803921568
## features.properties.sum.46 CL101730013                 0   63.9176470588235
## features.properties.sum.47 CL101730010                 0   7.74117647058823
## features.properties.sum.48 CL101730014                 0   144.137254901961
## features.properties.sum.49 CL101730016                 0   25.1607843137255
## features.properties.sum.50 CL101730017                 0   9.94901960784314
## features.properties.sum.51 CL101730015                 0   9.06666666666666
## features.properties.sum.52 CL101730018                 0   11.5607843137255
## features.properties.sum.53 CL101730004                 0   33.2745098039216
## features.properties.sum.54 CL101730003                 0   79.1529411764706
## features.properties.sum.55 CL101730006                 0   27.6901960784314
## features.properties.sum.56 CL101730007                 0   24.5843137254902
## features.properties.sum.57 CL101730005                 0   13.2509803921569
## features.properties.sum.58 CL101044001  3.98823529411765   45.4196078431373
## features.properties.sum.59 CL101044002  1.10588235294118   25.9176470588235
## features.properties.sum.60 CL101044003 0.831372549019608   96.3333333333333
## features.properties.sum.61 CL101044004                 0   31.4117647058824
## features.properties.sum.62 CL101044009                 0   13.6901960784314
## features.properties.sum.63 CL101044005  11.5725490196078   10.5450980392157
## features.properties.sum.64 CL101040008                 0   51.1137254901961
## features.properties.sum.65 CL101611001                 0   128.905882352941
## features.properties.sum.66 CL101610024                 0   53.6862745098039
## features.properties.sum.67 CL101610023                 0    256.78431372549
## features.properties.sum.68 CL101610022                 0   229.074509803922
## features.properties.sum.69 CL101610021                 0   180.870588235294
## features.properties.sum.70 CL101610020                 0   106.007843137255
## features.properties.sum.71 CL101610018                 0    145.83137254902
## features.properties.sum.72 CL101610017                 0   119.317647058824
## features.properties.sum.73 CL101610016                 0   63.3176470588235
## features.properties.sum.74 CL101040006                 0   53.2901960784314
## features.properties.sum.75 CL101040005                 0   55.3372549019608
## features.properties.sum.76 CL101040004  4.65490196078431   8.56862745098039
## features.properties.sum.77 CL101040003  1.01960784313725    6.4156862745098
## features.properties.sum.78 CL101610015                 0   9.92156862745098
## features.properties.sum.79 CL101610014                 0   59.9921568627451
## features.properties.sum.80 CL101610019                 0   253.152941176471
## features.properties.sum.81 CL101040002                 0   50.8196078431372
## features.properties.sum.82 CL101040007  7.82352941176471    41.756862745098
## features.properties.sum.83 CL101040001  34.6627450980392   14.8117647058824
## features.properties.sum.84 CL101040009  55.3882352941176   46.2588235294118
## features.properties.sum.85 CL101044024  410.019607843137 0.0549019607843137
## features.properties.sum.86 CL101070001  22.3176470588235   7.57647058823529
## features.properties.sum.87 CL101044020  30.8039215686274                  0
## features.properties.sum.88 CL101044022  24.9058823529412                  0
## features.properties.sum.89 CL101044021  15.3411764705882                  0
## features.properties.sum.90 CL101044023   11.956862745098                  0
##                              area_Km2(shp)              NA porcentaje_rojo
## features.properties.sum    0.0397270010711 0.0397270010711     100.0000000
## features.properties.sum.1   0.100876482862  0.100876482862      12.8071949
## features.properties.sum.2   0.109852809052  0.109852809052       3.8578053
## features.properties.sum.3  0.0873709651722 0.0873709651722       0.0000000
## features.properties.sum.4  0.0639482581884 0.0639482581884       1.3996240
## features.properties.sum.5  0.0390825631323 0.0390825631323       0.0000000
## features.properties.sum.6  0.0327605622869 0.0327605622869      27.1786135
## features.properties.sum.7  0.0382199791074 0.0382199791074       0.0000000
## features.properties.sum.8  0.0167628398049 0.0167628398049       0.0000000
## features.properties.sum.9   0.234736140147  0.234736140147       0.4733955
## features.properties.sum.10 0.0123072515574 0.0123072515574       0.0000000
## features.properties.sum.11 0.0490761861993 0.0490761861993       0.0000000
## features.properties.sum.12  0.181770894446  0.181770894446       0.2917210
## features.properties.sum.13  0.187241197168  0.187241197168       0.0000000
## features.properties.sum.14 0.0822324117868 0.0822324117868       0.0000000
## features.properties.sum.15 0.0811335051897 0.0811335051897       3.4438776
## features.properties.sum.16 0.0517317014976 0.0517317014976       0.0000000
## features.properties.sum.17 0.0342342838388 0.0342342838388      11.1635750
## features.properties.sum.18  0.025537434093  0.025537434093       0.0000000
## features.properties.sum.19   0.14556645363   0.14556645363       0.0000000
## features.properties.sum.20 0.0969989514597 0.0969989514597       0.0000000
## features.properties.sum.21 0.0133418051918 0.0133418051918       0.0000000
## features.properties.sum.22 0.0876500187874 0.0876500187874       0.0000000
## features.properties.sum.23 0.0722716222512 0.0722716222512       0.0000000
## features.properties.sum.24 0.0460606415757 0.0460606415757       0.0000000
## features.properties.sum.25 0.0863192881262 0.0863192881262       0.0000000
## features.properties.sum.26 0.0786716671606 0.0786716671606       0.0000000
## features.properties.sum.27 0.0996942476728 0.0996942476728       0.0000000
## features.properties.sum.28 0.0926249593749 0.0926249593749       0.0000000
## features.properties.sum.29  0.107310516435  0.107310516435       0.0000000
## features.properties.sum.30  0.125452215699  0.125452215699       0.0000000
## features.properties.sum.31 0.0370788979097 0.0370788979097       0.0000000
## features.properties.sum.32  0.121432954093  0.121432954093       0.0000000
## features.properties.sum.33 0.0156484103605 0.0156484103605       0.0000000
## features.properties.sum.34 0.0637414139919 0.0637414139919       0.0000000
## features.properties.sum.35  0.163304691787  0.163304691787       0.0000000
## features.properties.sum.36 0.0228456773695 0.0228456773695       0.0000000
## features.properties.sum.37 0.0161851633664 0.0161851633664       0.0000000
## features.properties.sum.38 0.0165129804368 0.0165129804368       0.0000000
## features.properties.sum.39 0.0198858464782 0.0198858464782       0.0000000
## features.properties.sum.40  0.121994691845  0.121994691845       0.0000000
## features.properties.sum.41 0.0391150993622 0.0391150993622       0.0000000
## features.properties.sum.42 0.0305984387632 0.0305984387632       0.0000000
## features.properties.sum.43 0.0225271588162 0.0225271588162       0.0000000
## features.properties.sum.44 0.0661863412709 0.0661863412709       0.0000000
## features.properties.sum.45 0.0641183674004 0.0641183674004       1.1481919
## features.properties.sum.46 0.0711972768362 0.0711972768362       0.0000000
## features.properties.sum.47 0.0127425943588 0.0127425943588       0.0000000
## features.properties.sum.48  0.130168840973  0.130168840973       0.0000000
## features.properties.sum.49 0.0425520003997 0.0425520003997       0.0000000
## features.properties.sum.50 0.0157886226892 0.0157886226892       0.0000000
## features.properties.sum.51 0.0155196563794 0.0155196563794       0.0000000
## features.properties.sum.52 0.0181153020415 0.0181153020415       0.0000000
## features.properties.sum.53 0.0346999506794 0.0346999506794       0.0000000
## features.properties.sum.54  0.071272988267  0.071272988267       0.0000000
## features.properties.sum.55 0.0311922006569 0.0311922006569       0.0000000
## features.properties.sum.56 0.0269613307305 0.0269613307305       0.0000000
## features.properties.sum.57 0.0182961256272 0.0182961256272       0.0000000
## features.properties.sum.58 0.0474476041426 0.0474476041426       8.0720692
## features.properties.sum.59  0.030925134966  0.030925134966       4.0922943
## features.properties.sum.60 0.0990839862833 0.0990839862833       0.8556322
## features.properties.sum.61 0.0316396276002 0.0316396276002       0.0000000
## features.properties.sum.62 0.0151640809662 0.0151640809662       0.0000000
## features.properties.sum.63 0.0261707721535 0.0261707721535      52.3226950
## features.properties.sum.64 0.0531454414338 0.0531454414338       0.0000000
## features.properties.sum.65  0.108143702435  0.108143702435       0.0000000
## features.properties.sum.66 0.0585930577914 0.0585930577914       0.0000000
## features.properties.sum.67  0.232784007153  0.232784007153       0.0000000
## features.properties.sum.68  0.189817275455  0.189817275455       0.0000000
## features.properties.sum.69  0.157189625326  0.157189625326       0.0000000
## features.properties.sum.70 0.0937245066846 0.0937245066846       0.0000000
## features.properties.sum.71  0.128313247654  0.128313247654       0.0000000
## features.properties.sum.72   0.11302848218   0.11302848218       0.0000000
## features.properties.sum.73 0.0630994631392 0.0630994631392       0.0000000
## features.properties.sum.74 0.0559590444049 0.0559590444049       0.0000000
## features.properties.sum.75 0.0578071886107 0.0578071886107       0.0000000
## features.properties.sum.76 0.0229483888793 0.0229483888793      35.2016607
## features.properties.sum.77 0.0126627103916 0.0126627103916      13.7130802
## features.properties.sum.78 0.0139952832236 0.0139952832236       0.0000000
## features.properties.sum.79 0.0548194450787 0.0548194450787       0.0000000
## features.properties.sum.80  0.224397707218  0.224397707218       0.0000000
## features.properties.sum.81 0.0476141336392 0.0476141336392       0.0000000
## features.properties.sum.82 0.0443321397721 0.0443321397721      15.7794827
## features.properties.sum.83 0.0490302735314 0.0490302735314      70.0618263
## features.properties.sum.84 0.0920502698629 0.0920502698629      54.4907407
## features.properties.sum.85  0.343230877798  0.343230877798      99.9866117
## features.properties.sum.86 0.0303463589963 0.0303463589963      74.6556474
## features.properties.sum.87 0.0304940656141 0.0304940656141     100.0000000
## features.properties.sum.88 0.0301132769341 0.0301132769341     100.0000000
## features.properties.sum.89 0.0145008622673 0.0145008622673     100.0000000
## features.properties.sum.90 0.0173113872952 0.0173113872952     100.0000000
##                            porcentaje_area_glaciar
## features.properties.sum               0.0397270011
## features.properties.sum.1             0.0129194478
## features.properties.sum.2             0.0042379074
## features.properties.sum.3             0.0000000000
## features.properties.sum.4             0.0008950352
## features.properties.sum.5             0.0000000000
## features.properties.sum.6             0.0089038666
## features.properties.sum.7             0.0000000000
## features.properties.sum.8             0.0000000000
## features.properties.sum.9             0.0011112303
## features.properties.sum.10            0.0000000000
## features.properties.sum.11            0.0000000000
## features.properties.sum.12            0.0005302638
## features.properties.sum.13            0.0000000000
## features.properties.sum.14            0.0000000000
## features.properties.sum.15            0.0027941386
## features.properties.sum.16            0.0000000000
## features.properties.sum.17            0.0038217700
## features.properties.sum.18            0.0000000000
## features.properties.sum.19            0.0000000000
## features.properties.sum.20            0.0000000000
## features.properties.sum.21            0.0000000000
## features.properties.sum.22            0.0000000000
## features.properties.sum.23            0.0000000000
## features.properties.sum.24            0.0000000000
## features.properties.sum.25            0.0000000000
## features.properties.sum.26            0.0000000000
## features.properties.sum.27            0.0000000000
## features.properties.sum.28            0.0000000000
## features.properties.sum.29            0.0000000000
## features.properties.sum.30            0.0000000000
## features.properties.sum.31            0.0000000000
## features.properties.sum.32            0.0000000000
## features.properties.sum.33            0.0000000000
## features.properties.sum.34            0.0000000000
## features.properties.sum.35            0.0000000000
## features.properties.sum.36            0.0000000000
## features.properties.sum.37            0.0000000000
## features.properties.sum.38            0.0000000000
## features.properties.sum.39            0.0000000000
## features.properties.sum.40            0.0000000000
## features.properties.sum.41            0.0000000000
## features.properties.sum.42            0.0000000000
## features.properties.sum.43            0.0000000000
## features.properties.sum.44            0.0000000000
## features.properties.sum.45            0.0007362019
## features.properties.sum.46            0.0000000000
## features.properties.sum.47            0.0000000000
## features.properties.sum.48            0.0000000000
## features.properties.sum.49            0.0000000000
## features.properties.sum.50            0.0000000000
## features.properties.sum.51            0.0000000000
## features.properties.sum.52            0.0000000000
## features.properties.sum.53            0.0000000000
## features.properties.sum.54            0.0000000000
## features.properties.sum.55            0.0000000000
## features.properties.sum.56            0.0000000000
## features.properties.sum.57            0.0000000000
## features.properties.sum.58            0.0038300034
## features.properties.sum.59            0.0012655475
## features.properties.sum.60            0.0008477945
## features.properties.sum.61            0.0000000000
## features.properties.sum.62            0.0000000000
## features.properties.sum.63            0.0136932533
## features.properties.sum.64            0.0000000000
## features.properties.sum.65            0.0000000000
## features.properties.sum.66            0.0000000000
## features.properties.sum.67            0.0000000000
## features.properties.sum.68            0.0000000000
## features.properties.sum.69            0.0000000000
## features.properties.sum.70            0.0000000000
## features.properties.sum.71            0.0000000000
## features.properties.sum.72            0.0000000000
## features.properties.sum.73            0.0000000000
## features.properties.sum.74            0.0000000000
## features.properties.sum.75            0.0000000000
## features.properties.sum.76            0.0080782140
## features.properties.sum.77            0.0017364476
## features.properties.sum.78            0.0000000000
## features.properties.sum.79            0.0000000000
## features.properties.sum.80            0.0000000000
## features.properties.sum.81            0.0000000000
## features.properties.sum.82            0.0069953823
## features.properties.sum.83            0.0343515051
## features.properties.sum.84            0.0501588739
## features.properties.sum.85            0.3431849251
## features.properties.sum.86            0.0226552708
## features.properties.sum.87            0.0304940656
## features.properties.sum.88            0.0301132769
## features.properties.sum.89            0.0145008623
## features.properties.sum.90            0.0173113873

Logica de la informacion temporal:

region_tara <- st_read("rrss/Region_01.shp", 
        quiet = TRUE) %>% 
        sf_as_ee()

glaciar_tara <- st_read("GlaciaresxRegion/Simple_TARAPACA.shp",
        quiet = TRUE) %>%
        sf_as_ee()

fn_cart <- function(a,m){

i_cuatrimestre <- switch(m, "01","05","09")
f_cuatrimestre <- switch(m, "04","08","12")


fecha_1 <- paste("20",a,"-",i_cuatrimestre,"-01", sep = "")
fecha_2 <- paste("20",a,"-",f_cuatrimestre,"-01", sep = "")
tara <- region_tara$geometry()


start <- ee$Date(fecha_1)
finish <- ee$Date(fecha_2)
cc <- 20

sentinel1 = ee$ImageCollection('COPERNICUS/S2')$filterDate(start, finish)$filterBounds(tara)$
  filter(ee$Filter$lt("CLOUDY_PIXEL_PERCENTAGE", cc))

first <- sentinel1$median()

vizParams <- list(
  bands = c("B8","B5" , "B3"),
  #  bands = c("B2", "B3"),
  min = 100,
  max = 1000,
  gamma = 2
)

getNDGI <- function(image) 
{
     image$normalizedDifference(c("B4", "B3")) > 0
}

ndgi <- getNDGI(first)

ndgiParams <- list(palette = c(
  "#d73027", "#f46d43", "#fdae61",
  "#fee08b", "#d9ef8b", "#a6d96a",
  "#66bd63", "#1a9850"
))



glaciar_tara <- glaciar_tara$geometry()

glaciares_adm_con_ngvi <- ndgi$clip(glaciar_tara)

tara_sin_glacia <- tara$difference(glaciar_tara, ee$ErrorMargin(1))

ndvi_tarapaca_sin_glaciares <- ndgi$clip(tara_sin_glacia)

limites_glaciares_con_rasters <- ndgi$clip(glaciar_tara)

# region_arica <- st_read("Regiones_separadas/Region_01.shp", quiet = TRUE)
# region_arica

muestras_in <- ndgi$sampleRegions(
    collection = ee$Feature(glaciar_tara),
    scale = 100,
    tileScale = 16,
    geometries = TRUE
  )

  muestras_in_gee_rojas = muestras_in$filter(ee$Filter$eq('nd', 0))



ndvi_tara_sin_glaciares <- ndgi$clip(tara_sin_glacia)


muestras_in <- ndgi$sampleRegions(
    collection = ee$Feature(ndvi_tara_sin_glaciares),
    scale = 10000,
    tileScale = 16,
    geometries = TRUE
  )

muestras_out_gee_verdes = muestras_in$filter(ee$Filter$eq('nd', 1))

tryCatch(
        {

union_total <- muestras_out_gee_verdes$merge(muestras_in_gee_rojas)

region_arica <- st_read("rrss/Region_01.shp", 
        quiet = TRUE) %>% 
        sf_as_ee()
arica <- region_arica$geometry()

start <- ee$Date(fecha_1)
finish <- ee$Date(fecha_2)
cc <- 20

sentinel1 = ee$ImageCollection('COPERNICUS/S2')$filterDate(start, finish)$filterBounds(arica)$
  filter(ee$Filter$lt("CLOUDY_PIXEL_PERCENTAGE", cc))

first <- sentinel1$median()

# Use these bands for prediction.
bands <- c("B2", "B3", "B4", "B8")

# This property of the table stores the land cover labels.
label <- "nd"

# Overlay the points on the imagery to get training.
training <- first$select(bands)$sampleRegions(
  collection = union_total,
  properties = list(label),
  scale = 10
)

# Train a CART classifier with default parameters.
trained <- ee$Classifier$smileRandomForest(10)$train(training, label, bands)

# Classify the image with the same bands used for training.
classified1 <- first$select(bands)$classify(trained)


periodo <- paste(fecha_1," / ",fecha_2," ", sep = "")
print(periodo)

###################################################################
mask_0 <- read_sf("GlaciaresxRegion/Simple_TARAPACA.shp")


reg_15_2 <- mask_0 %>% dplyr::select(COD_GLA,AREA_Km2)
varia <- reg_15_2 %>%
sf_as_ee()

region_0 <- varia$geometry()
sale_4 <- classified1$clip(region_0)

  puntos_r <- sale_4$select('classification')$eq(0)
  
  areas_rojas = puntos_r$reduceRegions(
  reducer= ee$Reducer$sum(),
  collection= varia,
  scale= 30
  )
  
  rojas <- areas_rojas$getInfo()
  
  rojas1 <- as.data.frame(rojas)
  
  
  
  puntos_v <- sale_4$select('classification')$eq(1)
  
  areas_verdes = puntos_v$reduceRegions(
  reducer= ee$Reducer$sum(),
  collection= varia,
  scale= 30
  )
  
  verdes <- areas_verdes$getInfo()
  
  verdes1 <- as.data.frame(verdes)
  
  
  ###########
  
  r_v <- rbind(colnames(rojas1),rojas1,verdes1)
  r_v <- r_v[1:4,]
  r_v <- r_v[-3,]
  r_v <- t(r_v)
  r_v <- as.data.frame(r_v)
  names(r_v)[1] <- "col1"
  names(r_v)[2] <- "col2"
  names(r_v)[3] <- "col3"
  r_v$col4 <- r_v$col3
  r_v$col5 <- r_v$col3
  
  ###########
  
  data1 <- filter(r_v, grepl("features.properties.COD_GLA",col1))
  data1 <- as.data.frame(data1[,2])
  
  data2 <- filter(r_v, grepl("features.properties.sum",col1))
  data2 <- as.data.frame(data2[,c(2,3)])
  
  data3 <- filter(r_v, grepl("features.properties.AREA_Km2",col1))
  data3 <- as.data.frame(data3[,c(2,3)])
  
  data4 <- cbind(data1,data2,data3)
  data4 <- data4[,-c(4)]
  
  colnames(data4) <- c("codigo_gla", "rojo", "verde", "area_Km2(shp)")
  
  nombre <- paste("TARAPACA_",fecha_1,"_",fecha_2,".xlsx", sep = "")
  
  data4$porcentaje_rojo <- (100*as.numeric(data4$rojo)) / (as.numeric(data4$rojo) + as.numeric(data4$verde))
  data4$porcentaje_area_glaciar <- (data4$porcentaje_rojo * as.numeric(data4$`area_Km2(shp)`)) / 100
  
  data4$fecha <- periodo
  
   write_xlsx(data4,nombre)
  # print(data4)
   ########################################

# box <- ee$Geometry$Rectangle(
#   coords = c(-70.28665 , -21.63061, -68.40454 , -18.9369),
#   proj = "EPSG:4326",
#   geodesic = FALSE
# )
# 
# img <- classified1
# 
# dem_raster2 <- img %>% ee_as_raster(
#       region = box,
#             scale = 10,)
# 
# writeRaster(dem_raster2,"dem_raster2.tif")

   ########################################

}, error = function(msg){
            print("No hay muestras rojas")
        })

write_xlsx(data4,nombre)

}

for (a in 20:20) {
for (m in 1:1) {

fn_cart(a,m)

}
}

## Warning: Ops.ee.image.Image will be deprecated in rgee v.1.1.0. Please install
## rgeeExtra (https://github.com/r-earthengine/rgeeExtra). Deeply sorry for the
## inconveniences.

## [1] "2020-01-01 / 2020-04-01 "

Vamos a exponer la data de dos cuatrimetres, el segundo y el tercero ambos para el año 2018

# segundo_cuatri_18 <- fn_cart(18,2)
# tercer_cuatri_18 <- fn_cart(18,3)

Analisis temporal

library(readxl)

## Warning: package 'readxl' was built under R version 4.0.5

# a <- read_xlsx("arica_2017-05-01_2017-08-01.xlsx")
# b <- read_xlsx("arica_2019-05-01_2019-08-01.xlsx")
# c <- read_xlsx("arica_2020-05-01_2020-08-01.xlsx")
# 
# abc <- rbind(a$codigo_gla,a$porcentaje_area_glaciar,b$porcentaje_area_glaciar,c$porcentaje_area_glaciar)
# 
# write_xlsx(abc,"2do_periodo_17_19_20.xlsx")
# 
# kbl(abc) %>%
#  kable_styling(bootstrap_options = c("striped", "hover")) %>%
#  kable_paper() %>%
#  scroll_box(width = "100%", height = "300px")

# tercer_cuatri_18 <- fn_cart(18,3)
# 
# union <- rbind(segundo_cuatri_18,tercer_cuatri_18)
# 
# kbl(union) %>%
# kable_styling(bootstrap_options = c("striped", "hover")) %>%
# kable_paper() %>%
# scroll_box(width = "100%", height = "300px")
# 
# write_xlsx(union,"union.xlsx")