1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
* Тема: Доверительные интервалы доли выживших на графике функции выживаемости.
* Ключевые слова: анализ выживаемости, графики, доверительные интервалы.
* Опубликован: 21.03.2008.
* Автор: Marta Garcia-Granero.
* Перевод: А. Балабанов.
* Размещение: http://www.spsstools.ru/Syntax/SurvivalAnalysis/Show95pcCIforFailurePointsOnSurvivalPlot.txt (.sps).

(Вопрос) Можно ли отобразить 95% доверительный интервал
 на графике функции выживаемости?

**************************
(Ответ) Для этого потребуется сначала подсчитать значения функции выживаемости,
 стандартные ошибки и границы доверительного интервала, а затем все это нарисовать
 на одном и том же графике.
* Автор решения: Marta Garcia-Granero, 06.18.2004.

* (1) Решение для одной группы наблюдений *.

* Пример данных *.
DATA LIST FREE /surtime(F8.0) status(F4.0).
BEGIN DATA
65 1 156 0 100 1 134 1 16 1 108 1 121 1 4 1
69 1 143 0 56 1 26 1 22 1 1 1 1 1 5 1 65 1
65 1 17 1 7 1 16 1 22 1 3 1 4 1 2 1 3 1 8 1
4 1 3 1 30 1 4 1 43 1 56 1
END DATA.
VALUE LABELS status 0 'Цензурированное' 1 'Событие'.

* Процедура KM с опцией /SAVE создаст 2 переменные: значения функции выживаемости (sur_1)
  и её стандартную ошибку (se_1). Для ориентира выведем также график функции*.
KM surtime /STATUS=status(1)
  /PRINT TABLE MEAN
  /PLOT SURVIVAL
  /SAVE SURVIVAL SE .

* Для цензурированных данных значения функции и её стандартная ошибка не
  вычисляются. Мы припишем им значения из ближайших предыдущих наблюдений *.
SORT CASES BY surtime(A).
IF (MISSING(sur_1)) sur_1 = LAG(sur_1) .
IF (MISSING(se_1)) se_1 = LAG(se_1) .

* Вычисляем нижние и верхние границы интервала, заменяя 
  границы, выходящие за пределы 0 или 1 на теоретический минимум/максимум *.
COMPUTE lower = MAX((sur_1-1.96*se_1),0) .
COMPUTE upper = MIN((sur_1+1.96*se_1),1) .
EXECUTE.

* Теперь рисуем доверительные уровни вместе с функцией выживаемости
  (сравните форму графика с формой графика из процедуры KV) *.
GRAPH
 /SCATTERPLOT(OVERLAY)=surtime WITH sur_1 lower upper.

* Далее следует "ручное" вмешательство (настройка графика).
* Потребуется: ввести интерполирующую линию, определив ей
  характер ступенчатой функции (Left Step), а также, по желанию,
  убрать маркеры.
  
* Кроме того, вы можете добавить заголовки к осям, поменять цвета и т.д.

* Можете также сохранить шаблон графика, чтобы в будущем меньше
  возиться с ручной настройкой. Использовать этот шаблон следует
  так (пример):
* GRAPH
   /SCATTERPLOT(OVERLAY)=surtime WITH sur_1 lower upper
   /TEMPLATE='C:\\temp\\surv.sgt'.


**************************
* (2) Решение для двух групп наблюдений. Несколько более сложное, т.к.
  процедура SCATTERPLOT(OVERLAY) не принимает группировку через ключевое поле BY *.

* Пример данных из учебника "Statistics at Square One" (Swinscow&Campell) *.
DATA LIST FREE/treatmen(F8.0) surtime(F8.1) status(F8.0).
BEGIN DATA
1 1 0 1 5 0 1 6 1 1 6 1 1 9 0 1 10 1 1 10 1 1 10 0 1 12 1
1 12 1 1 12 1 1 12 1 1 12 0 1 13 0 1 15 0 1 16 0 1 20 0
1 24 1 1 24 0 1 27 0 1 32 1 1 34 0 1 36 0 1 36 0 1 44 0
0 3 0 0 6 1 0 6 1 0 6 1 0 6 1 0 8 1 0 8 1 0 12 1 0 12 1
0 12 0 0 15 0 0 16 0 0 18 0 0 18 0 0 20 1 0 22 0 0 24 1
0 28 0 0 28 0 0 28 0 0 30 1 0 30 0 0 33 0 0 42 1
END DATA.

VAR LABELS treatmen 'Характер лечения'
 /surtime 'Время жизни (мес.)'
 /status 'Статус'.
VALUE LABELS treatmen 1'Gamma-Linoleic' 0'контрольная группа'
  /status 0 'Жив' 1 'Мёртв'.

* Получим оценки функции выживаемости и её стандартные ошибки для обеих групп *.
TEMPORARY.
SELECT IF (treatmen EQ 0).
KM surtime
  /STATUS=status(1)
  /PRINT TABLE MEAN
  /SAVE SURVIVAL SE .
TEMPORARY.
SELECT IF (treatmen EQ 1).
KM surtime
  /STATUS=status(1)
  /PRINT TABLE MEAN
  /SAVE SURVIVAL SE .

* Вычисление пределов 95%-го ДИ *.
DO REPEAT A=lower1 lower2 /B=upper1 upper2 /C=sur_1 sur_2 /D=se_1 se_2.
- COMPUTE A = C-1.96*D .
- COMPUTE B = C+1.96*D .
END REPEAT.
EXECUTE.

* Коррекция границ интервалов, выходящих за теоретические пределы функции (0-1) *.
DO IF (lower1 LT 0) AND (treatmen EQ 0).
- COMPUTE lower1=0.
END IF.
DO IF upper1 GT 1 AND (treatmen EQ 0).
- COMPUTE upper1=1.
END IF.
DO IF lower2 LT 0 AND (treatmen EQ 1).
- COMPUTE lower2=0.
END IF.
DO IF upper2 GT 1 AND (treatmen EQ 1).
- COMPUTE upper2=1.
END IF.
EXECUTE.

* Заполнение пропущенных значений для цензурированных данных *.
SORT CASES BY treatmen(A) surtime(A) .
DO IF ($casenum EQ 1) AND (MISSING(sur_1)).
- COMPUTE sur_1=1.
- COMPUTE lower1=1.
- COMPUTE upper1=1.
END IF.
DO IF MISSING(sur_1) AND (treatmen EQ 0).
- COMPUTE sur_1=LAG(sur_1).
- COMPUTE lower1=LAG(lower1).
- COMPUTE upper1=LAG(upper1).
END IF.
SORT CASES BY treatmen(D) surtime(A) .
DO IF ($casenum EQ 1) AND (MISSING(sur_2)).
- COMPUTE sur_2=1.
- COMPUTE lower2=1.
- COMPUTE upper2=1.
END IF.
DO IF MISSING(sur_2) AND (treatmen EQ 1).
- COMPUTE sur_2=LAG(sur_2).
- COMPUTE lower2=LAG(lower2).
- COMPUTE upper2=LAG(upper2).
END IF.
EXEC.


VARIABLE LABELS
 sur_1 'Функция1'
 sur_2 'Функция2'
 lower1 'НГ95%ДИ/ф-я1'
 upper1 'ВГ95%ДИ/ф-я1'
 lower2 'НГ95%ДИ/ф-я2'
 upper2 'ВГ95%ДИ/ф-я2'.
GRAPH
 /SCATTERPLOT(OVERLAY)=surtime WITH sur_1 lower1 upper1 sur_2 lower2 upper2
 /MISSING=VARIABLEWISE.

* Далее проводим аналогичные "ручные" настройки графика (после двойного щелчка мышкой на нём), 
  которые делали в Решении №1.