Квантильные оценки для подгрупп населения, где в некоторых подгруппах имеется только один случай с использованием пакетов srvyr и survey R - PullRequest
Новая
       50

Квантильные оценки для подгрупп населения, где в некоторых подгруппах имеется только один случай с использованием пакетов srvyr и survey R

1 голос
/ 06 апреля 2020

Я пытаюсь получить оценки 25-го процентиля непрерывной переменной для ряда подгрупп, где данные взяты из опроса, в котором используются веса выборки. Я делаю это в R с использованием пакетов survey и srvyr .

Эта проблема, с которой я сталкиваюсь, заключается в том, что в небольшом меньшинстве случаев подгруппа имеет только одну наблюдение и, следовательно, 25-й процентиль не имеет смысла. Это было бы хорошо, однако это приводит к ошибке, которая препятствует вычислению процентилей для тех подгрупп с достаточным количеством наблюдений. 

Error in approxfun(cum.w, xx[oo], method = method, f = f, yleft = min(xx),  : 
  need at least two non-NA values to interpolate

Код запускается, когда удаляющиеся группы удаляются, однако мне пришлось идентифицировать их вручную что далеко от идеала.

Есть ли способ достичь того же результата, но где для отдельных групп наблюдений выводится NA, или просто значение этого наблюдения, а не ошибка? В качестве альтернативы, есть ли удобный способ автоматического исключения таких групп из расчета?

Ниже приведен воспроизводимый пример, иллюстрирующий мою проблему с использованием набора данных apistrat из пакета survey . 

library(dplyr)
library(survey)
library(srvyr)
data(api)

#25th percentile of api00 by school type and whether school is year round  or not
apistrat %>% 
  as_survey(strata = stype, weights = pw) %>%
  group_by(yr.rnd, stype, .drop=TRUE) %>%
  summarise(survey_quantile(api00, 0.25, na.rm=T))

#Error in approxfun(cum.w, xx[oo], method = method, f = f, yleft = min(xx),  : 
#need at least two non-NA values to interpolate

apistrat %>% group_by(yr.rnd, stype) %>% tally() %>% filter(n==1)
#one group out of 6 has only a single api00 observation and therefore a quantile can't be interpolated

#Removing that one group means the code can now run as intended
apistrat %>% 
  as_survey(strata = stype, weights = pw) %>%
  filter(!(yr.rnd=="Yes"&stype=="H")) %>%
  group_by(yr.rnd, stype, .drop=TRUE) %>%
  summarise(survey_quantile(api00, 0.25, na.rm=T))

#Get the same error if you do it the 'survey' package way
dstrat <- svydesign(id=~1,strata=~stype,data=apistrat, fpc=~fpc)
svyby(~api99, ~stype+yr.rnd, dstrat, svyquantile, quantiles=0.25)

1 Ответ

1 голос
/ 09 апреля 2020

Один из обходных путей - заключить вызов в svyquantile(), используя tryCatch()

> svyq<-function( ...){tryCatch(svyquantile(...), error=function(e) matrix(NA,1,1))}
> svyby(~api99, ~stype+yr.rnd, dstrat, svyq, quantiles=0.25,keep.var=FALSE,na.rm=TRUE)
      stype yr.rnd statistic
E.No      E     No    560.50
H.No      H     No    532.75
M.No      M     No    509.00
E.Yes     E    Yes    456.00
H.Yes     H    Yes        NA
M.Yes     M    Yes    436.00

С квантилями и svyby, вам необходимо четко указать, хотите ли вы стандартных ошибок - код выше не Если вам нужны стандартные ошибки, вам потребуется ветвь error = tryCatch, чтобы вернуть фактический svyquantile объект с NA.

...