Как эффективно рассчитать взвешенную пропорцию и доверительный интервал?

Question

Мне интересно, существует ли эффективный способ расчета взвешенной пропорции, а также 95% -й ДИ для нескольких переменных, сгруппированных по другим переменным? Это мои данные выборки

      sex          hpvvac  wtmec4yr ohpv06 ohpv11 ohpv16 ohpv18 ohpv26
1    Male            <NA> 67814.750      0      0      0      0      0
2    Male No HPV vaccined 12641.213      0      0      0      0      0
3  Female No HPV vaccined 51039.316      0      0      0      0      0
4    Male    HPV vaccined 19676.654      0      0      0      0      0
5  Female No HPV vaccined 11778.582      0      0      0      0      0
6    Male No HPV vaccined  9124.663      0      0      0      0      0
7    Male No HPV vaccined 10034.331      0      0      1      1      0
8    Male No HPV vaccined 17836.982      0      0      1      0      0
9    Male No HPV vaccined 48500.992      0      0      0      0      0
10 Female No HPV vaccined 19340.266      0      0      0      0      0

structure(list(sex = structure(c(1L, 1L, 2L, 1L, 2L, 1L, 1L, 
1L, 1L, 2L), .Label = c("Male", "Female"), class = "factor"), 
    hpvvac = structure(c(NA, 1L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 
    1L), .Label = c("No HPV vaccined", "HPV vaccined"), class = "factor"), 
    wtmec4yr = c(67814.75, 12641.212890625, 51039.31640625, 19676.654296875, 
    11778.58203125, 9124.6630859375, 10034.3310546875, 17836.982421875, 
    48500.9921875, 19340.265625), ohpv06 = c(0, 0, 0, 0, 0, 0, 
    0, 0, 0, 0), ohpv11 = c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0), ohpv16 = c(0, 
    0, 0, 0, 0, 0, 1, 1, 0, 0), ohpv18 = c(0, 0, 0, 0, 0, 0, 
    1, 0, 0, 0), ohpv26 = c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0)), row.names = c(NA, 
10L), class = "data.frame")

wtmec4yr - мой вес для каждого наблюдения. ohpv* - двоичные переменные. Я хотел бы рассчитать долю значения 1 из группы ohpv* на sex и hpvvac, а также 95% -й ДИ. Я пытался сделать это вручную, но сомневаюсь, что это было неправильно, и это было неэффективно. Я также попробовал survey пакет, но я могу рассчитывать только для одной переменной каждый раз.

d.s <- svydesign(ids=~1, data=mydt, weights =~wtmec4yr)
a <- svyby(~ohpv06, ~hpvvac+sex,d.s,svymean, na.rm=F)
ftable(a)
confint(a)

Спасибо!

Answer 1

Вы можете указать формулу:

library(survey)
svyby(~ohpv06+ohpv11+ohpv16+ohpv18+ohpv26, ~hpvvac+sex,d.s,svymean, na.rm=F)

или если все ваши столбцы начинаются с ohpv:

FORM = paste("~",paste(grep("ohpv",colnames(mydt),value=TRUE),collapse=" + "))
svyby(as.formula(FORM), ~hpvvac+sex,d.s,svymean, na.rm=F)

Оба они дают вам следующее:

                                hpvvac    sex ohpv06 ohpv11
No HPV vaccined.Male   No HPV vaccined   Male      0      0
HPV vaccined.Male         HPV vaccined   Male      0      0
No HPV vaccined.Female No HPV vaccined Female      0      0
                          ohpv16   ohpv18 ohpv26 se.ohpv06 se.ohpv11
No HPV vaccined.Male   0.2840007 0.102247      0         0         0
HPV vaccined.Male      0.0000000 0.000000      0         0         0
No HPV vaccined.Female 0.0000000 0.000000      0         0         0
                       se.ohpv16 se.ohpv18 se.ohpv26
No HPV vaccined.Male   0.2211842  0.113473         0
HPV vaccined.Male      0.0000000  0.000000         0
No HPV vaccined.Female 0.0000000  0.000000         0