slr_change: Test whether the rate of sea level rise has changed

Description Usage Arguments Details Value See Also

View source: R/slr_change.R

Description

Fit a piecewise linear model to test whether a recent change in rate of sea level rise has occurred.

Usage

1
2
3
4
5
6
7
8
9
slr_change(
  .data,
  .sl,
  .dt,
  .span = 20L,
  .mode = c("year", "duration", "count"),
  t_fit = FALSE,
  retain_model = FALSE
)

Arguments

.data

Source data frame for data. Use NULL if no data frame is used and all data is passed from the enclosing environment.

.sl

The data variable (usually found in the data frame) showing sea level or mean sea level. Must be a named variable, not an expression.

.dt

Data variable containing corresponding midpoint dates for the period of averaging used to calculate .sl. Must be a named variable, not an expression. Midpoint dates for a given month can be approximated with as.Date(paste(year, month, 15, sep = '-').

.span

Span defining the "recent" period for fitting the model. integer or difftime object. Interpretation of .span depends on the value of the by_year. If by_year == TRUE, .span must be an integer, and is interpreted as the number of years included in the "recent" period for fitting the model. If by_year == FALSE, and .span is an integer, it is interpreted as the number of records to include in the 'recent' period.

.mode

one of c('year', 'duration', 'count') indicating whether the .span is expressed in years, time coordinates, or number of observations. If .mode == 'year', slopes are scaled to be expressed on a per year basis. Scaling requires the time coordinate to inherit from R Date or POSIXt objects. Rates based on Dates are converted to annual values by multiplying by 365.2422 days / year. Rates based on POSIXt objects are converted to annual values by multiplying slope estimates by multiplying by 31556926 seconds / year (with a warning about rounding error). No scaling is done for .mode == 'duration' or .mode == 'count'.

t_fit

Should the underlying generalized linear model's correlation structure be fit based on the time coordinate (.dt), or only on the sequence of observations in the data? Setting this to TRUE is safer if you are uncertain of the sequence of observations in the source data, or significant missing values in the data, but it significantly slows model fitting. Results tend to be very similar for complete or near complete data.

retain_model

Boolean. If TRUE, the piecewise linear model is returned via the 'model' slot of the return value.

Details

A frequent question raised in analysis of recent sea level records is whether the rate of sea level rise is increasing, as predicted by numerous climate models over the years. This is a simple question that is rather more complex to answer statistically than it at first appears.

This function checks for a change is rate of sea level rise by fitting a piecewise linear model and testing if a change in slope at a (user-defined) breakpoint improves the statistical fit enough to matter. A piecewise linear model fits the data with (here) two linear segments that join at a breakpoint or knot. The underlying model is a generalized linear model, with an AR1 autocorrelated error term.

The AR1 process, an "autoregressive process of order 1", is a simple assumption about correlation structures that reduces the risk of overstating the precision of trend estimates. Actual tidal data series tend to have more complex autocorrelation structure. If working with high frequency data, short-term structure should be removed before using this function, typically by calculating average sea level over convenient and meaningful time periods, like months or years. See additional comments on use of the AR1 model in analysis of sea level data in the help file for slr_slope().

This function is ripe for abuse. The problem is that rates of sea level rise year over year vary. By chance, some periods of the historic record will have higher rates of sea level rise than others. This makes it far too easy to "cherry pick" starting dates for a recent period of record, and claim a statistically meaningful change in rate of sea level rise. Use with appropriate humility.

For a more robust approach (with lower statistical power, but deeper justification), see the function all_sr_change()

Value

a list of three (optionally four) items, with the following components:

summary (Vector of numerical values.)

slope_ratio

The ratio of the slope for the "recent" period to the slope for the remainder of the historic record.

slope_old

The estimated slope of changes in sea level prior to the "recent" period.

slope_old_err

Standard error of slope_old.

slope_recent

Estimated slope during the "Recent" period.

slope_recent_err

Standard error of slope_recent.

details (Vector of numerical values.)

l_ratio

likelihood ratio comparing the piecewise linear model to a linear model.

p-value

the p-value (by ANOVA) comparing the two models.

delta-AIC

Change in AIC between the two models. A value less than about -2 is evidence that the piecewise model would outperform the strict linear model for "predicting" past sea levels.

sample

Sample size on which model is based, including both older and recent data.

recents

the number of samples in the "recent" time period.

settings (Vector of strings.)

mode

The value of the .mode argument specifying how to process the data.

span

The value of the .span argument, as a string, not a number. If .mode = 'year', the default, this will be the number of years in the 'recent' period. Otherwise, it will represent either a time period, (.mode = 'duration') or a number of samples (.mode = 'count').

cor_struct

Either 'Order-based' (for t_fit = FALSE) or 'Time-based', for t_fit = TRUE.

model The underlying piecewise GLS model, returned only if retain_model = TRUE

See Also

Other sea level rate functions: get_sl_trend(), slr_change_comp(), slr_slope()


ccb60/SLRSIM documentation built on Jan. 21, 2022, 1:31 a.m.