6

私はRで、学生が特定のコースで受け取ったレターグレードである応答変数を使用して作業しています。応答は序数であり、私の意見では、論理的に比例しているように見えます。私の理解では、multinom() の代わりに polr() を使用する前に、比例していることをテストする必要があります。

データのコースの 1 つについて、次のように比例性を「テスト」しました。

M1 <- logLik(polrModel)  #'log Lik.' -1748.180691 (df=8)
M2 <- logLik(multinomModel)  #'log Lik.' -1734.775727 (df=20)
G <- -2*(M1$1 - M2$2)) #I used a block bracket here in the real code
# 26.8099283
pchisq(G,12,lower.tail = FALSE) #DF is #of predictors
#0.008228890393     #THIS P-VAL TELLS ME TO REJECT PROPORTIONAL

比例オッズの仮定をテストする 2 つ目の方法として、2 つの vglm モデルも実行family=cumulative(parallel =TRUE)しましたfamily=cumulative(parallel =FALSE)。次にpchisq()、モデルの逸脱度の差と残差の自由度の差でテストを実行しました。

これらの方法のどちらかが立派ですか?そうでない場合は、比例オッズの仮定を受け入れるか拒否するかを決定するための実際のコーディングを支援したいと思います!

上記の 2 つのテストに加えて、各予測変数に対する累積確率を個別にグラフ化しました。これらの線を平行にしたいということを読みました。私が理解していないのは、polr()出力が各独立変数(係数)の単一の勾配であり、次に、作業している累積確率に応じて特定の切片であるということです(例:P(Y <= A)、P( Y<=B) など)。では、各方程式の勾配係数がすべて同じである場合、どのようにして直線が平行にならないのでしょうか?

Chris Bilder の YouTube クラスで基本的な知識を習得しました。彼はここの 42 分で平行グラフについて語っています。

どんな助けでも大歓迎です!ありがとうございました!

4

1 に答える 1

2

あなたのアプローチは本質的に正しいです。次のコードは、Fox の「R および S-PLUS の応用回帰へのコンパニオン」から着想を得たものです。第 5 章: 一般化線形モデルの適合。155 ~ 189 ページ。コードを使用する際は、本の章を引用してください。この章には、グラフ化に関するセクションもあります。

library(car)
library(nnet)
library(xlsx)
library(MASS)
options(warn=1)
options(digits = 3)
#
Trial <- read.xlsx("Trial.xls", "Sheet 1")
# Set up an out file structure
sink("Testing_adequacy_of_Prop_odds.txt")
# Trial$Outcome is assessed on a six point scale 0-5
schtyp_M_M.f <- factor(Trial$Outcome, labels = c("M0", "M1", "M2", "M3", "M4", "M5"))
#
cat("Multinomial logistic regression \n")
# Assign takes on a value of 1 (Treatment) or 0 (Control) 
mod.multinom <-multinom(schtyp_M_M.f~Assign, data = Trial)
print(summary(mod.multinom, cor=F, Wald=T))
x1<-logLik(mod.multinom)
cat("Degrees of freedom Multinomial logistic regression \n")
print(df_of_multinom_model <- attributes(x1)$df)
cat("Proportional odds logistic regression\n")
mod.polr <- polr(schtyp_M_M.f ~ Assign, data=Trial)
print(summary(mod.polr))
x2<-logLik(mod.polr)
cat("Degrees of freedom Proportional Odds Logistic Regression \n")
print(df_of_polr_model <- attributes(x2)$df)

cat("Answering the question: Is proportional odds model assumption violated\n")
cat("P value for difference in AIC between POLR and Multinomial Logit model\n")
# abs since the values could be negative. That is negative difference of degrees of freedom would produce p=NaN
print(1-pchisq(abs(mod.polr$deviance-mod.multinom$deviance),   abs(df_of_multinom_model-df_of_polr_model)))
sink()
于 2016-08-29T01:43:24.677 に答える