[R] svyby and svyratio in the Thomas Lumley's survey package

[Thomas Lumley]

On Thu, 23 Feb 2006, Smith, Phil wrote:

> Dear R-ers:
>
> I'm using Thomas Lumley's "survey" package.
>
> I'd like to compute survey ratio estimates (numerator=~utd , 
> denominator=~one) for each of serval domains using by=~factor(domain).
>
> I can't quite work out the syntax for the call to the "svyby" function. I try:
>
> svyby( numerator=~ utd, denominator=~one , by=~ factor(domain) , 
> design=nis , svyratio )
>
> and I get an error that says that "domain" is not found. (nis is the 
> design object)
>

Hmm. I don't get the same error message (perhaps a question of versions?), 
but there is a problem. svyby is not really set up to work with svyratio. 
I will try to fix this, but there are simple workarounds.

If you do
  data(api)
  dstrat<-svydesign(id=~1,strata=~stype, weights=~pw,
        data=apistrat, fpc=~fpc)

  result <- svyby(~api.stu,by=~factor(stype), denominator=~enroll,
        design=dstrat,svyratio)

(ie, don't name the numerator argument) then the computations are done as 
you want.  The result does not print correctly, because svyby() returns 
some extra attributes that print() doesn't know how to handle.

You can extract the results you need with coef() or just remove the 
extraneous attributes
> result
Error in unlist(x, recursive, use.names) :
         argument not a list
> coef(result)
   statistic.ratio statistic.var
E       0.8518163  4.945408e-05
H       0.8105702  0.0004193180
M       0.8356958  0.0003307829
> result$statistic.call<-NULL
> result
   factor(stype) statistic.ratio statistic.var         SE
E             E       0.8518163  4.945408e-05 0.00703236
H             H       0.8105702  0.0004193180 0.02047726
M             M       0.8356958  0.0003307829 0.01818744



Finally, the name of the denominator variable looks suspicious. If `one' 
is a variable whose values are all 1 and you just want domain means you 
can do the equivalent of
   svyby(~api00,~stype,design=dstrat,svymean)
to get the mean of api00 within each value of stype.  You don't need 
svyratio(). In fact, part of the test suite for the survey package is to 
check that svyratio, svymean, and svyglm give the same answers for domain 
means.

 	-thomas