Stat
LIBRERIA
library(terra)## terra version 1.0.10library(raster)## Loading required package: splibrary(sp)
library(psych)##
## Attaching package: 'psych'## The following objects are masked from 'package:terra':
##
## describe, rescalelibrary(corrplot)## corrplot 0.84 loadedlibrary(PerformanceAnalytics)## Loading required package: xts## Loading required package: zoo##
## Attaching package: 'zoo'## The following object is masked from 'package:terra':
##
## time<-## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric##
## Attaching package: 'PerformanceAnalytics'## The following object is masked from 'package:graphics':
##
## legendENTORNO Y DATOS
ruta <- getwd()
carpeta <- 'Datos'
file <- paste0(ruta,'/',carpeta,'/','LC08_L1TP_009057_20140820_20200911_02_T1_B', 1:7, ".tif")
landsat <-rast(file)
names(landsat) <- c('ultra-blue','blue','green','red','NIR','SWIR1','SWIR2')PROCESAMIENTO DE DATOS
xmx <- 374838
xmn <- 358682
ymx <- 467578
ymn <- 454847
ext(landsat)## SpatExtent : 273585, 501315, 362985, 595515 (xmin, xmax, ymin, ymax)e <- ext(xmn, xmx, ymn, ymx)
# crop landsat by the extent
landsatcrop <- crop(landsat, e)
Clandsat<- -0.1 + landsatcrop*2.0E-5Estadistica basica univariada
NIR <- as.list(Clandsat$NIR)
#promedio
mean(NIR)## [1] 0.2578855#varianza
var(NIR)## [1] 0.002819754#Desviacion estandar
sd(NIR)## [1] 0.05310136#oblicuidad, coeficiente de simetria
skew(NIR)## [1] -0.3566767#curtosis
kurtosi(NIR)## [1] 1.23192#histograma
hist(NIR,col = 'yellow',main = "NIR")
abline(v=mean(NIR),col='Blue')
abline(v=(mean(NIR)+sd(NIR)),col='red')
abline(v=(mean(NIR)-sd(NIR)),col='red')
legend(x ='topright',legend=c('mean = 0.2579','sd = 0.0531','var = 0.0028','skew = -0.3566767','Kurt = 1.23192'), fill=c('blue','red','white','white','white'))
#Estadistica datos espaciales
#######metodo numerico#######
l_max <- c()
l_min <- c()
l_mean <- c()
l_sd <- c()
l_skew <- c()
l_kurt <- c()
for (a in 1:7) {
b <- as.list(Clandsat[[a]])
x <- round(c(max(b)),4)
n <- round(c(min(b)),4)
m <- round(c(mean(b)),4)
s <- round(c(sd(b)),4)
w <- round(c(skew(b)),4)
k <- round(c(kurtosi(b)),4)
l_max <- c(l_max,x)
l_min <- c(l_min,n)
l_mean <- c(l_mean,m)
l_sd <- c(l_sd,s)
l_skew <- c(l_skew,w)
l_kurt <- c(l_kurt,k)
}
statL <- rbind(l_max[1:7],l_min[1:7],l_mean[1:7],l_sd[1:7],l_skew[1:7],l_kurt[1:7])
colnames(statL) <- c('ultrablue','blue','green','red','NIR','SWIR1','SWIR2')
rownames(statL) <- c('max','min','mean','sd','skew','kurt')
statL## ultrablue blue green red NIR SWIR1 SWIR2
## max 0.4690 0.5255 0.5799 0.6320 0.8783 1.2107 1.2107
## min 0.0899 0.0706 0.0511 0.0316 0.0115 -0.0064 -0.0031
## mean 0.1080 0.0926 0.0861 0.0750 0.2579 0.1820 0.0998
## sd 0.0080 0.0102 0.0121 0.0220 0.0531 0.0585 0.0453
## skew 4.8350 3.7420 3.2389 1.5064 -0.3567 0.8474 1.2831
## kurt 111.1462 77.1517 66.3253 10.1322 1.2319 3.2630 12.3001#######metodo cellstat########
writeRaster(landsatcrop, filename="brick.tif", overwrite=TRUE)
LSbrick <- brick(paste0(ruta,'/','brick.tif'))
LSbrick <- -0.1 + LSbrick*2.0E-5
LSb_max <- round(cellStats(LSbrick, 'max'),4)
LSb_min <- round(cellStats(LSbrick, 'min'),4)
LSb_mean <- round(cellStats(LSbrick, 'mean'),4)
LSb_sd <- round(cellStats(LSbrick, 'sd'),4)
LSb_skew <- round(cellStats(LSbrick, 'skew'),4)
LSb_kurtosi <- round(cellStats(LSbrick, 'kurtosi'),4)
stat <- rbind(LSb_max,LSb_min,LSb_mean,LSb_sd,LSb_skew,LSb_kurtosi)
stat## ultra.blue blue green red NIR SWIR1 SWIR2
## LSb_max 0.4690 0.5255 0.5799 0.6320 0.8783 1.2107 1.2107
## LSb_min 0.0899 0.0706 0.0511 0.0316 0.0115 -0.0064 -0.0031
## LSb_mean 0.1080 0.0926 0.0861 0.0750 0.2579 0.1820 0.0998
## LSb_sd 0.0080 0.0102 0.0121 0.0220 0.0531 0.0585 0.0453
## LSb_skew 4.8350 3.7420 3.2389 1.5064 -0.3567 0.8474 1.2831
## LSb_kurtosi 111.1462 77.1517 66.3253 10.1322 1.2319 3.2630 12.3001BRINK HISTOGRAMA
par(mfrow = c(1,2))
hist(LSbrick, c(1), maxpixels=1000000, plot=TRUE, main="Histogramas of B1 and B5", xlab= "reflectance", ylab="Frequency", col= "darkblue")
hist(LSbrick, c(5), maxpixels=1000000, plot=TRUE, main="Histogramas of B1 and B5", xlab= "reflectance", ylab="Frequency", col= "darkred")
#ANALISIS MULTIBANDA brick
par(mfrow=c(1,2))
# Create an RGB image from the raster stack
plotRGB(LSbrick, r = 4, g = 3, b = 2,main='L8 2013 RGB432', stretch='lin')
plotRGB(LSbrick, r = 5, g = 6, b = 4,main='L8 2013 RGB564', stretch='lin')
#CORRELACIONES
pairs(LSbrick[[1:3]], hist=TRUE, cor=TRUE, use="pairwise.complete.obs", maxpixels=100000)
pairs(LSbrick[[4:6]], hist=TRUE, cor=TRUE, use="pairwise.complete.obs", maxpixels=100000)
#Matriz de covarianza
covmat <- cov(as.matrix(LSbrick))
round(covmat*1000,4)#por mil para valores mas grandes## ultra.blue blue green red NIR SWIR1 SWIR2
## ultra.blue 0.0636 0.0805 0.0889 0.1557 -0.1470 0.2740 0.2627
## blue 0.0805 0.1042 0.1165 0.2089 -0.1968 0.3754 0.3530
## green 0.0889 0.1165 0.1453 0.2397 -0.1212 0.4125 0.3793
## red 0.1557 0.2089 0.2397 0.4833 -0.5078 1.0156 0.8833
## NIR -0.1470 -0.1968 -0.1212 -0.5078 2.8198 -0.7412 -0.9874
## SWIR1 0.2740 0.3754 0.4125 1.0156 -0.7412 3.4281 2.5271
## SWIR2 0.2627 0.3530 0.3793 0.8833 -0.9874 2.5271 2.0495Matriz de Correlacion
cmat <- cor(as.matrix(LSbrick))
round(cmat, 4)## ultra.blue blue green red NIR SWIR1 SWIR2
## ultra.blue 1.0000 0.9890 0.9245 0.8882 -0.3472 0.5868 0.7277
## blue 0.9890 1.0000 0.9467 0.9305 -0.3629 0.6280 0.7637
## green 0.9245 0.9467 1.0000 0.9042 -0.1894 0.5843 0.6949
## red 0.8882 0.9305 0.9042 1.0000 -0.4349 0.7890 0.8875
## NIR -0.3472 -0.3629 -0.1894 -0.4349 1.0000 -0.2384 -0.4107
## SWIR1 0.5868 0.6280 0.5843 0.7890 -0.2384 1.0000 0.9534
## SWIR2 0.7277 0.7637 0.6949 0.8875 -0.4107 0.9534 1.0000LA NAPA
ncmat <- round(cmat, 2)
corrplot(ncmat, type = "upper", order = "hclust", tl.col = "black", tl.srt = 45)
funcion chart.Correlation()[del paquete PerformanceAnalytics] para la matriz de correlacion
tmp1 <- as.matrix(sample(LSbrick,50000))
chart.Correlation(tmp1, histogram=TRUE, pch=19)
#reproducibilidad
sessionInfo()## R version 4.0.4 (2021-02-15)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows 10 x64 (build 17763)
##
## Matrix products: default
##
## locale:
## [1] LC_COLLATE=Spanish_Colombia.1252 LC_CTYPE=Spanish_Colombia.1252
## [3] LC_MONETARY=Spanish_Colombia.1252 LC_NUMERIC=C
## [5] LC_TIME=Spanish_Colombia.1252
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] PerformanceAnalytics_2.0.4 xts_0.12.1
## [3] zoo_1.8-8 corrplot_0.84
## [5] psych_2.0.12 raster_3.3-7
## [7] sp_1.4-2 terra_1.0-10
##
## loaded via a namespace (and not attached):
## [1] Rcpp_1.0.5 knitr_1.31 magrittr_1.5 mnormt_2.0.2
## [5] lattice_0.20-41 quadprog_1.5-8 rlang_0.4.10 highr_0.8
## [9] stringr_1.4.0 tools_4.0.4 rgdal_1.5-23 parallel_4.0.4
## [13] grid_4.0.4 tmvnsim_1.0-2 nlme_3.1-152 xfun_0.21
## [17] htmltools_0.5.1.1 yaml_2.2.1 digest_0.6.25 codetools_0.2-18
## [21] evaluate_0.14 rmarkdown_2.7 stringi_1.4.6 compiler_4.0.4