Poglavje 1: Abiotski in biotski parametri v pragozdu Krokar
Žan Kuralt
Temperatura tal
Najprej aktivirajmo vse potrebne pakete.
library(tidyverse)
library(ggpubr)
library(dabestr)
library(RColorBrewer)Analiza podatkov pridobljenih z avtomatskimi merilniki temperature.
# Uvoz podatkov
t_mlad <- read.csv(file = "~/duhturat/temperature/zan_mladovje_03_0.csv", stringsAsFactors = FALSE)[, -1]
t_mlad$DateTime <- as.POSIXct(t_mlad$DateTime, format = "%m/%d/%y %H:%M:%S %p")
t_mlad$stage <- "mladovje"
t_mlad$date <- as.Date(t_mlad$DateTime, format = "%m/%d/%y %H:%M:%S %p")
t_drog <- read.csv(file = "~/duhturat/temperature/zan_drogovnjak_04_0.csv", stringsAsFactors = FALSE)[, -1]
t_drog$DateTime <- as.POSIXct(t_drog$DateTime, format = "%m/%d/%y %H:%M:%S %p")
t_drog$stage <- "drogovnjak"
t_drog$date <- as.Date(t_drog$DateTime, format = "%m/%d/%y %H:%M:%S %p")
t_deb <- read.csv(file = "~/duhturat/temperature/zan_debeljak_05_0.csv", stringsAsFactors = FALSE)[, -1]
t_deb$DateTime <- as.POSIXct(t_deb$DateTime, format = "%m/%d/%y %H:%M:%S %p")
t_deb$stage <- "debeljak"
t_deb$date <- as.Date(t_deb$DateTime, format = "%m/%d/%y %H:%M:%S %p")
t_gap <- read.csv(file = "~/duhturat/temperature/zan_presvetlitev_01_0.csv", stringsAsFactors = FALSE)[, -1]
t_gap$DateTime <- as.POSIXct(t_gap$DateTime, format = "%m/%d/%y %H:%M:%S %p")
t_gap$stage <- "presvetlitev"
t_gap$date <- as.Date(t_gap$DateTime, format = "%m/%d/%y %H:%M:%S %p")
data <- rbind(t_mlad, t_drog, t_deb, t_gap)
data <- subset(data, data$date > as.Date("2018-10-18") & data$date < as.Date("2019-10-30"))
data$month <- months(data$date)
data$month <- factor(data$month, ordered = TRUE, levels = c("January", "February", "March", "April", "May", "June", "July", "August", "September", "October", "November", "December"))
# Povprečna letna temperatura
data %>%
summarise(mean = mean(Temp))## mean
## 1 8.027683#Povprečne temperature v različnih razvojnih fazah
data %>%
group_by(stage) %>%
summarise(mean = mean(Temp))## # A tibble: 4 x 2
## stage mean
## <chr> <dbl>
## 1 debeljak 8.04
## 2 drogovnjak 7.62
## 3 mladovje 7.90
## 4 presvetlitev 8.55temps <- ggplot(data) +
geom_smooth(aes(x = date, y = Temp, color = stage)) +
# geom_line(aes(x = date, y = Temp, color = stage), alpha = 0.2) +
theme_bw() +
theme(panel.grid = element_blank(), legend.position = "left", axis.title.x = element_blank()) +
labs(color = "Razvojna faza") +
scale_color_brewer(palette = "Set1") +
ylab("Temperatura [°C]")
dates <- unique(as.character(data$date))
res <- data.frame(stg = NA,
date = NA,
low = NA,
high = NA,
var = NA)
n <- 0
for(stage in unique(data$stage)){
stg <- stage
df <- data[data$stage == stg, ]
for(i in dates){
n <- n + 1
dat <- df[df$date == i, ]
low <- min(dat$Temp)
high <- max(dat$Temp)
range <- high - low
out <- c(stg, i, low, high, range)
res <- rbind(res, out)
}
}
res <- na.omit(res)
res$date <- as.Date(res$date)
res$var <- as.numeric(res$var)
res$stg <- as.factor(res$stg)
anomaly <- ggplot() +
geom_smooth(data = res, aes(x = date, y = var, color = stg), span = 0.4, se = TRUE) +
# geom_line(data = res, aes(x = date, y = var), alpha = 0.05) +
# ggtitle("Daily temperature anomaly") +
ylab("Temperaturna razlika [°C]") +
theme_bw() +
theme(panel.grid = element_blank(), legend.position = "none", axis.title.x = element_blank()) +
scale_color_brewer(palette = "Set1")
theme_set(theme_pubr())
ggarrange(temps, anomaly, ncol = 2, nrow = 1)
res$month <- months(res$date)
res$month <- factor(res$month, ordered = TRUE, levels = c("January", "February", "March", "April", "May", "June", "July", "August", "September", "October",
"November", "December"))
estimation.plot <- res %>%
dabest(x = stg, y = var,
idx = c("mladovje", "drogovnjak", "debeljak", "presvetlitev"),
paired = FALSE) %>%
mean_diff()
plot(estimation.plot, color.column = month, palette = "Set3")
LiDAR in Sentinel-2
Aktivirajmo pakete, ki jih bomo potrebovali za delo z podatki LiDAR.
library(lidR)
library(ForestTools)
library(rLiDAR)
library(rgl)
library(ggbiplot)
library(ggpubr)
library(foreach)
library(doParallel)
library(sf)
source("~/duhturat/lidar/functions.R")Abiotski parametri
Naklon, izpostavljenost, povprečno nadmorsko višino, razpon nadmorskih višin, TRI (angl. Terrain Ruggedness Index), TPI (angl. Topographic Position Index), grobost (angl. roughness) terena smo izračunali s pomočjo javno dostopnih LiDAR posnetkov območja raziskave. S funkcijo lidR::grid_terrain() smo izdelali digitalni model terena ločljivosti 1x1 m ter nato s funkcijo raster::terrain() izračunali prej omenjene parametre za vseh 77 vzorčnih kvadrantov.
# Uvoz oblaka točk, mreže vzorčnih kvadrantov ter območja pragozda Krokar.
krokar <- lidR::readLAScatalog("~/duhturat/lidar/krokar.las")
grid <- st_read("~/duhturat/lidar/sampling/grid_50_plus.shp")## Reading layer `grid_50_plus' from data source
## `C:\Users\zanku\OneDrive - Univerza v Ljubljani\Documents\duhturat\lidar\sampling\grid_50_plus.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 77 features and 5 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: 481880.4 ymin: 43547.11 xmax: 482880.5 ymax: 45247.13
## Projected CRS: MGI 1901 / Slovene National Gridkrok <- st_read("~/duhturat/lidar/sampling/krokar.shp")## Reading layer `krokar' from data source
## `C:\Users\zanku\OneDrive - Univerza v Ljubljani\Documents\duhturat\lidar\sampling\krokar.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 1 feature and 2 fields
## Geometry type: POLYGON
## Dimension: XYZ
## Bounding box: xmin: 481880.4 ymin: 43555.41 xmax: 482883.1 ymax: 45247.13
## z_range: zmin: 0 zmax: 0
## Projected CRS: MGI 1901 / Slovene National Gridplot(krokar)
plot(grid, add = TRUE)
plot(krok, add = TRUE)
# Razreži oblak točk na vzorčne kvadrante.
krokar_grid <- clip_roi(las = krokar, geometry = grid)
## Chunk 1 of 77 (1.3%): state <U+2713>
## Chunk 2 of 77 (2.6%): state <U+2713>
## Chunk 3 of 77 (3.9%): state <U+2713>
## Chunk 4 of 77 (5.2%): state <U+2713>
## Chunk 5 of 77 (6.5%): state <U+2713>
## Chunk 6 of 77 (7.8%): state <U+2713>
## Chunk 7 of 77 (9.1%): state <U+2713>
## Chunk 8 of 77 (10.4%): state <U+2713>
## Chunk 9 of 77 (11.7%): state <U+2713>
## Chunk 10 of 77 (13%): state <U+2713>
## Chunk 11 of 77 (14.3%): state <U+2713>
## Chunk 12 of 77 (15.6%): state <U+2713>
## Chunk 13 of 77 (16.9%): state <U+2713>
## Chunk 14 of 77 (18.2%): state <U+2713>
## Chunk 15 of 77 (19.5%): state <U+2713>
## Chunk 16 of 77 (20.8%): state <U+2713>
## Chunk 17 of 77 (22.1%): state <U+2713>
## Chunk 18 of 77 (23.4%): state <U+2713>
## Chunk 19 of 77 (24.7%): state <U+2713>
## Chunk 20 of 77 (26%): state <U+2713>
## Chunk 21 of 77 (27.3%): state <U+2713>
## Chunk 22 of 77 (28.6%): state <U+2713>
## Chunk 23 of 77 (29.9%): state <U+2713>
## Chunk 24 of 77 (31.2%): state <U+2713>
## Chunk 25 of 77 (32.5%): state <U+2713>
## Chunk 26 of 77 (33.8%): state <U+2713>
## Chunk 27 of 77 (35.1%): state <U+2713>
## Chunk 28 of 77 (36.4%): state <U+2713>
## Chunk 29 of 77 (37.7%): state <U+2713>
## Chunk 30 of 77 (39%): state <U+2713>
## Chunk 31 of 77 (40.3%): state <U+2713>
## Chunk 32 of 77 (41.6%): state <U+2713>
## Chunk 33 of 77 (42.9%): state <U+2713>
## Chunk 34 of 77 (44.2%): state <U+2713>
## Chunk 35 of 77 (45.5%): state <U+2713>
## Chunk 36 of 77 (46.8%): state <U+2713>
## Chunk 37 of 77 (48.1%): state <U+2713>
## Chunk 38 of 77 (49.4%): state <U+2713>
## Chunk 39 of 77 (50.6%): state <U+2713>
## Chunk 40 of 77 (51.9%): state <U+2713>
## Chunk 41 of 77 (53.2%): state <U+2713>
## Chunk 42 of 77 (54.5%): state <U+2713>
## Chunk 43 of 77 (55.8%): state <U+2713>
## Chunk 44 of 77 (57.1%): state <U+2713>
## Chunk 45 of 77 (58.4%): state <U+2713>
## Chunk 46 of 77 (59.7%): state <U+2713>
## Chunk 47 of 77 (61%): state <U+2713>
## Chunk 48 of 77 (62.3%): state <U+2713>
## Chunk 49 of 77 (63.6%): state <U+2713>
## Chunk 50 of 77 (64.9%): state <U+2713>
## Chunk 51 of 77 (66.2%): state <U+2713>
## Chunk 52 of 77 (67.5%): state <U+2713>
## Chunk 53 of 77 (68.8%): state <U+2713>
## Chunk 54 of 77 (70.1%): state <U+2713>
## Chunk 55 of 77 (71.4%): state <U+2713>
## Chunk 56 of 77 (72.7%): state <U+2713>
## Chunk 57 of 77 (74%): state <U+2713>
## Chunk 58 of 77 (75.3%): state <U+2713>
## Chunk 59 of 77 (76.6%): state <U+2713>
## Chunk 60 of 77 (77.9%): state <U+2713>
## Chunk 61 of 77 (79.2%): state <U+2713>
## Chunk 62 of 77 (80.5%): state <U+2713>
## Chunk 63 of 77 (81.8%): state <U+2713>
## Chunk 64 of 77 (83.1%): state <U+2713>
## Chunk 65 of 77 (84.4%): state <U+2713>
## Chunk 66 of 77 (85.7%): state <U+2713>
## Chunk 67 of 77 (87%): state <U+2713>
## Chunk 68 of 77 (88.3%): state <U+2713>
## Chunk 69 of 77 (89.6%): state <U+2713>
## Chunk 70 of 77 (90.9%): state <U+2713>
## Chunk 71 of 77 (92.2%): state <U+2713>
## Chunk 72 of 77 (93.5%): state <U+2713>
## Chunk 73 of 77 (94.8%): state <U+2713>
## Chunk 74 of 77 (96.1%): state <U+2713>
## Chunk 75 of 77 (97.4%): state <U+2713>
## Chunk 76 of 77 (98.7%): state <U+2713>
## Chunk 77 of 77 (100%): state <U+2713># Izračunaj parametre terena
source("~/duhturat/lidar/pozen_cluster.R")
# Hiter pregled podatkov in VIF analiza za ugotavljanje korelacije med spremenljivkami
summary(terrain_data)## altrange altmean tri_median tri_mean
## Min. : 20.27 Min. : 930.9 Min. :0.1600 Min. :0.1734
## 1st Qu.: 31.97 1st Qu.:1057.6 1st Qu.:0.2200 1st Qu.:0.2539
## Median : 55.33 Median :1107.8 Median :0.3013 Median :0.3461
## Mean : 58.09 Mean :1092.7 Mean :0.3621 Mean :0.4140
## 3rd Qu.: 76.35 3rd Qu.:1140.0 3rd Qu.:0.5125 3rd Qu.:0.5576
## Max. :134.98 Max. :1171.8 Max. :0.6388 Max. :0.8676
## tri_min tri_max tpi_median tpi_mean
## Min. :0.01000 Min. : 0.6887 Min. :-0.00500 Min. :-0.0041252
## 1st Qu.:0.02000 1st Qu.: 1.4525 1st Qu.:-0.00250 1st Qu.:-0.0007767
## Median :0.03750 Median : 3.9413 Median :-0.00125 Median : 0.0007661
## Mean :0.08016 Mean : 10.2766 Mean :-0.00181 Mean : 0.0009881
## 3rd Qu.:0.14500 3rd Qu.: 14.4413 3rd Qu.:-0.00125 3rd Qu.: 0.0019280
## Max. :0.27500 Max. :111.4912 Max. : 0.00500 Max. : 0.0134312
## tpi_min tpi_max roughness_median roughness_mean
## Min. :-43.4025 Min. : 0.5837 Min. :0.53 Min. :0.5635
## 1st Qu.: -7.6575 1st Qu.: 1.1100 1st Qu.:0.73 1st Qu.:0.8389
## Median : -2.1875 Median : 2.2462 Median :0.94 Median :1.1749
## Mean : -5.7106 Mean : 8.6046 Mean :1.20 Mean :1.3737
## 3rd Qu.: -0.9537 3rd Qu.: 10.5062 3rd Qu.:1.65 3rd Qu.:1.8680
## Max. : -0.2950 Max. :111.4912 Max. :2.34 Max. :2.9645
## roughness_min roughness_max slope_median slope_mean
## Min. :0.0000 Min. : 1.93 Min. :10.81 Min. :11.34
## 1st Qu.:0.0500 1st Qu.: 3.91 1st Qu.:15.43 1st Qu.:16.91
## Median :0.0900 Median : 8.69 Median :20.91 Median :21.48
## Mean :0.2174 Mean : 15.26 Mean :24.12 Mean :25.05
## 3rd Qu.:0.4000 3rd Qu.: 20.93 3rd Qu.:33.29 3rd Qu.:34.32
## Max. :0.8100 Max. :113.44 Max. :40.18 Max. :40.85
## slope_min slope_max aspect_median aspect_mean
## Min. : 0.0000 Min. :34.63 Min. : 28.20 Min. : 39.40
## 1st Qu.: 0.2265 1st Qu.:53.58 1st Qu.: 68.24 1st Qu.: 71.79
## Median : 0.7713 Median :70.18 Median : 90.28 Median :109.53
## Mean : 2.9880 Mean :67.14 Mean :110.79 Mean :117.14
## 3rd Qu.: 5.8541 3rd Qu.:82.47 3rd Qu.:146.65 3rd Qu.:154.88
## Max. :14.1503 Max. :88.00 Max. :279.01 Max. :269.04
## aspect_min aspect_max
## Min. : 0.0000 Min. :111.8
## 1st Qu.: 0.0000 1st Qu.:357.7
## Median : 0.1996 Median :359.8
## Mean : 2.6689 Mean :342.8
## 3rd Qu.: 2.7091 3rd Qu.:360.0
## Max. :23.2773 Max. :360.0HH::vif(terrain_data)## altrange altmean tri_median tri_mean
## 11.221140 3.672631 1038.419354 1130.732199
## tri_min tri_max tpi_median tpi_mean
## 10.232008 75.722942 4.898131 6.132053
## tpi_min tpi_max roughness_median roughness_mean
## 13.672421 72.669034 605.451555 1098.449131
## roughness_min roughness_max slope_median slope_mean
## 4.767894 66.987977 960.058451 566.178275
## slope_min slope_max aspect_median aspect_mean
## 5.841128 5.179202 9.066312 13.345459
## aspect_min aspect_max
## 2.437515 2.466287# Še grafično
terrain_data$id <- 1:nrow(terrain_data)
grid <- st_read("~/duhturat/lidar/sampling/grid_50_plus.shp")## Reading layer `grid_50_plus' from data source
## `C:\Users\zanku\OneDrive - Univerza v Ljubljani\Documents\duhturat\lidar\sampling\grid_50_plus.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 77 features and 5 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: 481880.4 ymin: 43547.11 xmax: 482880.5 ymax: 45247.13
## Projected CRS: MGI 1901 / Slovene National Gridgrid_abiotic <- merge(grid, terrain_data, by = "id")
plot(grid_abiotic[, -c(1:5, 28)], max.plot = 21)
Biotski parametri
Izračun lastnosti terena temelji na delu oblaka točk, ki je klasificiran kot tla (ground). V nalogi pa smo se ukvarjali tudi z delom oblaka točk, ki so klasificirane kot vegetacija. Tako smo s pomočjo funkcij, ki so del že omenjenega paketa lidR, zaznali posamezna drevesa ter izračunali njihove višine. Te so bile osnova za razdelitev dreves v tri razrede - razvojne faze - mladovje, drogovnjak, debeljak. Tako kot v primeru parametrov terena, smo tudi drevesa “šteli” v vsakem od 77 vzorčnih kvadrantov.
source("~/duhturat/lidar/biotski_parametri-Zans-ThinkPad.R")
#### DODAJ SENTINEL - MOGOČE TREE HULLS IN SENTINEL DATA COMBINED ####
# Hiter pregled podatkov in VIF analiza za ugotavljanje korelacije med spremenljivkami
summary(tree_data)## tree_count_dalponte mladovje drogovnjak debeljak
## Min. :167.0 Min. :0.000000 Min. :0.000000 Min. :0.3990
## 1st Qu.:303.0 1st Qu.:0.002494 1st Qu.:0.009036 1st Qu.:0.9331
## Median :332.0 Median :0.005682 Median :0.018450 Median :0.9700
## Mean :333.2 Mean :0.013889 Mean :0.053470 Mean :0.9326
## 3rd Qu.:381.0 3rd Qu.:0.016340 3rd Qu.:0.052083 3rd Qu.:0.9873
## Max. :484.0 Max. :0.229885 Max. :0.543269 Max. :1.0000
## tree_count_li mladovje drogovnjak debeljak
## Min. :150.0 Min. :0.000000 Min. :0.00000 Min. :0.2652
## 1st Qu.:278.0 1st Qu.:0.005587 1st Qu.:0.02299 1st Qu.:0.8764
## Median :302.0 Median :0.017241 Median :0.04319 Median :0.9270
## Mean :311.6 Mean :0.027950 Mean :0.08232 Mean :0.8897
## 3rd Qu.:353.0 3rd Qu.:0.035320 3rd Qu.:0.08815 3rd Qu.:0.9695
## Max. :453.0 Max. :0.274448 Max. :0.63258 Max. :0.9975
## mlad drog deb
## Min. :0.01051 Min. :0.01364 Min. :0.1945
## 1st Qu.:0.04056 1st Qu.:0.04583 1st Qu.:0.7722
## Median :0.06254 Median :0.08411 Median :0.8315
## Mean :0.07320 Mean :0.11389 Mean :0.8129
## 3rd Qu.:0.09071 3rd Qu.:0.12890 3rd Qu.:0.9049
## Max. :0.35781 Max. :0.60088 Max. :0.9610HH::vif(tree_data)## tree_count_dalponte mladovje drogovnjak debeljak
## 1.657001e+01 1.745251e+04 1.583279e+05 2.317082e+05
## tree_count_li mladovje drogovnjak debeljak
## 1.428227e+01 -2.000000e-01 -1.666667e-01 -1.428571e-01
## mlad drog deb
## Inf Inf Inftree_data$id <- 1:nrow(tree_data)
grid <- st_read("~/duhturat/lidar/sampling/grid_50_plus.shp")## Reading layer `grid_50_plus' from data source
## `C:\Users\zanku\OneDrive - Univerza v Ljubljani\Documents\duhturat\lidar\sampling\grid_50_plus.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 77 features and 5 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: 481880.4 ymin: 43547.11 xmax: 482880.5 ymax: 45247.13
## Projected CRS: MGI 1901 / Slovene National Gridgrid_biotic <- merge(grid, tree_data, by = "id")
plot(grid_biotic[, -c(1:5)])
Klasifikacija vzorčnih kvadrantov
Združimo abiotski in biotski del podatkov, dodamo sentinel ter izvedemo clustering.
library(stars)
library(corrplot)
grid <- cbind(grid_abiotic[, -c(1:5)], grid_biotic[, -c(1:5)])
# plot(grid[, -c(31,32)], max.plot = 30)
grid.wgs <- st_transform(grid, 4326)
# plot(grid.wgs[, -c(31,32)], max.plot = 30)
# dodamo še sentinel podatke
source("~/duhturat/sentinel/randomForest_classification-Zans-ThinkPad.R")
plot(rf_raster)
zones <- raster::extract(x = rf_raster, y = grid.wgs)
grid.wgs$sentinel_clusters <- unlist(lapply(zones, median))
colnames(grid.wgs)## [1] "altrange" "altmean" "tri_median"
## [4] "tri_mean" "tri_min" "tri_max"
## [7] "tpi_median" "tpi_mean" "tpi_min"
## [10] "tpi_max" "roughness_median" "roughness_mean"
## [13] "roughness_min" "roughness_max" "slope_median"
## [16] "slope_mean" "slope_min" "slope_max"
## [19] "aspect_median" "aspect_mean" "aspect_min"
## [22] "aspect_max" "tree_count_dalponte" "mladovje"
## [25] "drogovnjak" "debeljak" "tree_count_li"
## [28] "mladovje.1" "drogovnjak.1" "debeljak.1"
## [31] "mlad" "drog" "deb"
## [34] "geometry" "geometry.1" "sentinel_clusters"grid.wgs <- grid.wgs[, c(1:33, 36, 34, 35)]
colnames(grid.wgs)## [1] "altrange" "altmean" "tri_median"
## [4] "tri_mean" "tri_min" "tri_max"
## [7] "tpi_median" "tpi_mean" "tpi_min"
## [10] "tpi_max" "roughness_median" "roughness_mean"
## [13] "roughness_min" "roughness_max" "slope_median"
## [16] "slope_mean" "slope_min" "slope_max"
## [19] "aspect_median" "aspect_mean" "aspect_min"
## [22] "aspect_max" "tree_count_dalponte" "mladovje"
## [25] "drogovnjak" "debeljak" "tree_count_li"
## [28] "mladovje.1" "drogovnjak.1" "debeljak.1"
## [31] "mlad" "drog" "deb"
## [34] "sentinel_clusters" "geometry" "geometry.1"plot(grid.wgs, max.plot = 33)
data <- as.data.frame(grid.wgs)[, -c(35, 36)]
vif.res <- HH::vif(data)
# uncor <- vif.res < 4
cors <- cor(data)
corrplot(cors, method = 'square', order = 'FPC', type = 'lower', diag = FALSE)
colnames(data)## [1] "altrange" "altmean" "tri_median"
## [4] "tri_mean" "tri_min" "tri_max"
## [7] "tpi_median" "tpi_mean" "tpi_min"
## [10] "tpi_max" "roughness_median" "roughness_mean"
## [13] "roughness_min" "roughness_max" "slope_median"
## [16] "slope_mean" "slope_min" "slope_max"
## [19] "aspect_median" "aspect_mean" "aspect_min"
## [22] "aspect_max" "tree_count_dalponte" "mladovje"
## [25] "drogovnjak" "debeljak" "tree_count_li"
## [28] "mladovje.1" "drogovnjak.1" "debeljak.1"
## [31] "mlad" "drog" "deb"
## [34] "sentinel_clusters"# select non-correlated variables
# df <- data[, c(uncor)]
# df$deb <- data$deb
# df$tree_count <- data$tree_count_dalponte
df <- data[, c(1,2,3,7,19, 27, 33, 34)]
df <- as.data.frame(scale(df))
colnames(df)## [1] "altrange" "altmean" "tri_median"
## [4] "tpi_median" "aspect_median" "tree_count_li"
## [7] "deb" "sentinel_clusters"HH::vif(df)## altrange altmean tri_median tpi_median
## 3.120426 2.277998 3.651402 1.395247
## aspect_median tree_count_li deb sentinel_clusters
## 1.792459 1.359248 1.681289 1.339560d <- dist(df, method = "manhattan")
fit <- hclust(d, method = "ward.D2")
plot(fit)
knum <- 3
groups <- cutree(fit, k = knum)
rect.hclust(fit, k = knum, border = "red")
fil <- data.frame(id = 1:length(groups), group = factor(groups))
pr <- prcomp(df, center = TRUE, scale. = TRUE)
bplot <- ggbiplot::ggbiplot(pr, obs.scale = 1, groups = fil$group, var.scale = 1,
ellipse = TRUE, circle = TRUE, labels = rownames(fil)) +
theme_bw() +
scale_color_brewer(palette = "Dark2")
bplot
# pca3d::pca3d(pr, group = fil$group)
table(groups)## groups
## 1 2 3
## 42 32 3grid$group <- groups
centros <- data.frame(geosphere::centroid(as_Spatial(grid)))
centros$id <- 1:nrow(centros)
gridplot <- ggplot() +
geom_sf(data = krok, alpha = 0.3) +
geom_sf(data = grid, aes(fill = factor(group)), alpha = 0.8) +
geom_text(data = centros, aes(x = X1, y = X2, label = id)) +
theme_bw() +
scale_fill_brewer(palette = "Dark2")
gridplot
theme_set(theme_pubr())
ggarrange(bplot, gridplot, ncol = 2, nrow = 1)
# coords <- st_geometry(grid)
# out <- st_drop_geometry(grid)
# st_geometry(out) <- coords
# st_write(out, "grid_data.kml", delete_layer = TRUE, delete_dsn = TRUE)Prikaz kazalnikov NDVI in EVI za vse obhode.
ndvi.stars <- st_as_stars(ndvi)
# ggplot()+
# geom_stars(data = ndvi.stars) +
# coord_equal() +
# facet_wrap(~band) +
# theme_void() +
# scale_x_discrete(expand=c(0,0))+
# scale_y_discrete(expand=c(0,0)) +
# scale_fill_gradient(low = "#f7fcb9", high = "#006837")
evi1.stars <- st_as_stars(evi1)
# ggplot()+
# geom_stars(data = evi1.stars) +
# coord_equal() +
# facet_wrap(~band) +
# theme_void() +
# scale_x_discrete(expand=c(0,0))+
# scale_y_discrete(expand=c(0,0)) +
# scale_fill_gradient(low = "#f7fcb9", high = "#006837")
plot(evi1.stars)
plot(ndvi.stars)
plot(rf_raster) Izračunov povprečnega števila dreves in deleža debeljaka na podlagi obeh uporabljenih algoritmov (dalponte, li) ter izračun deleža debeljaka na podlagi modela višin krošenj (CHM)
mean(grid$tree_count_dalponte)## [1] 333.1948sd(grid$tree_count_dalponte)## [1] 72.20075sd(grid$tree_count_dalponte) / sqrt(nrow(grid))## [1] 8.228039mean(grid$tree_count_li)## [1] 311.6234sd(grid$tree_count_li)## [1] 66.06796sd(grid$tree_count_li) / sqrt(nrow(grid))## [1] 7.529142mean(grid$debeljak)## [1] 0.9325992sd(grid$debeljak)## [1] 0.1100366sd(grid$debeljak) / sqrt(nrow(grid))## [1] 0.01253983mean(grid$debeljak.1)## [1] 0.8896934sd(grid$debeljak.1)## [1] 0.1389572sd(grid$debeljak.1) / sqrt(nrow(grid))## [1] 0.01583564mean(grid$deb)## [1] 0.8129109sd(grid$deb)## [1] 0.1359584sd(grid$deb) / sqrt(nrow(grid))## [1] 0.01549389