# edps 589
# Fall 2018
# c.j.anderson
#
# PIRLS example ignoring missing data and unequal proportional sampling
#   design.  For correct analysis see paper on web-site by Anderson, Kim
#   and Keller
#
#*****************************************************************
#  Ordinal Model
#          P(Y_ij <= daily)             
#          P(Y_ij <= weekly)            
#          P(Y_ij <= monthly)           
#          
#*****************************************************************/
library(MASS)
library(ordinal)

pirls <- read.csv("D:\\Dropbox\\edps 589\\8 Multicategory_logit\\pirls_data.csv",sep=",")

# Response variable 
table(pirls$schnet)

#### Proportional odds
pirls$schnet <- as.factor(pirls$schnet)
summary(p1.a <- polr(schnet ~ short + girl + screenT, data=pirls,Hess=TRUE))

summary(p1.b <- clm(schnet ~ short + girl + screenT, data=pirls))
# Seems to be violating proportional odds assumption...
nominal_test(p1.b)

# Show that it matters-- fit separate binary logistic regression to each
# cummulative oddss ----> see lecture notes

pirls$schnet1 <- ifelse(pirls$schnet=="1",1,0)
summary(one <- glm(schnet1 ~ short + girl + screenT, data=pirls, family=binomial))

pirls$schnet2 <- ifelse(pirls$schnet=="2" | pirls$schnet=="1",1,0)
summary(two <- glm(schnet2 ~ short + girl + screenT, data=pirls, family=binomial))

pirls$schnet3 <- ifelse(pirls$schnet=="4",0,1)
summary(three <- glm(schnet3 ~ short + girl + screenT, data=pirls, family=binomial))

###########################################################
# Cumulative logits but relax proportional odds assumption#
###########################################################

summary(p2 <- vglm(schnet ~ short + girl + screenT, data=pirls,cumulative(parallel=FALSE)))

###########################################################
# Multinomial (Baseline)                                  #
###########################################################

summary(m <- vglm(schnet ~ short + girl + screenT, data=pirls,family=multinomial))