How are BRR scales used in the survey package for R?

Question

Does anyone know how to use the BRR scales in the Lumley survey package to measure variance, if your dataset already has a BRR, does it weigh it?

I work with PISA data and they already include 80 BRR replicas in their dataset. How can I get as.svrepdesign to use them instead of trying to create my own? I tried the following and got the following error:

dstrat <- svydesign(id=~uniqueID,strata=~strataVar, weights=~studentWeight, data=data, nest=TRUE) dstrat <- as.svrepdesign(dstrat, type="BRR") Error in brrweights(design$strata[, 1], design$cluster[, 1], ..., fay.rho = fay.rho, : Can't split with odd numbers of PSUs in a stratum 

Any help would be greatly appreciated, thanks.

Answer 1

no need to use as.svrepdesign() if you already have a data frame with replicating weights :) you can create a replicated weighted design directly from your data frame.

let's say that you have data with a main weight column named mainwgt and 80 replicated weight columns called repwgt1 via repwgt80 , you can use this -

 yoursurvey <- svrepdesign( weights = ~mainwgt , repweights = "repwgt[0-9]+" , type = "BRR", data = yourdata , combined.weights = TRUE ) 

- Thus, you do not need to specify the exact column numbers. then you can run regular browse commands, for example -

 svymean( ~variable , design = yoursurvey ) 

if you need another example, here is some sample code and an explanatory blog post using the current population survey.

Answer 2

I did not use the PISA data, I used the svprepdesign method last year using a microsphere for public use from the American Community Survey (US Census Bureau), which also came with 80 repetitive weights. They claim that they use the Fay method for this particular survey, so here you can build the svyrep object using this data:

 pums_p.rep<-svrepdesign(variables=pums_p[,2:7], repweights=pums_p[8:87], weights=pums_p[,1],combined.weights=TRUE, type="Fay",rho=(1-1/sqrt(4)),scale=1,rscales=1) attach(pums_p.rep) #CROSS - TABS #unweighted xtabs(~ is5to17youth + withinAMILimit) table(is5to17youth + withinAMILimit) #weighted, mean income by sex by race for select age groups svyby(~PINCP,~RAC1P+SEX,subset( pums_p.rep,AGEP > 25 & AGEP <35),na.rm = TRUE,svymean,vartype="se","cv") 

To make this work, I found an article article from A. Damico: Damico, A. (2009). Transitioning to R: Replicating SAS, Stata, and SUDAAN Analysis Techniques in Health Policy Data. The R Journal, 1(2), 37–44. Damico, A. (2009). Transitioning to R: Replicating SAS, Stata, and SUDAAN Analysis Techniques in Health Policy Data. The R Journal, 1(2), 37–44.