Land use land cover dynamics through time and their proximate drivers of change in a tropical mountain system: a case study in a highland landscape of northern Ecuador

Tropical mountain ecosystems are threatened by land use pressures, reducing the capacity of ecosystems to provide a large diversity of benefits to people and to be able to achieve them in the long term. The analysis of land use pressures is often superficial and very general, although they are characterized by numerous interactions and strong differences in their local dynamics. We used a variety of freely available geospatial and temporal data and methods to assess and explain patterns of land use land cover (LULC) change, focusing on native ecosystem dynamics, in a sensitive region of the northern Ecuadorian Andes. Our results demonstrate a dynamic and clear geographical pattern of distinct LULC transitions through time, explained by different combination of socio-economic factors, pressure variables and environmental parameters, from which ecological context variables, such as slope and elevation, were the main drivers of change in this landscape. We found that deforestation of remnant native forest and agricultural expansion still occur in higher elevations located, while land conversion toward anthropic environments were observed in lower elevations to the east of the studied territory. Our findings also reveal an unexpected stability trend of paramo and a successional recovery of previous agricultural land to the west and center of the territory which could be explained by agricultural land abandonment. However, the very low probability of persistence of montane forests in most of the studied landscape, highlights the risk that the remnant montane forests will be permanently lost in a few years, posing a greater threat to the already vulnerable biodiversity and limiting the capacity ecosystem service provisioning. The dynamic patterns through space and time and their explanatory drivers, found in our study, could help improve sustainably resource land management in vulnerable landscapes such as the tropical Andes in northern Ecuador.

9 164 groupings take into consideration the mean medium of relief surface roughness, type of forests 165 and the presence of urban areas [46]. Finally, the LULC classification for each year was layered 166 over both (1) the reclassified elevation map, and the (2) reclassified administrative map. Spatial 167 data assimilation, processing and overlaying analysis were conducted in the R environment [48].
168 Drivers of change 169 To understand what predictors could explain the LULC dynamics we tested a set of five group of 170 variables that have previously been reported as possible drivers of change [7,12] and were 171 described in the conceptual framework that interconnects: Driving forces -Pressure -State -172 Benefits / Impacts-Response adopted by the European Environment Agency [25]; a scheme also 173 used by the Ministry of the Environment of Ecuador as a tool to guide the formulation and 174 adjustment of policies to foster biodiversity conservation in Ecuador [49]. Within this approach, 175 we compiled a dataset of 13 variables of LULCC ranging from (1) socio-economic, (2) 176 topography, (3) anthropic pressures to natural ecosystems, (4), climate and (6) governance 177 decisions toward landscape development (Table 1). 178 In order to increase the number of units of analysis within parishes, all these variables were 179 obtained at the spatial resolution of census area [50]. After all the spatial data assimilation, 180 processing and visualisation necessary to obtain the drivers at the spatial unit of analysis, we 181 carried out a reduction dimension procedure using Principal Component Analysis [51] of the 182 drivers of change to summarize the distinct variables within the grouping drivers, all these 183 procedures were completed using R software (version 3.2.3) [48]. Correlated variables were 184 screened for the total variation explained by the first principal axes, and used to remove 185 correlated variables [52]. Coordinates of the principal components that accounted for more than 186 60% of the variation were then used as explanatory variables in a subsequent statistical model.  (Table 2), there was a 40% and 16% decrease of native forest and 227 paramo cover when compared the first and last periods of study (Table 2); but, by the last period 228 of study, areas of paramo still represent an important part (13%) of the study territory. Natural 229 water bodies (lakes and rivers) showed high persistence over time (Table 2).
230 Developed areas and floricultural crops were continuously increasing over time, and although 231 they were poorly represented in the first period of analysis (less than 0.4% in 1990), by 2014 232 they represented almost 5% of the study area (

Transitions of native ecosystems
253 In general, as expected the stability of native forests is decreasing through time in the entire 254 territory (Fig 3), with the exception of the western parish where the probability of remaining in 255 this LULC class increases through time probably due to agricultural land abandonment (Fig 3).
256 In contrast, areas located in the east tend to have lower values of stability through time and 257 higher probabilities to change to paramo and agricultural land; this pattern was more evident in 258 the last period evaluated (2008-2014) (Fig 4). Additionally, this trend is more evident along 259 elevation bands; where native forests located above 3300 m showed a lower probability of 260 remaining as forest along the years (Fig 4) and in the 2800-3300 altitudinal belt there is a high 261 probability of converting native to planted forests, especially in the center of the territory (Fig 4).  Transitions to anthropic environments 281 Developed areas demonstrate a differential trend through time in the study area (Fig 5). In the 283 period of evaluation (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008) and significantly increased again in the last time period (2008-284 2014), in contrast, the parishes located to the east exhibit a more stable probability of remaining 285 as developed areas through time probably due to the proximity of the larger towns (Fig 5). Since 286 the territory studied is in general a rural area, there is a dynamic trend towards converting 287 agricultural to urban areas, which follow a geographic pattern (Fig 5).   (Table 2), which also show a contrasting pattern in their response 340 variables, in such a way that when was prevalent an increase in agricultural areas it was a 341 decrease in shrub and herb extension (Fig 2); these models depicted the following grouping 342 drivers as significant parameters (p<0.001): pressure drivers (PC1) and a variable that describes 343 differences development of the distinct administrative areas within the study area (Table 2).  378 We estimated a paramo loss of 16% from 1990 to 2014 in Pedro Moncayo county (Table 1)