Moisture-driven fire early warning for the Bogota Cerros Orientales and the southern Amazon
Moisture-driven fire early warning for the Bogota Cerros Orientales and the southern Amazon
What this system is for
The operating purpose is early warning as Bogota and the Amazon begin to lose soil moisture. Fires are the visible event, but by the time a fire is visible the useful decision window has closed. What can be forecast at a useful lead time is not the fire itself but the moisture state that makes fire cheap to start and expensive to stop. This system therefore forecasts a moisture index, calibrates that index against observed fire so the forecast carries a fire meaning, and converts the result into dated probabilities that a preparedness office can act on.
The argument for doing it this way is economic as much as scientific. Prevention that is scheduled weeks ahead, when crews can be positioned, water points filled, access roads cleared and burn bans timed, costs a fraction of what suppression and post-fire repair cost once a season has gone bad. The information required is already free and public. The gap this work addresses is translation: turning open climate reanalysis and open fire archives into a dated, calibrated, honestly scored statement about the coming season.
Scope is fourteen study sites on four continents, weekly, from 1950 to the present, with the two target basins being the Bogota Cerros Orientales and a six-site southern Amazon composite. The twelve other sites serve three purposes: they supply candidate predictors, they provide out-of-region controls for the calibration, and they let the coherence analysis test whether the Bogota and Amazon signals are locally generated or part of a wider mode.
What was explored, in order
- Construction of a soil moisture index that behaves in a saturating tropical climate.
- Calibration of that index against fire counts and against burned area, site by site.
- Validation that the seasonal target actually separates bad fire seasons from ordinary ones.
- Seasonal forecast models for the target, with an honest cross-validated skill score.
- Wavelet coherence over all 91 site pairs, to find shared modes and their lead times.
- Direct coupling tests between the Amazon and Bogota, with three confounding guards.
- Northern hemisphere teleconnection, asking whether Canada, Europe and the western United States carry information about the southern seasons.
- A crossing model that converts a seasonal probability into dated probabilities.
- Long-term trends in moisture and in burned area, with an availability control.
Sources, records and their limits
2.1 Soil moisture
ERA5-Land volumetric soil water in layer 1, the top 0 to 7 cm, monthly and hourly aggregates reduced to a weekly site mean over each study polygon. The record runs from February 1950 to May 2026 in the monthly product and is resampled onto a contiguous weekly grid of 3981 steps running from 1950 week 1 to 2026 week 29. The 0 to 7 cm layer is chosen deliberately: it is the layer that dries fastest, that governs fine fuel moisture, and that responds within days rather than months. Deeper layers are smoother and would look more skilful while carrying less fire information.
2.2 Fire
Three fire products are used, each for a different job.
- MODIS active fire detections supply the weekly event counts used in the Poisson calibration. They are the only product with weekly resolution across the whole modern period, and they are the product with the worst small-fire sensitivity.
- FireCCI51 burned area, monthly, 2001 to 2020, supplies the burned area calibration and the Amazon seasonal response variable. It sees area, not events, and it is the product that behaves in the Amazon where the detection counts do not.
- Landsat burn scars derived per season from 1984 onward supply the long Cerros Orientales burned area series used in the trend section. This is the series that requires a scene availability control.
Supporting long records were catalogued for extension but are not in the operational chain: NBAC Canada 1972 to 2025, Geoscience Australia 1899 to 2024, MTBS United States 1984 to 2024, MapBiomas Fogo Brazil 1985 to 2024, FireCCILT11 1982 to 2018 with 1994 missing, UAECOB Bogota fire brigade polygons 2009 to 2025, and DesInventar Colombia 1914 to 2018 with roughly 9300 forest fire records.
2.3 Ocean state
CPC Oceanic Nino Index, three month running means of Nino 3.4 sea surface temperature anomaly, from the DJF season of 1950 through AMJ 2026. Used to form four seasonal windows and one tendency term, defined in section 4.
2.4 Sites
Fourteen polygons. Two target regions: Bogota_Cerros, and a southern Amazon composite. The composite has two definitions used throughout:
The five site core drops Guaviare, which sits at the northwestern edge and behaves more like the Colombian piedmont than like the southern arc. The remaining seven sites are Med_Europe, US_West, Canada_boreal, Siberia, Central_Africa, Indonesia and SE_Australia.
Every number quoted in this report was recomputed directly from the archived CSV outputs rather than copied from an intermediate summary. Where a figure differs from an earlier draft, this document supersedes it.
A normal score index that survives a saturating wet season
The obvious index is a Gaussian z-score of volumetric soil moisture against a climatological mean and standard deviation. In the tropics this fails. Wet season soil moisture in the Cerros and in the Amazon arc saturates against the soil porosity ceiling, so the distribution is strongly left-skewed with a hard right edge. The standard deviation is then dominated by the dry tail, and a wet week produces an absurd positive z-score while the dry tail is compressed. The early version of this analysis had exactly that defect.
The replacement is a rank-based normal score in the spirit of the SPI. For each site
and each ISO week, a climatological pool is built from the 1991 to 2020 WMO normal
period, pooling weeks k-2 through k+2 circularly across the
year boundary. That is 30 years times 5 weeks, so 150 samples per pool. The observation
is ranked within its pool using the Weibull plotting position, and the resulting
empirical percentile is passed through the inverse normal.
BASE0, BASE1 = 1991, 2020 # WMO normal
WIN = 2 # +/- weeks pooled into each week's climatology
idx = [((k - 1 + d) % 52) + 1 for d in range(-WIN, WIN + 1)] # circular
def normal_score(row):
p = pools.get((row.site, row.wk))
v = row.sm
if p is None or not np.isfinite(v):
return np.nan
n = len(p)
# Weibull plotting position, then probit, capped where the tail is
# unresolvable by a sample of this size.
q = (np.searchsorted(p, v, side='right') + 0.5) / (n + 1.0)
q = min(max(q, 0.5 / (n + 1.0)), 1 - 0.5 / (n + 1.0))
return float(np.clip(stats.norm.ppf(q), -3.5, 3.5))
Three properties matter operationally. The index is unitless and comparable across
sites and across weeks of the year, so a value of -1.0 in the Cerros in February means
the same percentile as -1.0 in Bolivia in September. The +/-2 week pooling
smooths the seasonal cycle without smoothing the anomaly, which a per-week climatology
built from only 30 samples could not do stably. And the clip at +/-3.5 is
an honest statement that a 150 sample pool cannot resolve the tail beyond roughly the
1 in 300 level; values are censored there rather than extrapolated.
Throughout this report the index is called the moisture index and is quoted in sigma
units. The column sm_z in data/weekly_sm_z_allsites.csv is
this quantity. A three week centred rolling mean of it, with
min_periods=2, is used wherever a smoothed index is needed, in particular
for the crossing rule in section 9.
3.1 Fixed reference levels
| Level | Value | Meaning |
|---|---|---|
| Watch line | -0.5 | Drier than roughly 69 per cent of weeks in the base period. |
| Action line, DRY | -1.0 | Drier than roughly 84 per cent of weeks. This is the fixed level used to count dry weeks. |
| Severe | -1.5 | Drier than roughly 93 per cent of weeks. |
These are properties of the normal score construction, not tuned parameters. The seasonal critical thresholds in section 4 are a separate and empirically derived quantity.
Eight ordinary least squares models with nested leave-one-season-out selection
4.1 Seasons and targets
Two fire seasons are targeted. Bogota runs December to March, defined as ISO weeks
wk >= 48 or wk <= 13 and indexed by the year in which the season
ends, so the season opening in December 2026 is fyear 2027. A season is admitted only
if at least 16 of its weeks are present. The southern Amazon runs September to
November, weeks 35 to 48, indexed by calendar year, admitted at 12 weeks or more.
Both give 76 usable cycles.
Four targets are forecast:
| Target | Definition | Critical |
|---|---|---|
bog_mean | Season mean moisture index, Cerros Orientales | -0.18183 |
bog_dry | Count of season weeks with index at or below -1.0 | 3 |
amz_mean | Season mean of the six site composite amz6 | -0.15423 |
amz_core | Season mean of the five site core amz5 | -0.19718 |
The critical thresholds are not chosen by hand. Each is the one third quantile
of the observed distribution of that target over the 76 cycles, so by construction
about one season in three is critical. For the count target the threshold is
ceil(quantile(bog_dry, 2/3)), which evaluates to 3. Defining the threshold
by quantile rather than by judgement is what makes the AUC and Brier scores in section
4.5 interpretable: the base rate is fixed and known.
4.2 Candidate predictors
Predictors are seasonal window means of the moisture index at various sites, plus five ocean index terms. All windows close strictly before the issue date, so no model sees its own season.
jj = wmean(D, yr, 18, 29, [...]) # June to mid July, suffix _jj
ya = wmean(D, yr, 1, 29, [...]) # year to date, suffix _ytd
au = wmean(D, yr, 30, 34, [...]) # August, suffix _aug
sn = wmean(D, yr, 35, 48, [...]) # September to November, suffix _son
oni_mam = mean ONI over Mar, Apr, May
oni_mjj = mean ONI over May, Jun, Jul
oni_jja = mean ONI over Jun, Jul, Aug
oni_ason = mean ONI over Aug, Sep, Oct, Nov
oni_tend = ONI(May) - ONI(February) # the tendency termThe tendency term deserves a note. It is not a level, it is the sign and speed of the ocean's move through the northern spring. A weak but rising Nino carries different information for the following southern dry season than a strong but decaying one, and the level alone cannot express that. It survives selection in five of the eight models.
Candidate pools are declared per basin and per issue date, and are deliberately small:
POOLS = {
('bogota', 'august'): ['oni_mam', 'oni_tend', 'Bogota_Cerros_jj',
'Guaviare_jj', 'amz6_jj', 'Canada_boreal_jj',
'Indonesia_jj', 'bog_ytd'],
('bogota', 'september'): ['oni_mjj', 'oni_tend', 'bog_aug', 'gua_aug',
'amz6_aug', 'Canada_boreal_jj', 'Indonesia_jj'],
('bogota', 'december'): ['oni_ason', 'oni_jja', 'bog_son', 'gua_son',
'amz6_son', 'amz5_son'],
('amazon', 'august'): ['oni_mam', 'oni_tend', 'amz6_jj', 'amz5_jj',
'Guaviare_jj', 'amz6_ytd', 'Bogota_Cerros_jj'],
('amazon', 'september'): ['oni_mjj', 'oni_tend', 'amz6_aug', 'amz6_jj',
'gua_aug'],
}
MAXP = 34.3 Estimator
Plain ordinary least squares with an intercept. No standardization, no regularization, no interaction terms. With 76 observations and at most three predictors this is the right complexity, and it keeps the coefficients directly readable in index units.
def fit(Xy, ycol, cols):
d = Xy.dropna(subset=cols + [ycol])
X = np.c_[np.ones(len(d)), d[cols].values]
beta, *_ = np.linalg.lstsq(X, d[ycol].values, rcond=None)
return beta, d4.4 Nested cross validation, and why the nesting matters
Forward selection from an eight member pool on 76 observations will find something that looks skilful whether or not skill exists. If selection is done once on the full record and the resulting model is then cross validated, the reported skill is contaminated: every held out point participated in choosing the predictors. The protocol here redoes the entire selection inside every fold.
def loo_select(Xy, ycol, pool, maxp=MAXP):
"""Nested leave one season out. Selection is redone for every held out
season, so the returned skill has no selection bias in it."""
d = Xy.dropna(subset=pool + [ycol]).reset_index(drop=True)
n = len(d); pred = np.full(n, np.nan); picks = []
for i in range(n):
tr = d.drop(index=i)
chosen = []; best_prev = np.inf
while len(chosen) < maxp:
cand, cbest = None, np.inf
for c in pool:
if c in chosen: continue
cols = chosen + [c]
e = []
for j in range(len(tr)): # inner LOO on training set only
t2 = tr.drop(index=tr.index[j])
beta, _ = fit(t2, ycol, cols)
e.append(tr[ycol].iloc[j] - predict(beta, tr.iloc[j], cols))
s = float(np.sqrt(np.mean(np.square(e))))
if s < cbest: cand, cbest = c, s
if cand is None or cbest > best_prev - 1e-4: break
chosen.append(cand); best_prev = cbest
if not chosen: chosen = [pool[0]]
beta, _ = fit(tr, ycol, chosen)
pred[i] = predict(beta, d.iloc[i], chosen)
picks.append(tuple(chosen))
return d, pred, picks
The outer loop holds out one season. The inner loop cross validates each candidate
addition on the training seasons only. A predictor is added only if it improves the
inner root mean square error by more than 1e-4, which is what stops the
models short of the three predictor cap when the third adds nothing. Because selection
is inside the fold, the by-product picks is a genuine measure of how
stable each model's structure is, reported in section 4.6.
Scoring, once the out of sample predictions exist:
The count target uses a half unit continuity correction because bog_dry
is an integer and the threshold of 3 means "3 or more". The AUC is the Mann-Whitney U
statistic normalised, tested one sided in the greater direction. The climatological
reference for skill is the standard deviation of the target itself, and the
climatological Brier reference is p(1-p) at the observed base rate.
4.5 The eight fitted models
Coefficients from the full record fit. Only the 1 August models are issued live today.
bog_mean (1 Aug) = +0.031 +0.368*Bogota_Cerros_jj -0.259*oni_tend -0.436*Canada_boreal_jj
bog_mean (1 Sep) = +0.134 -0.419*oni_mjj -0.417*Canada_boreal_jj -0.190*oni_tend
bog_mean (1 Dec) = +0.068 -0.482*oni_ason +0.479*bog_son +0.314*oni_jja
bog_dry (1 Aug) = +2.552 -1.631*Indonesia_jj +1.448*oni_tend -1.612*bog_ytd
bog_dry (1 Dec) = +2.174 +2.747*oni_ason -0.840*gua_son -1.617*oni_jja
amz_mean (1 Aug) = -0.022 +0.633*amz6_jj -0.098*oni_tend
amz_mean (1 Sep) = +0.042 +0.469*amz6_aug -0.214*oni_mjj +0.198*gua_aug
amz_core (1 Aug) = -0.043 +0.651*amz6_jj +0.136*Bogota_Cerros_jjThe signs are physically coherent. Wetter antecedent conditions at the target site predict a wetter season, a rising ocean index predicts a drier Bogota season, and a wet Canadian boreal spring predicts a dry Bogota season, which is the northern hemisphere link discussed in section 8. The Amazon models are dominated by their own antecedent moisture with the ocean playing only a supporting role, which is the clearest structural contrast between the two basins in the whole analysis.
Out of sample statistics, all with n = 76:
| Model | r | CV rmse | clim rmse | skill | events | AUC | Brier | Brier skill |
|---|---|---|---|---|---|---|---|---|
| bog_mean 1 Aug | +0.482 | 0.539 | 0.604 | +0.108 | 25/76 | 0.760 | 0.185 | +0.164 |
| bog_mean 1 Sep | +0.674 | 0.447 | 0.604 | +0.260 | 25/76 | 0.784 | 0.174 | +0.213 |
| bog_mean 1 Dec | +0.764 | 0.392 | 0.604 | +0.352 | 25/76 | 0.885 | 0.134 | +0.395 |
| bog_dry 1 Aug | +0.496 | 2.505 | 2.869 | +0.127 | 27/76 | 0.799 | 0.184 | +0.198 |
| bog_dry 1 Dec | +0.663 | 2.152 | 2.869 | +0.250 | 27/76 | 0.825 | 0.168 | +0.266 |
| amz_mean 1 Aug | +0.437 | 0.623 | 0.680 | +0.084 | 25/76 | 0.691 | 0.207 | +0.064 |
| amz_mean 1 Sep | +0.768 | 0.437 | 0.680 | +0.358 | 25/76 | 0.906 | 0.124 | +0.437 |
| amz_core 1 Aug | +0.438 | 0.656 | 0.717 | +0.086 | 25/76 | 0.730 | 0.192 | +0.132 |
Highlighted rows are the models issued live on 1 August 2026.
4.6 Predictor set stability
Because selection runs inside every fold, the distribution of selected sets across the 76 folds is a direct diagnostic. A model whose set changes every fold is fitting noise.
| Model | Modal set | Predictors | Runner up |
|---|---|---|---|
| bog_mean 1 Aug | 42/76 | Bogota_Cerros_jj, oni_tend, Canada_boreal_jj | 16/76 swaps oni_tend for Indonesia_jj |
| bog_mean 1 Sep | 75/76 | oni_mjj, Canada_boreal_jj, oni_tend | 1/76 same set, different order |
| bog_mean 1 Dec | 52/76 | oni_ason, bog_son, oni_jja | 24/76 swaps oni_jja for gua_son |
| bog_dry 1 Aug | 74/76 | Indonesia_jj, oni_tend, bog_ytd | 1/76 each, two others |
| bog_dry 1 Dec | 75/76 | oni_ason, gua_son, oni_jja | 1/76 swaps oni_jja for bog_son |
| amz_mean 1 Aug | 32/76 | amz6_jj, oni_tend | 25/76 selects amz6_jj alone |
| amz_mean 1 Sep | 75/76 | amz6_aug, oni_mjj, gua_aug | 1/76 swaps gua_aug for oni_tend |
| amz_core 1 Aug | 39/76 | amz6_jj, Bogota_Cerros_jj | 18/76 selects amz6_jj alone |
Five of the eight models select the same set in 52 or more folds out of 76, and two of
them in 75. The two 1 August Amazon models are the least stable, and the instability is
informative rather than alarming: in every fold amz6_jj is selected first,
and the disagreement is only over whether a second term earns its place. In 25 folds
out of 76 the model reduces to a single predictor. That is consistent with the low
continuous skill of those models and is one reason the 1 September reissue is a large
step up rather than a marginal one.
4.7 The live forecast issued 1 August 2026
| Target | Forecast | 80 % interval | Critical | P(critical) | AUC |
|---|---|---|---|---|---|
| bog_mean | -0.4291 | -1.1197 to +0.2615 | -0.18183 | 0.6769 | 0.760 |
| bog_dry | +4.5318 | +1.3216 to +7.7420 | 3 | 0.7914 | 0.799 |
| amz_mean | -0.3598 | -1.1581 to +0.4385 | -0.15423 | 0.6293 | 0.691 |
| amz_core | -0.2553 | -1.0956 to +0.5851 | -0.19718 | 0.5353 | 0.730 |
Predictor values entering the live forecast:
Bogota_Cerros_jj -0.009 Guaviare_jj +0.025 oni_mam +0.540
Canada_boreal_jj +0.248 amz6_jj -0.324 oni_tend +1.120
Indonesia_jj -0.513 bog_ytd +0.500 amz6_ytd -0.349
The dominant live signal is oni_tend at +1.120, a strongly rising ocean
index through the northern spring, together with a dry Indonesia and a dry southern
Amazon in June and July. Bogota's own antecedent moisture is close to neutral, which
is why the Bogota forecast rests more on remote predictors than on local memory this
cycle. The four probabilities are the input to the crossing model in section 9.
What one sigma of dryness is worth, in fires and in hectares
A moisture forecast is not a fire forecast until the index has been given a fire meaning. Two separate calibrations do this, one on event counts and one on burned area, and they are structurally different models because the two response variables have different statistics.
5.1 Poisson rate model on weekly detection counts
Weekly MODIS detections at a site are counts with a strong seasonal cycle. The model separates the cycle from the anomaly by putting the weekly climatology in as a fixed offset:
The offset coefficient is fixed at 1, not estimated. That is the whole point of
the specification: it forces the model to explain only the departure from the seasonal
norm, so b1 is a pure anomaly effect and cannot absorb seasonality. Fitting
is by hand-rolled iteratively reweighted least squares, 80 iterations, weight matrix
W = mu, ridge stabiliser 1e-8 * I, converging at
max|db| < 1e-9.
Detection counts are heavily overdispersed relative to Poisson, so standard errors are inflated by a quasi-Poisson dispersion estimate rather than taken from the naive Poisson information matrix:
Ignoring this would understate the standard errors by a factor of tens to nearly a
hundred at the large Amazon sites, where phi exceeds 8000. Off-season
weeks are excluded by keeping only weeks whose climatology is at least 10 per cent of
the site maximum, and a site is skipped if it has fewer than 200 site-weeks or fewer
than 50 total detections, then skipped again if fewer than 150 weeks survive the
off-season filter. 18928 site-weeks were matched in total.
Derived quantities, reported per unit of the index:
| Site | weeks | fires | slope | se | z | phi | rate ratio | double at | triple at |
|---|---|---|---|---|---|---|---|---|---|
| Bogota_Cerros | 541 | 1169 | -0.7562 | 0.1075 | -7.03 | 11.7 | 2.130 | -0.92 | -1.45 |
| Guaviare | 312 | 343937 | -0.8811 | 0.0521 | -16.90 | 671 | 2.414 | -0.79 | -1.25 |
| Amazon_arc | 631 | 13824370 | -0.1727 | 0.0222 | -7.78 | 6186 | 1.188 | -4.01 | -6.36 |
| Amazon_S | 631 | 19314792 | -0.2315 | 0.0220 | -10.53 | 8050 | 1.260 | -2.99 | -4.75 |
| Xingu | 659 | 1803206 | -0.1807 | 0.0302 | -5.99 | 1332 | 1.198 | -3.84 | -6.08 |
| Cerrado | 633 | 14502539 | -0.2413 | 0.0192 | -12.59 | 4595 | 1.273 | -2.87 | -4.55 |
| Bolivia | 605 | 1807319 | -0.5415 | 0.0284 | -19.10 | 1207 | 1.719 | -1.28 | -2.03 |
| Med_Europe | 972 | 3785487 | -0.3162 | 0.0214 | -14.78 | 1580 | 1.372 | -2.19 | -3.47 |
| US_West | 1149 | 3301627 | -0.2509 | 0.0192 | -13.05 | 912 | 1.285 | -2.76 | -4.38 |
| Canada_boreal | 485 | 2457494 | -0.7014 | 0.0668 | -10.51 | 10673 | 2.017 | -0.99 | -1.57 |
| Siberia | 511 | 5647767 | -0.6999 | 0.0433 | -16.16 | 8718 | 2.013 | -0.99 | -1.57 |
| Central_Africa | 975 | 61056025 | -0.1979 | 0.0127 | -15.59 | 7527 | 1.219 | -3.50 | -5.55 |
| Indonesia | 998 | 3915359 | -0.7902 | 0.0196 | -40.38 | 1177 | 2.204 | -0.88 | -1.39 |
| SE_Australia | 1251 | 2375240 | -0.4047 | 0.0283 | -14.33 | 1602 | 1.499 | -1.71 | -2.71 |
All fourteen slopes are negative and all fourteen are significant, the weakest at z = -5.99. This is the single most robust result in the project: drier means more fire, everywhere, on every continent, with no exceptions and no sign flips.
The magnitudes carry the operational content. Bogota's slope of -0.756 means the fire rate doubles at index -0.92, which is essentially the -1.0 action line, and triples at -1.45. Indonesia, Guaviare, Canada boreal and Siberia are similarly steep. The large Amazon interior sites are much shallower, requiring three to four sigma to double, which reflects that they are ignition limited and land use driven rather than purely moisture limited: dryness there permits fire that people then have to light.
Observed Bogota rate against the seasonal norm, by index band, as a direct check on the fitted curve:
| Index band | weeks | observed rate vs norm |
|---|---|---|
| -4.0 to -1.5 | 40 | 3.73× |
| -1.5 to -1.0 | 53 | 1.81× |
| -1.0 to -0.5 | 84 | 0.89× |
| -0.5 to 0.0 | 104 | 0.86× |
| 0.0 to +0.5 | 83 | 0.64× |
| +0.5 to +4.0 | 177 | 0.42× |
Monotone across all six bands, with the step change occurring between the -1.0 to -0.5 band and the -1.5 to -1.0 band. This is independent empirical support for placing the action line at -1.0 rather than at -0.5.
5.2 Burned area, a deliberately different model
Burned area is not a count and should not be modelled as one. FireCCI51 is monthly, so
weeks are mapped to months by
((week - 1) * 7 // 30.4 + 1).clip(1, 12). The monthly climatology is then
divided out to form a ratio against the seasonal norm, and a log-linear fit is taken on
that ratio with a floor to keep zeros finite:
slope, intercept = np.polyfit(zz, np.log(np.maximum(ratio, 1e-3)), 1)
mult_per_unit = np.exp(-slope)
rho, p = stats.spearmanr(zz, ratio)A Spearman correlation is reported alongside the parametric slope, because the log transform with a floor is a modelling choice and the rank statistic is not. Gates are fewer than 60 records or zero burned area, then fewer than 40 after filtering.
| Site | n | slope | multiplier per unit | Spearman rho | p |
|---|---|---|---|---|---|
| Bogota_Cerros | 60 | -1.4889 | 4.43× | -0.310 | 1.6e-2 |
| Guaviare | 60 | -1.0752 | 2.93× | -0.489 | 7.3e-5 |
| Amazon_arc | 80 | -0.3485 | 1.42× | -0.452 | 2.6e-5 |
| Amazon_S | 80 | -0.4247 | 1.53× | -0.481 | 6.2e-6 |
| Xingu | 80 | -0.5654 | 1.76× | -0.502 | 2.1e-6 |
| Cerrado | 80 | -0.5417 | 1.72× | -0.539 | 2.5e-7 |
| Bolivia | 80 | -0.8107 | 2.25× | -0.541 | 2.3e-7 |
| Med_Europe | 80 | -0.3149 | 1.37× | -0.386 | 4.1e-4 |
| US_West | 120 | -0.3809 | 1.46× | -0.320 | 3.7e-4 |
| Canada_boreal | 60 | -1.1251 | 3.08× | -0.477 | 1.2e-4 |
| Siberia | 100 | -0.7637 | 2.15× | -0.510 | 6.2e-8 |
| Central_Africa | 120 | -0.2465 | 1.28× | -0.491 | 1.2e-8 |
| Indonesia | 160 | -1.5779 | 4.84× | -0.718 | 1.2e-26 |
| SE_Australia | 180 | -0.8098 | 2.25× | -0.352 | 1.3e-6 |
Note also that Bogota's Spearman of -0.310 is the weakest of the fourteen even though its slope is the second steepest. The parametric fit is being carried by a small number of very dry months with very large ratios, and the rank statistic properly discounts that. Both numbers are reported for exactly this reason.
Does the seasonal target separate bad fire seasons from ordinary ones
Sections 4 and 5 establish that the target can be forecast and that the weekly index predicts weekly fire. Neither establishes that a season classified as critical by the one third quantile rule is actually a worse fire season. That is a separate test and it is the one that licenses the whole product.
Four tests are run per basin: a Spearman on the season index against season fire, a Spearman on the dry week count against season fire, a group mean comparison with a ratio, and a one sided Mann-Whitney of the critical group against the rest.
6.1 Bogota, 26 seasons with MODIS coverage, 2001 to 2026
| Test | Statistic | p |
|---|---|---|
| Spearman, season index vs detections | -0.544 | 0.00406 |
| Spearman, dry week count vs detections | +0.502 | 0.00893 |
| Mean detections, critical (11) vs not (15) | 64.6 vs 22.3, ratio 2.89× | — |
| Mann-Whitney, one sided | — | 0.00110 |
All four agree. A season the rule calls critical has on average 2.89 times the fire
detections of one it does not, and the separation is significant at p = 0.0011. The
dry week count test matters independently, because it is the count target
bog_dry that is forecast at 1 August, and it confirms the count carries
fire information in its own right rather than only as a proxy for the mean.
The season list is instructive. The eleven critical seasons are 2010, 2003, 2024, 2004, 2020, 2001, 2016, 2007, 2005, 2002 and 2013. Seven of the eleven exceed 47 detections. The two clear misses are 2024, which was called critical at index -0.850 with 10 dry weeks but recorded only 19 detections, and 2023, which was called not critical at +0.560 and recorded 86. Those two are honest failures of the rule and are not hidden here. The 2024 case is discussed in section 11 alongside the detection sensitivity problem, since 2024 was the largest burned area year on record in the Cerros at 48303 hectares while returning a low detection count.
6.2 Amazon, 25 years, 2001 to 2025, and why the response variable changed
| Test | Statistic | p | verdict |
|---|---|---|---|
| Spearman, index vs MODIS detections | -0.348 | 0.0878 | fails |
| Spearman, index vs FireCCI burned hectares (n = 20) | -0.814 | 1.29e-05 | passes strongly |
| Mean detections, critical (12) vs not (13) | 1497442 vs 1231219, ratio 1.22× | 0.101 | fails |
| Mean burned hectares, critical vs not | 28368321 vs 13739153, ratio 2.06× | — | passes |
Individual years support the reading. 2020 sat at index -1.192 with 26.98 million hectares burned, 2010 at -0.468 with 39.39 million, and 2007 at -0.302 with 45.86 million, against 2009 at +0.010 with only 4.92 million and 2001 at +0.393 with 10.70 million. The detection counts across those same years vary by less than a factor of three while the burned areas vary by nearly a factor of ten.
Shared modes and lead times across all 91 site pairs
Correlation over a 76 year record answers one question with one number. Wavelet coherence answers a better question: at which timescales do two sites share variance, when in the record did they share it, and which one leads. Everything in this section runs on the same weekly index used everywhere else, 3981 contiguous steps.
7.1 Transform
Continuous wavelet transform with a Morlet mother wavelet, following Torrence and Compo. Constants:
W0 = 6.0 # Morlet order
FOURIER = 4*pi / (W0 + sqrt(2 + W0**2)) # scale to Fourier period
COI_FACTOR= FOURIER / sqrt(2)
DOFMIN = 2 ; CDELTA = 0.776 ; GAMMA = 2.32 ; DJ0 = 0.60
Morlet 6 is the standard choice because it puts the time-frequency tradeoff near the
centre and because FOURIER is then close to 1.03, so scale and period are
nearly interchangeable and the plots can be read directly in weeks. The transform is
computed by FFT with zero padding to the next power of two, as in Torrence and Compo
section 3.h, which is what makes 91 pairs across 55 scales tractable.
Run parameters:
This gives 55 scales spanning 8.3 to 890 weeks, so from roughly two months to seventeen years. The lower bound is set by the weekly sampling and the upper bound by the record length: 890 weeks is about 4.5 cycles in 76 years, which is the fewest cycles worth reporting at all.
7.2 Coherence and phase
The smoothing follows Grinsted: a Gaussian in time applied in the frequency domain as
exp(-(s*k)^2 / 2), then a boxcar in scale of width
round(DJ0 / dj) forced to be odd. Both smoothings are required. Without
them the coherence of any two series is identically 1 everywhere, since a ratio of a
quantity to itself is being formed at every point.
Phase is converted to a lead in weeks:
def phase_to_lead(phase, periods):
return phase / (2*np.pi) * periods # positive: the FIRST series leadsThe sign convention was printed backwards in an early version and every lead in this report has been recomputed under the corrected convention.
7.3 Significance
Two different null models, because two different quantities are being tested.
- Power is tested against the analytic chi-square distribution for an AR(1) red noise background, with the AR(1) coefficient estimated per series. This is the standard Torrence and Compo test and is closed form.
- Coherence has no usable analytic null after the Grinsted smoothing, so it
is tested by Monte Carlo. For each pair, AR(1) surrogates are generated
matching both series' estimated coefficients, with
N_MC = 120realisations per alpha pair on a 0.05 grid of the coefficient, seed 20260726. The 95th percentile of surrogate coherence becomes the significance level for that pair.
Two consequences must be respected when reading the results. The Monte Carlo standard error on the 95 per cent level is roughly 0.01, so any pair whose coherence sits within about 0.01 of its level is a tie and must not be reported as significant or as not significant. And the level varies by pair, because it depends on both series' AR(1) coefficients, so coherence values are not directly comparable across pairs without their levels.
The cone of influence marks where zero padding contaminates the transform. Over this
record, 86.7 per cent of the time-scale plane lies inside the cone. Leads are reported
only from cells that are inside the cone and have R2 > 0.5, and only
where more than 20 such cells exist.
7.4 Self test
The wavelet code is hand written, so it is verified against synthetic series with known answers before any real result is trusted. All fourteen checks pass:
PASS recovers a 52-week period peak at 52.48 wk
PASS separates 26 and 104 week components maxima at 26.2, 105.0, 187.0, 296.8
PASS recovers a +0 week lead at period 52 measured -0.00 wk
PASS recovers a +4 week lead at period 52 measured +4.03 wk
PASS recovers a +8 week lead at period 52 measured +8.07 wk
PASS recovers a -6 week lead at period 52 measured -6.05 wk
PASS coherence of a series with itself is 1 min 1.0000
PASS coherence of independent noise is low mean R2 0.413
PASS AR(1) surrogate reproduces its coefficient estimated 0.706
PASS AR(1) surrogate has unit variance sd 1.003
PASS wavelet energy matches series variance ratio 0.900
PASS cone of influence is zero-width at the edges
PASS cone of influence is widest mid-series max at index 285
PASS about 5 percent of white-noise power exceeds 3.0% exceeds
the 95 percent level
all wavelet_core checks passedThe lead recovery checks are the important ones, since the lead is the operationally useful output and its sign convention had already been wrong once. Note that the independent noise coherence of 0.413 is not near zero, which is precisely why the Monte Carlo null is necessary: smoothed coherence between two unrelated red noise series is intrinsically high and an unwary reader would over-interpret a value of 0.4.
7.5 Bands and results
Coherence is averaged into four bands: sub-seasonal 8 to 26 weeks, annual 40 to 70 weeks, interannual 104 to 364 weeks, decadal 416 to 890 weeks. Across the 91 pairs, 16 clear the null in the sub-seasonal band, 35 in the annual band, 22 in the interannual band and 8 in the decadal band.
Bogota's thirteen pairs, sorted by annual coherence:
| Partner | annual | level | sig | interannual | level | sig | interann lead (wk) |
|---|---|---|---|---|---|---|---|
| Guaviare | 0.440 | 0.380 | yes | 0.412 | 0.401 | yes | +9.2 |
| Bolivia | 0.439 | 0.390 | yes | 0.401 | 0.402 | no | +35.5 |
| Amazon_S | 0.413 | 0.375 | yes | 0.406 | 0.404 | yes | +29.5 |
| Amazon_arc | 0.398 | 0.385 | yes | 0.396 | 0.424 | no | +18.4 |
| Xingu | 0.389 | 0.385 | yes | 0.451 | 0.424 | yes | +37.5 |
| Cerrado | 0.386 | 0.375 | yes | 0.416 | 0.404 | yes | +38.8 |
| SE_Australia | 0.376 | 0.375 | yes | 0.553 | 0.404 | yes | -21.5 |
| US_West | 0.371 | 0.375 | no | 0.372 | 0.404 | no | -63.4 |
| Canada_boreal | 0.339 | 0.383 | no | 0.376 | 0.421 | no | +2.4 |
| Siberia | 0.332 | 0.383 | no | 0.381 | 0.421 | no | +17.7 |
| Indonesia | 0.329 | 0.390 | no | 0.572 | 0.402 | yes | -14.8 |
| Central_Africa | 0.322 | 0.390 | no | 0.327 | 0.421 | no | +35.3 |
| Med_Europe | 0.317 | 0.389 | no | 0.387 | 0.418 | no | -14.2 |
Three readings. First, the six southern Amazon sites occupy six of the seven top positions in the annual band and all six are significant, which is the coherence evidence for the direct coupling tested in section 8. Second, distance decays the relationship in an orderly way: Guaviare at 479 km gives 0.440 while Xingu at 2973 km gives 0.389. Third, the two strongest interannual coherences with Bogota are not Amazonian at all: Indonesia at 0.572 and SE Australia at 0.553, both with negative leads, meaning Bogota lags them. Both are classic ENSO centres of action, and this is the coherence signature of a shared ocean driver rather than of any direct Amazon-to-Bogota transport at long timescales.
Two further diagnostics. The share of annual band mean variance attributable to the ENSO mode is 48.6 per cent at Indonesia, 31.5 per cent at Bogota, 24.2 per cent at SE Australia, and 2.6 per cent or less at every one of the six Amazon sites. And the internal Amazon coupling is far stronger than anything involving Bogota: Amazon_arc against Amazon_S gives 0.919 annual and 0.892 interannual. The southern Amazon is one coherent system that is largely not ENSO driven; Bogota is a partly ENSO driven system that is weakly attached to it.
Is the Amazon to Bogota link real, and can the north inform the south
8.1 Contemporaneous coupling with three guards
The wavelet result says Bogota and the southern Amazon share annual band variance. It does not say the link survives the obvious confounders. Three guards are applied to every pairwise Pearson correlation, all on season aligned data:
- Raw, with a Bartlett and Quenouille effective sample size correction
neff = n * (1 - a*b) / (1 + a*b)using both series' lag one autocorrelations, because weekly and seasonal series are autocorrelated and the nominal n overstates the degrees of freedom. - ENSO partialled, regressing the June to August ONI out of both sides before correlating. This is the guard that matters most, since section 7 shows a shared ocean mode exists.
- Linearly detrended, removing a linear trend from both sides, which prevents two independently drying series from producing a spurious correlation.
The season alignment regresses the Amazon September to November index of year Y against the Bogota December to March season opening in December of year Y, so the Amazon precedes Bogota by construction. A reverse direction placebo runs the same test with the time order flipped.
| Partner | r raw | p | r ENSO-partialled | p | r detrended | p | survives all three |
|---|---|---|---|---|---|---|---|
| Guaviare | 0.317 | 5.3e-3 | 0.288 | 0.0121 | 0.311 | 6.2e-3 | yes |
| Amazon_arc | 0.298 | 8.9e-3 | 0.173 | 0.138 | 0.294 | 1.0e-2 | no |
| Amazon_S | 0.281 | 1.4e-2 | 0.213 | 0.0668 | 0.274 | 1.7e-2 | no |
| Xingu | 0.166 | 0.151 | 0.105 | 0.368 | 0.169 | 0.144 | no |
| Cerrado | 0.209 | 0.0694 | 0.134 | 0.253 | 0.209 | 0.0697 | no |
| Bolivia | 0.254 | 2.7e-2 | 0.307 | 7.3e-3 | 0.244 | 3.4e-2 | yes |
| Amazon six site mean | 0.293 | 1.0e-2 | 0.233 | 0.0444 | 0.287 | 1.2e-2 | yes |
| ONI Jun-Aug | -0.671 | 3.3e-11 | — | — | -0.675 | 2.2e-11 | reference |
For context, the June to August ONI alone correlates with the Bogota December to March season at r = -0.671 over 76 seasons, and that correlation barely changes on detrending. The ocean is by a wide margin the strongest single predictor of a Bogota fire season, and no Amazon site comes close.
Lagged cross correlations were also computed over LAGS = arange(-52, 53)
weeks, with positive lag meaning the other site leads Bogota, and are archived in
data/coupling_lags.csv. Unlike the wavelet script, this analysis uses no
surrogates or bootstrap; its inference rests on the effective sample size correction
and on the three guards.
8.2 Northern hemisphere teleconnection
The question asked was whether Europe, the United States and Canada carry information about the southern fire seasons. The answer is yes for the Amazon, no for Bogota by the composite route, but yes for Bogota via Canada specifically.
| Link | raw | ENSO-partialled | detrended | adj R² before | adj R² after |
|---|---|---|---|---|---|
| Canada boreal → Bogota | -0.317 (p 0.0054) | -0.285 (p 0.0253) | -0.319 (p 0.0051) | +0.443 | +0.483 |
| Europe + US composite → Amazon | +0.333 | +0.372 | +0.335 | +0.026 | +0.149 |
| MTBS hectares → Amazon | significant | — | +0.081 (p 0.609) | — | rejected |
The Canada link is the one that reaches the operational models: Canada_boreal_jj
is selected in both the 1 August and 1 September Bogota mean models, with a negative
coefficient, so a wet Canadian boreal June and July forecasts a dry Bogota December to
March. It survives all three guards and lifts the adjusted R squared from +0.443 to
+0.483. For the Amazon the Europe plus US composite is a larger relative gain, taking
adjusted R squared from a near useless +0.026 to +0.149, though from a very low base.
The critical check on all of this is that none of the northern sites correlates significantly with the June to August ONI. They are therefore not the ocean signal in disguise, and the information they add is genuinely additional rather than a second copy of a predictor already in the model.
MTBS burned hectares were tested as a northern predictor and rejected. The raw correlation is significant, but it collapses to +0.081 at p = 0.609 on detrending. MTBS hectares have a strong upward trend, and so do several southern moisture series, so the raw correlation was trend against trend. This is exactly the failure mode the detrending guard exists to catch, and it is reported here rather than quietly dropped.
8.3 Moisture to fire lead time
The lag at which the weekly moisture index best predicts weekly fire is a direct measure of how moisture limited a site is.
| Site | best lag (wk) | r | p |
|---|---|---|---|
| Bogota_Cerros | 1 | -0.453 | 1.1e-15 |
| Guaviare | 2 | -0.114 | 0.0382 |
| Amazon_arc | 2 | -0.238 | 9.0e-4 |
| Amazon_S | 3 | -0.186 | 9.8e-3 |
| Xingu | 5 | -0.103 | 0.116 |
| Cerrado | 5 | -0.092 | 0.118 |
| Bolivia | 6 | -0.087 | 0.164 |
Bogota responds within one week and does so very strongly. The Amazon sites respond more slowly and more weakly on detections, with the southernmost and most seasonal sites showing the longest lags and the weakest significance. At the seasonal scale the same structure appears far more strongly in the burned area domain, where moisture leads fire by about five weeks in the moisture limited south, r = -0.82 at Bolivia and r = -0.85 for the southern basin composite, and collapses to zero lead at Guaviare, Xingu and Bogota. The apparent contradiction is a matter of which response is being measured: Bogota's fire counts react within a week, while the southern Amazon's burned area accumulates over more than a month after the moisture signal.
From a seasonal probability to a dated probability
A statement that there is a 68 per cent chance of a critical Bogota season is not yet actionable, because it carries no date. The crossing model converts the seasonal probability into a curve giving the probability that a sustained crossing has already happened by any given week.
9.1 The sustained crossing rule
A crossing is counted only when the three week centred smoothed index sits at or below a level for two consecutive weeks. Because the index is already smoothed over three weeks, this is effectively a four week dry spell. The rule replaced a naive first touch rule that generated artifact crossings from single dry weeks, in some cases months before any real drying. The three levels evaluated are -0.5, -1.0 and -1.5, with -1.0 as the operational action line.
9.2 Crossing climatology over 76 cycles
| Basin | level | critical seasons | median date | normal seasons |
|---|---|---|---|---|
| Bogota | -0.5 | 100 % | 9 Nov | 61 %, median 14 Dec |
| Bogota | -1.0 | 88 % | 7 Dec (IQR 9 Nov to 1 Feb) | 14 % |
| Bogota | -1.5 | 40 % | 11 Jan | 4 % |
| Amazon | -0.5 | 100 % | 8 Jun | 31 % |
| Amazon | -1.0 | 84 % | 6 Jul (IQR 25 May to 3 Aug) | 6 % |
| Amazon | -1.5 | 32 % | — | 2 % |
The -1.0 level is chosen as the action line because it discriminates most sharply: 88 per cent against 14 per cent for Bogota and 84 per cent against 6 per cent for the Amazon. The -0.5 level is crossed in every critical season but also in most normal ones, so it warns without discriminating. The -1.5 level discriminates cleanly in relative terms but fires in fewer than half of critical seasons, so it misses too much to be the primary trigger. It is retained as a severity escalation.
9.3 The blend
The dated curve is a two component mixture. Each component is an empirical cumulative distribution of crossing dates over the 76 cycle climatology, one conditional on the season being critical and one conditional on it not being, and the mixture weight is the calibrated seasonal probability from the 1 August forecast.
With P(critical) = 0.6769 for Bogota and 0.6293 for the Amazon. The
construction is deliberately transparent: no distributional assumption is made about
the crossing dates, the conditional curves are simply observed frequencies, and the
only model input is the single mixture weight. If the forecast is reissued on 1
September with a different probability, only the weight changes and the curve rescales.
Ordinal week mapping across the year boundary, which was an early source of an off by one error and is therefore written out explicitly:
(cw - 40) if cw >= 40 else (cw + 12) # forward, Bogota window
def ord_date(basin, o):
o = int(round(o))
if basin == 'Bogota':
return wk_date(2026, 40 + o) if o <= 12 else wk_date(2027, o - 12)
return wk_date(2026, 18 + o)9.4 Results for the 2026 to 2027 cycle
| Milestone | Amazon (Sep-Nov 2026) | Bogota (Dec 2026 - Mar 2027) |
|---|---|---|
| 25 per cent probability reached | 6 Jul 2026 (passed) | 23 Nov 2026 |
| 50 per cent probability reached | 7 Sep 2026 | 1 Feb 2027 |
| End of season probability | 55 % by late Nov | 61 to 64 % by 31 Mar |
| Both basins critical in the same cycle | 10 of 76 cycles, 13 per cent, with the later arrival by 30 November | |
9.5 Analog seasons
As an independent check on the blend, the ten nearest analog cycles in predictor space were identified for each basin. For Bogota: 2001, 2007, 1952, 1990, 1977, 1964, 2019, 2003, 1958 and 1973, of which 7 of 10 ended critical and 8 of 10 produced a sustained crossing, with a median crossing date of 23 November and an interquartile range of 2 November to 14 December. For the Amazon: 1968, 2008, 2000, 2023, 1957, 2006, 2018, 2011, 2017 and 2002, of which 4 of 10 ended critical.
The Bogota analogs at 7 of 10 sit slightly above the model's 0.677, and the Amazon analogs at 4 of 10 sit well below the model's 0.629. The Amazon disagreement is the larger and is not resolved here; it is reported as a source of uncertainty and is one reason the 1 September reissue matters more for the Amazon than for Bogota. An earlier median analog path calculation was found to be misleading and was removed rather than patched.
9.6 Youden J triggers for the watch windows
Each forecast reissue window has one variable worth watching. The threshold for each is set at the Youden J optimum on the 76 cycle record, reported with hit rate against false alarm rate and with the AUC of the underlying variable.
| Variable | Trigger | hits / false alarms (%) | AUC |
|---|---|---|---|
| bog_son | at or below +0.02 | 88 / 25 | 0.86 |
| amz6_aug | at or below -0.50 | 73 / 12 | 0.85 |
| gua_aug | at or below -0.32 | 77 / 16 | 0.84 |
| oni_ason | at or above +0.09 | 80 / 24 | 0.82 |
| oni_jja | at or above -0.18 | 92 / 41 | 0.78 |
| oni_mjj | at or above -0.23 | 92 / 51 | 0.71 |
| oni_tend | at or above +0.29 | 64 / 31 | 0.68 |
| gua_son | at or below -0.18 | 60 / 27 | 0.64 |
The oni_ason threshold was originally computed in the wrong direction,
producing an AUC of 0.18 and zero flagged seasons, which is the signature of a reversed
inequality. Recomputed as "at or above", it gives the +0.09 trigger with 80 per cent
hits at 24 per cent false alarms and AUC 0.82. The gua_son trigger at AUC
0.64 is weak and is listed for completeness rather than for use.
What has actually changed since 1950
10.1 Moisture trends
Theil-Sen slopes with Kendall tau significance, in index units per decade, over 77 years.
| Site | index per decade | tau | p | direction |
|---|---|---|---|---|
| Central_Africa | -0.278 | -0.470 | 1.4e-9 | drying, strong |
| Amazon_S | -0.096 | -0.267 | 5.8e-4 | drying |
| Cerrado | -0.093 | -0.241 | 2.0e-3 | drying |
| Bolivia | -0.090 | -0.252 | 1.2e-3 | drying |
| Xingu | -0.071 | -0.179 | 2.1e-2 | drying |
| Med_Europe | -0.062 | -0.197 | 1.1e-2 | drying |
| US_West | -0.055 | -0.183 | 1.9e-2 | drying |
| Amazon_arc | -0.054 | -0.143 | 0.066 | drying, not significant |
| Guaviare | -0.022 | -0.079 | 0.308 | flat |
| Bogota_Cerros | -0.011 | -0.027 | 0.732 | flat |
| Siberia | -0.028 | -0.075 | 0.338 | flat |
| Indonesia | +0.025 | +0.089 | 0.253 | flat |
| SE_Australia | +0.026 | +0.071 | 0.360 | flat |
| Canada_boreal | +0.072 | +0.210 | 6.9e-3 | wetting |
10.2 Dry week counts
The trend in the mean understates what has happened in the tail. Weeks per year with the index at or below -1.0, comparing 1980 to 1999 against 2015 to 2025:
| Site | 1980-1999 | 2015-2025 | change |
|---|---|---|---|
| Amazon_S | 2.0 | 20.0 | 10.0× |
| Amazon_arc | 3.1 | 19.1 | 6.2× |
| Xingu | 3.2 | 15.9 | 5.0× |
| Cerrado | 3.1 | 15.7 | 5.1× |
| Bolivia | 3.5 | 15.9 | 4.5× |
| Guaviare | 6.7 | 15.4 | 2.3× |
| Bogota_Cerros | 6.4 | 7.0 | 1.1× |
Severe weeks, index below -1.5, rose from 0.7 to 10.9 per year at Amazon_S alone. A tenfold increase in dry weeks against a moderate shift in the mean is the signature of a distribution whose left tail has extended rather than whose centre has moved, which matters because it is the tail that the fire calibration in section 5 says is dangerous. Bogota again shows almost no change, 6.4 to 7.0 weeks.
10.3 Cerros burned area, and the availability control
The Landsat derived Cerros burned area series runs from 1987 with 31 usable years. Decadal medians are 13901, 4764 and 11853 hectares, and the two largest years are 2024 at 48303 hectares and 1987 at 32955 hectares.
| Series | Theil-Sen slope (ha/yr) | 95 % CI | tau | p | n |
|---|---|---|---|---|---|
| Total burned area, raw | +142.65 | -83.2 to +434.9 | +0.166 | 0.198 | 31 |
| Largest season, raw | +5.43 | -225.7 to +220.0 | +0.011 | 0.947 | 31 |
Splitting into epochs shows no consistent direction either:
| Epoch | total slope | tau | p | largest season slope | p | year vs scenes r |
|---|---|---|---|---|---|---|
| 1987-1999 (n 6) | -1230.6 | -0.600 | 0.136 | -1469.4 | 0.017 | +0.455 |
| 2000-2012 (n 12) | +105.9 | +0.152 | 0.545 | -113.0 | 0.737 | -0.035 |
| 2013-2025 (n 13) | +14.9 | +0.026 | 0.952 | -144.1 | 0.435 | +0.324 |
For contrast, the MTBS United States west control series does show a clear and robust increase of +24733 hectares per year at p = 0.0015. That series has a consistent observing protocol across its whole length. The difference between the two cases is entirely about data provenance, and it is the reason the Cerros trend claim is withdrawn while the US one is not.
What did not work, and what that taught us
Six results in this project are negative or partly negative. They are collected here rather than scattered, because in an early warning system the failures constrain the operational envelope as much as the successes do.
1. Amazon detection counts do not respond to the seasonal index
Spearman -0.348 at p = 0.088 with a group ratio of only 1.22, against -0.814 at p = 1.3e-05 and a ratio of 2.06 for burned area on the same years. The interpretation in section 6.2 is that Amazon ignition is anthropogenic and roughly constant while moisture governs spread. The consequence is that the Amazon must be scored on burned area, which in turn ties the Amazon product to FireCCI51 and its 2001 to 2020 window, a real operational constraint.
2. The Cerros burned area increase does not survive an availability control
Detailed in section 10.3. The raw slope of +142.7 hectares per year is not significant, and the year coefficient turns significantly negative once Landsat scene counts enter the model. This is a straightforward observing system confound and the claim is withdrawn.
3. MTBS hectares are a trend against a trend
Significant raw, +0.081 at p = 0.609 detrended. Rejected as an Amazon predictor. The same detrending guard is what makes the surviving Canada boreal link credible, since it passes the identical test.
4. Four of six Amazon sites lose their Bogota link under an ENSO partial
Amazon_arc, Amazon_S, Xingu and Cerrado all fall below significance once the June to August ONI is regressed out. Only Guaviare, Bolivia and the six site mean survive all three guards. The Amazon to Bogota teleconnection is real at roughly r = 0.23 to 0.29 after the ocean is removed, which is a genuine but modest link and should not be described as strong.
5. Mosaic buffering does not appear at coarse resolution
The hypothesis that a more heterogeneous vegetation mosaic buffers fire spread returned a null at MODIS resolution. Investigation showed this to be a resolution limit rather than a real null: against the Bogota ground truth from fire brigade records, MODIS detected only 27 per cent of real fire-weeks. At 10 m resolution using AlphaEarth embeddings the relationship recovers cleanly: larger fires burn less heterogeneous vegetation, Spearman rho = -0.34 at p = 9e-10, tightening to rho = -0.39 on a clean pre-fire subset, with mean heterogeneity 0.016 in large fires against 0.022 in small ones at p = 2.5e-6. The lesson is that a coarse-resolution null on a landscape-structure question is uninformative, not negative.
6. Continuous skill at 1 August is low even where discrimination is good
Skill of 0.084 to 0.127 across the four live models. The system does not claim to forecast how dry the season will be four months ahead. It claims, and demonstrates at AUC 0.69 to 0.80, to forecast whether the season lands in the critical third. Products built on this system must be probability of exceedance products, not index value products, and the 1 September and 1 December reissues are where the point forecast becomes worth quoting.
Two anomalies worth recording
The 2024 Amazon season produced a burned area anomaly of +13 sigma including a single 2.78 million hectare megafire at 48 sigma, with the Indigenous Territory to frontier burned area ratio moving from 0.39 to 0.80. A moisture-driven system will not anticipate an event of that character from moisture alone. Separately, the Bogota 2024 season was called critical by the rule and recorded the largest Cerros burned area on record at 48303 hectares while returning only 19 MODIS detections, which is the detection sensitivity problem of item 5 appearing in the validation set.
What this system cannot do
- It forecasts moisture, not fire. The link to fire is a calibration, and calibration slopes are estimated over a period during which land use also changed. A change in ignition pressure would move the fire outcome without moving the forecast.
- Reanalysis is not observation. ERA5-Land soil moisture is a model product. Its skill in the Cerros, which are steep, small and heterogeneous relative to the grid, is not independently validated here. The most recent week is preliminary and subject to revision.
- Seventy six cycles is a small sample for a one-in-three event. Twenty five to twenty seven events support every AUC and Brier figure quoted. Confidence intervals on those figures are wide and are not reported, which is itself a limitation.
- The nested cross validation removes selection bias, not specification bias. The candidate pools, the window definitions, the three predictor cap and the OLS form were all chosen by the analyst with the full record visible.
- The crossing climatology is stationary by construction. Conditional crossing-date distributions are pooled over 76 years. Section 10 shows the southern Amazon has changed materially over that period, so the Amazon conditional curves probably run late.
- Bogota's fire response is measured against a detection product that misses most of its fires. The 27 per cent detection rate against brigade ground truth means the Bogota Poisson calibration is fitted to a biased sample of events, most likely biased toward larger ones.
- The Amazon analogs and the model disagree. Four of ten analogs critical against a model probability of 0.629. That gap is unresolved.
- Nothing here is a fire danger rating. Wind, fuel load, fuel continuity, slope, ignition source and suppression capacity are all absent. This is one input to a preparedness decision, not the decision.
Revision triggers
The forecast should be reissued, and this report revised, when any of the following occur: the 1 September Amazon reissue, which raises skill from 0.084 to 0.358 and AUC from 0.691 to 0.906; the 1 December Bogota reissue, which raises skill to 0.352 and AUC to 0.885; the August Amazon composite crossing the -0.50 Youden trigger; the ASON ocean index crossing +0.09; or the Bogota September to November index crossing +0.02.
Code, caches and outputs
| Script | Produces | Log |
|---|---|---|
ews_trends.py | the normal score index, trend statistics | ews_trends.log |
forecast_model.py | the eight seasonal models, live forecast | forecast_model.log |
fire_calibration.py | Poisson and burned area calibration | fire_calibration.log |
target_validation.py | the four validation tests per basin | target_validation.log |
wavelet_core.py, pairwise_wavelet.py | coherence over 91 pairs | pairwise_wavelet.log |
wavelet_selftest.py | the fourteen synthetic checks | wavelet_selftest.log |
coupling_analysis.py | guards, lags, placebo | coupling_analysis.log |
timing_analysis.py | crossing climatology, blended curve | — |
make_fig_timing.py | fig_timing.png, fig_timing.pdf | — |
Cross validation results are cached at data/cv_cache/<md5>.pkl, keyed
on repr((basin, issue, ycol, pool, MAXP)) together with a pandas hash of
the input data rounded to nine decimal places, so a cache entry cannot survive a change
in either the specification or the data. The wavelet Monte Carlo uses seed 20260726 and
is therefore exactly reproducible.
Principal data outputs: data/weekly_sm_z_allsites.csv,
data/forecast_live.csv, data/forecast_predictors_bogota.csv,
data/forecast_predictors_amazon.csv,
data/fire_calibration.csv, data/fire_calibration_area.csv,
data/coupling_guards.csv, data/coupling_lags.csv,
data/wavelet_pairs.csv, data/timing_curve.csv,
data/timing_summary.json, data/ews_trends_summary.json.
Companion documents: SITREP.html carries the operational statement, the crossing dates, the watch windows and the recommended actions. fire_risk_dashboard.html carries the lay-audience framing and the current-state charts. This report carries the methods, the full statistics and the negative results.
Note: Only a member of this blog may post a comment.