svymean_reg: Robust Generalized Regression Predictor (GREG) of the Mean and Total

Description

Generalized regression estimator (GREG) predictor of the mean and total, and robust GREG \(M\)-estimator predictor

Usage

svytotal_reg(object, totals, N = NULL, type, k = NULL, check.names = TRUE,
             keep_object = TRUE, ...)
svymean_reg(object, totals, N = NULL, type, k = NULL, check.names = TRUE,
            keep_object = TRUE, N_unknown = FALSE, ...)

Value

Object of class svystat_rob

Arguments

object: an object of class [svyreg_rob], e.g., result of the Huber regression \(M\)-estimator svyreg_huberM.
totals: [numeric] vector of population totals of the auxiliary variables.
N: [numeric] population size (see also N_unknown.
type: [character] type of predictor; see Details Section.
k: [numeric] robustness tuning constant of the psi-function used in the bias-correction term of the GREG. The definition of k depends on the type of predictor and is discussed in the Details Section.
check.names: [logical] if TRUE, the names of auxiliary are checked against the names of the independent variables of the fitted model object (default: TRUE).
keep_object: [logical] if TRUE, object is returned as an additional slot of the return value (default: TRUE).
N_unknown: [logical] if TRUE, it is assumed that the population size is unknown; thus, it is estimated (default: FALSE).
...: additional arguments (currently not used).

Details

Package survey must be attached to the search path in order to use the functions (see library or require).

The (robust) GREG predictor of the population total or mean is computed in two steps.

Step 1: Fit the regression model associated with the GREG predictor by one of the functions svyreg, svyreg_huberM, svyreg_huberGM, svyreg_tukeyM or svyreg_tukeyGM. The fitted model is called object.
Step 2: Based on the fitted model obtained in the first step, we predict the population total and mean, respectively, by the predictors svytotal_reg and svymean_reg, where object is the fitted regression model.

The following GREG predictors are available:

GREG (not robust, k = NULL)

The following non-robust GREG predictors are available:

type = "projective" ignores the bias correction term of the GREG predictor; see Särndal and Wright (1984).
type = "ADU" is the "standard" GREG, which is an asymptotically design unbiased (ADU) predictor; see Särndal et al.(1992, Chapter 6).

If the fitted regression model (object) does include a regression intercept, the predictor types "projective" and "ADU" are identical because the bias correction of the GREG is zero by design.

Robust GREG

The following robust GREG predictors are available:

type = "huber" and type = "tukey" are, respectively, the robust GREG predictors with Huber and Tukey bisquare (biweight) psi-function. The tuning constant must satisfy 0 < k <= Inf. We can use the Huber-type GREG predictor although the model has been fitted by the regression estimator with Tukey psi-function (and vice versa).
type = "BR" is the bias-corrected robust GREG predictor of Beaumont and Rivest (2009), which is inspired by the bias-corrected robust predictor of Chambers (1986). The tuning constant must satisfy 0 < k <= Inf.
type = "lee" is the bias-corrected predictor of Lee (1991; 1992). Tthe tuning constant k must satisfy 0 <= k <= 1.
type = "duchesne" is the bias-corrected, calibration-type estimator/ predictor of Duchesne (1999). The tuning constant k must be specified as a vector k = c(a, b), where a and b are the tuning constants of Duchesne's modified Huber psi-function (default values: a = 9 and b = 0.25).

Auxiliary data

Two types of auxiliary variables are distinguished: (1) population size \(N\) and (2) population totals of the auxiliary variables used in the regression model (i.e., non-constant explanatory variables).

The option N_unknown = TRUE can be used in the predictor of the population mean if \(N\) is unknown.

The names of the entries of totals are checked against the names of the regression fit (object), unless we specify check.names = FALSE.

Utility functions

The return value is an object of class svystat_rob. Thus, the utility functions summary, coef, SE, vcov, residuals, fitted, and robweights are available.

References

Beaumont, J.-F. and Rivest, L.-P. (2009). Dealing with outliers in survey data. In: Sample Surveys: Theory, Methods and Inference ed. by Pfeffermann, D. and Rao, C. R. Volume 29A of Handbook of Statistics, Amsterdam: Elsevier, Chap. 11, 247--280. tools:::Rd_expr_doi("10.1016/S0169-7161(08)00011-4")

Chambers, R. (1986). Outlier Robust Finite Population Estimation. Journal of the American Statistical Association 81, 1063--1069. tools:::Rd_expr_doi("10.1080/01621459.1986.10478374")

Duchesne, P. (1999). Robust calibration estimators, Survey Methodology 25, 43--56.

Gwet, J.-P. and Rivest, L.-P. (1992). Outlier Resistant Alternatives to the Ratio Estimator. Journal of the American Statistical Association 87, 1174--1182. tools:::Rd_expr_doi("10.1080/01621459.1992.10476275")

Lee, H. (1991). Model-Based Estimators That Are Robust to Outliers, in Proceedings of the 1991 Annual Research Conference, Bureau of the Census, 178--202. Washington, DC, Department of Commerce.

Lee, H. (1995). Outliers in business surveys. In: Business survey methods ed. by Cox, B. G., Binder, D. A., Chinnappa, B. N., Christianson, A., Colledge, M. J. and Kott, P. S. New York: John Wiley and Sons, Chap. 26, 503--526. tools:::Rd_expr_doi("10.1002/9781118150504.ch26")

Särndal, C.-E., Swensson, B. and Wretman, J. (1992). Model Assisted Survey Sampling, New York: Springer.

Särndal, C.-E. and Wright, R. L. (1984). Cosmetic Form of Estimators in Survey Sampling. Scandinavian Journal of Statistics 11, 146--156.

Examples

Run this code

head(workplace)

library(survey)
# Survey design for stratified simple random sampling without replacement
dn <- if (packageVersion("survey") >= "4.2") {
        # survey design with pre-calibrated weights
        svydesign(ids = ~ID, strata = ~strat, fpc = ~fpc, weights = ~weight,
                  data = workplace, calibrate.formula = ~-1 + strat)
    } else {
        # legacy mode
        svydesign(ids = ~ID, strata = ~strat, fpc = ~fpc, weights = ~weight,
                  data = workplace)
    }

# Robust regression M-estimator with Huber psi-function
reg <- svyreg_huberM(payroll ~ employment, dn, k = 3)

# Summary of the regression M-estimate
summary(reg)

# Diagnostic plots of the regression M-estimate (e.g., standardized
# residuals against fitted values)
plot(reg, which = 1L)

# Plot of the robustness weights of the regression M-estimate against
# its residuals
plot(residuals(reg), robweights(reg))

# ADU (asymptotically design unbiased) estimator
m <- svytotal_reg(reg, totals = 1001233, 90840, type = "ADU")
m

# Robust GREG estimator of the mean; the population means of the auxiliary
# variables are from a register
m <- svymean_reg(reg, totals = 1001233, 90840, type = "huber", k = 20)
m

# Summary of the robust GREG estimate
summary(m)

# Extract estimate
coef(m)

# Extract estimated standard error
SE(m)

# Approximation of the estimated mean square error
mse(m)

Run the code above in your browser using DataLab