Матрица Coxph X считается сингулярной без идеальной классификации?

Я пытаюсь провести самоконтролируемую серию случаев, используя пакет sccs из здесь. Это статистический метод, который берет «исходный» и «выставленный» периоды в течение времени, например, год жизни больного. Периоды воздействия могут представлять собой воздействие препарата, а измеряемый результат может быть побочным эффектом препарата, как, собственно, и было в моем случае.

Этот пакет по существу форматирует данные в интервалы исходного уровня и риска воздействия. Идентификатор пациента indivL (фактор), interval (целое число, количество дней), exposure статус (0/1), event статус (0/1). Затем он передает эти данные в survival::clogit как модель формы:

event ~ exposure + strata(indivL) + offset(log(interval))

Данные, подаваемые на clogit, представляют собой фрейм данных в форме:

   indivL event eventday  lower  upper interval age   exposure  indiv aevent astart   aend drugtype
 * <fct>  <dbl>    <int>  <dbl>  <dbl>    <dbl> <fct> <fct>    <dbl>  <dbl>  <dbl>  <dbl>     <dbl>
 1 1         0.    22361 22219. 22252.      34. 1     0           1. 22361. 22219. 22460.        0.
 2 1         0.    22361 22253. 22260.       8. 1     1           1. 22361. 22219. 22460.        0.
 3 1         1.    22361 22261. 22460.     200. 1     0           1. 22361. 22219. 22460.        0.
 4 2         0.    22401 22219. 22252.      34. 1     0           1. 22401. 22219. 22460.        0.
 5 2         0.    22401 22253. 22260.       8. 1     1           1. 22401. 22219. 22460.        0.
 6 2         1.    22401 22261. 22460.     200. 1     0           1. 22401. 22219. 22460.        0.
 7 3         0.    31071 30834. 30863.      30. 1     0           2. 31071. 30834. 31075.        0.
 8 3         0.    31071 30864. 30871.       8. 1     1           2. 31071. 30834. 31075.        0.
 9 3         1.    31071 30872. 31075.     204. 1     0           2. 31071. 30834. 31075.        0.
10 4         1.      261   207.   356.     150. 1     0           3.   261.   207.   425.        0.
# ... with 1,211,460 more rows

У меня есть моя модель, которая хорошо работает, чтобы дать мне результат при использовании вышеизложенного. Однако я хочу добавить другие независимые переменные. Это бинарные категориальные числа, и я пробовал их как целые числа 0/1, так и двухуровневые факторы. Одним из примеров может быть drugtype. В этом случае модель принимает вид:

event ~ exposure + drugtype + strata(indivL) + offset(log(interval))

Моя ошибка:

Warning message:
In coxph(formula = Surv(rep(1, 176241L), event) ~ exposure + drugtype +  :
  X matrix deemed to be singular; variable 2

Моя модель:

--SNIP--
coxph(formula = Surv(rep(1, 176241L), event) ~ exposure + drugtype + 
     strata(indivL) + offset(log(interval)), data = chopdat, method = "exact")

  n= 176049, number of events= 58602 
   (192 observations deleted due to missingness)

             coef exp(coef) se(coef)     z            Pr(>|z|)    
exposure1 0.70760   2.02912  0.01662 42.57 <0.0000000000000002 ***
drugtype       NA        NA  0.00000    NA                  NA    
--SNIP--

Как видите, ему не нравится drugtype, которая является бинарной переменной.

Оглядевшись, я наткнулся на несколько источников, которые предполагают, что проблема связана с «идеальной классификацией», то есть одна из моих переменных точно предсказывает наличие другой. Однако, используя xtabs(), я получаю:

> xtabs(~drugtype + event, data = chopdat)

         event
 drugtype      0      1
        0 778306 388279
        1  29344  14625

а также

> xtabs(~ exposure + event, data = chopdat)

        event
exposure      0      1
       0 427482 380101
       1 380788  23113

а также

> xtabs(~ drugtype + exposure, data = chopdat)

       drugtype
exposure      0      1
       0 777655  29308
       1 388930  14661

Предполагая, что существует хорошее распределение и нет идеальной классификации.

Может ли кто-нибудь указать мне в правильном направлении для получения дополнительной информации об этом? Я чувствую, что достиг предела своих возможностей с документацией и поиском других ответов на этот вопрос в StackOverflow.

Большое спасибо.


person JWC    schedule 11.06.2018    source источник
comment
А как насчет drugtype от events?   -  person Mike    schedule 11.06.2018
comment
Спасибо @Mike, я обновил вопрос этим xtabs()   -  person JWC    schedule 11.06.2018
comment
в чем разница между кадром данных dat и кадром данных chopdat, chopdat используется в модели, а dat используется в xtabs?   -  person Mike    schedule 11.06.2018
comment
Извините, @Mike, упущение с моей стороны - это один и тот же фреймворк данных. Я изменил вопрос соответственно.   -  person JWC    schedule 11.06.2018
comment
Когда вы включаете в модель все остальные переменные, может быть идеальное разделение, если вы просто запустите модель с типом наркотика, вы получите эту ошибку?   -  person Mike    schedule 11.06.2018
comment
@Mike: даже работая как event ~ drugtype, я получаю: X matrix deemed to be singular; variable 1   -  person JWC    schedule 12.06.2018
comment
Без данных я не могу воспроизвести ошибку/быть более полезным. Только из приведенного выше фрагмента я вижу, что drugtype имеет только значения 0 (я понимаю, что вы не показали весь фрейм данных). Также почему вы указываете компонент времени как вектор от 1 до длины вашего df?   -  person Mike    schedule 12.06.2018


Ответы (1)


Ладно, всем привет.

Большое спасибо @Mike за попытку помочь мне. Я нашел объяснение этому поведению, и это особенность метода моделирования SCCS, а не самого clogit.

Из Farrington, Whitaker and Weldeselassie Self-Controlled Series Studies, A Modeling with R:

Обратите внимание, что основной эффект ковариаты не включен в формулу модели, поскольку он не может быть оценен в модели SCCS, поскольку он выпадает из вероятности.

Таким образом, основной эффект не может быть оценен и поэтому возвращается как NA.

Извиняюсь за публикацию вопроса, на который я должен был найти ответ, но, надеюсь, это будет полезно всем, кто столкнется с этой «проблемой».

person JWC    schedule 18.06.2018