R : 계층 적 데이터에 대한 베이지안 로지스틱 회귀

Question

이것은 만족스러운 응답을 얻지 못한 stats.stackexchange의 재 게시입니다. 나는 두 개의 데이터 세트를 가지고 있는데 첫 번째는 학교에서, 두 번째는 각 학교에있는 학생들 에서 표준화 테스트에 실패한 사람을 나열합니다 (의도적 인 강조). 가짜 데이터 세트 (Tharen 덕분에)에 의해 생성 될 수있다 : 나는 (실패 = 1 | 학생 인종, 학교 수익) P를 추정하기 위해 노력하고

#random school data for 30 schools 
schools.num = 30 
schools.data = data.frame(school_id=seq(1,schools.num) 
         ,tot_white=sample(100:300,schools.num,TRUE) 
         ,tot_black=sample(100:300,schools.num,TRUE) 
         ,tot_asian=sample(100:300,schools.num,TRUE) 
         ,school_rev=sample(4e6:6e6,schools.num,TRUE) 
         ) 

#total students in each school 
schools.data$tot_students = schools.data$tot_white + schools.data$tot_black + schools.data$tot_asian 
#sum of all students all schools 
tot_students = sum(schools.data$tot_white, schools.data$tot_black, schools.data$tot_asian) 
#generate some random failing students 
fail.num = as.integer(tot_students * 0.05) 

students = data.frame(student_id=sample(seq(1:tot_students), fail.num, FALSE) 
         ,school_id=sample(1:schools.num, fail.num, TRUE) 
         ,race=sample(c('white', 'black', 'asian'), fail.num, TRUE) 
        ) 

. 학생 데이터 세트에서 다항식 이산 선택 모델을 실행하면 분명히 P (Race | Fail = 1)를 추정하게됩니다. 분명히 이것의 역함을 추정해야합니다. 모든 정보가 두 개의 데이터 세트 (P (Fail), P (Race), Revenue)에서 사용 가능하므로이 작업을 수행 할 수없는 이유는 없습니다. 그러나 실제로 R에서 구현하는 방법에 관해서는 난처한 상황에 있습니다. 어떤 포인터라도 크게 감사 할 것입니다. 감사.

Answer 1

단일 data.frame이 있으면 더 쉬울 것입니다.

library(reshape2) 
library(plyr) 
d1 <- ddply(
    students, 
    c("school_id", "race"), 
    summarize, 
    fail=length(student_id) 
) 
d2 <- with(schools.data, data.frame( 
    school_id = school_id, 
    white = tot_white, 
    black = tot_black, 
    asian = tot_asian, 
    school_rev = school_rev 
)) 
d2 <- melt(d2, 
    id.vars=c("school_id", "school_rev"), 
    variable.name="race", 
    value.name="total" 
) 
d <- merge(d1, d2, by=c("school_id", "race")) 
d$pass <- d$total - d$fail 

그런 다음 데이터

library(lattice) 
xyplot(d$fail/d$total ~ school_rev | race, data=d) 

보고 또는 당신이 원하는 무엇이든을 계산할 수 있습니다.

r <- glm(
    cbind(fail,pass) ~ race + school_rev, 
    data=d, 
    family=binomial() # Logistic regression (not bayesian) 
) 
summary(r) 

(EDIT)는 전달 된 것들에 실패한 학생, 만 집계 된 데이터에 대한 자세한 정보가있는 경우 다음과 같이 당신은 완벽한 데이터 집합을 다시 만들 수 있습니다.

# Unique student_id for the passed students 
d3 <- ddply(d, 
    c("school_id", "race"), 
    summarize, student_id=1:pass 
) 
d3$student_id <- - seq_len(nrow(d3)) 
# All students 
d3$result <- "pass" 
students$result <- "fail" 
d3 <- merge(# rather than rbind, in case there are more columns 
    d3, students, 
    by=c("student_id", "school_id", "race", "result"), 
    all=TRUE 
) 
# Students and schools in a single data.frame 
d3 <- merge(d3, schools.data, by="school_id", all=TRUE) 
# Check that the results did not change 
r <- glm(
    (result=="fail") ~ race + school_rev, 
    data=d3, 
    family=binomial() 
) 
summary(r) 

Answer 2

모든 학생에 대한 정보가 포함 된 데이터 세트가 필요합니다. 둘 다 실패하고 통과했습니다.

schools.num = 30 
schools.data = data.frame(school_id=seq(1,schools.num) 
          ,tot_white=sample(100:300,schools.num,TRUE) 
          ,tot_black=sample(100:300,schools.num,TRUE) 
          ,tot_asian=sample(100:300,schools.num,TRUE) 
          ,school_rev=sample(4e6:6e6,schools.num,TRUE) 
         ) 

library(plyr) 
fail_ratio <- 0.05 
dataset <- ddply(schools.data, .(school_id, school_rev), function(x){ 
    data.frame(Fail = rbinom(sum(x$tot_white, x$tot_asian, x$tot_black), size = 1, prob = fail_ratio), Race = c(rep("white", x$tot_white), rep("asian", x$tot_asian), rep("black", x$tot_black))) 
}) 
dataset$Race <- factor(dataset$Race) 

그런 다음 자주 사용하는 방법으로 lme4 패키지에 glmer()를 사용할 수 있습니다.

library(lme4) 
glmer(Fail ~ school_rev + Race + (1|school_id), data = dataset, family = binomial) 

베이지안 추정치가 필요한 경우 MCMCglmm 패키지를 살펴보십시오.