Может кто-нибудь сказать мне, почему R не использует весь data.frame для этого chisq.test?
Я не могу найти решение проблемы, с которой столкнулся, пытаясь создать свою собственную data.frame
и запустить количественный анализ (такой как chisq.test
) в теме.
Фон выглядит следующим образом: я суммировал полученные данные, касающиеся двух больниц. Оба измеряли одну и ту же категориальную переменную n раз. В данном случае это то, как часто бактерии, связанные со здравоохранением, были обнаружены в течение определенного периода наблюдения.
В таблице обобщенные данные выглядят следующим образом, где% указывает процентную долю всех измерений, выполненных за период времени.
n Hospital 1 (%) n Hospital 2 (%)
Healthcare associated bacteria 829 (59.4) 578 (57.6)
Community associated bacteria 473 (33.9) 372 (37.1)
Contaminants 94 (6.7) 53 (5.3)
Total 1396 (100.0) 1003 (100.0)
Теперь, глядя на проценты, видно, что пропорции очень похожи, и вы можете удивиться, почему, черт возьми, я хочу статистически сравнить две больницы. Но у меня есть другие данные, где пропорции разные, и поэтому целью этого вопроса является:
Как сравнить Больницу 1 с Больницей 2 в отношении измеренных категорий.
Поскольку данные представлены в обобщенном виде и в формате массива, я решил сделать data.frame
для каждой из отдельных переменных / категорий.
hosp1 <- rep(c("Yes", "No"), times=c(829,567))
hosp2 <- rep(c("Yes", "No"), times=c(578,425))
all <- cbind(hosp1, c(hosp2,rep(NA, length(hosp1)-length(hosp2))))
all <- data.frame(all)
names(all)[2]<-"hosp2"
summary(all)
Пока все хорошо, потому что резюме, кажется, выглядит правильно, чтобы иметь возможность теперь запустить chisq.test()
, Но теперь все становится странным.
with(all, chisq.test(hosp1, hosp2, correct=F))
Pearson's Chi-squared test
data: hosp1 and hosp2
X-squared = 286.3087, df = 1, p-value < 2.2e-16
Результаты, кажется, указывают, что есть существенная разница. Если вы кросс-таблицы данных, вы видите, что R суммирует их очень странным образом:
with(all, table(hosp1, hosp2))
No Yes
No 174 0
Yes 251 578
Поэтому, конечно, если данные будут обобщены таким образом, будет получен статистически значимый результат - потому что одна категория суммируется как не имеющая измерений вообще. С какой стати это происходит и что я могу сделать, чтобы исправить это? Наконец, вместо того, чтобы сделать отдельный data.frame
есть ли в любом случае очевидная функция для ее зацикливания? Я не могу придумать один.
Спасибо за вашу помощь!
ОБНОВЛЕНО НА ОСНОВЕ ЗАПРОСА ТЕЛАТЕМАЙЛА ДЛЯ RAW DATA.FRAME
dput(SO_Example_v1)
structure(list(Type = structure(c(3L, 1L, 2L), .Label = c("Community",
"Contaminant", "Healthcare"), class = "factor"), hosp1_WoundAssocType = c(464L,
285L, 24L), hosp1_BloodAssocType = c(73L, 40L, 26L), hosp1_UrineAssocType = c(75L,
37L, 18L), hosp1_RespAssocType = c(137L, 77L, 2L), hosp1_CathAssocType = c(80L,
34L, 24L), hosp2_WoundAssocType = c(171L, 115L, 17L), hosp2_BloodAssocType = c(127L,
62L, 12L), hosp2_UrineAssocType = c(50L, 29L, 6L), hosp2_RespAssocType = c(135L,
142L, 6L), hosp2_CathAssocType = c(95L, 24L, 12L)), .Names = c("Type",
"hosp1_WoundAssocType", "hosp1_BloodAssocType", "hosp1_UrineAssocType",
"hosp1_RespAssocType", "hosp1_CathAssocType", "hosp2_WoundAssocType",
"hosp2_BloodAssocType", "hosp2_UrineAssocType", "hosp2_RespAssocType",
"hosp2_CathAssocType"), class = "data.frame", row.names = c(NA,
-3L))
Пояснение: это data.frame
на самом деле является более сложным, чем то, что суммировано в приведенной выше таблице, поскольку оно также содержит информацию о том, где конкретные виды бактерий были культивированы (например, в ранах, посевах крови, катетерах и т. д.). Таким образом, таблица, которую я делаю, выглядит следующим образом:
All locations
n Hospital 1 (%) n Hospital 2 (%) p-val
Healthcare associated bacteria 829 (59.4) 578 (57.6) 0.39
Community associated bacteria 473 (33.9) 372 (37.1) ...
Contaminants 94 (6.7) 53 (5.3) ...
Total 1396 (100.0) 1003 (100.0) -
Если заголовок "Все места" будет впоследствии заменен на рану, кровь, мочу, катетер и т. Д.
1 ответ
Ответ на вопрос о том, как заставить работать p-значения, несколько прост. Вы можете получить два других p-значения, которые вы ищете, используя тот же синтаксис, что и @thelatemail, следующим образом:
#community (p = 0.1049)
chisq.test(cbind(c(473,923),c(372,631)),correct=FALSE)
#contaminants (p = 0.1443)
chisq.test(cbind(c(94,1302),c(53,950)),correct=FALSE)
Вы можете получить эти ответы более программно следующим образом:
out <- cbind(rowSums(SO_Example_v1[,2:6]),rowSums(SO_Example_v1[,7:11]))
chisq.test(rbind(out[1,],colSums(out[2:3,])),correct=FALSE)
chisq.test(rbind(out[2,],colSums(out[c(1,3),])),correct=FALSE)
chisq.test(rbind(out[3,],colSums(out[1:2,])),correct=FALSE)
Конечно, на данный момент мы выходим за рамки SO, но, возможно, более уместный вопрос, учитывая характер данных, заключается в том, есть ли разница между больницами в целом, на которую вы можете ответить (с точки зрения частых пациентов), используя тест хи-квадрат на основе всех трех типов:
chisq.test(out,correct=FALSE)