Repeat Oblique Photography Shows Terrain and Fire-Exposure Controls on Century-Scale Canopy Cover Change in the Alpine Treeline Ecotone

Alpine Treeline Ecotone (ATE), the typically gradual transition zone between closed canopy forest and alpine tundra vegetation in mountain regions, displays an elevational range that is generally constrained by thermal deficits. At landscape scales, precipitation and moisture regimes can suppress ATE elevation below thermal limits, causing variability in ATE position. Recent studies have investigated the relative effects of hydroclimatic variables on ATE position at multiple scales, but less attention has been given to interactions between hydroclimatic variables and disturbance agents, such as fire. Advances in monoplotting have enabled the extraction of canopy cover information from oblique photography. Using airborne lidar, and repeat photography from the Mountain Legacy Project, we observed canopy cover change in West Castle Watershed (Alberta, Canada; ~103 km2; 49.3◦ N, 114.4◦ W) over a 92-year period (1914–2006). Two wildfires, occurring 1934 and 1936, provided an opportunity to compare topographic patterns of mortality and succession in the ATE, while factoring by exposure to fire. Aspect was a strong predictor of mortality and succession. Fire-exposed areas accounted for 83.6% of all mortality, with 72.1% of mortality occurring on southand east-facing slope aspects. Succession was balanced between fire-exposed and unburned areas, with 62.0% of all succession occurring on northand east-facing slope aspects. The mean elevation increase in closed canopy forest (i.e., the lower boundary of ATE) on northand east-facing undisturbed slopes was estimated to be 0.44 m per year, or ~44 m per century. The observed retardation of treeline advance on south-facing slopes is likely due to moisture limitation.


Introduction
Alpine treeline ecotone (ATE), the transition zone between closed canopy forest and alpine tundra, occurs where tree growth is limited by thermal thresholds at high elevations [1,2]. At continental scales, the approximate elevation of ATE is predicted by the growth limitation hypothesis, which suggests an inverse relationship exists between latitude and ATE elevation [1] as a function of temperature in the growing season [2].
A meta-analysis [13] identifies 103 publications which describe long term records of treeline position (both ATE and latitudinal treeline) at 166 sites around the world, between 1900-2008.
Approximately 52% of these sites have some form of treeline advance. Two studies report treeline recession and both cases are associated with disturbances. Linear regression with climate records reveals that many of the advances in ATE position in [13] are associated with increases in atmospheric temperature.
Modulating factors, which suppress ATE below thermal limits regardless of rising atmospheric temperature, are often studied in isolation, and the potential for interaction between modulating factors is rarely considered. Our research investigated the interaction between two of these modulating factors, fire disturbance and moisture limitation, and their combined influence on ATE position in a spruce-fir forest in the Rocky Mountains. Using a novel monoplotting technique to extract quantitative spatial data from a historic photographic record, we investigated whether postfire moisture limitation caused topographic patterns of vegetative regeneration, thereby selectively suppressing ATE in dry, fire-exposed environments.

Fire and Topography in the ATE
In the Rocky Mountains, fires are less frequent in the ATE than at lower elevations, and fire return intervals (FRIs) are positively correlated with elevation [14,15]. While relatively infrequent, fire remains an important disturbance agent in the ATE [16,17], and the effects of fire on subalpine regions are potentially increasing with climate change [18].
Topographic factors like aspect and slope influence many of the processes that affect fire behavior in mountain landscapes. Topography can alter precipitation regimes, species composition, lighting strike frequency, and drying caused by both solar insolation and wind [14]. For example, in North America, south-facing mountain aspects receive more solar insolation and are generally drier than north-facing aspects. This impacts fire density (i.e., number of fires per unit area), with greater density observed on south-facing aspects [19][20][21]. However, the increase in fire density on warm aspects does not translate to an increase in area burned on south-facing aspects [22,23]. Thus, small fires are ignited more frequently on south-facing aspects, but a large fire will burn a range of topographic conditions, without aspect bias.
Patterns of moisture limitation also have the potential to impact postfire regrowth, potentially suppressing ATE elevation below levels that would be expected if ATE was only limited by atmospheric temperature. At high elevations in the Rocky Mountains, postfire areas are predominantly colonized by spruce or fir. The initial postfire colonization in these areas can take 20-100+ years [24][25][26][27], leaving alpine meadows that persist for a century or more. Postfire rates of succession in subalpine forest can be twice as high on mesic slopes as on xeric slopes [28], and growth of subalpine fir is limited on south-facing aspects [29,30]. If delayed regrowth can inhibit colonization for almost a century, then the present position of fire-exposed treelines on south-facing aspects may be temporarily suppressed, appearing lower than would be expected if ATE were solely limited by atmospheric temperature. While stand-age reconstructions are useful at identifying the processes involved in postfire regrowth, few direct observations of an aspect effect on postfire regrowth have been made using imagery.

ATE Observation Techniques
Determining the interaction between aspect and regrowth in postfire environments using imagery requires high resolution observations that extend more than a century. These are uncommon. A trade-off exists between spatial resolution and temporal extent in ATE observation methods [31]. Multitemporal remote sensing can monitor vegetative land cover change at spatial resolutions required to correlate with topographic variables (10-100 m) [17], but observations of ATE using passive optical satellite imagery, like Landsat TM [11,32,33] only extend back to the 1970-80s. Aerial photographs are similarly useful for treeline change monitoring [17,34], but aerial photography records only extend to the late 1940s and early 1950s for much of North America.
Conversely, ATE observation records with very long temporal extents often have poor spatial resolution, prohibiting topographic analysis. For example, palynological records extend thousands of years before present, with some exceeding 10,000 years [35], but these observations can only provide a spatial resolution equivalent to the pollen distribution distance (<100 m). Dendrochronological records have offered the best combination of high spatial resolution (<10 m) and high temporal extent (~1000 yr). Tree-ring records are often used in ATE studies [9,[36][37][38]. The spatial resolution and temporal extent of dendrochronological studies are both high, but the spatial extent is limited by onerous collection methods and the large sample size required to cover enough sites for topographic analysis [9,15].
Repeat photography can address the need for ATE observations that provide high spatial extent and resolution, and long temporal extent. This method is widely used to assess land cover change in the ATE. Patterns of ATE advance in eastern slope regions of the Canadian Rockies have been described using repeat photography [39][40][41]. High resolution photographs were used to identify shrub advance in the Swedish Scandes [42]. Repeat photography has also been used for a description of land cover change in the Ural Mountains [43]. While the temporal extent of repeat photographs, 100 years or more in cases, is useful for studies of change in the ATE, the method has been mostly restricted to qualitative analysis, given the previous inability to standardize spatial scale across oblique imagery. Quantitative spatial analysis of ATE has been attempted by draping a polygon fishnet over repeat photographs [44], in order to standardize the observation units in both historical and contemporary photographs, but spatial scale still varied between polygons, limiting the usefulness of this technique for spatial analysis.
Recent technological advances have enabled the fully quantitative spatial analysis of repeat photographs, using the WSL Monoplotting Tool (WSL-MT, note: WSL is the German acronym for "Eidgenössische Forschungsanstalt für Wald, Schnee und Landschaft", the Swiss Federal Institute for Forest, Snow and Landscape Research -see Supplementary Material for further information on WSL-MT), which produces georeferenced vector data from oblique photographs using tie points between oblique and aerial photographs, and high resolution topographic data [45]. The WSL-MT tool has been used to monitor land cover change by extracting spatial data from oblique repeat photographs [46] and to observe century-scale vegetation transitions in the Rocky Mountains of southern Alberta, Canada [47]. These experiments were conducted on photo pairs from the Mountain Legacy Project, a collection of over 120,000 historic survey images of the Canadian Rockies (1888-1958), of which over 6000 have 21 st -century repeat photographs available under creative commons license [48]. Previous research, validated using airborne lidar, established that fractional cover classes could be assessed from oblique imagery at a spatial resolution of 20 m, using the WSL-MT [49].

Research Goals
We used a technique to assess fractional cover from oblique photographs, described in [49], to observe postfire vegetative change in repeat photographs of a watershed in the Canadian Rockies. We sought to determine: (1) if there were observable differences in canopy cover change between fire-exposed and non-fire-exposed regions of the watershed, (2) if changes in canopy cover correlate to topographic patterns that are consistent with expected edaphic controls on regrowth in the region (e.g., less regrowth on south-facing aspects), and (3) if it was possible to detect the rate of ATE advance in non-fire-exposed areas, where ATE position is regulated by climate.

Study Area
The study was conducted over the West Castle Watershed (WCW), Alberta, Canada (49.3 • N, 114.4 • W). WCW is in the headwaters of the Oldman River, on the eastern slopes of the Canadian Rockies, with an area of~103 km 2 , and an elevation range of~1400-2600 m a.s.l. Aside from a small ski resort and village near the downstream end of the WCW (2.9% of the watershed area), and trails in the valley Remote Sens. 2020, 12, 1569 4 of 21 bottoms, the forested slopes demonstrate limited anthropogenic disturbance. Consequently, the WCW is an ideal study area for assessing recent natural shifts in the ATE of this part of the Canadian Rockies. Dominant species at ATE elevation in WCW are subalpine fir and Engelmann spruce.

Repeat Photographs
Repeat photographs of WCW were provided by the Mountain Legacy Project [48] (Figure 1). The historical images were taken in the summer of 1914 by Morrison Parsons Bridgeland. The repeat photographs were taken by the Mountain Legacy Project in July of 2006 (Table 1). Photograph pairs were selected to maximize the spatial coverage of the observed area of WCW. Using the same set of photographs, it was determined that 38.5 km 2 (37.7%) of the watershed was observable in these images [49] (Figure 2). Remote Sens. 2020, 12, x FOR PEER REVIEW 4 of 22

Repeat Photographs
Repeat photographs of WCW were provided by the Mountain Legacy Project [48] ( Figure 1). The historical images were taken in the summer of 1914 by Morrison Parsons Bridgeland. The repeat photographs were taken by the Mountain Legacy Project in July of 2006 (Table 1). Photograph pairs were selected to maximize the spatial coverage of the observed area of WCW. Using the same set of photographs, it was determined that 38.5 km 2 (37.7%) of the watershed was observable in these images [49] (Figure 2).  Table 1 for further details.   Table 1 for further details.  Figure 2. Visible region of WCW in each of the seven MLP photographs. In cases where the same region is visible in multiple photographs, the photograph used for analysis was determined by criteria described in [49].

Canopy Classification
Canopy cover was extracted from oblique photography with the WSL monoplotting tool [45], using a method described in [49]. The WSL-MT performs oblique to orthogonal vector transformation with a camera calibration, which generated by pairing control points between oblique photographs and a high resolution DEM, co-registered with aerial imagery. For the DEM requirement, a 1 m DEM was used, which was collected from an airborne lidar survey in October 2014. For the aerial imagery requirement, a SPOT 6 image of WCW was used (1.5 m, acquired July 31, 2014, © 2014 CNES -see Supplementary Material). In the absence of high resolution topographic and aerial data from 1914, the analysis relied on the camera calibration from the 2014 lidar DEM and 2006 SPOT image being applied to the 1914 oblique photographs. This approach was deemed reasonable, as no extreme changes in topography (e.g., rock slides, slumping) are apparent between the 1914 and 2006 datasets.
The monoplotting procedure results in a fishnet vector layer, where each cell represents a 20 × 20 m area (viewed orthogonally). This vector layer is draped over the oblique raster image, such that cells are scaled according to the perspective in the image (i.e., cells in the foreground appear larger, Figure 2. Visible region of WCW in each of the seven MLP photographs. In cases where the same region is visible in multiple photographs, the photograph used for analysis was determined by criteria described in [49].

Canopy Classification
Canopy cover was extracted from oblique photography with the WSL monoplotting tool [45], using a method described in [49]. The WSL-MT performs oblique to orthogonal vector transformation with a camera calibration, which generated by pairing control points between oblique photographs and a high resolution DEM, co-registered with aerial imagery. For the DEM requirement, a 1 m DEM was used, which was collected from an airborne lidar survey in October 2014. For the aerial imagery requirement, a SPOT 6 image of WCW was used (1.5 m, acquired July 31, 2014, © 2014 CNES -see Supplementary Material). In the absence of high resolution topographic and aerial data from 1914, the analysis relied on the camera calibration from the 2014 lidar DEM and 2006 SPOT image being applied to the 1914 oblique photographs. This approach was deemed reasonable, as no extreme changes in topography (e.g., rock slides, slumping) are apparent between the 1914 and 2006 datasets.
The monoplotting procedure results in a fishnet vector layer, where each cell represents a 20 × 20 m area (viewed orthogonally). This vector layer is draped over the oblique raster image, such that cells are scaled according to the perspective in the image (i.e., cells in the foreground appear larger, and cells in the background appear smaller) ( Figure 3). This method circumvents the problem of varying spatial scale in the oblique image, allowing quantitative spatial analysis to be performed on the resulting raster. and cells in the background appear smaller) ( Figure 3). This method circumvents the problem of varying spatial scale in the oblique image, allowing quantitative spatial analysis to be performed on the resulting raster. Calibrated fishnets were applied to each of the seven oblique images. Any fishnet grid cell that intersected with less than 75% of the viewshed (i.e., < 300 m 2 ) was omitted from the analysis (note: because these data were eventually rasterized, and each raster represented an area between 300-400 m 2 , canopy cover and change is reported by grid cell count rather than surface area).
Once in an oblique projection, grid cells were manually assigned into one of six canopy cover classes, based on observed canopy openness and texture: 1) No Cover-grid cells devoid of vegetation; 2) Low Vegetation-grid cells appear vegetated, but context and texture suggest shrubs or krummholz, not upright trees > 2 m; 3) Partial Canopy-trees are present, but ground is visible in > 50% of the grid cell; 4) Full Canopy-trees cover > 50% of a grid cell; 5) Snow-snow covers > 50% of a grid cell; 6) Structure-any anthropogenic structure or non-vegetated land cover is present in the grid cell (e.g., buildings, roads, and trails) ( Figure 4). Calibrated fishnets were applied to each of the seven oblique images. Any fishnet grid cell that intersected with less than 75% of the viewshed (i.e., < 300 m 2 ) was omitted from the analysis (note: because these data were eventually rasterized, and each raster represented an area between 300-400 m 2 , canopy cover and change is reported by grid cell count rather than surface area).
Once in an oblique projection, grid cells were manually assigned into one of six canopy cover classes, based on observed canopy openness and texture: (1) No Cover-grid cells devoid of vegetation; (2) Low Vegetation-grid cells appear vegetated, but context and texture suggest shrubs or krummholz, not upright trees > 2 m; (3) Partial Canopy-trees are present, but ground is visible in > 50% of the grid cell; (4) Full Canopy-trees cover > 50% of a grid cell; (5) Snow-snow covers > 50% of a grid cell; (6) Structure-any anthropogenic structure or non-vegetated land cover is present in the grid cell (e.g., buildings, roads, and trails) (  fires are provided for reference [50]. fires are provided for reference [50]. Given the focus on vegetative land cover change, further analysis focused on the first four, ordinal classes of vegetation (no cover through full canopy). Any grid cells having snow or structure in either 1914 or 2006 were omitted from analysis. The snow class was removed due to uncertainty of vegetation beneath the snow. The anthropogenic structure class was omitted because the processes causing vegetative land cover to change were distinct from the processes being studied. Of the total observed area, 7.2% was affected by the removal of snow cover, and 2.9% by the removal of anthropogenic structures.

Anthropogenic Disturbance
The Alberta Vegetation Inventory [51] was used to identify and omit areas of anthropogenic disturbance (see Supplementary Material). It was produced from interpretation of aerial photographs, and identifies several types of former and current disturbance, including former oil well pads, historical cut-blocks, settlements (including the Castle Mountain Ski Resort), and roads. These areas were aggregated into a single layer of anthropogenic disturbance so that topographic and other factors which affected land cover change could be considered independently of human activity.

Fire Disturbance
Two major historical fires occurred in the interval between the oblique photograph observations, one in 1934 and one in 1936. The extents of the historical fires were delineated in 2005, using a combination of visual evidence from a 1949 aerial photograph survey and descriptions from historical fire reports [50] (Figure 4). More detailed attributes of these fires, such as intensity and cause, are unknown, but previous research has demonstrated that fire intensity does not co-vary with topographic features such as aspect and slope in subalpine fir-spruce forests of the Rockies, as fires that reach these elevations are typically high intensity, crown replacing fires [22,23]. Thus, if we observe topographic differences in patterns of postfire regrowth, they are likely related environmental factors that affect regrowth, and not likely caused by topographic differences in fire intensity.
Fire suppression in early 20th-century was common practice, yet it is known that the Government of Alberta did not engage in fire suppression in the 1930s and 1940s, for reasons related to the Great Depression and the Second World War [52]. Fire suppression techniques did not affect the fire return interval of subalpine forests in the region [52]. Mean fire return intervals in the WCW range from 19.3 years near the mouth of the valley, to 70.8 years in the southern headwaters [53].

Topographic Analysis of Change
The four ordinal classes of canopy cover (1-no cover, 2-low vegetation, 3-partial canopy, 4-full canopy) were used to categorize vegetative change between 1914 and 2006. For each 20 x 20 m grid cell, the difference between the 2006 cover class and the 1914 cover class was used to generate a change class ( Figure 5). For example, if canopy cover in 2006 was no cover (1), and canopy cover in 1914 was full canopy (4), then the change class for that area was 1 − 4 = −3. The seven possible change classes were: −3 (high mortality), −2 (medium mortality), −1 (low mortality), 0 (no change), +1 (low succession), +2 (medium succession), +3 (high succession). Data were factored by fire exposure, with the 1934 and 1936 fire extents aggregated into a single layer. The distribution of topographic variables within each change cover class was contrasted between fire-exposed and non-fire-exposed areas in the time interval .
For aspect, we tested if the proportional distribution of aspects within each change class differed from expected aspect distribution in a way that suggested aspect influence on canopy cover. The 20 m DEM was sampled to establish the expected proportional distribution of aspect classes in the WCW. The observed proportional distribution of aspect within each change class was compared to the Data were factored by fire exposure, with the 1934 and 1936 fire extents aggregated into a single layer. The distribution of topographic variables within each change cover class was contrasted between fire-exposed and non-fire-exposed areas in the time interval .
For aspect, we tested if the proportional distribution of aspects within each change class differed from expected aspect distribution in a way that suggested aspect influence on canopy cover. The 20 m DEM was sampled to establish the expected proportional distribution of aspect classes in the WCW. The observed proportional distribution of aspect within each change class was compared to the expected distribution with a chi-square goodness-of-fit test. To test the significance to each comparison, a Monte Carlo analysis using 100 repetitions was used to build a sampling distribution.
For elevation and slope, we tested if the distribution of these variables differed in between comparisons of change class. One-way ANOVA tested whether any of the seven classes' mean values differed from the others. For both elevation and slope comparisons, Welch's t-test was used on each of the 20 pairwise combinations of seven change classes to establish which classes differed from each other. Note that Welch's t-test corrects for both unequal variance between distributions and unequal sample size [54]. A Bonferroni correction was applied to the t-tests to account for multiple comparisons. All topographic analyses were completed solely outside of the anthropogenic disturbance area.
Finally, the rate of upslope ATE advance in undisturbed areas was estimated along five transects of canopy change cells. Along each transect, the boundary of closed canopy forest was identified in the 1914 and 2006 scenes. The closed canopy forest boundary is the lower extent of ATE (sometimes referred to as timberline). Given that the upper boundary (species limit) could vary within a 20 m grid cell, and may have been obstructed by snow, the closed canopy forest boundary was used as a proxy for total ATE change. Transects were selected where a distinct timberline (lower ATE boundary) could be observed across at least five grid cells (100 m), in areas without orographic boundaries to treeline advance (i.e., climatic treelines), and where there was no exposure to fire between observations. Five of these transects were identified. Elevation was sampled from each of the 20 m grid cells, and the mean elevations of the closed canopy fronts in the 1914 and 2006 scenes were compared.

Aspect
Observed distribution of aspect significantly differed from expected distribution, with the majority of mortality seen on fire-exposed, warm aspects ( Figure 6). In total, fire affected 63.2% of the observed area in WCW and 83.6% of all mortality occurred in fire-exposed areas. A strong interaction between aspect and fire was observed, with 72.0% of all mortality occurring on SW, S, SE, and E aspects. There was a high occurrence of succession classes in N and NE aspects, in both fire-exposed and non-fire-exposed areas. Succession was also notably high in fire-exposed SW aspects, and non-fire-exposed NE aspects. In total, north-and east-facing slope aspects accounted for 62.0% of all succession. Mortality is largely restricted to fire exposed areas, and is predominantly in the SW, SE, and E aspects.
The expected proportion distribution of aspect in WCW varied, as a function of the orientation of the valley; SW, E, and NE aspect had higher occurrence, while S, SE, and NW aspects had lower occurrence (Figure 7). We compared this expected proportion of aspects for each of the seven change classes (8 aspects x 7 change classes = 56 comparisons). The results of the chi-squared test show that significant differences between expected and observed aspect proportion occurred frequently; 46 out of 56 (82.1%) of the change class/aspect cases showed significant differences ( Table 2). The greatest departures from expected proportions were seen in high occurrences of mortality on S and SE aspects, Mortality is largely restricted to fire exposed areas, and is predominantly in the SW, SE, and E aspects.
The expected proportion distribution of aspect in WCW varied, as a function of the orientation of the valley; SW, E, and NE aspect had higher occurrence, while S, SE, and NW aspects had lower occurrence (Figure 7). We compared this expected proportion of aspects for each of the seven change classes (8 aspects x 7 change classes = 56 comparisons). The results of the chi-squared test show that significant differences between expected and observed aspect proportion occurred frequently; 46 out of 56 (82.1%) of the change class/aspect cases showed significant differences ( Table 2). The greatest departures from expected proportions were seen in high occurrences of mortality on S and SE aspects, and high occurrences of succession on N and NE aspects. Mortality is largely restricted to fire exposed areas, and is predominantly in the SW, SE, and E aspects.
The expected proportion distribution of aspect in WCW varied, as a function of the orientation of the valley; SW, E, and NE aspect had higher occurrence, while S, SE, and NW aspects had lower occurrence (Figure 7). We compared this expected proportion of aspects for each of the seven change classes (8 aspects x 7 change classes = 56 comparisons). The results of the chi-squared test show that significant differences between expected and observed aspect proportion occurred frequently; 46 out of 56 (82.1%) of the change class/aspect cases showed significant differences ( Table 2). The greatest departures from expected proportions were seen in high occurrences of mortality on S and SE aspects, and high occurrences of succession on N and NE aspects.  . Expected proportion of aspects within the WCW (in black) compared to the observed proportion of aspect for each change class. For example, roughly 0.1 (or~10%) of the land in WCW is on a SE aspect (black bar). If there were no aspect effect, we would expect the same proportion of high mortality cases to be on SE aspects. Instead, > 0.4 (or >40%) of high mortality occurs on SE aspects (red bar).

Elevation
The distribution of elevation within change classes shows that the greatest canopy changes were seen in high elevation. High mortality and high succession classes had the highest mean elevation, while the lowest mean elevation was in the no change class (Figure 8). This result suggests that the degree of change, for both mortality and succession, covaried with elevation. To test this hypothesis, change classes were pooled by absolute value (e.g., −3 combined with 3, −2 combined with 2), and a Spearman rank test was performed. The results showed that degree of canopy change was significantly correlated with elevation (r = 0.174, p < 0.001), though the correlation was weak.  Figure 8. Elevation distributions in WCW, normalized with a probability density function. For each category, the hypsometric distribution of elevation within that category is displayed. The further right a curve extends from each baseline, the more likely it is to find values at that elevation within each category. Across all WCW, in black, elevation is evenly distributed, with elevation above 2200 m a.s.l. being less common. Fire, in orange, was more likely below 1600 m a.s.l., and high mortality and high succession were more likely at ATE elevations, between 1800-2100 m a.s.l. Corresponding t-tests can be seen in Table 3.  Figure 8. Elevation distributions in WCW, normalized with a probability density function. For each category, the hypsometric distribution of elevation within that category is displayed. The further right a curve extends from each baseline, the more likely it is to find values at that elevation within each category. Across all WCW, in black, elevation is evenly distributed, with elevation above 2200 m a.s.l. being less common. Fire, in orange, was more likely below 1600 m a.s.l., and high mortality and high succession were more likely at ATE elevations, between 1800-2100 m a.s.l. Corresponding t-tests can be seen in Table 3.
Between change class differences in elevation distribution were generally significant. One-way ANOVA showed that the mean elevation of at least one of the change classes was significantly different from the remaining groups (F = 743.16, df = 6, p < 0.001). The result of the pairwise t-tests showed that, of the 21 possible pairwise combinations, 20 showed significant differences (Table 3).
Mean elevation of fire-exposed areas was lower than mean elevation of non-fire-exposed areas in every change class (Figure 9). This is as expected, given the decreasing probability of fire with increasing elevation [15]. While factoring by fire showed between group differences in elevations, the significant correlation between degree of canopy change and elevation was preserved in both fire-exposed cases (r = 0.340, p < 0.001), and non-fire-exposed cases (r = 0.110, p < 0.001).
Remote Sens. 2020, 12, x FOR PEER REVIEW 14 of 22 Between change class differences in elevation distribution were generally significant. One-way ANOVA showed that the mean elevation of at least one of the change classes was significantly different from the remaining groups (F = 743.16, df = 6, p < 0.001). The result of the pairwise t-tests showed that, of the 21 possible pairwise combinations, 20 showed significant differences (Table 3).
Mean elevation of fire-exposed areas was lower than mean elevation of non-fire-exposed areas in every change class (Figure 9). This is as expected, given the decreasing probability of fire with increasing elevation [15]. While factoring by fire showed between group differences in elevations, the significant correlation between degree of canopy change and elevation was preserved in both fireexposed cases (r = 0.340, p < 0.001), and non-fire-exposed cases (r = 0.110, p < 0.001). Figure 9. Distribution of elevation within change class, factored by fire exposure. In each canopy change class, median elevation is higher in areas that were not exposed to fire.

Slope
Significant differences in slope distribution were observed among change classes. The highest mean slope values were in the high mortality and high succession classes. The lowest mean slope value was in the no change category ( Figure 10). As with elevation, the degree of change correlated positively with slope. Change classes were pooled by absolute value, and a Spearman rank test was performed. A weak correlation was found between slope and the degree of canopy change (r = 0.126, p < 0.001). Figure 9. Distribution of elevation within change class, factored by fire exposure. In each canopy change class, median elevation is higher in areas that were not exposed to fire.

Slope
Significant differences in slope distribution were observed among change classes. The highest mean slope values were in the high mortality and high succession classes. The lowest mean slope value was in the no change category ( Figure 10). As with elevation, the degree of change correlated positively with slope. Change classes were pooled by absolute value, and a Spearman rank test was performed. A weak correlation was found between slope and the degree of canopy change (r = 0.126, p < 0.001).   Table 4.  Table 4.
The differences in slope distribution between change classes were generally significant. The slope distribution of at least one of the change classes was significantly different from the remaining groups (One-way ANOVA, F = 213.69, df = 6, p < 0.001). The result of the pairwise t-tests showed that, of the 20 possible pairwise combinations, 17 had significant differences (Table 4).
In a majority of change classes (−3 through 1), mean slope was lower in fire-exposed area than in non-fire-exposed areas ( Figure 11). The hypsometry of this watershed shows a positive correlation between elevation and slope. Given that fires predominantly occurred at lower elevations, lower slope in fire-exposed areas was expected. In medium and maximum succession classes, fire-exposed areas had slightly higher mean slope than non-fire-exposed areas. Degree of canopy change and slope were significantly correlated in both fire-exposed cases (r = 0.244, p < 0.001), and non-fire-exposed cases (r = 0.027, p < 0.001).
Remote Sens. 2020, 12, x FOR PEER REVIEW 16 of 22 between elevation and slope. Given that fires predominantly occurred at lower elevations, lower slope in fire-exposed areas was expected. In medium and maximum succession classes, fire-exposed areas had slightly higher mean slope than non-fire-exposed areas. Degree of canopy change and slope were significantly correlated in both fire-exposed cases (r = 0.244, p < 0.001), and non-fire-exposed cases (r = 0.027, p < 0.001). Figure 11. Distribution of slope within change class, factored by fire exposure. Median slope was generally higher in regions that were not exposed to fire.

Upslope Advance
Of the five closed canopy transects used to determine upslope elevation advance, all occurred in areas where Engelmann spruce was the dominant species in 2006 (Table 5). In each of these transects, mean elevation was higher in 2006 than in 1914. The mean elevation of closed canopy forest in these unburned areas increased by 40.5 m over 92 years, which translates to an annual increase of 0.44 m.

Discussion
Our study demonstrated an example of topographic patterns of postfire colonization that conforms to the edaphic process model of postfire colonization, using direct quantitative observation Figure 11. Distribution of slope within change class, factored by fire exposure. Median slope was generally higher in regions that were not exposed to fire.

Upslope Advance
Of the five closed canopy transects used to determine upslope elevation advance, all occurred in areas where Engelmann spruce was the dominant species in 2006 (Table 5). In each of these transects, mean elevation was higher in 2006 than in 1914. The mean elevation of closed canopy forest in these unburned areas increased by 40.5 m over 92 years, which translates to an annual increase of 0.44 m.

Discussion
Our study demonstrated an example of topographic patterns of postfire colonization that conforms to the edaphic process model of postfire colonization, using direct quantitative observation with monoplotting on repeat photography. In fir-spruce forests, postfire regrowth in high elevation regions, like the ATE, likely has edaphic controls. As a result, colonization in xeric areas, like south-facing slopes, can lag several decades or even centuries behind regrowth in mesic areas, like north-facing slopes. These edaphic controls on colonization were generally suggested based on research using stand level age-reconstructions, but given the long gap in colonization between mesic and xeric colonization, few direct observations of this effect have been made in imagery. By increasing the temporal range of observation using monoplotting on historic repeat photography, the topographic differences in postfire colonization are now visible in imagery.

Aspect
Greater than expected proportions of succession occurred on N and NW aspects, while higher than expected proportion of mortality occurred on S and SE aspects. Correlation between aspect and canopy cover change was associated with fire-exposure, suggesting soil moisture stress levels may have been a significant barrier to postfire recruitment.
It is known that seedling recruitment of subalpine fir and Engelmann spruce is inhibited on warm-aspects in the ATE, primarily as a function of soil moisture stress caused by high temperature and solar radiation [55]. This process has been demonstrated experimentally: heating experiments in trees at ATE elevations show that moisture limitation both inhibits conifer seedling recruitment [56] and reduces carbon gain [57]. In a postfire environment, insolation would increase due to reduced interception and increased gap fraction in the burned canopy. This increased insolation could have caused aspect specific mortality on warm aspects in WCW, particularly the SW aspect. Postfire erosion may also contribute to mortality. A study in a watershed adjacent to WCW [58] shows that total suspended sediment levels are eight times higher in fire-exposed areas, an effect that is exaggerated by steep terrain at high elevations. High levels of erosion in fire-exposed, subalpine areas may suppress the moisture retention capacity of the soil and exacerbate other obstacles to recruitment.
On the north-and east-facing slopes in the WCW, environmental factors may have mitigated soil moisture stress. Reduced energy inputs and wind exposure diminish the moisture losses to snowpack sublimation, soil and vegetation evapotranspiration. For example, it is known that the end of winter mean snowpack depth on north-facing slopes can be~40% greater than south-facing slopes in the southern Canadian Rockies [59]. North-facing snow packs also persist longer into the growing season, thus ensuring melt water is available to support processes of germination and plant maintenance. Aspect specific patterns of snow accumulation can affect the functional length of the growing season, in turn altering ATE position [9]. Additionally, aspect can predict postfire recruitment of subalpine fir [17], with cool, north-facing aspects having rates of establishment twice as high as on south-facing aspects. Interactions among solar radiation (aspect), fire history, and elevation (in addition to anthropogenic activity) predict broad patterns of succession in a large region of the Rocky Mountains of southern Alberta [47]; this includes the encroachment on closed canopy forest on alpine meadow, or ATE advance. Together, these findings suggest postfire colonization is highly dependent on soil moisture, and that this limitation is exacerbated at high elevations.
Aspect effects unrelated to fire disturbance are reported in a number of studies of ATE dynamics, globally [9,[60][61][62]. However, the proposed mechanisms of these aspect effects are often site-specific, and vary by study region. These include differences in lapse rates among aspects [61] and diurnal cycles of cloud cover causing preferential aspects in equatorial rainforests [60].
The suggestion that soil moisture limitation, caused by warm and dry conditions, prevented postfire recruitment is consistent with known climatic history of the region. The early 20th-century is noted as being a historically warm period in the eastern slopes of the Canadian Rockies, with dendroclimatic records showing the periods from 1917-1936 and 1939-1958 as having the highest temperature anomalies since the 1400s [63]. In fact, one of these reconstructions listed 1936, the year of the Pass Creek fire in WCW, as the second warmest year in the~1000 year record [63]. This observation suggests the possibility that abnormally high atmospheric temperatures may have a two-fold effect on subalpine forests, first causing the desiccated conditions required for a burn at high elevation, then limiting expected regrowth.

Elevation
High degree change was almost entirely restricted to elevations greater than 1800 m, which coincide with the lower extent of ATE at the study site. This finding is consistent with the expectation that, as a transitional boundary that occurs as a function of elevation, canopy cover at ATE would be subject to a high level of variability [13,64,65]. Fire-exposed areas had lower mean elevations in all change classes, but the majority of canopy cover change was seen in these areas. Elevations greater than 1800 m accounted for a proportionally small region of fire-exposed areas, and yet the highest degree of change was seen at those elevations. The positive correlation between change degree and elevation is maintained in fire exposed areas suggesting that postfire mortality and re-establishment are potentially affected by factors that covary with elevation.
In the Canadian Rockies, there is evidence of both upslope advance and downslope recession of ATE elevation over the past 300 years. Our estimate of ATE advance in areas undisturbed by fire is consistent with other observations in the region; reported a subalpine fir advance rate of 0.28-0.62 m yr −1 along an upslope transect [36].

Slope
Treelines in high slope areas are known to be subject to mechanical disturbances by avalanches [11,66], or enhanced soil erosion following disturbance [17,67], which introduces the prospect that treelines in WCW may be orographically regulated [6]. Cyclical patterns of disturbance and regrowth may explain the high degree of both mortality and succession seen in high slope areas. However, in a nearby site in Glacier National Park (GNP), avalanche-regulated treelines were stable over a >70 year repeat photograph record, showing no cycle of disturbance and advance [40]. If the observed slope pattern in WCW is a result of cyclical processes, then this discrepancy between WCW and GNP warrants further investigation.
The slope of fire-exposed areas was generally lower than non-fire-exposed areas. High terrain slopes can accelerate the spread of fire [68], but this effect was outweighed by the extent of the WCW fire that occurred on the low slopes of the valley bottom. In two cases, medium and high succession classes, slopes were greater in fire-exposed areas than in non-fire-exposed areas, which agrees with research showing that slope was one of the factors that predicted re-establishment of subalpine fir following disturbance by fire [17]. In that study, 75% of establishment occurs on slopes with angles between 40 • -60 • . The authors reason that this intermediate slope range is steep enough to avoid competition by shrubs during re-establishment, but shallow enough to avoid frequent disturbance by rockslide and avalanche. That explanation is inconsistent with the present finding, as areas with slopes above 30 • saw a high proportion of both mortality and succession.

Conclusions
The WSL-MT monoplotting technique provided the ability to perform quantitative spatial analysis of mountain landcover change on oblique, historic repeat photography from the MLP. Canopy cover estimates from these images agree with canopy cover measurements using airborne lidar [49]. The method provides the ability to observe ATE at a spatial resolution and temporal extent which is sufficient to identify topographic correlation to change, which is an important step in understanding interactions between modulating factors. The MLP dataset is ideal for multi-scale investigation of change in the ATE; the spatial extent of MLP photographs ranges from individual meadows (~100 m), to valleys (~10 km), at sites which can be several hundred kilometers apart [47].
In ATE in the WCW, patterns of postfire colonization related to edaphic controls have hindered regrowth and locally suppressed ATE on south-facing aspects, while considerable postfire colonization was seen on north-facing aspects. This observation appears to directly support hypothesis of edaphic controls on colonization postfire. While fire is not a direct control on ATE position, it is in this case a modulating factor, and acts to locally suppress ATE below the thermal limits. While much has been learned about these edaphic controls using stand reconstructions, these can be limited in their understanding of microsite-processes and larger regional trends. Repeat photography provides the opportunity to observe the patterns at a high spatial resolution, over broad regions and diverse landscapes. Future research in the region would benefit from investigating the interactions between fire, moisture limitation, and other modulating factors, using the multi-scale approach to understanding ATE dynamics [4,5].
In unburned areas on north-facing slopes, where moisture is not limiting, the ATE has systematically advanced. These observations of successful regeneration and ATE advance are consistent with increasing alpine isotherms and energy availability, as is expected under many climatic change scenarios, provided sufficient moisture or nutrient availability. More research is needed to further refine these hypotheses, particularly with respect to snowpack quantities and persistence at and above ATE.