Мультиномиальная регрессия (разные результаты, тот же набор данных, R vs SPSS). nnet пакет мультином функция

Недавно мне пришлось работать с R и SPSS для анализа набора данных с помощью структуры полиномиальной регрессии. Мы опросили некоторых участников (10–12 лет) и спросили, какая «профессиональная сфера» им нравится больше всего, а затем спросили, как часто они выходят в Интернет. Таким образом, результатом является «Категориальная переменная»: профессиональная сфера - «военный», «Я не знаю» и «Другая профессия»; а независимая переменная также является категориальной переменной (как часто вы заходите в Интернет ( «У меня нет доступа», «1-3 часа в день», «3-5 часов в день»).

Я запускал модель с использованием R (с пакетом nnet, с помощью функции multinom), а другой статистик запускал с помощью SPSS. Все ссылочные категории определены правильно.

Теперь, когда мы сравниваем результаты, они не совпадают для второй категории моей независимой переменной. Первый в порядке.

Взгляните на весь код:

library(tidyverse)
library(stargazer)
library(nnet)

ds <- ds %>% mutate(internet = factor(internet))
ds <- ds %>% mutate(internet = relevel(internet, ref = "I dont have internet access"))

ds <- ds %>% mutate(field = factor(field))
ds <- ds %>% mutate(fielf = relevel(field, ref = "I dont know"))

mod <- multinom(field ~ internet, data = ds, maxit=1000, reltol = 1.0e-9)
stargazer(mod, type = 'text')

и результат SPSS  Результаты SPSS

Для ясности, когда независимая переменная имеет только две категории (например, пол, мужской и женский), и R, и SPSS согласны с его результатами.

Результаты SPSS 2

После огромных усилий, пытающихся понять несоответствие между обоими результатами, я прочитал, что оценка nnet может иметь некоторые проблемы (проблема оптимизации?), и что расхождение результатов не так странно, как я думал вначале ..

Может кто-нибудь объяснить мне, что здесь происходит? Что мне не хватает ?! Я предполагаю, что SPSS и R должны иметь одинаковые результаты, если мы используем одну и ту же модель.

Спасибо

Это DS, который я использую в этом примере:

ds <- structure(list(sex = structure(c(2L, 1L, 2L, 1L, 2L, 1L, 2L, 
                                             2L, 1L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 1L, 
                                             2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
                                             2L, 1L, 2L, 1L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 
                                             1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 
                                             2L, 1L), .Label = c("male", "female"), class = "factor"), internet = structure(c(3L, 
                                                                                                                              3L, 2L, 3L, 2L, 3L, 3L, 3L, 2L, 2L, 3L, 3L, 3L, 2L, 2L, 3L, 2L, 
                                                                                                                              2L, 2L, 2L, 3L, 3L, 2L, 2L, 3L, 1L, 3L, 2L, 2L, 2L, 3L, 3L, 3L, 
                                                                                                                              2L, 3L, 3L, 2L, 3L, 3L, 3L, 3L, 3L, 1L, 2L, 2L, 3L, 3L, 2L, 2L, 
                                                                                                                              2L, 3L, 3L, 3L, 3L, 3L, 1L, 2L, 3L, 1L, 2L, 2L, 2L, 3L, 3L, 2L, 
                                                                                                                              2L, 1L, 3L, 2L, 2L, 3L, 2L, 2L), .Label = c("I dont have internet access", 
                                                                                                                                                                          "1-3 hours/day", "3-5 hours/day"), class = "factor"), field = structure(c(1L, 
                                                                                                                                                                                                                                                    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
                                                                                                                                                                                                                                                    1L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
                                                                                                                                                                                                                                                    1L, 1L, 1L, 2L, 1L, 1L, 1L, 3L, 1L, 1L, 1L, 3L, 1L, 1L, 1L, 1L, 
                                                                                                                                                                                                                                                    1L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 
                                                                                                                                                                                                                                                    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L), .Label = c("I dont know", "Military", 
                                                                                                                                                                                                                                                                                                "Other profession"), class = "factor")), class = "data.frame", row.names = c(NA, 
                                                                                                                                                                                                                                                                                                                                                                             -73L))

person Luis    schedule 24.04.2019    source источник


Ответы (1)


В качестве альтернативы вы можете использовать mlogit, который больше похож на результаты SPSS . Значения SPSS должны быть действительными, поскольку Stata дает аналогичные результаты (-14.88 (982.95), 11.58 (982.95), 11.44 (982.95)). Остальные отклонения могут быть связаны с нелепым значением «другой профессии».

library(mlogit)
ml.dat <- mlogit.data(ds, choice="field", shape="wide")
ml <- mlogit(field ~ 1 | internet, data=ml.dat)

Податливость

texreg::screenreg(ml)
=========================================================
                                                Model 1  
---------------------------------------------------------
Military:(intercept)                               -0.41 
                                                   (0.91)
Other profession:(intercept)                      -16.89 
                                                (2690.89)
Military:factor(internet)1-3 hours/day             -1.50 
                                                   (1.06)
Other profession:factor(internet)1-3 hours/day     13.60 
                                                (2690.89)
Military:factor(internet)3-5 hours/day             -1.64 
                                                   (1.06)
Other profession:factor(internet)3-5 hours/day     13.46 
                                                (2690.89)
---------------------------------------------------------
AIC                                                85.49 
Log Likelihood                                    -36.74 
Num. obs.                                          73    
=========================================================
*** p < 0.001, ** p < 0.01, * p < 0.05
person jay.sf    schedule 24.04.2019
comment
Спасибо! Сейчас результаты близки, но я все еще пытаюсь выяснить, что здесь происходит. Я вижу, что информационные критерии R и SPSS очень сильно различаются (R AIC = 85,49 против SPSS AIC = 24). Между прочим, пакет texreg создает красивую таблицу, очень похожую на Stargazer! Я не знал этого пакета, спасибо! - person Luis; 24.04.2019