前往小程序,Get更优阅读体验!
立即前往
首页
学习
活动
专区
工具
TVP
发布
社区首页 >专栏 >MADlib——基于SQL的数据挖掘解决方案(17)——回归之Cox比例风险回归

MADlib——基于SQL的数据挖掘解决方案(17)——回归之Cox比例风险回归

作者头像
用户1148526
发布2019-05-25 19:36:20
1.1K0
发布2019-05-25 19:36:20
举报
文章被收录于专栏:Hadoop数据仓库

版权声明:本文为博主原创文章,未经博主允许不得转载。 https://cloud.tencent.com/developer/article/1433073

一、Cox比例风险回归简介

代码语言:txt
复制
    Cox比例风险回归模型(Cox’s proportional hazards regression model),简称Cox回归模型,由英国统计学家D.R.Cox于1972年提出,主要用于肿瘤和其它慢性病的预后分析,也可用于队列研究的病因探索。

1. 基本概念

  • 生存函数:又称累计生存率,简称生存率,常用S(t)表示。代表被观察对象的生存时间大于t时刻的概率,实际中t时刻的取值估算公式为:S(t) ≈ 生存时间长于t的患者数 / 患者总数
  • 死亡概率函数:简称为死亡概率,常用F(t)表示。代表一个被观察对象从开始观察到时间t为止的死亡概率,它与生存函数的关系为:F(t) = 1 - S(t)。
  • 死亡密度函数:简称为密度函数,用f(t)表示。代表所有被观察对象在t时刻的瞬时死亡率。
  • 风险函数:即生存时间已达到t的一群被观察对象在t时刻的瞬时死亡率,用h(t)表示。代表已存活到时间t的每个观察对象从t到t+∆t这一非常小的区间内死亡的概率极限,它与生存函数、死亡密度函数的关系为:h(t) = f(t) / S(t)。

2. Cox回归模型结构

代码语言:txt
复制
    Cox回归模型不直接考察生存函数与协变量(影响因素)的关系,而是用风险函数作为因变量。设有n名病人(i=1,2,...,n),第i名病人的生存时间为 
代码语言:txt
复制
                         ,同时该病人具有一组协变量 

,则模型为:

:在时间t处的风险函数。

:基准风险函数,为所有协变量取零时t时刻的风险函数,即没有协变量下的风险函数。这是模型中的非参数部分,因此Cox回归是一种半参数分析方法。

:协变量。

:根据观察值估算出的回归系数。

:预后指数。大于0时表示该病人对应的危险度大于平均水平;等于0时为达到平均水平;小于0时表示该病人的危险度小于平均水平。

  • 回归系数

时,协变量的取值越大,风险函数

的值越大,表示病人死亡的风险越大。

  • 回归系数

时,表示协变量对风险函数

没有影响。

  • 回归系数

时,协变量的取值越大,风险函数

的值越小,表示病人死亡的风险越小。

3. Cox回归模型的前提条件

代码语言:txt
复制
    Cox回归模型必须满足比例风险假设(Proportional Hazards Assumption,PHA):

(1)任何两个个体的风险函数之比,即风险比(HazardRatio,HR)保持一个恒定的比例,而与时间t无关。

(2)模型中协变量的效应不随时间改变而改变。

代码语言:txt
复制
    检查某协变量是否满足PHA,最简单的方法是观察该变量分组的生存曲线。若生存曲线交叉,表示不满足PHA,此时可采用分层比例风险模型。

4. 参数估算与假设检验

代码语言:txt
复制
    Cox回归的参数估计同逻辑回归分析一样采用最大似然估计法。其基本思想是先建立偏似然函数或对数偏似然函数,求偏似然函数或对数偏似然函数达到极大时参数的取值,即为参数的最大似然估计值。
代码语言:txt
复制
     假设检验的方法有时协变量法、线性相关检验法、加权残差Score法等。这三种检验法有较高的准确率,且三种方法的检验效能相近。MADlib的Cox模型PHA检验函数使用线性相关检验法实现。

5. Cox模型的注意事项

  • 研究的协变量在被研究对象中的分布要适中,否则会给回归参数的估计带来困难。
  • Cox模型应用较灵活,被观察对象进入研究队列的早晚、时间长短可以不一致,但如果研究的变量随时间而变化,可以采用时依协变量模型进行分析。
  • 样本含量不宜过小,一般应在40例以上,经验上的做法是样本含量为协变量个数的5-20倍。
  • 尽管可以分析删失数据,但在观察时,要尽量避免观察对象的失访,过多失访容易造成结果的偏倚。
  • Cox模型对异常值较为敏感,所以在进行模型拟合时要注意拟合优度的检验。

二、MADlib中Cox比例风险回归相关函数

1. 训练函数

(1) 语法

代码语言:javascript
复制
coxph_train( source_table,  
             output_table,  
             dependent_variable,  
             independent_variable,  
             right_censoring_status,  
             strata,  
             optimizer_params  
           ) 

(2) 参数

参数名称

数据类型

描述

source_table

VARCHAR

包含训练数据的源表名。

output_table

VARCHAR

保存模型的输出表名,主输出表列和概要输出表列分别如表2、表3所示。

dependent_variable

VARCHAR

因变量名称,指死亡时间,不需要对数据进行预排序。

independent_variable

VARCHAR

自变量名称数组。

right_censoring_status(可选)

VARCHAR

缺省值为TRUE。表示右删失状态的字符串,如果无删失则为TRUE,否则为FALSE。该参数可以包含是右删失状态的列名,或者是一个可以对每个观察值进行评估的布尔表达式,如‘column_name < 10’。

strata(可选)

VARCHAR

缺省值为NULL,不做任何分层,可以是逗号分隔的列名。

optimizer_params(可选)

VARCHAR

缺省值为NULL,此时使用缺省的优化参数:max_iter=100, optimizer=newton, tolerance=1e-8, array_agg_size=10000000, sample_size=1000000。应该是一个逗号分隔的‘key=value’字符串,参数含义如下: l max_iter:最大迭代数。如果迭代次数达到此数值则停止计算,这通常意味着无法收敛。 l optimizer:优化方法,当前仅支持‘newton’。 l tolerance:停止条件。当连续两次迭代的对数似然值之差小于此参数,计算已经收敛并停止。 l array_agg_size:为了加速计算,将原始数据表切分成多个数据片,每片数据聚合成一个大行。计算处理时,整个大行被导入内存以提高运算速度。此参数控制一个大行中包含多少数据,参数值越大速度越快,但由于PostgreSQL数据库的限制,一个大行的大小不能超过1G。 l sample_size:为了将数据切分成大致相等的片,先进行数据采样,然后利用采样数据找出断点。该参数值越大,产生的断点越精确。

表1 coxph_train函数参数说明

列名

数据类型

描述

Coef

FLOAT8[]

回归系数向量。

loglikelihood

FLOAT8

极大似然估计的对数似然值。

std_err

FLOAT8[]

回归系数标准差向量。

stats

FLOAT8[]

回归系数统计向量。

p_values

FLOAT8[]

回归系数p值向量。

hessian

FLOAT8[]

Hessian矩阵。

num_iterations

INTEGER

优化器执行的迭代次数。

表2 coxph_train函数主输出表列说明

代码语言:txt
复制
    训练函数在产生输出表的同时,还会创建一个名为<output\_table>\_summary的概要表,具有以下列:

列名

数据类型

描述

method

VARCHAR

'coxph',描述模型的字符串。

source_table

VARCHAR

源表名。

out_table

VARCHAR

模型表名。

dependent_varname

VARCHAR

因变量表达式。

independent_varname

VARCHAR

自变量表达式。

right_censoring_status

VARCHAR

右删失状态。

Strata

VARCHAR

分层列。

num_processed

INTEGER

处理行数。

num_missing_rows_skipped

INTEGER

跳过行数。

表3 coxph_train函数概要输出表列说明

2. 比例风险假设检验函数

代码语言:txt
复制
    cox\_zph()函数检验Cox回归的比例风险假设,它通过计算coxph\_train()输出模型中残差与时间的相关性验证比例风险假设。

(1) 语法

代码语言:javascript
复制
cox_zph(cox_model_table, output_table)

(2) 参数

  • cox_model_table:TEXT类型,包含Cox模型的表名,应该是coxph_train()函数的输出表。
  • output_table:TEXT类型,保存测试结果的表名,主输出表列和残差输出表列分别如表4、表5所示。

列名

数据类型

描述

covariate

TEXT

协变量。

rho

FLOAT8[]

生存时间与Schoenfeld残差相关系数向量。

chi_square

FLOAT8[]

相关分析卡方检验统计量。

p_value

FLOAT8[]

卡方统计双尾p值。

表4 cox_zph函数输出表列说明

代码语言:txt
复制
    检验函数还创建一个保存残差值的输出表,名为<output\_table>\_residual,具有以下列:

列名

数据类型

描述

<dep_column_name>

FLOAT8

原始源表中存在的时间值(因变量)。

residual

FLOAT8[]

原始协变量与coxph_train模型中的期望协变量之差。

scaled_residual

FLOAT8[]

由系数的方差来衡量的残差值。

表5 cox_zph函数残差输出表列说明

3. 预测函数

(1) 语法

代码语言:javascript
复制
coxph_predict(model_table,  
              source_table,  
              id_col_name,  
              output_table,  
              pred_type,  
              reference) 

(2) 参数

参数名称

数据类型

描述

model_table

TEXT

包含cox模型的表名。

source_table

TEXT

包含预测数据的表名。

id_col_name

TEXT

id列名。

output_table

TEXT

存储预测结果的输出表名,输出表具有以下列: l id:TEXT类型,id列。 l predicted_result:FLOAT8类型,基于预测类型参数值的预测结果。

pred_type(可选)

TEXT

预测类型。缺省为‘linear_predictors’,可以是‘linear_predictors’、‘risk’或‘terms’。 l linear_predictors:计算自变量和系数的点积,也即预后指数。 l risk:预后指数预测结果的指数值。 l terms:每个协变量与它们相应的系数值的乘积,不求和。

reference(可选)

TEXT

用于中心预测。可以为‘strata’、‘overall’,缺省值为‘strata’。

表6 coxph_predict函数参数说明

代码语言:txt
复制
    注:Cox回归模型的因变量是风险函数,因此与其它模型的预测函数不同,它不直接返回生存时间的预测值。

三、示例

1. 查看联机帮助

代码语言:javascript
复制
select madlib.coxph_train();  
-- 训练函数用法  
select madlib.coxph_train('usage');  
-- 示例  
select madlib.coxph_train('example'); 

2. 建立测试数据表并装载原始数据

代码语言:javascript
复制
drop table if exists sample_data;  
create table sample_data (  
    id integer not null,    -- 逻辑主键  
    grp double precision,   -- 分组(治疗方法)  
    wbc double precision,   -- 白细胞值  
    timedeath integer,      -- 生存天数  
    status boolean          -- 结局  
);  
copy sample_data from stdin with delimiter '|';  
  1 |   0 | 1.45 |        35 | t  
  2 |   0 | 1.47 |        34 | t  
  3 |   0 |  2.2 |        32 | t  
  4 |   0 | 1.78 |        25 | t  
  5 |   0 | 2.57 |        23 | t  
  6 |   0 | 2.32 |        22 | t  
  7 |   0 | 2.01 |        20 | t  
  8 |   0 | 2.05 |        19 | t  
  9 |   0 | 2.16 |        17 | t  
 10 |   0 |  3.6 |        16 | t  
 11 |   1 |  2.3 |        15 | t  
 12 |   0 | 2.88 |        13 | t  
 13 |   1 |  1.5 |        12 | t  
 14 |   0 |  2.6 |        11 | t  
 15 |   0 |  2.7 |        10 | t  
 16 |   0 |  2.8 |         9 | t  
 17 |   1 | 2.32 |         8 | t  
 18 |   0 | 4.43 |         7 | t  
 19 |   0 | 2.31 |         6 | t  
 20 |   1 | 3.49 |         5 | t  
 21 |   1 | 2.42 |         4 | t  
 22 |   1 | 4.01 |         3 | t  
 23 |   1 | 4.91 |         2 | t  
 24 |   1 |    5 |         1 | t  
\.
代码语言:txt
复制
    说明:在这个假设的生存分析案例中,将24名患者分为两组(如模拟两种治疗方法)进行观察。协变量有两个,分组与白细胞值,样本量是协变量个数的12倍。因变量为生存天数。所有患者结局已知,不存在删失情况。

3. 训练模型

代码语言:javascript
复制
drop table if exists sample_cox, sample_cox_summary;  
select madlib.coxph_train( 'sample_data',  
                           'sample_cox',  
                           'timedeath',  
                           'array[grp,wbc]',  
                           'status'  
                         ); 

4. 查看回归结果

代码语言:javascript
复制
\x on  
select * from sample_cox;
代码语言:txt
复制
    结果:
代码语言:javascript
复制
-[ RECORD 1 ]--+----------------------------------------------------------------------------
coef           | {2.54407073265254,1.67172094779487}
loglikelihood  | -37.8532498733
std_err        | {0.677180599294897,0.387195514577534}
z_stats        | {3.7568570855419,4.31751114064138}
p_values       | {0.000172060691513886,1.5779844638453e-05}
hessian        | {{2.78043065745617,-2.25848560642414},{-2.25848560642414,8.50472838284472}}
num_iterations | 5
代码语言:txt
复制
    经过5次迭代后求得两个协变量的回归系数。z检验得到的p值很小,说明样本间的差异由抽样误差所致的概率极小,具有显著统计意义。

5. 检验所得模型的PHA

代码语言:javascript
复制
drop table if exists sample_zph_output, sample_zph_output_residual;  
select madlib.cox_zph( 'sample_cox',  
                       'sample_zph_output'  
                     );  
  
select * from sample_zph_output; 
代码语言:txt
复制
    结果:
代码语言:javascript
复制
-[ RECORD 1 ]-----------------------------------------  
covariate  | array[grp,wbc]  
rho        | {0.00237308357328641,0.037560056884043}  
chi_square | {0.000100675718191977,0.0232317400546175}  
p_value    | {0.991994376850758,0.878855984657948} 
代码语言:txt
复制
    madlib.cox\_zph函数使用卡方线性相关法检验比例风险假设是否成立。它的原理很简单:如果数据满足比例风险假设,则Schoenfeld残差和生存时间的秩次之间无相关性。计算步骤按以下3步进行:①用未删失数据计算每个协变量Schoenfeld残差;②将未删失的生存时间排序,并以新变量(协变量残差)记录秩次1、2、3...,如出现相同生存时间(结点),则以平均秩次记录。③进行假设检验,检验Schoenfeld残差和生存时间秩次间的相关性(原假设为H0:r=0,无相关性)。若原假设被拒绝,则表明数据不满足比例风险假设。
代码语言:txt
复制
    从本例的检验p值结果看,协变量对应的双尾p值接近于1,说明应该接受原假设,模型满足比例风险假设。

6. 用模型进行预测

代码语言:txt
复制
    本例使用源数据表演示预测。

(1)条目预测

代码语言:javascript
复制
\x off  
drop table if exists sample_pred_terms;  
select madlib.coxph_predict('sample_cox',  
                            'sample_data',  
                            'id',  
                            'sample_pred_terms',  
                            'terms');  
  
select id, terms, sum(c) linear_predictors  
  from (select id, terms, unnest(terms) c  
          from sample_pred_terms) t  
 group by id, terms  
 order by linear_predictors;
代码语言:txt
复制
    结果:
代码语言:javascript
复制
 id |                  terms                   | linear_predictors    
----+------------------------------------------+--------------------  
  1 | {-0.848023577550848,-2.12308560369949}   |  -2.97110918125034  
  2 | {-0.848023577550848,-2.08965118474359}   |  -2.93767476229444  
  4 | {-0.848023577550848,-1.57141769092718}   |  -2.41944126847803  
  7 | {-0.848023577550848,-1.18692187293436}   |  -2.03494545048521  
  8 | {-0.848023577550848,-1.12005303502256}   |  -1.96807661257341  
  9 | {-0.848023577550848,-0.936163730765128}  |  -1.78418730831598  
  3 | {-0.848023577550848,-0.869294892853333}  |  -1.71731847040418  
 19 | {-0.848023577550848,-0.685405588595898}  |  -1.53342916614675  
  6 | {-0.848023577550848,-0.668688379117949}  |   -1.5167119566688  
  5 | {-0.848023577550848,-0.250758142169231}  |  -1.09878171972008  
 14 | {-0.848023577550848,-0.200606513735385}  |  -1.04863009128623  
 15 | {-0.848023577550848,-0.0334344189558975} | -0.881457996506746  
 16 | {-0.848023577550848,0.133737675823589}   | -0.714285901727259  
 12 | {-0.848023577550848,0.267475351647179}   | -0.580548225903669  
 13 | {1.6960471551017,-2.03949955630974}      | -0.343452401208047  
 10 | {-0.848023577550848,1.47111443405949}    |  0.623090856508639  
 11 | {1.6960471551017,-0.702122798073847}     |   0.99392435702785  
 17 | {1.6960471551017,-0.668688379117949}     |   1.02735877598375  
 21 | {1.6960471551017,-0.501516284338462}     |   1.19453087076323  
 18 | {-0.848023577550848,2.85864282072923}    |   2.01061924317838  
 20 | {1.6960471551017,1.28722512980205}       |   2.98327228490375  
 22 | {1.6960471551017,2.15652002265538}       |   3.85256717775708  
 23 | {1.6960471551017,3.66106887567077}       |   5.35711603077247  
 24 | {1.6960471551017,3.81152376097231}       |     5.507570916074  
(24 rows)
代码语言:txt
复制
    条目预测的输出是一个协变量数组,每个元素值是对应协变量与回归系数的点积。

(2)预后指数预测

代码语言:javascript
复制
\x off  
drop table if exists sample_pred;  
select madlib.coxph_predict('sample_cox',  
                            'sample_data',  
                            'id',  
                            'sample_pred');  
  
select id, linear_predictors, exp(linear_predictors) risk  
  from sample_pred order by linear_predictors; 
代码语言:txt
复制
    结果:
代码语言:javascript
复制
 id | linear_predictors  |        risk          
----+--------------------+--------------------  
  1 |  -2.97110918125034 | 0.0512464372091503  
  2 |  -2.93767476229444 |  0.052988797150351  
  4 |  -2.41944126847803 | 0.0889713146524469  
  7 |  -2.03494545048521 |  0.130687611255372  
  8 |  -1.96807661257341 |  0.139725343898993  
  9 |  -1.78418730831598 |  0.167933483703554  
  3 |  -1.71731847040418 |  0.179546963459175  
 19 |  -1.53342916614675 |  0.215794402223054  
  6 |   -1.5167119566688 |  0.219432204682558  
  5 |  -1.09878171972008 |   0.33327686110022  
 14 |  -1.04863009128623 |  0.350417460388247  
 15 | -0.881457996506746 |  0.414178600294289  
 16 | -0.714285901727259 |  0.489541567796517  
 12 | -0.580548225903669 |  0.559591499901242  
 13 | -0.343452401208048 |  0.709317242984782  
 10 |  0.623090856508639 |   1.86468261037507  
 11 |   0.99392435702785 |   2.70181658768287  
 17 |   1.02735877598375 |   2.79367735395697  
 21 |   1.19453087076323 |   3.30200833843654  
 18 |   2.01061924317838 |   7.46794038691976  
 20 |   2.98327228490375 |   19.7523463121038  
 22 |   3.85256717775708 |   47.1138577624205  
 23 |   5.35711603077247 |   212.112338046552  
 24 |     5.507570916074 |   246.551504798816  
(24 rows) 
代码语言:txt
复制
    可以看到,1的风险度最小,24的风险度最大,并且预后指数预测的结果,就是对条目预测结果数组的元素求和而得。

(3)风险预测

代码语言:javascript
复制
drop table if exists sample_pred_risk;  
select madlib.coxph_predict('sample_cox',  
                            'sample_data',  
                            'id',  
                            'sample_pred_risk',  
                            'risk');  
  
select * from sample_pred_risk order by risk;
代码语言:txt
复制
    结果:
代码语言:javascript
复制
 id |        risk          
----+--------------------  
  1 | 0.0512464372091503  
  2 |  0.052988797150351  
  4 | 0.0889713146524469  
  7 |  0.130687611255372  
  8 |  0.139725343898993  
  9 |  0.167933483703554  
  3 |  0.179546963459175  
 19 |  0.215794402223054  
  6 |  0.219432204682558  
  5 |   0.33327686110022  
 14 |  0.350417460388247  
 15 |  0.414178600294289  
 16 |  0.489541567796517  
 12 |  0.559591499901242  
 13 |  0.709317242984782  
 10 |   1.86468261037507  
 11 |   2.70181658768287  
 17 |   2.79367735395697  
 21 |   3.30200833843654  
 18 |   7.46794038691976  
 20 |   19.7523463121038  
 22 |   47.1138577624205  
 23 |   212.112338046552  
 24 |   246.551504798816  
(24 rows)  
代码语言:txt
复制
    可以看到,风险预测就是预后指数结果的指数值。
本文参与 腾讯云自媒体同步曝光计划,分享自作者个人站点/博客。
原始发表:2018年02月06日,如有侵权请联系 cloudcommunity@tencent.com 删除

本文分享自 作者个人站点/博客 前往查看

如有侵权,请联系 cloudcommunity@tencent.com 删除。

本文参与 腾讯云自媒体同步曝光计划  ,欢迎热爱写作的你一起参与!

评论
登录后参与评论
0 条评论
热度
最新
推荐阅读
目录
  • 一、Cox比例风险回归简介
    • 1. 基本概念
      • 2. Cox回归模型结构
        • 3. Cox回归模型的前提条件
          • 4. 参数估算与假设检验
            • 5. Cox模型的注意事项
            • 二、MADlib中Cox比例风险回归相关函数
              • 1. 训练函数
                • 2. 比例风险假设检验函数
                  • 3. 预测函数
                  • 三、示例
                    • 1. 查看联机帮助
                      • 2. 建立测试数据表并装载原始数据
                        • 3. 训练模型
                          • 4. 查看回归结果
                            • 5. 检验所得模型的PHA
                              • 6. 用模型进行预测
                              相关产品与服务
                              云数据库 PostgreSQL
                              腾讯云数据库 PostgreSQL(TencentDB for PostgreSQL,云 API 使用 postgres 作为简称)能够让您在云端轻松设置、操作和扩展目前功能最强大的开源数据库 PostgreSQL。腾讯云将负责绝大部分处理复杂而耗时的管理工作,如 PostgreSQL 软件安装、存储管理、高可用复制、以及为灾难恢复而进行的数据备份,让您更专注于业务程序开发。
                              领券
                              问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档