The Iterative Process.

An exploratory statistical analysis of Canadian Wildfires

---

Canadian wildfires have become a subject of increasingly confident claims. Climate change is blamed for worsening fires; forestry practices blamed for changing fuel loads; herbicides, particularly Glyphosate, have became a recurring subject of debate; and question about firefighting capacity inevitably follow every catastrophic fire season.

Some claims are well established. Others are plausible hypotheses. Some are difficult to distinguish from one another using datasets that are actually available.

This article is an attempt to interrogate those questions statistically.

It is not intended to produce a grand unified model of Canadian wildfires, nor to establish causal effects for every variable included in the analysis. The data is provincial and depending on the specification, annual or monthly. This data cannot capture every relevant feature of a wildfire landscape, and some variables are necessarily imperfect proxies for the things that we would ideally like to measure.

The interesting question isn’t whether a regression produces a statistically significant coefficient. It is what that coefficient actually represents, whether the model is capable of distinguishing competing explanations, and what happens when we change the way the problem is specified. This is not so much an article on Canadian wildfires using statistics, but rather an article about statistics, statistical modelling and analysis using Canadian wildfires as the case study.

As such, in this article, I’ll walk through the modelling process, rather than simply present the final results. Some variables that initially seemed useful turned out to introduce endogeneity. Annual aggregation obscured relationships that became much clearer at monthly granularity. In one particularly interesting case, a statistically significant relationship between glyphosate use and wildfire outcomes turned out to be a good demonstration of why statistical significance doesn’t tell us what the variable actually means.

The objective is less to discover a previously unknown theory of Canadian wildfire behaviour, and more to see what the data can reasonably tell us, what it cannot, and how easily an apparently convincing statistical result can be made to say more than the underlying data can justify.

---

Pre-emptive caveat

Before getting into it, It’s worth underscoring that while these are deliberately simple models using relatively accessible data, they are very early, preliminary steps in an iterative exploratory modelling process, rather than toy models deliberately constructed to highlight statistical pitfalls and traps.

You start simple and test the relatively straightforward options first. You run the models, interpret the results, and determine what works, what doesn’t, what is missing, and what would be required to address those gaps. You don’t skip straight to a massively complex spatial model with tens of millions of observations before establishing what the simpler models can tell you. With that being said…

Conceptual framework

The primary objective of this exercise is to quantify the structural, ecological, and climatic drivers of Canadian wildfire regimes. Specifically, it aims to isolate the factors governing ignition frequency (fire counts) and fire extent (area burned) across both natural and human-caused ignition sources to inform evidence-based forestry and suppression policy.

The core conceptual framework rests on a fundamental premise: ignitions do not equal extent. The variables that dictate whether a fire ignites (e.g. fine fuel moisture, atmospheric instability, human access) are distinct from those that govern whether an ignition escapes initial attack and expands into a mega-fire. Differentiating between natural and human ignition regimes is equally critical.

• Natural (lightning) regimes: driven by convective atmospheric instability and landscape-scale fuel dynamics across both managed and unmanaged wilderness.

• Human regimes: strictly tethered to human infrastructure, industrial operations, and recreational access. Humans cause fires where humans exist, not in isolated wilderness.

Using a 35-year panel dataset (1990-2024) across Canadian provinces and territories, I estimated Two-Way Fixed Effects (TWFE) regression models across four distinct target variables.

1. Natural area burned | area burned by naturally ignited fires.

2. Human area burned | area burned by human-ignited fires.

3. Natural fire counts | total number of nature-caused ignitions.

4. Human fire counts | total number of human-caused ignitions.

The two-way fixed effects model

1. Provincial Fixed Effects: accounts for permanent baseline differences between provinces (e.g. mountainous terrain in British Columbia as opposed to flatlands in the prairies).

2. Year Fixed Effects: accounts for nationwide annual weather events and climate trends that affect all provinces simultaneously.

Econometric specification

The TWFE framework absorbs unobserved heterogeneity through two dimensions:

\(\text{Outcome}_{it} = \beta_0 + \mathbf{X}_{it}\boldsymbol{\beta} + \gamma_i + \delta_t + \varepsilon_{it}\)

Provincial fixed effects Controls for time-invariant baseline spatial differences between jurisdictions (e.g. mountainous topography in British Columbia versus flat prairie terrain or lake-dense Boreal Shield).

Year Fixed Effects: Controls for nationwide annual macroeconomic, climate, or synoptic weather shocks that affect all jurisdictions simultaneously.

Variable selection

• log_burned_natural: Log hectares burned by natural (lightning) ignitions.

• log_burned_human: Log hectares burned by human ignitions.

• log_num_natural: Log count of nature-caused fire ignitions.

• log_num_human: Log count of human-caused fire ignitions.

• Mean Temperature (mean_temp): captures macro thermal conditions. Higher temperatures increase vapour pressure, drying out forest fuels, while promoting unstable convective atmospheric conditions for lightning.

• Total Precipitation (total_precip): captures fuel moisture saturation, this is expected to suppress fire ignition, though the impact on large infernos is limited once fuels cross critical dryness thresholds.

• Mean Wind Speed (mean_wind_speed): contributes to fire spread.

• Current Insect Defoliation (log_insect_defol): insect outbreaks kill leaves and needles, creating abundant dry fuel first in the canopy, and later on the ground.

• 1-Year Lagged Insect Defoliation (log_insect_defol_lag1): attempts to capture multi-year canopy dynamics. By year two, dead leaves and needles should be on the ground, breaking the canopy ladder, but theoretically increasing spread on the ground.

• Ground Herbicides (log_herb_ground): ground-level vegetation management targets below-canopy brush, shrubs, and invasive growth. This includes glyphosate usage to address present discourse.

• Aerial Herbicides (log_herb_helicopter, log_herb_fixed_wing): typically deployed for dense, remote areas which carry heavy fuel loads. Includes glyphosate usage to address present discourse.

• Wood Supply / Harvest Volume (log_wood_supply): proxies active commercial harvesting, introduces access roads, and proxies human access.

• Lagged Targets (log_burned_*_lag1, log_num_*_lag1): Large fires alter fuel loads, open landscapes, and leave standing dead wood that influences ignition and spread in the following season. For human fires, lag terms capture structural persistence in human recreation and industrial access.

Omitted variables and Endogeneity constraints.

Suppression Mobilization Data: inter-provincial Mutual Aid Resource Sharing (MARS) statistics (aircraft, personnel, pumps) were omitted due to severe simultaneous endogeneity. Resource mobilization is a direct function of fire extent; high fire activity triggers increased resource sharing, creating a false positive correlation between suppression deployment and total area burned.

Budgetary metrics: Annual provincial suppression and forestry budgets (~470 budget documents across 13 jurisdictions and Parks Canada) were excluded. Budget-to-actual variance severely penalizes jurisdiction facing anomalous, black-swan fire seasons (ex. Quebec in 2023, Alberta in 2016). Instead, silvicultural treatments, commercial harvest volumes, and fleet assets serve as cleaner proxies for active preventative management.

---

Annual Model Results

The annual TWFE model yields high overall explanatory power but reveals stark divergences between ignition frequency and total burn acreage:

Standardized coefficients are calculated after within-province z-scoring of the dependent and independent variables. They therefore represent standardized within-province associations and are not conventional globally standardized coefficients.

What the models appear to show us

• Natural Area Burned: The natural-area model has a high R-squared of 0.8765, conventionally interpreted as the model “explaining” approximately 88% of the observed variance. The quotation marks are deliberate: R-squared describes how much variation is accounted for by the fitted model, not whether the included variables actually explain that variation causally.
  
  Mean temperature has a positive and statistically significant relationship with naturally burned area (raw coefficient = 0.2459, standardized β = 0.1269, p = 0.0248).The direction is consistent with the expected relationship between warmer conditions and wildfire activity, although perhaps tempting the coefficient should not be read as meaning that a 1°C increase produces a ~24.6% increase in burned area (despite what it suggests on the surface), the model does not actually show this; the dependent variable is log-transformed, and the model is a fixed-effects panel specification, so the coefficient represents a conditional association within the model rather than a straightforward percentage prediction. While it can be interpreted mathematically, it should not be interpreted as a causal effect.
  
  Annual precipitation is not statistically distinguishable from zero (β = −0.0496, p = 0.7032). This does not demonstrate that precipitation has no effect on wildfire behaviour; rather annual aggregation may obscure relationships that may occur on shorter timescales.
  
  The strongest statistically significant relationship in this model is the lagged burned-area variable (β = 0.21, p = 0.0013). Provinces with greater burned area in the previous year tend to have greater burned area in the current year, conditional on the other variables and fixed effects. There are several potential explanations for this; regional climatic conditions, landscape characteristics, changes to fuel structure, etc; but the annual model cannot distinguish between them.

• Human-caused area burned: The human area model has an R-squared of 0.7629. Neither mean temperature (p = 0.4207) nor precipitation (p = 0.2563) reaches conventional statistical significance thresholds. None of the forestry treatment variables are statistically significant, either.
  
  The lagged human-burn variable has a positive coefficient (β = 0.1072, p = 0.0922), providing weak evidence at the 10% level of temporal persistence, but not at the conventional 5% threshold.
  
  The important result is therefore partly what the model does not find: at annual provincial resolution, human-caused burned area is not clearly explained by the broad climatic variables included in the model.

• Natural fire counts: here, the R-squared is 0.9077, making it one of the highest-fitting models in the exercise. Mean temperature has a positive relationship with with natural fire counts (β = 0.2167, p = 0.0024), while precipitation has a negative one (β = −0.1997, p = 0.0313). Both are statistically significant.
  
  Current year insect defoliation is also positively associated with natural fire counts (β = 0.1322, p = 0.0182), while lagged defoliation has a negative relationship (β = −0.1197, p = 0.0084). This is an interesting pattern, but the model does not establish the mechanism behind it. We can interpret it as insect damage and its effects on fuels change over time, but the coefficients alone cannot tell us whether that is actually what is taking place.
  
  Two herbicide columns are also statistically significant. Helicopter-applied herbicide has a positive relationship with natural fire counts (β = 0.1106, p = 0.0432), and ground-applied herbicide has a negative one (−0.0953, p = 0.0312).
  
  These coefficients are particularly useful as an illustration of why statistical significance needs to be interpreted cautiously. The model has identified relationships, but it does not tell us that herbicide application itself causes either more or fewer fires. The variables may be associated with particular types of forestry activity, geography, access, land management, or other characteristics that are not fully represented in the data.

• Human-caused fire counts: the R-squared here is 0.9233, the highest fit of the four models. The dominant relationship is temporal persistence: the previous year’s human-fire counts has a standardized coefficient of β = 0.4963 (p < 0.001). Provinces with more human-caused fires in one year, tend to have more in the following years. As above, human-caused fires occur where there are humans. The relationship may be proxying for human presence.
  
  Precipitation is negatively associated with human fire counts (β = −0.2218, p = 0.0172). Mean temperature, wind speed, insect defoliation, forestry production, and the herbicide variables are not statistically significant in this specification.
  
  The absence of a significant temperature relationship is noteworthy, in particular when compared with the natural fire model. It is consistent with the basic observation that human ignitions require human activity: their geographic and temporal distribution is constrained by where people, infrastructure, recreation, and industrial activity occur. Weather may affect whether those ignitions grow or persist, but the annual count of human-caused ignitions is not simply a function of broad provincial temperature trends.

---

Caveats and limitations: the next steps

We will perhaps have more, better insight be adapting the model for more granular data, using monthly, instead of annual averages. With monthly aggregates we can gleam more nuance that specifically targets the wildfire season, and annual averages hide temperature and precipitation impacts. For example, a 2C average temperature increase in spring may suggest a potentially earlier thaw and extended fire season, the same average increase in July is significant as well, but these are diluted over annual aggregates. Similarly, increased snowfall during winter may dilute reduced precipitation in relevant months.

Unfortunately not all of our dataset is available as monthly breakdowns. I was able to find fire number and area data by province, as well as monthly climate data. There are two options to proceed using hybrid data:

1. Average out annual data we already have, to a monthly basis. This is to be avoided because it’s not real data. OR

2. Adapt the extant model to accommodate one-to-many relationships between monthly and annual data.

Option two is preferred as it does not invent data that isn’t in our dataset.

As well, a curious relationship emerged between ground-applied herbicides and natural fire counts, in this case, aside from Quebec, natural fires are far more common than human-caused, and ground-applied herbicides are typically applied in managed areas where there is road access. Because out dataset is not spatially-aware, it will not capture these distinctions.

---

Standardized coefficients in the monthly models are calculated using the conventional pooled-sample standardization, β* = β(σX/σY). Unlike the annual models, the monthly regressions do not z-score variables separately within each province.

---

Shifts at monthly granularity

Moving from annual to monthly observations, even partially, produces an immediate and superficially undesirable result: R-squared falls substantially. The four monthly models range from 0.35 to 0.61, compared with 0.76 to 0.92 in the annual models.

This is not necessarily a deterioration in the model. A higher R-squared is not inherently evidence of a better model, particularly when comparing different levels of aggregation. The annual models capture broad provincial and year-level patterns that are less prominent at monthly resolution. The monthly models instead allow us to examine whether relationships that appear important in the annual data visible when the data is disaggregated.

There are significant changes from the annual models:

1. Weather becomes much more prominent. In the natural area model, mean temperature has a standardized coefficient of β = 0.4883 (p < 0.001), while precipitation has a negative coefficient of β = −0.1307 (p < 0.001). Both are statistically significant at monthly resolution, whereas precipitation was not significant in the annual model. This is consistent with the possibility that annual aggregation obscures short-term weather relationships; a particularly warm or dry period during fire season can have a substantial relationship with fire activity even when those conditions are diluted by the rest of the year.
   
   Temperature is also strongly associated with natural fire counts β = 0.5214, p < 0.001), while precipitation is negatively associated with them (β = −0.1581, p < 0.001). The monthly specification therefore produces a much clearer relationship between contemporaneous weather conditions and wildfire activity than the annual specification.

2. The intertemporal relationship changes substantially. In the annual models, last year’s burned area strongly predicts this year’s burned area (β = 0.21, p = 0.0013). In the monthly model, the corresponding lag loses statistical significance (β = 0.0348, p = 0.2860). This suggests that the persistence observed at annual resolution may reflect broader landscape or climatic conditions rather than a simple short-term carryover effect. A large fire during one period, for example, may alter fuel availability in the area immediately afterward, while broader regional conditions continue to influence fire activity over longer periods. The model cannot distinguish between these mechanisms, but the change between temporal resolutions is itself informative.

3. Human-caused burned area becomes much more weather sensitive. Unlike the annual model, the monthly human-area model finds statistically significant relationships for temperature (β = 0.2894, p < 0.001), precipitation (β = −0.2480, p < 0.001), and wind speed (β = 0.1805, p = 0.0094).
   
   This does not establish that any of these variables “govern” human-caused fires. Rather, it suggests that weather conditions have a much clearer association with the amount of area burned by human-caused fires when examined at fire season timescales. One plausible interpretation is that weather affects whether human caused ignitions develop into larger fires, even if the underlying occurrence of those ignitions is primarily determined by human activity and geography.

4. Forestry variables produce additional relationships at monthly granularity. In the natural area model, biological treatment (log_bt) is statistically significant (β = 0.1026, p = 0.007), as is helicopter-applied herbicide (β = 0.1780, p < 0.001). In the human-area model, fixed-wing herbicide is negatively associated with burned area (β = −0.1543, p = 0.0038), while helicopter herbicide is positively associated with it (β = 0.1295, p = 0.0085)..
   
   These results are useful primarily as examples of the interpretive problem explored throughout this article. They identify statistical relationships, but the model cannot establish that the treatments themselves cause these changes. The variables may instead be associated with particular forestry practices, locations, access conditions, fuel characteristics, or other factors not represented in the dataset. The changing coefficients between specifications are therefore more interesting than any simplistic claim that a particular treatment increases or decreases wildfire activity.

5. Human fire counts show a similar weather relationship while retaining strong temporal persistence. Temperature (β = 0.2664, p < 0.001), precipitation (β = −0.1910, p < 0.001), and wind speed (β = 0.1252, p = 0.0234) are all statistically significant. The lagged human-fire count is also strongly associated with current counts (β = 0.3693, p < 0.001). This is consistent with persistent differences in human activity, infrastructure, recreation, and industrial presence between jurisdictions, although the model cannot distinguish among those explanations.

In brief:

Annual granularity

• Higher R-squared

• Broad provincial and year-level patterns

• Weather relationships can be diluted by aggregation

• Stronger annual persistence

Monthly granularity

• Lower R-squared

• Stronger contemporaneous weather relationships.

• Some relationships appear or disappear with temporal resolution

• Forestry variables produce additional relationships requiring careful interpretation.

The important result is not that the monthly models are “better” because they have lower or higher R-squared. They answer a somewhat different question. Disaggregation exposes relationships that annual averages may obscure, while also demonstrating that relationships appearing strong at one level of granularity may weaken at another.

That is arguably the more important finding: the level at which we choose to measure a phenomenon can materially change the relationships the model appears to find.

When a one-month lag isn’t a one-month lag

There is another, considerably more mundane problem hidden in the monthly specification.

The monthly dataset was restricted to the active fire season, May through September before autoregressive lagged variables were calculated. The model then used a one-period shift(1) within each province to construct the lagged dependent variables.

At first glance, this appears to create a one-month lag as one would expect, but there’s a catch: It doesn’t.

June’s lag is May, July’s is June, and so on, so forth, but for May, the lagged value is the previous September. The final relationship is not a one-month lag at all, rather it’s the end of the previous active fire season.

This is not a particularly exotic statistical problem. It is a mundane data-processing mistake. The code did exactly what it was instructed to do: shift(1) returned the previous observation. The previous observation is not what we’d expect, and it skews the results and fundamentally changes the question being posed.

This may nevertheless yield interesting results, September-to-May persistence could reflect overwintering fuel conditions, landscape state, or other process spanning the offseason. But again that’s an entirely different hypothesis from short-term monthly persistence.

The important point is that nothing in the regression output would necessarily tell us that this had happened. The variable has a perfectly plausible name, produces perfectly plausible numbers, and can produce a statistically significant coefficient. Without inspecting how it was constructed, it would be very easy to interpret it as a one-month lag simply because it was labelled lag1.

This is precisely why iterative statistical analysis requires interrogation of the data-generating process, rather than treating the regression table as the end of the analysis. The model can be mathematically correct, but nevertheless answer a subtly different question from what its author intended.

The implication is not that the monthly model is invalid, overwinter conditions can be useful. The implication is that the result needs to be interpreted according to the variable that was ultimately constructed, rather than what the author intended. Correcting the lag would be an obvious next iteration, but comparing the two specifications can be informative: one measures short term persistence, the other incidentally measures persistence across the fire-season boundary.

The regression cannot tell us that we constructed the variable incorrectly. We have to ask the regression what variable we actually gave it. Mistakes are part of the iterative process, they sometimes produce useful results, if we’re aware of them.

---

The Glyphosate question and spatial proxy confounding

To Evaluate ongoing policy claims regarding silvicultural herbicide bans and wildfire suppression, 010_flyphosate_quebec_null_hypothesis.py isolated glyphosate application (log_glyphosate) as a separate variable, across a 2001-2024 panel, using Quebec as an empirical benchmark following its 2001 ban of public forest land.

The raw data includes herbicide use by type and method of application, the dataset aggregates them into each other for totals. The obvious step is to separate out glyphosate usage into its own variable, but as we’ll see, there are numerous issues and pitfalls.

econometric pitfalls: spatial omitted variable bias

Running the baseline panel yields a statistically significant negative coefficient for glyphosate (Beta = -0.1452, p - 0.0194), suggesting that higher glyphosate use correlates with lower burn acreage. Reiterating, on the surface, it looks like the model is telling us that glyphosate actually suppresses fire spread. However, this is a textbook example of proxy confounding,

The glyphosate variable is inadvertently acting as a proxy for several other factors; namely forestry management and other preventative measure, absorbing their explanatory power.

1. The managed tenure trap: Silvicultural glyphosate is applied exclusively within managed forestry zones. As such, log_glyphosate acts as a spatial proxy for intense forestry infrastructure, road networks, clear site access, proximity to initial-attack ground crews, often well withing scooper or water-bomber range, and active suppression monitoring.

2. Macro vs. Micro distortion: Provincial-level models compare managed forestry zones against millions of hectares of unmanaged northern wilderness (with no, or limited access and no herbicide treatment). The model misattributes suppression access to chemical application.

3. Quebec empirical evidence: Arguments that Quebec’s low average fire incidence is mono-causally linked to its 2001 glyphosate ban fail under historical analysis. Post-ban years such as 2002, 2013, and 2023 saw massive burned area, while several pre-ban years recorded exceptionally low fire counts. Single-variable causality cannot hold across historical volatility. Similarly, pre-ban years such as 1991, 1995, and 1996 had fewer fires than 2023, but nowhere near as much area burned.

Even if we were to account for zoning (wilderness vs managed zones) we’d still likely see log_glyphosate acting as a false positive when using Quebec as a null hypothesis simply because in Quebec, human-caused fires, which typically occur in managed areas are far more common than natural ignitions, while the trend is the inverse in other regions.

The core issue is that simply separating out glyphosate application into a separate variable, isn’t sufficient to actually isolate its effects. A far more complex model and dataset would be required, perhaps one employing spatial grid data, identifying each grid tile as manage or not, in addition to ascertaining how to isolate glyphosate and its impacts from all other factors, all at a grid-level. This adds a whole lot of complexity for little explanatory power, but is perhaps worthwhile in a model dedicated to determining glyphosate’s impact on forestry.

---

Brief structural and operational view of provincial assets

This section is a brief glimpse at figures which seem intuitive and may seem more explanatory than they actually are. Forest area, water bomber fleets, number of fires and extent for 2025, and coverage per bomber in square kilometres. This section serves to demonstrate that the recurring traps and interpretative hurdles mentioned throughout this article aren’t exclusive to statistical modelling and analysis, they apply to simple data tables as well.

At first glance, one would assume that Quebec, which has the worst area coverage of any jurisdiction would perhaps have worse results than Manitoba or NWT which have much better coverage density. However, both have substantially more area burned on significantly fewer fires than Quebec does.

Simply looking at something like coverage obscures that in Quebec human-caused ignitions (in managed areas) are far more frequent than natural ignitions, which facilitates rapid response, with the province’s scooper fleet taking advantage of lake density in the affected areas for fire suppression at the micro level before they spread. The higher spatial load becomes a vulnerability in times of black swan weather events (such as 2023) where multiple simultaneous remote fires quickly overwhelm suppression capacity.

The Northwest Territories, however, has among the highest area burned despite being at the low end of ignition frequency. Because NWT is sparsely populated and comprised mainly of uninhabited wilderness, suppression efforts are likely prioritised in remote inhabited communities and infrastructure, while wilderness fires are left to burn. The spatial burden is similar to other eastern provinces, but the operational logic is wildly different.

Otherwise intuitive seeming metrics like km2/bomber don’t tell us this story. It’s not very useful on its own.

---

Conclusions

The purpose of this exercise was not to construct a definitive model of Canadian wildfire behaviour. It was to examine what can be learned from relatively simple models, how those models can potentially mislead, and how the process of interrogating their results can inform the next stage of research and analysis.

Several of the results illustrate why statistical modelling is an iterative process rather than a search for the single correct specification.

The first lesson is that model fit is not the same as model usefulness. The annual models produced substantially higher R-squared values than their monthly counterparts. That does not make them inherently better models. Annual aggregates can absorb broad structural differences and produce impressive explanatory power while simultaneously obscuring relationships that occur in more nuanced, more granular datasets. The monthly models explained less variance, but in several cases produced relationships that were more interpretable once the analysis was restricted to the fire season. A lower R-squared can thus represent a more useful model if it better corresponds to the relationship being investigated.

The second lesson is that statistical significance does not explain a variable. The glyphosate exercise provides a useful example. the model produced a statistically significant relationship, the kind that could drive click-bait headlines, but that result did not establish the mechanism suggested by the variable’s name. Instead, it raised the possibility that glyphosate usage was correlated with a collection of unmeasured characteristics of managed forestry areas. The important finding was therefore not that the coefficient was positive or negative, but that the available data could not cleanly distinguish the competing explanations.

This is also why apparently intuitive variables require scrutiny. A measure such as waterbomber area coverage can appear to be an obvious indicator of suppression capacity, yet the statistic itself does not contain information about operational doctrine, geography, water access, deployment priorities, response times, or the distribution of fires. The variable may still contain useful information, but its interpretation depends on what it actually measures.

The appropriate response to these problems is not necessarily to discard the model. It is to ask what the result tells us about the next iteration.

If a variable appears to be acting as a proxy, what additional variables would allow the underlying mechanisms to be separated? If annual aggregation produces an ambiguous relationship, would monthly or weekly resolve it? If provincial averages conceal important spatial variation, would grid-level data be necessary? If a suppression metric appears counterintuitive, what operational information would be required before it could be interpreted meaningfully?

This exercise provides a practical argument for starting small.

It is possible to construct increasingly complex and elaborate models involving hundred of thousands of observations, extensive environmental covariates, detailed operational records, and large collections of administrative data. Such models may be necessary when the objective is to isolate a specific causal mechanism. They are not necessarily required for discovering which questions require that level of complexity in the first place. A simpler model can serve as a diagnostic instrument to identify relationships worth investigating, demonstrate where aggregation hides information, expose variables that may be unintended proxies, and offer insight into what hidden variables are likely to matter.

Simpler models also raise questions in terms of what to investigate and what could be interesting to look into and how. For instance, in this exercise, the question arises as to why in Quebec, human-caused fires are far more common than naturally ignited ones, and what the operational implications are. It also reveals what modelling glyphosate effects may require. They can serve as an indicator not only of how the model can be refined, but of what avenues of research to pursue.

Moreover it exposes an analytical trap. While these simple models paint a picture that is in line with mainstream research; that climate factors are the main driver for increased wildfire frequency, extent, and severity, it doesn’t establish that link nor statistically demonstrate it. Alignment with consensus is not a substitute for methodological rigour. Much more granular and extensive grid data, along with substantially more intricate and complex modelling is required to make that determination. Failure to make this distinction drives the common view that statistics can be cajoled and manipulated into painting any picture; a simple narrative focusing only on variance coverage and statistical significance can for example, result in pushing that glyphosate usage actively contributing to reduced fire incidence and extent, but we don’t actually know that, not from these models, at least.

Statistical analysis is therefore not simply a process of selecting variables, running regressions, and accepting coefficients that meet a significance threshold. It is an iterative process of specification, observation, criticism, and refinement. The model is tested against the data, the results are tested against the conceptual framework, and discrepancies between the two become information for the next iteration.

Canadian wildfire data happens to provide an unusually good case study for this process because the underlying system is sufficiently complex to generate all of these problems at once: temporal aggregation, spatial heterogeneity, proxy variables, omitted variables, endogenous responses, operational confounding, and substantial differences between what a variable measures and what we would ideally like it to measure.

The objective, ultimately, is not to find the model that produces the most impressive statistics. It is to understand what those statistics actually tell us.

“Technological advance is an inherently iterative process.
One does not simply take sand from the beach and produce a Dataprobe.
We use crude tools to fashion better tools, and then our better
tools to fashion more precise tools, and so on. Each minor refinement
is a step in the process, and all of the steps must be taken.”

— Chairman Sheng-ji Yang, “Looking God in the Eye”

Replicating the exercise

Project files are available on GitHub

003_timeseries_regressions_no_mobilised.py runs the annual regression.
008_direct_monthly_adaptation.py runs the monthly hybrid regression.
010_glyphosate_quebec_null_hypothesis.py runs the glyphosate regression.

This codebase represents an iterative exploratory process. While the primary regression scripts documented in the article are fully functional and reproducible, several helper scripts, legacy data-transformation passes, and feature-engineering experiments remain in the repo. Some of these served as temporary stepping stones wand I don’t remember what they do.

Data acquisition

Number of fires by cause: http://nfdp.ccfm.org/en/data/fires.php#:~:text=3.2.3%20Area%20burned%20by%20fire%20size%20class.,1990.%201991.%201992.%201993.%201994.%201995.%201996.

Area burned by cause: http://nfdp.ccfm.org/en/data/fires.php#:~:text=3.2.3%20Area%20burned%20by%20fire%20size%20class.,1990.%201991.%201992.%201993.%201994.%201995.%201996.

Number of fires by size: http://nfdp.ccfm.org/en/data/fires.php#:~:text=3.2.3%20Area%20burned%20by%20fire%20size%20class.,1990.%201991.%201992.%201993.%201994.%201995.%201996.

Area burned by size: http://nfdp.ccfm.org/en/data/fires.php#:~:text=3.2.3%20Area%20burned%20by%20fire%20size%20class.,1990.%201991.%201992.%201993.%201994.%201995.%201996.

mobilization stats: https://ciffc.ca/mobilization-stats/

Climate data: data is collected in 10,000 entry increments from https://api.weather.gc.ca/collections/ahccd-annual/items?f=csv&datetime=1990-01-01/2024-12-31&limit=10000 and https://api.weather.gc.ca/collections/ahccd-monthly/items?f=csv&datetime=1990-01-01/2024-12-31&limit=10000&offset=300000 . Stitched into the first file, removing the first line of headers with “awk ‘FNR>1’ file2 file3 … fileN » consolidated.csv”.

bilingual AHCCD headers were mapped to English headers, monitoring station data was grouped and aggregated by province and year. arithmetic mean was computed across all reporting station in a province for each metric, skipping NaN entries to avoid skewing. Derived station_count by calculating the number of unique active stations identifiers (station_id) reporting in which province for a given year. Final results were filtered to retain only data between 1990 and 2024.

Because data is only available until 2024, a python script was used to loop through the station IDs and hit the climate-daily endpoint on the weather.qc.ca API to fill gaps.