Блог

Главная - Вопросы по программированию - Существует ли пакет r для коррекции для многократного тестирования в GLM, в частности, для относительного риска или разницы в рисках?

Существует ли пакет r для коррекции для многократного тестирования в GLM, в частности, для относительного риска или разницы в рисках?

#r #glm #confidence-interval

Вопрос:

В своем исследовании я рассматриваю половые различия в эффективности биологических препаратов (=тип лекарства) у пациентов с аксиальным спондилоартритом. Я рассматриваю 7 различных исходов (оценки активности заболевания) за 6, 12 и 24 месяца. Я использую GLM с биномиальным семейством для расчета относительного риска и разницы в рисках. Я хочу расширить доверительные интервалы, так как я хочу внести коррективы для многократного тестирования.

Как бы я это сделал?

Вот небольшой репрезентативный набор данных и код, который я запустил:

 dat <- structure(list(gender = c(0, 1, 0, 0, 1, 1, 1, 1, 0, 0, 1, 0, 
1, 0, 1, 1, 1, 0, 0, 1), age_std = structure(c(-0.467384378359854, 
-0.623520001679169, 0.703632796535007, -0.935791248317799, -1.16999468329677, 
1.32817528981227, -1.16999468329677, 0.235225926577062, -0.857723436658141, 
-0.389316566700197, -0.467384378359854, -0.31124875504054, -1.24806249495643, 
-1.09192687163711, -0.623520001679169, -1.40419811827574, -1.79453717657403, 
-1.71646936491437, -1.87260498823369, -1.09192687163711), .Dim = c(20L, 
1L)), bio_drug_start_year_centered = c(-6, -2, -1, -1, -1, -9, 
-1, -1, -1, -1, -1, -1, 3, 2, 0, 2, 0, 0, 1, 0), country = structure(c(1L, 
1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 1L, 1L, 1L, 
1L, 1L, 1L), .Label = c("CH", "CZ", "DK", "ES", "IS", "IT", "NL", 
"NO", "PT", "RO", "SE", "SF", "SI", "TR", "UK"), class = "factor"), 
    bio_drug_start_year = c(2007, 2011, 2012, 2012, 2012, 2004, 
    2012, 2012, 2012, 2012, 2012, 2012, 2016, 2015, 2013, 2015, 
    2013, 2013, 2014, 2013), asdas_crp_cii_6month = c(1, 1, 0, 
    1, 0, 0, 1, 0, 1, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 0)), na.action = structure(c(`114174` = 2755L, 
`116484` = 2770L, `231118` = 3050L), class = "omit"), row.names = c("463", 
"7729", "7756", "8306", "8324", "128", "8440", "8450", "8663", 
"8809", "8840", "8857", "9020", "9033", "9101", "9324", "9377", 
"9523", "9702", "9718"), class = "data.frame")
    
rr = rd = numeric()
for( i in 1: 1e3){
  logmod = glm( asdas_crp_cii_6month ~ gender   country   age_std   bio_drug_start_year_centered, 
                family = 'binomial', data = dat[sample(1:nrow(dat),nrow(dat),replace=T),])
  summary(logmod)
  dat.1 = dat.0 = dat
  dat.1$gender = 1
  dat.0$gender = 0
  p1 = predict(logmod, newdata = dat.1, 'response' )
  p0 = predict(logmod, newdata = dat.0, 'response' )
  rr[i] = mean(p1)/mean(p0)
  rd[i] = mean(p1)-mean(p0)
}
twosidep<-function(data, test = 0){
  p1<-sum(data>test)/length(data)
  p2<-sum(data<test)/length(data)
  p<-min(p1,p2)*2
  return(p)
}
m = rbind(c(mean(rr), quantile(rr, c(0.025, 0.975)), twosidep(rr, test = 1) ),
          c(mean(rd), quantile(rd, c(0.025, 0.975)), twosidep(rd) ))
rownames(m) = c('rr', 'rd')
colnames(m) = c('Estimate', 'Lower CI', 'Upper CI', 'p-value')
m

Вот вывод кода

      Estimate   Lower CI  Upper CI p-value
rr  0.6337688  0.1612215 1.1929914    0.25
rd -0.3050054 -0.5731631 0.1239957    0.25

Вопрос: Как я могу увеличить эти доверительные интервалы с помощью метода коррекции Бонферрони?

Комментарии:

1. Вместо того, чтобы пытаться найти функцию, вам сначала нужно определить, существует ли статистический метод, подходящий для ваших данных, чтобы выполнить то, что вы хотите. Сначала вам следует обратиться за статистической помощью в Перекрестную проверку . Как только вы узнаете, какой метод подходит, вам будет легче найти, как выполнить это действие в R. Вам нужно точно знать, как вы хотите исправить для многократного тестирования.

2. Я должен был добавить это в вопрос. Я хочу использовать коррекцию бонферрони.

Метки: , ,