Introduction

Landforms represent the long-term development of geologic and geomorphologic processes (Bolongaro-Crevenna et al. 2005; Oldroyd and Grapes 2008; Kleman et al. 2016). They tend to reflect the interaction of climate, tectonics, erosion and deposition (Castelltort et al. 2015; Zhang et al. 2016; Marshall et al. 2017). An increasing amount of the research (Szabó et al. 2010; Hooke 2012; Ellis et al. 2013; Goudie and Viles 2016; Tarolli and Sofia 2016; Tarolli 2016; Brown et al. 2017; Migoń and Latocha 2018; Goudie 2018; Tarolli et al. 2019) has pointed out that human activity has played a pivotal role as geomorphic forcing. For instance, agriculture is susceptible to accelerate soil erosion (Tóth 2010; Curebal et al. 2015; Borrelli et al. 2017), dams and reservoirs engineering interrupt the continuity of sediment transport in rivers system (Tessler et al. 2016; Wang et al. 2016; Poeppl et al. 2017), road network construction is associated with slope stability of roadcut and other geological risks (Csima 2010; Sidle and Ziegler 2012; Penna et al. 2014; Ramos-scharrón 2018).

With this literature, the concept of surface reshaping from both abiotic and biotic forces has emerged (Ellis 2004; Steiger and Corenblit 2012; Pietrasiak et al. 2014; Tarolli and Sofia 2016). As suggested by Dietrich and Perron (2006), small-scale biotic processes can influence the form of landscapes and create a distinctive topography, but this has yet to be investigated for human-made landforms. Identifying natural and anthropogenic features and further distinguishing the landform signatures still poses a significant challenge for the geomorphological community (Tarolli et al. 2019). Thanks to the progress in remote sensing techniques and open-access datasets, the recognition of large-scale geomorphic signatures is now possible at various scales (Evans 1980; Nagel et al. 2014; Sofia et al. 2014a; Tarolli 2014; Byun and Seong 2015; Jordan et al. 2016). However, an explicit characterisation, from a morphological point of view, of natural and anthropogenic surfaces and for different landscape contexts, is still missing. This study showcases how high-resolution topographic data can offer the basis to (1) characterise specific signatures with land covers on the basis of an objective geomorphometric analysis; (2) demonstrate that the anthropogenic and natural land covers show a statistically different underlining morphology and (3) understand (where present) the degree of anthropogenic impact due to the various land covers.

Since a significant concern is how natural systems are being modified or transformed by anthropogenic land uses, one crucial issue is how the different land surface should be disaggregated for modelling and further analysis, and if any generic relationships can be identified between land uses and morphological transformations to the landscape. Geosciences must advance towards empirical and theoretical frameworks that integrate the natural and sociocultural forces that are now among the leading shapers of Earth’s surface processes (Tarolli et al. 2019) to understand the causes and consequences of these transformations and contribute to building a sustainable future. This work offers an example of such an empirical framework, providing a diagnostic tool to infer objectively morphological differences within various landscapes. Processes happen in three dimensions and observing the topographic differences among land covers offer a basis to potentially infer differences in the processes happening in these landscapes.

Study area

This study investigates four study areas of 10 × 10 km in northeastern Italy, representing different landscapes, from floodplains to mountains (Fig. 1): the Veneto floodplain (floodplain, Fig. 1a), Colli Euganei (plain to hilly, Fig. 1b), Venetian Prealps (hilly, Fig. 1c) and Trentino (Alpine mountains, Fig. 1d). These sites were selected because they share close geographic locations and distributions of land covers, but they differ in landforms.

Fig. 1
figure 1

Considered four study areas: a floodplain, b plain to hilly, c hilly, d mountain

The elevations in Veneto floodplain (Fig. 1a) range from − 2 to 10 m a.s.l., with an average of 4 m a.s.l. Seventy-five percent of the area has a height lower than 5 m a.s.l. This area is characterised by a higher level of anthropogenic pressure, especially agricultural landscapes, due to urbanisation and industrialisation. The area is intensively drained for reclamation and irrigation purposes through a dense network of channels and ditches (Fig. 2a). The plain to hilly area (Fig. 1b) has an elevation ranging from 0.4 to 601 m a.s.l. (average 112 m a.s.l.). Seventy-eight percent of the area has a height lower than 200m a.s.l. These hills are of volcanic origins and rise between 300 and 600 m. Vineyard cultivation is widespread in this area (Fig. 2b). The elevation in hilly area (Fig. 1c) ranges between 88 and 889 m a.s.l. (average 251 m a.s.l.). Ninety-five percent of the area is concentrated on the height between 100 and 500 m a.s.l. As in the Euganei, vineyard is also a typical characteristic of this area (Fig. 2c). The fourth area (Figs. 1d and 2d) is an Alpine landscape. The elevation ranges from 541 to 2488 m a.s.l. with a mean value of 1577 m a.s.l. Eighty percent of the area has the height from 1000 to 2000 m a.s.l.

Fig. 2
figure 2

The field overview of the study areas: a floodplain, b plain to hilly, c hilly, d mountain (photo in 2a by P. Claudio; photo in 2b by M. Luca; photo in 2c by B. Eros, photo in 2d from www.abfotografia.it)

Method

Data

Light detection and ranging–derived digital terrain models (DTMs) at 2-m resolution (Fig. 3) are available thanks to public authorities in Italy (Italian Ministry for Environment, Land and Sea; Treviso Province; Trentino Alto-Adige Autonomous Region). The datasets refer to the year’s range 2010–2012.

Fig. 3
figure 3

LiDAR DTMs of the four study areas: a floodplain, b plain to hilly, c hilly, d mountain

Information about land cover is available through the Corine-Land-Cover database (CLC, Coordination of Information on the Environment Land Cover) classification, as also reported by the local authorities. The considered CLC data come from an updated version of the Urban Atlas (European Environment Agency 2012) provided by the local government (Regione del Veneto 2012). The original Urban Atlas is mainly based on the combination of statistical image classification and visual interpretation of very high resolution (VHR) satellite imagery. Multispectral SPOT 5 & 6 and Formosat-2 pan-sharpened imagery with a 2-to 2.5-m spatial resolution is used as input data. The built-up classes are combined with density information on the level of sealed soil derived from the high-resolution layer imperviousness to provide more detail in the density of the urban fabric (European Environment Agency 2012). The updated version was enriched by the local government (Regione del Veneto) with functional information (road network, services, utilities…) using ancillary data sources such as regional cartography, forest inventories, road network graphs, aerial photographs and ground surveys.

For the purpose of this work, we focused on artificial surfaces, agriculture and forest (level 1 of the CLC classification). However, due to the large-scale cultivation of vineyard in the plain to hilly and hilly areas, which we expect to have a significant impact on the morphology of the surfaces, we defined vineyard as an independent classification from agriculture. As well, we considered grass as an independent land cover because it may occur naturally or as the result of human activity (pastures, park and recreational sites), and this allows us to understand better the associated anthropogenic impact on land covers. The land cover classification can be seen from Fig. 4.

Fig. 4
figure 4

The land cover classification (Corine Land Cover related to 2012) in the study areas: a floodplain, b plain to hilly, c hilly, d mountain. The black rectangle areas are the case studies in Fig. 12

Geomorphometric parameters

To make quantitative measurements of landscape properties, we considered three geomorphometric parameters: slope and mean curvature proposed by Evans (1980) and the Spc developed by Sofia et al. (2014a).

Slope and curvature

Evans (1980) describes the DTM surface is approximated to a bivariate quadratic function in the form of:

$$ Z=a{x}^2+b{y}^2+ cxy+ dx+ ey+f $$
(1)

where x, y and z are local coordinates, and a to f are quadratic coefficients.

From such a surface, it is possible to compute the first (slope, Eq. (2)) and second (curvature, Eq. (5)) derivative. Slope (Fig. 5) is calculated as:

$$ \mathrm{Slope}=\arctan \sqrt{d^2+{e}^2} $$
(2)
Fig. 5
figure 5

Slope maps for the four study areas: a floodplain, b plain to hilly, c hilly, d mountain. Slope is colour-coded according to 1 to 2 times intervals of standard deviation σ from the mean μ

where d and e are coefficients from Eq. s(1).

Curvature is the second derivative of the surface, also referred to the change rate of slope gradient or direction (Wilson and Gallant 2000), and it emphasises convex and concave elements in the landscape. Evans (1980) proposes two measure of curvature, maximum and minimum, and Wood (1996) testifies that only the resolution of the DTM and the neighbouring cells relevant to these parameters and further defined as

$$ {\mathrm{curvature}}_{\mathrm{max}}=k\times g\left(-a-b+\sqrt{{\left(a-b\right)}^2+{c}^2}\right) $$
(3)
$$ {\mathrm{curvature}}_{\mathrm{min}}=k\times g\left(-a-b-\sqrt{{\left(a-b\right)}^2+{c}^2}\right) $$
(4)

where a, b and c are quadratic coefficients from Eq. (1), g is the grid resolution of the DTM (2 m) and k is the size of the moving window.

From the Eq. (3) and (4), mean curvature (Fig. 6) can be defined as:

$$ {\mathrm{curvature}}_{\mathrm{mean}}=\frac{{\mathrm{curvature}}_{\mathrm{max}}+\kern0.5em {\mathrm{curvature}}_{\mathrm{min}}\ }{2} $$
(5)
Fig. 6
figure 6

Mean curvature maps for the four study areas: a floodplain, b plain to hilly, c hilly, d mountain. Mean curvature is colour-coded according to 1 to 2 times intervals of standard deviation σ from the mean μ

Surface peak curvature

The Spc is inversely correlated with anthropogenic pressure (Chen et al. 2015; Sofia et al. 2016; Xiang et al. 2018). Surface morphology (slope) of regions presenting anthropogenic structures tends to be well organised (low Spc) and, in general, self-similar at a long distance. The basis for the evaluation of the Spc is the Slope Local Length of Auto-Correlation (SLLAC). This index quantifies the local self-similarities of slope (Sofia et al. 2014a). It is based on the (demonstrated) assumption that natural areas present low correlations within a neighbourhood because they are inherently irregular, while artificial surfaces to satisfy human needs for mobility and machine access tend to display a higher level of self-similarity with surroundings (Sofia et al. 2014a; Xiang et al. 2019). Describing the algorithm in detail is beyond the scope of this study: the authors refer to Sofia et al. (2014a) for a complete description of the procedure and to other examples of applications (Chen et al. 2015; Sofia et al. 2016; Tarolli and Sofia 2016; Xiang et al. 2018; Xiang et al. 2019).

Briefly, the steps to obtain the Spc (Fig. 7) are as follows:

  1. 1)

    Evaluate correlation between a moving window (W) and a patch (T) centred at the centre of the moving window. The implemented algorithm computes a normalised cross-correlation between a template and the patch, in the spatial frequency domain, and reports a standardised value that ranges between 0 (no correlation) and 1 (perfect correlation). The larger the absolute values, the stronger of the correlation.

Fig. 7
figure 7

Spc maps for the four study areas: a floodplain, b plain to hilly, c hilly, d mountain. Spc is colour-coded according to 1 to 2 times intervals of standard deviation σ from the mean μ

$$ {\mathrm{Corr}}_{\left(i,j\right)}=\frac{\sum_{u,v}\left({W}_{\left(i+u,j+v\right)}-{\overline{W}}_{i,j}\right)\left({T}_{u,v}-\overline{T}\right)}{{\left({\sum}_{u.v}{\left({W}_{\left(i+u,j+v\right)}-{\overline{W}}_{i,j}\right)}^2{\sum}_{u,v}{\left({T}_{u,v}-\overline{T}\right)}^2\right)}^{0.5}} $$
(6)
  1. 2)

    Evaluate the correlation length (L) thresholding at 37% (ISO 2013; Whitehouse 2011), the maximum correlation value (Eq. (6)). The length of correlation is the length of the longest line passing through the central pixel and connecting two boundary pixels on the extracted area connected to the central pixel (SLLAC map in Additional file 1).

  2. 3)

    Evaluate the Spc (surface peak curvature) of the SLLAC map defined as:

$$ \mathrm{Spc}=-\frac{1}{2n}{\sum}_{i=1}^n\left[\left(\frac{\partial^2z\left(x,y\right)}{\partial^2x}\right)+\left(\frac{\partial^2z\left(x,y\right)}{\partial^2y}\right)\right], $$
(7)

for every peak (pixel higher than its eight nearest neighbours). Where z stands for SLLAC value, x and y represent the cell spacing, n is the number of considered peaks.

Please refer to the supplement to infer statistic values (mean, median, STD, MAD, skewness…) of each geomorphometric parameters within each land cover.

Statistical analysis

We expect that the topographic signature of anthropogenic activities may be more subtle than the presence of a specific landform and that it would likely be a signature on the frequency of occurrence of the various degrees of the investigated landscape properties (slope, curvature, Spc). That is, the frequency distributions of these measurements would be very different, even though all observed landform types would be found in both natural or anthropogenically modified landscapes. Therefore, we observed the probability density function (PDF) of the considered landscape parameters to (1) investigate statistical differences in geomorphological surfaces between land covers under different landforms contexts and (2) explore the specific topographic signatures of land uses. For this work, the PDFs are a probability density estimate for the sampled data. The estimate is based on a normal kernel function and is evaluated at equally spaced points that cover the range of the sampled data. The distance between points is chosen automatically, based on the range of values. This means that it can be very narrow (< 0.001) for landscape parameter with small magnitude. In these cases, the PDFs can reach values much greater than 1, but their integral over any interval is always less or equal to 1.

After statistically ensuring that the datasets did not present a normal distribution and they exhibit heteroscedasticity, we decided to consider a Kruskal-Wallis test (McKight and Najab 2010) to evaluate whether there were significant differences between landscape properties underneath a specific land cover, across multiple landscapes, and we set a p value threshold of 0.05 for significance. The null hypothesis for this test is that the data for each group are statistically equal.

To investigate the similarities in PDFs between land covers, we applied the two-sample Kolmogorov-Smirnov test, which specifies the equality of probability distribution between two samples (Wilcox 2005; Razali and Wah 2011). One thousand points within each land cover were randomly selected and tested ten times to ensure the robustness of the results.

Results and discussion

Signatures recognition with different land covers

The PDF of slope (Fig. 8) exhibits very different appearances throughout the investigated landscapes. The central tendency moves from lower value to higher value, and the PDF itself tends to be more dispersed, as we increase landscape elevation. Even though the slope PDFs are always skewed, those in steeper topography present (as expected) a much longer tail than that of more gentle landscapes (i.e. floodplain). Taking a closer view of land cover distributions, the forest distribution in hilly and mountain areas present lower skewness respect to that in the floodplain. As well, most of the land covers in the mountain site show lower asymmetry. This could be an underlining symptom that humans activities are less marked in the mountains rather than in floodplains. Sofia et al. (2017) showed, for the Veneto region, different trends in anthropogenic expansion depending on the topographic location, highlighting a significant pressure in floodplains rather than in high mountains. Other works also highlighted how anthropogenic processes in the Alps are not fundamentally different from the processes in the floodplains, but they occur with a time lag and on a smaller scale (Perlik et al. 2001). Consequentially, the human signature on morphology might be less marked (Sofia et al. 2016). At the same time, the anthropogenic signatures on morphology in the Alpine environments reflect the fact that activities are generally shaped through valley bottoms and ridges, and by limits due to the slope and the steepness of the terrain (Forman et al. 2003).

Fig. 8
figure 8

The PDFs of the slope with different land covers in the four study areas: a floodplain, b plain to hilly, c hilly, d mountain. The vineyard in plain to hilly (b) and hilly (c) areas present bimodal curve

A further interesting result is the striking similarity between the PDFs of vineyards in the plain to hilly (b) and the hilly sites (c): for these landscapes, the PDFs present a double peak.

The first peak falls at a range of 0–3°, while the second peak falls around 3–10°. It is possible to note that in this landscape (Fig. 9), vineyards are constructed over terraces: the terrace walls present slope with the highest values (the second peak); on the other hand, the slope with lower values (the first peak) is related to the terrace benches. This peculiar double peak in terraced landscapes was also showed by (Tarolli and Sofia 2016), for terraces related to urbanisation over hillslopes.

Fig. 9
figure 9

Overview of a typical vineyard in the plain to hilly landscape (from (a), the blue color represents the first slope peak value which is the terrace bench, and the green shows the second slope peak value which represents the terrace walls; (b) and (c) presents the corresponding elevation and satellite image)

Focusing on mean curvature (Fig. 10), all landforms present a distribution that peaks around 0. The extreme values on the positive side are related to divergent-convex landforms, and they are generally associated with the dominance of hillslopes. The presence of extreme values on the negative side is related to convergent-concave landforms associated generally with fluvial-dominated erosion (Tarolli et al. 2012; Evans 2013). As shown in Fig. 10, the long tails of extreme value are related to artificial land covers in the floodplain and to forests in the plain to hilly and hilly areas.

Fig. 10
figure 10

The PDFs of the mean curvature with different land covers in the four areas: a floodplain, b plain to hilly, c hilly, d mountain

To better identify the reason behind these long tails, we used boxplots to detect the positive and negative outliers (Fig. 11) The idea behind this is that convex features can be identified as curvature values above the upper bound, and on the contrary, concave features can be identified as value below the lower bound (Sofia et al. 2014b).

Fig. 11
figure 11

The boxplot of mean curvature and the negative/positive threshold in the floodplain

As we can see from the satellite image (Fig. 12), the negative outliers of curvature in the floodplain are mostly related to channel networks, while the positive outliers are related to scarps, levees and small banks around them. Some noise in the curvature map is given by the footprint of the urban area, where negative and positive outliers can be found around buildings. For the hilly area, the outliers on the positive side are related to ridges, while the negative ones are channelized valleys, where forests are mostly present.

Fig. 12
figure 12

The positive and negative outliers identified from mean curvature in the floodplain (a) and hilly area (c), the left side are the DTMs, and the right side (b) and (d) are the satellite images. The location of the study areas is marked in black rectangle in Fig. 4

Statistic test on the morphology of different land covers under various landforms

The results of Kruskal-Wallis test (Table 1) show that significant differences (p < 0.05) exist among land covers regarding their slope.

Table 1 Kruskal-Wallis significance tests of difference in slope between different land covers in each study area. Source means the origin of variance, and it includes the variance between groups (different land covers) and error (geomorphometric value in the same land cover); SS means sum of square; MS represent standard deviation, and it could reflect the degree of dispersion of a dataset; chi sq. is the H-statistic of the Kruskal-Wallis test, which is approximately chi-square distributed. The Prob. > chi sq. is the p value

Looking at the p values for the Kolmogorov-Smirnov test (Table 2), it is possible to highlight that (1) there are significant differences in slope between any pair of land covers in each case study areas, except for grass and vineyard in hilly sites (marked with asterisk symbol); (2) among all the p value, the most remarkable difference always relates to the forest and any of the other land covers.

Table 2 The significance level of two-sample Kolmogorov-Smirnov tests of pairwise differences in slope between land covers within topographic types

We addressed the same analysis on mean curvature. However, the results (Table 3) show no significant difference among land covers (p value > 0.05).

Table 3 Kruskal-Wallis significance tests of mean curvature in each study area based on different land covers

As a trial test, we randomly sampled 10,000 points from the maps: with this enlarged dataset, the results show a significant difference among various land covers as a group or pairwise (Tables 4 and 5). Crops and vineyards give exceptions to this in the plain to the hilly area and also crop and forest in the mountain area (p value > 0.05).

Table 4 Kruskal-Wallis significance tests of mean curvature in each study area based on different land covers with 10,000 sample
Table 5 The significant level of two-sample Kolmogorov-Smirnov tests of pairwise differences in mean curvature between land covers within topographic types with 10,000 sample

When observing the Spc and the result of the Kruskal-Wallis test (Table 6), it is confirmed that the different land covers present different topography signatures (p value < 0.05). The exceptional case which doesn't show the significant difference makred with asterisk symbol.  

Table 6 Kruskal-Wallis significance tests of Spc in each study area based on different land covers

As it is possible to infer from Table 7: (1) pairwise differences exist in all land covers except for grass and vineyards(marked with asterisk symbol). (2) Forest differentiates from other land covers on the Spc, and this is evident for all landforms’ context considered.

Table 7 The significant level of two-sample Kolmogorov-Smirnov test of pairwise differences in Spc between land covers within topographic types

The distinct anthropogenic impact analysis

The Spc is mathematically related to the percentage of anthropogenic activity (Chen et al. 2015; Xiang et al. 2018; Xiang et al. 2019). As a consequence, it is a proxy to illustrate the extent of human impact on morphology for each land cover in different landforms (Fig. 13). The most recognisable topographic signature is that given by the forests in the floodplain (Fig. 13a). Besides, the peak of forest distribution (higher due to the small and uniform surface) falls within values of Spc on the range of those obtained in literature for highly anthropogenic surfaces. The forest distribution in both plain to hilly (Fig. 13b) and hilly areas (Fig. 13c) present a similar trend, but the peak in the hilly area tends to be to the right side (where Spc values are referred to be more ‘natural’ if compared with the mentioned literature). When moving to the mountain (Fig. 13d), the forest shape appears less skewed, which might indicate that lower human interference is present.

Fig. 13
figure 13

The PDFs of the Spc. with different land covers in the four areas. a floodplain, b plain to hilly, c hilly, d mountain

In Fig. 14a and b, the map of Spc and forest as seen from the satellite is shown in detail. The forest here is closed to the channel and mostly on the levees. This also gives a reasonable interpretation of the lower value of Spc due to the human alteration on the forest in the floodplain. By contrast, we highlight a small area with different Spc values (Fig. 14c and e, marked with different colours) and the LiDAR-derived shaded relief map (Fig. 14d) in a different area where the transition from plain to hilly is evident. The forests with relatively lower Spc value are located in areas surrounded by agricultural terraces and other anthropogenic surfaces. On the other side, the forests with higher Spc value are distributed on the top of small hills and tend to be more natural. Forest not only presents a remarkable difference with other land covers but also shows a different morphology based on the degree of anthropogenic disturbance.

Fig. 14
figure 14

Forest with Spc value in floodplain (a) and 3D overview (b) of the surrounding from Google Earth (the red outline extracted from forest); forest with Spc value in plain to hilly (c), LiDAR (d) and 3D overview (e); the outline with different colour denote the different Spc values

Xiang et al. (2019) highlighted how morphological differences under forest cover emerged by considering natural forests or artificial plantations, with higher Spc for natural forests. In actuality, the forest (mixed of shrubs and medium trees) in the floodplain area considered in this study have been altered in their structure and distribution, thus appearing as small patches surrounded by agricultural and urban areas, in lands highly disturbed by human activities. Forests in this floodplain are also managed to adopt peculiar forestry technique to preserve and maintain the vegetation, through new plantations near the ancient wood (Bellio and Pividori 2009). Our results confirm that, from a morphological perspective, the described forest for the floodplain is related to a surface affected by anthropogenic pressure, while a lower anthropogenic disturbance might be present on forests on hilly places, with the existence of more natural forests. Some land covers do not exhibit apparent differences in specific landforms. For example, grass and vineyards present some similarities in hilly areas (Table 2). As it is shown on Fig. 15, the floodplain grass (left side) has lower values of Spc, which implies that anthropogenic disturbance might be relatively significant in this environment. Taking a closer view (Fig. 15d), this area is related to a sports field and a park. By contrast, the grass in the hilly area (right side) shows a relatively higher Spc value than that of the floodplain. This could imply that anthropogenic modifications in this grassland still exist, even though there are less marked than in floodplain, and from the satellite image (Fig. 15h), the grass in hilly can be identified as pastures.

Fig. 15
figure 15

Grass with Spc value in floodplain (a) and hilly (d), LiDAR-derived shaded relief map in floodplain (b) and hilly (e), and satellite image in floodplain (c) and hilly (f); panels d and f are the 3D view from Google map in these two study areas

The similarities in morphology (Spc) between grass and vineyard in hilly areas can be explained by the fact these pastures are intensively farmed to maximise forage production. Some of them have been abandoned (‘prati vegri’) or converted to vineyards. Therefore, the human-morphological signature is still evident and appears in the Spc.

Conclusions

The primary goal of this work was to investigate the statistical differences of surface morphologies of anthropogenic and natural land covers, testifying that human activities alter the landforms from a statistical point of view. The work highlights how, if we were to make quantitative measurements of landscape properties on a landscape without human interference, and compare them to measurements of a landscape where humans activities are preponderant; the frequency distributions of these measurements would be very different, even though all observed landform types would be found in both realities. Possibly, if we had to model and describe a landscape without humans, it would look much different from the one we are used to, but it would not exhibit much different landforms. Rather, the subtle differences would lie in the frequency distributions of specific landform properties. This work also confirms the possibility to recognise with a pure geomorphometric analysis the signatures of anthropogenic activities within a specific landscape context and further demonstrate how people use of the land changes the Earth surface in three dimensions. Different utilisations of the same land cover show a different extent of anthropogenic impacts, underlining the opportunity for future analyses of the ‘magnitude’ and ‘type’ of human forcing on Earth. The results provide robust evidence of the human activity’s impact on some terrestrial surfaces, fostering therefore future studies linking the relationship between humans, land use, and geomorphological alterations. Our study offers a new insight to understand the present geomorphology coupling the function of human activities and poses a challenge for future research of the geomorphic and human systems in a world increasingly affected by anthropogenic activities.