8.2 Exercices
Dans cette section, vous mettrez en pratique certains concepts vus dans la section leçon de ce module. Bien que la réponse à chaque question soit disponible, il est très important de tenter d’y répondre par vous même!
Question 1
Identifier le point d’élévation maximale sur une carte du parc du Mont-Orford.
Réponse
Isoler le polygone du parc du Mont-Orford du shapefile parcs
en utilisant la fonction subset
:
subset(parcs, TRQ_NM_== "Parc national du Mont-Orford") parc_orford <-
Filtrer le raster dem_lcc
pour ne retenir que les cellules situées à l’intérieure du polygone délimitant les limites du parc du Mont-Orford
mask(dem_lcc, parc_orford)
dem_orford <-mapview(dem_orford)
Déterminer l’indice de la cellule ou des cellules d’élévation maximale en utilisant les fonctions getValues()
et which.max()
.
which.max(getValues(dem_orford)) imax_orford<-
Trouver les coordonnées correspondant à cet indice en utilisant la fonction xyFromCell
.
xyFromCell(dem_orford,imax_orford)
coordmax_orford <- coordmax_orford
x y
[1,] -293815 154024
Transformer ces coordonnées en une donnée vectorielle de type POINT
en utilisant la fonction st_as_sf()
.
#Créer d'abord un data.frame à partir de coordmax_orford
as.data.frame(coordmax_orford) #ou simplement data.frame(coordmax_orford)
coordmax_orford_df <-
#Créer un point
st_as_sf(coordmax_orford_df,
pt_max_orford <-coords = c("x", "y"),
crs = st_crs(dem_orford)) #attribuer le même SCR
Visualiser la carte du parc du Mont-Orford ainsi que le point d’élévalion maximale.
mapview(dem_orford) + mapview(pt_max_orford, col.regions="red")
Question 2
Reclassifier le raster dem_orford
selon quatre classes correspondant aux valeurs comprises entre le zéro et le \(1^{er}\) quantile, le \(1^{er}\) et le \(2^{ième}\), le \(2^{ième}\) et \(3^{ième}\), et le \(3^{ième}\) et \(4^{ième}\).
Réponse
Trouver d’abord les quantiles des valeurs d’élévation dans le parc du Mont-Orford.
quantile(getValues(dem_orford), na.rm = TRUE)
quantile_orford <- quantile_orford
0% 25% 50% 75% 100%
270.1 338.6 394.2 536.6 816.6
25 <- as.numeric(quantile_orford["25%"]) # as.numeric pour conserver seulement les chiffres
q_50 <- as.numeric(quantile_orford["50%"])
q_75 <- as.numeric(quantile_orford["75%"])
q_100 <- as.numeric(quantile_orford["100%"]) q_
En se servant de ces valeurs, construire une matrice qui détermine la nouvelle classification.
matrix(c(0, q_25, 1, #1ere classe
classes_orford <-25, q_50, 2, #2e classe
q_50, q_75, 3, #3e classe
q_75, q_100, 4), #4e classe
q_nrow = 4, ncol = 3, byrow = TRUE)
# donner des titres aux colonnes de cette matrice
colnames(classes_orford) <- c("Limite_min", "Limite_max", "Classement_quantile")
Utiliser la fonction reclassify()
pour assigner les nouvelles classes au raster dem_orford
.
reclassify(dem_orford, classes_orford, rigth = FALSE) dem_orford_quantile <-
Confirmer visuellement la reclassification en utilisant la fonction mapview()
mapview(dem_orford_quantile)
Noter que nous aurions pu vouloir attribuer une nouvelle classification de type catégorique au raster dem_orford
. Par exemple,
c("Faible", "Intermediare", "Forte", "Très forte")
Categorie_elevation <-
# Ou encore
c("[0%, 25%[", "[25%, 50%[", "[50%, 75%[", "[75%, 100%]") Categorie_quantile <-
Dans ce cas, il faut ajouter ces catégories au raster dem_orford
en utilisant la fonction levels()
.
levels(dem_orford_quantile) <- data.frame(ID=1:4, quantile = Categorie_quantile)
dem_orford_quantile
class : RasterLayer
dimensions : 178, 304, 54112 (nrow, ncol, ncell)
resolution : 521, 741 (x, y)
extent : -336797, -178413, 110676, 242574 (xmin, xmax, ymin, ymax)
crs : +proj=lcc +lat_0=44 +lon_0=-68.5 +lat_1=60 +lat_2=46 +x_0=0 +y_0=0 +datum=NAD83 +units=m +no_defs
source : memory
names : DEM
values : 1, 4 (min, max)
attributes :
ID quantile
1 [0%, 25%[
2 [25%, 50%[
3 [50%, 75%[
4 [75%, 100%]
La classification est maintenant traitée comme une variable catégorique.
mapview(dem_orford_quantile)
Question 3
À partir du raster dem_lcc
, créer un seul raster composé de zones tampons circulaires, de 10 km de rayon, autour des points d’élévation maximale de chaque parc national de la région de Sherbrooke.
Réponse
Répéter les étapes de la Question 1 pour créer un shapefile correspondant au point d’élévation maximale de chaque parc.
Pour le Parc du Mont-Orford:
# Isoler les limites du parc
subset(parcs, TRQ_NM_== "Parc national du Mont-Orford")
parc_orford <-# Trouver dem à l'intérieur du parc
mask(dem_lcc, parc_orford)
dem_orford <-# Trouver l'indice de la cellule d'élévation maximale
which.max(getValues(dem_orford))
imax_orford<-# Trouver les coordonnées de cette cellule
xyFromCell(dem_orford,imax_orford)
coordmax_orford <-# Transformer en data.frame
as.data.frame(coordmax_orford)
coordmax_orford_df <-# Transformer en point
st_as_sf(coordmax_orford_df,
pt_max_orford <-coords = c("x", "y"),
crs = st_crs(dem_orford)) #attribuer le même SCR
Pour le Parc national du Mont-Mégantic:
subset(parcs, TRQ_NM_== "Parc national du Mont-Mégantic")
parc_megantic <- mask(dem_lcc, parc_megantic)
dem_megantic <- which.max(getValues(dem_megantic))
imax_megantic <- xyFromCell(dem_megantic,imax_megantic)
coordmax_megantic <- as.data.frame(coordmax_megantic)
coordmax_megantic_df <- st_as_sf(coordmax_megantic_df,
pt_max_megantic <-coords = c("x", "y"),
crs = st_crs(dem_orford))
Pour le Parc national de la Yamaska:
subset(parcs, TRQ_NM_== "Parc national de la Yamaska")
parc_yamaska <- mask(dem_lcc, parc_yamaska)
dem_yamaska <- which.max(getValues(dem_yamaska))
imax_yamaska <- xyFromCell(dem_yamaska,imax_yamaska)
coordmax_yamaska <- as.data.frame(coordmax_yamaska)
coordmax_yamaska_df <- st_as_sf(coordmax_yamaska_df,
pt_max_yamaska <-coords = c("x", "y"),
crs = st_crs(dem_orford))
Pour le Parc national de Frontenac:
subset(parcs, TRQ_NM_== "Parc national de Frontenac")
parc_frontenac <- mask(dem_lcc, parc_frontenac)
dem_frontenac <- which.max(getValues(dem_frontenac))
imax_frontenac <- xyFromCell(dem_frontenac,imax_frontenac)
coordmax_frontenac <- as.data.frame(coordmax_frontenac)
coordmax_frontenac_df <- st_as_sf(coordmax_frontenac_df,
pt_max_frontenac <-coords = c("x", "y"),
crs = st_crs(dem_orford))
Confirmer vos calculs en visualisant les points trouvés.
mapview(pt_max_orford, col.regions = "red") +
mapview(pt_max_megantic, col.regions = "blue") +
mapview(pt_max_yamaska, col.regions = "green") +
mapview(pt_max_frontenac, col.regions = "yellow")
mapview(pt_max_orford, col.regions = "red") +
m <- mapview(pt_max_megantic, col.regions = "blue") +
mapview(pt_max_yamaska, col.regions = "green") +
mapview(pt_max_frontenac, col.regions = "yellow")
@map m
Calculer une zone tampon de 10 km de rayon autour de chaque point en utilisant la fonction st_buffer()
.
st_buffer(pt_max_orford, dist = 10e3) #nous savons que le CRS est d'unité métrique
tampon_orford <- st_buffer(pt_max_megantic, dist = 10e3)
tampon_megantic <- st_buffer(pt_max_yamaska, dist = 10e3)
tampon_yamaska <- st_buffer(pt_max_frontenac, dist = 10e3) tampon_frontenac <-
Filtrer le raster dem_lcc
pour ne retenir que les cellules comprises à l’intérieur des zones tampons grâce à la fonction mask()
.
mask(dem_lcc, tampon_orford)
dem_tampon_orford <- mask(dem_lcc, tampon_megantic)
dem_tampon_megantic <- mask(dem_lcc, tampon_yamaska)
dem_tampon_yasmaka <- mask(dem_lcc, tampon_frontenac) dem_tampon_frontenac <-
Combiner les quatre raster en un seul en utilisant la fonction merge()
merge(dem_tampon_orford,
dem_tampon_max <-merge(dem_tampon_megantic,
merge(dem_tampon_yasmaka, dem_tampon_frontenac)))
Noter que la fonction merge()
s’utilise sur deux rasters à la fois. Il faut donc l’appliquer successivement pour combiner quatre rasters.
Visualiser le raster final contenant les quatre zones tampons:
mapview(dem_tampon_max)