[R] Can you get the DEFT from svyratio?

Chris Webb iknowchris at gmail.com
Wed Dec 7 23:03:38 CET 2016


To Dr. Lumley or anyone who may know the answer,

I am trying to obtain ratio estimates from Levy and Lemeshow's Sampling of
Populations 4th ed. page 281. The results in the book are from STATA.
According to the STATA output, the DEFT is 0.830749

I can recreate all of the results except for DEFT. For svytotal and svymean
I can use the option deff="replace" to obtain DEFT results (by taking the
square root), but I get an error when using this option with svyratio. The
problem can be my poor understanding of how to calculate DEFT, but perhaps
it's not implemented for svyratio?


R code to that fails:

library(survey)

# Creating the dataset
df_tbl_10_1 <-
  data.frame(
    center = c(1,1,1,2,2,2,3,3,3,4,4,4,5,5,5),
    nurse = c(rep(c(1,2,3),5)),
    seen = c(58,44,18,42,63,10,13,18,37,16,32,10,25,23,23),
    referred = c(5,6,6,3,19,2,12,6,30,5,14,4,17,9,14)
  )
df_tbl_10_2 <- df_tbl_10_1[c(2,3,4,6,10,11),]

# Defining the cluster sampling design
svy_tbl_10_2 <-
  svydesign(id=~center + nurse,
            data=df_tbl_10_2,
            fpc= ~M + Nbar)

# Ratio estimates
svyratio(~referred, ~seen, svy_tbl_10_2)
confint(svyratio(~referred, ~seen, svy_tbl_10_2), df=degf(svy_tbl_10_2))

# DEFF
deff(svyratio(~referred, ~seen, svy_tbl_10_2, deff=TRUE))

# DEFT (fails)
sqrt(deff(svyratio(~referred, ~seen, svy_tbl_10_2, deff="replace")))

Fail message:
Error in if (deff) deffs <- matrix(ncol = nd, nrow = nn) : argument is not
interpretable as logical


For other individuals, I have included code that will calculate DEFF and
DEFT for svytotal on page 280, This code doesn't fail.

svytotal(~referred, svy_tbl_10_2)
confint(svytotal(~referred, svy_tbl_10_2), df=degf(svy_tbl_10_2))
deff(svytotal(~referred, svy_tbl_10_2, deff=TRUE))
sqrt(deff(svytotal(~referred, svy_tbl_10_2, deff="replace")))



To recap: can you get the DEFT from svyratio?

Sincerely,
Chris Webb

	[[alternative HTML version deleted]]



More information about the R-help mailing list