Регрессия в цикле, но получим "Ошибка в коэффиценте (сводка (подгонка))[2, 4]: индекс за пределами"

Question

При запуске следующего кода в наборе данных большого размера я получаю сообщение об ошибке

Ошибка в coef(итоговая (подходящая))[2, 4]: индекс за пределами

Вектор beta для которого p-значения моделей логистической регрессии сохраняются, имеет длину 19481. Если я перебираю различные независимые переменные модели регрессии до 100 раз, я не получаю эту ошибку.

Кто-нибудь может дать мне подсказку, почему мой код не работает гладко?

beta = rep(0, 19481)
for (i in 25:19505) {
  fit = glm(mdr.mdr ~ an.mdr[,i], family=binomial)
  beta[i-24] = coef(summary(fit))[2,4]
  }

Answer 1

Как ошибка приходит при попытке извлечь [2,4] элемент таблицы коэффициентов, то есть значение р наклона, я уверен, что у вас есть NA оценка для склона.

Это означает, что для некоторых i Ваша модель имеет недостаток ранга, и нет информации для оценки наклона.

Обратите внимание, что, coef(summary(fit)) упадет NA оцените, так что в этой ситуации ваша таблица коэффициентов имеет только одну строку вместо двух строк (что объясняет ошибку "вне границы"). См. Таблицу коэффициентов, не имеющую рядов NA в подгонке с недостаточным рангом; как их вставить?

Я предлагаю следующее:

beta = rep(NA, 19481)
for (i in 25:19505) {
  fit = glm(mdr.mdr ~ an.mdr[,i], family = binomial)
  slope <- coef(fit)[2]
  if (!is.na(slope)) beta[i-24] = coef(summary(fit))[2,4]
  }

Еще один потенциальный сбой этого цикла - "нет завершенных случаев", то есть sum(complete.cases(mdr.mdr, an.mdr[, i])) дает вам 0. Если это произойдет, вы можете захотеть:

beta = rep(NA, 19481)
for (i in 25:19505) {
  if (sum(complete.cases(mdr.mdr, an.mdr[, i])) > 0) {
    fit = glm(mdr.mdr ~ an.mdr[,i], family = binomial)
    slope <- coef(fit)[2]
    if (!is.na(slope)) beta[i-24] = coef(summary(fit))[2,4]
    }
  }

Answer 2

Запустите цикл и, когда он остановится, наберите i в консоль и нажмите ввод. Это скажет вам, в какой итерации цикла цикл завершился неудачно. Затем осмотрите an.mdr[,i] за что-то неожиданное

Попробуй это:

  beta = rep(0, 19481) 
  for (i in 25:19505) {
      fit = glm(mdr.mdr ~ an.mdr[,i], family=binomial) 
      if (is.na(coef(summary(fit))[2,4]) { beta[i-24] = NA }
      if (!is.na(coef(summary(fit))[2,4]) { beta[i-24] = coef(summary(fit))[2,4] }
  }