Changepoint Methods in Climatology

Perhaps the chief use of changepoint methods in climatology is to homogenize series. Climate homogenizers seek to adjust series for any man-made effects by subtracting/adding estimates of the $\Delta_{i}$’s to the series. This rids the series of the man-made effects induced by gauge changes, station relocations, or other measurement change methods, leaving a record that can be attributed only to natural forcings. Climatologists prefer to leave changepoint effects from natural causes in the series, deeming it part of the ‘true record’. Adjusting a series for an excessive number of changepoints risks removing too much variation in the record, making other fluctuations seem more significant than they truly are.

To homogenize series, climatologists often compare a series for the location under study, called the target series, to a series of the same quantity collected at a nearby recording station that experiences similar weather, called a reference series. Often, many reference series are used to make conclusions, sometimes as many as 40 in practice. If $\{X_{t}\}$ denotes the target series and $\{Y_{t}\}$ the reference series, climatologists often analyze the difference $\{X_{t}-Y_{t}\}$ for changepoints. If the reference series matches the target series well, it should have similar fluctuations (weather) as the target series. The subtraction should hence eliminate most naturally occurring fluctuations. However, if the target series has a shift induced my a man-made changepoint (say a station relocation or a gauge change), then it will also be a changepoint in the target minus reference series.

Unfortunately, by making target minus reference comparisons, changepoints in either the target or reference series become changepoints in the target minus reference series. This in essence doubles the number of changepoints in the series, complicating matters for the homogenizer. Modern homogeneity methods combat this issue by examining many target minus reference series: if the target minus reference series flags the same changepoint time for ten distinct station references, then evidence suggests that the changepoint time is attributable to the target series. There exist stations where multiple good reference series do not exist; the North Atlantic Basin tropical cyclone analysis in the next subsection has no references series whatever.

Another use of changepoint methods lies with stationarity assessments. A common null hypothesis for climatologists is that some aspect of climate is static (stationary). Changepoint methods can be used to assess this hypothesis: if a multiple changepoint procedure flags one or more changepoints, then the record is deemed non-stationary. The form of the mean $\mu_{t}=E[X_{t}]$ is not important here; linear trends, sinusoidal components, and even hinge shaped structures have arisen in practice. In general, a multiple mean shift changepoint procedure will flag changepoints that attempt to follow the mean of the series. For example, if the series has an increasing linear trend, then a mean shift changepoint routine will attempt to approximate this by flagging multiple changepoints that result in an increasing staircase.

As an example of a stationarity check, we consider the North Atlantic Ocean Basin’s tropical cyclone counts for the 53-year period 1970-2022. These are plotted in Figure 2 against a mean shift changepoint fit that is described below.

Annual tropical cyclone counts in the North Atlantic Basin are well modeled by a Poisson marginal distribution. Hence, a reasonable model posits that the cyclone count for year $t$, denoted by $X_{t}$, has a Poisson distribution with parameter $\lambda_{t}$. Since $\mu_{t}=\lambda_{t}$ in this case, we conduct a multiple changepoint analysis with

The year-to-year storm counts are taken to be statistically independent; indeed, there is very little predictive power in a tropical cyclone count forecast even 12 months in advance (annual conferences are held each year in May, just prior to the start of hurricane season, where various modelers release their annual count forecast).

The multiple changepoint fit estimates one changepoint in 1995. The estimated mean of the series is plotted in Figure 2 against the storm counts. Clearly, we have entered into an era of enhanced tropical cyclone activity circa the mid 1990s, a concerning aspect for Eastern Coastal residents.

Binary segmentation methods essentially turn an AMOC procedure into a multiple changepoint technique. The process is simple: one examines the entire series for the existence of a single changepoint. If none is found, the procedure is stopped and no chagepoints are declared and the entire series is declared changepoint free. Should a changepoint be flagged, the series is split about the estimated changepoint time and the two series subsegments (before and after the changepoint time) are scrutinized for additional single changepoints. The procedure is repeated iteratively, stopping when all series subsegments are declared changepoint free. Binary segmentation methods can euse any AMOC test — it need not be a CUSUM procedure.

Binary segmentation calculations proceed very rapidly in comparison to some other multiple changepoint techniques. In fact, the calculations are almost instantaneous, even for long series. Unfortunately, binary segmentation is one of the worst performing multiple changepoint procedures, often being fooled when multiple mean shifts move the series in different directions. Simulations reinforce this point.

Attempted fixes to binary segmentation include wild binary segmentation (WBS) and its related procedures. WBS examines randomly selected subintervals for a single changepoint with a CUSUM AMOC test and then combines the conclusions from all random subintervals to estimate the multiple changepoint configuration. With its current settings, WBS is an aggressive procedure that often overestimates the number of changepoints, especially when few changepoints in truth exist. Research on these methods continues today.

Another way to estimate a multiple changepoint configuration uses likelihoods that are penalized to prevent fitting an excessive number of changepoints. Changepoint configurations were estimated in this manner earlier for the Tuscaloosa temperature data and the North Atlantic tropical cyclone counts.

Suppose that $m$ changepoints occur at the times $\tau_{1},\ldots,\tau_{m}$. The statistical likelihood of such a model will be denoted by $L(m;\tau_{1},\ldots,\tau_{m})$. For the Tuscaloosa data, this likelihood was based on a Gaussian time series model with first order autoregressive errors; for the North Atlantic tropical cyclone counts, the likelihood is based on independent Poisson counts.

As more changepoints are added to the model, the resulting likelihood becomes bigger, or equivalently, $-2\ln(L(m;\tau_{1},\ldots,\tau_{m}))$ becomes smaller. However, once the optimal model is reached, adding additional changepoints to the fit does little to improve its likelihood. A penalized likelihood procedure adds a penalty to the negative log likelihood for having $m$ changepoints at the times $\tau_{1},\ldots,\tau_{m}$. We denote the penalty by $P(m;\tau_{1},\ldots,\tau_{m})$. The penalized likelihood estimated changepoint configuration is the model that minimizes

Many different penalty types have been proposed. One that seems to work well in mean shift problems is the BIC penalty $P(m;\tau_{1},\ldots,\tau_{m})=m\ln(N)$. The classical AIC penalty $P(m;\tau_{1},\ldots,\tau_{m})=2m$ tends to estimate too many changepoints. There exist more exotic penalties that depend on the times of the changepoints (one such is a minimum description length).

While penalized likelihood methods perform better than binary segmentation in simulations, they are much more computationally intensive. In particular, the objective function above is highly non-convex, making locating its minimum unwieldy. An exhaustive search over all changepoint configurations requires calculations to fit $2^{N-1}$ distinct models (this tally is not $2^{N}$ since time 1 cannot be a changepoint). This exceeds one billion model fits when $N\geq 32$. To combat this issue, researchers have devised genetic algorithms to optimize the penalized likelihood. Genetic algorithms are intelligent random walk searches that are unlikely to fit models at suboptimal changepoint configurations. These techniques rely on random generations of changepoint configurations that evolve according to elements of genetic fitness and mutation. While many "off the shelf" genetic algorithms can now reliably handle the problem, they often take hours of computing time.