Как понять вывод из функции R Polr (упорядоченная логистическая регрессия)?

Я новичок в R, заказал логистическую регрессию, и polr.

В разделе «Примеры» в нижней части страницы справки для polr (которая соответствует логистической или пробит-регрессионной модели для упорядоченного факторного отклика) показано

options(contrasts = c("contr.treatment", "contr.poly"))
house.plr <- polr(Sat ~ Infl + Type + Cont, weights = Freq, data = housing)
pr <- profile(house.plr)
plot(pr)
pairs(pr)

Какая информация prсодержит? Страница справки в профиле является общей и не дает указаний для polr.
Что plot(pr)показывает? Я вижу шесть графиков. Каждый из них имеет ось X, которая является числовой, хотя метка является переменной индикатора (выглядит как входная переменная, которая является индикатором для порядкового значения). Тогда ось Y представляет собой «тау», что совершенно необъяснимо.
Что pairs(pr)показывает? Это выглядит как график для каждой пары входных переменных, но опять же я не вижу объяснения осям X или Y.
Как понять, хорошо ли подойдет модель? summary(house.plr)показывает остаточное отклонение 3479,149 и AIC (информационный критерий Акаике?) 3495,149. Это хорошо? В том случае, если они полезны только в качестве относительных показателей (то есть для сравнения с другой моделью), что является хорошим абсолютным показателем? Распределение остаточного отклонения приблизительно по критерию хи-квадрат? Можно ли использовать «% правильно предсказанный» на исходных данных или какую-либо перекрестную проверку? Какой самый простой способ сделать это?
Как можно применять и интерпретировать anovaэту модель? В документах говорится: «Существуют методы для стандартных функций подбора моделей, включая предикт, сводку, vcov, anova». Тем не менее, anova(house.plr)результаты работы вanova is not implemented for a single "polr" object
Как интерпретировать значения t для каждого коэффициента? В отличие от некоторых моделей, здесь нет значений P.

Я понимаю, что это много вопросов, но для меня имеет смысл задавать один пакет («как я могу использовать эту вещь?»), А не 7 разных вопросов. Любая информация приветствуется.

r logistic dfrankow
источник

@dfrankow Довольно грубая и, конечно, очень частичная помощь для ваших первых двух вопросов, но methods("profile")даст вам (в данном случае S3) методы, связанные с profileобъектом R , тогда вы увидите, что есть специальный метод для polrрезультатов, который вы можете просмотреть в режиме онлайн, набрав getAnywhere("profile.polr")в приглашении R.

хл

Благодарность! Исходный код хорош. Объяснение было бы еще лучше. :)

dfrankow

Кто-то указал мне на «Современную прикладную статистику с S» от Venables и Ripley. В разделе 7.3 приведен «Пример таблицы частот в четырех направлениях», который широко охватывает эту модель дома. Чтение ..

dfrankow

На самом деле раздел является «моделью пропорциональных шансов»

dfrankow

Ответы:

Я бы посоветовал вам взглянуть на книги по анализу категориальных данных (см. Alan Agresti, Категориальный анализ данных, 2002), чтобы лучше объяснить и понять упорядоченную логистическую регрессию . На все вопросы, которые вы задаете, в основном отвечают несколько глав в таких книгах. Если вы заинтересованы только в Rродственных примерах Расширение линейных моделей в R по юлианскому Faraway (CRC Press, 2008) является отличным справочником.

Прежде чем я отвечу на ваши вопросы, упорядоченная логистическая регрессия - это случай полиномиальных логит-моделей, в которых упорядочены категории. Предположим, что у нас есть упорядоченные категории и для индивидуального с порядковым ответом , $J$ $i$ $Y_i$ для . С упорядоченным ответом часто легче работать с кумулятивными вероятностями, $p_{ij}=P(Yi=j)$ $j=1,..., J$ . Кумулятивные вероятности растут и инвариантны к объединению смежных категорий. Кроме того, , поэтому нам нужны тольковероятностимодели . $\gamma_{ij}=P(Y_i \le j)$ $\gamma_{iJ}=1$ $J–1$

Теперь мы хотим связать s с ковариатами . В вашем случае, имеет 3 упорядоченные уровни: , , . Более разумно обращаться с ними как с заказанными, а не неупорядоченными. Остальные переменные - ваши ковариаты. Конкретная модель, которую вы рассматриваете, является моделью пропорциональных шансов и математически эквивалентна: $\gamma_{ij}$ $x$ Satlowmediumhigh

logit γ_{j} (x_{i}) = θ_{j} - β^{T} x_{i}, j = 1 \dots J - 1

$\mbox{logit } \gamma_j(x_i) = \theta_j - \beta^T x_i, j = 1 \ldots J-1$

where γ_{j} (x_{i}) = P (Y_{i} \leq j | x_{i})

$\mbox{where }\gamma_j(x_i)=P(Y_i \le j | x_i)$

Это так называется, потому что относительные шансы для сравнении и : $Y \le j$ $x_1$ $x_2$

(\frac{γ_{j} (x_{1})}{1 - γ_{j} (x_{1})}) / (\frac{γ_{j} (x_{2})}{1 - γ_{j} (x_{2})}) = \exp (- β^{T} (x_{1} - x_{2}))

$\left(\frac {\gamma_j(x_1)}{1-\gamma_j(x_1)}\right) / \left(\frac {\gamma_j(x_2)}{1-\gamma_j(x_2)}\right)=\exp(-\beta^T (x_1-x_2))$

Обратите внимание, что приведенное выше выражение не зависит от . Конечно, предположение о пропорциональных шансах действительно необходимо проверить для данного набора данных. $j$

Теперь я отвечу на некоторые (1, 2, 4) вопросы.

Как понять, хорошо ли подошла модель? Summary (house.plr) показывает остаточное отклонение 3479,149 и AIC (информационный критерий Акаике?) 3495,149. Это хорошо? В том случае, если они полезны только в качестве относительных показателей (то есть для сравнения с другой моделью), что является хорошим абсолютным показателем? Распределение остаточного отклонения приблизительно по критерию хи-квадрат? Можно ли использовать «% правильно предсказанный» на исходных данных или какую-либо перекрестную проверку? Какой самый простой способ сделать это?

Подходящая модель polrявляется особенной glm, поэтому все допущения, которые имеют место для традиционного glmудержания здесь. Если вы позаботитесь о параметрах правильно, вы можете выяснить распределение. В частности, чтобы проверить, является ли модель хорошей или нет, вы можете захотеть выполнить тест на пригодность , который проверяет следующее значение null (обратите внимание, что это неуловимо, в основном вы хотите отклонить значение null, но здесь вы не хотите отклоните это, чтобы получить хорошую подгонку):

H_{o} : current model is good enough

$H_o: \mbox{ current model is good enough }$

Вы бы использовали для этого критерий хи-квадрат . Значение р получается как:

1-pchisq(deviance(house.plr),df.residual(house.plr))

В большинстве случаев вы надеетесь получить значение p больше 0,05, чтобы не отклонять ноль, чтобы сделать вывод о том, что модель подходит (философская корректность здесь игнорируется).

AIC должен быть высоким для хорошей подгонки, в то время как вы не хотите иметь большое количество параметров. stepAICхороший способ проверить это.

Да, вы определенно можете использовать перекрестную проверку, чтобы убедиться в правильности прогнозов. Смотрите predictфункцию (опция:) type = "probs"в ?polr. Все, что вам нужно, это ковариаты.

Какую информацию содержит PR? Страница справки в профиле является общей и не дает указаний для опроса

Как указано @chl и др., prСодержит всю информацию, необходимую для получения КИ и другую информацию, связанную с вероятностью polr fit. Все glms подгоняются с использованием итеративно взвешенного метода наименьших квадратов для логарифмической вероятности. В этой оптимизации вы получаете много информации (см. Ссылки), которая понадобится для расчета матрицы дисперсии ковариации, CI, t-значения и т. Д. Она включает в себя все это.

Как интерпретировать значения t для каждого коэффициента? В отличие от некоторых моделей> подходит, здесь нет значений P.

В отличие от нормальной линейной (специальной glm) модели, другие glms не имеют хорошего t-распределения для коэффициентов регрессии. Поэтому все, что вы можете получить, - это оценки параметров и их ковариационную матрицу асимптотической дисперсии с использованием теории максимального правдоподобия. Следовательно:

Variance (\hat{β}) = (X^{T} W X)^{- 1} \hat{ϕ}

$\text{Variance}(\hat \beta) = (X^T W X)^{-1}\hat \phi$

Оценка, разделенная на стандартную ошибку, - это то, что BDR и WV называют t-значением (здесь я предполагаю MASSсоглашение). Это эквивалентно t-значению из нормальной линейной регрессии, но не следует t-распределению. Используя CLT, он асимптотически нормально распределен. Но они предпочитают не использовать этот прибл (я думаю), следовательно, нет р-значений. (Я надеюсь, что я не ошибаюсь, и если я ошибаюсь, я надеюсь, что BDR нет на этом форуме. Я также надеюсь, что кто-то исправит меня, если я ошибаюсь.)

suncoolsu
источник

Я добавлю больше.

Suncoolsu

Спасибо за это. Я прочитал это несколько раз. Осталось много вопросов. 1. Как функционально в R можно проверить предположение о пропорциональных шансах? 2. Вы уверены, что тест хи-квадрат подходит? В этом примере он возвращает 0, что означает .. дерьмовое соответствие? Но некоторые из значений t довольно высоки (InflHigh 10.1, InflMedium 5.4, ContHigh 3.7). 3. Что показывают графики или пары?

dfrankow

Спасибо за ваш обширный ответ suncoolsu. Я в похожей ситуации и у меня есть пара вопросов. 1. Я также получаю 0 для каждой модели, используя ваше тестовое уравнение хи-кв. 2. На странице википедии на AIC сказано, что «предпочтительная модель - это модель с минимальным значением AIC», но вы сказали: «AIC должен быть высоким для хорошей подгонки». Я пытаюсь согласовать эти счета.

Сэм Свифт

@dfrankow и @Sam Swift. Извините, я был немного занят написанием статей. Хорошо - если вы получите значение p = 0, это означает, что модель НЕ подходит, так как проверка на соответствие не проходит. Что касается проблемы AIC, мы с Википедией используем для этого другое соглашение. Я использую тот, который используется BDR и WV. (ср. Расширение линейных моделей в R, д-р Джулиан

Фарауэй

Есть несколько специальных вопросов для значений 0/1 p и интерпретации AIC, которые могут оказаться полезными: stats.stackexchange.com/questions/15223/… stats.stackexchange.com/questions/81427/…

Скотт,

Мне очень понравился разговор здесь, однако я чувствую, что ответы не правильно адресовали все (очень хорошие) компоненты к поставленному вами вопросу. Вторая половина страницы примера для polr- это все о профилировании. Хороший технический справочник здесь - «Венераблс» и Рипли, которые обсуждают профилирование и что он делает. Это критический метод, когда вы выходите за пределы зоны комфорта, когда модели с экспоненциальным семейством соответствуют полной вероятности (обычные GLM).

$k-1$ $k$ lmernls, polrи glm.nb.

Страница справки ?profile.glmдолжна быть polrполезной, поскольку объекты по сути являются GLM (плюс категориальные пороги). Наконец, вы можете использовать пиковый исходный код, если он вам нужен getS3method('profile', 'polr'). Я getS3methodчасто использую эту функцию, потому что, хотя R, по-видимому, настаивает на том, что многие методы должны быть скрыты, можно удивительно много узнать о реализации и методах, просмотрев код.

• Какую информацию содержит PR? Страница справки в профиле является общей и не дает указаний для polr.

prэто profile.polr, profileобъект (унаследованный класс profile). Там есть запись для каждого ковариата. Профилировщик зацикливается на каждом ковариате, пересчитывает оптимальное соответствие модели с этим ковариатом, фиксированным до некоторой немного другой величины. Выходные данные показывают фиксированное значение ковариаты, измеренное как масштабированное различие «z-показателя» от его оценочного значения и результирующих фиксированных эффектов в других ковариатах. Например, если вы посмотрите pr$InflMedium, вы заметите, что, когда «z» равно 0, другие фиксированные эффекты такие же, как и в исходной подгонке.

• Что показывает сюжет (пр)? Я вижу шесть графиков. Каждый из них имеет ось X, которая является числовой, хотя метка является переменной индикатора (выглядит как входная переменная, которая является индикатором для порядкового значения). Тогда ось Y представляет собой «тау», что совершенно необъяснимо.

Опять ?plot.profileдает описание. График грубо показывает, как коэффициенты регрессии ковариально. тау - это масштабированная разница, а z-значение раньше, поэтому значение 0 дает оптимальные коэффициенты подгонки, обозначенные галочкой. Вы бы не сказали, что эта подгонка так хорошо себя ведет, но эти «линии» на самом деле являются сплайнами. Если бы вероятность очень нерегулярно вела себя при оптимальной подгонке, вы бы увидели странное и непредсказуемое поведение в сюжете. Это позволит вам оценить выход, используя более надежную оценку ошибок (начальная загрузка / складной нож), рассчитать CI, использовать method='profile'для перекодирования переменных или выполнить другую диагностику.

• Что показывает пара (pr)? Это выглядит как график для каждой пары входных переменных, но опять же я не вижу объяснения осям X или Y.

В файле справки говорится: «Метод пар показывает для каждой пары параметров x и y две кривые, пересекающиеся при оценке максимального правдоподобия, которые дают локусы точек, в которых касательные к контурам вероятности двумерного профиля становятся вертикальными и горизонтальный, соответственно. В случае точно двумерной вероятности нормального профиля эти две кривые были бы прямыми линиями, дающими условные средние y | x и x | y, и контуры были бы точно эллиптическими. " По сути, они снова помогают визуализировать эллипсы доверия. Неортогональные оси указывают на очень коваризуемые меры, такие как InfMedium и InfHigh, которые интуитивно очень связаны. Опять же, нерегулярные вероятности могут привести к тому, что образы здесь довольно непонятные.

• Как понять, хорошо ли подойдет модель? Summary (house.plr) показывает остаточное отклонение 3479,149 и AIC (информационный критерий Акаике?) 3495,149. Это хорошо? В том случае, если они полезны только в качестве относительных показателей (то есть для сравнения с другой моделью), что является хорошим абсолютным показателем? Распределение остаточного отклонения приблизительно по критерию хи-квадрат? Можно ли использовать «% правильно предсказанный» на исходных данных или какую-либо перекрестную проверку? Какой самый простой способ сделать это?

Одно предположение, которое хорошо оценить, это предположение о пропорциональных шансах. Это в некоторой степени отражается в глобальном тесте (который оценивает polr по насыщенной логлинеарной модели). Ограничение здесь состоит в том, что с большими данными глобальные тесты всегда терпят неудачу. В результате, использование графики и проверка оценок (бета-версий) и точности (SE) для логлинейной модели и подгонки polr является хорошей идеей. Если они массово не согласны, возможно, что-то не так.

При упорядоченных результатах трудно определить процентное согласие. Как вы выберете классификатор, основанный на модели, и если вы это сделаете, то как вы будете считать плохую производительность из плохого классификатора. modeэто плохой выбор. Если у меня есть 10 категорийных логитов и мой прогноз всегда только на одну категорию, возможно, это не плохо. Кроме того, моя модель может правильно предсказать 40% вероятности ответа 0, а также 20% вероятности 8, 9, 10. Итак, если я наблюдаю 9, это хорошо или плохо? Если вы должны измерить соглашение, используйте взвешенную каппу или даже MSE. Логлинейная модель всегда будет давать наилучшее согласие. Это не то, что делает POLR.

• Как применить и интерпретировать анову на этой модели? В документах говорится: «Существуют методы для стандартных функций подбора моделей, включая предикт, сводку, vcov, anova». Однако выполнение anova (house.plr) приводит к тому, что anova не реализована для одного объекта «polr»

Вы можете тестировать вложенные модели в пакете waldtestи lrtestв lmtestпакете R. Это эквивалентно ANOVA. Интерпретация точно такая же, как с GLM.

• Как интерпретировать значения t для каждого коэффициента? В отличие от некоторых моделей, здесь нет значений P.

Опять же, в отличие от линейных моделей, модель POLR может иметь проблемы с нерегулярной вероятностью, поэтому вывод, основанный на гессиане, может быть очень нестабильным. Это аналогично подгонке смешанных моделей, см., Например, файл справки confint.merModдля пакета lme4. Здесь оценки, сделанные с помощью профилирования, показывают, что ковариация ведется хорошо. Программисты сделали бы это по умолчанию, за исключением того, что профилирование может быть очень сложным в вычислительном отношении, и, таким образом, они оставляют это на ваше усмотрение. Если вы должны увидеть вывод на основе Wald, используйте coeftest(house.plr)из lrtestпакета.

Adamo
источник

Чтобы «проверить» (т. Е. Оценить) предположение о пропорциональных шансах в R, вы можете использовать residuals.lrm () в пакете дизайна Фрэнка Харрелла-младшего. Если вы наберете? Residuals.lrm, есть быстрый для повторения пример того, как Фрэнк Харрелл рекомендует оценивать предположение о пропорциональных коэффициентах (т. Е. Визуально, а не с помощью кнопочного теста). Проектируйте оценки упорядоченных логистических регрессий, используя lrm (), который вы можете заменить на polr () из MASS.

Более формальный пример того, как визуально проверить предположение о пропорциональных шансах в R, см. В статье: «Модели регрессионного анализа ординальных реакций в экологии» Автор (ы): Антуан Гисан и Фрэнк Э. Харрелл Источник: Journal of Vegetation Science, Vol. 11, № 5 (октябрь 2000 г.), с. 617-626

mBrewster
источник

Я искренне ценю ваш ответ. Однако целью StackExchange является предоставление ответов, а не ссылок. Статистики кажутся особенно склонными к этой проблеме ссылки. Можете ли вы добавить какие-либо подробности о том, как использовать residuals.lrm? Например, пример команды и пример интерпретации графика для примера house.plr?

Дфранков

Обновление с сайта автора: «Пакет Design теперь устарел. Пользователи R должны вместо этого использовать пакет rms». Марк, твой ответ был очень полезным для меня.

Тал Галили