Transforming between estimates of mean and total (complex survey design)

Question

I want to calculate the relative fraction of schools classified as E(lementary), M(iddle), and H(igh) from a complex survey. I do this using the survey package in R:

library('survey')
data(api)
a = svydesign(id=~dnum, weights=~pw, data=apiclus1)
r1 = svymean(~stype,a)
r2 = svytotal(~stype, a)

My naive expectation was that the results of r1 and r2 are just transformations of each other. And to some extent, this is true:

#Are the relative fractions implied by svytotal equal to the results of svymean
all(r2/sum(r2) == r1) #TRUE

#can you recover svytotal results from svymean and N
all.equal(as.numeric(r1) * sum(1/a$prob) , as.numeric(r2)) #TRUE

However, the standard errors (extracted via SE(r1) and SE(r2)) do not seem to have the same transformation. My in-depth knowledge of survey statistics (and statistics more generally) is more conceptual rather than mathematical, so I'm having trouble reasoning out why the following condition is not (even close to) true:

all(SE(r1)*sum(1/a$prob) == SE(r2))

I've skimmed through Lumley's "Complex Surveys" and Lohr's "Sampling: Design and Analysis" and haven't been able to find the section that provides insight (or I didn't understand it) into why this might be the case. Most of what I have gathered is that (at least in simple cases) they are transformations of each other as a function of N (e.g. sum of the inverse selection probabilities (N-hat?)).

So in short: Is there a way to translate between an estimate of the mean (and standard error) and an estimate of the total (and standard error) or are they fundamentally different calculations?

Answer 1

The uncertainty in the mean is not just a rescaling of the total because there is uncertainty in the denominator as well. You get the mean by dividing the estimated total by the estimated population size, so there's a fairly complicated formula involving variances and covariances for the SE of the mean.

A simple way to reproduce the calculation is to use svyratio(), with the numerator being the variable whose total you are calculating and the denominator being the number 1 for each record.

> svyratio(~I(stype=="E"), ~I(!is.na(stype)), design=a)
Ratio estimator: svyratio.survey.design2(~I(stype == "E"), ~I(!is.na(stype)), 
    design = a)
Ratios=
                I(!is.na(stype))
I(stype == "E")        0.7868852
SEs=
                I(!is.na(stype))
I(stype == "E")       0.04680257
> svymean(~stype,a)
           mean     SE
stypeE 0.786885 0.0468
stypeH 0.076503 0.0271
stypeM 0.136612 0.0299

You might think that the population size would be known. That's often true, but not always. For example, in cluster sampling you would need to know the sizes of the clusters you didn't sample.

For that reason, software typically uses the estimated population size to rescale the total to the mean. Also, the standard error is often smaller that way, because the estimated total tends to be positively correlated with the estimated population size and some of the uncertainty cancels out.