How to calculate medians svyquantiles using svyby from survey package?

Question

I am reading the book Complex Survey and trying to follow the book codes, perhaps because the book is a little old now the codes are not working as intended, I don't know for sure.

I am trying to calculate the median of subpopulations using the svyby() function from the survey package.

The data used can be found in here.

library(survey)
library(haven)
chis_adult <- read_dta("Adult.dta")

The following line was not in the book, but for whatever reason I couldn't set the design and setting the chis_adult as data.frame made the svrepdesign work, not sure if this is relevant and has anything to do with the error:

chis_adult <- as.data.frame(chis_adult)

chis <- svrepdesign(variables = chis_adult[,1:418],
                    repweights = chis_adult[,420:499],
                    weights = chis_adult[,419], combined.weights = T,
                    type= "other", scale = 1, rscales =1)

Now, I am trying to calculate the medians for the subpopulations using the following codes:

medians <- svyby(~bmi_p, ~srsex+racehpr, svyquantile, design = chis, 
                 quantiles = 0.5, covmat = T)

Which gives me the error:

Error in (function (cond)  : 
  error in evaluating the argument 'x' in selecting a method for function 't': Can't merge the outer name `ci` with a vector of length > 1.
Please supply a `.name_spec` specification.

I am trying to set options such as na.rm = T, na.rm.by = T, na.rm.all = T, none of them seems to work. Any idea what could be happening?

Answer 1

I don't have a definitive answer, but the error seems to be due to the labelling of the variables in the Stata file, which is causing issues in the svyby function.

Try reading the file again using foreign::read.dta

library(foreign)
chis_adult <- read.dta("Adult.dta")

chis <- svrepdesign(variables = chis_adult[,1:418],
                    repweights = chis_adult[,420:499],
                    weights = chis_adult[,419], combined.weights = T,
                    type= "other", scale = 1, rscales =1)

The following command now returns a different error:

svyby(~bmi_p, ~srsex+racehpr, 
      FUN=svyquantile, 
      design = chis, 
      quantiles = 0.5, 
      covmat = T)
#Error in attr(replicates, "scale") <- design$scale : 
#  attempt to set an attribute on NULL

Which I managed to prevent by setting covmat to FALSE (the default).

svyby(~bmi_p, ~srsex+racehpr, 
       FUN=svyquantile, 
       design = chis, 
       quantiles = 0.5)
                                       srsex                        racehpr bmi_p   se.bmi_p
MALE.LATINO                             MALE                         LATINO 27.44 0.12057562
FEMALE.LATINO                         FEMALE                         LATINO 26.26 0.17583981
MALE.PACIFIC ISLANDER                   MALE               PACIFIC ISLANDER 28.87 0.70448510
FEMALE.PACIFIC ISLANDER               FEMALE               PACIFIC ISLANDER 27.26 0.67903816
MALE.AMERICAN INDIAN/ALASKAN NATIVE     MALE AMERICAN INDIAN/ALASKAN NATIVE 28.37 0.61292679
FEMALE.AMERICAN INDIAN/ALASKAN NATIVE FEMALE AMERICAN INDIAN/ALASKAN NATIVE 25.09 0.49486275
MALE.ASIAN                              MALE                          ASIAN 24.37 0.13564793
FEMALE.ASIAN                          FEMALE                          ASIAN 22.14 0.11806357
MALE.AFRICAN AMERICAN                   MALE               AFRICAN AMERICAN 27.41 0.37428760
FEMALE.AFRICAN AMERICAN               FEMALE               AFRICAN AMERICAN 27.20 0.22356736
MALE.WHITE                              MALE                          WHITE 26.39 0.06028781
FEMALE.WHITE                          FEMALE                          WHITE 24.28 0.04521598
MALE.OTHER SINGLE/MULTIPLE RACE         MALE     OTHER SINGLE/MULTIPLE RACE 26.63 0.40443079
FEMALE.OTHER SINGLE/MULTIPLE RACE     FEMALE     OTHER SINGLE/MULTIPLE RACE 25.69 0.36675145