Version: 1.1.6
Date: 2023-10-10
Title: Non-Parametric Trend Tests and Change-Point Detection
Depends: R (≥ 3.0)
Description: The analysis of environmental data often requires the detection of trends and change-points. This package includes tests for trend detection (Cox-Stuart Trend Test, Mann-Kendall Trend Test, (correlated) Hirsch-Slack Test, partial Mann-Kendall Trend Test, multivariate (multisite) Mann-Kendall Trend Test, (Seasonal) Sen's slope, partial Pearson and Spearman correlation trend test), change-point detection (Lanzante's test procedures, Pettitt's test, Buishand Range Test, Buishand U Test, Standard Normal Homogeinity Test), detection of non-randomness (Wallis-Moore Phase Frequency Test, Bartels rank von Neumann's ratio test, Wald-Wolfowitz Test) and the two sample Robust Rank-Order Distributional Test.
Imports: extraDistr (≥ 1.8.0)
Suggests: strucchange, Kendall, psych, datasets
Classification/MSC-2010: 62F03, 62G10, 62M10, 62P12
License: GPL-3
Encoding: UTF-8
ZipData: yes
NeedsCompilation: yes
RoxygenNote: 7.2.2
Packaged: 2023-10-10 07:57:22 UTC; thorsten
Author: Thorsten Pohlert ORCID iD [aut, cre]
Maintainer: Thorsten Pohlert <thorsten.pohlert@gmx.de>
Repository: CRAN
Date/Publication: 2023-10-10 08:50:02 UTC

Simulated data of Page (1955) as test-example for change-point detection

Description

Simulated data of Page (1955) as test-example for change-point detection taken from Table 1 of Pettitt (1979)

Usage

data(PagesData)

Format

a vector that contains 40 elements

Details

According to the publication of Pettitt (1979), the series comprise a significant p = 0.014 change-point at i = 17. The function pettitt.test computes the same U statistics as given by Pettitt (1979) in Table1, row 4.

References

Page, E. S. (1954), A test for a change in a parameter occuring at an unknown point. Biometrika 41, 100–114.

Pettitt, A. N., (1979). A non-parametric approach to the change point problem. Journal of the Royal Statistical Society Series C, Applied Statistics 28, 126–135.

See Also

pettitt.test

Examples

data(PagesData)
pettitt.test(PagesData)

Bartels Test for Randomness

Description

Performes a rank version of von Neumann's ratio test as proposed by Bartels. The null hypothesis of randomness is tested against the alternative hypothesis

Usage

bartels.test(x)

Arguments

x

a vector of class "numeric" or a time series object of class "ts"

Details

In this function, the test is implemented as given by Bartels (1982), where the ranks r_1, \ldots, r_n of the X_i, \ldots, X_n are used for the statistic:

T = \frac{\sum_{i=1}^n (r_i - r_{i+1})^2}{\sum_{i=1}^n (r_i - \bar{r})^2}

As proposed by Bartels (1982), the p-value is calculated for sample sizes in the range of (10 \le n < 100) with the non-standard beta distribution for the range 0 \le x \le 4 with parameters:

a = b = \frac{5 n \left( n + 1\right) \left(n - 1\right)^2} {2 \left(n - 2\right) \left(5n^2 - 2n - 9\right)} - \frac{1}{2}

For sample sizes n \ge 100 a normal approximation with N(2, 20/(5n + 7)) is used for p-value calculation.

Value

A list with class "htest"

data.name

character string that denotes the input data

p.value

the p-value

statistic

the test statistic

alternative

the alternative hypothesis

method

character string that denotes the test

Note

The current function is for complete observations only.

References

R. Bartels (1982), The Rank Version of von Neumann's Ratio Test for Randomness, Journal of the American Statistical Association 77, 40–46.

See Also

ww.test, wm.test

Examples

# Example from Schoenwiese (1992, p. 113)
## Number of frost days in April at Munich from 1957 to 1968
## 
frost <- ts(data=c(9,12,4,3,0,4,2,1,4,2,9,7), start=1957)
bartels.test(frost)

## Example from Sachs (1997, p. 486)
x <- c(5,6,2,3,5,6,4,3,7,8,9,7,5,3,4,7,3,5,6,7,8,9)
bartels.test(x)

## Example from Bartels (1982, p. 43)
x <- c(4, 7, 16, 14, 12, 3, 9, 13, 15, 10, 6, 5, 8, 2, 1, 11, 18, 17)
bartels.test(x)
 

Buishand Range Test for Change-Point Detection

Description

Performes the Buishand range test for change-point detection of a normal variate.

Usage

br.test(x, m = 20000)

Arguments

x

a vector of class "numeric" or a time series object of class "ts"

m

numeric, number of Monte-Carlo replicates, defaults to 20000

Details

Let X denote a normal random variate, then the following model with a single shift (change-point) can be proposed:

x_i = \left\{ \begin{array}{lcl} \mu + \epsilon_i, & \qquad & i = 1, \ldots, m \\ \mu + \Delta + \epsilon_i & \qquad & i = m + 1, \ldots, n \\ \end{array} \right.

with \epsilon \approx N(0,\sigma). The null hypothesis \Delta = 0 is tested against the alternative \Delta \ne 0.

In the Buishand range test, the rescaled adjusted partial sums are calculated as

S_k = \sum_{i=1}^k \left(x_i - \hat{x}\right) \qquad (1 \le i \le n)

The test statistic is calculated as:

Rb = \left(\max S_k - \min S_k\right) / \sigma

.

The p.value is estimated with a Monte Carlo simulation using m replicates.

Critical values based on m = 19999 Monte Carlo simulations are tabulated for Rb / \sqrt{n} by Buishand (1982).

Value

A list with class "htest" and "cptest"

data.name

character string that denotes the input data

p.value

the p-value

statistic

the test statistic

null.value

the null hypothesis

estimates

the time of the probable change point

alternative

the alternative hypothesis

method

character string that denotes the test

data

numeric vector of Sk for plotting

Note

The current function is for complete observations only.

References

T. A. Buishand (1982), Some Methods for Testing the Homogeneity of Rainfall Records, Journal of Hydrology 58, 11–27.

G. Verstraeten, J. Poesen, G. Demaree, C. Salles (2006), Long-term (105 years) variability in rain erosivity as derived from 10-min rainfall depth data for Ukkel (Brussels, Belgium): Implications for assessing soil erosion rates. Journal of Geophysical Research 111, D22109.

See Also

efp sctest.efp

Examples

data(Nile)
(out <- br.test(Nile))
plot(out)


data(PagesData) ; br.test(PagesData)
 

Buishand U Test for Change-Point Detection

Description

Performes the Buishand U test for change-point detection of a normal variate.

Usage

bu.test(x, m = 20000)

Arguments

x

a vector of class "numeric" or a time series object of class "ts"

m

numeric, number of Monte-Carlo replicates, defaults to 20000

Details

Let X denote a normal random variate, then the following model with a single shift (change-point) can be proposed:

x_i = \left\{ \begin{array}{lcl} \mu + \epsilon_i, & \qquad & i = 1, \ldots, m \\ \mu + \Delta + \epsilon_i & \qquad & i = m + 1, \ldots, n \\ \end{array} \right.

with \epsilon \approx N(0,\sigma). The null hypothesis \Delta = 0 is tested against the alternative \Delta \ne 0.

In the Buishand U test, the rescaled adjusted partial sums are calculated as

S_k = \sum_{i=1}^k \left(x_i - \bar{x}\right) \qquad (1 \le i \le n)

The sample standard deviation is

D_x = \sqrt{n^{-1} \sum_{i=1}^n \left(x_i - \bar{x}\right)}

The test statistic is calculated as:

U = \left[n \left(n + 1 \right) \right]^{-1} \sum_{k=1}^{n-1} \left(S_k / D_x \right)^2

.

The p.value is estimated with a Monte Carlo simulation using m replicates.

Critical values based on m = 19999 Monte Carlo simulations are tabulated for U by Buishand (1982, 1984).

Value

A list with class "htest" and "cptest"

data.name

character string that denotes the input data

p.value

the p-value

statistic

the test statistic

null.value

the null hypothesis

estimates

the time of the probable change point

alternative

the alternative hypothesis

method

character string that denotes the test

data

numeric vector of Sk for plotting

Note

The current function is for complete observations only.

References

T. A. Buishand (1982), Some Methods for Testing the Homogeneity of Rainfall Records, Journal of Hydrology 58, 11–27.

T. A. Buishand (1984), Tests for Detecting a Shift in the Mean of Hydrological Time Series, Journal of Hydrology 73, 51–69.

See Also

efp sctest.efp

Examples

data(Nile)
(out <- bu.test(Nile))
plot(out)

data(PagesData)
bu.test(PagesData)
 

Cox and Stuart Trend Test

Description

Performes the non-parametric Cox and Stuart trend test

Usage

cs.test(x)

Arguments

x

a vector or a time series object of class "ts"

Details

First, the series is devided by three. It is compared, whether the data of the first third of the series are larger or smaller than the data of the last third of the series. The test statistic of the Cox-Stuart trend test for n > 30 is calculated as:

z = \frac{|S - \frac{n}{6}|}{\sqrt{\frac{n}{12}}}

where S denotes the maximum of the number of signs, i.e. + or -, respectively. The z-statistic is normally distributed. For n \le 30 a continuity correction of -0.5 is included in the denominator.

Value

An object of class "htest"

method

a character string indicating the chosen test

data.name

a character string giving the name(s) of the data

statistic

the Cox-Stuart z-value

alternative

a character string describing the alternative hypothesis

p.value

the p-value for the test

Note

NA values are omitted. Many ties in the series will lead to reject H0 in the present test.

References

L. Sachs (1997), Angewandte Statistik. Berlin: Springer.

C.-D. Schoenwiese (1992), Praktische Statistik. Berlin: Gebr. Borntraeger.

D. R. Cox and A. Stuart (1955), Quick sign tests for trend in location and dispersion. Biometrika 42, 80-95.

See Also

mk.test

Examples

## Example from Schoenwiese (1992, p. 114)
## Number of frost days in April at Munich from 1957 to 1968
## z = -0.5, Accept H0
frost <- ts(data=c(9,12,4,3,0,4,2,1,4,2,9,7), start=1957)
cs.test(frost)

## Example from Sachs (1997, p. 486-487)
## z ~ 2.1, Reject H0 on a level of p = 0.0357
x <- c(5,6,2,3,5,6,4,3,7,8,9,7,5,3,4,7,3,5,6,7,8,9)
cs.test(x)

cs.test(Nile)


Correlated Seasonal Mann-Kendall Test

Description

Performs a Seasonal Mann-Kendall test under the presence of correlated seasons.

Usage

csmk.test(x, alternative = c("two.sided", "greater", "less"))

Arguments

x

a time series object with class ts comprising >= 2 seasons; NA values are not allowed

alternative

the alternative hypothesis, defaults to two.sided

Details

The Mann-Kendall scores are first computed for each season seperately. The variance - covariance matrix is computed according to Libiseller and Grimvall (2002). Finally the corrected Z-statistics for the entire series is calculated as follows, whereas a continuity correction is employed for n \le 10:

z = \frac{\mathbf{1}^T \mathbf{S}} {\sqrt{\mathbf{1}^T \mathbf{\Gamma}~\mathbf{1}}}

where

z denotes the quantile of the normal distribution, 1 indicates a vector with all elements equal to one, \mathbf{S} is the vector of Mann-Kendall scores for each season and \mathbf{\Gamma} denotes the variance - covariance matrix.

Value

An object with class "htest"

data.name

character string that denotes the input data

p.value

the p-value for the entire series

statistic

the z quantile of the standard normal distribution for the entire series

null.value

the null hypothesis

estimates

the estimates S and varS for the entire series

alternative

the alternative hypothesis

method

character string that denotes the test

cov

the variance - covariance matrix

Note

Ties are not corrected. Current Version is for complete observations only.

References

Hipel, K.W. and McLeod, A.I. (1994), Time Series Modelling of Water Resources and Environmental Systems. New York: Elsevier Science.

Libiseller, C. and Grimvall, A., (2002), Performance of partial Mann-Kendall tests for trend detection in the presence of covariates. Environmetrics 13, 71–84, doi:10.1002/env.507.

See Also

cor, cor.test, mk.test, smk.test

Examples

csmk.test(nottem)

Monthly concentration of particle bound HCB, River Rhine

Description

Time series of monthly concentration of particle bound Hexachlorobenzene (HCB) in \mug/kg at six different monitoring sites at the River Rhine, 1995.1-2006.12

Usage

data(hcb)

Format

a time series object of class "mts"

Details

NO DATA values in the series were filled with estimated values using linear interpolation (see approx.

The Rhine Kilometer (RKM) is in increasing order from source to mouth of the River Rhine.

Source

International Commission for the Protection of the River Rhine

References

T. Pohlert, G. Hillebrand, V. Breitung (2011), Trends of persistent organic pollutants in the suspended matter of the River Rhine, Hydrological Processes 25, 3803–3817. doi:10.1002/hyp.8110

Examples

data(hcb)
plot(hcb)
mult.mk.test(hcb)

Lanzante's Test for Change Point Detection

Description

Performes a non-parametric test after Lanzante in order to test for a shift in the central tendency of a time series. The null hypothesis, no shift, is tested against the alternative, shift.

Usage

lanzante.test(x, method = c("wilcox.test", "rrod.test"))

Arguments

x

a vector of class "numeric" or a time series object of class "ts"

method

the test method. Defaults to "wilcox.test".

Details

Let X denote a continuous random variable, then the following model with a single shift (change-point) can be proposed:

x_i = \left\{ \begin{array}{lcl} \theta + \epsilon_i, & \qquad & i = 1, \ldots, m \\ \theta + \Delta + \epsilon_i & \qquad & i = m + 1, \ldots, n \\ \end{array} \right.

with \theta(\epsilon) = 0. The null hypothesis, H:\Delta = 0 is tested against the alternative A:\Delta \ne 0.

First, the data are transformed into increasing ranks and for each time-step the adjusted rank sum is computed:

U_k = 2 \sum_{i=1}^k r_i - k \left(n + 1\right) \qquad k = 1, \ldots, n

The probable change point is located at the absolute maximum of the statistic:

m = k(\max |U_k|)

.

For method = "wilcox.test" the Wilcoxon-Mann-Whitney two-sample test is performed, using m to split the series. Otherwise, the robust rank-order distributional test (rrod.test is performed.

Value

A list with class "htest" and "cptest".

References

Lanzante, J. R. (1996), Resistant, robust and non-parametric techniques for the analysis of climate data: Theory and examples, including applications to historical radiosonde station data, Int. J. Clim., 16, 1197–1226.

See Also

pettitt.test

Examples

data(maxau) ; plot(maxau[,"s"])
s.res <- lanzante.test(maxau[,"s"])
n <- s.res$nobs
i <- s.res$estimate
s.1 <- mean(maxau[1:i,"s"])
s.2 <- mean(maxau[(i+1):n,"s"])
s <- ts(c(rep(s.1,i), rep(s.2,(n-i))))
tsp(s) <- tsp(maxau[,"s"])
lines(s, lty=2)
print(s.res)


data(PagesData) ; lanzante.test(PagesData)

Annual suspended sediment concentration and flow data, River Rhine

Description

Annual time series of average suspended sediment concentration (s) in mg/l and average discharge (Q) in m^3 / s at the River Rhine, 1965.1-2009.1

Usage

data(maxau)

Format

a time series object of class "mts"

Source

Bundesanstalt für Gewässerkunde, Koblenz, Deutschland (Federal Institute of Hydrology, Koblenz, Germany)

Examples

data(maxau)
plot(maxau)

Mann-Kendall Trend Test

Description

Performs the Mann-Kendall Trend Test

Usage

mk.test(x, alternative = c("two.sided", "greater", "less"), continuity = TRUE)

Arguments

x

a vector of class "numeric" or a time series object of class "ts"

alternative

the alternative hypothesis, defaults to two.sided

continuity

logical, indicates whether a continuity correction should be applied, defaults to TRUE.

Details

The null hypothesis is that the data come from a population with independent realizations and are identically distributed. For the two sided test, the alternative hypothesis is that the data follow a monotonic trend. The Mann-Kendall test statistic is calculated according to:

S = \sum_{k = 1}^{n-1} \sum_{j = k + 1}^n \mathrm{sgn}\left(x_j - x_k\right)

with \mathrm{sgn} the signum function (see sign).

The mean of S is \mu = 0. The variance including the correction term for ties is

\sigma^2 = \left\{n \left(n-1\right)\left(2n+5\right) - \sum_{j=1}^p t_j\left(t_j - 1\right)\left(2t_j+5\right) \right\} / 18

where p is the number of the tied groups in the data set and t_j is the number of data points in the j-th tied group. The statistic S is approximately normally distributed, with

z = S / \sigma

If continuity = TRUE then a continuity correction will be employed:

z = \mathrm{sgn}(S) ~ \left(|S| - 1\right) / \sigma

The statistic S is closely related to Kendall's \tau:

\tau = S / D

where

D = \left[\frac{1}{2}n\left(n-1\right)- \frac{1}{2}\sum_{j=1}^p t_j\left(t_j - 1\right)\right]^{1/2} \left[\frac{1}{2}n\left(n-1\right) \right]^{1/2}

Value

A list with class "htest"

data.name

character string that denotes the input data

p.value

the p-value

statistic

the z quantile of the standard normal distribution

null.value

the null hypothesis

estimates

the estimates S, varS and tau

alternative

the alternative hypothesis

method

character string that denotes the test

Note

Current Version is for complete observations only.

References

Hipel, K.W. and McLeod, A.I. (1994), Time Series Modelling of Water Resources and Environmental Systems. New York: Elsevier Science.

Libiseller, C. and Grimvall, A., (2002), Performance of partial Mann-Kendall tests for trend detection in the presence of covariates. Environmetrics 13, 71–84, doi:10.1002/env.507.

See Also

cor.test, MannKendall, partial.mk.test, sens.slope

Examples

data(Nile)
mk.test(Nile, continuity = TRUE)

##
n <- length(Nile)
cor.test(x=(1:n),y=Nile, meth="kendall", continuity = TRUE)


Multivariate (Multisite) Mann-Kendall Test

Description

Performs a Multivariate (Multisite) Mann-Kendall test.

Usage

mult.mk.test(x, alternative = c("two.sided", "greater", "less"))

Arguments

x

a time series object of class "ts"

alternative

the alternative hypothesis, defaults to two.sided

Details

The Mann-Kendall scores are first computed for each variate (side) seperately.

S = \sum_{k = 1}^{n-1} \sum_{j = k + 1}^n \mathrm{sgn}\left(x_j - x_k\right)

with \mathrm{sgn} the signum function (see sign).

The variance - covariance matrix is computed according to Libiseller and Grimvall (2002).

\Gamma_{xy} = \frac{1}{3} \left[K + 4 \sum_{j=1}^n R_{jx} R_{jy} - n \left(n + 1 \right) \left(n + 1 \right) \right]

with

K = \sum_{1 \le i < j \le n} \mathrm{sgn} \left\{ \left( x_j - x_i \right) \left( y_j - y_i \right) \right\}

and

R_{jx} = \left\{ n + 1 + \sum_{i=1}^n \mathrm{sgn} \left( x_j - x_i \right) \right\} / 2

Finally, the corrected z-statistics for the entire series is calculated as follows, whereas a continuity correction is employed for n \le 10:

z = \frac{\sum_{i=1}^d S_i}{\sqrt{\sum_{j=1}^d \sum_{i=1}^d \Gamma_{ij}}}

where

z denotes the quantile of the normal distribution S is the vector of Mann-Kendall scores for each variate (site) 1 \le i \le d and \Gamma denotes symmetric variance - covariance matrix.

Value

An object with class "htest"

data.name

character string that denotes the input data

p.value

the p-value for the entire series

statistic

the z quantile of the standard normal distribution for the entire series

null.value

the null hypothesis

estimates

the estimates S and varS for the entire series

alternative

the alternative hypothesis

method

character string that denotes the test

cov

the variance - covariance matrix

Note

Ties are not corrected. Current Version is for complete observations only.

References

Hipel, K.W. and McLeod, A.I. (1994), Time Series Modelling of Water Resources and Environmental Systems. New York: Elsevier Science.

Lettenmeier, D.P. (1988), Multivariate nonparametric tests for trend in water quality. Water Resources Bulletin 24, 505–512.

Libiseller, C. and Grimvall, A. (2002), Performance of partial Mann-Kendall tests for trend detection in the presence of covariates. Environmetrics 13, 71–84, doi:10.1002/env.507.

See Also

cor, cor.test, mk.test, smk.test

Examples

data(hcb)
mult.mk.test(hcb)


Partial Correlation Trend Test

Description

Performs a partial correlation trend test with either Pearson's or Spearman's correlation coefficients (r(tx.z)).

Usage

partial.cor.trend.test(x, z, method = c("pearson", "spearman"))

Arguments

x

a "vector" or "ts" object that contains the variable, which is tested for trend (i.e. correlated with time)

z

a "vector" or "ts" object that contains the co-variate, which will be partialled out

method

a character string indicating which correlation coefficient is to be computed. One of "pearson" (default) or "spearman", can be abbreviated.

Details

This function performs a partial correlation trend test using either the "pearson" correlation coefficient, or the "spearman" rank correlation coefficient (Hipel and McLoed (1994), p. 882). The partial correlation coefficient for the response variable "x" with time "t", when the effect of the explanatory variable "z" is partialled out, is defined as:

r_{tx.z} = \frac{r_{tx} - r_{tz}~r_{xz}} {\sqrt{1 - r_{tz}^2} ~ \sqrt{1-r_{xz}^2}}

The H0: r_{tx.z} = 0 (i.e. no trend for "x", when effect of "z" is partialled out) is tested against the alternate Hypothesis, that there is a trend for "x", when the effect of "z" is partialled out.

The partial correlation coefficient is tested for significance with the student t distribution on df = n - 2 degree of freedom.

Value

An object of class "htest"

method

a character string indicating the chosen test

data.name

a character string giving the name(s) of the data

statistic

the value of the test statistic

estimate

the partial correlation coefficient r(tx.z)

parameter

the degrees of freedom of the test statistic in the case that it follows a t distribution

alternative

a character string describing the alternative hypothesis

p.value

the p-value of the test

null.value

The value of the null hypothesis

Note

Current Version is for complete observations only.

References

Hipel, K.W. and McLeod, A.I. (1994), Time Series Modelling of Water Resources and Environmental Systems. New York: Elsevier Science.

Bahrenberg, G., Giese, E. and Nipper, J., (1992): Statistische Methoden in der Geographie, Band 2 Multivariate Statistik, Teubner, Stuttgart.

See Also

cor, cor.test, partial.r, partial.mk.test,

Examples

data(maxau)
a <- tsp(maxau) ; tt <- a[1]:a[2]
s <- maxau[,"s"] ; Q <- maxau[,"Q"]
maxau.df <- data.frame(Year = tt, s =s, Q = Q)
plot(maxau.df)

partial.cor.trend.test(s,Q, method="pearson")
partial.cor.trend.test(s,Q, method="spearman")


Partial Mann-Kendall Trend Test

Description

Performs a partial Mann-Kendall Trend Test

Usage

partial.mk.test(x, y, alternative = c("two.sided", "greater", "less"))

Arguments

x

a "vector" or "ts" object that contains the variable, which is tested for trend (i.e. correlated with time)

y

a "vector" or "ts" object that contains the variable, which effect on "x" is partialled out

alternative

character, the alternative method; defaults to "two.sided"

Details

According to Libiseller and Grimvall (2002), the test statistic for x with its covariate y is

z = \frac{S_x - r_{xy} S_y}{\left[ \left( 1 - r_{xy}^2 \right) n \left(n - 1 \right) \left(2 n + 5 \right) / 18 \right]^{0.5}}

where the correlation r is calculated as:

r_{xy} = \frac{\sigma_{xy}} {n \left(n - 1\right) \left(2 n + 5 \right) / 18}

The conditional covariance between x and y is

\sigma_{xy} = \frac{1}{3} \left[K + 4 \sum_{j=1}^n R_{jx} R_{jy} - n \left(n + 1 \right) \left(n + 1 \right) \right]

with

K = \sum_{1 \le i < j \le n} \mathrm{sgn} \left\{ \left( x_j - x_i \right) \left( y_j - y_i \right) \right\}

and

R_{jx} = \left\{ n + 1 + \sum_{i=1}^n \mathrm{sgn} \left( x_j - x_i \right) \right\} / 2

Value

A list with class "htest"

method

a character string indicating the chosen test

data.name

a character string giving the name(s) of the data

statistic

the value of the test statistic

estimate

the Mann-Kendall score S, the variance varS and the correlation between x and y

alternative

a character string describing the alternative hypothesis

p.value

the p-value of the test

null.value

the null hypothesis

Note

Current Version is for complete observations only. The test statistic is not corrected for ties.

References

Libiseller, C. and Grimvall, A., (2002). Performance of partial Mann-Kendall tests for trend detection in the presence of covariates. Environmetrics 13, 71–84, doi:10.1002/env.507.

See Also

partial.cor.trend.test,

Examples

data(maxau)
s <- maxau[,"s"]; Q <- maxau[,"Q"]
partial.mk.test(s,Q)


Pettitt's Test for Change-Point Detection

Description

Performes a non-parametric test after Pettitt in order to test for a shift in the central tendency of a time series. The H0-hypothesis, no change, is tested against the HA-Hypothesis, change.

Usage

pettitt.test(x)

Arguments

x

a vector of class "numeric" or a time series object of class "ts"

Details

In this function, the test is implemented as given by Verstraeten et. al. (2006), where the ranks r_1, \ldots, r_n of the X_i, \ldots, X_n are used for the statistic:

U_k = 2 \sum_{i=1}^k r_i - k \left(n + 1\right) \qquad k = 1, \ldots, n

The test statistic is the maximum of the absolute value of the vector:

\hat{U} = \max |U_k|

.

The probable change-point K is located where \hat{U} has its maximum. The approximate probability for a two-sided test is calculated according to

p = 2 \exp^{-6K^2 / (T^3 + T^2)}

Value

A list with class "htest" and "cptest"

Note

The current function is for complete observations only. The approximate probability is good for p \le 0.5.

References

CHR (ed., 2010), Das Abflussregime des Rheins und seiner Nebenfluesse im 20. Jahrhundert, Report no I-22 of the CHR, p. 172.

Pettitt, A. N. (1979), A non-parametric approach to the change point problem. Journal of the Royal Statistical Society Series C, Applied Statistics 28, 126-135.

G. Verstraeten, J. Poesen, G. Demaree, C. Salles (2006), Long-term (105 years) variability in rain erosivity as derived from 10-min rainfall depth data for Ukkel (Brussels, Belgium): Implications for assessing soil erosion rates. Journal of Geophysical Research 111, D22109.

See Also

efp sctest.efp

Examples

data(maxau) ; plot(maxau[,"s"])
s.res <- pettitt.test(maxau[,"s"])
n <- s.res$nobs
i <- s.res$estimate
s.1 <- mean(maxau[1:i,"s"])
s.2 <- mean(maxau[(i+1):n,"s"])
s <- ts(c(rep(s.1,i), rep(s.2,(n-i))))
tsp(s) <- tsp(maxau[,"s"])
lines(s, lty=2)
print(s.res)


data(PagesData) ; pettitt.test(PagesData)
 

Plotting cptest-objects

Description

Plotting method for objects inheriting from class "cptest"

Usage

## S3 method for class 'cptest'
plot(x, ...)

Arguments

x

an object of class "cptest"

...

further arguments, currently ignored

Examples

data(Nile)
(out <- br.test(Nile))
par(mfrow=c(2,1))
plot(Nile) ; plot(out)


Robust Rank-Order Distributional Test

Description

Performs Fligner-Pollicello robust rank-order distributional test for location.

Usage

rrod.test(x, ...)

## Default S3 method:
rrod.test(x, y, alternative = c("two.sided", "less", "greater"), ...)

## S3 method for class 'formula'
rrod.test(formula, data, subset, na.action, ...)

Arguments

x

a vector of data values.

...

further arguments to be passed to or from methods.

y

an optional numeric vector of data values.

alternative

the alternative hypothesis. Defaults to "two.sided".

formula

a formula of the form response ~ group where response gives the data values and group a vector or factor of the corresponding groups.

data

an optional matrix or data frame (or similar: see model.frame) containing the variables in the formula formula. By default the variables are taken from environment(formula).

subset

an optional vector specifying a subset of observations to be used.

na.action

a function which indicates what should happen when the data contain NAs. Defaults to getOption("na.action").

Details

The non-parametric RROD two-sample test can be used to test for differences in location, whereas it does not assume variance homogeneity.

Let X and Y denote two samples with sizes n_x and n_y of a continuous variable.First, the combined sample is transformed into ranks in increasing order. Let S_{xi} and S_{yj} denote the counts of Y (X) values having a lower rank than x_i (y_j). The mean counts are:

\bar{S}_x = \sum_{i=1}^{n_x} S_{xi} / n_x

\bar{S}_y = \sum_{j=1}^{n_y} S_{yj} / n_y

The variances are:

s^2_{Sx} = \sum_{i=1}^{n_x} \left( S_{xi} - \bar{S}_x \right)^2

s^2_{Sy} = \sum_{j=1}^{n_y} \left( S_{yj} - \bar{S}_y \right)^2

The test statistic is:

z = \frac{1}{2}~ \frac{n_x \bar{S}_x - n_y \bar{S}_y} {\left( \bar{S}_x \bar{S}_y + s^2_{Sx} + s^2_{Sy} \right)^{1/2}}

The two samples have significantly different location parameters, if |z| > z_{1-\alpha/2}. The function calculates the p-values of the null hypothesis for the selected alternative than can be "two.sided", "greater" or "less".

Value

A list with class "htest".

References

Fligner, M. A., Pollicello, G. E. III. (1981), Robust Rank Procedures for the Behrens-Fisher Problem, Journal of the American Statistical Association, 76, 162–168.

Lanzante, J. R. (1996), Resistant, robust and non-parametric techniques for the analysis of climate data: Theory and examples, including applications to historical radiosonde station data, Int. J. Clim., 16, 1197–1226.

Siegel, S. and Castellan, N. (1988), Nonparametric Statistics For The Behavioural Sciences, New York: McCraw-Hill.

See Also

wilcox.test

Examples

## Two-sample test.
## Hollander & Wolfe (1973), 69f.
## Permeability constants of the human chorioamnion (a placental
##  membrane) at term (x) and between 12 to 26 weeks gestational
##  age (y).  The alternative of interest is greater permeability
##  of the human chorioamnion for the term pregnancy.
x <- c(0.80, 0.83, 1.89, 1.04, 1.45, 1.38, 1.91, 1.64, 0.73, 1.46)
y <- c(1.15, 0.88, 0.90, 0.74, 1.21)
rrod.test(x, y, alternative = "g")

## Formula interface.
boxplot(Ozone ~ Month, data = airquality)
rrod.test(Ozone ~ Month, data = airquality,
            subset = Month %in% c(5, 8)) 

Seasonal Sen's Slope

Description

Computes seasonal Sen's slope for linear rate of change

Usage

sea.sens.slope(x)

Arguments

x

a time series object of class "ts"

Details

Acccording to Hirsch et al. (1982) the seasonal Sen's slope is calculated as follows:

d_{ijk} = \frac{x_{ij} - x_{ik}}{j - k}

for each (x_{ij}, x_{ik}) pair i = 1, \ldots, m, where (1 \leq k < j \leq n_i and n_i is the number of known values in the i-th season. The seasonal slope estimator is the median of the d_{ijk} values.

Value

numeric, Seasonal Sen's slope.

Note

Current Version is for complete observations only.

References

Hipel, K.W. and McLeod, A.I. (1994), Time Series Modelling of Water Resources and Environmental Systems. New York: Elsevier Science.

Hirsch, R., J. Slack, R. Smith (1982), T echniques of Trend Analysis for Monthly Water Quality Data. Water Resources Research 18, 107-121.

Sen, P.K. (1968), Estimates of the regression coefficient based on Kendall's tau, Journal of the American Statistical Association 63, 1379–1389.

See Also

smk.test,

Examples

sea.sens.slope(nottem)


Sen's slope

Description

Computes Sen's slope for linear rate of change and corresponding confidence intervalls

Usage

sens.slope(x, conf.level = 0.95)

Arguments

x

numeric vector or a time series object of class "ts"

conf.level

numeric, the level of significance

Details

This test computes both the slope (i.e. linear rate of change) and confidence levels according to Sen's method. First, a set of linear slopes is calculated as follows:

d_{k} = \frac{x_j - x_i}{j - i}

for \left(1 \le i < j \le n \right), where d is the slope, x denotes the variable, n is the number of data, and i, j are indices.

Sen's slope is then calculated as the median from all slopes: b_{Sen} = \textnormal{median}(d_k).

This function also computes the upper and lower confidence limits for sens slope.

Value

A list of class "htest".

estimates

numeric, Sen's slope

data.name

character string that denotes the input data

p.value

the p-value

statistic

the z quantile of the standard normal distribution

null.value

the null hypothesis

conf.int

upper and lower confidence limit

alternative

the alternative hypothesis

method

character string that denotes the test

Note

Current Version is for complete observations only.

References

Hipel, K.W. and McLeod, A.I. (1994), Time Series Modelling of Water Resources and Environmental Systems. New York: Elsevier Science.

Sen, P.K. (1968), Estimates of the regression coefficient based on Kendall's tau, Journal of the American Statistical Association 63, 1379–1389.

Examples

data(maxau)
sens.slope(maxau[,"s"])
mk.test(maxau[,"s"])


Seasonal Mann-Kendall Trend Test

Description

Performs a Seasonal Mann-Kendall Trend Test (Hirsch-Slack Test)

Usage

smk.test(x, alternative = c("two.sided", "greater", "less"), continuity = TRUE)

Arguments

x

a time series object with class ts comprising >= 2 seasons; NA values are not allowed

alternative

the alternative hypothesis, defaults to two.sided

continuity

logical, indicates, whether a continuity correction should be done; defaults to TRUE

Details

The Mann-Kendall statistic for the $g$-th season is calculated as:

S_g = \sum_{i = 1}^{n-1} \sum_{j = i + 1}^n \mathrm{sgn}\left(x_{jg} - x_{ig}\right), \qquad (1 \le g \le m)

with \mathrm{sgn} the signum function (see sign).

The mean of S_g is \mu_g = 0. The variance including the correction term for ties is

\sigma_g^2 = \left\{n \left(n-1\right)\left(2n+5\right) - \sum_{j=1}^p t_{jg}\left(t_{jg} - 1\right)\left(2t_{jg}+5\right) \right\} / 18 ~~ (1 \le g \le m)

The seasonal Mann-Kendall statistic for the entire series is calculated according to

\begin{array}{ll} \hat{S} = \sum_{g = 1}^m S_g & \hat{\sigma}_g^2 = \sum_{g = 1}^m \sigma_g^2 \end{array}

The statistic S_g is approximately normally distributed, with

z_g = S_g / \sigma_g

If continuity = TRUE then a continuity correction will be employed:

z = \mathrm{sgn}(S_g) ~ \left(|S_g| - 1\right) / \sigma_g

Value

An object with class "htest" and "smktest"

data.name

character string that denotes the input data

p.value

the p-value for the entire series

statistic

the z quantile of the standard normal distribution for the entire series

null.value

the null hypothesis

estimates

the estimates S and varS for the entire series

alternative

the alternative hypothesis

method

character string that denotes the test

Sg

numeric vector that contains S scores for each season

varSg

numeric vector that contains varS for each season

pvalg

numeric vector that contains p-values for each season

Zg

numeric vector that contains z-quantiles for each season

References

Hipel, K.W. and McLeod, A.I. (1994), Time Series Modelling of Water Resources and Environmental Systems. New York: Elsevier Science.

Libiseller, C. and Grimvall, A. (2002), Performance of partial Mann-Kendall tests for trend detection in the presence of covariates. Environmetrics 13, 71–84, doi:10.1002/env.507.

R. Hirsch, J. Slack, R. Smith (1982), Techniques of Trend Analysis for Monthly Water Quality Data, Water Resources Research 18, 107–121.

Examples

res <- smk.test(nottem)
## print method
res
## summary method
summary(res)


Standard Normal Homogeinity Test (SNHT) for Change-Point Detection

Description

Performes the Standard Normal Homogeinity Test (SNHT) for change-point detection of a normal variate.

Usage

snh.test(x, m = 20000)

Arguments

x

a vector of class "numeric" or a time series object of class "ts"

m

numeric, number of Monte-Carlo replicates, defaults to 20000

Details

Let X denote a normal random variate, then the following model with a single shift (change-point) can be proposed:

x_i = \left\{ \begin{array}{lcl} \mu + \epsilon_i, & \qquad & i = 1, \ldots, m \\ \mu + \Delta + \epsilon_i & \qquad & i = m + 1, \ldots, n \\ \end{array} \right.

with \epsilon \approx N(0,\sigma). The null hypothesis \Delta = 0 is tested against the alternative \Delta \ne 0.

The test statistic for the SNHT test is calculated as follows:

T_k = k z_1^2 + \left(n - k\right) z_2^2 \qquad (1 \le k < n)

where

\begin{array}{l l} z_1 = \frac{1}{k} \sum_{i=1}^k \frac{x_i - \bar{x}}{\sigma} & z_2 = \frac{1}{n-k} \sum_{i=k+1}^n \frac{x_i - \bar{x}}{\sigma}. \\ \end{array}

The critical value is:

T = \max T_k.

The p.value is estimated with a Monte Carlo simulation using m replicates.

Critical values based on m = 1,000,000 Monte Carlo simulations are tabulated for T by Khaliq and Ouarda (2007).

Value

A list with class "htest" and "cptest"

data.name

character string that denotes the input data

p.value

the p-value

statistic

the test statistic

null.value

the null hypothesis

estimates

the time of the probable change point

alternative

the alternative hypothesis

method

character string that denotes the test

data

numeric vector of Tk for plotting

Note

The current function is for complete observations only.

References

H. Alexandersson (1986), A homogeneity test applied to precipitation data, Journal of Climatology 6, 661–675.

M. N. Khaliq, T. B. M. J. Ouarda (2007), On the critical values of the standard normal homogeneity test (SNHT), International Journal of Climatology 27, 681–687.

G. Verstraeten, J. Poesen, G. Demaree, C. Salles (2006), Long-term (105 years) variability in rain erosivity as derived from 10-min rainfall depth data for Ukkel (Brussels, Belgium): Implications for assessing soil erosion rates. Journal of Geophysical Research 111, D22109.

See Also

efp sctest.efp

Examples

data(Nile)
(out <- snh.test(Nile))
plot(out)

data(PagesData) ; snh.test(PagesData)
 

Object summaries

Description

Generic function "summary" for objects of class smktest.

Usage

## S3 method for class 'smktest'
summary(object, ...)

Arguments

object

an object of class smktest

...

further arguments, currently ignored


Wallis and Moore Phase-Frequency Test

Description

Performes the non-parametric Wallis and Moore phase-frequency test for testing the H0-hypothesis, whether the series comprises random data, against the HA-Hypothesis, that the series is significantly different from randomness (two-sided test).

Usage

wm.test(x)

Arguments

x

a vector or a time series object of class "ts"

Details

The test statistic of the phase-frequency test for n > 30 is calculated as:

z = \frac{| h - \frac{2 n - 7}{3}|}{\sqrt{\frac{16 n - 29}{90}}}

where h denotes the number of phases, whereas the first and the last phase is not accounted. The z-statistic is normally distributed. For n \le 30 a continuity correction of -0.5 is included in the denominator.

Value

An object of class "htest"

method

a character string indicating the chosen test

data.name

a character string giving the name(s) of the data

statistic

the Wallis and Moore z-value

alternative

a character string describing the alternative hypothesis

p.value

the p-value for the test

Note

NA values are omitted. Many ties in the series will lead to reject H0 in the present test.

References

L. Sachs (1997), Angewandte Statistik. Berlin: Springer.

C.-D. Schoenwiese (1992), Praktische Statistik. Berlin: Gebr. Borntraeger.

W. A. Wallis and G. H. Moore (1941): A significance test for time series and other ordered observations. Tech. Rep. 1. National Bureau of Economic Research. New York.

See Also

mk.test

Examples

## Example from Schoenwiese (1992, p. 113)
## Number of frost days in April at Munich from 1957 to 1968
## z = -0.124, Accept H0
frost <- ts(data=c(9,12,4,3,0,4,2,1,4,2,9,7), start=1957)
wm.test(frost)

## Example from Sachs (1997, p. 486)
## z = 2.56, Reject H0 on a level of p < 0.05
x <- c(5,6,2,3,5,6,4,3,7,8,9,7,5,3,4,7,3,5,6,7,8,9)
wm.test(x)

wm.test(nottem)


Wald-Wolfowitz Test for Independence and Stationarity

Description

Performes the non-parametric Wald-Wolfowitz test for independence and stationarity.

Usage

ww.test(x)

Arguments

x

a vector or a time series object of class "ts"

Details

Let x_1, x_2, ..., x_n denote the sampled data, then the test statistic of the Wald-Wolfowitz test is calculated as:

R = \sum_{i=1}^{n-1} x_i x_{i+1} + x_1 x_n

The expected value of R is:

E(R) = \frac{s_1^2 - s_2}{n - 1}

The expected variance is:

V(R) = \frac{s_2^2 - s_4}{n - 1} - E(R)^2 + \frac{s_1^4 - 4 s_1^2 s_2 + 4 s_1 s_3 + s_2^2 - 2 s_4}{(n - 1) (n - 2)}

with:

s_t = \sum_{i=1}^{n} x_i^t, ~~ t = 1, 2, 3, 4

For n > 10 the test statistic is normally distributed, with:

z = \frac{R - E(R)}{\sqrt{V(R)}}

ww.test calculates p-values from the standard normal distribution for the two-sided case.

Value

An object of class "htest"

method

a character string indicating the chosen test

data.name

a character string giving the name(s) of the data

statistic

the Wald-Wolfowitz z-value

alternative

a character string describing the alternative hypothesis

p.value

the p-value for the test

Note

NA values are omitted.

References

R. K. Rai, A. Upadhyay, C. S. P. Ojha and L. M. Lye (2013), Statistical analysis of hydro-climatic variables. In: R. Y. Surampalli, T. C. Zhang, C. S. P. Ojha, B. R. Gurjar, R. D. Tyagi and C. M. Kao (ed. 2013), Climate change modelling, mitigation, and adaptation. Reston, VA: ASCE. doi = 10.1061/9780784412718.

A. Wald and J. Wolfowitz (1943), An exact test for randomness in the non-parametric case based on serial correlation. Annual Mathematical Statistics 14, 378–388.

WMO (2009), Guide to Hydrological Practices. Volume II, Management of Water Resources and Application of Hydrological Practices, WMO-No. 168.

Examples

ww.test(nottem)
ww.test(Nile)

set.seed(200)
x <- rnorm(100)
ww.test(x)