logo
搜索
我的图书馆

发文章

发文工具

其他工具

留言交流
搜索
分享

R生存分析AFT

 勤悦轩 2018-07-08

 γ = 1/scale =1/0.902
 α = exp(−(Intercept)γ)=exp(-(7.111)*γ)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
> library(survival)
> myfit=survreg(Surv(futime, fustat)~1 , ovarian, dist="weibull",scale=0)
> summary(myfit)
Call:
survreg(formula = Surv(futime, fustat) ~ 1, data = ovarian, dist = "weibull",
    scale = 0)
             Value Std. Error      z         p
(Intercept)  7.111      0.293 24.292 2.36e-130
Log(scale)  -0.103      0.254 -0.405  6.86e-01
Scale= 0.902
Weibull distribution
Loglik(model)= -98   Loglik(intercept only)= -98
Number of Newton-Raphson Iterations: 5
n= 26

 

画生存函数图

1
2
3
4
5
6
7
8
9
10
11
d<- seq(0, 2000, length.out=10000)
h<-1-pweibull(d,shape=1/0.902,scale=exp(7.111))
df<-data.frame(t=d,s=h)
library(ggplot2)
ggplot(df,aes(x=t,y=s))+
  geom_line(colour="green")+
  ggtitle("s(t) \n 生存函数")

 

 

 

1. Surv

Description

创建一个生存对象，通常用作模型公式中的响应变量。 参数匹配是此功能的特殊功能，请参阅下面的详细信息。

Surv(time, time2, event,
    type=c('right', 'left', 'interval', 'counting', 'interval2', 'mstate'),
    origin=0)
is.Surv(x)

Arguments
time
对于右删失数据，这是一个跟踪时间。对于区间数据，第一个参数是区间的开始时间。 
event 
状态指示，通常，0=活着，1=死亡。其他选择是TRUE/FALSE (TRUE = 死亡) or 1/2 (2=死亡)。对于区间删失数据，状态指示，0=右删失, 1=事件时间, 2=左删失, 3=区间删失.
右删失(Right Censoring)：只知道实际寿命大于某数；
左删失(Left Censoring)：只知道实际寿命小于某数；
区间删失(Interval Censoring)：只知道实际寿命在一个时间区间内。
time2 
区间删失区间的结束时间或仅对过程数据进行计数。
type 
指定删失类型。 "right", "left", "counting", "interval", "interval2" or "mstate".
如果event变量是一个因子，假定type="mstate"。如果没有指定参数time2，type="right"；如果指定参数time2，type="counting" 

 

Surv使用示例

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
> str(lung)
'data.frame':   228 obs. of  10 variables:
 $ inst     : num  3 3 3 5 1 12 7 11 1 7 ...
 $ time     : num  306 455 1010 210 883 ...
 $ status   : num  2 2 1 2 2 1 2 2 2 2 ...
 $ age      : num  74 68 56 57 60 74 68 71 53 61 ...
 $ sex      : num  1 1 1 1 1 1 2 2 1 1 ...
 $ ph.ecog  : num  1 0 0 1 0 1 2 2 1 2 ...
 $ ph.karno : num  90 90 90 90 100 50 70 60 70 70 ...
 $ pat.karno: num  100 90 90 60 90 80 60 80 80 70 ...
 $ meal.cal : num  1175 1225 NA 1150 NA ...
 $ wt.loss  : num  NA 15 15 11 0 0 10 1 16 34 ...
> with(lung, Surv(time, status))
  [1]  306   455  1010+  210   883  1022+  310   361   218
 [10]  166   170   654   728    71   567   144   613   707
 [19]   61    88   301    81   624   371   394   520   574
 [28]  118   390    12   473    26   533   107    53   122
 [37]  814   965+   93   731   460   153   433   145   583
 [46]   95   303   519   643   765   735   189    53   246
 [55]  689    65     5   132   687   345   444   223   175
 [64]   60   163    65   208   821+  428   230   840+  305
 [73]   11   132   226   426   705   363    11   176   791
 [82]   95   196+  167   806+  284   641   147   740+  163
 [91]  655   239    88   245   588+   30   179   310   477
[100]  166   559+  450   364   107   177   156   529+   11
[109]  429   351    15   181   283   201   524    13   212
[118]  524   288   363   442   199   550    54   558   207
[127]   92    60   551+  543+  293   202   353   511+  267
[136]  511+  371   387   457   337   201   404+  222    62
[145]  458+  356+  353   163    31   340   229   444+  315+
[154]  182   156   329   364+  291   179   376+  384+  268
[163]  292+  142   413+  266+  194   320   181   285   301+
[172]  348   197   382+  303+  296+  180   186   145   269+
[181]  300+  284+  350   272+  292+  332+  285   259+  110
[190]  286   270    81   131   225+  269   225+  243+  279+
[199]  276+  135    79    59   240+  202+  235+  105   224+
[208]  239   237+  173+  252+  221+  185+   92+   13   222+
[217]  192+  183   211+  175+  197+  203+  116   188+  191+
[226]  105+  174+  177+
> str(heart)
'data.frame':   172 obs. of  8 variables:
 $ start     : num  0 0 0 1 0 36 0 0 0 51 ...
 $ stop      : num  50 6 1 16 36 39 18 3 51 675 ...
 $ event     : num  1 1 0 1 0 1 1 1 0 1 ...
 $ age       : num  -17.16 3.84 6.3 6.3 -7.74 ...
 $ year      : num  0.123 0.255 0.266 0.266 0.49 ...
 $ surgery   : num  0 0 0 0 0 0 0 0 0 0 ...
 $ transplant: Factor w/ 2 levels "0","1": 1 1 1 2 1 2 1 1 1 2 ...
 $ id        : num  1 2 3 3 4 4 5 6 7 7 ...
> Surv(heart$start, heart$stop, heart$event)
  [1] (  0.0,  50.0]  (  0.0,   6.0]  (  0.0,   1.0+]
  [4] (  1.0,  16.0]  (  0.0,  36.0+] ( 36.0,  39.0]
  [7] (  0.0,  18.0]  (  0.0,   3.0]  (  0.0,  51.0+]
 [10] ( 51.0, 675.0]  (  0.0,  40.0]  (  0.0,  85.0]
 [13] (  0.0,  12.0+] ( 12.0,  58.0]  (  0.0,  26.0+]
 [16] ( 26.0, 153.0]  (  0.0,   8.0]  (  0.0,  17.0+]
 [19] ( 17.0,  81.0]  (  0.0,  37.0+] ( 37.0,1387.0]
 [22] (  0.0,   1.0]  (  0.0,  28.0+] ( 28.0, 308.0]
 [25] (  0.0,  36.0]  (  0.0,  20.0+] ( 20.0,  43.0]
 [28] (  0.0,  37.0]  (  0.0,  18.0+] ( 18.0,  28.0]
 [31] (  0.0,   8.0+] (  8.0,1032.0]  (  0.0,  12.0+]
 [34] ( 12.0,  51.0]  (  0.0,   3.0+] (  3.0, 733.0]
 [37] (  0.0,  83.0+] ( 83.0, 219.0]  (  0.0,  25.0+]
 [40] ( 25.0,1800.0+] (  0.0,1401.0+] (  0.0, 263.0]
 [43] (  0.0,  71.0+] ( 71.0,  72.0]  (  0.0,  35.0]
 [46] (  0.0,  16.0+] ( 16.0, 852.0]  (  0.0,  16.0]
 [49] (  0.0,  17.0+] ( 17.0,  77.0]  (  0.0,  51.0+]
 [52] ( 51.0,1587.0+] (  0.0,  23.0+] ( 23.0,1572.0+]
 [55] (  0.0,  12.0]  (  0.0,  46.0+] ( 46.0, 100.0]
 [58] (  0.0,  19.0+] ( 19.0,  66.0]  (  0.0,   4.5+]
 [61] (  4.5,   5.0]  (  0.0,   2.0+] (  2.0,  53.0]
 [64] (  0.0,  41.0+] ( 41.0,1408.0+] (  0.0,  58.0+]
 [67] ( 58.0,1322.0+] (  0.0,   3.0]  (  0.0,   2.0]
 [70] (  0.0,  40.0]  (  0.0,   1.0+] (  1.0,  45.0]
 [73] (  0.0,   2.0+] (  2.0, 996.0]  (  0.0,  21.0+]
 [76] ( 21.0,  72.0]  (  0.0,   9.0]  (  0.0,  36.0+]
 [79] ( 36.0,1142.0+] (  0.0,  83.0+] ( 83.0, 980.0]
 [82] (  0.0,  32.0+] ( 32.0, 285.0]  (  0.0, 102.0]
 [85] (  0.0,  41.0+] ( 41.0, 188.0]  (  0.0,   3.0]
 [88] (  0.0,  10.0+] ( 10.0,  61.0]  (  0.0,  67.0+]
 [91] ( 67.0, 942.0+] (  0.0, 149.0]  (  0.0,  21.0+]
 [94] ( 21.0, 343.0]  (  0.0,  78.0+] ( 78.0, 916.0+]
 [97] (  0.0,   3.0+] (  3.0,  68.0]  (  0.0,   2.0]
[100] (  0.0,  69.0]  (  0.0,  27.0+] ( 27.0, 842.0+]
[103] (  0.0,  33.0+] ( 33.0, 584.0]  (  0.0,  12.0+]
[106] ( 12.0,  78.0]  (  0.0,  32.0]  (  0.0,  57.0+]
[109] ( 57.0, 285.0]  (  0.0,   3.0+] (  3.0,  68.0]
[112] (  0.0,  10.0+] ( 10.0, 670.0+] (  0.0,   5.0+]
[115] (  5.0,  30.0]  (  0.0,  31.0+] ( 31.0, 620.0+]
[118] (  0.0,   4.0+] (  4.0, 596.0+] (  0.0,  27.0+]
[121] ( 27.0,  90.0]  (  0.0,   5.0+] (  5.0,  17.0]
[124] (  0.0,   2.0]  (  0.0,  46.0+] ( 46.0, 545.0+]
[127] (  0.0,  21.0]  (  0.0, 210.0+] (210.0, 515.0+]
[130] (  0.0,  67.0+] ( 67.0,  96.0]  (  0.0,  26.0+]
[133] ( 26.0, 482.0+] (  0.0,   6.0+] (  6.0, 445.0+]
[136] (  0.0, 428.0+] (  0.0,  32.0+] ( 32.0,  80.0]
[139] (  0.0,  37.0+] ( 37.0, 334.0]  (  0.0,   5.0]
[142] (  0.0,   8.0+] (  8.0, 397.0+] (  0.0,  60.0+]
[145] ( 60.0, 110.0]  (  0.0,  31.0+] ( 31.0, 370.0+]
[148] (  0.0, 139.0+] (139.0, 207.0]  (  0.0, 160.0+]
[151] (160.0, 186.0]  (  0.0, 340.0]  (  0.0, 310.0+]
[154] (310.0, 340.0+] (  0.0,  28.0+] ( 28.0, 265.0+]
[157] (  0.0,   4.0+] (  4.0, 165.0]  (  0.0,   2.0+]
[160] (  2.0,  16.0]  (  0.0,  13.0+] ( 13.0, 180.0+]
[163] (  0.0,  21.0+] ( 21.0, 131.0+] (  0.0,  96.0+]
[166] ( 96.0, 109.0+] (  0.0,  21.0]  (  0.0,  38.0+]
[169] ( 38.0,  39.0+] (  0.0,  31.0+] (  0.0,  11.0+]
[172] (  0.0,   6.0]

 

2.survreg

拟合参数生存回归模型。 这些是用于时间变量的任意变换的位置尺度模型; 最常见的情况使用对数变换，建立加速失效时间模型。

survreg(formula, data, weights, subset,
        na.action, dist="weibull", init=NULL, scale=0,
        control,parms=NULL,model=FALSE, x=FALSE,
        y=TRUE, robust=FALSE, score=FALSE, ...)
  
dist
y变量的假设分布。"weibull", "exponential", "gaussian", "logistic","lognormal" and "loglogistic"。

scale 
可选的固定值。如果设置<=0，scale将被估计

#   survreg's scale  =    1/(rweibull shape)
#   survreg's intercept = log(rweibull scale) 
#   survreg结果中输出的scale与“rweibull scale”不同

 

control 
控制值列表，参考survreg.control() 

 

survreg.control
survreg.control(maxiter=30, rel.tolerance=1e-09,
toler.chol=1e-10, iter.max, debug=0, outer.max=10)

 

maxiter 
最大迭代次数
rel.tolerance 
“相对容忍度”来声明收敛
toler.chol 
Tolerance to declare Cholesky decomposition singular
iter.max 
与maxiter相同
debug 
打印调试信息
outer.max 
用于选择惩罚参数的外部迭代的最大数目

 

convergence tolerance

/x(k+1)/=/x(k)/*Tr+Ta
其中：
k 迭代次数；
x(k+1) x的k次迭代计算值；
x(k) x的k次迭代初始值；
Tr 相对误差
Ta绝对误差
// 表示绝对值
因此，从上述公式我们可以得到一个重要结论：设定计算迭代误差的时候，要全面权衡物理量的绝对值大小，同时要衡量收敛迭代值的相对误差，这样才能正确满意的设定自己需要的计算误差！
不能简单的以为绝对误差和相对误差越小越好，也要杜绝tolerance设定时的随意性（按照公式进行合理的设定是一个优秀的过程工程师的基本素质）

 

survreg使用示例

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
> library(survival)
> str(ovarian)
'data.frame':   26 obs. of  6 variables:
 $ futime  : num  59 115 156 421 431 448 464 475 477 563 ...
 $ fustat  : num  1 1 1 0 1 0 1 1 0 1 ...
 $ age     : num  72.3 74.5 66.5 53.4 50.3 ...
 $ resid.ds: num  2 2 2 2 2 1 2 2 2 1 ...
 $ rx      : num  1 1 1 2 1 1 2 2 1 2 ...
 $ ecog.ps : num  1 1 2 1 1 2 2 2 1 2 ...
  
> # 拟合指数模型，以下2个模型是一样的
> survreg(Surv(futime, fustat) ~ ecog.ps + rx, ovarian, dist='weibull',
+         scale=1)
Call:
survreg(formula = Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian,
    dist = "weibull", scale = 1)
Coefficients:
(Intercept)     ecog.ps          rx
  6.9618376  -0.4331347   0.5815027
Scale fixed at 1 # 注意这里
Loglik(model)= -97.2   Loglik(intercept only)= -98
    Chisq= 1.67 on 2 degrees of freedom, p= 0.43
n= 26
> survreg(Surv(futime, fustat) ~ ecog.ps + rx, ovarian,
+         dist="exponential")
Call:
survreg(formula = Surv(futime, fustat) ~ ecog.ps + rx, data = ovarian,
    dist = "exponential")
Coefficients:
(Intercept)     ecog.ps          rx
  6.9618376  -0.4331347   0.5815027
Scale fixed at 1 # 注意这里
Loglik(model)= -97.2   Loglik(intercept only)= -98
    Chisq= 1.67 on 2 degrees of freedom, p= 0.43
n= 26
> ##########################################
>
> str(lung)
'data.frame':   228 obs. of  10 variables:
 $ inst     : num  3 3 3 5 1 12 7 11 1 7 ...
 $ time     : num  306 455 1010 210 883 ...
 $ status   : num  2 2 1 2 2 1 2 2 2 2 ...
 $ age      : num  74 68 56 57 60 74 68 71 53 61 ...
 $ sex      : num  1 1 1 1 1 1 2 2 1 1 ...
 $ ph.ecog  : num  1 0 0 1 0 1 2 2 1 2 ...
 $ ph.karno : num  90 90 90 90 100 50 70 60 70 70 ...
 $ pat.karno: num  100 90 90 60 90 80 60 80 80 70 ...
 $ meal.cal : num  1175 1225 NA 1150 NA ...
 $ wt.loss  : num  NA 15 15 11 0 0 10 1 16 34 ...
> # weibull分布,如果设置scale<=0，模型将使用最大似然估计，估计最优的scale
> myfit<-survreg(Surv(time, status) ~ ph.ecog + age + sex, data=lung,dist = "weibull")
> myfit
Call:
survreg(formula = Surv(time, status) ~ ph.ecog + age + sex, data = lung,
    dist = "weibull")
Coefficients:
 (Intercept)      ph.ecog          age          sex
 6.273435252 -0.339638098 -0.007475439  0.401090541
Scale= 0.731109
Loglik(model)= -1132.4   Loglik(intercept only)= -1147.4
    Chisq= 29.98 on 3 degrees of freedom, p= 1.4e-06
n=227 (1 observation deleted due to missingness)
> summary(myfit)
Call:
survreg(formula = Surv(time, status) ~ ph.ecog + age + sex, data = lung,
    dist = "weibull")
               Value Std. Error     z        p
(Intercept)  6.27344    0.45358 13.83 1.66e-43
ph.ecog     -0.33964    0.08348 -4.07 4.73e-05
age         -0.00748    0.00676 -1.11 2.69e-01
sex          0.40109    0.12373  3.24 1.19e-03
Log(scale)  -0.31319    0.06135 -5.11 3.30e-07
Scale= 0.731
Weibull distribution
Loglik(model)= -1132.4   Loglik(intercept only)= -1147.4
    Chisq= 29.98 on 3 degrees of freedom, p= 1.4e-06
Number of Newton-Raphson Iterations: 5
n=227 (1 observation deleted due to missingness)
  
 #   survreg's scale  =    1/(rweibull shape)
 #   survreg's intercept = log(rweibull scale)
 # survreg结果中输出的scale与“rweibull scale”不同
  
##############################################
> # weibull分布，以sex对数据分组，产生2个不同的scale
> myfit1<-survreg(Surv(time, status) ~ ph.ecog + age + strata(sex), data=lung,dist = "weibull")
> summary(myfit1)
Call:
survreg(formula = Surv(time, status) ~ ph.ecog + age + strata(sex),
    data = lung, dist = "weibull")
               Value Std. Error      z        p
(Intercept)  6.73235    0.42396 15.880 8.75e-57
ph.ecog     -0.32443    0.08649 -3.751 1.76e-04
age         -0.00581    0.00693 -0.838 4.02e-01
sex=1       -0.24408    0.07920 -3.082 2.06e-03
sex=2       -0.42345    0.10669 -3.969 7.22e-05
Scale:
sex=1 sex=2
0.783 0.655
Weibull distribution
Loglik(model)= -1137.3   Loglik(intercept only)= -1146.2
    Chisq= 17.8 on 2 degrees of freedom, p= 0.00014
Number of Newton-Raphson Iterations: 5
n=227 (1 observation deleted due to missingness)
################################################
# 具有高斯误差的线性回归，以及左删失数据
> survreg(Surv(durable, durable>0, type='left') ~ age + quant,
+         data=tobin, dist='gaussian', scale = 0)
Call:
survreg(formula = Surv(durable, durable > 0, type = "left") ~
    age + quant, data = tobin, dist = "gaussian", scale = 0)
Coefficients:
(Intercept)         age       quant
15.14486636 -0.12905928 -0.04554166
Scale= 5.57254
Loglik(model)= -28.9   Loglik(intercept only)= -29.5
    Chisq= 1.1 on 2 degrees of freedom, p= 0.58
n= 20

 

3.survdiff
测试在两条或更多条生存曲线之间是否存在差异
survdiff(formula, data, subset, na.action, rho=0)
rho = 0 表示log-rank or Mantel-Haenszel检验；
rho = 1 表示Wilcoxon检验

log-rank检验(时序检验)--生存时间分布近似weibull分布或者属于比例风险模型时效率较高
似然比检验(likelihood ratio test)--生存时间分布近似指数分布时效率较高
wilcoxon检验--生存时间分布近似对数正态分布效率较高

survdiff使用示例

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
> survdiff(Surv(futime, fustat) ~ rx,data=ovarian)
Call:
survdiff(formula = Surv(futime, fustat) ~ rx, data = ovarian)
      N Observed Expected (O-E)^2/E (O-E)^2/V
rx=1 13        7     5.23     0.596      1.06
rx=2 13        5     6.77     0.461      1.06
 Chisq= 1.1  on 1 degrees of freedom, p= 0.303
> # 用strata来控制协变量的影响
> survdiff(Surv(time, status) ~ pat.karno + strata(inst), data=lung)
Call:
survdiff(formula = Surv(time, status) ~ pat.karno + strata(inst),
    data = lung)
n=224, 4 observations deleted due to missingness.
               N Observed Expected (O-E)^2/E (O-E)^2/V
pat.karno=30   2        1    0.692   0.13720   0.15752
pat.karno=40   2        1    1.099   0.00889   0.00973
pat.karno=50   4        4    1.166   6.88314   7.45359
pat.karno=60  30       27   16.298   7.02790   9.57333
pat.karno=70  41       31   26.358   0.81742   1.14774
pat.karno=80  50       38   41.938   0.36978   0.60032
pat.karno=90  60       38   47.242   1.80800   3.23078
pat.karno=100 35       21   26.207   1.03451   1.44067
 Chisq= 21.4  on 7 degrees of freedom, p= 0.00326
> survdiff(Surv(time, status) ~ pat.karno + sex, data=lung)
Call:
survdiff(formula = Surv(time, status) ~ pat.karno + sex, data = lung)
n=225, 3 observations deleted due to missingness.
                      N Observed Expected (O-E)^2/E (O-E)^2/V
pat.karno=30, sex=1   1        1    0.246  2.31e+00  2.33e+00
pat.karno=30, sex=2   1        0    0.412  4.12e-01  4.16e-01
pat.karno=40, sex=1   1        1    0.132  5.68e+00  5.72e+00
pat.karno=40, sex=2   1        0    1.204  1.20e+00  1.22e+00
pat.karno=50, sex=1   2        2    0.111  3.23e+01  3.25e+01
pat.karno=50, sex=2   2        2    0.968  1.10e+00  1.11e+00
pat.karno=60, sex=1  18       17    9.173  6.68e+00  7.14e+00
pat.karno=60, sex=2  12       10    6.064  2.56e+00  2.68e+00
pat.karno=70, sex=1  30       23   16.737  2.34e+00  2.65e+00
pat.karno=70, sex=2  11        8    9.527  2.45e-01  2.62e-01
pat.karno=80, sex=1  32       26   24.570  8.32e-02  9.91e-02
pat.karno=80, sex=2  19       13   16.311  6.72e-01  7.55e-01
pat.karno=90, sex=1  31       25   24.992  2.66e-06  3.21e-06
pat.karno=90, sex=2  29       13   24.420  5.34e+00  6.33e+00
pat.karno=100, sex=1 21       15   14.002  7.11e-02  7.91e-02
pat.karno=100, sex=2 14        6   13.131  3.87e+00  4.25e+00
 Chisq= 66.1  on 15 degrees of freedom, p= 2.17e-08

 

anova
计算一个或多个拟合模型的方差（或偏差）表的分析

1
2
3
4
5
6
7
8
9
10
11
12
mod1<-survreg(Surv(time, status) ~ ph.ecog + age + sex, data=lung,dist = "weibull")
mod2<-survreg(Surv(time, status) ~ ph.ecog + age + sex+ph.karno*pat.karno, data=lung,dist = "weibull")
    
anova(mod1,mod2,test="Chisq")
                                       Terms Resid. Df    -2*LL Test Df Deviance     Pr(>Chi)
1                        ph.ecog + age + sex       222 2264.877      NA       NA           NA
2 ph.ecog + age + sex + ph.karno * pat.karno       215 2209.372    =  7 55.50527 1.183848e-09
   
library(MASS)
# 简洁模型更好
stepAIC(mod1)
stepAIC(mod2)

 

    本站是提供个人知识管理的网络存储空间，所有内容均由用户发布，不代表本站观点。请注意甄别内容中的联系方式、诱导购买等信息，谨防诈骗。如发现有害或侵权内容，请点击一键举报。
    转藏 分享 献花（0）

    来自： 勤悦轩 > 《网络统计学习》

    举报/认领

    0条评论

    发表

    请遵守用户 评论公约

    类似文章 更多