내 상황에 대한 약간의 배경은 다음과 같습니다. 내 데이터는 포식자가 성공적으로 먹은 먹이 수를 나타냅니다. 각 시험에서 먹이의 수가 제한되어 있으므로 (25 개의 이용 가능), 이용 가능한 먹이의 수를 나타내는 “샘플”열 (각 시험에서 25)과 성공 횟수 ( “Count”)가 있습니다 ( 먹이를 몇 마리 먹었는지) 비율 데이터 (578 페이지)를 기준으로 R 책의 예를 바탕으로 분석을 수행했습니다. 설명 변수는 온도 (내가 요인으로 취급 한 4 단계)와 포식자의 성별 (분명히 남성 또는 여성)입니다. 그래서 나는이 모델로 끝납니다 :
model <- glm(y ~ Temperature+Sex+Temperature*Sex data=predator, family=quasibinomial) 편차 분석 표를 얻은 후 온도와 성별 (상호 작용은 아님)이 먹이 소비에 중요한 영향을 미치는 것으로 나타났습니다. 이제 내 문제 : 나는 어떤 온도가 다른지 알아야합니다. 즉, 4 개의 온도를 서로 비교해야합니다. 선형 모델이 있다면 TukeyHSD 기능을 사용하지만 GLM을 사용함에 따라 사용할 수 없습니다. MASS 패키지를 살펴보고 대비 행렬을 설정하려고했지만 어떤 이유로 작동하지 않습니다. 어떤 제안이나 참조?
내 모델에서 얻은 요약은 다음과 같습니다.
y <- cbind(data$Count, data$Sample-data$Count)
model <- glm(y ~ Temperature+Sex+Temperature*Sex data=predator, family=quasibinomial)
> summary(model)
# Call:
# glm(formula = y ~ Temperature + Sex + Temperature * Sex, family=quasibinomial, data=data)
# Deviance Residuals:
# Min 1Q Median 3Q Max
# -3.7926 -1.4308 -0.3098 0.9438 3.6831
# Coefficients:
# Estimate Std. Error t value Pr(>|t|)
# (Intercept) -1.6094 0.2672 -6.024 3.86e-08 ***
# Temperature8 0.3438 0.3594 0.957 0.3414
# Temperature11 -1.0296 0.4803 -2.144 0.0348 *
# Temperature15 -1.2669 0.5174 -2.449 0.0163 *
# SexMale 0.3822 0.3577 1.069 0.2882
# Temperature8:SexMale -0.2152 0.4884 -0.441 0.6606
# Temperature11:SexMale 0.4136 0.6093 0.679 0.4990
# Temperature15:SexMale 0.4370 0.6503 0.672 0.5033
# ---
# Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# (Dispersion parameter for quasibinomial family taken to be 2.97372)
# Null deviance: 384.54 on 95 degrees of freedom
# Residual deviance: 289.45 on 88 degrees of freedom
# AIC: NA
# Number of Fisher Scoring iterations: 5답변
Anne, 나는 그러한 다중 비교를 일반적으로 수행하는 방법을 간단히 설명 할 것이다. 특정 상황에서 이것이 작동하지 않는 이유는 모르겠습니다. 죄송 해요.
그러나 일반적으로 multcomp패키지와 함수로 할 수 있습니다 glht. 예를 들면 다음과 같습니다.
mydata <- read.csv("http://www.ats.ucla.edu/stat/data/binary.csv")
mydata$rank <- factor(mydata$rank)
my.mod <- glm(admit~gre+gpa*rank, data=mydata, family=quasibinomial)
summary(my.mod)
#
# Coefficients:
# Estimate Std. Error t value Pr(>|t|)
# (Intercept) -4.985768 2.498395 -1.996 0.0467 *
# gre 0.002287 0.001110 2.060 0.0400 *
# gpa 1.089088 0.731319 1.489 0.1372
# rank2 0.503294 2.982966 0.169 0.8661
# rank3 0.450796 3.266665 0.138 0.8903
# rank4 -1.508472 4.202000 -0.359 0.7198
# gpa:rank2 -0.342951 0.864575 -0.397 0.6918
# gpa:rank3 -0.515245 0.935922 -0.551 0.5823
# gpa:rank4 -0.009246 1.220757 -0.008 0.9940
# ---
# Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1rankTukey의 HSD 사용 간의 쌍별 비교를 계산하려면 다음과 같이하십시오.
library(multcomp)
summary(glht(my.mod, mcp(rank="Tukey")))
#
# Simultaneous Tests for General Linear Hypotheses
#
# Multiple Comparisons of Means: Tukey Contrasts
#
# Fit: glm(formula = admit ~ gre + gpa * rank, family = quasibinomial, data = mydata)
#
# Linear Hypotheses:
# Estimate Std. Error z value Pr(>|z|)
# 2 - 1 == 0 0.5033 2.9830 0.169 0.998
# 3 - 1 == 0 0.4508 3.2667 0.138 0.999
# 4 - 1 == 0 -1.5085 4.2020 -0.359 0.984
# 3 - 2 == 0 -0.0525 2.6880 -0.020 1.000
# 4 - 2 == 0 -2.0118 3.7540 -0.536 0.949
# 4 - 3 == 0 -1.9593 3.9972 -0.490 0.960
# (Adjusted p values reported -- single-step method)
#
# Warning message:
# In mcp2matrix(model, linfct = linfct) :
# covariate interactions found -- default contrast might be inappropriate참고 : 주석에서 @gung이 언급했듯이 가능하면 온도를 범주 변수가 아닌 연속으로 포함해야합니다. 교호 작용에 대해 : 교호 작용 항이 모형 적합을 유의하게 개선하는지 여부를 확인하기 위해 우도 비율 검정을 수행 할 수 있습니다. 귀하의 경우 코드는 다음과 같습니다.
# Original model
model <- glm(y ~ Temperature+Sex+Temperature*Sex, data=predator, family=quasibinomial)
# Model without an interaction
model2 <- glm(y ~ Temperature+Sex data=predator, family=quasibinomial)
# Likelihood ratio test
anova(model, model2, test="LRT")이 테스트가 중요하지 않으면 모델에서 상호 작용을 제거 할 수 있습니다. 아마도 효과 glht가 있을까요?