R множественная логистическая регрессия (пакет mlogit)

Я пытаюсь оценить полиномиальную логистическую регрессию с помощью пакета mlogit, но у меня возникла ошибка. У меня есть 205109 наблюдений, и ниже я попытался ввести код и ошибку.

> Data <- read.csv("NWEScsv.csv",header=T)
> colnames(Data)
 [1] "NWSE"  "W.S"   "W.E"   "W.N"   "W.W"   "U.S"   "U.E"   "U.N"   "U.W"  
[10] "W.Not" "U.Not" "BW"    "BU" 
> movement = mlogit.data(Data, shape ="wide", varying = 2:11, choice = "NWSE")
> head(movement,15)
       NWSE     BW BU alt        W           U chid
1.E   FALSE 194.42  0   E 198.3434 2.47404e-10    1
1.N   FALSE 194.42  0   N 194.4160 1.41319e-10    1
1.Not  TRUE 194.42  0 Not 194.4200 0.00000e+00    1
1.S   FALSE 194.42  0   S 212.7249 2.08726e-10    1
1.W   FALSE 194.42  0   W 198.3434 1.37143e-10    1
2.E    TRUE 257.24  0   E 202.3502 8.05144e-10    2
2.N   FALSE 257.24  0   N 200.3368 4.59906e-10    2
> result.movement = mlogit(NWSE ~ 0 |BW+BU | W+U, movement)
 solve.default(H, g[!fixed]) error: 
  Lapack routine dgesv: system is exactly singular: U[20,20] = 0 

Как решить проблему?

========

Извините за длинный вопрос. Насколько мне известно, когда мы пытаемся оценить множественную логистическую регрессию, нам необходимо разделить независимые переменные на три типа. (например: U_ {ij} = α_ {j} + β * X_ {ij} + γ_ {j} Z_ {i} + δ_ {j} W_ {ij} + ε_ {ij}, X_ {ij}, Z_ {i}, W_ {ij} - три типа) В пакете mlogit нам нужно разделить переменные с помощью "|" и уравнение примет вид X_ {ij} | Z_ {i} | W_ {ij}.

Независимо от того, сколько раз я пытался, я не мог получить результаты, когда помещал переменные в W_ {ij}. Ниже я помещаю часть структуры данных.

> dput(head(movement, n = 25))
structure(list(NWSE = c(`1.E` = FALSE, `1.N` = FALSE, `1.Not` = TRUE, 
`1.S` = FALSE, `1.W` = FALSE, `2.E` = TRUE, `2.N` = FALSE, `2.Not` = FALSE, 
`2.S` = FALSE, `2.W` = FALSE, `3.E` = TRUE, `3.N` = FALSE, `3.Not` = FALSE, 
`3.S` = FALSE, `3.W` = FALSE, `4.E` = TRUE, `4.N` = FALSE, `4.Not` = FALSE,
`4.S` = FALSE, `4.W` = FALSE, `5.E` = TRUE, `5.N` = FALSE, `5.Not` = FALSE, 
`5.S` = FALSE, `5.W` = FALSE), BW = c(`1.E` = 194.42, `1.N` = 194.42,
`1.Not` = 194.42, `1.S` = 194.42, `1.W` = 194.42, `2.E` = 257.24, 
`2.N` = 257.24, `2.Not` = 257.24, `2.S` = 257.24, `2.W` = 257.24, 
`3.E` = 262.43, `3.N` = 262.43, `3.Not` = 262.43, `3.S` = 262.43, 
`3.W` = 262.43, `4.E` = 183.09, `4.N` = 183.09, `4.Not` = 183.09, 
`4.S` = 183.09, `4.W` = 183.09, `5.E` = 311.06, `5.N` = 311.06, 
`5.Not` = 311.06, `5.S` = 311.06, `5.W` = 311.06), BU = c(`1.E` = 0, 
`1.N` = 0, `1.Not` = 0, `1.S` = 0, `1.W` = 0, `2.E` = 0, `2.N` = 0, 
`2.Not` = 0, `2.S` = 0, `2.W` = 0, `3.E` = 0, `3.N` = 0, `3.Not` = 0, 
`3.S` = 0, `3.W` = 0, `4.E` = 0, `4.N` = 0, `4.Not` = 0, `4.S` = 0, 
`4.W` = 0, `5.E` = 0, `5.N` = 0, `5.Not` = 0, `5.S` = 0, `5.W` = 0
), alt = c(`1.E` = "E", `1.N` = "N", `1.Not` = "Not", `1.S` = "S", 
`1.W` = "W", `2.E` = "E", `2.N` = "N", `2.Not` = "Not", `2.S` = "S", 
`2.W` = "W", `3.E` = "E", `3.N` = "N", `3.Not` = "Not", `3.S` = "S", 
`3.W` = "W", `4.E` = "E", `4.N` = "N", `4.Not` = "Not", `4.S` = "S", 
`4.W` = "W", `5.E` = "E", `5.N` = "N", `5.Not` = "Not", `5.S` = "S", 
`5.W` = "W"), W = c(`1.E` = 198.3434254, `1.N` = 194.4159624, 
`1.Not` = 194.42, `1.S` = 212.7249464, `1.W` = 198.3434254, `2.E` = 
202.3502284, 
`2.N` = 200.33681, `2.Not` = 257.24, `2.S` = 219.2033856, `2.W` = 
204.383882, 
`3.E` = 208.5127103, `3.N` = 206.4379742, `3.Not` = 262.43, `3.S` = 
225.8791225, 
`3.W` = 210.6082979, `4.E` = 198.3434254, `4.N` = 194.4159624, 
`4.Not` = 183.09, `4.S` = 212.7249464, `4.W` = 198.3434254, `5.E` = 
 244.6919323, 
`5.N` = 239.8467074, `5.Not` = 311.06, `5.S` = 262.4340992, `5.W` = 
244.6919323
), U = c(`1.E` = 2.47e-10, `1.N` = 1.41e-10, `1.Not` = 0, `1.S` = 2.09e- 
  10, 
`1.W` = 1.37e-10, `2.E` = 8.05e-10, `2.N` = 4.6e-10, `2.Not` = 0, 
`2.S` = 6.73e-10, `2.W` = 4.46e-10, `3.E` = 6.53e-10, `3.N` = 3.73e-10, 
`3.Not` = 0, `3.S` = 5.51e-10, `3.W` = 3.62e-10, `4.E` = 2.47e-10, 
`4.N` = 1.41e-10, `4.Not` = 0, `4.S` = 2.09e-10, `4.W` = 1.37e-10, 
`5.E` = 2.65e-10, `5.N` = 1.52e-10, `5.Not` = 0, `5.S` = 2.24e-10, 
`5.W` = 1.47e-10), chid = c(`1.E` = 1L, `1.N` = 1L, `1.Not` = 1L, 
`1.S` = 1L, `1.W` = 1L, `2.E` = 2L, `2.N` = 2L, `2.Not` = 2L, 
`2.S` = 2L, `2.W` = 2L, `3.E` = 3L, `3.N` = 3L, `3.Not` = 3L, 
`3.S` = 3L, `3.W` = 3L, `4.E` = 4L, `4.N` = 4L, `4.Not` = 4L, 
`4.S` = 4L, `4.W` = 4L, `5.E` = 5L, `5.N` = 5L, `5.Not` = 5L, 
`5.S` = 5L, `5.W` = 5L)), reshapeLong = list(varying = structure(list(
W = c("W.S", "W.E", "W.N", "W.W", "W.Not"), U = c("U.S", 
"U.E", "U.N", "U.W", "U.Not")), v.names = c("W", "U"), times = c("S", 
"E", "N", "W", "Not")), v.names = c("W", "U"), idvar = "chid", 
timevar = "alt"), index = structure(list(chid = structure(c(1L, 
1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 
4L, 4L, 4L, 5L, 5L, 5L, 5L, 5L), .Label = c("1", "2", "3", "4", 
"5"), class = "factor"), alt = structure(c(1L, 2L, 3L, 4L, 5L, 
1L, 2L, 3L, 4L, 5L, 1L, 2L, 3L, 4L, 5L, 1L, 2L, 3L, 4L, 5L, 1L, 
2L, 3L, 4L, 5L), .Label = c("E", "N", "Not", "S", "W"), class = "factor")), 
class = "data.frame", row.names = c("1.E", 
"1.N", "1.Not", "1.S", "1.W", "2.E", "2.N", "2.Not", "2.S", "2.W", 
"3.E", "3.N", "3.Not", "3.S", "3.W", "4.E", "4.N", "4.Not", "4.S", 
"4.W", "5.E", "5.N", "5.Not", "5.S", "5.W")), choice = "NWSE", row.names = 
c("1.E", 
"1.N", "1.Not", "1.S", "1.W", "2.E", "2.N", "2.Not", "2.S", "2.W", 
"3.E", "3.N", "3.Not", "3.S", "3.W", "4.E", "4.N", "4.Not", "4.S", 
"4.W", "5.E", "5.N", "5.Not", "5.S", "5.W"), class = c("mlogit.data", 
"data.frame"))

person SEIKA UTA    schedule 25.11.2019    source источник
comment
У вас идеальное разделение или только одно наблюдение на категорию?   -  person MDEWITT    schedule 25.11.2019
comment
Я добавил главу данных. BW и BU зависят только от человека, а W и U зависят от человека и его выбора.   -  person SEIKA UTA    schedule 25.11.2019
comment
Вы можете сделать table(movement$BU) и посмотреть, есть ли вариации. Другой вариант - lapply( movement[,sapply(movement,is.numeric)] , sd). Это проверит, есть ли у вас какие-либо вариации в ваших столбцах. Это может быть причиной проблем с посадкой.   -  person MDEWITT    schedule 26.11.2019
comment
@MDEWITT Спасибо за ваш комментарий. Я попытался увидеть вариацию. ›Таблица (движение $ BU) 0 0,01 0,02 0,03 0,04 0,05 0,06 0,07 0,08 0,09 0,1 0,11 0,12 0,13 526520 155455 64380 35180 24465 19310 18560 13700 14055 9865 9765 8350 8295 7230› lapply (движение [, sapply (движение, числовое)] , sd) $ BW [1] 93.83399 $ BU [1] 0.1472054 $ W [1] 89.31311 $ U [1] 0.1319652 $ chid [1] 59209.61 Но как я могу оценить результаты?   -  person SEIKA UTA    schedule 26.11.2019
comment
гипотеза заключалась в том, что у вас не было достаточно вариаций в ваших предикторах, чтобы обеспечить стабильную подгонку из-за коллинеарности (например, ваша матрица предикторов не имела полного ранга, имела линейные комбинации друг друга или не менялась и действовала как вторая перехват и, таким образом, был коллинеарен с перехватом). Похоже, что это так. Не могли бы вы попробовать caret::nearZeroVar(movement)?   -  person MDEWITT    schedule 26.11.2019
comment
@MDEWITT Я получил результат ниже ›caret :: nearZeroVar (движение) целое число (0)   -  person SEIKA UTA    schedule 27.11.2019
comment
попробуйте caret::findLinearCombos(movement)   -  person MDEWITT    schedule 27.11.2019
comment
@ MDEWITTqr.default (объект) Ошибка в qr.default (.swts * attr (rhs, gradient)): NA / NaN / Inf в вызове внешней функции (аргумент 1)   -  person SEIKA UTA    schedule 27.11.2019
comment
Я думаю, что у меня есть все возможное, чтобы устранить неполадки без дополнительных данных (например, dput(head(movement, n = 25)) и поместите это в сообщение, если можете). Я думаю, у вас есть несколько коллинеарных предикторов, поэтому вы получаете сингулярную матрицу, которую нельзя инвертировать.   -  person MDEWITT    schedule 27.11.2019
comment
@MDEWITT Я отредактировал свой вопрос и добавил результаты, например. dput (голова (движение, n = 25)).   -  person SEIKA UTA    schedule 27.11.2019


Ответы (1)


Таким образом, похоже, что ваши данные не имеют полного ранга, например высокая коллинеарность. Мы можем проверить это, исследуя ранг данных по сравнению с количеством предикторов, если ранг не равен количеству предикторов, у вас будет неуникальное соответствие, которое является ошибкой, которую вы получаете.

Например:

library(dplyr)
movement %>% 
    select_if(is.numeric) %>% 
    {
        cat("num predictors", ncol(.), "\n")
        qr(.) %>% 
            .$rank
    } 

Урожайность:

num predictors 5 
[1] 4

Вы также можете увидеть это, посмотрев на корреляцию предикторов:

movement %>% 
    select_if(is.numeric) %>% 
    cor()

Который дает:

           BW BU          W          U       chid
BW   1.0000000 NA  0.7360761  0.2079872  0.4765022
BU          NA  1         NA         NA         NA
W    0.7360761 NA  1.0000000 -0.2686001  0.4933215
U    0.2079872 NA -0.2686001  1.0000000 -0.1964596
chid 0.4765022 NA  0.4933215 -0.1964596  1.0000000

Это указывает, например, что BW и W предоставляют в основном одинаковую информацию. Я не знаю вашего предмета, но вам, вероятно, не следует включать оба этих предиктора в вашу модель. BU также не имеет никаких вариаций, так что это вызовет ошибки.

person MDEWITT    schedule 27.11.2019