I want to test the hypothesis that the coefficient for tobacco use by age is greater in women than men using the couples data (CR) from the Demographic and Health Surveys (DHS), which has a complex survey design. Tobacco use (yes/no) and age data are sampled from couples in a household. Each row therefore includes columns:

household_id, age_female, age_male, tobacco_use_female, tobacco_use_male, ...
     12,          25,        32,           0,                1, ...
     13, ...

The survey design would be specified:

design <- svydesign(ids = ~v021, strata = ~v023, weights = ~I(v005/1e6), nest = T, data=ALCR51)

However, to fit the logistic regression I want, I need to convert to long format:

household_id, age, sex, tobacco_use, ...
12,            25,  F,     0, ...
12,            32,  M,     1, ...
13, ...
13, ...

The regression model would be specified like so, with the interaction coefficient the term of interest:

m <- svyglm(tobacco_use ~ age * sex, family = binomial, design)

But by converting to long format, observations are now clustered within households, thus requiring clustered robust errors or something equivalent. I know how to use the sandwich and lmtest packages with a standard glm model, but not with svyglm. I'm a beginner using complex survey designs in general, and the survey package in particular. Any guidance appreciated!


1 Answer 1


svyglm always returns standard errors account for the clustering in the design - see the Notes section of ?svyglm.

  • $\begingroup$ Thanks. My worry is that by converting the data to long format, I've introduced household level clustering that is not accounted for in the design specification. $\endgroup$
    – Ed Hagen
    Commented Jul 25, 2023 at 15:16
  • 1
    $\begingroup$ @EdHagen you need to refit the long data as a new design object where ...strata = ~household + v023, ...nest=T. $\endgroup$
    – AdamO
    Commented Jul 25, 2023 at 18:50
  • $\begingroup$ oh right... I didn't see that. Yes what AdamO is saying. The id argument is "specifying cluster ids from largest level to smallest level, ~0 or ~1 is a formula for no clusters.". So it'll be ~ hhid (you don't need to specify the person id if they are the observation unit in this data frame). Then keep strata = ~ v023 $\endgroup$
    – Alex J
    Commented Jul 25, 2023 at 23:02
  • 1
    $\begingroup$ You don't need to specify both stages of the id, since you're already approximating the design as sampling with replacement, you can specify just the coarsest id, ~v021. It might be clearer to write ~v021+hhid but it won't affect the results $\endgroup$ Commented Jul 25, 2023 at 23:31
  • 1
    $\begingroup$ Yes, that's correct. The 'with-replacement' approximation requires the weights, the first-stage PSUs and the first-stage strata. You might have to fix up the weights, though: the ones you have presumably add to the number of households and you now have a variable number of repeats of each weight. $\endgroup$ Commented Jul 26, 2023 at 0:20

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.