**Contextualizing the 2019–2020 Kangaroo Island Bushfires: Quantifying Landscape-Level Influences on Past Severity and Recovery with Landsat and Google Earth Engine**

#### **Mitchell T. Bonney 1,\*, Yuhong He <sup>1</sup> and Soe W. Myint <sup>2</sup>**


Received: 27 October 2020; Accepted: 30 November 2020; Published: 2 December 2020

**Abstract:** The 2019–2020 Kangaroo Island bushfires in South Australia burned almost half of the island. To understand how to avoid future severe 'mega-fires' and how vegetation may recover from 2019–2020, we can utilize information from the bulk of historical fires in an area. Landsat time-series of vegetation change provide this opportunity, but there has been little analysis of large numbers of fires to build a landscape-level understanding and quantify drivers in an Australian context. In this study, we built a yearly cloud-free surface reflectance normalized burn ratio (NBR) time-series (1988–2020) using all available summer Landsat images over Kangaroo Island. Data were collected in Google Earth Engine and fitted with LandTrendr. Burn severity and post-fire recovery were quantified for 47 fires, with a new recovery metric facilitating comparison where fire frequency is high. Variables representing the current burn, fire history, vegetation structure, and topography were related to severity and yearly recovery with random forest and bivariate analysis. Results show that the 2019–2020 bushfires were the most widespread and severe, followed by 2007–2008. Vegetation recovers quickly, with NBR stabilizing ten years post-fire on average. Severity is most influenced by fire frequency, vegetation capacity and land use with more severe burns in nature conservation areas with dense vegetation and a history of frequent fires. Influence on recovery varied with time since fire, with initial (year 1–3) faster recovery observed in areas with less surviving vegetation. Later (year 6–10) recovery was most influenced by a variable representing burn year and further investigation indicates that precipitation increases in later post-fire years likely facilitated faster recovery. The relative abundance of eucalypt woodlands also has a positive influence on recovery in middle and later years. These results provide valuable information to land managers on Kangaroo Island and in similar environments, who should consider adjusting practices to limit future mega-fire risk and potential ecosystem shifts if severe fires become more frequent with climate change.

**Keywords:** Landsat; time-series; Google Earth Engine; LandTrendr; NBR; random forest; fire history

#### **1. Introduction**

The 2019–2020 Australia bushfires were among the largest and most devastating on record. It is estimated that close to 200,000 km2 of land burned across the country, with projected losses of nearly three billion animals and potentially dozens of rare species [1,2]. There was also a major impact on society, with thousands of houses destroyed, 33 people killed and billions in economic losses. Kangaroo Island, off the coast of South Australia, was heavily impacted with about half of its land area burned. The 2020 fires have especially impacted the island's animal populations, with thousands of kangaroos,

koalas, sheep, and other species killed and endangered species, including the dunnart and glossy black cockatoo, being pushed close to extinction [3,4].

Bushfires are becoming increasingly common in Australia as the climate warms [5]. Temperatures in Australia since 2010 are over 1 ◦C warmer than as recently as the mid-1970s [6] and combined with a period of drought creates a high fire risk scenario. The 2019–2020 fires occurred in the context of Australia's hottest and driest year on record, and a drought that had been ongoing for two years [2]. Despite increased mega-fire occurrence, including the 2009 Black Saturday fires that led to 173 deaths in Victoria, there is still policy debate about the effectiveness and implementation of forest management techniques such as fuel-reduction burning [7]. The occurrence of another devasting fire season a decade later has undoubtedly revived the debate over best practices in terms of land management policies and reducing fire risk. As fires are expected to become larger and more common [2], a better understanding of fire and recovery processes across the landscape is needed. If land managers knew landscape-level drivers of burn severity and recovery from past fires, they could then make more informed decisions about how and where to allocate limited resources across their landscape to mitigate future burns. Quantification of severity, recovery, and their drivers are also important at broad scales in terms of grasping their impact on forest carbon budgets [8].

A compelling way to obtain a landscape-level understanding of fire processes is to build a comprehensive dataset of fire history (e.g., when, where, severity, recovery) for an area. Remote sensing provides this opportunity, with the freely available and 30 m spatial resolution Landsat platform being used to observe burn severity and recovery since at least the 1990s [9,10]. These early studies primarily focused on observing single or small groups of fires with at least two images, but not enough to achieve a long-term yearly resolution of post-fire recovery. However, since the opening of the Landsat archive in 2008, the remote sensing community can support more dense Landsat time-series analysis [11]. More recent studies have made use of this to build longer time-series at near-annual temporal resolutions, but still often focus on a single fire [12]. Many of these types of studies use dozens of images downloaded from USGS to first calculate the difference between pre- and post-fire image spectral values with a given vegetation index and then quantify trends from subsequent post-fire images.

In the last few years, the development of segmentation algorithms such as LandTrendr that fit noisy yearly data to simplified trajectories and their release on Google Earth Engine (GEE), a cloud-based big data processing platform, has enabled more advanced analysis [13–15]. GEE allows for time-series to be built with hundreds or thousands of images since they do not have to be downloaded and can instead be accessed in the cloud. Studies have begun to make use of these platforms and algorithms to examine fires by building spectrally consistent surface reflectance time-series and quantifying severity and recovery as large negative trajectories followed by gradual positive trajectories [16,17]. However, the use of Landsat time-series to analyze the bulk of recent fires across a region and build a general understanding of fire processes related to these changes has focused on boreal and alpine forests [18–21].

Despite the long history of fire-focused remote sensing, there are still disagreements between fire studies that warrant more exploration. The normalized difference vegetation index (NDVI) is commonly used to quantify burn severity and recovery [12,19,22]. However, it is the normalized burn ratio (NBR) that has proven more effective in many different forested environments [16,20,23,24] because the shortwave infrared (SWIR) band is more sensitive to forest structure and moisture [25]. In Australian eucalypt forests, Hislop et al. [26] concluded that NBR is more reliable than NDVI and six other indices for estimating severity and recovery after using Glass's delta to compare pre- and post-fire data. Studies have reported spectral recovery periods differing from as low as three years to multiple decades with variations linked to location (e.g., vegetation type, terrain, climate etc.), the fire itself (e.g., severity) and the vegetation index used [12,16,23,27,28]. Hislop et al. [26] found that NDVI recovered to pre-fire levels after 3–5 years, while NBR took 8–10 years to recover.

There has been a wide range of studies concerned with identifying drivers behind differences in post-fire recovery. Burn severity, often quantified by dNBR or ΔNDVI, is commonly used as a predictor of recovery rate, with some studies indicating that more severely burned areas recover faster [16,23,29] and other studies linking more severe burns with slower recovery [24,27,30] or showing no significant relationship [12]. Other variables, such as topography [12,22,31], climate [16,19,24] and vegetation type [16,28] are used, with differing relationships to recovery. In south-east Australia, vegetation type analysis has shown differing spectral response to fire for forested and shrubland classes of various vegetation densities [32]. There has been less study on the influence of the historical fire regime, such as time since previous fire or longer-term frequency, on recovery [31,33]. On Kangaroo Island, where fires are common, fire frequency is an important variable controlling the distribution and survival of many fire-adapted vegetation groups, with recommended intervals from as low as six years to multiple decades [34,35]. Within the long-term Landsat record, there is thus an opportunity to monitor recovery in areas impacted by regular fires, with studies often focusing on single burn pixels and removing those burned multiple times from analysis [28]. Furthermore, there has been little study of how variables influencing recovery themselves change through time, with studies conducting comparisons for a single selected post-fire year [16]. By comparing at a yearly interval, we can build a full temporal understanding of influences on recovery from initial to late recovery stages. Although there has been extensive work on identifying drivers of recovery, less consideration has been given to influences on severity with fuel amount and time since last fire considered [36–39].

In this study, we use Landsat time-series, fire history spatial data, and GEE to understand landscape-level influences on burn severity and post-fire recovery in eucalypt-dominated forests across a region influenced by frequent fires. This study is based on remote sensing: vegetation, severity, and recovery through time are inferred from Landsat-based proxies of these variables and the LandTrendr algorithm in GEE [13–15]. We build a methodology that is useful for areas that undergo multiple burns within the Landsat record and use these methods to best understand the fire dynamics of the region. The primary objectives are to (1) monitor the fire history (i.e., extent, severity, recovery) of Kangaroo Island forested areas, especially in the context of the 2019–2020 bushfires; and (2) quantify the influence of landscape-level variability (e.g., in current burn context, fire history, vegetation structure and topography) on severity and yearly recovery. Identifying how influences on recovery change through time after fires is an important aspect. We use these objectives to gain information on why the 2019–2020 bushfires were so severe and how vegetation may recover from these bushfires in the coming years based on the historical data.

#### **2. Materials and Methods**

#### *2.1. Study Area*

Our focus is on Kangaroo Island (Figure 1). This island, located in South Australia southwest of Adelaide, consists of a mix of mostly eucalypt forest (nature conservation, residual cover, plantations) and agricultural (grazing pastures, cropping) landscapes [40,41]. Dominant tree species vary across the island, with *Eucalyptus diversifolia* being common along the coast, *E. remota* common in the dense western forests of Flinders Chase National Park (NP), *E. cladocalyx* common in riparian landscapes and *E. baxteri* fragmented across Flinders Chase NP and in central agricultural settings [35,42]. Species grow as either trees with a single stem or as mallees which have multiple stems and may take on a shrubby growth form [42]. Woodlands grow denser and are more widespread in the conserved western and extreme south of the island. However, residual native cover remains along riparian corridors, ridges and road networks in the more developed central and eastern parts of the island in notable quantities [41]. Growing development in the last half century has put some of these remnant and forest edge landscapes under threat [40,41]. The warmest months are December to March, and the wettest months are June to August [6]. In the summer, temperatures reach an average high of about 24 ◦C (19 ◦C annual average) and have been steadily warming in recent decades (Figure 2A). Annual

precipitation fluctuates, with an average of 567 mm since 1988, and is consistently highest in the forested west and central parts (Figure 2B).

**Figure 1.** Kangaroo Island, Australia. Land use classes based on the Australian Land Use and Management Classification V8 [43]. Population centres and other points of interest shown. Flinders Chase NP (and the associated Ravine Des Casoars Wilderness Protection Area) shown with a dashed black outline [44].

**Figure 2.** Mean maximum temperatures (**A**) and total annual rainfall (**B**) across Kangaroo Island from 1988-2020. Red triangles indicate the two most widespread fire seasons on record (2008, 2020). Climate data collected from on-island weather stations [6] and merged into regional clusters. West: stations within and near Flinders Chase NP (4 stations); Central: nearby Parndana and Vivonne Bay (9); Northeast: nearby Kingscote (5); Southeast: nearby American River and Cape Willoughby (5). Yearly is the average of these clusters.

Kangaroo Island's climate and landscape make it vulnerable to regular bushfires, especially during drier periods (Figure 2B) and in forested areas. Spatial records of fires have been kept since the early 1930s, although data collection is sparse until the mid-1950s [45]. The two most widespread fire seasons occurred in summer 2007–2008 and 2019–2020 (below as 2008 and 2020) where 950 km2 and 2100 km<sup>2</sup> burned, respectively. Both fires occurred after multiple years of below average rainfall. The 2008 fires were initially ignited by lightning strikes and burned for two weeks, although it was

mostly contained within Flinders Chase NP [46]. The 2020 fire season was devastating for Australia overall, but was especially widespread on the island. The fires were not contained to the park and burned throughout the agricultural countryside, with population centres across the island being evacuated [3]. Although the island's eucalypt vegetation is adapted to fire, the severity and close temporal proximity of these mega-fires over the same areas may impact their ability to recover. In the years before these fires, local analysis indicated that much of the area of all vegetation groups besides *E. remota* have not burned in over 40 years and recommend the introduction of prescribed burning [35].

#### *2.2. Building Landsat Fire History Time-Series*

Using GEE, a yearly time-series (1988–2020) of cloud-free surface reflectance medoid data from Landsat Thematic Mapper (TM), Enhanced Thematic Mapper Plus (ETM+), and Operational Land Imager (OLI) was constructed over Kangaroo Island. ETM+ imagery includes SLC-on and SLC-off data. Clouds and cloud shadows were automatically removed from input images with the CFMASK algorithm [47]. Georeferenced images were converted to surface reflectance using the LEDAPS algorithm for TM and ETM+ and LaSRC for OLI [48,49]. To minimize between-sensor discrepancies, OLI bands were converted to ETM+ spectral equivalents using the regression equations reported in Roy et al. [50]. When multiple cloud-free values were available for a given pixel in a given year, medoid compositing reduced them by finding the value that is closest to the center (i.e., median) for both the X (i.e., time of year) and Y (i.e., reflectance) axes. [51]. These processing steps were completed in GEE using the buildSRcollection code within LandTrendr to assure a consistent time-series and limit atmospheric and sensor influences to the highest extent possible [14].

For each year, satellite images were eligible for the time-series if collected between December 25 and March 25. This approximately corresponds to the Australian summer used in other Landsat time-series studies in the region [26,52] while also maximizing the amount of cloud-free data. However, three years in the study period (1988–2020) are heavily impacted by missing data (i.e., 1990, 1994, 1997), with 1997 having no available images (Figure 3). Some other years (e.g., 1989, 1992, 1998, 2011, 2012) also have small areas of missing data. Data availability for each year was checked in GEE using the buildClearPixelCountCollection code within LandTrendr [14].

To infer vegetation abundance in connection to moisture, the surface reflectance time-series was converted to the normalized burn ratio (*NBR*), which is calculated as follows:

$$NBR = \frac{NIR - SWIR}{NIR + SWIR} \tag{1}$$

where *NIR* represents the near-infrared Landsat band, and *SWIR* represents a shortwave infrared band [53]. In Landsat TM and ETM+, *NIR* is band 4 and *SWIR* is band 7. In Landsat OLI, *NIR* is band 5, and *SWIR* is band 7. We decided to use *NBR* because it has best quantified fire damage and recovery in other eucalypt forests [26]. A related index called the differenced *NBR (dNBR)* is used to quantify burn severity in remote sensing. It is defined as the difference between NBR in the year before a fire (*preNBR*) and the NBR in the year after the same fire (*postNBR*) based on available imagery (Figure 4A).

$$dNBR = preNBR - postNBR \tag{2}$$

If a *dNBR* year is missing imagery, LandTrendr (discussed below) will fit the decline over two years when filling in the missing value (Figure A1). In this case, *preNBR* is defined as the *NBR* value before the decline and *postNBR* is the lowest *NBR* value at the end of the decline. In this study, all *NBR*-based indices are multiplied by 1000. Based on visual interpretation of burns, the focus on forested land cover and the need to minimize more gradual declines being classified as burns, we set 100 *dNBR* as the 'burned' threshold.

**Figure 3.** Number of cloud free images available for every pixel in each summer (December 25–March 25) across Kangaroo Island.

**Figure 4.** Conceptual model outlining LandTrendr and the quantification of burn severity (*dNBR*) and recovery *(%RtoMV*). (**A**) LandTrendr fitted and medoid *NBR* values for an example pixel that experienced fires in the 1992 and 2001 seasons. (**B**–**D**) Comparison of methods to quantify post-fire recovery for a multi-burn pixel. (**B**) Fitted *NBR* slope. (**C**) *NBR* recovery (%) to *preNBR* [16]. (**D**) *NBR* recovery (%) to the maximum amount of vegetation for that pixel in the time-series (*mNBR*).

The *NBR* time-series was input into the GEE implementation of LandTrendr, which is an algorithm that creates a time-series of spectral-temporal segments similar to a piecewise linear function [14]. LandTrendr functions on the pixel-level by first removing ephemeral spikes due to non-real change (e.g., un-masked clouds) and fitting a simple linear regression to the time-series [15]. Real disturbances are differentiated from non-real change based on how quickly the spectral index recovers to the pre-change values, with non-real changes returning to pre-change values after a single or very small number of years (the threshold of which is defined by recoveryThreshold). Breakpoint years are then identified based on their deviation from the regression line, with the number of final breakpoint years (and therefore the number of segments) defined by the maxSegments and vertexCountOvershoot parameters. Straight line segments are used to connect breakpoint years, with successively simpler models created from the maxSegments number down to one segment (i.e., simple linear regression). The simplest model with the best fit is chosen as the final fitted time-series for the pixel in question (Figure 4A).

LandTrendr provides a default set of parameters for general use [14] (Table 1). However, for more accurate quantification of spectral-temporal trends, it is recommended to find a set of parameters that best fits the process being observed over a study area [54]. To do this, we selected 100 random pixels across Kangaroo Island and built their fitted time-series using different variations of LandTrendr parameters. The LandTrendr developers provide a web-app that allows for efficient testing of different parameter sets (https://emaprlab.users.earthengine.app/view/lt-gee-pixel-time-series). For each of these pixels, we recorded any year that experienced a *dNBR* value considered a burn (i.e., 100+) in the fitted time-series. We then cross-referenced this information with a spatial dataset of fire polygons created by the Government of South Australia [45], who record historical burn scars using digitization from aerial photography, *NBR* calculation from satellite imagery and other sources and methods. A successful match occurred when the pixel in question experienced a recorded fire, and *dNBR* was at least 100 in the same year. Note that not all *dNBR* values above 100 are due to fires; plantation harvest, and other human and natural disturbances may also lead to high *dNBR* values.

**Table 1.** LandTrendr parameters. Default: generalized LandTrendr parameters. Selected: LandTrendr parameters that best characterize fire history on Kangaroo Island. Fire Identification: percent of sample pixels in a recorded fire polygon that exhibits an *NBR* loss of at least 100 during the correct year.


For our random sample of pixels, we determined which parameter set most accurately identified recorded fires (Table 1). *NBR* more consistently identified fires compared to other indices, such as the NDVI and Tasseled Cap Angle. The default parameters correctly identified 70% of recorded fires overall and 80% in forested land cover, which is where fires are expected to be more easily observed by Landsat and is the focus of this study. However, the most accurate LandTrendr parameters correctly identified 81% of recorded fires overall and 93% of forest fires. In this parameter set, maxSegments was increased from 6 to 8 to best fit the time-series of some multi-burn pixels. The recoveryThreshold was also increased from 0.25 to 0.75 because post-fire spectral vegetation indices recover quickly in this area, sometimes incorrectly activating this filter (which is meant to remove erroneous segments) at the default value. This parameter set was selected for use in subsequent analysis. A full description of all LandTrendr parameters is available in Kennedy et al. [14,15].

With the selected parameters, the yearly surface reflectance time-series was converted to a LandTrendr fitted time-series across Kangaroo Island in GEE. Pixel-level information was extracted from the fitted time-series. However, since our focus is on monitoring fire history, we used the spatial dataset of recorded fires for subsequent analysis [45]. From this dataset, we defined our unit of analysis as a fire area extent over one km2 that burned between 1988 and 2020. There are 47 fires that fit these criteria (Table 2).

Another key aspect of this study focuses on monitoring post-fire recovery. LandTrendr is well suited to quantify this using the fitted *NBR* slope of the post-fire segment (Figure 4B). A more useful method, which quantifies percent *NBR* recovery for years post fire in terms of how much a burned pixel has recovered to pre-fire vegetation amounts (*%RtoPfV*), allows for comparison of different fires [16] (Figure 4C). It is calculated for a given year *X* as:

$$\% \text{RtoP} \text{f}V = \frac{(\text{Year} \, \text{X} \, \text{NBR} - \text{postNBR})}{(\text{preNBR} - \text{postNBR})} \ast 100 \tag{3}$$

This method is robust in many scenarios but may be misleading in areas that experience additional fires before the area in question has fully recovered from the previous fire. For example, in Figure 4C, *%RtoPfV* quantifies recovery for two fires that burn the same pixel nine years apart. Since Fire 2 occurred when the pixel was still recovering from Fire 1, it more quickly reaches *preNBR* and thus indicates a faster recovery rate even though their *NBR* slopes are very similar.


**Table 2.** Recorded fires on Kangaroo Island [45] from 1988 to 2020 with extent polygons larger than 1 km2. Vegetation species determined by the most common National Vegetation Information System (NVIS) 5.1 [42] Level 5–6 vegetation species overlapping burned pixels for each fire. Fires listed in chronological order.

Therefore, we modified *%RtoPfV* to better account for multi-burn situations (Figure 4D). Rather than recovery to pre-fire vegetation, we instead quantify the percent recovery to the inferred maximum vegetation capacity for a given pixel (*%RtoMV*):

$$\% \text{Rto} MV = \frac{(Y \text{arr } X \text{ NBR} - \text{postNBR})}{(m \text{NBR} - \text{postNBR})} \text{\*100} \tag{4}$$

where *mNBR* is defined as the highest fitted *NBR* value in the time-series. With this approach, multiple recoveries in the same pixel are gauged against the same value, which allows for a more direct comparison of post-fire recovery. In Figure 4D, *%RtoMV* accurately quantifies that the recovery rate for both fires is very similar despite Fire 2 occurring amid recovery from Fire 1. Kangaroo Island experiences regular fires, with parts of Flinders Chase NP having burned at least nine times since 1930 (Figure 5), meaning that *%RtoMV* is more robust in this setting.

**Figure 5.** Number of times burned since 1930. Based on the size of fire history polygons [45] to be consistent across the entire record. Polygons from 1950 and onwards were used in the creation of fire history variables (BLX, SLB).

Burn severity (i.e., *dNBR*) and post-fire recovery (i.e., *%RtoMV*) were calculated at the pixel-level and averaged for each fire since that is the unit of analysis for this landscape-level study. Eligible fire polygons were used to extract *NBR* values from the fitted time-series. Only pixels with a *dNBR* value of at least 100 within these polygons were used in subsequent analysis. For *%RtoMV* calculation, the *postNBR* year is considered year 0, and %RtoMV is quantified for the following years. The recovery trajectory for each pixel ends when either the time-series ends (i.e., 2020), when another recorded fire occurs (i.e., spatial-temporal overlap with a fire polygon), or when some other notable damage occurs to vegetation (i.e., fitted *dNBR* of at least 50). When less than 25% of pixels are left recovering from a fire, the average recovery trajectory also ends.

#### *2.3. Assessing Landscape-Level Influences*

To quantify the influence of landscape-level variability on burn severity and recovery, 20 variables were built for comparison (Table 3). These variables mostly fall into four categories: the current burn (i.e., burn extent, burn year, time of year, *postNBR*, *dNBR*), fire history (i.e., years since last burn, longer-term fire frequency), vegetation structure (i.e., *preNBR*, *mNBR*, land use, vegetation group, species) and topography (i.e., elevation, slope, aspect, topographic position, topographic wetness). One additional variable assesses the influence of missing data in the Landsat time-series (i.e., part of burn area impacted by missing data in a *dNBR* year). These variables were built from the Landsat time-series data, the fire polygon data [45], Australian land use data [43], NVIS [42] and SRTM 1 Arc-Second Global DEM data [55]. Fire history variables are based on the spatial-temporal overlap

of fire polygons (e.g., Figure 5). Final land use classes and vegetation groups were chosen based on how common they were across all fires. Information was extracted at the pixel and landscape level. If acquired at the pixel-level, mean values were calculated for each fire.

**Table 3.** Descriptions of variables tested for influence on burn severity (*dNBR*) and recovery (*%RtoMV* year 1–10). Pixel-level variables are based on burned pixels (*dNBR* ≥ 100) only and averaged for each fire to match landscape-level analysis.


Random forest (RF) modeling was used to quantify the influence of these landscape-level variables on burn severity and recovery. RF models create many high variance decision trees that are built using bootstrap sampling (i.e., with replacement) from the original dataset [58]. Each RF model included 500 decision trees (ntree). Due to the relatively small sample size of the analysis, we ran each RF model iteration 30 times with the goal of finding the smallest set of variables that explained the highest average variation in *dNBR* or *%RtoMV*. Separate average RF runs were conducted for each year of *%RtoMV* from years 1 to 10. Only variables that improved overall model performance and were not highly correlated to each other were included in the final model runs. Modeling was conducted using the "randomForest" package in R [59].

Similar to Bright et al. [16], we evaluated overall model performance with percent root mean square error (*%RMSE*) and percent variance explained (*%VE*):

$$\%RMSE = \frac{\sqrt{MSE}}{\overline{X}} \ast 100\tag{5}$$

$$\%VE = \frac{1 - MSE}{var(X)} \ast 100\tag{6}$$

where *X* is either *dNBR* or *%RtoMV*, and *MSE* is the mean square error. *%VE* explained is the pseudo R<sup>2</sup> for RF models. We also evaluated the individual influence of each explanatory variable with the model improvement ratio (MIR). MIR is a standardized measure of variable importance (i.e., influence on *dNBR* or *%RtoMV*) between 0 and 1, where 0 is least important, and 1 is most important [60].

Further comparisons between those variables and both *dNBR* and *%RtoMV* were completed at the fire level. Scatterplots were created to relate these variables to *dNBR*, with strength (R2) and statistical significance (*p* < 0.05) tested. Accompanying boxplots were produced by splitting average *dNBR* values into three groups based on natural breaks in each explanatory variable and testing for significant differences with the nonparametric Mann-Whitney (MW) and Kruskal-Wallis (KW) tests. Time-series plots relating the variables and *%RtoMV* from years 1 to 10 were also generated. Average *%RtoMV* trajectories were also split into three groups, and significant differences were tested for each year with the same tests.

#### **3. Results**

#### *3.1. Fire History*

Kangaroo Island experiences regular fires. The four most widespread fire seasons by area burned were 1954 (17.4% of the island burned), 1959 (16.4%), 2008 (19.5%) and 2020 (45.9%) (Figure 6A).

**Figure 6.** Percent of Kangaroo Island burned in each fire season from 1929–30 to 2019–20 (**A**), based on the size of fire history polygons [45] to be consistent across the entire record. Polygons from 1950 and onwards were used in the creation of fire history variables (BLX, SLB). Prescribed burns are planned fires on government land and are minor compared to bushfires. Inset graph (**B**) with all recorded burns since 2000 shows that prescribed burns never burned more than 0.2% of Kangaroo Island in a single season.

These fires have historically been most common in the densely forested west, south, and extreme east parts (Figure 7), with some parts of Flinders Chase NP burning every ten years on average (Figure 5). Landsat 30 m data have been collected over the island since 1988, where 47 fires over one km2 have been recorded (Table 2). These fires have varied in location, size (BExt: 0.75 <sup>±</sup> 0.90), and burn severity (*dNBR* 484 ± 173) over time (Figures 7 and A2). Fires in densely forested areas, especially Flinders Chase NP, are generally larger and burn more severely. Fires in residual forest cover, fragmented across agricultural areas, are usually smaller and burn less severely. The only fire season where large swaths of the agricultural landscape burned was 2020. Broadly, the 1990s fires burned less severely than more recent fires. The two largest fires, occurring in 2008 and 2020, both followed a roughly decade-long period of relatively little burning. This is especially true for the 2020 fires, where only 65 km2 (1.5% of the island) burned between 2009 and 2019. These two fires were also the most severe, with 2008 having a *dNBR* of 797 ± 210 and 2020 having a *dNBR* of 820 ± 294.

**Figure 7.** Burn severity (dNBR) for all recorded fires (>1 km2) during the Landsat record (1988–1989 to 2019–20). Split into six periods to reduce overlap in multi-burn areas and highlight differences in burn extent through time. See Figure A2 for violin plot distributions of burn severity for each fire. See Figure A3 for a closer view of two smaller fires.

Fires have burned across Kangaroo Island in many different landscapes. However, when isolating burned pixels based on our *dNBR* analysis, it is clear that the vast majority of pixels burned in forested and shrubby areas (e.g., Figure A3). Based on land use (Figure 8A), most fires burned either primarily on nature conservation land (i.e., Flinders Chase NP, other large forested areas) or on residual native cover (i.e., fragmented forest patches located near agricultural land). One fire burned primarily on a plantation forest. Giving all fires equal weight, only 5% of burning occurred on land use classified as grazing modified pastures or cropping, and closer inspection of some of these areas indicates possible land use change (e.g., pasture to residual cover) between data date and fire date. Based on vegetation group (Figure 8B), most fires burned either primarily in denser eucalypt woodlands or mallee woodlands, which are also mostly eucalypts in a different growth form. Some fires burned primarily in other woodlands and shrublands, including *Casuarina* and *Acacia* dominated areas, but most of this class includes residual native cover that was not classified by NVIS [42] and plantation forests. Again, the influence of non-forest and shrubland is minimal, with cleared land covering similar areas as

pastures and cropping land use. Eucalypts dominate the species composition of burned areas, with 42 of 47 fires consisting primarily of *Eucalyptus* spp. (Table 2).

**Figure 8.** Land use and vegetation group composition of burned pixels (dNBR ≥ 100) for each fire. (**A**) Land use classes derived from the Australian Land Use and Management Classification V8 [43]. (**B**) Vegetation groups based on NVIS Major Vegetation Groups [42]. Visual inspection against ArcGIS Basemap imagery was used to confirm accurate land use and vegetation group classification. Fires listed chronologically as in Table 2.

The average recovery trajectories of the 40 eligible fires (i.e., not including recent fires with less than two years of recovery data), as quantified by *%RtoMV*, indicate that burned vegetation usually recovers quickly on Kangaroo Island (Figure 9).

In the first year, *NBR* recovers 22.3 ± 8.2% to its inferred maximum vegetation capacity on average. By year three, *NBR* recovers over half (51.3 ± 12.5%) of its capacity. Average *%RtoMV* follows a logarithmic pattern, reaching 80.3 ± 10.6% in year nine and stabilizing around 85% after ten years. Most fires follow a similar recovery pattern. The fastest recovering fire reached 97 *%RtoMV* in year six before the area was burned again in 2020. This fire (P5\_3) occurred in a sloped riparian mallee *E. remota* woodland adjacent to pastures. The slowest recovering fire (P2\_1) reached 48 *%RtoMV* in year 10, before the area was burned again in 2006. This fire burned in a mallee *E. diversifolia* woodland and shrubland area that is transitioning from a pasture land use. It also had the lowest burn severity (*dNBR*: 190 ± 66) of any recorded fire.

#### *3.2. Severity: Landscape-Level Influences*

RF modeling identifies seven landscape-level variables that together most influence burn severity (*dNBR*) (Figure 10). From most to least important, these variables are BL50, *mNBR*, percent nature conservation land use (%NC), BExt, BYr, TPI65, and MD (Table 3). Longer-term fire frequency, inferred by BL50, and vegetation capacity, inferred by *mNBR*, are practically identical in terms of their strong influence on *dNBR*. BL50 and *mNBR* have an average MIR of 0.964 ± 0.04 and 0.960 ± 0.06 respectively. Both BL50 and *mNBR* outcompeted similar variables, SLB and *preNBR*, that were also correlated to *dNBR* but not as strongly. %NC also has a strong influence on *dNBR*, with an MIR of 0.91 ± 0.07. A RF model of just BL50, *mNBR* and %NC can explain over 50% of the variation in *dNBR*. After these variables, there is a clear drop in influence, with BExt having an MIR of 0.44 ± 0.05, BYr 0.26 ± 0.04,

TPI65 0.20 ± 0.04 and MD 0.11 ± 0.02. However, each of these variables slightly increases model performance compared to them not being included. Overall, well over half (58.7 ± 1.3%) of the variance in *dNBR* is explained by these seven variables, and %RMSE is 22.9 ± 0.4%.

**Figure 9.** Average post-fire recovery (*%RtoMV*) trajectories for all recorded fires (≥1 km2) during the Landsat record (1988–2020). Only burned pixels (*dNBR* ≥ 100) analyzed. Red: overall average recovery. Green: Google Earth image six years after the fastest recovering fire. Blue: eight years after slowest recovering fire. Black line indicates fire polygon area for both fires [45]. Note that a fire can only reach 100 *%RtoMV* if every burned pixel comes to its *mNBR* value in the same post-fire year.

**Figure 10.** RF outputs for burn severity (*dNBR*) over 30 runs (Mean ± SD). Variable importance quantified by MIR. See Table 3 for full variable names and descriptions. %NC is the percentage of burned pixels in each fire comprised of nature conservation land use (Figures 1 and 8A).

The strength and direction of the relationships between each of these variables and *dNBR* vary at the fire scale. BL50 has a significant positive relationship with *dNBR*, with 30% variance explained (Figure 11A). This relationship is significant throughout the range of BL50; fires in areas that have burned more in the past burn more severely. When there are fires, areas that have not burned in the last 50 years burn less severely (*dNBR*: 367 ± 129) than areas that had burned once or only partially (*dNBR*: 497 ± 138), which in turn burn less severely than areas that have burned more than once (*dNBR*: 666 ± 115). *mNBR* has a significant positive relationship with *dNBR*, with 36% variance explained (Figure 11B). This relationship is most strong at lower levels of severity, with 68% variance in *dNBR* explained by *mNBR* when *mNBR* is below 450. Fires in areas with a greater vegetation capacity (*mNBR* > 600) burn significantly more severely (*dNBR*: 579 ± 139) than fires in low vegetation areas

(*mNBR* < 450; *dNBR*: 348 ± 98), but this mostly levels off past 450 *mNBR* where there is no significant difference from 600 *mNBR*. %NC has a significant positive relationship with *dNBR*, with 23% variance explained (Figure 11C). %NC forms three distinct groups, where fires in areas with little to no nature conservation land use burn significantly less severely (*dNBR*: 327 ± 144) than the two groups with higher %NC (524 ± 117 and 496 ± 167). Interestingly, the middle group of %NC, representing fires with nature conservation cover between 16% and 64%, has a strong relationship with *dNBR* by itself (R2 = 0.74\*\*\*). This may be because fires that contain an increasingly mixed land use of nature conservation and mostly residual native cover are generally larger fires that spill over into both land use settings. Relatedly, BExt has a significant positive relationship with *dNBR*, with 35% variance explained (Figure 11D). Within the relationship, there is a significant increase between fires that burn an area less than 10 km<sup>2</sup> (BExt < 1; *dNBR*: 426 <sup>±</sup> 152) and larger fires (621 <sup>±</sup> 136), especially the most extensive and severe fires from 2008 and 2020. Although BYr, TPI65 and MD marginally improve RF model performance, they do not have a significant relationship with *dNBR* by themselves.

#### *3.3. Recovery: Landscape-Level Influences*

RF modeling identifies six landscape-level variables that together most influence post-fire recovery (*%RtoMV*) (Figure 12). In the initial years 1–3, surviving vegetation (*postNBR*) is most important. This importance relative to the other variables decreases over time, especially in year 3 where eucalypt woodland percentage (EucPer), dNBR, SLB, BL40 and BYr all increase in MIR. In year 4 SLB is the most important and in year 5 EucPer is the most important. In these middle years, all included variables have a stronger relative influence on *%RtoMV* with no variable having a MIR below 0.4. In the later years 6–10, BYr is the most important. This relative influence on *%RtoMV* becomes stronger towards year 10, with all other variables having a MIR below 0.5 by year 8. Generally, *postNBR* linearly decreases in importance, while BYr follows a sigmoid shape of increasing importance from years 1–10. The other variables (EucPer, dNBR, SLB, BL40) are less important in the initial (i.e., 1–2) and later (i.e., 8–10) years but peak in importance during the middle (i.e., 4–5) years as influence on *%RtoMV* transitions from *postNBR* to BYr. Overall, model performance increases from years 1–10 with %VE increasing from 41.3 ± 1.5% to 69.9 ± 0.9% and %RMSE remaining mostly stable from 7.6 ± 0.1% to 7.0 ± 0.1%.

The relationship between each of these variables and *%RtoMV* varies in strength, direction, and through time. Burned areas with little vegetation remaining (i.e., *postNBR* < 0) recover faster than areas with more vegetation remaining (i.e., *postNBR* > 0), and these differences are most significant in post-fire year 1–3 (Figure 13A). In year 1, *postNBR* < −200 fires have recovered almost double (*%RtoMV* 28.3 ± 7.4) compared to *postNBR* > 0 fires (15.7 ± 6.6%). Higher severity (*dNBR* ≥ 400) burn areas recover faster than lower severity (*dNBR* < 400) burns, and these differences are most significant in years 3–10 (Figure 13B). In year 3, when dNBR reaches its highest influence on *%RtoMV*, higher severity burns have recovered 55.7 ± 9.4% while lower severity burns have recovered 42.8 ± 13.2.%. Areas that burn where historical fires have not occurred (i.e., BL40 = 0) recover somewhat faster than areas with frequent fires (i.e., BL40 > 1), and this difference is significant in years 4–7 (Figure 13C). Areas that burn where there was no recent previous fire (i.e., SLB 30+) recover faster than areas that had already burned recently (i.e., SLB < 30), and these differences are most significant in years 4–6 (Figure 13D). In year 4, when SLB is the most influential variable, SLB < 30 fires have recovered 53.9 ± 15.0% while SLB 30+ fires have recovered 66.0 ± 11.1%. Areas that have a higher percentage of eucalypt woodlands (i.e., EucPer > 40%) recover faster than areas that have lower amounts of this vegetation group and these differences are most significant in years 3–10 (Figure 13E). In year 5, when EucPer is the most influential variable, EucPer > 40% fires have recovered 75.7 ± 8.2% while EucPer < 10% fires have recovered 59.5 ± 13.9%. This difference is among the largest and most significant of any variable at any year. More recent fires (i.e., BYr > 2000) have recovered faster than older fires (i.e., BYr < 2000), and this difference is most significant in years 4–10 (Figure 13F). The difference between 1990s fires and 2000s fires in years 8–10 represents the most significant difference between any two groups of fires in terms of *%RtoMV*, with 1990s fires at 70.8 ± 9.1% and 2000s fires at 87.9 ± 5.0% in year 10.

**Figure 11.** Relationship between variables that improved RF model performance and *dNBR*. Only variables that have a significant relationship with *dNBR* shown. From most to least important: (**A**) BL50 (**B**) *mNBR*, (**C**) %NC, (**D**) BExt. See Table 3 for full variable names and descriptions. Each scatterplot includes a boxplot with fires split into groups based on differences within each variable. MW between group test significance is shown above each boxplot, and KW overall test significance is shown below.

**Figure 12.** RF outputs for fire recovery (*%RtoMV* years 1–10) over 30 runs for each year. (**A**) Overall model results over time (Mean ± SD). (**B**) Variable importance over time as quantified by MIR. SD bars are not included here to make the figure legible. See Table 3 for full variable names and descriptions. EucPer is the percentage of burned pixels in each fire comprised of the eucalypt woodlands vegetation group (Figure 8B).

**Figure 13.** Relationship between variables that improved model performance and *%RtoMV* for the first 10 years post-fire. (**A**) *postNBR*, (**B**) *dNBR*, (**C**) BL40, (**D**) SLB, (**E**) EucPer (% Eycalypt woodlands vegetation group), (**F**) BYr. See Table 3 for full variable names and descriptions. Each time-series shows fires split into groups based on differences within each variable. MW between group test significance is shown in a table next to each plot and is based on numbers in each legend. KW overall test significance is shown in the far-right column.

#### **4. Discussion**

#### *4.1. Fire History in the Context of the 2019–2020 Bushfires*

Remote sensing time-series results confirm that the 2019–2020 bushfires were the most widespread and devastating to occur on Kangaroo Island since the late-1980s and likely since at least the early 1930s. They were the largest in extent, burned the most severely, and affected large areas that had not experienced burns in the modern record (Figures 5–7). That the second most widespread and severe fire season was only 12 years prior should be concerning if we hope to maintain the island's unique natural environment. Both the 2020 and 2008 fires occurred in the context of increasingly hot temperatures across the island (Figure 2A). Furthermore, both fire seasons occurred after at least three years of below-average annual precipitation (Figure 2B), which undoubtedly led to dry and fire-prone environments.

Although local and national authorities do not have direct control of temperature and precipitation, they can mitigate damage through the implementation of fuel-reduction burning, firebreaks, and other forest management techniques [7]. The relative lack of fire, especially in the decade before 2020, likely contributed to fuel build-up and higher burn severity (Figures 6 and 7). Although a local project on the island recommended fire intervals of 17–40 years [35], scientific studies on eucalypt forests

elsewhere in southern Australia have found that shorter intervals of 6–10 years, through prescribed burns, are needed to effectively reduce fire risk [34,61]. Prescribed burns have been recorded on the island since the early 2000s, but only at small scales (32.6 km2 recorded; Figure 6B) and as part of trial projects [62]. More involved forest management practices may be required to prevent future mega-fires. This is especially true since the four most common tree species on Kangaroo Island (*E. diversifolia*, *E. remota*, *E. cladocalyx* and *E. baxteri*) are fast growing (leading to large fuel buildup), highly flammable (from combustible oils and leaves and dense litter) and have the capacity to withstand and recover from fires (with lignotubers or epicormic buds). Prescribed burns can be used to limit fuel buildup, allow for controlled burning of flammable material and approximate a natural part of the life cycle for these fire-adapted species.

As in another Landsat time-series study of forest recovery in south-east Australia [26], post-fire NBR for Kangaroo Island fires recovers and stabilizes after about ten years (Figure 9). This is notably faster than in the western United States [16,28], Spain [12], boreal Canada [28] and elsewhere. These studies come from many different environments and employ different methods but our results indicate that Kangaroo Island eucalypt forests are relatively well adapted to fire and can quickly assume the regeneration process [35]. Depending on the specific eucalypt species impacted, this may occur with resprouting from lignotubers that protect buds from fire damage or from epicormic buds that sprout post-fire. However, even fire-tolerant eucalypts may not be able to survive higher severity fires, as shown from a study in south-eastern Australia [63]. These fires produce greater seedling regeneration but also a higher likelihood of eucalypt mortality, with the potential to eventually transition eucalypt forests to more open woodland or savannah environments. NBR recovery from the large 2008 fires (93 ± 6% after 10 years) indicates that eucalypt forests may be able to recover in the same manner after the 2020 fires. However, the unprecedented extent and severity of the 2020 fires should still be concerning because they may kill so many trees that recovery from re-sprouting becomes inadequate at the landscape-level, which means that eucalypt forests may not be able to recover their pre-fire density [63]. To end, it should be noted that we did not see any significant impact of specific eucalypt species (Table 2) on results at this scale, but did see significantly higher severity in nature conservation land use (which is mostly mallee and eucalypt woodlands) and faster recovery in the eucalypt woodland vegetation group.

#### *4.2. Influences on Burn Severity*

Included variables explain well over half of the variation in burn severity (i.e., *dNBR*). Severity is influenced mostly by fire frequency (i.e., BL50), vegetation capacity (i.e., *mNBR*) and land use (i.e., %NC) (Figure 10), with nature conservation areas of higher fuel loads and more historical fires burning more severely (Figure 11). Burn extent also has moderately strong and significant relationships with severity when compared individually, but RF results show that in a multivariate scenario, much of those relationships are explained by *mNBR*, BL50 and %NC. The influence of vegetation capacity on severity is not surprising, as ground-level and remote sensing analysis have come to similar conclusions elsewhere [37]. On Kangaroo Island, Peace and Mills [46] explain that Flinders Chase NP has both dense eucalypt canopies (esp. *E. remota*) and understory vegetation, especially in riparian areas, which led to high fuel loads and thus the severity of the 2008 fires.

Flinders Chase NP and other nature conservation land use areas (Figure 1) overlap with the mallee woodlands vegetation group (Figure 8). In this vegetation group, eucalypt species grow in multiple stems with a shorter and shrubby growth form [42]. The abundance of mallees in these areas is not surprising, as they burn often (Figure 5) are highly flammable and this type of growth form is common when lignotuber eucalypts are resprouting [64]. It is also possible that frequent burning leads to an increase in low canopy biomass for combustible mallees and understory shrub species, which would provide amble fuel for higher severity fires. Furthermore, the lack of modification and interference in nature conservation land use (esp. wilderness areas), means fuel build up is less managed [43]. It is thus likely that the importance of %NC in influencing dNBR is related in part to both fire frequency (i.e., BL50) and increased biomass fuel (i.e., *mNBR*).

In terms of fire history, studies in western North America have found that burn severity is lower when time since last fire (i.e., SLB) is lower, depressing severity even 20–30 years later [38,39]. We do not see this on Kangaroo Island, where *dNBR* is above average in areas that burned less than 10 years previously and have experienced multiple fires in the last 50 years. This difference can likely be explained by the relatively quick recovery of vegetation observed (Figure 9). However, our results do broadly correspond with Estes et al. [36], who found that severity was highest 30–50 years after the previous fire. The five highest severity fires on the island occurred in areas that burned between 1 and 3 times in the last 50 years and 24–36 years after the previous fire, with the somewhat shorter interval potentially related to faster regrowth of fire-adapted eucalypts. These studies do not consider a longer timeline in terms of the number of fires over multiple decades (i.e., BL50), which was a more important influencer than SLB.

The ecological reasoning for higher *mNBR* leading to more severe burns is clear and related to increased fuel loads in these areas, as discussed previously. This type of connection for BL50 is less clear, but it may be that areas subject to more frequent fires contain species (i.e., *Eucalyptus* spp.) that are better adapted to them and even require higher severity burns (although not too high) for most productive bud resprouting, seedling regeneration and subsequent regrowth [63]. Whatever the reasoning, these results provide valuable information to help mitigate future mega-fires like that of 2020. A key requirement for mitigation should be the limitation of fuel build-up through active pre-fire management methods [7,37]. However, it may not be feasible to manage all densely vegetated areas. Therefore, we recommended forest managers be especially conscious of fuel build up in those densely vegetated areas that have also burned many times in the past and to prioritize these types of settings if resources are limited. Furthermore, care should be given to limiting fire expansion from nature conservation areas (where managers may have less control on things like fuel buildup) to other settings. In 2020, bushfires were able to make the jump from Flinders Chase NP to forest plantations and residual native cover areas in the central part of the island, where corridors of remnant eucalypt forest likely facilitated fire expansion into large swaths of agricultural land.

#### *4.3. Influences on Post-Fire Recovery*

Included variables explain much of the variation in post-fire recovery (i.e., *%RtoMV*), with %VE increasing in later years. Recovery can be split into three periods with different landscape-level influences (Figures 12 and 13). In the initial years (i.e., 1–3), surviving vegetation (i.e., *postNBR*) most dictates recovery, with areas that have less remaining vegetation recovering more quickly. *postNBR* is significantly related to burn severity (R<sup>2</sup> = 0.45 \*\*\*\*), with higher severity burns usually leaving less remaining vegetation. Here, we find that higher *dNBR* burns recover more quickly than lower *dNBR* burns, supporting those studies that have come to similar conclusions [16,23,29]. However, *postNBR* outcompetes *dNBR* for importance in our RF analysis, especially in the initial years, indicating that more attention should be paid to recovery in different post-fire environments separate from severity.

In the middle years (i.e., 4–5), time since last burn and eucalypt woodlands percentage is most influential, but this appears to be a transitional period. After this, variables related to the past situation of the burn (i.e., *postNBR*, *dNBR*), previous fire history (i.e., SLB, BL40), and vegetation structure (i.e., EucPer) decline in relative influence on recovery. However, the strong relationship between eucalypt woodlands abundance and recovery rate for most post-fire years (Figure 13E) deserves more discussion. Eucalypt woodlands are common in the most densely vegetated interior areas of Flinders Chase NP and along riparian corridors here and across Kangaroo Island. We have previously established that eucalypts burn severely and that severe burns are associated with faster recovery. Eucalypt woodlands have a denser canopy than mallee woodlands [42] and their abundance in riparian areas provide both more opportunity for survival to recovery and better access to moisture [65], allowing for a relatively fast transition from burn scar to recovering mallees to denser eucalypt woodlands. That eucalypt

woodlands are common in riparian corridors may explain why certain topographic variables like topographic position and wetness did not influence results.

It is in the later years (i.e., 6–10) that BYr becomes the clearly dominant influencer with 1990s fires recovering significantly slower than more recent fires. This indicates that the current (relative to burn year) context for recovering areas takes precedent over past disturbances. Bright et al. [16] found, nine years post-fire, that precipitation and temperature variables explained much of the variation in recovery. On Kangaroo Island, it also appears that climate change between fires in different time periods is responsible for variations in later year recovery. In the western portion, where most fires have occurred, summer temperatures 6–10 years after the faster recovering 2000s fires were 1.6 ◦C warmer on average, and annual temperatures 1◦ C warmer, than the later recovery years for 1990s fires (Figure 2A). Total annual precipitation 6–10 years after the 2000s fires (693 mm) averaged 76 mm more annually than the same period for 1990s fires (617 mm), a 12% increase (Figure 2B). It appears that BYr may represent a proxy of differences in climate through time, becoming most influential in later recovery years. Malak and Pausas [29] found that precipitation begins to directly influence post-fire vegetation only after seven years in a Spanish dryland setting. Eucalypt forests in this region are also moisture limited, meaning that short-term precipitation increase in the 2010s has likely contributed to the faster recovery of more recent fires, even in the broader context of precipitation decrease in the region due to global climate change [65].

As Kangaroo Island recovers from the 2019–2020 bushfires in the coming decade, historical recovery trends and drivers can provide some insights into how this may occur, with the caveat that their larger extent and higher severity are unprecedented in the modern record (Figures 6 and 7). Based on our findings, land managers could expect recovery to be fastest in the first few years post-fire, especially where it burned more severely, and there is little surviving vegetation. Severely burned eucalypt woodlands growing along riparian corridors are prime candidates for fast recovery based on our results. Parts of the severely burned Flinders Chase NP is already going through early recovery stages with bud resprouting, epicormic growth and fast-growing eucalypts and shrub species beginning to green up just weeks after the burn [66]. The impact of previous burn history should also be considered, with areas that have experienced more frequent and recent fires (i.e., higher BL40 and lower SLB) recovering slower than less burned areas. This should be especially concerning when you consider that climate change in south-east Australia is expected to further increase fire frequency in the coming decades [5], which may lead to heavily burned areas failing to recover their historical tree canopy density or transitioning to a new type of ecosystem [31]. In the long term, it is climate, especially precipitation, which appears to most influence the recovery rate.

#### *4.4. Limitations*

This study made use of Landsat imagery and other freely available data sources to monitor the recent fire history of Kangaroo Island. As such, we did not conduct field work, but did collect information on vegetation types and species composition [42]. As our and other results show, vegetation type in Australia does play a role in remotely sensed burn severity [32] and recovery [22]. However, it is not a static variable with factors like fire frequency potentially leading to lasting vegetation type change [31]. Even though vegetation type changes are likely minor over our 33-year study period, we should note that differences between fire dates and data dates for land use and vegetation groups may have some impact.

Climate was not directly incorporated into our analysis, which instead focused on landscape-level variables related to the current burn, fire history, vegetation structure, and topography. We did not anticipate that climate would play a major role since the analyzed fires occurred on an island with a mostly consistent climate. However, considering these fires occurred in the context of a changing climate (Figure 2), it becomes clear there is an influence on the probability of a severe fire occurring and on long-term post-fire recovery. A more robust examination of climatic influences on fire history in the region, especially in the context of a warming and drying climate, would be beneficial.

Topographical variables had only a minor influence on burn severity (i.e., TPI65) and no influence on post-fire recovery. Other studies, mostly concerned with recovery, have found topographical influence from strong [31] to weak [16] to none [12]. The lack of topographical influence in this study may be due to the landscape-scale of the analysis. When variables, including topography, are averaged at the fire level, the spatial variability within the burn area is lost. Having the unit of analysis being the fire itself is advantageous for building broader-scale understanding and monitoring large numbers of fires, but future research may consider combining this full fire history approach with finer-scale analysis within each fire. Recent studies have considered the spatial variability within multiple fires in their analysis, but with a smaller number of fires than analyzed here [16,24]. Future research may consider splitting individual fires into vegetation group, topography, and other subsets, while still analyzing a large number of fires to build more robust conclusions.

Missing data (i.e., MD), especially in the early years of analysis (Figure 3), was also a concern. If no Landsat imagery is available over a fire the year that it burns, LandTrendr will fit the decline in NBR based on the next available year (Figure A1). Some vegetation will have regrown in the intervening period, leading to the expectation that these fires will be artificially less severe than others even if adjusting *dNBR* years account for this. Despite this, MD only has a minor influence on severity in RF analysis and a comparison between those fires where part of the burn area had missing data, and other fires shows no significant difference in *dNBR* (Figure A4A). More vegetation may have also regrown depending on the fire date, compared to Landsat image collection date. Fires dated before December 25 (i.e., shorter time until Landsat date) were compared to fires dated after March 25 (i.e., longer time) and there is a minor difference in *dNBR* that is not significant (Figure 4B). A combined approach, highlighting fires that have either missing data or are dated after March 25 similarly show a minor, but not significant difference in *dNBR* (Figure A4C). MD does not influence post-fire recovery and fire date had no influence on either severity or recovery RF models. In future, near-real time (as opposed to yearly) fitting algorithms for Landsat such as Stochastic Continuous Change Detection can be used to limit the time gap between fire disturbance and post-fire detection [67].

#### **5. Conclusions**

Landsat images through time, supplemented with other spatial datasets available at no cost, provide useful observations about recent fire history and associated drivers. On Kangaroo Island, we found a recent (1988–200) history of frequent fires, with the 2008 and 2020 seasons being most widespread and severe. The fire-adapted eucalypt-dominant vegetation on the island shows a remarkable capacity to regrow, with %RtoMV stabilizing in only ten years on average. However, there is cause for concern, as those areas that have most frequently burned are also recovering more slowly. As climate change leads to warmer and dryer conditions, more frequent fires may curtail the recovery of these forests and potentially trigger transitions to more open woodlands or entirely new ecosystems. More generally, we found that burn severity was most influenced by variables representing fire frequency, vegetation capacity and land use (i.e., nature conservation). Land managers interested in reducing severity should consider fuel-reduction practices where allowed, especially in areas that have been subject to more frequent fires. Furthermore, we found varying influences on yearly post-fire recovery, with variables representing surviving vegetation important in early years and an inferred proxy of post-fire climate important in later years after a transitional period where the time since the previous burn and eucalypt woodland percentage was important. In other eucalypt environments, we might expect those most severely burned areas with less remaining vegetation to have faster initial recovery and in later years for growth to be guided by precipitation in moisture limited settings.

Results from this study suggest that even without access to plot-level information or other fine-scale data, a comprehensive understanding of fire history and influences on severity and recovery processes can be obtained. This type of analysis can be conducted efficiently, freely, and remotely, especially at present with travel restrictions due to COVID-19. The cloud computing and big data capabilities of GEE and the trajectory fitting approach of LandTrendr prove to be an effective combination for

building landscape-level understandings of environmental change, especially when those changes are recent, and there is a need for timely analysis.

**Author Contributions:** Conceptualization, M.T.B., Y.H. and S.W.M.; methodology, M.T.B.; software, M.T.B.; validation, M.T.B.; formal analysis, M.T.B.; investigation, M.T.B. and Y.H.; resources, M.T.B. and Y.H.; data curation, M.T.B.; writing—original draft preparation, M.T.B.; writing—review and editing, M.T.B., Y.H. and S.W.M.; visualization, M.T.B.; supervision, Y.H.; project administration, Y.H.; funding acquisition, M.T.B. and Y.H. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the Natural Sciences and Engineering Research Council of Canada (NSERC) through the NSERC Discovery Program awarded to Dr. Yuhong He. It was also supported through the NSERC Postgraduate Scholarship, Centre for Urban Environments (University of Toronto Mississauga), and Department of Geography, Geomatics and Environment funding to Mitchell Bonney.

**Acknowledgments:** Analysis was conducted at the Remote Sensing and Spatial Ecosystem Modeling laboratory in Mississauga. We thank lab members for their advice and support throughout the project. We thank Justin Braaten for his correspondence related to LandTrendr and GEE. We thank Benjamin Bright for his correspondence related to NBR trajectories.

**Conflicts of Interest:** The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

#### **Appendix A**

**Figure A1.** Example pixel that experienced a fire in December 1996. Black dots are medoid surface reflectance values, and the red line is the LandTrendr fitted model as in Figure 4. Vertical dashed black lines indicate years with no available imagery. The *postNBR* year 1997 is missing, so LandTrendr fits a two-year decline to 1998. In this example, *preNBR* is 1996 and *postNBR* is now 1998.

**Figure A2.** Violin plots showing the distribution of burn severity (*dNBR*) around the mean for all recorded fires. Split into six periods, as in Figure 7.

**Figure A3.** Closer view of two select fires (P3\_7, P2\_0), with *dNBR* set to 75% transparency. P3\_7 burned in February 2001 at the northern edge of Flinders Chase NP. P2\_0 burned in December 1994 just west of Vivonne Bay. In both cases (and for other fires that burn in these types of landscapes), *dNBR* values above 100 almost exclusively occur in treed and dense shrubby land cover. Grasslands, agriculture and other non-woody land cover, although within the fire polygons [45] have very limited influence. Note that the ArcGIS basemap imagery over these areas are dated to 2019 (P3\_7) and 2018 (P2\_0) so there may have been landscape change in intervening period.

**Figure A4.** Influence of Landsat image time-series on *dNBR*. (**A**) Missing data: No = little to no influence of missing data on *dNBR* years. Yes = Notable portion of burn area impacted by missing data in *dNBR* years. (**B**) Fire dates [45]: < Dec 25 = fire dated before December 25 in fire season. > Mar 25 = fire dated after March 25 in fire season. Between = fire dated between December 25 and March 25 in fire season. Because image collection dates were December 25–March 25, the fires dated before these dates get a *postNBR* value soon after the burn while the fires dated after these dates get a post *NBR* value in the next year. The *postNBR* lag for fires dated during the image collection window depends on the medoid image that was selected over that area. (**C**) Combination: Yes = fires that either were impacted by missing data or are dated after March 25. These fires have the longest lag between burn date and *postNBR* date. No = all other fires.

#### **References**


**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).
