HTvar: Variance of the Horvitz-Thompson estimator

Description Usage Arguments Details Examples

View source: R/HTvar.R

Description

Compute or estimate the variance of the Horvitz-Thompson total estimator by the Horvitz-Thompson or Sen-Yates-Grundy variance estimators.

Usage

1
HTvar(y, pikl, sample = TRUE, method = "HT")

Arguments

y

numeric vector representing the variable of interest

pikl

matrix of second-order (joint) inclusion probabilities; the diagonal must contain the first-order inclusion probabilities.

sample

logical value indicating whether sample or population values are provided. If sample=TRUE, the function returns a sample estimate of the variance, while if sample=FALSE, the Variance is computed over all population units. Default is TRUE.

method

string, indicating if the Horvitz-Thompson ("HT") or the Sen-Yates-Grundy ("SYG") estimator should be computed.

Details

The Horvitz-Thompson variance is defined as

∑_U ∑_U [ π(ij) - π(i)π(j) ] y(i) y(j) / [ π(i)π(j) ]

which is estimated by

∑_s ∑_s [ π(ij) - π(i)π(j) ] y(i) y(j) / [ π(i)π(j)π(ij) ]

The Sen-Yates-Grundy variance is obtained from the Horvitz-Thompson variance by conditioning on the sample size n, and is therefore only appliable to fixed size sampling designs:

∑__U∑_{j > i} [ π(i)π(j) - π(ij) ] [ y(i)/π(i) - y(j)/π(j) ]^2

Its estimator is

∑__U∑_{j > i} [ π(i)π(j) - π(ij) ] [ y(i)/π(i) - y(j)/π(j) ]^2 / π(ij)

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
### Generate population data ---
N <- 500; n <- 50

set.seed(0)
x <- rgamma(500, scale=10, shape=5)
y <- abs( 2*x + 3.7*sqrt(x) * rnorm(N) )

pik  <- n * x/sum(x)
pikl <- jip_approx(pik, method='Hajek')

### Dummy sample ---
s   <- sample(N, n)


### Compute Variance ---
HTvar(y=y, pikl=pikl, sample=FALSE, method="HT")
HTvar(y=y, pikl=pikl, sample=FALSE, method="SYG")


### Estimate Variance ---
#' HTvar(y=y[s], pikl=pikl[s,s], sample=TRUE, method="HT")
#' HTvar(y=y[s], pikl=pikl[s,s], sample=TRUE, method="SYG")

rhobis/jipApprox documentation built on Sept. 26, 2018, 5:20 p.m.