Package 'trend'

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 $\mu$g/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"

• we first column, series of station Weil (RKM 164.3)

• ka second column, series of station Karlsruhe-Iffezheim (RKM 333.9)

• mz third column, series of station Mainz (RKM 498.5)

• ko fourth column, series of station Koblenz (RKM 590.3)

• bh fith column, series of station Bad Honnef(RKM 645.8)

• bi sixth column, series of station Bimmen (RKM 865.0)

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"

• s. first column, suspended sediment concentration

• Q. second column, average discharge

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)

---

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)

---

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)

---