Resplandy et al. Part 2: Regression in the presence of trend and scale systematic errors

Originally a guest post on

In a recent article I set out why I thought that the trend in ΔAPOClimate was overstated, and its uncertainty greatly understated, in the Resplandy et al. ocean heat uptake study. In this article I expand on the brief explanation of the points made about “trend errors” and “scale systematic errors” given in my original article, as these are key issues involved in estimating the trend in ΔAPOClimate and its uncertainty.

I will illustrate the trend error point using a 26-year time series that is a slightly modified version (‘pseudo-ΔAPOAtmD‘) of Resplandy et al.’s ΔAPOAtmD time series and associated uncertainties. The pseudo-ΔAPOAtmD best-estimate time series increases linearly by 0.27 per meg each year, starting at zero in 1991. In each year, the 1-σ uncertainty (error standard deviation) in its value is 50% of its best-estimate value (with σ set to 10-6 in 1991 to avoid divide-by-zero errors). I will assume errors have a Normal distribution, with zero mean. Figure 1 illustrates this situation.

Figure 1. Estimated pseudo-ΔAPOAtmD best values and their uncertainties. The black line goes through all years’ best estimates, marked with small black crosses. The pink area represents the ± 1-σ uncertainty limits in continuous time, and the red bars show that uncertainty for each year’s data point.

The case when errors are uncorrelated

If errors are uncorrelated then the usual conditions for ordinary least squares (OLS) to give valid estimates are satisfied. Obviously, regressing the best fit values will give a perfect fit with, correctly, a trend of 0.27 per meg yr−1, since they increase linearly with time. The trend and the uncertainty in it can be estimated as follows:

  1. a) for (say) 10,000 sample realizations of the pseudo-ΔAPOAtmD time series, add random, independent, errors drawn from each year’s error distribution to the pseudo-ΔAPOAtmD best estimate time series values;
  2. b) find the trend for each sample realization using OLS regression;
  3. c) take the mean of the 10,000 trends as the trend (best-)estimate and their standard deviation as the 1-σ uncertainty of that estimate.

On doing so, I obtained a trend estimate of 0.271 with uncertainty of ± 0.056. (Note: with only 10,000 samples the figure in the third decimal place is unreliable.).

I then repeated the calculation, with the same set of cases (sample realizations) but using weighted linear regression (WLS). The weight for each year were set at 1/σ2, where σ is the error standard deviation for that year, as is usual.

The resulting trend estimate was 0.270, with uncertainty of ± 0.033.

In both the OLS and WLS cases, the mean of the trend error estimates from each of the 10,000 regressions was close to the standard deviation of the 10,000 trend estimate, as one would expect.

So, if errors are uncorrelated between years but their estimated uncertainty varies between years, both OLS and WLS on average provide unbiased trend estimates, but with WLS the variation in the trend estimate from its true value between possible sample realisations is smaller. WLS in effect uses the available information more efficiently and produces more precise estimates.

The case when errors are perfectly correlated (trend errors)

Now consider the case where the uncertainty arises from the trend in the data being uncertain. For any given realisation of the pseudo-ΔAPOAtmD time series, there is then only a single uncertain error. That error scales up with time, pro rata with the magnitude of the time difference from the baseline year. The error distributions for each year are the same as previously, but they are now perfectly correlated across all years rather than being independent. For any given realisation, all data points will therefore lie on a straight line – a version of the black line rotated about its zero starting point in 1991, and most likely lying within the area shaded pink.

Repeating steps a) to c), but this time with each year’s error perfectly correlated with other years’ errors, I obtained the following results.

With OLS regression, a trend estimate of 0.268 with uncertainty of ± 0.135 per meg yr−1.

With WLS regression, a trend estimate of 0.268 with uncertainty of ± 0.135 per meg yr−1.

In both cases the trend error estimate for each of the 10,000 realisations was zero, as although the fit differed between realisations, it was a perfect fit in every case. And in both cases the trend uncertainty was the same as that in the original data – 50% of 0.27 per meg yr−1.

The key point is that if the data has a trend error, regression cannot reduce it at all. Weighting data values to reflect their absolute precision does not help, but it does no harm either in this case.

Which case applies to the actual ΔAPOAtmD values?

Resplandy et al. state the results of their model simulations of the fertilisation effects of combined, N, Fe and P aerosol deposition as: “The overall impact on ΔAPOAtmD is +0.27 per meg yr−1 over 27 years of simulation (1980–2007), which we extrapolate to our 1991–2016 period”, and that “Uncertainties at the 1σ level on ΔAPOAtmD are assumed to be ±50%”. That corresponds exactly to the way my pseudo-ΔAPOAtmD best-estimate time series and uncertainties were calculated. It is quite clear that the ΔAPOAtmD best-estimate time series should represent an exact linear trend of +0.27 per meg yr−1 with its errors being perfectly correlated, non-independent trend errors.

Scale systematic error case

Consider now the case where the best-estimate data, while broadly trending, has not been derived from a linear trend. That applies to the ‘scale systematic’ error component of ΔAPOOBS.

ΔAPOOBS is calculated from measured changes in the atmospheric (δO2/N2) ratio and CO2 concentration (XCO2, in p.p.m.) as:

ΔAPOOBS = (δO2/N2) + 1.1 ⤬ (XCO2 − 350) / XO2

where XO2 (= 0.2094) is the reference atmospheric mole fraction needed to convert XCO2 from p.p.m. to per meg units.

Resplandy et al. Extended Data Table 3 states that, in addition to Corrosion, Leakage and Desorption errors of respectively ± 0.3, ± 0.2, and ± 0.1 per meg yr−1, and random errors in each year of ± 2 per meg (± 4 before July 1992), there is a scale systematic error of 2% on (δO2/N2).

I downloaded the (δO2/N2) data for the three monitoring stations used, combined it using the weights given in the paper that Resplandy et al. cited, and deducted the 1991 value so as to match Resplandy et al.’s baseline.

The OLS linear trend of the resulting data was −19.87 per meg yr−1. Using WLS regression, it was −16.05 per meg yr−1.

I then drew samples from a normal  distribution with unit mean and a standard deviation of 0.02 and, for each sample multiplied all years’ (δO2/N2) data values by the sample value.

When using OLS regression, the trend estimate was −19.87 with uncertainty of ± 0.396 per meg yr−1

When using WLS regression, the trend estimate was −16.05 with uncertainty of ± 0.320 per meg yr−1

So in both cases sampling gives the same trend estimate as regression on the original data points, but WLS estimates a substantially shallower (less negative) trend than OLS. Note that, as with a trend error, the regression trend estimate uncertainty is the same as for the original data, here 2% (of the trend estimate).

Figure 2, which shows the data values (black circles, joined by thin black line), and linear regression fits using OLS (cyan line) and WLS (red line), helps provide insight into the difference between the OLS and WLS regression fits.

Figure 2. Average (δO2/N2) annual mean data for from La Jolla, Alert, and Cape Grim monitoring stations weighted as to 0.25, 0.25 and 0.5 respectively (black circles, joined by thin black line),the OLS linear fit to them (cyan line), and the WLS linear fit to them (red line).

The reason why the WLS fit has a shallower line is that earlier years are weighted much more highly than the later years when determining the WLS fit, as the scale systematic error is applied to larger (δO2/N2) values in the later years. Indeed, the near-infinite weighting given to the 1991 data point forces the WLS fit through zero in 1991.

It can be seen that there is little evidence of random errors being significant in any year, nor greater in the later years; the black line is fairly smooth throughout the period. But the trend appears to become more negative over time, so a linear fit is strongly affected by the weightings given in different years.  Confirming that data errors are trivial, a quadratic fit is extremely close (adjusted R2 0.9998, versus 0.9950 for a linear fit). Moreover, a quadratic fit shows no evidence of fit errors being greater in later years – they are slightly higher in the first half than in the second half of the period.

There is no justification for using WLS in a case like this. Because the errors are perfectly correlated between years, the conditions for WLS to be valid are not satisfied. That the WLS result cannot possibly be valid can be seen as follows. If the method is valid, it should give the same trend whatever year is used as a baseline for the data. So, rebase all the data to a zero baseline in 2016, which implies small errors in later years and large errors in early years. That will result in the WLS weights being high in later years and low in early years. The WLS fit will then pass through the 2016 data point, and its slope will be close to that exhibited by later years’ data and, as a result, much steeper – indeed, steeper than the slope of the OLS fit (which is unchanged by the rebasing).

A known method for validly regressing when the data errors are perfectly correlated or nearly so is to transform the data by taking first differences, and then regress. The regression intercept term then corresponds to a linear trend in the original data and the slope coefficient to a quadratic trend term. When one does that using OLS, the mean trend over the full period is very close to the trend from OLS linear regression on the original data. And WLS regression on the transformed data gives a mean trend whose magnitude is only 2.5% below that for OLS regression on the transformed data, as compared with 19.2% below the OLS trend when the original data is regressed.

What this implies for trend estimation in Resplandy et al.

I have demonstrated two important points. One is that regression cannot reduce trend or scale systematic errors in the data, because they are perfectly correlated across all data points. The second is that WLS regression is liable to produce seriously inaccurate trend estimates where a scale systematic error is involved.

With linear regression models, the estimated trend of the sum of several components equals the sum of the trends of the individual components. Moreover, provided there is no correlation between errors in the different components, the trend estimation error for their sum equals the sum in quadrature of the trend estimation errors for the individual components. That enables one to place a theoretical lower limit on the ΔAPOClimate trend uncertainty, by adding in quadrature identified trend and scale systematic errors.

The largest contribution to error in ΔAPOClimate comes from its ΔAPOOBS component. That in turn contains a scale systematic error, of (at ± 1σ) 2% of a trend of circa −19.8 per meg yr−1, or –0.396 per meg yr−1, and three trend errors, or 0.3, 0.2 and 0.1 per meg yr−1. Since these four errors sources are all independent, they can be added in quadrature. So can the scale error of 0.135 per meg yr−1 in ΔAPOAtmD. I will leave errors in the remaining components,  ΔAPOFF and ΔAPOCant, out of account; although it is highly likely that they also contain trend or scale systematic errors their magnitude is small and when added in quadrature they would not add significantly to the total.

Adding all the quantified error sources in quadrature gives an irreducible minimum ΔAPOClimate trend uncertainty of ± 0.56 per meg yr−1. Even if the much lower, clearly incorrect, WLS estimate of the (δO2/N2) trend of −16.0 were substituted for the OLS-based trend, the irreducible trend uncertainty would still be 0.51 per meg yr−1.  Resplandy et al.’s ± 0.15 ΔAPOClimate trend uncertainty estimate is completely infeasible!

 

Nicholas Lewis                                                                       7 November 2018

By | 2018-11-24T18:28:08+00:00 November 24th, 2018|General climate science, Posts originally elsewhere|0 Comments

About the Author:

Leave A Comment

This site uses Akismet to reduce spam. Learn how your comment data is processed.