# R/Hartung.R In punitroots: Unit Roots in Panels of Time Series

```Hartung <- function(p, lambda=rep(1,length(p)), kappa=0.2, alpha=0.10)
#
# This function applies the modified inverse normal method for the combination of dependent p-values.
#
# Arguments:
# p:         vector of p-values.
# lambda:    vector of weights. It must be of the same length of p.
# kappa:     adjustment parameter. It can be either a positive value (0.2 is the default value) or "formula". If
#            k == "formula", then it is computed as in Hartung, p. 853.
# alpha:     level for the 1-alpha confidence interval for rho (0.10 is the default).
#
# Value:
# 	     a list of class ``htest'' containing
#	     statistic: the Ht test statistic
#	     parameter: the number of combined tests (p-values)
#	     p.value: the combined test p-value
#	     conf.int: the confidence interval for the estimated correlation
#	     estimate: the estimated correlation
#	     null.value: the specified hypothesized value under the null
#	     alternative: a character string describing the alternative hypothesis
#	     method: a character srting indicating the type of combination test (only Z-test is implemented)
#	     data.name: a character string giving the name of the data set
#
# Reference:
# Hartung, J. (1999): "A note on combining dependent tests of significance",
#                     Biometrical Journal, 41(7), 849--855.
#
# Author:        Claudio Lupi
# This version:  August 22, 2010.
#
{
t       <- qnorm(p)
n       <- length(p)
avt     <- sum(t)/n
q       <- sum((t - avt)^2)/(n-1)                          # Hartung, eqn. (2.2)
rhohat  <- 1 - q
rhostar <- max(-1/(n-1), rhohat)                           # Hartung, p. 851
if (kappa=="formula") kappa <- (1 + 1/(n-1) - rhostar)/10  # Hartung, p. 853
if (kappa=="formula2") kappa <- (1 + 1/(n-1) - rhostar)/5  # Hartung, p. 853

# Hartung inverse normal corrected. See eqn. (2.4)
Ht <- sum(lambda*t)/sqrt(sum(lambda^2)+((sum(lambda))^2-sum(lambda^2))*(rhostar+kappa*sqrt(2/(n-1))*(1-rhostar)))
lower <- 1 - (n-1)/qchisq(alpha/2, (n-1)) * q
upper <- 1 - (n-1)/qchisq((1-alpha/2), (n-1)) * q          # Hartung, eqn. (2.3)

output <- list(statistic=c("Ht"=Ht),
parameter=c("N"=n),
p.value=pnorm(Ht),
conf.int=c(lower, upper),
estimate=c("average estimated correlation"=as.vector(rhohat)),
null.value=("Ht"=0),
alternative="less",
method="modified inverse normal combination",
data.name=deparse(substitute(p)))

class(output) <- "htest"

return(output)
}
```

## Try the punitroots package in your browser

Any scripts or data that you put into this service are public.

punitroots documentation built on May 2, 2019, 5:16 p.m.