Introduction à la géostatistique avec R

Guillaume Larocque (glaroc@gmail.com), 27 février 2020. Université de Montréal.

Présentation sur Prezi

Présentation simultanée

• Applied spatial data analysis with R, Bivand, Pebesma et Gómez-Rubio. - Bon survol des outils d'analyse spatiale dans R et un bon chapitre sur la géostatistique.
• Geostatistics for Natural Resources Evaluation, Pierre Goovaerts. - Très bon livre pour une connaissance approfondie de la géostatistique.
• Georges Matheron - La théorie des variables régionalisées et ses applications (PDF). - Text original décrivant les bases de la géostatistique.
• Analyses spatiales sous R - Document préparé par Nicolas Casajus pour un cours offert à l'UQAR à l'hiver 2013.

• Geostatistical Tools for Modeling and Interpreting Ecological Spatial Dependence
• GIS and geostatistics: Essential partners for spatial analysis
• Spatial autocorrelation and statistical tests in ecology
• Spatial covariance in plant communities: intergrating ordination, geostatistics and variance testing
• Coregionalization analysis with a drift for multi-scale assessment of spatial relationships between ecological variables

• Site web de Gstat.
• Gstat package - Gstat R package.
• geoR package geoR - Autre package pour la géostatistique dans R.

Cliquez ici pour télécharger le fichier zip contenant les scripts et autres fichiers qui seront utilisés pour l'atelier. Veuillez extraire ce fichier dans un dossier facile d'accès sur votre ordinateur.

Vous pouvez numériser le périmètre de votre aire d'étude à partir d'une couche Google Satellite en suivant les étapes suivantes:

1. Ouvrez QGIS.
2. Activez l'extension QuickMapServices.
3. Ajoutez une couche 'Google Satellite' en utilisant l'extension QuickMapServices (>web>QuickMapServices). Notez que le SCR du canevas est maintenant WGS 84 / Pseudo mercator. Déplacez cette couche en dessous du canevas en cliquant-glissant le nom de la couche vers le bas dans la table des couches à gauche. Si vous avez d'autres couches vecteur ou raster pouvant vous aider à définir les limites de la zone d'étude, ajoutez les au canevas an cliquant sur les icônes appropriés (ajouter couche vecteurs ou raster).
4. Numérisez la zone définissant votre aire d'étude. Vous devez d'abord ajouter une nouvelle couche shapefile de type polygone sous Couche>Nouveau>Nouvelle couche shapefile. Spécifiez le système de référence approprié selon l'échelle de vos données et donnez à ce fichier shapefile un nom approprié. Il est toujours préférable de travailler avec un système de référence avec des unités en mètres comme UTM, MTM ou Lambert, et non pas en degrés comme latitude/longitude. Activez le mode édition (­>Couche­­­­>­Basculer en mode édition)­­­­­­­­­­ et numérisez le polygone en cliquant sur l'icône et spécifiez '1' comme identifiant (ID).

1. Trouvez la fonction Rasterize sous Raster>Conversion. Spécifiez le fichier polygone que vous venez de créer à l'étape précédente comme fichier entrant, et spécifiez un nom de fichier avec une extension .tif comme fichier sortant. Choisissez 'Résolution exprimée en unité de la carte par pixel'. Vous devez ensuite spécifier la taille désirée de chaque pixel sur la carte en unités du système de référence du fichier polygone en entrée. Par exemple, si le fichier polygone est en coordonnées UTM, cette résolution est en mètre. Assurez-vous de choisir une résolution raisonnable vous permettant d'obtenir un fichier raster avec un nombre de pixels gérable (par.ex. quelques centaines de lignes et de colonnes, maximum 2000-3000). Par exemple, si votre aire d'étude fait environ 2km x 3km, vous pouvez spécifier 5 (mètres) comme résolution de façon à avoir un fichier raster d'environ 400 pixels x 600 pixels.
2. Maintenant, ce fichier tif est enregistré sur votre ordinateur. Dans votre aire d'étude, tous les pixels ont une valeur de 1, alors que les pixels à l'extérieur ont une valeur de 0.

1. Vous pouvez maintenant importer le masque dans R:

mask<-readGDAL(file.choose())

et choisissez le fichier raster en format tif.

ensuite, vous pouvez faire

mask@data[mask@data==0]=NA

pour changer tous les zéros en valeurs NA (nulles).

Votre masque est maintenant prêt à être utilisé!