Analyse de séries temporelles de NDVI

Cours réseaux informatiques, tutoriel & guide de travaux pratiques en pdf.

Comblement des valeurs manquantes

Une fois l’image ‘nettoyée’, il reste la question du traitement des données manquantes issues du masque. Selon les sites et les saisons, celles-ci peuvent être importantes et il peut être nécessaire de les combler. En résumé, il existe trois approches pour le remplissage des données manquantes (Weiss et al. 2014).
La première est basée sur l’information spatiale et consiste à utiliser une seule image source (ex. Rossi et al. 1994). Les vides peuvent être reconstruit par ‘krigeage’ ou par interpolation traditionnelle en utilisant d’autres régions de la même image. Différents algorithmes (recherche des plus proches voisins, algorithme bilinéaire, interpolation bilinéaire, interpolation bicubique, etc.) peuvent être utilisés à cet effet, bien qu’en général, cette approche produit des résultats insatisfaisants (Jabar et al. 2014).
La seconde approche est basée sur l’information temporelle et requiert plusieurs images sources. Une analyse de Fourier peut être réalisée sur plusieurs images afin de modéliser une série temporelle d’observations pixelliques permettant de combler les données manquantes (ex. Roerink et al. 2000). Aussi, d’autres algorithmes (fonction gaussienne, filtre de Savitzky-Golay, etc.) permettent de lisser des séries temporelles selon plusieurs paramètres (ex. Kan-dasamy et al. 2013). D’autres méthodes consistent simplement à remplacer une donnée man-quante par une donnée existante provenant d’autres images (ex. Surya & Simon 2015) ou d’une seule image de référence (ex. Jin et al. 2013).
Finalement l’approche mixte combine l’information spatiale et temporelle. Par exemple Kang et al. (2005) ont développé un algorithme utilisant une interpolation temporelle lorsque l’information spatiale n’est pas suffisamment précise pour combler les vides.
Plusieurs facteurs dépendant du site, de l’image, de la rapidité de calcul souhaitée, de l’analyse à effectuer, etc. sont à prendre en compte pour choisir l’une ou l’autre méthode. Par exemple, l’utilisation d’images multi-temporelles doit être considérée avec précaution pour le calcul du NDVI, compte tenu des variations inter-annuelles de ce dernier (ex. Gamon et al. 2013; Olusegun & Adeyewa 2013). Aussi, une méthode simple d’interpolation spatiale pourrait s’avérer suffisante pour une image où les nuages sont de tailles réduites.
Voici une liste non-exhaustive des différents programmes/fonctions disponibles associés à la problématique des données manquantes dans le cadre d’analyses d’images satellites.
– r.fillnulls (GRASS) comble les données manquantes en utilisant une interpolation spline des pixels avoisinants.
– Le script Python gdal_fillnodata.py de la librairie GDAL comble les données man-quantes par interpolation des pixels avoisinants.
– t.rast.gapfill (GRASS) utilise une interpolation linéaire entre plusieurs couches d’une série temporelle et comble les données manquantes de chaque couche.
– La fonction TSGFlinear du package ‘greenbrown’ (R) comble les données manquantes d’une série temporelle par interpolation linéaire.
– La fonction sgolayfilt du package ‘signal’ (R) permet de lisser les données en exécutant un filtre de Savitzy-Golay sur une série temporelle.
– La fonction approxNA du package ‘raster’ (R) estime la valeurs des données manquantes par interpolation d’une série temporelle.
Le choix de la méthode permettant de combler les données manquantes dépend de la nature des inputs (une seule image ou série temporelle), de leur qualité et des objets d’analyses postérieures (constants ou fluctuants) et de la praticité de la méthode. Par exemple, si la taille des patchs de données manquantes est importante, une interpolation spatiale ne serait pas recommandée. De même, si l’objet d’analyse fluctue en fonction des saisons (ex. NDVI) un lissage des données sur une série temporelle aurait tendance à sous-estimer les variations intra-annuelles. Dans le cadre d’analyses de NDVI sur des séries temporelles, tel que présenté dans le projet LiMES, la fonction R approxNA du package ‘raster’ représente un solution adaptée pour les raisons mentionnées ci-dessus. En effet, elle effectue une interpolation linéaire entre les données les plus proches d’un point de vue temporel et permet de conserver le caractère oscillatoire de l’indice si les images ne sont pas trop éloignées les unes des autres. La fiabilité des images issus de l’interpolation dépend de la quantité de données existantes pour chaque pixel. Quoiqu’il en soit l’interpolation linéaire peut potentiellement provoquer une sous-estimation des variations (Figure 25).

Script

Description

Un script Python intégrant certaines fonctions R offre à l’utilisateur la possibilité d’obtenir une série d’output complémentaires ainsi qu’une fiche de résumé en format PDF lui permettant d’appréhender l’évolution du NDVI au cours du temps sous plusieurs angles. Les output offerts à l’utilisateur incluent entre autre une représentation de la déviation standard des valeurs de NDVI permettant de mesurer la variabilité intra-annuelle de la végétation (ex. Hall-Beyer 2012), une représentation des coefficients de pente significatives issus de régressions linéaires effectuées sur chaque pixel (ex. Hartemink et al. 2008) ainsi qu’un graphique en aires illustrant l’évolution de la proportion des différentes classes de valeurs de NDVI au cours du temps. Toutes les analyses sont effectuées une première fois sans traiter les données manquantes, et une deuxième fois en leur assignant une valeur de NDVI via une interpolation linéaire tenant compte de la distance temporelle entre chaque scène. Ceci offre à l’utilisateur une certaine transparence et lui permet de définir, le cas échéant, la prudence avec laquelle les résultats complets (sans données manquantes) doivent être interprétés. La liste exhaustive des output issus de ce script est la suivante :
1. Un fichier ENVI multibande des valeurs de NDVI correspondant à chaque scène (par défaut)
2. Un raster (TIF) et une image (PNG) représentant les valeurs de déviation standard pour chaque pixel (argument optionnel -v)
3. Un raster (TIF) et une image (PNG) représentant les valeurs significatives de coefficient de pente (p-value < 0.05) issus des régressions linéaires réalisées sur chaque pixel (argument optionnel -s).
4. Un raster (TIF) et une image (PNG) représentant les différences de valeurs de NDVI entre la dernière et la première scène calculées sur chaque pixel (argument optionnel -d).
5. Un graphique (PNG) illustrant l’évolution de la distribution des valeurs de NDVI au cours du temps associé à un un fichier texte (TXT) contenant les données utilisées pour générer le graphique (argument optionnel -p).
6. Un graphique en aires (PNG) illustrant l’évolution de la proportion de différentes classes de valeurs de NDVI au cours du temps associé à un un fichier texte (TXT) contenant les données utilisées pour générer le graphique (argument optionnel -a).
7. Un raster (TIF) et une image (PNG) où la valeur de chaque pixel correspond au pourcentage du nombre de fois où celui-ci a été masqué, ceci permettant d’évaluer la fiabilité des données pour différentes zones des rasters en output (argument optionnel -m).
8. Un fichier PDF résumant l’analyse effectuée (données, résultats, etc.) et incluant les output précédents (par défaut).

Structure et processus

Le script est configuré de manière à être exécuté depuis le dossier LiMES contenant les différents scripts (LiMES_batcher, LiMES_grasser, etc.) et requiert l’exécution préalable du script LiMES_grasser afin d’obtenir pour chaque scène les réflectances des bandes 4 et 5 ainsi que les masques de nuages issus de la fonction Fmask. L’argument principal à fournir lors de l’exécution du script est le nom du site sur lequel l’analyse du NDVI va être réalisée. Par défaut, le script identifie les couches à traiter en fonction des suffixes T4, T5 et Fmask. Si ces derniers diffèrent, une modification du code permet de laisser à l’utilisateur la possibilité de définir les suffixes correspondant aux différents rasters (lignes de code commentées). Une fois les différentes scènes identifiées et les différentes variables établies, une boucle s’effectue sur les différentes scènes, laquelle permet d’appliquer les masques et d’effectuer les calculs de NDVI. Ceci réalisé, un RasterStack (collection de rasters ayant la même extension et la même résolution) est créé. Les différents output optionnels sont créés en fonction des choix de l’utilisateur d’abord une première fois sans traiter les données manquantes, puis une deuxième fois en ayant réalisé une interpolation linéaire afin de combler ces dernières. Le script a été pensé de manière à automatiser au maximum les processus et le format des output, notamment le format des graphiques où les titres, les labels et les valeurs de l’axe des abscisses sont obtenus par le nom des différentes scènes (eux-même contenant l’année et le jour julien de celle-ci). Tous les output ainsi que la fiche de résumé en format PDF sont exportés. Le script complet figure dans l’Annexe E.1.
Le module Python rpy2 importé dans le script principal (ndvi.py) permet d’intégrer des fonctions R dans le langage Python. Il permet aussi de compiler des fonctions R pures (pos-sible grâce à la fonction Signature Translated Anonymous Package). Ceci est d’autant plus utile que cette dernière permet de profiter des subtilités de R et de contourner certains obstacles liés aux différences de syntaxe. C’est pour cette raison qu’un script auxiliaire (RFunctions.py) a été élaboré contenant différentes fonctions codées en langage R (voir Annexe E.2). Aussi, le module Python pylatex importé dans le script principal (ndvi.py) permet de créer une interface entre Python et LATEXet ainsi de générer des fichiers PDF. Le code permettant cela est écrit dans un script auxiliaire (ndviPDF.py) figurant dans l’Annexe E.3. Les ordinagrammes structurel (Figure 26) et fonctionnel (Figure 27) résument le processus d’automatisation des analyses.

Exemple

Exécution

Le script ndvi.py a été exécuté sur une série de 8 scènes du site geneva_2 aux couvertures nuageuses variables (période avril 2013 – novembre 2015) (La liste des différentes scènes utilisées figure dans la Section 1 de la fiche de résumé de l’analyse, Annexe E.4). Tous les arguments possibles ont été ajoutés de manière à obtenir une vision d’ensemble des différents output possibles. Le temps total nécessaire à l’exécution du script a été de 3’648 secondes (1.4 million de pixels) et les temps d’exécution de chaque analyse figurent dans l’Annexe E.4 4.

Déviation standard des valeurs de NDVI

L’argument optionnel (-v) permet d’obtenir un raster représentant le degré de variation du NDVI au cours du temps. Une déviation standard des valeurs de NDVI pour chaque pixel est calculée sur la base des 8 scènes. Ceci permet entre autre d’identifier les zones où la végétation n’est pas stable, comme notamment les zones agricoles (Figure 28).
(a) Données manquantes non-traitées. (b) Données manquantes comblées par inter-polation linéaire.
Figure 28 – Déviation standard des valeurs de NDVI pour le site de geneva_2.
Le raster identifie les pixels où le NDVI varie, notamment de manière saisonnière (ex. les champs de culture). De même, les zones qui ne présentent pas de variation (zones urbaines) sont clairement identifiées. Il est intéressant de voir que le Lac Léman présente un certain degré de variation, probablement dû au cycle de croissance des algues. Il est à noter que la face orientée nord-ouest du Salève présente des valeurs de déviation standard relativement élevées. Ceci est probablement causé par l’ombrage, lequel oscille et fausse les valeurs de NDVI.