Thank you for your positive feedback and the constructive comments. You will find a response to each part below, including the changes we would propose based on your comments.

INTRODUCTION

We will reorganise the introduction as you suggested, and we hope it will be clearer.

At the scale of the French Pyrenees, vectorised maps of forest cover were only available for the three dates we mentioned. An ongoing project aims to vectorise forest land use from the 1950s map and will be used in a near future in the French Alps, but is not yet available for the French Pyrenees. Concerning the most recent forest cover map, the next edition of the French forest map (BD Forêt® v3) is actually scheduled but for 2027 only.

METHODOLOGY

Thank you for this comment. We will improve the readability of this section in the next version of the manuscript, as you suggested.

Of the 25 m, 50 m and 100 m elevational bands that we tested, the 100 m band showed the most consistent distribution of forest cover according to elevation with the least background noise (Figure 2). Consequently, the 100 m elevation band was the most appropriate for estimating the elevation of the limit at a canopy cover threshold in the next step of the method. We already specified this point in the supplementary material, but we will mention it in the main text in the revised version.

In our manuscript, we distinguished between temporal and spatial analyses. To assess the temporal variations in temperature, only data from the Pic-du-Midi weather station allowed us to go back to the 19th century. The other weather stations in the Pyrenees had shorter time series and were located at lower elevations. To assess the spatial variations, we focused on recent data and used the Aurelhy database to extract climatic data. We tested the derived variables “mean temperature in the warmest month of summer”, “mean temperature in the coldest month of winter” and “total annual precipitation” for the period 1981-2010 and in the 300 m above the forest line (in 1851 or in 1993) in the linear models, but they were rarely retained in the most parsimonious models. During preliminary analyses, we also tested temperature changes between the periods 1961-1985 and 1986-2010 in each municipality, but these were not retained in the models either. Furthermore, when we compared temperature changes between the periods 1961-1990 and 1981-2010 from the homogenised series (seven stations in the study area) with those from the climatic models (Aurelhy, SAFRAN, CHELSA v2.1 and CHELSAcruts), we found no correlation. This shows that the spatial models are unable to accurately characterise temperature change in the study area of the French Pyrenees. We will clarify this in the new version.

Based on the data we have, we cannot assess the proportion of livestock grazing in the mountains or whether farming is extensive or intensive. However, we can distinguish between sheep and cattle. After investigation, we observed a decrease in the number of sheep and a slight increase in the number of cattle. This detailed temporal trend could thus be added to the next version of the manuscript. However, including change in livestock density for cattle and sheep rather than change in total livestock density did not significantly improve the models explaining spatial variations. Therefore, we will keep the previous variables and models.

RESULTS, MAPS, TABLES AND FIGURES

Thank you for all these comments, which will be used to revise the tables and figures in the next version of the preprint.

DISCUSSION

In the preprint, we argued that the higher rate of change observed in this study compared to previous studies is probably related to differences in the study area or methodology when working in the same area. We didn’t address the potential deficiencies of the État-Major map, as it is considered to be the most accurate document for working at the scale of an entire massif, and a good quality map for assessing forest cover in the 1850s (Dupouey et al., 2007). Moreover, the État-Major map was based on the Napoleonic cadastre where available, which enhanced its accuracy (Rochel et al., 2017). Indeed, this was probably the case in the French Pyrenees, where the cadastral data predated the État-Major map. Verifying the correspondence between the Napoleonic cadastre and the État-Major map would be valuable for completing the historical data, but this is a complex task for such a regional study area. Nevertheless, we agree that some of the discrepancy may be due to differences between the maps or potential deficiencies in the État-Major map. However, we believe that the methods developed to vectorise the historical map and estimate forest-line elevations minimised these discrepancies. Firstly, the standardisation of the maps that we have carried out should account for most of the differences between the products. Secondly, as we have already discussed, the definition has limitations, as the forest on the État-Major map probably falls somewhere between the IGN definitions of “forest” and “closed forest”. Thirdly, as for the État-Major map itself, a method of standardising vectorisation and georeferencing across France has been proposed and implemented in the French Pyrenees, thereby reducing discrepancies between regions, particularly for vectorisation and georeferencing (Favre et al., 2017). Indeed, the État-Major map was drawn by different people, so there may be inconsistencies between map tiles in the cartographic figures (Dupouey et al., 2007; Thomas et al., 2017). Furthermore, forest boundary inaccuracies are more prevalent in mountainous areas and can be related to vectorisation and georeferencing (Thomas et al., 2017). Last but not least, using a forest cover threshold to determine the position of the forest line reduces the risk of error for forest patches located at the top of municipalities, for example where scree is found today.  We will develop these considerations in the next version of the manuscript, as you suggested.

References:

Dupouey, J.-L., Bachacou, J., Cosserat-Mangeot, R., Aberdam, S., Vallauri, D., Chappart, G., and Corvisier-de-Villèle, M.-A.: Vers la réalisation d’une carte géoréférencée des forêts anciennes de France, Le Monde des Cartes, Paris: Comité français de cartographie, 85–98, 2007.

Favre, C., Grel, A., Granier, E., Cosserat-Mangeot, R., Bachacou, J., Leroy, N., and Dupouey, J.-L.: Digitalisation des cartes anciennes. Manuel pour la vectorisation de l’usage des sols et le géoréférencement des minutes 1:40 000 de la carte d’Etat-Major. Version 13.3, INRA, 79, 2017.

Rochel, X., Abadie, J., Avon, C., Bergès, L., Chauchard, S., Defever, S., Grel, A., Jeanmonod, J., Leroy, N., and Dupouey, J.-L.: Quelles sources cartographiques pour la définition des usages anciens du sol en France ?, Rev. For. Fr., 353, https://doi.org/10.4267/2042/67866, 2017.

Thomas, M., Bec, R., Abadie, J., Avon, C., Bergès, L., Grel, A., and Dupouey, J.-L.: Changements à long terme des paysages forestiers dans cinq parcs nationaux métropolitains et le futur parc national des forêts de Champagne et Bourgogne, Rev. For. Fr., 387, https://doi.org/10.4267/2042/67868, 2017.

SPECIFIC COMMENTS

Thank you for these relevant comments. All of the following answers will be incorporated in the next version of the manuscript.

L36: These studies have been conducted in Europe, more precisely in the Alps and Pyrenees, for these examples.

L41-44: I propose changing the second sentence to make it clearer as follows: “Recent upward shifts in the forest line have also been reported at large spatial scales thanks to aerial photographs (Améztegui et al., 2016; Gehrig-Fasel et al., 2007). However, only a limited response of the forest line, or a lag in the response, has been documented in several locations across Europe despite increasing temperatures (Gehrig-Fasel et al., 2007; Körner and Hiltbrunner, 2024; Paulsen et al., 2000).”

L49-50: Across Europe, historical maps can cover large areas, as the État-Major map (Kaim et al., 2016). However, the use of historical maps was usually restricted to small areas in comparison to entire massifs (eg. Egarter Vigl et al., 2016; Mainieri et al., 2020; Mietkiewicz et al., 2017; Tasser et al., 2007).

References:

Egarter Vigl, L., Schirpke, U., Tasser, E., and Tappeiner, U.: Linking long-term landscape dynamics to the multiple interactions among ecosystem services in the European Alps, Landscape Ecol, 31, 1903–1918, https://doi.org/10.1007/s10980-016-0389-3, 2016.

Kaim, D., Kozak, J., Kolecka, N., Ziółkowska, E., Ostafin, K., Ostapowicz, K., Gimmi, U., Munteanu, C., and Radeloff, V. C.: Broad scale forest cover reconstruction from historical topographic maps, Applied Geography, 67, 39–48, https://doi.org/10.1016/j.apgeog.2015.12.003, 2016.

Mainieri, R., Favillier, A., Lopez-Saez, J., Eckert, N., Zgheib, T., Morel, P., Saulnier, M., Peiry, J.-L., Stoffel, M., and Corona, C.: Impacts of land-cover changes on snow avalanche activity in the French Alps, Anthropocene, 30, 100244, https://doi.org/10.1016/j.ancene.2020.100244, 2020.

Mietkiewicz, N., Kulakowski, D., Rogan, J., and Bebi, P.: Long‐term change in sub‐alpine forest cover, tree line and species composition in the Swiss Alps, J Vegetation Science, 28, 951–964, https://doi.org/10.1111/jvs.12561, 2017.

Tasser, E., Walde, J., Tappeiner, U., Teutsch, A., and Noggler, W.: Land-use changes and natural reforestation in the Eastern Central Alps, Agriculture, Ecosystems & Environment, 118, 115–129, https://doi.org/10.1016/j.agee.2006.05.004, 2007.

L111: “The Pyrenees range stretches over 300 km between the Atlantic Ocean and the Mediterranean Sea, at the French-Spanish border, and is almost 100 km wide in its central part.”

L116-118: The values are for the French Pyrenees. This will be indicated. “The eastern region is under Mediterranean influence and is characterised by lower precipitation and a warmer average temperature than the western region, under oceanic influence: 1060 vs 2298 mm of average annual precipitation and 5.9 vs 5.3°C between 1958 and 2008, respectively for Cerdagne (eastern French region) and the Pays-Basque massifs (western French region) (Maris et al., 2009).”

L143-144: We did not have information on forest cover in the other countries. We believe that such information would have been valuable without altering the results. Indeed, I could not find a suitable forest map in Spain dating from the 19th century (Ruiz Del Castillo Y Navascues et al., 2006). Furthermore, in most cases, the boundary follows the ridge. Therefore, the probability that a small patch extend in the neighbouring countries is reduced, as the alpine belt is often reached. Moreover, in the case where the forest continued at higher elevation in the neighbouring country, the forest most likely reached the highest elevation in the municipality, this municipality was thus excluded from the analyses.

Reference:

Ruiz Del Castillo Y Navascues, J., López Leiva, C., García Viñas, J. I., Villares Muyo, J. M., Tostado Rivera, P., and García Rodríguez, C.: The Forest Map of Spain 1:200,000. Methodology and analysis of general results, For. syst., 15, S24–S39, https://doi.org/10.5424/srf/200615S1-00979, 2006.

L167: “For each municipality in the study area, we rasterised the forest present in each vector map (Fig. 2a).”

L222: The Daubrée forest inventory only reports forest area for each municipality but this historical source is not at all a land-use map. Therefore, we only use it to estimate changes in forest cover for the entire study area. Nevertheless, we will properly cite this source of information in M&M.

L247: Including livestock density data from 1838 allows us to visualise changes in the temporal trend in the study area since the 1850s: an increase before 1850, followed by a decrease. Changes in livestock density between 1852 and 2000 are used in the models for this period to assess spatial variations in the study area. Here again, the distinction between temporal and spatial analysis will be clarified.

L255: The models were built with the lme function of the R package nlme (Pinheiro et al., 1999). The response variables were (1) forest-line shift velocity (between 1851 and 1993, between 1993 and 2010 and between 1851 and 2010), (2) the differences in forest-line shift velocity between the two periods (1851-1993 and 1993-2010), (3) the differences in shift velocity between the forest line and closed forest line, therefore accounting for municipality pairing. These three models were built using only the “arrondissement” as a random effect, with no fixed effects, to test if the differences were significant. While they are similar to t-tests, they allowed us to account for the “arrondissement” effect, i.e. spatial autocorrelation. This was probably unclear in the preprint, but we will clarify this in the next version. For the last test (4) the forest-line shift velocity between 1851 and 1993 and between 1993 and 2010 were the response variables and the effect of the spatial group was tested (west vs east, separated at an RGF93 longitude of 520 km).

Reference:

Pinheiro, J., Bates, D., and R Core Team: nlme: Linear and Nonlinear Mixed Effects Models, https://doi.org/10.32614/CRAN.package.nlme, 1999.

L284: I chose the threshold of 1.75 to simplify the selection of the most parsimonious models, by comparing fewer extracted models while retaining those with only small variations in the adjusted R² and with at least two variables. We will consider using a threshold of 2 in the revised manuscript and adjusting the results accordingly.

L355: A lower mean summer water balance means a lower quantity of water available for plants, after accounting for evaporation and transpiration, i.e. drier conditions. We will specify this in the next version of the manuscript.

L400-408: The data I used are monthly mean temperatures. Therefore, we cannot check the daily temperatures to ascertain whether the conditions described by Paulsen and Körner (2014) are met: a growing season of at least 94 days, with a daily mean temperature of at least 0.9 °C and a mean temperature of 6.4 °C over all these days. Still, we can calculate the mean temperature for the summer months (July, August, September), when the mean temperature usually exceeds 6.4 °C on average. This produces a similar trend to that of the annual mean temperature: an increase in growing season temperatures of 0.75°C between 1851 and 1993, resulting in a theoretical shift in the forest line of 136 m (and a velocity of 1.0 m.yr^-1); and an increase of 0.70°C between 1993 and 2010, resulting in a theoretical shift of 126 m (and a velocity of 7.4 m.yr^‑1).

Reference:

Paulsen, J. and Körner, C.: A climate-based model to predict potential treeline position around the globe, Alp Botany, 124, 1–12, https://doi.org/10.1007/s00035-014-0124-0, 2014.

Thank you for your review and valuable comments. We will emphasise the novelty of the methodology in the discussion section of the next version of the manuscript.

Please, find a detailed answer to each point below.

We are aware that the two periods under comparison are quite different. As you pointed out, the long period between the 1850s and the 1990s probably includes periods of stagnation and upward shifts, and the upward shift during such a long period is unlike to be linear. Despite the limiting availability of dates, we found it valuable to compare the periods in order to identify recent variations in forest-line dynamics trends. Moreover, calculating shift velocities enabled us to make this comparison despite the different lengths of the two periods. Furthermore, previous studies have found evidence of an upward shift in forest line in the Spanish Pyrenees early after the 1850s (Camarero and Gutiérrez, 2004). As land was abandoned early in the French Pyrenees, at least on the eastern side, the forest line could have started to shift upward early too. However, studies in the Alps found later forest-line upward shifts, starting in the 1950s, and this pattern may also be observed in the French Pyrenees. Adding a supplementary date would be very valuable and may be possible in the near future with the preparation of a digital version of the 1950 French forest map. Conversely, the recent period between the 1990s and 2010 is very short. Nevertheless, we were able to capture forest-line dynamics and spatial variations. Therefore, we believe that these two periods are still relevant to this study, and that the comparison is valuable when all of this is taken into account. We will address these issues and the limitations of this comparison in the next version of the manuscript.

Reference:

Camarero, J. J. and Gutiérrez, E.: Pace and pattern of recent treeline dynamics: response of ecotones to climatic variability in the Spanish Pyrenees, Climatic Change, 63, 181–200, https://doi.org/10.1023/B:CLIM.0000018507.71343.46, 2004.

The original État-Major map has a scale of 1:40,000. BD Forêt® V1 and V2 were created from the interpretation of aerial photographs at scales ranging from 1:17,000 to 1:25,000. We will mention this in the next version of the manuscript. As you underlined, the minimum mapping unit differs between the three maps. That’s why, in a first step, we have downgraded the accuracy of the État-Major map and BD Forêt® v2 at the level of the BD Forêt® v1. Thus, these differences should not impact the estimations of forest-line dynamics. We will clarify this point in the next version.

The Daubrée forest inventory data only included forest area, with no map or precise spatial indications. This prevented us from mapping forest cover and estimating forest-line elevation at this date. However, we believe that using the forest area information for comparison is still interesting. Figure 3 therefore shows four points for forest area (light green triangles) and three for forest-line elevation (dark green points). We will clarify the presentation of this additional data source and the figure 3 in the next version of the manuscript.

During the preliminary analyses, we examined the maximum forest elevation per municipality, which should correspond to your proposal for the position of the local maximum, albeit at a larger scale. However, the maximum forest elevation was more susceptible to potential mapping errors. Moreover, although present, the differences between the two limits (5% forest cover in the highest elevation band of 100 m and maximum forest limit) were still quite small in most municipalities thus insufficient to exclude the impact of human activity. Therefore, we don’t think this additional analysis will be relevant in the case of the French Pyrenees.

Potential treeline elevation modelling relies on climatic data. While it is possible to do this for recent times, spatially consistent climatic reconstructions going back to the 1850s are lacking. This prevents us from estimating the deviation from the potential treeline elevation at the beginning of the study period. As Ameztégui et al. (2016) previously proposed using forest-line elevation as a proxy for human impact, and this allowed quantification even in the 1850s, we have chosen to use a similar indicator in our study. Furthermore, we believe that estimating the potential treeline elevation would be more relevant at a finer scale than for the entire municipality, although we agree on the usefulness of this approach for answering these questions.

Reference:

Améztegui, A., Coll, L., Brotons, L., and Ninot, J. M.: Land‐use legacies rather than climate change are driving the recent upward shift of the mountain tree line in the Pyrenees, Global Ecology and Biogeography, 25, 263–273, https://doi.org/10.1111/geb.12407, 2016.

These relevant references will be added to the next version of the preprint.

Answers to minor comments

L62: We used the term 'potential treeline' to clarify the target object, considering the possibility of a 'realised treeline'. This is because the treeline is not considered a potential line as a consensus in the literature.

L201: We chose 300 m because it is small enough to be representative of conditions where forest-line shifts should occur, but large enough to include enough points for a meaningful average.

L265: We also checked the VIF, which resulted in a similar selection of variables. However, comparing the correlations enabled us to test several variable combinations in the preliminary analysis. We therefore continued with this approach for the final set selection.

Chapter 4.2: This was not so evident to us from the literature, as the effects of climate and land-use changes probably interact. Evidence of the impact of both climate change (Hagedorn et al., 2014; Harsch et al., 2009) and land-use change (Améztegui et al., 2016; Gehrig-Fasel et al., 2007) on forest-line dynamics has previously been observed. Moreover, whether or not it is a coincidence, the similar pace of forest-line shift and temperature increase does not allow us to reject the hypothesis that climate change is a driver. Conversely, the fact that the forest line is shifting more rapidly than the temperature is warming suggests that climate change is not the only driver. We will clarify this idea in the revised version of the manuscript.

References:

Améztegui, A., Coll, L., Brotons, L., and Ninot, J. M.: Land‐use legacies rather than climate change are driving the recent upward shift of the mountain tree line in the Pyrenees, Global Ecology and Biogeography, 25, 263–273, https://doi.org/10.1111/geb.12407, 2016.

Gehrig-Fasel, J., Guisan, A., and Zimmermann, N. E.: Tree line shifts in the Swiss Alps: climate change or land abandonment?, Journal of vegetation science, 18, 571–582, 2007.

Hagedorn, F., Shiyatov, S. G., Mazepa, V. S., Devi, N. M., Grigor’ev, A. A., Bartysh, A. A., Fomin, V. V., Kapralov, D. S., Terent’ev, M., Bugman, H., Rigling, A., and Moiseev, P. A.: Treeline advances along the Urals mountain range – driven by improved winter conditions?, Global Change Biology, 20, 3530–3543, https://doi.org/10.1111/gcb.12613, 2014.

Harsch, M. A., Hulme, P. E., McGlone, M. S., and Duncan, R. P.: Are treelines advancing? A global meta‐analysis of treeline response to climate warming, Ecology Letters, 12, 1040–1049, https://doi.org/10.1111/j.1461-0248.2009.01355.x, 2009.