# Radiometric landscape: a new conceptual framework and operational approach for landscape characterisation and mapping

Louise Lemettais, Samuel Alleaume, Sandra Luque, Anne-Élisabeth Laques, Yonas Alim, Laurent Demagistri, Agnès Bégué

<figure><img src="image_1.png" /></figure>

## To cite this version:

Louise Lemettais, Samuel Alleaume, Sandra Luque, Anne-Élisabeth Laques, Yonas Alim, et al. Radiometric landscape: a new conceptual framework and operational approach for landscape characterisation and mapping. Geo-spatial Information Science, inPress, pp.1-23. 10.1080/10095020.2024.2314558 . ird-04514513

HAL Id: ird-04514513

https://ird.hal.science/ird-04514513

Submitted on 26 Mar 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>
# Radiometric landscape: a new conceptual framework and operational approach for landscape characterisation and mapping

Louise Lemettais, Samuel Alleaume, Sandra Luque, Anne-Élisabeth Laques, Yonas Alim, Laurent Demagistri, Agnès Bégué

<figure><img src="image_1.png" /></figure>

## To cite this version:

Louise Lemettais, Samuel Alleaume, Sandra Luque, Anne-Élisabeth Laques, Yonas Alim, et al. Radiometric landscape: a new conceptual framework and operational approach for landscape characterisation and mapping. Geo-spatial Information Science, In press, pp.1-23. 10.1080/10095020.2024.2314558 . ird-04514513

HAL Id: ird-04514513

https://ird.hal.science/ird-04514513v1

Submitted on 26 Mar 2024

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_3.png" /></figure>

# Radiometric landscape: a new conceptual framework and operational approach for landscape characterisation and mapping

Louise Lemettais, Samuel Alleaume, Sandra Luque, Anne-Élisabeth Laques, Yonas Alim, Laurent Demagistri & Agnès Bégué

**To cite this article:** Louise Lemettais, Samuel Alleaume, Sandra Luque, Anne-Élisabeth Laques, Yonas Alim, Laurent Demagistri & Agnès Bégué (19 Mar 2024): Radiometric landscape: a new conceptual framework and operational approach for landscape characterisation and mapping, *Geo-spatial Information Science*, DOI: 10.1080/10095020.2024.2314558

**To link to this article:** https://doi.org/10.1080/10095020.2024.2314558

<figure><img src="image_4.png" /></figure>

© 2024 Wuhan University. Published by Informa UK Limited, trading as Taylor & Francis Group.

<figure><img src="image_5.png" /></figure>

Published online: 19 Mar 2024.

<figure><img src="image_6.png" /></figure>

Submit your article to this journal

<figure><img src="image_7.png" /></figure>

Article views: 168

<figure><img src="image_8.png" /></figure>

View related articles

<figure><img src="image_9.png" /></figure>

View Crossmark data

CrossMark
<figure><img src="image_3.png" /></figure>

# Radiometric landscape: a new conceptual framework and operational approach for landscape characterisation and mapping

Louise Lemettais, Samuel Alleaume, Sandra Luque, Anne-Élisabeth Laques, Yonas Alim, Laurent Demagistri & Agnès Bégué

**To cite this article:** Louise Lemettais, Samuel Alleaume, Sandra Luque, Anne-Élisabeth Laques, Yonas Alim, Laurent Demagistri & Agnès Bégué (19 Mar 2024): Radiometric landscape: a new conceptual framework and operational approach for landscape characterisation and mapping, *Geo-spatial Information Science*, DOI: 10.1080/10095020.2024.2314558

**To link to this article:** https://doi.org/10.1080/10095020.2024.2314558

<figure><img src="image_4.png" /></figure>

© 2024 Wuhan University. Published by Informa UK Limited, trading as Taylor & Francis Group.

<figure><img src="image_5.png" /></figure>

Published online: 19 Mar 2024.

<figure><img src="image_6.png" /></figure>

Submit your article to this journal

<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /></figure>

Article views: 168

<figure><img src="image_9.png" /></figure>

View related articles

<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /></figure>

View Crossmark data

<figure><img src="image_12.png" /></figure>

CrossMark
<figure><img src="image_10.png" /></figure>

# Radiometric landscape: a new conceptual framework and operational approach for landscape characterisation and mapping

Louise Lemettais IDa, Samuel Alleaume IDb,e, Sandra Luque Id,b,e, Anne-Élisabeth Laques Id,a,d, Yonas Alimc,e, Laurent Demagistri Id,a and Agnès Bégué Id,c,e

aUMR ESPACE-DEV, Univ Montpellier, IRD, Univ Antilles, Univ Guyane, Univ Réunion, Montpellier, France; bINRAE, UMR TETIS, Montpellier, France; cCIRAD, UMR TETIS, Montpellier, France; dLMI PAYSAGES, IRD, CNRE, Université d'Antananarivo, Antananarivo, Madagascar; eUMR TETIS, AgroParisTech CIRAD CNRS INRAE Univ Montpellier, Montpellier, France

## ABSTRACT

Landscape mapping has the potential to address some of the most pressing research issues of our time, including climate change, sustainable development, and human well-being. In this paper, we propose an original method that lays the foundations for landscape mapping and overcomes some of the major limitations of existing biophysical methods. Based on the assumption that the primary components of the landscape can be extracted directly from the radiometric information of satellite image time series, this paper presents a new approach to landscape characterization and mapping based solely on remote sensing data. The approach relies on a conceptual model, which links the description, characteristics, structure and functions of the landscape to a set of Remote Sensing-based Essential Landscape Variables (RS-ELVs). The RS-ELVs are then processed according to geographic object-based image analysis (GEOBIA) approach to produce a radiometric landscape map. The model and the remote sensing data processing chain are tested on a case study in central Madagascar (about 13 000 km²) composed of contrasting landscapes resulting from different climatic conditions and agricultural practices. The RS-ELVs are extracted from MODIS image time series for the temporal and spectral variables, and from MODIS and Sentinel-2 images for the texture variables. The parameterization of the segmentation and clustering algorithms is determined by statistical optimization. The final result is a radiometric landscape map in six classes. The landscape classes are then characterized using an independent set of remote sensing variables, a global land cover map and ground observations. The approach successfully identifies and delineates the gradient and major landscape types of the complex region of central Madagascar, confirming our initial hypothesis. The production of such radiometric landscape maps opens the way for integrated territorial development, including the planning and protection of the living environment and human well-being, and the implementation of sectoral policies.

## ARTICLE HISTORY

Received 31 July 2023
Accepted 31 January 2024

## KEYWORDS

Remote sensing; MODIS;
Sentinel-2; essential
variables; satellite image
time series; Madagascar

## 1. Introduction

Due to the multidisciplinary nature of landscape studies and their almost infinite variability, the simplification of landscapes into spatial units is a complex task, but nevertheless important for communication in management and research (e.g. Bunce et al. 1996; Hazeu et al. 2011; Turner and Gardner 2015). Landscape mapping is a cornerstone of the integrated territorial development, including the planning and protection of the living environment and human well-being (e.g. Reed, Ros-Tonen, and Sunderland 2020), and the implementation of sectoral policies such as the green and blue infrastructures (Bourget and le Dû-Blayo 2010; Sun, Wei, and Han 2022). Landscape characterization and mapping is also important in biogeography and ecology (e.g. Cullum et al. 2016) to guide the setting of conservation priorities, or to better understand and manage biodiversity (e.g. Préau et al.

2022). More generally, landscape mapping is very valuable for identifying and mapping key natural resources, such as water, minerals, and forests to inform decisions on how to manage these resources in a sustainable manner (cf. the emerging field of “landscape agronomy”; Rizzo, Marraccini, and Lardon 2022). Enhancing landscape mapping by dividing it into distinct units significantly enhances the representation of the complex interplay between biophysical and human factors across the land. Not only does this characterization better capture the heterogeneity and spatial variability of these factors, but it also offers substantial utility in model processing. By assigning specific attribute values to each delineated unit, biophysical models gain a heightened ability to simulate and analyze the intricate interactions and dynamics of various processes existing within the landscape (Opdam et al. 2018). For instance,

<figure><img src="image_11.png" /></figure>
<figure><img src="image_13.png" /></figure>

# Radiometric landscape: a new conceptual framework and operational approach for landscape characterisation and mapping

Louise Lemettais IDª, Samuel Alleaume IDᵇ,ᵉ, Sandra Luque IDᵇ,ᵉ, Anne-Élisabeth Laques IDª,ᵈ, Yonas Alimᶜ,ᵉ, Laurent Demagistri IDª and Agnès Bégué IDᶜ,ᵉ

ªUMR ESPACE-DEV, Univ Montpellier, IRD, Univ Antilles, Univ Guyane, Univ Réunion, Montpellier, France; ᵇINRAE, UMR TETIS, Montpellier, France; ᶜCIRAD, UMR TETIS, Montpellier, France; ᵈLMI PAYSAGES, IRD, CNRE, Université d'Antananarivo, Antananarivo, Madagascar; ᵉUMR TETIS, AgroParisTech CIRAD CNRS INRAE Univ Montpellier, Montpellier, France

## ABSTRACT

Landscape mapping has the potential to address some of the most pressing research issues of our time, including climate change, sustainable development, and human well-being. In this paper, we propose an original method that lays the foundations for landscape mapping and overcomes some of the major limitations of existing biophysical methods. Based on the assumption that the primary components of the landscape can be extracted directly from the radiometric information of satellite image time series, this paper presents a new approach to landscape characterization and mapping based solely on remote sensing data. The approach relies on a conceptual model, which links the description, characteristics, structure and functions of the landscape to a set of Remote Sensing-based Essential Landscape Variables (RS-ELVs). The RS-ELVs are then processed according to geographic object-based image analysis (GEOBIA) approach to produce a radiometric landscape map. The model and the remote sensing data processing chain are tested on a case study in central Madagascar (about 13 000 km²) composed of contrasting landscapes resulting from different climatic conditions and agricultural practices. The RS-ELVs are extracted from MODIS image time series for the temporal and spectral variables, and from MODIS and Sentinel-2 images for the texture variables. The parameterization of the segmentation and clustering algorithms is determined by statistical optimization. The final result is a radiometric landscape map in six classes. The landscape classes are then characterized using an independent set of remote sensing variables, a global land cover map and ground observations. The approach successfully identifies and delineates the gradient and major landscape types of the complex region of central Madagascar, confirming our initial hypothesis. The production of such radiometric landscape maps opens the way for integrated territorial development, including the planning and protection of the living environment and human well-being, and the implementation of sectoral policies.

## ARTICLE HISTORY

Received 31 July 2023
Accepted 31 January 2024

## KEYWORDS

Remote sensing; MODIS;
Sentinel-2; essential
variables; satellite image
time series; Madagascar

## 1. Introduction

Due to the multidisciplinary nature of landscape studies and their almost infinite variability, the simplification of landscapes into spatial units is a complex task, but nevertheless important for communication in management and research (e.g. Bunce et al. 1996; Hazeu et al. 2011; Turner and Gardner 2015). Landscape mapping is a cornerstone of the integrated territorial development, including the planning and protection of the living environment and human well-being (e.g. Reed, Ros-Tonen, and Sunderland 2020), and the implementation of sectoral policies such as the green and blue infrastructures (Bourget and le Dû-Blayo 2010; Sun, Wei, and Han 2022). Landscape characterization and mapping is also important in biogeography and ecology (e.g. Cullum et al. 2016) to guide the setting of conservation priorities, or to better understand and manage biodiversity (e.g. Préau et al.

2022). More generally, landscape mapping is very valuable for identifying and mapping key natural resources, such as water, minerals, and forests to inform decisions on how to manage these resources in a sustainable manner (cf. the emerging field of “landscape agronomy”; Rizzo, Marraccini, and Lardon 2022). Enhancing landscape mapping by dividing it into distinct units significantly enhances the representation of the complex interplay between biophysical and human factors across the land. Not only does this characterization better capture the heterogeneity and spatial variability of these factors, but it also offers substantial utility in model processing. By assigning specific attribute values to each delineated unit, biophysical models gain a heightened ability to simulate and analyze the intricate interactions and dynamics of various processes existing within the landscape (Opdam et al. 2018). For instance,

<figure><img src="image_14.png" /></figure>
<figure><img src="image_12.png" /></figure>

employing this approach greatly benefits simulations such as crop models (Leenhardt et al. 2010). Furthermore, the application of land stratification techniques plays a crucial role in refining satellite image classification accuracy by mitigating the inherent heterogeneity of class signatures, as exemplified in studies by Bellón et al. (2018) and Cano et al. (2017).

The review proposed by Simensen, Halvorsen, and Erikstad (2018) provides a good overview of the methods used for landscape characterization and mapping (LCM). The review proposes three groups of methodological approaches: 1. A holistic approach, generally based on visual landscape (from ground photographs or satellite images) and surveys of human perception (e.g. Tenerelli, Püffel, and Luque 2017; Tveit, Ode, and Fry 2006; Yazici 2018); 2. A holistic-biophysical methodology, based on an upscaling/bottom-up approach of a priori selection (based on expertise) of geoenvironmental and land-use-related characteristics of the landscape (e.g. Koç and Yılmaz 2020; Mücher et al. 2010; Silva et al. 2020); 3. Biophysical methods, in which the landscapes are characterized using multivariate statistical analyses or hierarchical clustering techniques of a large number of physical landscape attributes (e.g. Cabral et al. 2018; García-Llamas et al. 2016). The frequency of use of these three LCM methods tracks the development of advanced statistical analysis methods, geographical information systems (GIS), and the accessibility of open environmental databases (Simensen, Halvorsen, and Erikstad 2018). More recently, artificial intelligence has been added to the list of the methods, as shown in van Strien and Grêt-Regamey (2022) who generated a landscape typology for Switzerland by applying an unsupervised deep learning method to environmental-demographic data and satellite images. Hence, all these approaches have limitations. For example, holistic methods are non-reproducible, non-quantitative, and difficult to implement on a regional scale. Biophysical methods depend on the quality of the input data (such as the accuracy of the thematic maps), on the correspondence between the spatial resolution of the data and the scale at which the landscapes are intended to be characterized (García-Llamas et al. 2016), and, above all, on the subjectivity of the map producer in the initial selection of variables (Alcántara Manzanares and Muñoz Álvarez 2015) and the choice of statistical methods (Simensen, Halvorsen, and Erikstad 2018).

Earth observation has a great potential for landscape mapping, as it allows the detection of features and changes in land cover, land use and land conditions. The variety of images available makes it possible to work in a multiscalar mode (from local to continental scale) according to the chosen spatial resolution (from meters to kilometers). Most of the recent satellite systems acquire multispectral images at a frequency that allows to capture the seasonal

variations of the land conditions. Finally, satellite data are objective (pixel values are physical measures), repetitive and available over long periods of time. Since the early 2000s, the field of remote sensing has witnessed the emergence of various mapping methods. Supervised and unsupervised classifications of satellite data are commonplace for thematic mapping, including vegetation cover mapping (Ouchra, Belangour, and Erraissi 2022). These advantages make Earth observation the main source of data used to identify landscapes, as confirmed by the literature review of Simensen, Halvorsen, and Erikstad (2018) which showed that remote sensing-derived digital elevation model, land cover/land use and vegetation maps are the most frequently used variables (96%, 83% and 81% respectively) in the landscape characterization and mapping studies.

Nevertheless, according to Newton et al. (2009), applications of remote sensing in landscape ecology must progress beyond the simplistic approach of thematic mapping, such as the land cover/land use mapping, and the derivation of two-dimensional pattern metrics (Luque et al. 2018). In line with this recommendation, several applications (Alleaume et al. 2018; Ferreira et al. 2020; M. Lang et al. 2018) make use of continuous indices derived from remote sensing as proxies for essential biodiversity variables (EBV). Bisquert, Bégué, and Deshayes (2015) used spectral and textural variables derived from MODIS NDVI satellite image time series and an object-based approach to divide the French territory into 300 land units. In particular, they worked in the Pyrenean region, where the resulting land units proved to be consistent in terms of land cover (Bisquert et al. 2017). More recently, Bellón et al. (2017) used a similar approach in Tocantins State in Brazil, and obtained land units classified into different cropping systems to generate agricultural land units of the State. In another field, based on the hypothesis that the features visible in satellite images can be used to estimate the beauty of the landscape, Levering, Marcos, and Tuia (2021) produced a landscape scenicness map of the United Kingdom by applying a deep learning algorithm trained on a landscape beauty estimator using Sentinel-2 images.

Despite the many possibilities and methods available for landscape mapping, there is evidence that the potential of Earth observation data is still under-utilized, mainly due to a lack of awareness or technical knowledge. This paper aims to contribute to the awareness of the potential of Earth observation by providing an example of methods and tools that can be used to map landscapes and to show the potential of remote sensing applications. In line with the aforementioned landscape approaches, we propose here an operational innovative approach to landscape characterization and
<figure><img src="image_15.png" /></figure>

employing this approach greatly benefits simulations such as crop models (Leenhardt et al. 2010). Furthermore, the application of land stratification techniques plays a crucial role in refining satellite image classification accuracy by mitigating the inherent heterogeneity of class signatures, as exemplified in studies by Bellón et al. (2018) and Cano et al. (2017).

The review proposed by Simensen, Halvorsen, and Erikstad (2018) provides a good overview of the methods used for landscape characterization and mapping (LCM). The review proposes three groups of methodological approaches: 1. A holistic approach, generally based on visual landscape (from ground photographs or satellite images) and surveys of human perception (e.g. Tenerelli, Püffel, and Luque 2017; Tveit, Ode, and Fry 2006; Yazici 2018); 2. A holistic-biophysical methodology, based on an upscaling/bottom-up approach of a priori selection (based on expertise) of geoecological and land-use-related characteristics of the landscape (e.g. Koç and Yılmaz 2020; Mücher et al. 2010; Silva et al. 2020); 3. Biophysical methods, in which the landscapes are characterized using multivariate statistical analyses or hierarchical clustering techniques of a large number of physical landscape attributes (e.g. Cabral et al. 2018; García-Llamas et al. 2016). The frequency of use of these three LCM methods tracks the development of advanced statistical analysis methods, geographical information systems (GIS), and the accessibility of open environmental databases (Simensen, Halvorsen, and Erikstad 2018). More recently, artificial intelligence has been added to the list of the methods, as shown in van Strien and Grêt-Regamey (2022) who generated a landscape typology for Switzerland by applying an unsupervised deep learning method to environmental-demographic data and satellite images. Hence, all these approaches have limitations. For example, holistic methods are non-reproducible, non-quantitative, and difficult to implement on a regional scale. Biophysical methods depend on the quality of the input data (such as the accuracy of the thematic maps), on the correspondence between the spatial resolution of the data and the scale at which the landscapes are intended to be characterized (García-Llamas et al. 2016), and, above all, on the subjectivity of the map producer in the initial selection of variables (Alcántara Manzanares and Muñoz Álvarez 2015) and the choice of statistical methods (Simensen, Halvorsen, and Erikstad 2018).

Earth observation has a great potential for landscape mapping, as it allows the detection of features and changes in land cover, land use and land conditions. The variety of images available makes it possible to work in a multiscalar mode (from local to continental scale) according to the chosen spatial resolution (from meters to kilometers). Most of the recent satellite systems acquire multispectral images at a frequency that allows to capture the seasonal

variations of the land conditions. Finally, satellite data are objective (pixel values are physical measures), repetitive and available over long periods of time. Since the early 2000s, the field of remote sensing has witnessed the emergence of various mapping methods. Supervised and unsupervised classifications of satellite data are commonplace for thematic mapping, including vegetation cover mapping (Ouchra, Belangour, and Erraissi 2022). These advantages make Earth observation the main source of data used to identify landscapes, as confirmed by the literature review of Simensen, Halvorsen, and Erikstad (2018) which showed that remote sensing-derived digital elevation model, land cover/land use and vegetation maps are the most frequently used variables (96%, 83% and 81% respectively) in the landscape characterization and mapping studies.

Nevertheless, according to Newton et al. (2009), applications of remote sensing in landscape ecology must progress beyond the simplistic approach of thematic mapping, such as the land cover/land use mapping, and the derivation of two-dimensional pattern metrics (Luque et al. 2018). In line with this recommendation, several applications (Alleaume et al. 2018; Ferreira et al. 2020; M. Lang et al. 2018) make use of continuous indices derived from remote sensing as proxies for essential biodiversity variables (EBV). Bisquert, Bégué, and Deshayes (2015) used spectral and textural variables derived from MODIS NDVI satellite image time series and an object-based approach to divide the French territory into 300 land units. In particular, they worked in the Pyrenean region, where the resulting land units proved to be consistent in terms of land cover (Bisquert et al. 2017). More recently, Bellón et al. (2017) used a similar approach in Tocantins State in Brazil, and obtained land units classified into different cropping systems to generate agricultural land units of the State. In another field, based on the hypothesis that the features visible in satellite images can be used to estimate the beauty of the landscape, Levering, Marcos, and Tuia (2021) produced a landscape scenicness map of the United Kingdom by applying a deep learning algorithm trained on a landscape beauty estimator using Sentinel-2 images.

Despite the many possibilities and methods available for landscape mapping, there is evidence that the potential of Earth observation data is still under-utilized, mainly due to a lack of awareness or technical knowledge. This paper aims to contribute to the awareness of the potential of Earth observation by providing an example of methods and tools that can be used to map landscapes and to show the potential of remote sensing applications. In line with the aforementioned landscape approaches, we propose here an operational innovative approach to landscape characterization and
<figure><img src="image_13.png" /></figure>

mapping based solely on remote sensing data that are available worldwide and free of charge. The aim is to provide a transposable approach that could be reproducible at different scales and in different landscape contexts, in order to overcome the main limitations associated with the biophysical approaches, i.e. the quality of the thematic maps and the dependence on the observer. Our approach is based on the assumption that the primary components of the landscape can be extracted directly from the radiometric information of satellite images, especially from time series, and processed according to geographic object-based image analysis (GEOBIA) concepts. We first present the underlying conceptual model, and then we apply it by developing a remote sensing processing chain. The approach is tested on a case study in the highlands of Madagascar with the aim of producing a landscape map of a significant area (about 13,000 km²). The case serves as a demonstrator and can be replicable in other regions.

## 2. The conceptual approach and workflow

This section presents the conceptual basis of radiometric landscapes and the resulting workflow for processing the remote sensing data. We first present the landscape terminology used, and then link the descriptive variables of a landscape with variables extracted from satellite imagery to produce a list of Remote Sensing-based Essential Landscape Variables (RS-ELVs). These RS-ELVs are then processed using a GEOBIA approach to produce a radiometric landscape map. RS-ELVs computation, segmentation and clustering are part of the workflow presented below, along with the characterization and labeling of the resulting landscape classes.

### 2.1. The terminology

In order to clarify different issues and concepts, we first devote a sub-section to the terminology used. In this paper landscapes are considered to be recognizable, though often heterogeneous, parts of the Earth's surface that display a characteristic arrangement of elements resulting from the long- and medium-term interactions between natural abiotic, biotic and anthropogenic processes (e.g. Antrop 2005; Hazeu et al. 2011; Wu 2013). Landscapes correspond to an area of land that can be perceived by the eye. These landscapes are called physiognomic landscapes, and are composed of nested spatial components (landscape element, feature, unit and type; Figure 1).

A landscape element is defined as “a natural or human-induced object, category or characteristic, including ecosystem type, which is observable at landscape scale” (Erikstad, Andre Uttakleiv, and Halvorsen 2015). A landscape element (field, tree, house, etc.) is usually identifiable in aerial photographs; it is often 1 m to 10 m wide (Forman and Godron 1986; Turner and Gardner 2015), which corresponds to a pixel or a group of pixels.

The landscape feature corresponds to a system formed by elements of the landscape and the material and immaterial relationships that link them (cultivated lowlands, plots cultivated in terraces, villages, etc.). The landscape features represent, in a way, the pieces of the puzzle that need to be put together to form landscape units. They can be obtained by grouping adjacent pixels with similar characteristics (image segmentation), and are therefore referred to as segments or objects.

The landscape unit is characterized by a set of similar landscape features that are spatially adjacent. The unit is distinguished from the neighboring units by a difference in the presence, organization or shape

<figure><img src="image_14.png" /><figcaption>**Figure 1.** Landscape equivalent terminology in geography and remote sensing at different scales, and graphical representation. From left to right: landscape elements (fields, village, ...) visible on a UAV photograph (credit: S. Alleaume); a landscape feature, composed of landscape elements; a landscape unit made of adjacent similar landscape features; a landscape type made of similar landscape units.</figcaption></figure>
<figure><img src="image_16.png" /></figure>

mapping based solely on remote sensing data that are available worldwide and free of charge. The aim is to provide a transposable approach that could be reproducible at different scales and in different landscape contexts, in order to overcome the main limitations associated with the biophysical approaches, i.e. the quality of the thematic maps and the dependence on the observer. Our approach is based on the assumption that the primary components of the landscape can be extracted directly from the radiometric information of satellite images, especially from time series, and processed according to geographic object-based image analysis (GEOBIA) concepts. We first present the underlying conceptual model, and then we apply it by developing a remote sensing processing chain. The approach is tested on a case study in the highlands of Madagascar with the aim of producing a landscape map of a significant area (about 13,000 km²). The case serves as a demonstrator and can be replicable in other regions.

## 2. The conceptual approach and workflow

This section presents the conceptual basis of radiometric landscapes and the resulting workflow for processing the remote sensing data. We first present the landscape terminology used, and then link the descriptive variables of a landscape with variables extracted from satellite imagery to produce a list of Remote Sensing-based Essential Landscape Variables (RS-ELVs). These RS-ELVs are then processed using a GEOBIA approach to produce a radiometric landscape map. RS-ELVs computation, segmentation and clustering are part of the workflow presented below, along with the characterization and labeling of the resulting landscape classes.

### 2.1. The terminology

In order to clarify different issues and concepts, we first devote a sub-section to the terminology used. In this paper landscapes are considered to be recognizable, though often heterogeneous, parts of the Earth's surface that display a characteristic arrangement of elements resulting from the long- and medium-term interactions between natural abiotic, biotic and anthropogenic processes (e.g. Antrop 2005; Hazeu et al. 2011; Wu 2013). Landscapes correspond to an area of land that can be perceived by the eye. These landscapes are called physiognomic landscapes, and are composed of nested spatial components (landscape element, feature, unit and type; Figure 1).

A landscape element is defined as “a natural or human-induced object, category or characteristic, including ecosystem type, which is observable at landscape scale” (Erikstad, Andre Uttakleiv, and Halvorsen 2015). A landscape element (field, tree, house, etc.) is usually identifiable in aerial photographs; it is often 1 m to 10 m wide (Forman and Godron 1986; Turner and Gardner 2015), which corresponds to a pixel or a group of pixels.

The landscape feature corresponds to a system formed by elements of the landscape and the material and immaterial relationships that link them (cultivated lowlands, plots cultivated in terraces, villages, etc.). The landscape features represent, in a way, the pieces of the puzzle that need to be put together to form landscape units. They can be obtained by grouping adjacent pixels with similar characteristics (image segmentation), and are therefore referred to as segments or objects.

The landscape unit is characterized by a set of similar landscape features that are spatially adjacent. The unit is distinguished from the neighboring units by a difference in the presence, organization or shape

<figure><img src="image_17.png" /><figcaption>**Figure 1.** Landscape equivalent terminology in geography and remote sensing at different scales, and graphical representation. From left to right: landscape elements (fields, village, ...) visible on a UAV photograph (credit: S. Alleaume); a landscape feature, composed of landscape elements; a landscape unit made of adjacent similar landscape features; a landscape type made of similar landscape units.</figcaption></figure>
<figure><img src="image_15.png" /></figure>

of its characteristics (Raymond et al. 2015). The landscape unit can vary in size down to less than a few kilometers in diameter (Simensen, Halvorsen, and Erikstad 2018). Here, it is obtained by unsupervised classification (clustering) of the landscape segments, based on a set of variables, and corresponds to one spatial unit (one polygon).

The landscape type is representative of portions of homogeneous space that are coherent at the biophysical and socioeconomic levels (nature, arrangement and frequency of the constituent elements) at a regional scale. In general, the term “landscape type” is used when the landscape resembles a model, an archetype to which it is necessary to refer in order to produce landscape maps. The landscape type is obtained by clustering, and composed of polygon(s) (landscape unit(s)) that are not necessarily spatially contiguous. Once labeled, the cluster is referred to as a class.

## 2.2. The radiometric landscape concept

A physiognomic landscape can be described by a set of descriptors – geometry, colors, grain size, complexity – in three dimensions, and without reference to terrestrial processes or to the visual perception of the observer. Based on the model of Antrop and Van Eetvelde (2017), Karasov, Külvik, and Burdun (2021) classified physiognomic landscape descriptors into five

categories: point, line/polygon, surface, color and texture. In this section, we show that each of these landscape descriptors can be assessed from satellite imagery (Table 1; Figure 2) using a set of variables hereafter referred to as remote sensing-essential landscape variables (RS-ELVs).

The spatial information (shape and size of the landscape elements, distance between the elements, etc.) is intrinsic to satellite images. Mathematical morphology, originally developed for binary images, can be used on greyscale images to suppress noise, enhance the image and detect edges. This type of information is not retained in our approach, because the morphological indices are generally used for urban landscape analysis using very high spatial resolution (VHSR) images. At the scale at which we intend to produce landscape maps (territorial, regional), the spatial information is partly contained in the texture of the image (see textural information below).

The three-dimensional information (landform) is an important component of landscape, and is very often included in the biophysical approaches. The Digital Elevation Models (DEMs) are more and more produced from remote sensing data as stereopairs of optical satellite images, radar interferometry or lidar data. In our approach, we assume that the spectral response of the vegetation (vegetation type, condition and seasonality) is strongly dependent on

<table><caption>**Table 1.** Correspondence between the landscape descriptors (see text), the Remote Sensing-Essential Landscape Variables (RS-ELVs), and the landscape variables; the RS-ELVs used in this study are highlighted in bold.</caption><thead><tr><th>Landscape descriptors</th><th>Information type</th><th>Remote Sensing-Essential Landscape Variables (RS-ELVs)</th><th>Landscape traits, structure and functions</th></tr></thead><tbody><tr><td>Point, Line/Polygon</td><td>Spatial</td><td>Mathematical morphology indices</td><td>Shape and size of the landscape elements, distance between the elements, etc.</td></tr><tr><td>Surface</td><td>Three-dimensional</td><td>Altitude, Slope, Aspect</td><td>Landform, exposure, etc.</td></tr><tr><td>Color</td><td>Spectral</td><td>Spectral bands, <strong>Spectral indices</strong></td><td>Land cover, biomass, vegetation conditions, soil type, etc.</td></tr><tr><td>Texture</td><td>Textural</td><td><strong>Textural indices</strong></td><td>Landscape structure and composition, landform, hydrographic network, etc.</td></tr><tr><td></td><td>Temporal</td><td>Temporal indices, <strong>Principal Components</strong></td><td>Phenology, climate, cropping calendar and practices, etc.</td></tr></tbody></table>

<figure><img src="image_16.png" /><figcaption>**Figure 2.** Example of seasonal changes in colours and texture observed in a pair of Sentinel-2 sub-images (true-colour composites) acquired on 11 April 2019 at the end of the rainy season (on the left), and 28 October 2019 at the end of the dry season (on the right).</figcaption></figure>
<figure><img src="image_18.png" /></figure>

of its characteristics (Raymond et al. 2015). The landscape unit can vary in size down to less than a few kilometers in diameter (Simensen, Halvorsen, and Erikstad 2018). Here, it is obtained by unsupervised classification (clustering) of the landscape segments, based on a set of variables, and corresponds to one spatial unit (one polygon).

The landscape type is representative of portions of homogeneous space that are coherent at the biophysical and socioeconomic levels (nature, arrangement and frequency of the constituent elements) at a regional scale. In general, the term “landscape type” is used when the landscape resembles a model, an archetype to which it is necessary to refer in order to produce landscape maps. The landscape type is obtained by clustering, and composed of polygon(s) (landscape unit(s)) that are not necessarily spatially contiguous. Once labeled, the cluster is referred to as a class.

## 2.2. The radiometric landscape concept

A physiognomic landscape can be described by a set of descriptors – geometry, colors, grain size, complexity – in three dimensions, and without reference to terrestrial processes or to the visual perception of the observer. Based on the model of Antrop and Van Eetvelde (2017), Karasov, Külvik, and Burdun (2021) classified physiognomic landscape descriptors into five

categories: point, line/polygon, surface, color and texture. In this section, we show that each of these landscape descriptors can be assessed from satellite imagery (Table 1; Figure 2) using a set of variables hereafter referred to as remote sensing-essential landscape variables (RS-ELVs).

The spatial information (shape and size of the landscape elements, distance between the elements, etc.) is intrinsic to satellite images. Mathematical morphology, originally developed for binary images, can be used on greyscale images to suppress noise, enhance the image and detect edges. This type of information is not retained in our approach, because the morphological indices are generally used for urban landscape analysis using very high spatial resolution (VHSR) images. At the scale at which we intend to produce landscape maps (territorial, regional), the spatial information is partly contained in the texture of the image (see textural information below).

The three-dimensional information (landform) is an important component of landscape, and is very often included in the biophysical approaches. The Digital Elevation Models (DEMs) are more and more produced from remote sensing data as stereopairs of optical satellite images, radar interferometry or lidar data. In our approach, we assume that the spectral response of the vegetation (vegetation type, condition and seasonality) is strongly dependent on

<table><caption>**Table 1.** Correspondence between the landscape descriptors (see text), the Remote Sensing-Essential Landscape Variables (RS-ELVs), and the landscape variables; the RS-ELVs used in this study are highlighted in bold.</caption><thead><tr><th>Landscape descriptors</th><th>Information type</th><th>Remote Sensing-Essential Landscape Variables (RS-ELVs)</th><th>Landscape traits, structure and functions</th></tr></thead><tbody><tr><td>Point, Line/Polygon</td><td>Spatial</td><td>Mathematical morphology indices</td><td>Shape and size of the landscape elements, distance between the elements, etc.</td></tr><tr><td>Surface</td><td>Three-dimensional</td><td>Altitude, Slope, Aspect</td><td>Landform, exposure, etc.</td></tr><tr><td>Color</td><td>Spectral</td><td>Spectral bands, <strong>Spectral indices</strong></td><td>Land cover, biomass, vegetation conditions, soil type, etc.</td></tr><tr><td>Texture</td><td>Textural</td><td><strong>Textural indices</strong></td><td>Landscape structure and composition, landform, hydrographic network, etc.</td></tr><tr><td></td><td>Temporal</td><td>Temporal indices, <strong>Principal Components</strong></td><td>Phenology, climate, cropping calendar and practices, etc.</td></tr></tbody></table>

<figure><img src="image_19.png" /><figcaption>**Figure 2.** Example of seasonal changes in colours and texture observed in a pair of Sentinel-2 sub-images (true-colour composites) acquired on 11 April 2019 at the end of the rainy season (on the left), and 28 October 2019 at the end of the dry season (on the right).</figcaption></figure>
<figure><img src="image_17.png" /></figure>

its three-dimensional environment, and thus that the landform does not need to be explicitly taken into account as it is already included in the image radiometry (see spectral and temporal information below).

The spectral information (land cover type, biophysical state of the landscape elements, etc.) contained in multispectral imagery provides access to different landscape properties: the colour of the components (visible bands), the vegetation density and structure (visible and near-infrared bands), the vegetation biomass (cumulative spectral indices), and the inland water and vegetation water status (short-wave infrared band).

The textural information (landscape structure and composition, landform, road and hydrographic networks, etc.) depends on the size and spatial arrangement of the landscape elements, and on the spatial resolution of the image. The most common texture metrics are the Haralick metrics based on the Grey-Level Co-occurrence Matrix (Haralick, Shanmugam, and Dinstein 1973), but other approaches exist such as wavelet texture analysis or Fourier-based textural ordination (Couteron, Barbier, and Gautier 2006).

The temporal information (vegetation phenology, climate conditions, cropping practices) contained in dense satellite image time series (SITS) depends on the satellite system, the local atmospheric conditions and the quality of the preprocessing. As far as vegetation is concerned, phenological metrics can be calculated directly from SITS (start-of-season, length of season, maximum VI, integrated VI, etc.), but generally speaking, the Principal Component Analysis (PCA) applied to SITS is the simplest way to extract and synthetize the seasonal spectral variation of the surface (Bellón et al. 2017; Bisquert, Bégué, and Deshayes 2015).

As our approach to landscape characterization and mapping is based solely on remote sensing ELVs, it is hereafter referred to as Radiometric Landscape.

## 2.3. Overview of the workflow

To produce the radiometric landscape maps, we used multi-sensor imagery (MODIS and Sentinel-2) to take advantage of both the temporal information (daily) provided by MODIS satellites, and the spatial information (10 m) provided by Sentinel-2 satellites. These image datasets were chosen because they are available worldwide, free of charge, and they complement each other. The core of the processing chain is GEOBIA (GEographic Object-Based Image Analysis) unsupervised classification that responds to the need of image partitioning into objects according to a conceptual model (Blaschke 2010).

The workflow to produce the radiometric landscape map is summarized in Figure 3. It consists of three steps for mapping, and one step for characterization:

(1) Selection of satellite data and calculation of the essential landscape variables or RS-ELVs (seasonal, spectral and textural indices);

(2) Segmentation of the satellite data set, based on the RS-ELVs, to produce landscape features;

(3) Clustering of the segments into radiometric landscape types.

Each of these steps was evaluated using statistical indices in order to select the optimal dataset (sensors, bands, spectral indices, acquisition dates) for step 1 (ELVs calculation), and the optimal input variables, algorithm and parametrization for step 2 (segmentation) and step 3 (clustering).

(1) Characterization of the clusters using remote sensing variables and additional information (land cover, local expertise).

GEOBIA is based on the hypothesis that the partitioning of an image into objects is related to the way

<figure><img src="image_18.png" /><figcaption>Figure 3. Overview of the 4-step workflow to produce and characterize radiometric landscape maps.</figcaption></figure>
<figure><img src="image_20.png" /></figure>

its three-dimensional environment, and thus that the landform does not need to be explicitly taken into account as it is already included in the image radiometry (see spectral and temporal information below).

The spectral information (land cover type, biophysical state of the landscape elements, etc.) contained in multispectral imagery provides access to different landscape properties: the colour of the components (visible bands), the vegetation density and structure (visible and near-infrared bands), the vegetation biomass (cumulative spectral indices), and the inland water and vegetation water status (short-wave infrared band).

The textural information (landscape structure and composition, landform, road and hydrographic networks, etc.) depends on the size and spatial arrangement of the landscape elements, and on the spatial resolution of the image. The most common texture metrics are the Haralick metrics based on the Grey-Level Co-occurrence Matrix (Haralick, Shanmugam, and Dinstein 1973), but other approaches exist such as wavelet texture analysis or Fourier-based textural ordination (Couteron, Barbier, and Gautier 2006).

The temporal information (vegetation phenology, climate conditions, cropping practices) contained in dense satellite image time series (SITS) depends on the satellite system, the local atmospheric conditions and the quality of the preprocessing. As far as vegetation is concerned, phenological metrics can be calculated directly from SITS (start-of-season, length of season, maximum VI, integrated VI, etc.), but generally speaking, the Principal Component Analysis (PCA) applied to SITS is the simplest way to extract and synthetize the seasonal spectral variation of the surface (Bellón et al. 2017; Bisquert, Bégué, and Deshayes 2015).

As our approach to landscape characterization and mapping is based solely on remote sensing ELVs, it is hereafter referred to as Radiometric Landscape.

## 2.3. Overview of the workflow

To produce the radiometric landscape maps, we used multi-sensor imagery (MODIS and Sentinel-2) to take advantage of both the temporal information (daily) provided by MODIS satellites, and the spatial information (10 m) provided by Sentinel-2 satellites. These image datasets were chosen because they are available worldwide, free of charge, and they complement each other. The core of the processing chain is GEOBIA (GEographic Object-Based Image Analysis) unsupervised classification that responds to the need of image partitioning into objects according to a conceptual model (Blaschke 2010).

The workflow to produce the radiometric landscape map is summarized in Figure 3. It consists of three steps for mapping, and one step for characterization:

(1) Selection of satellite data and calculation of the essential landscape variables or RS-ELVs (seasonal, spectral and textural indices);

(2) Segmentation of the satellite data set, based on the RS-ELVs, to produce landscape features;

(3) Clustering of the segments into radiometric landscape types.

Each of these steps was evaluated using statistical indices in order to select the optimal dataset (sensors, bands, spectral indices, acquisition dates) for step 1 (ELVs calculation), and the optimal input variables, algorithm and parametrization for step 2 (segmentation) and step 3 (clustering).

(1) Characterization of the clusters using remote sensing variables and additional information (land cover, local expertise).

GEOBIA is based on the hypothesis that the partitioning of an image into objects is related to the way

<figure><img src="image_21.png" /><figcaption>Figure 3. Overview of the 4-step workflow to produce and characterize radiometric landscape maps.</figcaption></figure>
<figure><img src="image_19.png" /></figure>

people conceptually organize the landscape in order to understand it (Hay and Castilla 2006). GEOBIA approach has been widely used since the 2000s driven by the need to invent new data processing solutions to deal with the arrival of Earth observation data acquired at very high spatial resolution (e.g. Kalinicheva, Sublime, and Trocan 2020; Ma et al. 2017). Since then, several studies have shown that GEOBIA can be applied to other types of remote sensing data (e.g. Bégué et al. 2015; Cano et al. 2017), provided that they are H-resolution images. The concept of H-resolution (and its counterpart the L-resolution) refers to the relative size of the object compared to the image spatial resolution (Blaschke et al. 2014; Strahler, Woodcock, and Smith 1986). The H-resolution model applies when scene objects are much larger than the spatial resolution of the image (several pixels can represent a single object), typically in landscape mapping (see Figure 1).

Many authors agree that the segmentation process is the most fundamental step in GEOBIA as the segmentation quality significantly affects the classification results (e.g. Blaschke 2010; Clinton et al. 2010; Kavzoglu and Tonbul 2018). The choice of the variables used (the image spectral information, but also the textural, spatial and contextual information can be used in GEOBIA), and the segmentation parameters are particularly sensitive. Therefore, in this work, careful attention has been paid to the segmentation stage, testing different datasets (according to our conceptual model) and parameters (with statistical criteria) in order to

make the segmentation step as automatic and objective as possible.

# 3. Material and methods

## 3.1. Study area

The study was carried out in the Haute Matsiatra region of central Madagascar (Figure 4). The study area covers 13,000 km² and is composed, from west to east, of humid temperate highlands (altitude between 1000 and 1400 m), an ecological corridor with the Ranomafana Park (between 900 and 1200 m), and a humid tropical zone (between 100 and 800 m). The area is strongly influenced by the currents of the Indian Ocean creating contrasting environmental conditions between the eastern facade, mostly cloudy with 2000–2500 mm rainfall per year, and the vast western plateau, sunny with 1200–1600 mm per year (Ollivier et al. 2023).

The landscapes of the western highlands are characterized by red ferralitic soils and consist of areas of cultivated land, large savannahs, as well as reforested areas (eucalyptus, pins, acacias) and some patches of mature forest as one approaches the corridor. The corridor consists almost entirely of protected primary (mature) forests. East of the corridor, the hillside is covered by a dense formation dominated by bamboo and ravenalas (a plant formation locally called Savoka). A few cultivated plots and rice fields occupy the lower slopes and bottomlands.

<figure><img src="image_20.png" /><figcaption>**Figure 4.** (a) Location map of the study area in Madagascar (red square); (b) Digital elevation Model of the study area (source: SRTM); (c) Study area limits (red square) over a true-color composition of Landsat/SPOT image (source: google earth), with main roads (white lines; source World Bank).</figcaption></figure>
<figure><img src="image_22.png" /></figure>

people conceptually organize the landscape in order to understand it (Hay and Castilla 2006). GEOBIA approach has been widely used since the 2000s driven by the need to invent new data processing solutions to deal with the arrival of Earth observation data acquired at very high spatial resolution (e.g. Kalinicheva, Sublime, and Trocan 2020; Ma et al. 2017). Since then, several studies have shown that GEOBIA can be applied to other types of remote sensing data (e.g. Bégué et al. 2015; Cano et al. 2017), provided that they are H-resolution images. The concept of H-resolution (and its counterpart the L-resolution) refers to the relative size of the object compared to the image spatial resolution (Blaschke et al. 2014; Strahler, Woodcock, and Smith 1986). The H-resolution model applies when scene objects are much larger than the spatial resolution of the image (several pixels can represent a single object), typically in landscape mapping (see Figure 1).

Many authors agree that the segmentation process is the most fundamental step in GEOBIA as the segmentation quality significantly affects the classification results (e.g. Blaschke 2010; Clinton et al. 2010; Kavzoglu and Tonbul 2018). The choice of the variables used (the image spectral information, but also the textural, spatial and contextual information can be used in GEOBIA), and the segmentation parameters are particularly sensitive. Therefore, in this work, careful attention has been paid to the segmentation stage, testing different datasets (according to our conceptual model) and parameters (with statistical criteria) in order to

make the segmentation step as automatic and objective as possible.

# 3. Material and methods

## 3.1. Study area

The study was carried out in the Haute Matsiatra region of central Madagascar (Figure 4). The study area covers 13,000 km² and is composed, from west to east, of humid temperate highlands (altitude between 1000 and 1400 m), an ecological corridor with the Ranomafana Park (between 900 and 1200 m), and a humid tropical zone (between 100 and 800 m). The area is strongly influenced by the currents of the Indian Ocean creating contrasting environmental conditions between the eastern facade, mostly cloudy with 2000–2500 mm rainfall per year, and the vast western plateau, sunny with 1200–1600 mm per year (Ollivier et al. 2023).

The landscapes of the western highlands are characterized by red ferralitic soils and consist of areas of cultivated land, large savannahs, as well as reforested areas (eucalyptus, pins, acacias) and some patches of mature forest as one approaches the corridor. The corridor consists almost entirely of protected primary (mature) forests. East of the corridor, the hillside is covered by a dense formation dominated by bamboo and ravenalas (a plant formation locally called Savoka). A few cultivated plots and rice fields occupy the lower slopes and bottomlands.

<figure><img src="image_23.png" /><figcaption>**Figure 4.** (a) Location map of the study area in Madagascar (red square); (b) Digital elevation Model of the study area (source: SRTM); (c) Study area limits (red square) over a true-color composition of Landsat/SPOT image (source: google earth), with main roads (white lines; source World Bank).</figcaption></figure>
<figure><img src="image_21.png" /></figure>

## 3.2. Ground observations

A field campaign was carried out in the study area between 29 April and 3 May 2022 to identify and label the main Malagasy landscapes. The approach was developed by the LMI PAYSAGES working group in Madagascar. The approach consists of four steps: (i) analysis of Google Earth images and road maps combined with knowledge of the terrain to propose an initial landscape zoning and organize the field trip; (ii) field trip to the core of the identified zones and data collection for each stop using a standardized observation grid; (iii) processing of the collected data to produce a landscape typology; (iv) preparation of a document describing each identified landscape type.

The approach was complemented by aerial drone photography. The protocol consists of flying over the landscapes at a height of about 100 meters and at a distance of more than 150 meters from the road, in order to reach the core of the landscape unit. Panoramic 360° shots are then taken using a Anafi drone, Parrot. Back in the laboratory, the shots were assembled and transformed into a panorama using the Hugin software (freeware, GNU General public license, http://hugin.sourceforge.net/). A total of 22 panoramas were produced, and deposited on a dataverse platform (Alleaume 2022).

## 3.3. Satellites data and pre-processing

### 3.3.1. MODIS NDVI time series

The MOD13Q1 V6 product (Terra Vegetation Indices 16-Day Global 250 m) provides vegetation indices computed from atmospherically corrected bidirectional surface reflectances masked for water, clouds, heavy aerosols, and cloud shadows. The product is available on NASA's AppEEARS platform (https://appeears.earthdatacloud.nasa.gov). For regional landscape studies, the Normalised Difference Vegetation Index (NDVI) has proven to perform well (e.g. Soudani et al. 2008; Vintrou et al. 2012), and thus a total of 115 NDVI images were downloaded (23 composite images per year) for the 2016–2020 period.

Each annual NDVI series was smoothed for noise, mainly caused by cloud contamination and atmospheric variability, using the Savitzky-Golay filter (J. Chen et al. 2004) parametrized with a half-window of 2 and a polynomial of degree 2.

### 3.3.2. Sentinel-2

The Sentinel-2 Level 2A products (tile T38KQB) were downloaded from the Theia platform (https://cataloque.theia-land.fr). This product consists of surface reflectance data corrected for atmospheric effects, and a mask of clouds and their shadows obtained

with the MAJA processing chain (Lonjou et al. 2016). To account for the seasonal dynamics of the landscape (Kalinicheva, Sublime, and Trocan 2020), we selected “cloud-free” images acquired on two dates: 11 April 2019, characteristic of the end of the rainy season with significant vegetation development, and 28 October 2019 characteristic of the dry season. To correct for residual cloud (10% and 3% of the study area for the October and April images respectively), a temporal gap-filling technique was applied to replace the value of the masked cloud pixels with the interpolated value of pixels from the images acquired on the closest dates (Inglada 2016).

Two vegetation indices were then calculated from the gap-filled images: the NDVI at 10 m resolution (from B4 and B8 bands; Rouse et al. 1973), and the NDWI (Normalised Difference Water Index) at 20 m resolution (from B8a and B12 bands; Gao 1996). These two spectral indices are complementary, with NDVI being particularly sensitive to the spatial and temporal variations in vegetation types and conditions, while NDWI is particularly sensitive to surface water conditions (e.g. Dupuy et al. 2020) for irrigated field mapping.

## 3.4. Calculation of the remote sensing ELVs

### 3.4.1. Spectral and temporal ELVs

Both spectral and temporal ELVs are based on an annual MODIS NDVI time series, hereafter referred to as the MODIS reference year. The MODIS reference year is designed to represent the mean seasonal variation of the surface, while attenuating any local or sudden surface changes due to particular climatic conditions. The reference year consists of 23 NDVI images (16-day composites), each calculated as the mean of the corresponding composites over a 5-year period (from November 2016 to October 2020, to include 5 full climatic years). The time period of 5 years was considered to be a good compromise, long enough to reduce the effect of annual climate variability and short enough to avoid large changes in land cover and land use.

The spectral ELV is calculated as the mean annual NDVI calculated over the MODIS reference year. The temporal ELVs are obtained by applying a Principal Component Analysis (PCA) transformation to the 23 NDVI values of the MODIS reference year. PCA is particularly suitable for capturing seasonal variations and identifying their relationship with different factors, whether climatic or related to human activities (Deng et al. 2008). In this study, PCA was performed using the FactoMineR package (Lê, Josse, and Husson 2008) of the R software and the rasterPCA tool of the RStoolbox package (Leutner et al. 2022).
<figure><img src="image_24.png" /></figure>

## 3.2. Ground observations

A field campaign was carried out in the study area between 29 April and 3 May 2022 to identify and label the main Malagasy landscapes. The approach was developed by the LMI PAYSAGES working group in Madagascar. The approach consists of four steps: (i) analysis of Google Earth images and road maps combined with knowledge of the terrain to propose an initial landscape zoning and organize the field trip; (ii) field trip to the core of the identified zones and data collection for each stop using a standardized observation grid; (iii) processing of the collected data to produce a landscape typology; (iv) preparation of a document describing each identified landscape type.

The approach was complemented by aerial drone photography. The protocol consists of flying over the landscapes at a height of about 100 meters and at a distance of more than 150 meters from the road, in order to reach the core of the landscape unit. Panoramic 360° shots are then taken using a Anafi drone, Parrot. Back in the laboratory, the shots were assembled and transformed into a panorama using the Hugin software (freeware, GNU General public license, http://hugin.sourceforge.net/). A total of 22 panoramas were produced, and deposited on a dataverse platform (Alleaume 2022).

## 3.3. Satellites data and pre-processing

### 3.3.1. MODIS NDVI time series

The MOD13Q1 V6 product (Terra Vegetation Indices 16-Day Global 250 m) provides vegetation indices computed from atmospherically corrected bidirectional surface reflectances masked for water, clouds, heavy aerosols, and cloud shadows. The product is available on NASA's AppEEARS platform (https://appeears.earthdatacloud.nasa.gov). For regional landscape studies, the Normalised Difference Vegetation Index (NDVI) has proven to perform well (e.g. Soudani et al. 2008; Vintrou et al. 2012), and thus a total of 115 NDVI images were downloaded (23 composite images per year) for the 2016–2020 period.

Each annual NDVI series was smoothed for noise, mainly caused by cloud contamination and atmospheric variability, using the Savitzky-Golay filter (J. Chen et al. 2004) parametrized with a half-window of 2 and a polynomial of degree 2.

### 3.3.2. Sentinel-2

The Sentinel-2 Level 2A products (tile T38KQB) were downloaded from the Theia platform (https://catalogue.theia-land.fr). This product consists of surface reflectance data corrected for atmospheric effects, and a mask of clouds and their shadows obtained

with the MAJA processing chain (Lonjou et al. 2016). To account for the seasonal dynamics of the landscape (Kalinicheva, Sublime, and Trocan 2020), we selected “cloud-free” images acquired on two dates: 11 April 2019, characteristic of the end of the rainy season with significant vegetation development, and 28 October 2019 characteristic of the dry season. To correct for residual cloud (10% and 3% of the study area for the October and April images respectively), a temporal gap-filling technique was applied to replace the value of the masked cloud pixels with the interpolated value of pixels from the images acquired on the closest dates (Inglada 2016).

Two vegetation indices were then calculated from the gap-filled images: the NDVI at 10 m resolution (from B4 and B8 bands; Rouse et al. 1973), and the NDWI (Normalised Difference Water Index) at 20 m resolution (from B8a and B12 bands; Gao 1996). These two spectral indices are complementary, with NDVI being particularly sensitive to the spatial and temporal variations in vegetation types and conditions, while NDWI is particularly sensitive to surface water conditions (e.g. Dupuy et al. 2020) for irrigated field mapping.

## 3.4. Calculation of the remote sensing ELVs

### 3.4.1. Spectral and temporal ELVs

Both spectral and temporal ELVs are based on an annual MODIS NDVI time series, hereafter referred to as the MODIS reference year. The MODIS reference year is designed to represent the mean seasonal variation of the surface, while attenuating any local or sudden surface changes due to particular climatic conditions. The reference year consists of 23 NDVI images (16-day composites), each calculated as the mean of the corresponding composites over a 5-year period (from November 2016 to October 2020, to include 5 full climatic years). The time period of 5 years was considered to be a good compromise, long enough to reduce the effect of annual climate variability and short enough to avoid large changes in land cover and land use.

The spectral ELV is calculated as the mean annual NDVI calculated over the MODIS reference year. The temporal ELVs are obtained by applying a Principal Component Analysis (PCA) transformation to the 23 NDVI values of the MODIS reference year. PCA is particularly suitable for capturing seasonal variations and identifying their relationship with different factors, whether climatic or related to human activities (Deng et al. 2008). In this study, PCA was performed using the FactoMineR package (Lê, Josse, and Husson 2008) of the R software and the rasterPCA tool of the RStoolbox package (Leutner et al. 2022).
<figure><img src="image_22.png" /></figure>

a)

<figure><img src="image_23.png" /></figure>

b)

<figure><img src="image_24.png" /><figcaption>**Figure 5.** Semi-variogram calculated on (a) the MODIS mean annual NDVI image (2016–2020); (b) The NDVI Sentinel-2 image (October 2019).</figcaption></figure>

### 3.4.2. Textural ELVs

Texture images are widely used for vegetation and land use mapping (e.g. Hudak and Wessman 1998; Kupidura 2019; Thierion et al. 2014). The textural ELVs used in this study were obtained using the gray level co-occurrence matrix (GLCM). Eight GLCM indices (energy, entropy, correlation, homogeneity, inertia, cluster shade, cluster prominence, Haralick's correlation; Haralick, Shanmugam, and Dinstein 1973) were calculated from the MODIS mean annual NDVI and the Sentinel-2 seasonal NDVI, using the Orfeo ToolBox (Grizonnet et al. 2017). To run the algorithm, two parameters have to be set up: the window size (windows radius), and the direction of the shift (offset).

The optimal window size is a trade-off between the largest variance and the smallest window size. The results of the semi-variograms (Figure 5) – using the R software package usdm (Naimi et al. 2014) – show that the optimal variance distance calculated on the MODIS mean annual NDVI image is about 1500 m (i.e. an optimal window size of 6 pixels), while on the Sentinel-2 images it is about 200 m (i.e. an optimal window size of 20 pixels), which corresponds to the order of magnitude of the distance between two crests in the study area. The directional displacement parameter is set to (1;1).

As the Haralick texture dataset contains highly correlated indices, a Principal Component Analysis (PCA) was applied to reduce the dimensionality and increase the interpretability of the data. Based on statistical criteria (a low correlation between indices and an high contribution to the axis), two indices were retained as textural ELVs (Appendix 1): the entropy, a measure of the disorder or complexity of an image, and the Haralick's correlation (H-correlation), a measure of the likelihood of finding two pixels of similar intensity separated by a given distance. The texture images were then resampled at 250 m resolution.

### 3.5. The data processing chain

The data processing chain consists of RS-ELVs image segmentation and objects clustering (Figure 3). In order to select the optimal RS-ELV dataset for each step, a large number of datasets combining different sensors and information types were tested. For the sake of clarity, only a few key results are presented here. Since the segmentation step is the most critical step in the GEOBIA approach (Chen et al. 2018), we present segmentation results obtained with three ELV combinations (Table 2) representative of a monosensor (A and B) and a multisensor (C) datasets. Then, only the “best” segmentation obtained is used for the clustering and landscape characterization steps.

#### 3.5.1. Image segmentation

For the segmentation step, we present the results obtained for 3 datasets composed of different RS-ELV types and sensors (Table 2):

* Dataset A (MODIS-based ELVs): MODIS mean annual NDVI + MODIS PC2, PC3, PC4 (4 variables);

* Dataset B (MODIS-based ELVs): Dataset A + 2 MODIS Haralick's indices (6 variables);

* Dataset C (MODIS-based and S2-based ELVs): Dataset A + 2 NDVI Sentinel-2 Haralick's indices on two dates (8 variables).

These datasets allow investigating the benefit of adding textural ELVs to spectral and temporal ELVs, when using MODIS data (dataset A vs. dataset B) and when using Sentinel-2 data (dataset B vs. dataset C).

The segmentations were performed using the Generic Region Merging algorithm of the Orfeo Tool Box (OTB) library (https://www.orfeo-toolbox.org/; Grizonnet et al. 2017). It is known to be a very effective image segmentation method for the analysis of remote
<figure><img src="image_25.png" /></figure>

a)

<figure><img src="image_26.png" /></figure>

b)

<figure><img src="image_27.png" /><figcaption>**Figure 5.** Semi-variogram calculated on (a) the MODIS mean annual NDVI image (2016–2020); (b) The NDVI Sentinel-2 image (October 2019).</figcaption></figure>

### 3.4.2. Textural ELVs

Texture images are widely used for vegetation and land use mapping (e.g. Hudak and Wessman 1998; Kupidura 2019; Thierion et al. 2014). The textural ELVs used in this study were obtained using the gray level co-occurrence matrix (GLCM). Eight GLCM indices (energy, entropy, correlation, homogeneity, inertia, cluster shade, cluster prominence, Haralick's correlation; Haralick, Shanmugam, and Dinstein 1973) were calculated from the MODIS mean annual NDVI and the Sentinel-2 seasonal NDVI, using the Orfeo ToolBox (Grizonnet et al. 2017). To run the algorithm, two parameters have to be set up: the window size (windows radius), and the direction of the shift (offset).

The optimal window size is a trade-off between the largest variance and the smallest window size. The results of the semi-variograms (Figure 5) – using the R software package usdm (Naimi et al. 2014) – show that the optimal variance distance calculated on the MODIS mean annual NDVI image is about 1500 m (i.e. an optimal window size of 6 pixels), while on the Sentinel-2 images it is about 200 m (i.e. an optimal window size of 20 pixels), which corresponds to the order of magnitude of the distance between two crests in the study area. The directional displacement parameter is set to (1;1).

As the Haralick texture dataset contains highly correlated indices, a Principal Component Analysis (PCA) was applied to reduce the dimensionality and increase the interpretability of the data. Based on statistical criteria (a low correlation between indices and an high contribution to the axis), two indices were retained as textural ELVs (Appendix 1): the entropy, a measure of the disorder or complexity of an image, and the Haralick's correlation (H-correlation), a measure of the likelihood of finding two pixels of similar intensity separated by a given distance. The texture images were then resampled at 250 m resolution.

### 3.5. The data processing chain

The data processing chain consists of RS-ELVs image segmentation and objects clustering (Figure 3). In order to select the optimal RS-ELV dataset for each step, a large number of datasets combining different sensors and information types were tested. For the sake of clarity, only a few key results are presented here. Since the segmentation step is the most critical step in the GEOBIA approach (Chen et al. 2018), we present segmentation results obtained with three ELV combinations (Table 2) representative of a monosensor (A and B) and a multisensor (C) datasets. Then, only the “best” segmentation obtained is used for the clustering and landscape characterization steps.

#### 3.5.1. Image segmentation

For the segmentation step, we present the results obtained for 3 datasets composed of different RS-ELV types and sensors (Table 2):

* Dataset A (MODIS-based ELVs): MODIS mean annual NDVI + MODIS PC2, PC3, PC4 (4 variables);

* Dataset B (MODIS-based ELVs): Dataset A + 2 MODIS Haralick's indices (6 variables);

* Dataset C (MODIS-based and S2-based ELVs): Dataset A + 2 NDVI Sentinel-2 Haralick's indices on two dates (8 variables).

These datasets allow investigating the benefit of adding textural ELVs to spectral and temporal ELVs, when using MODIS data (dataset A vs. dataset B) and when using Sentinel-2 data (dataset B vs. dataset C).

The segmentations were performed using the Generic Region Merging algorithm of the Orfeo Tool Box (OTB) library (https://www.orfeo-toolbox.org/; Grizonnet et al. 2017). It is known to be a very effective image segmentation method for the analysis of remote
<figure><img src="image_25.png" /><figcaption>**Table 2.** Synthesis of the RS-ELV datasets used for the segmentation (A, B, and C; see text), the clustering and the landscape characterization steps.</figcaption></figure>

<table><thead><tr><th rowspan="2">Source</th><th rowspan="2">RS-ELV</th><th colspan="3">Segmentation</th><th rowspan="2">Clustering</th><th rowspan="2">Characterisation</th></tr><tr><th>A</th><th>B</th><th>C</th></tr></thead><tbody><tr><td rowspan="4">MODIS</td><td>Annual NDVI (MODIS PC1)</td><td>X</td><td>X</td><td>X</td><td>X</td><td>X</td></tr><tr><td>MODIS PC2, PC3, PC4</td><td>X</td><td>X</td><td>X</td><td>X</td><td>X</td></tr><tr><td>Annual NDVI Entropy</td><td></td><td>X</td><td></td><td></td><td></td></tr><tr><td>Annual NDVI H-correlation</td><td></td><td>X</td><td></td><td></td><td></td></tr><tr><td rowspan="8">Sentinel-2</td><td>April NDVI Entropy</td><td></td><td></td><td>X</td><td></td><td>X</td></tr><tr><td>April NDVI Correlation</td><td></td><td></td><td>X</td><td></td><td>X</td></tr><tr><td>October NDVI Entropy</td><td></td><td></td><td>X</td><td></td><td>X</td></tr><tr><td>October NDVI Correlation</td><td></td><td></td><td>X</td><td></td><td>X</td></tr><tr><td>April NDVI</td><td></td><td></td><td></td><td></td><td>X</td></tr><tr><td>October NDVI</td><td></td><td></td><td></td><td></td><td>X</td></tr><tr><td>April NDWI</td><td></td><td></td><td></td><td></td><td>X</td></tr><tr><td>October NDWI</td><td></td><td></td><td></td><td></td><td>X</td></tr></tbody></table>

sensing images (e.g. Meinel and Neubert 2004). The algorithm starts with a pixel forming an object and merges the neighboring pixels until the Baatz and Schape homogeneity criterion is reached (Baatz and Schape 2000). Two parameters are defined to perform the segmentation: the scale size, which determines the maximum allowed heterogeneity of the objects and thus the size of the segments, and the relative weight of the spectral and spatial homogeneity criteria. In this study, 25 levels of segmentation were tested (scales from 100 to 2500, in steps of 100), and only the spectral criterion was considered as no hypothesis was made about the shape of the landscape unit.

The statistical method chosen to evaluate the segmentations is derived from the method of Johnson and Xie (2011) and applied by Bisquert, Bégué, and Deshayes (2015) for landscape mapping. This method is adapted to multispectral imagery and allows the analysis of the disparity between neighboring segments. It is based on variance as a measure of homogeneity within segments (V), and on Moran's autocorrelation index (M) as a measure of similarity between segments (Espindola et al. 2006). V and M are calculated for each variable, averaged, and normalized for each segmentation scale (see Appendix 2 for the equations). The global score is then calculated as the sum of the normalized V and M values. Two versions of the global score are presented below: the original J score from (Johnson and Xie 2011), and a modified version, the JB score, from Böck, Immitzer, and Atzberger (2017) who proposed to normalize V and M in a way that mitigates the J score instability to the segmentation scales. In either case, the optimal segmentation was identified as the one with the lowest average global score.

### 3.5.2. Clustering

Unsupervised approach is largely used in land mapping from imagery, because it requires no training dataset. Among the traditional method, the k-means

clustering algorithm is widely employed for the analysis of remote sensing images because it is easy to apply and widely available in image processing and statistical software packages (e.g. Duda and Canty 2002; Xie, Sha, and Yu 2008).

The landscape units were obtained through a k-means clustering of the data objects. The segmented objects are the one produced in the previous step, and the variables are the mean values of the MODIS annual NDVI, MODIS PC2, PC3 and PC4 (Table 2). In the study, a range of cluster numbers between 2 and 15 was tested, and the optimal number of clusters was determined using the Elbow criterion method. For the clustering, we used the Scikit-learn package (Pedregosa et al. 2011), a Python module that integrates state-of-the-art machine learning algorithms (https://scikit-learn.org/stable/modules/clustering).

### 3.5.3. Landscape characterisation

Landscape characterization is a way of labeling landscape classes into landscape types (Figure 1). Each landscape class was characterized by a set of MODIS and Sentinel-2 ELVs (Table 2), complementary to the variables used for the data processing, and by the land cover composition. For the latter, we used the ESA WorldCover project 2020 global land cover map produced from Sentinel-2 and Sentinel-1 data at 10 m spatial resolution (Zanaga et al. 2021).

## 4. Results

### 4.1. The remote sensing ELVs

#### 4.1.1. The spectral and temporal ELVs

A PCA was performed on the MODIS reference year. The first 4 principal components (Figure 6) concentrate most of the spatial and phenological information, with 99.66% of the total variance explained.

The explained variance of PC1 is 95.33% with an equal contribution of each variable to -0.28 approximately (Figure 6c) and a magnitude of the
<figure><img src="image_28.png" /><figcaption>**Table 2.** Synthesis of the RS-ELV datasets used for the segmentation (A, B, and C; see text), the clustering and the landscape characterization steps.</figcaption></figure>

<table><thead><tr><th rowspan="2">Source</th><th rowspan="2">RS-ELV</th><th colspan="3">Segmentation</th><th rowspan="2">Clustering</th><th rowspan="2">Characterisation</th></tr><tr><th>A</th><th>B</th><th>C</th></tr></thead><tbody><tr><td rowspan="4">MODIS</td><td>Annual NDVI (MODIS PC1)</td><td>X</td><td>X</td><td>X</td><td>X</td><td>X</td></tr><tr><td>MODIS PC2, PC3, PC4</td><td>X</td><td>X</td><td>X</td><td>X</td><td>X</td></tr><tr><td>Annual NDVI Entropy</td><td></td><td>X</td><td></td><td></td><td></td></tr><tr><td>Annual NDVI H-correlation</td><td></td><td>X</td><td></td><td></td><td></td></tr><tr><td rowspan="8">Sentinel-2</td><td>April NDVI Entropy</td><td></td><td></td><td>X</td><td></td><td>X</td></tr><tr><td>April NDVI Correlation</td><td></td><td></td><td>X</td><td></td><td>X</td></tr><tr><td>October NDVI Entropy</td><td></td><td></td><td>X</td><td></td><td>X</td></tr><tr><td>October NDVI Correlation</td><td></td><td></td><td>X</td><td></td><td>X</td></tr><tr><td>April NDVI</td><td></td><td></td><td></td><td></td><td>X</td></tr><tr><td>October NDVI</td><td></td><td></td><td></td><td></td><td>X</td></tr><tr><td>April NDWI</td><td></td><td></td><td></td><td></td><td>X</td></tr><tr><td>October NDWI</td><td></td><td></td><td></td><td></td><td>X</td></tr></tbody></table>

sensing images (e.g. Meinel and Neubert 2004). The algorithm starts with a pixel forming an object and merges the neighboring pixels until the Baatz and Schape homogeneity criterion is reached (Baatz and Schape 2000). Two parameters are defined to perform the segmentation: the scale size, which determines the maximum allowed heterogeneity of the objects and thus the size of the segments, and the relative weight of the spectral and spatial homogeneity criteria. In this study, 25 levels of segmentation were tested (scales from 100 to 2500, in steps of 100), and only the spectral criterion was considered as no hypothesis was made about the shape of the landscape unit.

The statistical method chosen to evaluate the segmentations is derived from the method of Johnson and Xie (2011) and applied by Bisquert, Bégué, and Deshayes (2015) for landscape mapping. This method is adapted to multispectral imagery and allows the analysis of the disparity between neighboring segments. It is based on variance as a measure of homogeneity within segments (V), and on Moran's autocorrelation index (M) as a measure of similarity between segments (Espindola et al. 2006). V and M are calculated for each variable, averaged, and normalized for each segmentation scale (see Appendix 2 for the equations). The global score is then calculated as the sum of the normalized V and M values. Two versions of the global score are presented below: the original J score from (Johnson and Xie 2011), and a modified version, the JB score, from Böck, Immitzer, and Atzberger (2017) who proposed to normalize V and M in a way that mitigates the J score instability to the segmentation scales. In either case, the optimal segmentation was identified as the one with the lowest average global score.

### 3.5.2. Clustering

Unsupervised approach is largely used in land mapping from imagery, because it requires no training dataset. Among the traditional method, the k-means

clustering algorithm is widely employed for the analysis of remote sensing images because it is easy to apply and widely available in image processing and statistical software packages (e.g. Duda and Canty 2002; Xie, Sha, and Yu 2008).

The landscape units were obtained through a k-means clustering of the data objects. The segmented objects are the one produced in the previous step, and the variables are the mean values of the MODIS annual NDVI, MODIS PC2, PC3 and PC4 (Table 2). In the study, a range of cluster numbers between 2 and 15 was tested, and the optimal number of clusters was determined using the Elbow criterion method. For the clustering, we used the Scikit-learn package (Pedregosa et al. 2011), a Python module that integrates state-of-the-art machine learning algorithms (https://scikit-learn.org/stable/modules/clustering).

### 3.5.3. Landscape characterisation

Landscape characterization is a way of labeling landscape classes into landscape types (Figure 1). Each landscape class was characterized by a set of MODIS and Sentinel-2 ELVs (Table 2), complementary to the variables used for the data processing, and by the land cover composition. For the latter, we used the ESA WorldCover project 2020 global land cover map produced from Sentinel-2 and Sentinel-1 data at 10 m spatial resolution (Zanaga et al. 2021).

## 4. Results

### 4.1. The remote sensing ELVs

#### 4.1.1. The spectral and temporal ELVs

A PCA was performed on the MODIS reference year. The first 4 principal components (Figure 6) concentrate most of the spatial and phenological information, with 99.66% of the total variance explained.

The explained variance of PC1 is 95.33% with an equal contribution of each variable to -0.28 approximately (Figure 6c) and a magnitude of the
<figure><img src="image_26.png" /></figure>

<figure><img src="image_27.png" /><figcaption>**Figure 6.** Results of the MODIS NDVI PCA calculated over the reference year (2016–2020): (a) Images of the principal components 1 to 4; (b) Colored composition of the first three principal components (PC1 in red, PC2 in green, PC3 in Blue); (c) Eigenvector magnitudes of the first four principal components and mean monthly rainfall (source: GPCP data) over a climatic year.</figcaption></figure>

first eigenvector of 1. That means that each variable is highly correlated to this axis which carries the most important part of the image signal, that is the spatial information (see high PC1 values in the forest corridor; Figure 6a). It is therefore expected that the temporal information is carried by the remaining principal components. PC2 and PC3 (3.25% and 0.89% of variance, respectively) show indeed a strong seasonal cycle, in relation with the rainfall temporal pattern (Figure 6c), with a one period repetition in phase opposition regarding PC2, and at a 4-month time shift regarding PC3. PC2 and PC3 images (Figure 6a) show a strong east-west gradient of the vegetation phenology, with high values in the western lowlands for PC3 possibly linked to irrigated agriculture. PC4 presents a 2-period temporal signal. Even though carrying only 0.19% of the information, PC4 captures isolated changes and episodes that deviate from the dominant seasonal variation (Bellón et al. 2017). These latter could be related to different natural or anthropogenic (such

as double cropping in the lowlands) factors or combinations of factors.

### 4.1.2. The textural ELVs

Figure 7(a, d) show the entropy and H-correlation indices calculated from MODIS annual NDVI. As expected these two indices provide different information, with a longitudinal alternance of high (west of corridor with mosaic of land cover) and low (forest corridor with a “simple” texture) values of entropy, and an east-west gradient for H-correlation. Both indices tend to emphasize the contours of certain landscape elements such as the forest corridor or the hydrographic network (this is particularly evident for H-Correlation).

Looking at the Sentinel-2 texture images, we observe that the textural ELVs have a strong seasonal dimension, with higher texture index values in October (Figure 7c, f) than in April (Figure 7b, e). At the end of the rainy season (April), the vegetation, both natural and cultivated, is well developed and the
<figure><img src="image_29.png" /></figure>

<figure><img src="image_30.png" /><figcaption>**Figure 6.** Results of the MODIS NDVI PCA calculated over the reference year (2016–2020): (a) Images of the principal components 1 to 4; (b) Colored composition of the first three principal components (PC1 in red, PC2 in green, PC3 in Blue); (c) Eigenvector magnitudes of the first four principal components and mean monthly rainfall (source: GPCP data) over a climatic year.</figcaption></figure>

first eigenvector of 1. That means that each variable is highly correlated to this axis which carries the most important part of the image signal, that is the spatial information (see high PC1 values in the forest corridor; Figure 6a). It is therefore expected that the temporal information is carried by the remaining principal components. PC2 and PC3 (3.25% and 0.89% of variance, respectively) show indeed a strong seasonal cycle, in relation with the rainfall temporal pattern (Figure 6c), with a one period repetition in phase opposition regarding PC2, and at a 4-month time shift regarding PC3. PC2 and PC3 images (Figure 6a) show a strong east-west gradient of the vegetation phenology, with high values in the western lowlands for PC3 possibly linked to irrigated agriculture. PC4 presents a 2-period temporal signal. Even though carrying only 0.19% of the information, PC4 captures isolated changes and episodes that deviate from the dominant seasonal variation (Bellón et al. 2017). These latter could be related to different natural or anthropogenic (such

as double cropping in the lowlands) factors or combinations of factors.

### 4.1.2. The textural ELVs

Figure 7(a, d) show the entropy and H-correlation indices calculated from MODIS annual NDVI. As expected these two indices provide different information, with a longitudinal alternance of high (west of corridor with mosaic of land cover) and low (forest corridor with a “simple” texture) values of entropy, and an east-west gradient for H-correlation. Both indices tend to emphasize the contours of certain landscape elements such as the forest corridor or the hydrographic network (this is particularly evident for H-Correlation).

Looking at the Sentinel-2 texture images, we observe that the textural ELVs have a strong seasonal dimension, with higher texture index values in October (Figure 7c, f) than in April (Figure 7b, e). At the end of the rainy season (April), the vegetation, both natural and cultivated, is well developed and the
<figure><img src="image_28.png" /></figure>

<figure><img src="image_29.png" /><figcaption>**Figure 7.** Entropy (up line) and H-correlation (bottom line) images calculated for MODIS annual NDVI (a and d, respectively), NDVI Sentinel-2 April (b and e, respectively) and NDVI Sentinel-2 October (c and f, respectively).</figcaption></figure>

<figure><img src="image_30.png" /><figcaption>**Figure 8.** Zoom-in (50x50 km) of segmentations A, B and C, obtained with a scale parameter of 900: (a) MODIS spectral and temporal indices; (b) MODIS spectral, temporal and textural indices; (c) MODIS spectral and temporal indices, and Sentinel-2 textural indices. See Table 2 for ELVs details. The background image is Bing Aerial.</figcaption></figure>

land surface is essentially green, resulting in low values of texture indices, except for non-vegetated areas such as open waters or irrigated lowlands. At the end of the dry season (October), the contrast between vegetation cover types is exacerbated by the diversity of the ecosystems present (herbaceous and woody vegetation) and by human activity (irrigated and rainfed cropping). This high contrast in terms of vegetation conditions leads to high values of texture indices especially in the western part of the study area, which has drier conditions than the eastern part.

## 4.2. Segmentation and clustering results

### 4.2.1. Results obtained with different input variables

Figure 8 shows three examples of segmentation obtained with different ELVs datasets (Table 2), at scale 900. The results show different image partitioning depending on the dataset, with segments for segmentations A and C (Figure 8a, c) that appear to correctly follow the contours of the forest corridor. On the other hand, segmentation B (Figure 8b) shows artifacts, such as duplication
<figure><img src="image_31.png" /></figure>

<figure><img src="image_32.png" /><figcaption>**Figure 7.** Entropy (up line) and H-correlation (bottom line) images calculated for MODIS annual NDVI (a and d, respectively), NDVI Sentinel-2 April (b and e, respectively) and NDVI Sentinel-2 October (c and f, respectively).</figcaption></figure>

<figure><img src="image_33.png" /><figcaption>**Figure 8.** Zoom-in (50x50 km) of segmentations A, B and C, obtained with a scale parameter of 900: (a) MODIS spectral and temporal indices; (b) MODIS spectral, temporal and textural indices; (c) MODIS spectral and temporal indices, and Sentinel-2 textural indices. See Table 2 for ELVs details. The background image is Bing Aerial.</figcaption></figure>

land surface is essentially green, resulting in low values of texture indices, except for non-vegetated areas such as open waters or irrigated lowlands. At the end of the dry season (October), the contrast between vegetation cover types is exacerbated by the diversity of the ecosystems present (herbaceous and woody vegetation) and by human activity (irrigated and rainfed cropping). This high contrast in terms of vegetation conditions leads to high values of texture indices especially in the western part of the study area, which has drier conditions than the eastern part.

## 4.2. Segmentation and clustering results

### 4.2.1. Results obtained with different input variables

Figure 8 shows three examples of segmentation obtained with different ELVs datasets (Table 2), at scale 900. The results show different image partitioning depending on the dataset, with segments for segmentations A and C (Figure 8a, c) that appear to correctly follow the contours of the forest corridor. On the other hand, segmentation B (Figure 8b) shows artifacts, such as duplication
<figure><img src="image_31.png" /></figure>

<figure><img src="image_32.png" /><figcaption>**Figure 9.** Results for the global scores J (a) and JB (b) calculated for the set of test segmentations (scale between 100 and 2500, and datasets A, B and C detailed in Table 2); the blue line corresponds to the MODIS spectral and temporal ELVs, the orange line corresponds to the MODIS spectral, temporal and textural ELVs; the green line is for the MODIS spectral and temporal, and sentinel-2 textural ELVs. c) mean number of segments obtained for each scale parameter (logarithmic y-axis).</figcaption></figure>

of the forest corridor segmentation contours, which are clearly due to the MODIS texture indices (Figure 7d).

### 4.2.2. Statistical selection of the scale parameter

The best combination of ELVs and the optimal scale parameter for delineating homogeneous landscape structure units were identified through the global scores analysis. The results of the J and JB scores obtained for the three segmentation datasets are shown in Figure 9a, b, respectively. Low values of J scores are obtained in the ranges 900–1100, 100–600, 900–1800 for the segmentations with datasets A, B and C respectively, with high variations between the datasets, especially for dataset B. The range of variation is lower for JB between, and the minimum values appear less distinct. For the rest of the analysis, the segmentation performed with dataset A and a scale of 900 was chosen, as it gives the lowest J and JB scores, while maintaining a sufficient number of segments for the clustering step (18 segments; Figure 9c).

This choice is supported by the examples in Figure 10, where the scale parameter 100 produces an over-segmentation, the parameter 2000 an under-segmentation, while the parameter 900 produces a more balanced segmentation.

### 4.2.3. Clustering results

The elbow test shows that for segmentation of dataset A and scale 900, the first decrease in the slope of inertia occurs at 6 clusters (Figure 11). The final landscape map in 6 classes is shown in Figure 12.

## 4.3. Radiometric landscape characterisation

### 4.3.1. Landscape characterisation with remote sensing variables

The six classes of the radiometric landscape map were first evaluated using a radar chart (Figure 13) representing the normalized mean values of a selection of MODIS and Sentinel-2 ELVs (the mean and standard deviation values for each class are given in Appendix 3).

A closer look at Figure 13 and the ELVs statistics (Appendix 3) shows that the classes 5, 3 and 4 are characterized by high spectral index values (>0.67 NDVI), and low texture indices throughout the year. They correspond to rather homogeneous permanent vegetated areas; among them, class 4 stands out with a very high vegetation (>0.8) and water content (>0.6) index values throughout the year, corresponding to a dense tropical forest. The spectral and spatial proximity of classes 3 and 4 (Figure 12), with its elongated north-south shape
<figure><img src="image_34.png" /></figure>

<figure><img src="image_35.png" /><figcaption>**Figure 9.** Results for the global scores J (a) and JB (b) calculated for the set of test segmentations (scale between 100 and 2500, and datasets A, B and C detailed in Table 2); the blue line corresponds to the MODIS spectral and temporal ELVs, the orange line corresponds to the MODIS spectral, temporal and textural ELVs; the green line is for the MODIS spectral and temporal, and sentinel-2 textural ELVs. c) mean number of segments obtained for each scale parameter (logarithmic y-axis).</figcaption></figure>

of the forest corridor segmentation contours, which are clearly due to the MODIS texture indices (Figure 7d).

### 4.2.2. Statistical selection of the scale parameter

The best combination of ELVs and the optimal scale parameter for delineating homogeneous landscape structure units were identified through the global scores analysis. The results of the J and JB scores obtained for the three segmentation datasets are shown in Figure 9a, b, respectively. Low values of J scores are obtained in the ranges 900–1100, 100–600, 900–1800 for the segmentations with datasets A, B and C respectively, with high variations between the datasets, especially for dataset B. The range of variation is lower for JB between, and the minimum values appear less distinct. For the rest of the analysis, the segmentation performed with dataset A and a scale of 900 was chosen, as it gives the lowest J and JB scores, while maintaining a sufficient number of segments for the clustering step (18 segments; Figure 9c).

This choice is supported by the examples in Figure 10, where the scale parameter 100 produces an over-segmentation, the parameter 2000 an under-segmentation, while the parameter 900 produces a more balanced segmentation.

### 4.2.3. Clustering results

The elbow test shows that for segmentation of dataset A and scale 900, the first decrease in the slope of inertia occurs at 6 clusters (Figure 11). The final landscape map in 6 classes is shown in Figure 12.

## 4.3. Radiometric landscape characterisation

### 4.3.1. Landscape characterisation with remote sensing variables

The six classes of the radiometric landscape map were first evaluated using a radar chart (Figure 13) representing the normalized mean values of a selection of MODIS and Sentinel-2 ELVs (the mean and standard deviation values for each class are given in Appendix 3).

A closer look at Figure 13 and the ELVs statistics (Appendix 3) shows that the classes 5, 3 and 4 are characterized by high spectral index values (>0.67 NDVI), and low texture indices throughout the year. They correspond to rather homogeneous permanent vegetated areas; among them, class 4 stands out with a very high vegetation (>0.8) and water content (>0.6) index values throughout the year, corresponding to a dense tropical forest. The spectral and spatial proximity of classes 3 and 4 (Figure 12), with its elongated north-south shape
<figure><img src="image_33.png" /></figure>

<figure><img src="image_34.png" /><figcaption>**Figure 10.** Zoom-in (50x50 km) of segmentation obtained with different scale parameters: (a) 100, (b) 900, and (c) 2000. The background image is Bing Aerial.</figcaption></figure>

<figure><img src="image_35.png" /><figcaption>**Figure 11.** Cluster inertia plotted against the number of classes. The vertical line indicates the number of clusters used for the final clustering.</figcaption></figure>

following the local topography, indicates that class 3 is potentially a degraded forest.

In contrast, classes 1 and 6 are characterized by high texture values in all seasons, and low to medium vegetation indices (MODIS NDVI mean between 0.47 and 0.67), and could thus be interpreted as mosaic landscapes of crop and natural vegetation. However class 1 is characterized by a high contrast in seasonal NDWI and a very low value in October (0.045) indicating dry vegetation conditions; this is consistent with the geographic location of the class 1 in the western part of the study area, in the highlands, where the annual rainfall is moderate, in contrast to class 6 which is in a very rainy area.

Class 2 has a particular signature, with medium vegetation development during the year and high texture values especially at the end of the dry season; this

<figure><img src="image_36.png" /><figcaption>**Figure 12.** The final radiometric landscape map in six classes (the grey lines represent the roads network).</figcaption></figure>
<figure><img src="image_36.png" /></figure>

<figure><img src="image_37.png" /><figcaption>**Figure 10.** Zoom-in (50x50 km) of segmentation obtained with different scale parameters: (a) 100, (b) 900, and (c) 2000. The background image is Bing Aerial.</figcaption></figure>

<figure><img src="image_38.png" /><figcaption>**Figure 11.** Cluster inertia plotted against the number of classes. The vertical line indicates the number of clusters used for the final clustering.</figcaption></figure>

following the local topography, indicates that class 3 is potentially a degraded forest.

In contrast, classes 1 and 6 are characterized by high texture values in all seasons, and low to medium vegetation indices (MODIS NDVI mean between 0.47 and 0.67), and could thus be interpreted as mosaic landscapes of crop and natural vegetation. However class 1 is characterized by a high contrast in seasonal NDWI and a very low value in October (0.045) indicating dry vegetation conditions; this is consistent with the geographic location of the class 1 in the western part of the study area, in the highlands, where the annual rainfall is moderate, in contrast to class 6 which is in a very rainy area.

Class 2 has a particular signature, with medium vegetation development during the year and high texture values especially at the end of the dry season; this

<figure><img src="image_39.png" /><figcaption>**Figure 12.** The final radiometric landscape map in six classes (the grey lines represent the roads network).</figcaption></figure>
<figure><img src="image_37.png" /></figure>

<figure><img src="image_38.png" /><figcaption>**Figure 13.** Radar chart composed of the normalised mean values of MODIS and Sentinel-2 ELVs, for the six classes of the radiometric landscape map.</figcaption></figure>

<figure><img src="image_39.png" /><figcaption>**Figure 14.** Landscape type characterisation in terms of land cover : (a) ESA World cover map (2020) with cluster boundaries; (b) Frequency distribution (%) of the ESA global cover land cover classes (2020), by landscape class (labels of the LULC classes representing less than 1% are not shown).</figcaption></figure>

particularity could be related to the presence of irrigated agriculture (rice in the lowlands) which increases the spatial heterogeneity.

### **4.3.2. Landscape class characterisation with land cover map**

When the radiometric landscape map is overlaid with the ESA WorldCover project 2020 global land use map, we observe that the landscape classes are well differentiated in terms of land use and land cover (LULC) composition and proportion (Figure 14). Ranked by landscape complexity, we have:

* Landscape type with one dominant LULC class (>90%): class 4 is composed almost exclusively of tree cover (90.2%) and corresponds to the forest corridor;

* Landscape type with two dominant LULC classes (total >90%): class 1, in the western part of the study area, is composed of about ⅔ of grassland and ⅓ of cropland, and corresponds to the highlands;

* Landscape type composed of three LULC classes: class 5, located in the eastern part of the study area, is dominated by shrubland (58.4%), while class 3, adjacent to the forest corridor, has a large proportion of tree cover (50.7%). Both landscape classes are equally composed of grassland (around 20%);

* Landscape type composed of mixed LULC classes : mosaic of four land cover types (forest, cropland, grassland and shrubland) for classes 2 and 6.

Finally, it is interesting to note the exclusive presence of herbaceous wetlands in landscape class 6, along the hydrographic network.

### **4.3.3. Landscape class characterisation with ground observations**

Based on field observations, we were able to associate the panoramic photographs to five of the six classes of the radiometric landscape map (Figure 15). As the south-eastern part of the study
<figure><img src="image_40.png" /></figure>

<figure><img src="image_41.png" /><figcaption>**Figure 13.** Radar chart composed of the normalised mean values of MODIS and Sentinel-2 ELVs, for the six classes of the radiometric landscape map.</figcaption></figure>

<figure><img src="image_42.png" /><figcaption>**Figure 14.** Landscape type characterisation in terms of land cover : (a) ESA World cover map (2020) with cluster boundaries; (b) Frequency distribution (%) of the ESA global cover land cover classes (2020), by landscape class (labels of the LULC classes representing less than 1% are not shown).</figcaption></figure>

particularity could be related to the presence of irrigated agriculture (rice in the lowlands) which increases the spatial heterogeneity.

### **4.3.2. Landscape class characterisation with land cover map**

When the radiometric landscape map is overlaid with the ESA WorldCover project 2020 global land use map, we observe that the landscape classes are well differentiated in terms of land use and land cover (LULC) composition and proportion (Figure 14). Ranked by landscape complexity, we have:

* Landscape type with one dominant LULC class (>90%): class 4 is composed almost exclusively of tree cover (90.2%) and corresponds to the forest corridor;

* Landscape type with two dominant LULC classes (total >90%): class 1, in the western part of the study area, is composed of about ⅔ of grassland and ⅓ of cropland, and corresponds to the highlands;

* Landscape type composed of three LULC classes: class 5, located in the eastern part of the study area, is dominated by shrubland (58.4%), while class 3, adjacent to the forest corridor, has a large proportion of tree cover (50.7%). Both landscape classes are equally composed of grassland (around 20%);

* Landscape type composed of mixed LULC classes : mosaic of four land cover types (forest, cropland, grassland and shrubland) for classes 2 and 6.

Finally, it is interesting to note the exclusive presence of herbaceous wetlands in landscape class 6, along the hydrographic network.

### **4.3.3. Landscape class characterisation with ground observations**

Based on field observations, we were able to associate the panoramic photographs to five of the six classes of the radiometric landscape map (Figure 15). As the south-eastern part of the study
<figure><img src="image_40.png" /></figure>

<figure><img src="image_41.png" /><figcaption>**Figure 15.** Illustration of the landscape classes (classes 1 to 5, from top to bottom), with panoramic photographs acquired by drone (Alleaume 2022).</figcaption></figure>

area could not be visited during the field trip, the landscape class 6 is not illustrated and discussed in this section.

* Landscape class 1 is an open landscape. It takes the form of large areas of grazed grass-land on hilltops, ridges and slopes. There are a few root crops (mainly cassava and ground-nuts) and patches of eucalyptus. The lowlands are narrow to medium and are cultivated with rice. The settlements and villages associated with this landscape are scattered and sparse.

* Landscape class 2 is characterized by large patches of cultivated land, eucalyptus and savan-nah. The lowlands are narrow to very wide, and are mainly cultivated with rice. The hills are reg-ular, sometimes with cliffs or granite boulders. The terraces at the foot of the slopes are planted with rice or other crops, and sometimes extend up into the secondary valleys. On the slopes, crops and grasslands are present. Eucalyptus patches are more or less wide and numerous, often at the top of the hills. Villages range from

small to large, and are sometimes densely meshed.

* Landscape class 3 is characterized by large areas of eucalyptus with some patches of mature and mixed forest (pine, acacia associated with mature forest species). The eucalyptus plantations are of different ages, resulting in a heterogeneous woodland cover. The lowlands are more or less narrow, and are cultivated with rice and other staple crops. Few villages are visible in this land-scape, but in areas of greater human pressure, terraces can rise up the slopes and the built-up area becomes denser.

* Landscape class 4 corresponds to the mature forest zone of the corridor. It is characterized by a large area of mature forest, with few openings made of wetlands, rock walls with varying vegetation cover.

* Landscape class 5 corresponds to a large area covered by a mixture of crops, scrubland, exotic woodland, Ravenala (traveler's palm), banana trees, small patches of mature, mixed and euca-lyptus forests. The blurred boundaries between landscape components contribute to an
<figure><img src="image_43.png" /></figure>

<figure><img src="image_44.png" /><figcaption>**Figure 15.** Illustration of the landscape classes (classes 1 to 5, from top to bottom), with panoramic photographs acquired by drone (Alleaume 2022).</figcaption></figure>

area could not be visited during the field trip, the landscape class 6 is not illustrated and discussed in this section.

* Landscape class 1 is an open landscape. It takes the form of large areas of grazed grassland on hilltops, ridges and slopes. There are a few root crops (mainly cassava and groundnuts) and patches of eucalyptus. The lowlands are narrow to medium and are cultivated with rice. The settlements and villages associated with this landscape are scattered and sparse.

* Landscape class 2 is characterized by large patches of cultivated land, eucalyptus and savannah. The lowlands are narrow to very wide, and are mainly cultivated with rice. The hills are regular, sometimes with cliffs or granite boulders. The terraces at the foot of the slopes are planted with rice or other crops, and sometimes extend up into the secondary valleys. On the slopes, crops and grasslands are present. Eucalyptus patches are more or less wide and numerous, often at the top of the hills. Villages range from

small to large, and are sometimes densely meshed.

* Landscape class 3 is characterized by large areas of eucalyptus with some patches of mature and mixed forest (pine, acacia associated with mature forest species). The eucalyptus plantations are of different ages, resulting in a heterogeneous woodland cover. The lowlands are more or less narrow, and are cultivated with rice and other staple crops. Few villages are visible in this landscape, but in areas of greater human pressure, terraces can rise up the slopes and the built-up area becomes denser.

* Landscape class 4 corresponds to the mature forest zone of the corridor. It is characterized by a large area of mature forest, with few openings made of wetlands, rock walls with varying vegetation cover.

* Landscape class 5 corresponds to a large area covered by a mixture of crops, scrubland, exotic woodland, Ravenala (traveler's palm), banana trees, small patches of mature, mixed and eucalyptus forests. The blurred boundaries between landscape components contribute to an
<figure><img src="image_42.png" /></figure>

impression of great heterogeneity. The hills are generally of the same height. However, this type of landscape is also present at the edges of the great cliffs that follow one another between the corridor and the east coast. The upper parts of the large cliffs are often occupied by patches of mature forest.

# 5. Discussion

## 5.1. The radiometric landscape: a new concept and method

In this paper we have proposed an original concept that lays the foundation for an innovative approach to landscape mapping based solely on Earth Observation data. We first established a list of Remote Sensing-Landscape Essential Variables (RS-ELVs), based on physiognomic landscape modeling, and applied GEOBIA approach to produce the so-called radiometric landscape map.

We are aware that the objective of characterizing and mapping the landscape using remote sensing data is an ambitious one. Indeed, it combines two difficulties related to the object of study, the landscape, and the approach adopted, GEOBIA. On the one hand, the notion of landscape is multifaceted and the landscape depends on the observer's point of view, which raises the problem of evaluating the quality of landscape maps, among others. On the other hand, the GEOBIA approach raises conceptual and implementation issues (Hay and Castilla 2006) such as the segmentation, an ill-posed problem in the sense that it has no unique solution, and unsupervised classification, with an optimal number of clusters difficult to define. The potential and the main limitations of our approach are discussed hereafter.

### 5.1.1. The RS-ELVs

For the temporal and spectral ELVs, we used MODIS SITS which proved to be a data source offering a good compromise between temporal and spatial resolution. In agreement with other authors (e.g. Knight et al. 2006), who used MODIS NDVI SITS 250 m to include plant phenology for land cover characterization, our results confirm that the temporal information contained in the principal components of the MODIS NDVI SITS is essential for landscape characterization and mapping – as they appear to be directly correlated with temporal cycles such as vegetation phenology, rainfall, and agricultural practices – and are important RS-ELVs. The ELVs approach developed here provides a multi-temporal baseline that both enhances the seasonal variation of the landscape, and reduces the sensitivity to local and regional climatic events

(and thus the sensitivity to the data image quality). These temporal and spectral RS-ELVs could theoretically be extracted from other satellite image time series, but currently there are few alternatives as the Sentinel-2 acquisition frequency (5 days) is too loose to obtain acceptable SITS in tropical areas due to the atmospheric conditions.

For the textural ELVs, results showed that MODIS data were not suitable due to texture artifacts, and we therefore recommended the use of Sentinel-2 texture indices calculated at different seasons instead. However, we found that the joint use of temporal variables (MODIS) and S2 texture variables does not provide significant added value compared to the use of temporal variables alone.

Despite the initial assumption that the spectral and temporal RS-ELVs indirectly include the three-dimensional environment, the inconsistency observed for landscape class 2 (composed of landscape units located on both the western and eastern sides of the forest corridor) could mean that the RS-ELVs we propose are not discriminating enough. One way of correcting this type of problem would be to use a DEM, as originally considered (Table 1). Topographic data derived from Earth observation are available worldwide. Therefore we recommend, for regions with high relief, to include topographic variables in the landscape mapping process.

### 5.1.2. The algorithms

The approach is based on GEOBIA with segmentation and clustering steps. As shown in this paper, GEOBIA allows to overcome the pixel view and provide image objects that make sense, opening a new dimension in rule-based automated image analysis (S. Lang et al. 2019). However, GEOBIA implementation is not straightforward, and different technical options are available for each step.

In this paper, we have used algorithms (region merging and k-means) that are known to be robust and adapted to our goal. We are aware of the inherent instability of some of these algorithms, and their sensitivity to parameterization. The k-means algorithm is known to be highly dependent on the initial cluster centers (Peña, Lozano, and Larrañaga 1999), which can lead to different clustering results (Li and Wu 2012), but it is rather simple to implement and efficient. Regarding the sensitivity of the algorithms to parameterization, a statistical optimization of the parameters was proposed (number of classes with an unsupervised method that maximizes the homogeneity within segments and heterogeneity between adjacent segments, and number of clusters with Elbow's method) to guide the landscape mapping process.
<figure><img src="image_45.png" /></figure>

impression of great heterogeneity. The hills are generally of the same height. However, this type of landscape is also present at the edges of the great cliffs that follow one another between the corridor and the east coast. The upper parts of the large cliffs are often occupied by patches of mature forest.

# 5. Discussion

## 5.1. The radiometric landscape: a new concept and method

In this paper we have proposed an original concept that lays the foundation for an innovative approach to landscape mapping based solely on Earth Observation data. We first established a list of Remote Sensing-Landscape Essential Variables (RS-ELVs), based on physiognomic landscape modeling, and applied GEOBIA approach to produce the so-called radiometric landscape map.

We are aware that the objective of characterizing and mapping the landscape using remote sensing data is an ambitious one. Indeed, it combines two difficulties related to the object of study, the landscape, and the approach adopted, GEOBIA. On the one hand, the notion of landscape is multifaceted and the landscape depends on the observer's point of view, which raises the problem of evaluating the quality of landscape maps, among others. On the other hand, the GEOBIA approach raises conceptual and implementation issues (Hay and Castilla 2006) such as the segmentation, an ill-posed problem in the sense that it has no unique solution, and unsupervised classification, with an optimal number of clusters difficult to define. The potential and the main limitations of our approach are discussed hereafter.

### 5.1.1. The RS-ELVs

For the temporal and spectral ELVs, we used MODIS SITS which proved to be a data source offering a good compromise between temporal and spatial resolution. In agreement with other authors (e.g. Knight et al. 2006), who used MODIS NDVI SITS 250 m to include plant phenology for land cover characterization, our results confirm that the temporal information contained in the principal components of the MODIS NDVI SITS is essential for landscape characterization and mapping – as they appear to be directly correlated with temporal cycles such as vegetation phenology, rainfall, and agricultural practices – and are important RS-ELVs. The ELVs approach developed here provides a multi-temporal baseline that both enhances the seasonal variation of the landscape, and reduces the sensitivity to local and regional climatic events

(and thus the sensitivity to the data image quality). These temporal and spectral RS-ELVs could theoretically be extracted from other satellite image time series, but currently there are few alternatives as the Sentinel-2 acquisition frequency (5 days) is too loose to obtain acceptable SITS in tropical areas due to the atmospheric conditions.

For the textural ELVs, results showed that MODIS data were not suitable due to texture artifacts, and we therefore recommended the use of Sentinel-2 texture indices calculated at different seasons instead. However, we found that the joint use of temporal variables (MODIS) and S2 texture variables does not provide significant added value compared to the use of temporal variables alone.

Despite the initial assumption that the spectral and temporal RS-ELVs indirectly include the three-dimensional environment, the inconsistency observed for landscape class 2 (composed of landscape units located on both the western and eastern sides of the forest corridor) could mean that the RS-ELVs we propose are not discriminating enough. One way of correcting this type of problem would be to use a DEM, as originally considered (Table 1). Topographic data derived from Earth observation are available worldwide. Therefore we recommend, for regions with high relief, to include topographic variables in the landscape mapping process.

### 5.1.2. The algorithms

The approach is based on GEOBIA with segmentation and clustering steps. As shown in this paper, GEOBIA allows to overcome the pixel view and provide image objects that make sense, opening a new dimension in rule-based automated image analysis (S. Lang et al. 2019). However, GEOBIA implementation is not straightforward, and different technical options are available for each step.

In this paper, we have used algorithms (region merging and k-means) that are known to be robust and adapted to our goal. We are aware of the inherent instability of some of these algorithms, and their sensitivity to parameterization. The k-means algorithm is known to be highly dependent on the initial cluster centers (Peña, Lozano, and Larrañaga 1999), which can lead to different clustering results (Li and Wu 2012), but it is rather simple to implement and efficient. Regarding the sensitivity of the algorithms to parameterization, a statistical optimization of the parameters was proposed (number of classes with an unsupervised method that maximizes the homogeneity within segments and heterogeneity between adjacent segments, and number of clusters with Elbow's method) to guide the landscape mapping process.
<figure><img src="image_43.png" /></figure>

## 5.2. The thematic consistency of the radiometric landscape map

Drawing boundaries in a continuous environment is by nature artificial, and the resulting landscape map may appear unsatisfactory. However, Simensen, Halvorsen, and Erikstad (2018) emphasize the need to focus on the accuracy and reliability of the landscape maps. If a landscape map cannot be validated in the classical sense of the term, it can at least be assessed by qualifying its global thematic consistency. The consistency of the final radiometric landscape map was checked using three different and independent sources of data : a dataset composed of spectral and textural RS-ELVs derived from two seasonal Sentinel-2 images, a land cover map, and ground observations. These three sources of data show consistency in the discrimination and characterization of the six landscape types in the study area.

From north to south, the forest corridor follows the highest relief of the island, dividing Madagascar's landscape into two groups. On its eastern side, which receives the most rainfall, the vegetation is tropical and humid. The western slopes of the plateaus are much drier. The climatic organization (east-west), combined with the distribution of urban centers along the RN7, explains the longitudinal organization of the local landscapes. The radiometric landscape map clearly shows this spatial structure in strips on either side of the forest corridor. The forest corridor (landscape class 4) is characterized by high vegetation indices throughout the year, low texture values, and consists almost entirely of tree cover. The strip to the east of the corridor consists mainly of dense secondary vegetation (shrubland) with a few cultivated plots and low texture values (class 5). In the coastal region of the Indian Ocean, where the population density is high, appears a new landscape (class 6) with high texture values, typical of mosaics of crop and natural vegetation. To the west, on the high plateaux, the savannah landscape has a low vegetation productivity and a high spatial heterogeneity marked by cultivation, depending on the distance from the north-south road axis and the presence of towns (classes 1 and 2). These results show that our approach can identify and delineate the gradient and major regions of a complex region in a way that satisfies experts and is consistent with existing land cover products.

## 6. Conclusions

In this paper we have proposed an ambitious and original method that lays the foundations for a new approach to landscape mapping that overcomes some of the major limitations of the existing biophysical methods. The method is based on an original concept that does not depend on measured (climate data) or interpreted products (such as land cover maps), but only on remote sensing data.

The method is statistically sound, and can be replicated at different scales and locations, at no cost.

Applied to a site in central Madagascar, the approach resulted in a radiometric landscape map that highlights typical characteristic areas, such as ecological corridors, contrasting landscapes from east to west as a result of different climatic conditions and agricultural practices. This map is in agreement with the ground observations. However, it is important to keep in mind that there is no single landscape map, but rather a variety of maps depending on the objective being pursued. The method proposed in this paper, like all existing methods, cannot take into account all the dimensions of the landscape without important trade-offs. Despite all the possible ways still to be explored, we believe that the method produces a generic map that has the merit of being as objective as possible and not dependent on existing thematic layers.

We hope that the method presented in this paper serves as a first step toward a generic framework for producing radiometric landscapes, and that it could be applied and tested in different parts of the world. The data sets, the satellite image pre-processing, the segmentation and clustering algorithms all have an impact on the final map. Further research is thus needed to stabilize the automatic Earth observation data processing chain, to analyze the impact of adding other RS-ELVs to the process, and to evaluate the approach in different application domains (biodiversity conservation, land planning, agro-ecological zoning, etc.) using ground data and expertise.

## Acknowledgements

The authors would like to thank Samuel Razanaka (CNRE), Solofo Rakotondraompiana (IOGA), Ibrahim Rijasoa Ravonjimalala (CNRE), Avisoatolona Andrianarivo (université d'Antananarivo), and Eric Delaitre (IRD) for their constructive contribution during the project.

## Disclosure statement

No potential conflict of interest was reported by the author(s).

## Funding

This work was supported by the French research institutes IRD, CIRAD and INRAE; the Centre National d'Études Spatiales (TOSCA CES-PAYSAGE project); and the International Joint Laboratory based in Madagascar (LMI Paysages).

## Notes on contributors

Louise Lemettais received the double master's degree in risk geography from the University of Grenoble, France, in 2020,
<figure><img src="image_46.png" /></figure>

## 5.2. The thematic consistency of the radiometric landscape map

Drawing boundaries in a continuous environment is by nature artificial, and the resulting landscape map may appear unsatisfactory. However, Simensen, Halvorsen, and Erikstad (2018) emphasize the need to focus on the accuracy and reliability of the landscape maps. If a landscape map cannot be validated in the classical sense of the term, it can at least be assessed by qualifying its global thematic consistency. The consistency of the final radiometric landscape map was checked using three different and independent sources of data : a dataset composed of spectral and textural RS-ELVs derived from two seasonal Sentinel-2 images, a land cover map, and ground observations. These three sources of data show consistency in the discrimination and characterization of the six landscape types in the study area.

From north to south, the forest corridor follows the highest relief of the island, dividing Madagascar's landscape into two groups. On its eastern side, which receives the most rainfall, the vegetation is tropical and humid. The western slopes of the plateaus are much drier. The climatic organization (east-west), combined with the distribution of urban centers along the RN7, explains the longitudinal organization of the local landscapes. The radiometric landscape map clearly shows this spatial structure in strips on either side of the forest corridor. The forest corridor (landscape class 4) is characterized by high vegetation indices throughout the year, low texture values, and consists almost entirely of tree cover. The strip to the east of the corridor consists mainly of dense secondary vegetation (shrubland) with a few cultivated plots and low texture values (class 5). In the coastal region of the Indian Ocean, where the population density is high, appears a new landscape (class 6) with high texture values, typical of mosaics of crop and natural vegetation. To the west, on the high plateaux, the savannah landscape has a low vegetation productivity and a high spatial heterogeneity marked by cultivation, depending on the distance from the north-south road axis and the presence of towns (classes 1 and 2). These results show that our approach can identify and delineate the gradient and major regions of a complex region in a way that satisfies experts and is consistent with existing land cover products.

## 6. Conclusions

In this paper we have proposed an ambitious and original method that lays the foundations for a new approach to landscape mapping that overcomes some of the major limitations of the existing biophysical methods. The method is based on an original concept that does not depend on measured (climate data) or interpreted products (such as land cover maps), but only on remote sensing data.

The method is statistically sound, and can be replicated at different scales and locations, at no cost.

Applied to a site in central Madagascar, the approach resulted in a radiometric landscape map that highlights typical characteristic areas, such as ecological corridors, contrasting landscapes from east to west as a result of different climatic conditions and agricultural practices. This map is in agreement with the ground observations. However, it is important to keep in mind that there is no single landscape map, but rather a variety of maps depending on the objective being pursued. The method proposed in this paper, like all existing methods, cannot take into account all the dimensions of the landscape without important trade-offs. Despite all the possible ways still to be explored, we believe that the method produces a generic map that has the merit of being as objective as possible and not dependent on existing thematic layers.

We hope that the method presented in this paper serves as a first step toward a generic framework for producing radiometric landscapes, and that it could be applied and tested in different parts of the world. The data sets, the satellite image pre-processing, the segmentation and clustering algorithms all have an impact on the final map. Further research is thus needed to stabilize the automatic Earth observation data processing chain, to analyze the impact of adding other RS-ELVs to the process, and to evaluate the approach in different application domains (biodiversity conservation, land planning, agro-ecological zoning, etc.) using ground data and expertise.

## Acknowledgements

The authors would like to thank Samuel Razanaka (CNRE), Solofo Rakotondraompiana (IOGA), Ibrahim Rijasoa Ravonjimalala (CNRE), Avisoatolona Andrianarivo (université d'Antananarivo), and Eric Delaitre (IRD) for their constructive contribution during the project.

## Disclosure statement

No potential conflict of interest was reported by the author(s).

## Funding

This work was supported by the French research institutes IRD, CIRAD and INRAE; the Centre National d'Études Spatiales (TOSCA CES-PAYSAGE project); and the International Joint Laboratory based in Madagascar (LMI Paysages).

## Notes on contributors

Louise Lemettais received the double master's degree in risk geography from the University of Grenoble, France, in 2020,
<figure><img src="image_44.png" /></figure>

and in geomatic sciences for the environment and planning from the University of Toulouse, France, in 2021. She is currently pursuing a PhD with CNES and the University of Lille at the Laboratoire d'océanologie et de géosciences (LOG), France. Her research focuses both on the development of satellite methods for assessing vegetation dynamics in relation to climate change, and on the development of methods for mapping landscapes using satellite time series.

**Samuel Alleaume**, after studying ecology, focused his skills on geography and remote sensing, obtaining a Master's degree in Geography at the Université de Rennes 2 in 1994, followed by a Master's degree in Geography in 1999 at the Université du Québec à Montréal (Canada). Subsequently, he was employed as a study engineer in various research organizations : the University of Virginia (USA), CNRS, CEMAGREF. He now works at INRAE in Montpellier as a research engineer. His main research activities involve developing and implementing biodiversity mapping and assessment methodologies based on remote sensing and modeling.

**Sandra Luque** Sandra Luque (PhD), Landscape Ecologist, graduated from Rutgers University (USA) and Habilitation (HDR) in France. She is currently Research Director Exceptional Class at INRAE. Former NASA EOS Fellow, Global Change Program (USA), she currently represents CNES (France) for CEOS (https://ceos.org/) in the Ecosystems Task Force. She works on temporal and spatial change, especially in forested landscapes, biodiversity indicators, habitat modeling and has recently published on conservation values, sustainability science and ecosystem services.

**Anne-Élisabeth Laques** (Ph.D.) is Director of Research in geography at IRD, a French public multidisciplinary research institute. She is currently based in Madagascar, where she co-directs the LMI (Laboratoire Mixte International) on behalf of the IRD. Her research primarily focuses on landscape analysis and the utilization of satellite imagery to construct socio-environmental indicators. Her objective is to generate cartographic assessments that aid in modeling future scenarios, aimed at assessing the spatial effects of public policies. Additionally, she fosters collaborative research programs with partners from Brazil and Madagascar.

**Yonas Alim** is a computer scientist, who graduated from the Master's program in BDMA (Big Data Management and Analytics) at the University of Tours, Blois campus in 2022.

**Laurent Demagistri** obtained a Master's degree in Computer Vision from the University of Nice, France, in 1992. He was initially trained in applied mathematics. He started his career at a private company, where he worked on research projects that used satellite imagery and participated in the development of a commercial image processing software. Since 2002, he has been working at IRD, a French public multidisciplinary research institute that collaborates with countries in the global south. He is currently a Research Engineer focusing on methodological research for environmental characterization using satellite imagery workflows. More specifically in recent projects on urban areas.

**Agnès Bégué** received the Ph.D. degree in physics from the University Paris VII, France, in 1991. She has initial training in agricultural engineering. After Postdoctoral Research experiences at the University of Maryland (USA) and at CNES (France), she got a permanent position in 1995 with

CIRAD, a French scientific organization specialized in development-oriented agricultural research for the tropics and subtropics, Montpellier, France. Her research interests include the development of satellite data-based methods to characterize land use, in particular the African agricultural systems, through landscape analysis.

## ORCID

Louise Lemettais http://orcid.org/0000-0002-9510-507X

Samuel Alleaume http://orcid.org/0000-0002-9200-8338

Sandra Luque http://orcid.org/0000-0002-4002-3974

Anne-Élisabeth Laques http://orcid.org/0000-0003-4469-270X

Laurent Demagistri http://orcid.org/0000-0003-1910-200X

Agnès Bégué http://orcid.org/0000-0002-9289-1052

## Data availability statement

Data that support the results and analyses presented in this paper are freely available online. Sentinel-2 data are available via the Theia platform (https://catalogue.theia land.fr); MODIS products can be downloaded from the NASA's AppEEARS platform (https://appeears.earthdata cloud.nasa.gov); The ESA WorldCover 2020 map is available at https://worldcover2020.esa.int/download; UAV photographs panoramas for the 2022 campaign are available at https://data.inrae.fr/citation?persistentId=doi:10. 57745/DBGRDI. The landscape maps produced in this study are available from the corresponding author on request.

## References

Alcántara Manzanares, J., and J. M. Muñoz Álvarez. 2015. "Landscape Classification of Huelva (Spain): An Objective Method of Identification and Characterization." *Estudios Geográficos* 76 (279): 447-471. https://doi.org/10.3989/estgeogr.201516.

Alleaume, S. 2022. *Panoramas Paysages Madagascar Acquis Par Drone*. Portail Data INRAE. https://doi.org/10.57745/DBGRDI.

Alleaume, S., P. Dusseux, V. Thierion, L. Commagnac, S. Laventure, M. Lang, J. Féret, L. Hubert-Moy, S. Luque, and P. Vihervaara. 2018. "A Generic Remote Sensing Approach to Derive Operational Essential Biodiversity Variables (EBVs) for Conservation Planning." Edited by Petteri Vihervaara. *Methods in Ecology and Evolution* 9 (8): 1822-1836. https://doi.org/10.1111/2041-210X.13033.

Antrop, M. 2005. "Why Landscapes of the Past are Important for the Future." *Landscape and Urban Planning* 70 (1-2): 21-34. https://doi.org/10.1016/j.land urbplan.2003.10.002.

Antrop, M., and V. Van Eetvelde. 2017. "Approaches in Landscape Research." In *Landscape Perspectives, by Marc Antrop and Veerle Van Eetvelde*, Vol. 23, 61-80. Landscape Series. Dordrecht: Springer Netherlands. https://doi.org/10.1007/978-94-024-1183-6_4.

Baatz, M., and A. Schape. 2000. "Multiresolution Segmentation: An Optimization Approach for High Quality Multi-Scale Image Segmentation." In *Angewandte Geographische Informations-Verarbeitung*
<figure><img src="image_47.png" /></figure>

and in geomatic sciences for the environment and planning from the University of Toulouse, France, in 2021. She is currently pursuing a PhD with CNES and the University of Lille at the Laboratoire d'océanologie et de géosciences (LOG), France. Her research focuses both on the development of satellite methods for assessing vegetation dynamics in relation to climate change, and on the development of methods for mapping landscapes using satellite time series.

Samuel Alleaume, after studying ecology, focused his skills on geography and remote sensing, obtaining a Master’s degree in Geography at the Université de Rennes 2 in 1994, followed by a Master’s degree in Geography in 1999 at the Université du Québec à Montréal (Canada). Subsequently, he was employed as a study engineer in various research organizations : the University of Virginia (USA), CNRS, CEMAGREF. He now works at INRAE in Montpellier as a research engineer. His main research activities involve developing and implementing biodiversity mapping and assessment methodologies based on remote sensing and modeling.

Sandra Luque Sandra Luque (PhD), Landscape Ecologist, graduated from Rutgers University (USA) and Habilitation (HDR) in France. She is currently Research Director Exceptional Class at INRAE. Former NASA EOS Fellow, Global Change Program (USA), she currently represents CNES (France) for CEOS (https://ceos.org/) in the Ecosystems Task Force. She works on temporal and spatial change, especially in forested landscapes, biodiversity indicators, habitat modeling and has recently published on conservation values, sustainability science and ecosystem services.

Anne-Élisabeth Laques (Ph.D.) is Director of Research in geography at IRD, a French public multidisciplinary research institute. She is currently based in Madagascar, where she co-directs the LMI (Laboratoire Mixte International) on behalf of the IRD. Her research primarily focuses on landscape analysis and the utilization of satellite imagery to construct socio-environmental indicators. Her objective is to generate cartographic assessments that aid in modeling future scenarios, aimed at assessing the spatial effects of public policies. Additionally, she fosters collaborative research programs with partners from Brazil and Madagascar.

Yonas Alim is a computer scientist, who graduated from the Master's program in BDMA (Big Data Management and Analytics) at the University of Tours, Blois campus in 2022.

Laurent Demagistri obtained a Master's degree in Computer Vision from the University of Nice, France, in 1992. He was initially trained in applied mathematics. He started his career at a private company, where he worked on research projects that used satellite imagery and participated in the development of a commercial image processing software. Since 2002, he has been working at IRD, a French public multidisciplinary research institute that collaborates with countries in the global south. He is currently a Research Engineer focusing on methodological research for environmental characterization using satellite imagery workflows. More specifically in recent projects on urban areas.

Agnès Bégué received the Ph.D. degree in physics from the University Paris VII, France, in 1991. She has initial training in agricultural engineering. After Postdoctoral Research experiences at the University of Maryland (USA) and at CNES (France), she got a permanent position in 1995 with

CIRAD, a French scientific organization specialized in development-oriented agricultural research for the tropics and subtropics, Montpellier, France. Her research interests include the development of satellite data-based methods to characterize land use, in particular the African agricultural systems, through landscape analysis.

## ORCID

Louise Lemettais http://orcid.org/0000-0002-9510-507X

Samuel Alleaume http://orcid.org/0000-0002-9200-8338

Sandra Luque http://orcid.org/0000-0002-4002-3974

Anne-Élisabeth Laques http://orcid.org/0000-0003-4469-270X

Laurent Demagistri http://orcid.org/0000-0003-1910-200X

Agnès Bégué http://orcid.org/0000-0002-9289-1052

## Data availability statement

Data that support the results and analyses presented in this paper are freely available online. Sentinel-2 data are available via the Theia platform (https://catalogue.theia land.fr); MODIS products can be download from the NASA's AppEEARS platform (https://appeears.earthdata cloud.nasa.gov); The ESA WorldCover 2020 map is avail- able at https://worldcover2020.esa.int/download; UAV photographs panoramas for the 2022 campaign are avail- able at https://data.inrae.fr/citation?persistentId=doi:10. 57745/DBGRDI. The landscape maps produced in this study are available from the corresponding author on request.

## References

Alcántara Manzanares, J., and J. M. Muñoz Álvarez. 2015. "Landscape Classification of Huelva (Spain): An Objective Method of Identification and Characterization." *Estudios Geográficos* 76 (279): 447-471. https://doi.org/10.3989/estgeogr.201516.

Alleaume, S. 2022. *Panoramas Paysages Madagascar Acquis Par Drone*. Portail Data INRAE. https://doi.org/10.57745/DBGRDI.

Alleaume, S., P. Dusseux, V. Thierion, L. Commagnac, S. Laventure, M. Lang, J. Féret, L. Hubert-Moy, S. Luque, and P. Vihervaara. 2018. "A Generic Remote Sensing Approach to Derive Operational Essential Biodiversity Variables (EBVs) for Conservation Planning." Edited by Petteri Vihervaara. *Methods in Ecology and Evolution* 9 (8): 1822-1836. https://doi.org/10.1111/2041-210X.13033.

Antrop, M. 2005. "Why Landscapes of the Past are Important for the Future." *Landscape and Urban Planning* 70 (1-2): 21-34. https://doi.org/10.1016/j.land urbplan.2003.10.002.

Antrop, M., and V. Van Eetvelde. 2017. "Approaches in Landscape Research." In Landscape Perspectives, by *Marc Antrop and Veerle Van Eetvelde*, Vol. 23, 61-80. Landscape Series. Dordrecht: Springer Netherlands. https://doi.org/10.1007/978-94-024-1183-6_4.

Baatz, M., and A. Schape. 2000. "Multiresolution Segmentation: An Optimization Approach for High Quality Multi-Scale Image Segmentation." In *Angewandte Geographische Informations-Verarbeitung*
<figure><img src="image_45.png" /></figure>

XII, edited by Strobl, J., Blaschke, T., Griesebner, G., 12-23. Karlsruhe, Germany: Wichmann Verlag.

Bégué, A., D. Arvor, C. Lelong, E. Vintrou, and M. Simoes. 2015. "Agricultural Systems Studies Using Remote Sensing." In *Remote Sensing Handbook, Volume 2: Land Resources Monitoring, Modeling, and Mapping with Remote Sensing, Edition*, edited by Prasad Thenkabail, Francis and Taylor Group, Publisher: CRC Press.

Bellón, B., A. Bégué, D. Lo Seen, C. De Almeida, and M. Simões. 2017. "A Remote Sensing Approach for Regional-Scale Mapping of Agricultural Land-Use Systems Based on NDVI Time Series." *Remote Sensing* 9 (6): 600. https://doi.org/10.3390/rs9060600.

Bellón, B., A. Bégué, D. Lo Seen, V. Lebourgeois, B. Antônio Evangelista, M. Simões, and R. Peçanha Demonte Ferraz. 2018. "Improved Regional-Scale Brazilian Cropping Systems' Mapping Based on a Semi-Automatic Object-Based Clustering Approach." *International Journal of Applied Earth Observation and Geoinformation* 68: 127-138. https://doi.org/10.1016/j.jag.2018.01.019.

Bisquert, M., A. Bégué, and M. Deshayes. 2015. "Object-Based Delineation of Homogeneous Landscape Units at Regional Scale Based on MODIS Time Series." *International Journal of Applied Earth Observation and Geoinformation* 37 (May): 72-82. https://doi.org/10.1016/j.jag.2014.10.004.

Bisquert, M., A. Bégué, M. Deshayes, and D. Ducrot. 2017. "Environmental Evaluation of MODIS-Derived Land Units." *GIScience & Remote Sensing* 54 (1): 64-77. https://doi.org/10.1080/15481603.2016.1256861.

Blaschke, T. 2010. "Object Based Image Analysis for Remote Sensing." *ISPRS Journal of Photogrammetry and Remote Sensing* 65 (1): 2-16. https://doi.org/10.1016/j.isprsjprs.2009.06.004.

Blaschke, T., G. J. Hay, M. Kelly, S. Lang, P. Hofmann, E. Addink, R. Queiroz Feitosa, et al. 2014. "Geographic Object-Based Image Analysis – Towards a New Paradigm." *ISPRS Journal of Photogrammetry and Remote Sensing* 87 (January): 180-191. https://doi.org/10.1016/j.isprsjprs.2013.09.014.

Böck, S., M. Immitzer, and C. Atzberger. 2017. "On the Objectivity of the Objective Function – Problems with Unsupervised Segmentation Evaluation Based on Global Score and a Possible Remedy." *Remote Sensing* 9 (8): 769. https://doi.org/10.3390/rs9080769.

Bourget, É., and L. le Dû-Blayo. 2010. "Définition d'unités Paysagères Par Télédétection En Bretagne : Méthodes et Critiques." *Norois*, (216): 69-83. https://doi.org/10.4000/norois.3399.

Bunce, R. G. H., C. J. Barr, R. T. Clarke, D. C. Howard, and A. M. J. Lane. 1996. "Land Classification for Strategic Ecological Survey." *Journal of Environmental Management* 47 (1): 37-60. https://doi.org/10.1006/jema.1996.0034.

Cabral, A. I. R., C. Saito, H. Pereira, and A. Elisabeth Laques. 2018. "Deforestation Pattern Dynamics in Protected Areas of the Brazilian Legal Amazon Using Remote Sensing Data." *Applied Geography* 100 (November): 101-115. https://doi.org/10.1016/j.apgeog.2018.10.003.

Cano, E., J.-P. Denux, M. Bisquert, L. Hubert-Moy, and V. Chéret. 2017. "Improved Forest-Cover Mapping Based on MODIS Time Series and Landscape Stratification." *International Journal of Remote Sensing* 38 (7): 1865-1888. https://doi.org/10.1080/01431161.2017.1280635.

Chen, J., J. Per, M. Tamura, G. Zhihui, B. Matsushita, and L. Eklundh. 2004. "A Simple Method for Reconstructing

a High-Quality NDVI Time-Series Data Set Based on the Savitzky-Golay Filter." *Remote Sensing of Environment* 91 (3-4): 332-344. https://doi.org/10.1016/j.rse.2004.03.014.

Chen, Y., D. Lu, E. Moran, M. Batistella, L. Dutra, I. Sanches, R. Silva, J.-F. Huang, A. Luiz, and M. A. F. de Oliveira. 2018. "Mapping Croplands, Cropping Patterns, and Crop Types Using MODIS Time-Series Data." *International Journal of Applied Earth Observation and Geoinformation* 69: 133-147. https://doi.org/10.1016/j.jag.2018.03.005.

Clinton, N., A. Holt, J. Scarborough, L. Yan, and P. Gong. 2010. "Accuracy Assessment Measures for Object-Based Image Segmentation Goodness." *Photogrammetric Engineering & Remote Sensing* 76 (3): 289-299. https://doi.org/10.14358/PERS.76.3.289.

Couteron, P., N. Barbier, and D. Gautier. 2006. "Textural Ordination Based on Fourier Spectral Decomposition: A Method to Analyze and Compare Landscape Patterns." *Landscape Ecology* 21 (4): 555-567. https://doi.org/10.1007/s10980-005-2166-6.

Cullum, C., K. H. Rogers, G. Brierley, and E. T. F. Witkowski. 2016. "Ecological Classification and Mapping for Landscape Management and Science: Foundations for the Description of Patterns and Processes." *Progress in Physical Geography: Earth and Environment* 40 (1): 38-65. https://doi.org/10.1177/0309133315611573.

Deng, J. S., K. Wang, Y. H. Deng, and G. J. Qi. 2008. "PCA-Based Land-Use Change Detection and Analysis Using Multitemporal and Multisensor Satellite Data." *International Journal of Remote Sensing* 29 (16): 4823-4838. https://doi.org/10.1080/01431160801950162.

Duda, T., and M. Canty. 2002. "Unsupervised Classification of Satellite Imagery: Choosing a Good Algorithm." *International Journal of Remote Sensing* 23 (11): 2193-2212. https://doi.org/10.1080/01431160110078467.

Dupuy, S., L. Defrise, V. Lebourgeois, R. Gaetano, P. Burnod, and J.-P. Tonneau. 2020. "Analyzing Urban Agriculture's Contribution to a Southern City's Resilience Through Land Cover Mapping: The Case of Antananarivo, Capital of Madagascar." *Remote Sensing* 12 (12): 1962. https://doi.org/10.3390/rs12121962.

Erikstad, L., L. Andre Uttakleiv, and R. Halvorsen. 2015. "Characterisation and Mapping of Landscape Types, a Case Study from Norway". *Belgeo*, (3). https://doi.org/10.4000/belgeo.17412.

Espindola, G. M., G. Camara, I. A. Reis, L. S. Bins, and A. M. Monteiro. 2006. "Parameter Selection for Region-Growing Image Segmentation Algorithms Using Spatial Autocorrelation." *International Journal of Remote Sensing* 27 (14): 3035-3040. https://doi.org/10.1080/01431160600617194.

Ferreira, T. R., B. B. Da Silva, M. S. B. De Moura, A. Verhoef, and R. L. B. Nóbrega. 2020. "The Use of Remote Sensing for Reliable Estimation of Net Radiation and Its Components: A Case Study for Contrasting Land Covers in an Agricultural Hotspot of the Brazilian Semiarid Region." *Agricultural and Forest Meteorology* 291 (September): 108052. https://doi.org/10.1016/j.agrfor met.2020.108052.

Forman, R. T. T., and M. Godron. 1986. *Landscape Ecology*. New York: John Wiley&Sons Inc.

Gao, B.-C. 1996. "NDWI – A Normalized Difference Water Index for Remote Sensing of Vegetation Liquid Water from Space." *Remote Sensing of Environment* 58 (3): 257-266. https://doi.org/10.1016/S0034-4257(96)00067-3.
<figure><img src="image_48.png" /></figure>

XII, edited by Strobl, J., Blaschke, T., Griesebner, G., 12-23. Karlsruhe, Germany: Wichmann Verlag.

Bégué, A., D. Arvor, C. Lelong, E. Vintrou, and M. Simoes. 2015. "Agricultural Systems Studies Using Remote Sensing." In *Remote Sensing Handbook, Volume 2: Land Resources Monitoring, Modeling, and Mapping with Remote Sensing, Edition*, edited by Prasad Thenkabail, Francis and Taylor Group, Publisher: CRC Press.

Bellón, B., A. Bégué, D. Lo Seen, C. De Almeida, and M. Simões. 2017. "A Remote Sensing Approach for Regional-Scale Mapping of Agricultural Land-Use Systems Based on NDVI Time Series." *Remote Sensing* 9 (6): 600. https://doi.org/10.3390/rs9060600.

Bellón, B., A. Bégué, D. Lo Seen, V. Lebourgeois, B. Antônio Evangelista, M. Simões, and R. Peçanha Demonte Ferraz. 2018. "Improved Regional-Scale Brazilian Cropping Systems' Mapping Based on a Semi-Automatic Object-Based Clustering Approach." *International Journal of Applied Earth Observation and Geoinformation* 68: 127-138. https://doi.org/10.1016/j.jag.2018.01.019.

Bisquert, M., A. Bégué, and M. Deshayes. 2015. "Object-Based Delineation of Homogeneous Landscape Units at Regional Scale Based on MODIS Time Series." *International Journal of Applied Earth Observation and Geoinformation* 37 (May): 72-82. https://doi.org/10.1016/j.jag.2014.10.004.

Bisquert, M., A. Bégué, M. Deshayes, and D. Ducrot. 2017. "Environmental Evaluation of MODIS-Derived Land Units." *GIScience & Remote Sensing* 54 (1): 64-77. https://doi.org/10.1080/15481603.2016.1256861.

Blaschke, T. 2010. "Object Based Image Analysis for Remote Sensing." *ISPRS Journal of Photogrammetry and Remote Sensing* 65 (1): 2-16. https://doi.org/10.1016/j.isprsjprs.2009.06.004.

Blaschke, T., G. J. Hay, M. Kelly, S. Lang, P. Hofmann, E. Addink, R. Queiroz Feitosa, et al. 2014. "Geographic Object-Based Image Analysis - Towards a New Paradigm." *ISPRS Journal of Photogrammetry and Remote Sensing* 87 (January): 180-191. https://doi.org/10.1016/j.isprsjprs.2013.09.014.

Böck, S., M. Immitzer, and C. Atzberger. 2017. "On the Objectivity of the Objective Function – Problems with Unsupervised Segmentation Evaluation Based on Global Score and a Possible Remedy." *Remote Sensing* 9 (8): 769. https://doi.org/10.3390/rs9080769.

Bourget, É., and L. le Dû-Blayo. 2010. "Définition d'unités Paysagères Par Télédétection En Bretagne : Méthodes et Critiques." *Norois*, (216): 69-83. https://doi.org/10.4000/norois.3399.

Bunce, R. G. H., C. J. Barr, R. T. Clarke, D. C. Howard, and A. M. J. Lane. 1996. "Land Classification for Strategic Ecological Survey." *Journal of Environmental Management* 47 (1): 37-60. https://doi.org/10.1006/jema.1996.0034.

Cabral, A. I. R., C. Saito, H. Pereira, and A. Elisabeth Laques. 2018. "Deforestation Pattern Dynamics in Protected Areas of the Brazilian Legal Amazon Using Remote Sensing Data." *Applied Geography* 100 (November): 101-115. https://doi.org/10.1016/j.apgeog.2018.10.003.

Cano, E., J.-P. Denux, M. Bisquert, L. Hubert-Moy, and V. Chéret. 2017. "Improved Forest-Cover Mapping Based on MODIS Time Series and Landscape Stratification." *International Journal of Remote Sensing* 38 (7): 1865-1888. https://doi.org/10.1080/01431161.2017.1280635.

Chen, J., J. Per, M. Tamura, G. Zhihui, B. Matsushita, and L. Eklundh. 2004. "A Simple Method for Reconstructing

a High-Quality NDVI Time-Series Data Set Based on the Savitzky-Golay Filter." *Remote Sensing of Environment* 91 (3-4): 332-344. https://doi.org/10.1016/j.rse.2004.03.014.

Chen, Y., D. Lu, E. Moran, M. Batistella, L. Dutra, I. Sanches, R. Silva, J.-F. Huang, A. Luiz, and M. A. F. de Oliveira. 2018. "Mapping Croplands, Cropping Patterns, and Crop Types Using MODIS Time-Series Data." *International Journal of Applied Earth Observation and Geoinformation* 69: 133-147. https://doi.org/10.1016/j.jag.2018.03.005.

Clinton, N., A. Holt, J. Scarborough, L. Yan, and P. Gong. 2010. "Accuracy Assessment Measures for Object-Based Image Segmentation Goodness." *Photogrammetric Engineering & Remote Sensing* 76 (3): 289-299. https://doi.org/10.14358/PERS.76.3.289.

Couteron, P., N. Barbier, and D. Gautier. 2006. "Textural Ordination Based on Fourier Spectral Decomposition: A Method to Analyze and Compare Landscape Patterns." *Landscape Ecology* 21 (4): 555-567. https://doi.org/10.1007/s10980-005-2166-6.

Cullum, C., K. H. Rogers, G. Brierley, and E. T. F. Witkowski. 2016. "Ecological Classification and Mapping for Landscape Management and Science: Foundations for the Description of Patterns and Processes." *Progress in Physical Geography: Earth and Environment* 40 (1): 38-65. https://doi.org/10.1177/0309133315611573.

Deng, J. S., K. Wang, Y. H. Deng, and G. J. Qi. 2008. "PCA-Based Land-Use Change Detection and Analysis Using Multitemporal and Multisensor Satellite Data." *International Journal of Remote Sensing* 29 (16): 4823-4838. https://doi.org/10.1080/01431160801950162.

Duda, T., and M. Canty. 2002. "Unsupervised Classification of Satellite Imagery: Choosing a Good Algorithm." *International Journal of Remote Sensing* 23 (11): 2193-2212. https://doi.org/10.1080/01431160110078467.

Dupuy, S., L. Defrise, V. Lebourgeois, R. Gaetano, P. Burnod, and J.-P. Tonneau. 2020. "Analyzing Urban Agriculture's Contribution to a Southern City's Resilience Through Land Cover Mapping: The Case of Antananarivo, Capital of Madagascar." *Remote Sensing* 12 (12): 1962. https://doi.org/10.3390/rs12121962.

Erikstad, L., L. Andre Uttakleiv, and R. Halvorsen. 2015. "Characterisation and Mapping of Landscape Types, a Case Study from Norway". *Belgeo*, (3). https://doi.org/10.4000/belgeo.17412.

Espindola, G. M., G. Camara, I. A. Reis, L. S. Bins, and A. M. Monteiro. 2006. "Parameter Selection for Region-Growing Image Segmentation Algorithms Using Spatial Autocorrelation." *International Journal of Remote Sensing* 27 (14): 3035-3040. https://doi.org/10.1080/01431160600617194.

Ferreira, T. R., B. B. Da Silva, M. S. B. De Moura, A. Verhoef, and R. L. B. Nóbrega. 2020. "The Use of Remote Sensing for Reliable Estimation of Net Radiation and Its Components: A Case Study for Contrasting Land Covers in an Agricultural Hotspot of the Brazilian Semiarid Region." *Agricultural and Forest Meteorology* 291 (September): 108052. https://doi.org/10.1016/j.agrfor met.2020.108052.

Forman, R. T. T., and M. Godron. 1986. *Landscape Ecology*. New York: John Wiley&Sons Inc.

Gao, B.-C. 1996. "NDWI - A Normalized Difference Water Index for Remote Sensing of Vegetation Liquid Water from Space." *Remote Sensing of Environment* 58 (3): 257-266. https://doi.org/10.1016/S0034-4257(96)00067-3.
<figure><img src="image_46.png" /></figure>

García-Llamas, P., L. Calvo, J. Manuel Álvarez-Martínez, and S. Suárez-Seoane. 2016. "Using Remote Sensing Products to Classify Landscape. A Multi-Spatial Resolution Approach." *International Journal of Applied Earth Observation and Geoinformation* 50 (August): 95-105. https://doi.org/10.1016/j.jag.2016.03.010.

Grizonnet, M., J. Michel, V. Poughon, J. Inglada, M. Savinaud, and R. Cresson. 2017. "Orfeo ToolBox: Open Source Processing of Remote Sensing Images." *Open Geospatial Data, Software and Standards* 2 (1): 15. https://doi.org/10.1186/s40965-017-0031-6.

Haralick, R. M., K. Shanmugam, and I. Dinstein. 1973. "Textural Features for Image Classification." *IEEE Transactions on Systems, Man, and Cybernetics SMC-3* (6): 610-621. https://doi.org/10.1109/TSMC.1973.4309314.

Hay, G. J., and G. Castilla. 2006. "Object-Based Image Analysis: Strengths, Weaknesses, Opportunities and Threats (SWOT)." In *Commission VI, WG VI/4. file://C:/Users/begue/Downloads/Object-Based_Image_Analysis-Strengths_Weaknesses_O.pdf*.

Hazeu, G. W., M. J. Metzger, C. A. Mücher, M. Perez-Soba, C. H. Renetzeder, and E. Andersen. 2011. "European Environmental Stratifications and Typologies: An Overview." *Agriculture, Ecosystems and Environment* 142 (1-2): 29-39. https://doi.org/10.1016/j.agee.2010.01.009.

Hudak, A. T., and C. A. Wessman. 1998. "Textural Analysis of Historical Aerial Photography to Characterize Woody Plant Encroachment in South African Savanna." *Remote Sensing of Environment* 66 (3): 317-330. https://doi.org/10.1016/S0034-4257(98)00078-9.

Inglada, J. 2016. "OTB Gapfilling, a Temporal Gapfilling for Image Time Series Library." *Zenodo*. https://doi.org/10.5281/ZENODO.45572.

Johnson, B., and Z. Xie. 2011. "Unsupervised Image Segmentation Evaluation and Refinement Using a Multi-Scale Approach." *ISPRS Journal of Photogrammetry and Remote Sensing* 66 (4): 473-483. https://doi.org/10.1016/j.isprsjprs.2011.02.006.

Kalinicheva, E., J. Sublime, and M. Trocan. 2020. "Unsupervised Satellite Image Time Series Clustering Using Object-Based Approaches and 3D Convolutional Autoencoder." *Remote Sensing* 12 (11): 1816. https://doi.org/10.3390/rs12111816.

Karasov, O., M. Külvik, and I. Burdun. 2021. "Deconstructing Landscape Pattern: Applications of Remote Sensing to Physiognomic Landscape Mapping." *Geo Journal* 86 (1): 529-555. https://doi.org/10.1007/s10708-019-10058-6.

Kavzoglu, T., and H. Tonbul. 2018. "An Experimental Comparison of Multi-Resolution Segmentation, SLIC and K-Means Clustering for Object-Based Classification of VHR Imagery." *International Journal of Remote Sensing* 39 (18): 6020-6036. https://doi.org/10.1080/01431161.2018.1506592.

Knight, J. F., R. S. Lunetta, J. Ediriwickrema, and S. Khorram. 2006. "Regional Scale Land Cover Characterization Using MODIS-NDVI 250 M Multi-Temporal Imagery: A Phenology-Based Approach." *GIScience & Remote Sensing* 43 (1): 1-23. https://doi.org/10.2747/1548-1603.43.1.1.

Koç, A., and S. Yılmaz. 2020. "Landscape Character Analysis and Assessment at the Lower Basin-Scale." *Applied Geography* 125 (December): 102359. https://doi.org/10.1016/j.apgeog.2020.102359.

Kupidura, P. 2019. "The Comparison of Different Methods of Texture Analysis for Their Efficacy for Land Use Classification in Satellite Imagery." *Remote Sensing* 11 (10): 1233. https://doi.org/10.3390/rs11101233.

Lang, M., S. Alleaume, S. Luque, N. Baghdadi, and J.-B. Féret. 2018. "Monitoring and Characterizing Heterogeneous Mediterranean Landscapes with Continuous Textural Indices Based on VHSR Imagery." *Remote Sensing* 10 (6): 868. https://doi.org/10.3390/rs10060868.

Lang, S., G. J. Hay, A. Baraldi, D. Tiede, and T. Blaschke. 2019. "Geobia Achievements and Spatial Opportunities in the Era of Big Earth Observation Data." *ISPRS International Journal of Geo-Information* 8 (11): 474. https://doi.org/10.3390/ijgi8110474.

Lè, S., J. Josse, and F. Husson. 2008. "FactoMineR: An R Package for Multivariate Analysis." *Journal of Statistical Software* 25 (1). https://doi.org/10.18637/jss.v025.i01.

Leenhardt, D., Angevin, F., Biarnès, A., Colbach, N., and Mignolet, C.l. 2010. "Describing and Locating Cropping Systems on a Regional Scale - A Review." *Agronomy for Sustainable Development* 30: 131-138. https://doi.org/10.1051/agro/2009002.

Leutner, R., N. H. Benjamin, J. Schwalb-Willmann, and R. J. Hijmans. 2022. "Package 'RStoolbox.' Tools for Remote Sensing Data Analysis. https://cran.r-project.org/web/packages/RStoolbox/RStoolbox.pdf.

Levering, A., D. Marcos, and D. Tuia. 2021. "On the Relation Between Landscape Beauty and Land Cover: A Case Study in the U.K. at Sentinel-2 Resolution with Interpretable AI." *ISPRS Journal of Photogrammetry and Remote Sensing* 177 (July): 194-203. https://doi.org/10.1016/j.isprsjprs.2021.04.020.

Li, Y., and H. Wu. 2012. "A Clustering Method Based on K-Means Algorithm." *Physics Procedia* 25: 1104-1109. https://doi.org/10.1016/j.phpro.2012.03.206.

Lonjou, V., C. Desjardins, O. Hagolle, B. Petrucci, T. Tremas, M. Dejus, A. Makarau, and S. Auer. 2016. "MACCS-ATCOR Joint Algorithm (MAJA)." In edited by A. Comerón, E. I. Kassianov, and K. Schäfer, 1000107. Edinburgh, United Kingdom. https://doi.org/10.1117/12.2240935.

Luque, S., N. Pettorelli, P. Vihervaara, M. Wegmann, and J. Vamosi. 2018. "Improving Biodiversity Monitoring Using Satellite Remote Sensing to Provide Solutions Towards the 2020 Conservation Targets." Edited by Jana Vamosi. *Methods in Ecology and Evolution* 9 (8): 1784-1786. https://doi.org/10.1111/2041-210X.13057.

Ma, L., M. Li, X. Ma, L. Cheng, P. Du, and Y. Liu. 2017. "A Review of Supervised Object-Based Land-Cover Image Classification." *Isprsjprs Journal of Photogrammetry & Remote Sensing* 130 (August): 277-293. https://doi.org/10.1016/j.isprsjprs.2017.06.001.

Meinel, G., and M. Neubert. 2004. "A Comparison of Segmentation Programs for High Resolution Remote Sensing Data." In *Geo-Imagery Bridging Continents. XXth ISPRS Congress, Istanbul, Turkey, 12-23 July (ISPRS International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences; vol. XXXV, Part B4)*, 1097-1102.

Mücher, C. A., J. A. Klijn, D. M. Wascher, and J. H. J. Schaminée. 2010. "A New European Landscape Classification (LANMAP): A Transparent, Flexible and User-Oriented Methodology to Distinguish Landscapes." *Ecological Indicators* 10 (1): 87-103. https://doi.org/10.1016/j.ecolind.2009.03.018.
<figure><img src="image_49.png" /></figure>

García-Llamas, P., L. Calvo, J. Manuel Álvarez-Martínez, and S. Suárez-Seoane. 2016. "Using Remote Sensing Products to Classify Landscape. A Multi-Spatial Resolution Approach." *International Journal of Applied Earth Observation and Geoinformation* 50 (August): 95-105. https://doi.org/10.1016/j.jag.2016.03.010.

Grizonnet, M., J. Michel, V. Poughon, J. Inglada, M. Savinaud, and R. Cresson. 2017. "Orfeo ToolBox: Open Source Processing of Remote Sensing Images." *Open Geospatial Data, Software and Standards* 2 (1): 15. https://doi.org/10.1186/s40965-017-0031-6.

Haralick, R. M., K. Shanmugam, and I. Dinstein. 1973. "Textural Features for Image Classification." *IEEE Transactions on Systems, Man, and Cybernetics SMC-3* (6): 610-621. https://doi.org/10.1109/TSMC.1973.4309314.

Hay, G. J., and G. Castilla. 2006. "Object-Based Image Analysis: Strengths, Weaknesses, Opportunities and Threats (SWOT)." In *Commission VI, WG VI/4. file://C:/Users/begue/Downloads/Object-Based_Image_Analysis-Strengths_Weaknesses_O.pdf*.

Hazeu, G. W., M. J. Metzger, C. A. Mücher, M. Perez-Soba, C. H. Renetzeder, and E. Andersen. 2011. "European Environmental Stratifications and Typologies: An Overview." *Agriculture, Ecosystems and Environment* 142 (1-2): 29-39. https://doi.org/10.1016/j.agee.2010.01.009.

Hudak, A. T., and C. A. Wessman. 1998. "Textural Analysis of Historical Aerial Photography to Characterize Woody Plant Encroachment in South African Savanna." *Remote Sensing of Environment* 66 (3): 317-330. https://doi.org/10.1016/S0034-4257(98)00078-9.

Inglada, J. 2016. "OTB Gapfilling, a Temporal Gapfilling for Image Time Series Library." *Zenodo*. https://doi.org/10.5281/ZENODO.45572.

Johnson, B., and Z. Xie. 2011. "Unsupervised Image Segmentation Evaluation and Refinement Using a Multi-Scale Approach." *ISPRS Journal of Photogrammetry and Remote Sensing* 66 (4): 473-483. https://doi.org/10.1016/j.isprsjprs.2011.02.006.

Kalinicheva, E., J. Sublime, and M. Trocan. 2020. "Unsupervised Satellite Image Time Series Clustering Using Object-Based Approaches and 3D Convolutional Autoencoder." *Remote Sensing* 12 (11): 1816. https://doi.org/10.3390/rs12111816.

Karasov, O., M. Külvik, and I. Burdun. 2021. "Deconstructing Landscape Pattern: Applications of Remote Sensing to Physiognomic Landscape Mapping." *Geo Journal* 86 (1): 529-555. https://doi.org/10.1007/s10708-019-10058-6.

Kavzoglu, T., and H. Tonbul. 2018. "An Experimental Comparison of Multi-Resolution Segmentation, SLIC and K-Means Clustering for Object-Based Classification of VHR Imagery." *International Journal of Remote Sensing* 39 (18): 6020-6036. https://doi.org/10.1080/01431161.2018.1506592.

Knight, J. F., R. S. Lunetta, J. Ediriwickrema, and S. Khorram. 2006. "Regional Scale Land Cover Characterization Using MODIS-NDVI 250 M Multi-Temporal Imagery: A Phenology-Based Approach." *GIScience & Remote Sensing* 43 (1): 1-23. https://doi.org/10.2747/1548-1603.43.1.1.

Koç, A., and S. Yılmaz. 2020. "Landscape Character Analysis and Assessment at the Lower Basin-Scale." *Applied Geography* 125 (December): 102359. https://doi.org/10.1016/j.apgeog.2020.102359.

Kupidura, P. 2019. "The Comparison of Different Methods of Texture Analysis for Their Efficacy for Land Use Classification in Satellite Imagery." *Remote Sensing* 11 (10): 1233. https://doi.org/10.3390/rs11101233.

Lang, M., S. Alleaume, S. Luque, N. Baghdadi, and J.-B. Féret. 2018. "Monitoring and Characterizing Heterogeneous Mediterranean Landscapes with Continuous Textural Indices Based on VHSR Imagery." *Remote Sensing* 10 (6): 868. https://doi.org/10.3390/rs10060868.

Lang, S., G. J. Hay, A. Baraldi, D. Tiede, and T. Blaschke. 2019. "Geobia Achievements and Spatial Opportunities in the Era of Big Earth Observation Data." *ISPRS International Journal of Geo-Information* 8 (11): 474. https://doi.org/10.3390/ijgi8110474.

Lè, S., J. Josse, and F. Husson. 2008. "FactoMineR: An R Package for Multivariate Analysis." *Journal of Statistical Software* 25 (1). https://doi.org/10.18637/jss.v025.i01.

Leenhardt, D., Angevin, F., Biarnès, A., Colbach, N., and Mignolet, C.l. 2010. "Describing and Locating Cropping Systems on a Regional Scale - A Review." *Agronomy for Sustainable Development* 30: 131-138. https://doi.org/10.1051/agro/2009002.

Leutner, R., N. H. Benjamin, J. Schwalb-Willmann, and R. J. Hijmans. 2022. "Package 'RStoolbox.' Tools for Remote Sensing Data Analysis. https://cran.r-project.org/web/packages/RStoolbox/RStoolbox.pdf.

Levering, A., D. Marcos, and D. Tuia. 2021. "On the Relation Between Landscape Beauty and Land Cover: A Case Study in the U.K. at Sentinel-2 Resolution with Interpretable AI." *ISPRS Journal of Photogrammetry and Remote Sensing* 177 (July): 194-203. https://doi.org/10.1016/j.isprsjprs.2021.04.020.

Li, Y., and H. Wu. 2012. "A Clustering Method Based on K-Means Algorithm." *Physics Procedia* 25: 1104-1109. https://doi.org/10.1016/j.phpro.2012.03.206.

Lonjou, V., C. Desjardins, O. Hagolle, B. Petrucci, T. Tremas, M. Dejus, A. Makarau, and S. Auer. 2016. "MACCS-ATCOR Joint Algorithm (MAJA)." In edited by A. Comerón, E. I. Kassianov, and K. Schäfer, 1000107. Edinburgh, United Kingdom. https://doi.org/10.1117/12.2240935.

Luque, S., N. Pettorelli, P. Vihervaara, M. Wegmann, and J. Vamosi. 2018. "Improving Biodiversity Monitoring Using Satellite Remote Sensing to Provide Solutions Towards the 2020 Conservation Targets." Edited by Jana Vamosi. *Methods in Ecology and Evolution* 9 (8): 1784-1786. https://doi.org/10.1111/2041-210X.13057.

Ma, L., M. Li, X. Ma, L. Cheng, P. Du, and Y. Liu. 2017. "A Review of Supervised Object-Based Land-Cover Image Classification." *Isprsjprs Journal of Photogrammetry & Remote Sensing* 130 (August): 277-293. https://doi.org/10.1016/j.isprsjprs.2017.06.001.

Meinel, G., and M. Neubert. 2004. "A Comparison of Segmentation Programs for High Resolution Remote Sensing Data." In *Geo-Imagery Bridging Continents. XXth ISPRS Congress, Istanbul, Turkey, 12-23 July (ISPRS International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences; vol. XXXV, Part B4)*, 1097-1102.

Mücher, C. A., J. A. Klijn, D. M. Wascher, and J. H. J. Schaminée. 2010. "A New European Landscape Classification (LANMAP): A Transparent, Flexible and User-Oriented Methodology to Distinguish Landscapes." *Ecological Indicators* 10 (1): 87-103. https://doi.org/10.1016/j.ecolind.2009.03.018.
<figure><img src="image_47.png" /></figure>

Naimi, B., N. A. S. Hamm, T. A. Groen, A. K. Skidmore, and A. G. Toxopeus. 2014. "Where is Positional Uncertainty a Problem for Species Distribution Modelling?" *Holarctic Ecology* 37 (2): 191–203. https://doi.org/10.1111/j.1600-0587.2013.00205.x.

Newton, A. C., R. A. Hill, C. Echeverría, D. Golicher, J. M. Rey Benayas, L. Cayuela, and S. A. Hinsley. 2009. "Remote Sensing and the Future of Landscape Ecology." *Progress in Physical Geography: Earth and Environment* 33 (4): 528–546. https://doi.org/10.1177/0309133309346882.

Ollivier, C. C., S. D. Carrière, T. Heath, A. Olioso, Z. Rabefitia, H. Rakoto, L. Oudin, and F. Satgé. 2023. "Ensemble Precipitation Estimates Based on an Assessment of 21 Gridded Precipitation Datasets to Improve Precipitation Estimations Across Madagascar." *Journal of Hydrology: Regional Studies* 47 (June): 101400. https://doi.org/10.1016/j.ejrh.2023.101400.

Opdam, P., S. Luque, J. Nassauer, P. H. Verburg, and J. Wu. 2018. "How Can Landscape Ecology Contribute to Sustainability Science?" *Landscape Ecology* 33 (1): 1–7. https://doi.org/10.1007/s10980-018-0610-7.

Ouchra, H., A. Belangour, and A. Erraissi. 2022. "Machine Learning for Satellite Image Classification: A Comprehensive Review." In *2022 International Conference on Data Analytics for Business and Industry (ICDABI)*, 1–5. Sakhir, Bahrain: IEEE. https://doi.org/10.1109/ICDABI56818.2022.10041606.

Pedregosa, F., G. Varoquaux, A. Gramfort, V. Michel, B. Thierion, O. Grisel, M. Blondel, et al. 2011. "Scikit-Learn: Machine Learning in Python." *Journal of Machine Learning Research*. https://www.jmlr.org/papers/volume12/pedregosa11a/pedregosa11a.pdf.

Peña, J. M., J. A. Lozano, and P. Larrañaga. 1999. "An Empirical Comparison of Four Initialization Methods for the K-Means Algorithm." *Pattern Recognition Letters* 20 (10): 1027–1040. https://doi.org/10.1016/S0167-8655(99)00069-0.

Préau, C., N. Dubos, M. Lenormand, P. Denelle, M. Le Louarn, S. Alleaume, and S. Luque. 2022. "Dispersal-Based Species Pools as Sources of Connectivity Area Mismatches." *Landscape Ecology* 37 (3): 729–743. https://doi.org/10.1007/s10980-021-01371-y.

Raymond, R., Y. Luginbüh, J.-F. Seguin, and Q. Cedelle. 2015. *Landscape Atlases Landscape Identification, Characterisation and Assessment Method*. Paris (FR): Ministère de l'Ecologie, du Développement Durable et de l'Energie. https://www.ecologie.gouv.fr/sites/default/files/Landscape%20Atlases%2C%20Landscape%20identification%2C%20characterisation%20and%20assessment%20method.pdf.

Reed, J., M. Ros-Tonen, and T. C. H. Sunderland. 2020. *Operationalizing Integrated Landscape Approaches in the Tropics*. Center for International Forestry Research (CIFOR). https://doi.org/10.17528/cifor/007800.

Rizzo, D., E. Marraccini, and S. Lardon, eds. 2022. *Landscape Agronomy: Advances and Challenges of a Territorial Approach to Agricultural Issues*. Cham: Springer International Publishing. https://doi.org/10.1007/978-3-031-05263-7.

Rouse, J. W., R. H. Hass, J. A. Schell, and D. W. Deering. 1973. "Monitoring Vegetation Systems in the Great Plains with ERTS." In *Third ERTS Symposium, Goddard Space Flight Center*, 309–317. Washington, DC: NASA.

Silva, R. G. P., S. Araújo Zagallo, A.-E. Laques, and C. Hiroo Saito. 2020. "Landscape Signature as an Integrative View of Landscape Metrics: A Case Study in Brazil-French Guiana Border." *Landscape Online* 85 (November): 1–18. https://doi.org/10.3097/LO.202085.

Simensen, T., R. Halvorsen, and L. Erikstad. 2018. "Methods for Landscape Characterisation and Mapping: A Systematic Review." *Land Use Policy* 75 (June): 557–569. https://doi.org/10.1016/j.landusepol.2018.04.022.

Soudani, K., G. le Maire, E. Dufrêne, C. François, N. Delpierre, E. Ulrich, and S. Cecchini. 2008. "Evaluation of the Onset of Green-Up in Temperate Deciduous Broadleaf Forests Derived from Moderate Resolution Imaging Spectroradiometer (MODIS) Data." *Remote Sensing of Environment* 112 (5): 2643–2655. https://doi.org/10.1016/j.rse.2007.12.004.

Strahler, A. H., C. E. Woodcock, and J. A. Smith. 1986. "On the Nature of Models in Remote Sensing." *Remote Sensing of Environment* 20 (2): 121–139. https://doi.org/10.1016/0034-4257(86)90018-0.

Sun, H., J. Wei, and Q. Han. 2022. "Assessing Land-Use Change and Landscape Connectivity Under Multiple Green Infrastructure Conservation Scenarios." *Ecological Indicators* 142 (September): 109236. https://doi.org/10.1016/j.ecolind.2022.109236.

Tenerelli, P., C. Püffel, and S. Luque. 2017. "Spatial Assessment of Aesthetic Services in a Complex Mountain Region: Combining Visual Landscape Properties with Crowdsourced Geographic Information." *Landscape Ecology* 32 (5): 1097–1115. https://doi.org/10.1007/s10980-017-0498-7.

Thierion, V., S. Alleaume, C. Jacqueminet, C. Vigneau, K. Michel, and S. Luque. 2014. "The Potential of Pléiades Imagery for Vegetation Mapping: An Example of Grasslands and Pastoral Environments." *Revue Française de Photogrammétrie et de Télédétection*, (208): 105–110. https://doi.org/10.52638/rfpt.2014.124.

Turner, M. G., and R. H. Gardner. 2015. *Landscape Ecology in Theory and Practice: Pattern and Process*. New York, NY: Springer New York. https://doi.org/10.1007/978-1-4939-2794-4.

Tveit, M., Å. Ode, and G. Fry. 2006. "Key Concepts in a Framework for Analysing Visual Landscape Character." *Landscape Research* 31 (3): 229–255. https://doi.org/10.1080/01426390600783269.

van Strien, M. J., and A. Grêt-Regamey. 2022. "Unsupervised Deep Learning of Landscape Typologies from Remote Sensing Images and Other Continuous Spatial Data." *Environmental Modelling & Software* 155 (September): 105462. https://doi.org/10.1016/j.envsoft.2022.105462.

Vintrou, E., M. Soumaré, S. Bernard, A. Bégué, C. Baron, and D. Lo Seen. 2012. "Mapping Fragmented Agricultural Systems in the Sudano-Sahelian Environments of Africa Using Random Forest and Ensemble Metrics of Coarse Resolution MODIS Imagery." *Photogrammetric Engineering & Remote Sensing* 78 (8): 839–848. https://doi.org/10.14358/PERS.78.8.839.

Wu, J. 2013. "Landscape Sustainability Science: Ecosystem Services and Human Well-Being in Changing Landscapes." *Landscape Ecology* 28 (6): 999–1023. https://doi.org/10.1007/s10980-013-9894-9.

Xie, Y., Z. Sha, and M. Yu. 2008. "Remote Sensing Imagery in Vegetation Mapping: A Review." *Journal of Plant Ecology* 1 (1): 9–23. https://doi.org/10.1093/jpe/rtm005.

Yazici, K. 2018. "Evaluation of Visual Landscape Quality in the Wetlands North of Sivas (Turkey)." *Applied Ecology and Environmental Research* 16 (4): 4183–4197. https://doi.org/10.15666/aeer/1604_41834197.

Zanaga, D., R. Van De Kerchove, W. De Keersmaecker, N. Souverijns, C. Brockmann, R. Quast, J. Wevers, et al. 2021. "ESA WorldCover 10 M 2020 V100." *Zenodo*. https://doi.org/10.5281/ZENODO.5571936.
<figure><img src="image_50.png" /></figure>

Naimi, B., N. A. S. Hamm, T. A. Groen, A. K. Skidmore, and A. G. Toxopeus. 2014. "Where is Positional Uncertainty a Problem for Species Distribution Modelling?" *Holarctic Ecology* 37 (2): 191–203. https://doi.org/10.1111/j.1600-0587.2013.00205.x.

Newton, A. C., R. A. Hill, C. Echeverría, D. Golicher, J. M. Rey Benayas, L. Cayuela, and S. A. Hinsley. 2009. "Remote Sensing and the Future of Landscape Ecology." *Progress in Physical Geography: Earth and Environment* 33 (4): 528–546. https://doi.org/10.1177/0309133309346882.

Ollivier, C. C., S. D. Carrière, T. Heath, A. Olioso, Z. Rabefitia, H. Rakoto, L. Oudin, and F. Satgé. 2023. "Ensemble Precipitation Estimates Based on an Assessment of 21 Gridded Precipitation Datasets to Improve Precipitation Estimations Across Madagascar." *Journal of Hydrology: Regional Studies* 47 (June): 101400. https://doi.org/10.1016/j.ejrh.2023.101400.

Opdam, P., S. Luque, J. Nassauer, P. H. Verburg, and J. Wu. 2018. "How Can Landscape Ecology Contribute to Sustainability Science?" *Landscape Ecology* 33 (1): 1–7. https://doi.org/10.1007/s10980-018-0610-7.

Ouchra, H., A. Belangour, and A. Erraissi. 2022. "Machine Learning for Satellite Image Classification: A Comprehensive Review." In *2022 International Conference on Data Analytics for Business and Industry (ICDABI)*, 1–5. Sakhir, Bahrain: IEEE. https://doi.org/10.1109/ICDABI56818.2022.10041606.

Pedregosa, F., G. Varoquaux, A. Gramfort, V. Michel, B. Thierion, O. Grisel, M. Blondel, et al. 2011. "Scikit-Learn: Machine Learning in Python." *Journal of Machine Learning Research*. https://www.jmlr.org/papers/volume12/pedregosa11a/pedregosa11a.pdf.

Peña, J. M., J. A. Lozano, and P. Larrañaga. 1999. "An Empirical Comparison of Four Initialization Methods for the K-Means Algorithm." *Pattern Recognition Letters* 20 (10): 1027–1040. https://doi.org/10.1016/S0167-8655(99)00069-0.

Préau, C., N. Dubos, M. Lenormand, P. Denelle, M. Le Louarn, S. Alleaume, and S. Luque. 2022. "Dispersal-Based Species Pools as Sources of Connectivity Area Mismatches." *Landscape Ecology* 37 (3): 729–743. https://doi.org/10.1007/s10980-021-01371-y.

Raymond, R., Y. Luginbüh, J.-F. Seguin, and Q. Cedelle. 2015. *Landscape Atlases Landscape Identification, Characterisation and Assessment Method*. Paris (FR): Ministère de l'Ecologie, du Développement Durable et de l'Energie. https://www.ecologie.gouv.fr/sites/default/files/Landscape%20Atlases%2C%20Landscape%20identification%2C%20characterisation%20and%20assessment%20method.pdf.

Reed, J., M. Ros-Tonen, and T. C. H. Sunderland. 2020. *Operationalizing Integrated Landscape Approaches in the Tropics*. Center for International Forestry Research (CIFOR). https://doi.org/10.17528/cifor/007800.

Rizzo, D., E. Marraccini, and S. Lardon, eds. 2022. *Landscape Agronomy: Advances and Challenges of a Territorial Approach to Agricultural Issues*. Cham: Springer International Publishing. https://doi.org/10.1007/978-3-031-05263-7.

Rouse, J. W., R. H. Hass, J. A. Schell, and D. W. Deering. 1973. "Monitoring Vegetation Systems in the Great Plains with ERTS." In *Third ERTS Symposium, Goddard Space Flight Center*, 309–317. Washington, DC: NASA.

Silva, R. G. P., S. Araújo Zagallo, A.-E. Laques, and C. Hiroo Saito. 2020. "Landscape Signature as an Integrative View of Landscape Metrics: A Case Study in Brazil-French Guiana Border." *Landscape Online* 85 (November): 1–18. https://doi.org/10.3097/LO.202085.

Simensen, T., R. Halvorsen, and L. Erikstad. 2018. "Methods for Landscape Characterisation and Mapping: A Systematic Review." *Land Use Policy* 75 (June): 557–569. https://doi.org/10.1016/j.landusepol.2018.04.022.

Soudani, K., G. le Maire, E. Dufrêne, C. François, N. Delpierre, E. Ulrich, and S. Cecchini. 2008. "Evaluation of the Onset of Green-Up in Temperate Deciduous Broadleaf Forests Derived from Moderate Resolution Imaging Spectroradiometer (MODIS) Data." *Remote Sensing of Environment* 112 (5): 2643–2655. https://doi.org/10.1016/j.rse.2007.12.004.

Strahler, A. H., C. E. Woodcock, and J. A. Smith. 1986. "On the Nature of Models in Remote Sensing." *Remote Sensing of Environment* 20 (2): 121–139. https://doi.org/10.1016/0034-4257(86)90018-0.

Sun, H., J. Wei, and Q. Han. 2022. "Assessing Land-Use Change and Landscape Connectivity Under Multiple Green Infrastructure Conservation Scenarios." *Ecological Indicators* 142 (September): 109236. https://doi.org/10.1016/j.ecolind.2022.109236.

Tenerelli, P., C. Püffel, and S. Luque. 2017. "Spatial Assessment of Aesthetic Services in a Complex Mountain Region: Combining Visual Landscape Properties with Crowdsourced Geographic Information." *Landscape Ecology* 32 (5): 1097–1115. https://doi.org/10.1007/s10980-017-0498-7.

Thierion, V., S. Alleaume, C. Jacqueminet, C. Vigneau, K. Michel, and S. Luque. 2014. "The Potential of Pléiades Imagery for Vegetation Mapping: An Example of Grasslands and Pastoral Environments." *Revue Française de Photogrammétrie et de Télédétection*, (208): 105–110. https://doi.org/10.52638/rfpt.2014.124.

Turner, M. G., and R. H. Gardner. 2015. *Landscape Ecology in Theory and Practice: Pattern and Process*. New York, NY: Springer New York. https://doi.org/10.1007/978-1-4939-2794-4.

Tveit, M., Å. Ode, and G. Fry. 2006. "Key Concepts in a Framework for Analysing Visual Landscape Character." *Landscape Research* 31 (3): 229–255. https://doi.org/10.1080/01426390600783269.

van Strien, M. J., and A. Grêt-Regamey. 2022. "Unsupervised Deep Learning of Landscape Typologies from Remote Sensing Images and Other Continuous Spatial Data." *Environmental Modelling & Software* 155 (September): 105462. https://doi.org/10.1016/j.envsoft.2022.105462.

Vintrou, E., M. Soumaré, S. Bernard, A. Bégué, C. Baron, and D. Lo Seen. 2012. "Mapping Fragmented Agricultural Systems in the Sudano-Sahelian Environments of Africa Using Random Forest and Ensemble Metrics of Coarse Resolution MODIS Imagery." *Photogrammetric Engineering & Remote Sensing* 78 (8): 839–848. https://doi.org/10.14358/PERS.78.8.839.

Wu, J. 2013. "Landscape Sustainability Science: Ecosystem Services and Human Well-Being in Changing Landscapes." *Landscape Ecology* 28 (6): 999–1023. https://doi.org/10.1007/s10980-013-9894-9.

Xie, Y., Z. Sha, and M. Yu. 2008. "Remote Sensing Imagery in Vegetation Mapping: A Review." *Journal of Plant Ecology* 1 (1): 9–23. https://doi.org/10.1093/jpe/rtm005.

Yazici, K. 2018. "Evaluation of Visual Landscape Quality in the Wetlands North of Sivas (Turkey)." *Applied Ecology and Environmental Research* 16 (4): 4183–4197. https://doi.org/10.15666/aeer/1604_41834197.

Zanaga, D., R. Van De Kerchove, W. De Keersmaecker, N. Souverijns, C. Brockmann, R. Quast, J. Wevers, et al. 2021. "ESA WorldCover 10 M 2020 V100." *Zenodo*. https://doi.org/10.5281/ZENODO.5571936.
# Appendices

## Appendix 1.

Example of the PCA correlation circle of the Haralick MODIS texture indices calculated on the mean MODIS NDVI image (2016–2020); the redder the indices, the stronger their contribution to the dimensions axis.

<figure><img src="image_48.png" /></figure>

## Appendix 2.

### The segmentation Global Scores

Johnson and Xie (2011) proposed an unsupervised method to evaluate the quality of multispectral image segmentations. The method relies on the weighted variance wV (Equation 1) as a measure of intra-segment homogeneity, and the Moran's autocorrelation index M (Equation 2) as a measure of inter-segment disparity (Espindola et al. 2006).

$$ wV = \frac{\sum_{i=1}^{n} a_i * v_i}{\sum_{i=1}^{n} a_i} \qquad (1) $$

where n is the total number of segments, a is the area of the segment i and v the variance of the segment i.

$$ M = \frac{n \sum_{i=1}^{n} \sum_{j=1}^{n} \omega_{ij} (y_i - \bar{y})(y_j - \bar{y})}{\sum_{j=1}^{n} (y_i - \bar{y})^2 \left( \sum_{i \neq j} \sum \omega_{ij} \right)} \qquad (2) $$

where n is the total number of regions, w<sub>ij</sub> is a measure of the spatial proximity of the segments, y<sub>i</sub> is the mean value of segment i, and y is the mean value of the image. Each weight w<sub>ij</sub> is a measure of the spatial proximity of segments i and j (if the segments are adjacent, w<sub>ij</sub> is equal to 1; otherwise, it is equal to 0).

For each segmentation scale, wV and M are calculated for each variable (or band), and normalized (wV<sub>norm</sub> and M<sub>norm</sub>) in order to balance their relative importance. Two versions of wV<sub>norm</sub> and M<sub>norm</sub> are proposed: a min-max normalization (Equation 3), and the normalization proposed by Böck et al. (2017) (Equation 4) which was shown to make the global score less dependent on the tested segmentations (stabilization of the optimal value).
# Appendices

## Appendix 1.

Example of the PCA correlation circle of the Haralick MODIS texture indices calculated on the mean MODIS NDVI image (2016–2020); the redder the indices, the stronger their contribution to the dimensions axis.

<figure><img src="image_51.png" /></figure>

## Appendix 2.

### The segmentation Global Scores

Johnson and Xie (2011) proposed an unsupervised method to evaluate the quality of multispectral image segmentations. The method relies on the weighted variance $wV$ (Equation 1) as a measure of intra-segment homogeneity, and the Moran's autocorrelation index $M$ (Equation 2) as a measure of inter-segment disparity (Espindola et al. 2006).

$$wV = \frac{\sum_{i=1}^{n} a_i * v_i}{\sum_{i=1}^{n} a_i} \qquad (1)$$

where $n$ is the total number of segments, $a$ is the area of the segment $i$ and $v$ the variance of the segment $i$.

$$M = \frac{n \sum_{i=1}^{n} \sum_{j=1}^{n} \omega_{ij} (y_i - \bar{y})(y_j - \bar{y})}{\sum_{j=1}^{n} (y_i - \bar{y})^2 \left( \sum_{i \neq j} \sum \omega_{ij} \right)} \qquad (2)$$

where $n$ is the total number of regions, $w_{ij}$ is a measure of the spatial proximity of the segments, $y_i$ is the mean value of segment $i$, and $y$ is the mean value of the image. Each weight $w_{ij}$ is a measure of the spatial proximity of segments $i$ and $j$ (if the segments are adjacent, $w_{ij}$ is equal to 1; otherwise, it is equal to 0).

For each segmentation scale, $wV$ and $M$ are calculated for each variable (or band), and normalized ($wV_{norm}$ and $M_{norm}$) in order to balance their relative importance. Two versions of $wV_{norm}$ and $M_{norm}$ are proposed: a min-max normalization (Equation 3), and the normalization proposed by Böck et al. (2017) (Equation 4) which was shown to make the global score less dependent on the tested segmentations (stabilization of the optimal value).
<figure><img src="image_49.png" /></figure>

$$wV_{norm} = \frac{(V - V_{min})}{(V_{max} - V_{min})} \quad \& \quad M_{norm} = \frac{(M - M_{min})}{(M_{max} - M_{min})} \qquad (3)$$

$$wV_{norm} = \frac{V}{V'} \quad \& \quad M_{norm} = \frac{(M + 1)}{2} \qquad (4)$$

with V' the variance of the image.

To assess the optimum segmentation scale, the evaluation global score J is provided (Equation 5).

$$J = wV_{norm} + M_{norm} \qquad (5)$$

For multiband images, Johnson and Xie (2011) recommended averaging the global scores calculated for each band individually.

## Appendix 3.

Mean (and standard deviation) of RS-ELV values for the six classes of the radiometric landscape map. Mean NDVI corresponds to the mean NDVI of the MODIS time series (2016-2020). Seasonal (Oct and April) NDVI and NDWI, and NDVI-textural indices were calculated from Sentinel-2 images acquired in 2019 (see text for details).

<table><thead><tr><th rowspan="2">Landscape class</th><th>MODIS</th><th>S2</th><th>S2</th><th>S2</th><th>S2</th><th>S2</th><th>S2</th><th>S2</th><th>S2</th></tr><tr><th>NDVI mean</th><th>NDVI April</th><th>NDVI Oct</th><th>NDWI April</th><th>NDWI Oct</th><th>Entropy April</th><th>H-correl April</th><th>Entropy Oct</th><th>H-correl Oct</th></tr></thead><tbody><tr><td>5</td><td>0.742<br>(0.033)</td><td>0.720<br>(0.088)</td><td>0.670<br>(0.123)</td><td>0.608<br>(0.107)</td><td>0.419<br>(0.198)</td><td>0.134<br>(0.271)</td><td>149.97<br>(11.16)</td><td>0.539<br>(0.397)</td><td>166.15<br>(22.01)</td></tr><tr><td>1</td><td>0.465<br>(0.060)</td><td>0.602<br>(0.099)</td><td>0.478<br>(0.098)</td><td>0.356<br>(0.145)</td><td>0.045<br>(0.179)</td><td>0.782<br>(0.455)</td><td>177.42<br>(29.38)</td><td>1.451<br>(0.372)</td><td>197.21<br>(55.90)</td></tr><tr><td>2</td><td>0.602<br>(0.068)</td><td>0.671<br>(0.091)</td><td>0.569<br>(0.112)</td><td>0.444<br>(0.151)</td><td>0.187<br>(0.197)</td><td>0.359<br>(0.357)</td><td>158.49<br>(20.67)</td><td>1.196<br>(0.511)</td><td>209.34<br>(41.09)</td></tr><tr><td>3</td><td>0.726<br>(0.060)</td><td>0.729<br>(0.078)</td><td>0.696<br>(0.114)</td><td>0.565<br>(0.148)</td><td>0.409<br>(0.210)</td><td>0.136<br>(0.219)</td><td>150.04<br>(13.01)</td><td>0.445<br>(0.414)</td><td>162.76<br>(23.13)</td></tr><tr><td>4</td><td>0.835<br>(0.029)</td><td>0.803<br>(0.046)</td><td>0.817<br>(0.076)</td><td>0.698<br>(0.070)</td><td>0.611<br>(0.142)</td><td>0.157<br>(0.238)</td><td>152.11<br>(15.69)</td><td>0.639<br>(0.470)</td><td>182.48<br>(45.48)</td></tr><tr><td>6</td><td>0.667<br>(0.058)</td><td>0.640<br>(0.183)</td><td>0.578<br>(0.169)</td><td>0.535<br>(0.161)</td><td>0.370<br>(0.191)</td><td>0.455<br>(0.567)</td><td>165.04<br>(34.14)</td><td>1.038<br>(0.579)</td><td>200.44<br>(49.37)</td></tr></tbody></table>
<figure><img src="image_52.png" /></figure>

$$wV_{norm} = \frac{(V - V_{min})}{(V_{max} - V_{min})} \quad \& \quad M_{norm} = \frac{(M - M_{min})}{(M_{max} - M_{min})} \qquad (3)$$

$$wV_{norm} = \frac{V}{V'} \quad \& \quad M_{norm} = \frac{(M + 1)}{2} \qquad (4)$$

with V' the variance of the image.

To assess the optimum segmentation scale, the evaluation global score J is provided (Equation 5).

$$J = wV_{norm} + M_{norm} \qquad (5)$$

For multiband images, Johnson and Xie (2011) recommended averaging the global scores calculated for each band individually.

## Appendix 3.

Mean (and standard deviation) of RS-ELV values for the six classes of the radiometric landscape map. Mean NDVI corresponds to the mean NDVI of the MODIS time series (2016-2020). Seasonal (Oct and April) NDVI and NDWI, and NDVI-textural indices were calculated from Sentinel-2 images acquired in 2019 (see text for details).

<table><thead><tr><th rowspan="2">Landscape class</th><th>MODIS</th><th>S2</th><th>S2</th><th>S2</th><th>S2</th><th>S2</th><th>S2</th><th>S2</th><th>S2</th></tr><tr><th>NDVI</th><th>NDVI</th><th>NDVI</th><th>NDWI</th><th>NDWI</th><th>Entropy</th><th>H-correl</th><th>Entropy</th><th>H-correl</th></tr></thead><tbody><tr><td>5</td><td>0.742<br>(0.033)</td><td>0.720<br>(0.088)</td><td>0.670<br>(0.123)</td><td>0.608<br>(0.107)</td><td>0.419<br>(0.198)</td><td>0.134<br>(0.271)</td><td>149.97<br>(11.16)</td><td>0.539<br>(0.397)</td><td>166.15<br>(22.01)</td></tr><tr><td>1</td><td>0.465<br>(0.060)</td><td>0.602<br>(0.099)</td><td>0.478<br>(0.098)</td><td>0.356<br>(0.145)</td><td>0.045<br>(0.179)</td><td>0.782<br>(0.455)</td><td>177.42<br>(29.38)</td><td>1.451<br>(0.372)</td><td>197.21<br>(55.90)</td></tr><tr><td>2</td><td>0.602<br>(0.068)</td><td>0.671<br>(0.091)</td><td>0.569<br>(0.112)</td><td>0.444<br>(0.151)</td><td>0.187<br>(0.197)</td><td>0.359<br>(0.357)</td><td>158.49<br>(20.67)</td><td>1.196<br>(0.511)</td><td>209.34<br>(41.09)</td></tr><tr><td>3</td><td>0.726<br>(0.060)</td><td>0.729<br>(0.078)</td><td>0.696<br>(0.114)</td><td>0.565<br>(0.148)</td><td>0.409<br>(0.210)</td><td>0.136<br>(0.219)</td><td>150.04<br>(13.01)</td><td>0.445<br>(0.414)</td><td>162.76<br>(23.13)</td></tr><tr><td>4</td><td>0.835<br>(0.029)</td><td>0.803<br>(0.046)</td><td>0.817<br>(0.076)</td><td>0.698<br>(0.070)</td><td>0.611<br>(0.142)</td><td>0.157<br>(0.238)</td><td>152.11<br>(15.69)</td><td>0.639<br>(0.470)</td><td>182.48<br>(45.48)</td></tr><tr><td>6</td><td>0.667<br>(0.058)</td><td>0.640<br>(0.183)</td><td>0.578<br>(0.169)</td><td>0.535<br>(0.161)</td><td>0.370<br>(0.191)</td><td>0.455<br>(0.567)</td><td>165.04<br>(34.14)</td><td>1.038<br>(0.579)</td><td>200.44<br>(49.37)</td></tr></tbody></table>