stata生存分析：stset命令详解

发布时间：2020-12-21 阅读 4943
分享到微信 分享到微博 收藏文章

 

Stata连享会   主页 || 视频 || 推文 || 知乎

温馨提示： 定期 清理浏览器缓存，可以获得最佳浏览体验。

New！ lianxh 命令发布了：
随时搜索连享会推文、Stata 资源，安装命令如下：
  . ssc install lianxh
使用详情参见帮助文件 (有惊喜)：
  . help lianxh

 

✌ 课程详情： https://gitee.com/arlionn/Course   |   lianxh.cn

 

✨ 课程主页： https://gitee.com/arlionn/Course

 

 

作者： 郝立（复旦大学）
邮箱： haolihayley@outlook.com

 

目录

 

许多新手在学习生存分析时，都苦恼于数据结构的设定问题。今天，我们就和大家分享一下如何处理生存分析数据结构以及 stset 命令各种选项的使用。

 

1. 生存分析数据结构处理

生存分析数据结构既可以用截面数据，也可以用面板数据。本推文中以因变量为初婚风险为例给大家示范截面数据与面板数据分别应该如何处理。因变量是个体进入初婚的风险率，是指风险时间里个体发生初婚事件的概率。

1.1 截面数据处理

首先，下图为本例初始截面数据结构。

. use data1,clear
. dataex birthy gender marriage marriageyear

---------- copy starting from the next line -------------
[CODE]
* Example generated by -dataex-. 
* To install: ssc install dataex
clear
input int birthy byte(gender marriage) int marriageyear
1989 0 1 2009
1986 0 1 2007
1973 1 1 1999
1984 0 1 2015
1985 0 1 2009
1972 1 1 2000
1972 0 1 1995
1963 1 1 1983
1975 1 1 2008
1988 1 0    .
1965 0 0    .
end
[/CODE]
---------- copy up to and including the previous line ---

Listed 11 out of 11 observations

其次，人口学中一般把 50 岁时仍没有结婚的人口界定为终身不婚群体，初婚年龄在 50 岁以上的概率非常小，因此50岁截止时绝大多数样本的初婚年龄已经完全可以观察到，因此可以删除 50 岁以上的样本。

第三，我们需要设定风险开始时间和风险结束时间。本推文中风险集的起始风险时间设定为 15 岁，也即假设一个人从 15 岁开始进入初婚风险期。风险的终止时间即样本发生初婚事件时的时间，直至调查截止时间（2017 年）仍未进入婚姻的样本为删失样本。对于事件发生者，风险时间是从起始风险年龄到初婚年龄期间所有的时间。对于事件未发生者，风险时间则是从起始风险年龄到风险终止年龄之间的所有时间。

第四，我们要产生一个风险持续时间变量（如：duration）

第五，用 stset 命令设置生存分析数据结构，这样我们的截面数据结构就处理好了。

. gen age=2017-birthy 
. drop if age>50 //删除年龄大于50的样本
. gen id=_n //给每个样本产生一个id（此步可以省略）
. order id birthy
. rename marriage event //是否初婚（此步可以省略）
. tab event  //定义事件
      event |      Freq.     Percent        Cum.
------------+-----------------------------------
          0 |          1       11.11       11.11
          1 |          8       88.89      100.00
------------+-----------------------------------
      Total |          9      100.00

. gen date0=birthy+16 //进入风险年份

. sum date0
Variable | Obs      Mean  Std. Dev.    Min    Max
---------+---------------------------------------
   date0 |   9  1996.444  7.264832    1988   2005

*-将结婚年份为缺失值的替换为调查年份2017年，也即已结婚的样本，结婚年份即为风险结束年份；
*-未结婚的样本，调查年份即为风险结束年份
. replace marriageyear=2017 if marriageyear==. 
(1 real change made)

//产生一个风险持续时间变量，以风险结束时间减去风险开始时间
. gen duration=marriageyear-date0 

. stset duration, failure(event)

     failure event:  event != 0 & event < .
obs. time interval:  (0, duration]
 exit on or before:  failure

------------------------------------------------------
  9  total observations
  0  exclusions
------------------------------------------------------
  9  observations remaining, representing
  8  failures in single-record/single-failure data
 91  total analysis time at risk and under observation
                                 at risk from t =    0
                      earliest observed entry t =    0
                           last observed exit t =   17

. 
end of do-file

下图是经过以上步骤处理之后的数据结构

. dataex id birthy gender event marriageyear ///
  age date0 duration _st _d _t _t0

---------- copy starting from the next line -------------
[CODE]
* Example generated by -dataex-.
* To install: ssc install dataex
clear
input float id int birthy byte(gender event) ///
      int marriageyear float(age date0 duration) ///
	  byte(_st _d _t _t0)
1 1989 0 1 2009 28 2005  4 1 1  4 0
2 1986 0 1 2007 31 2002  5 1 1  5 0
3 1973 1 1 1999 44 1989 10 1 1 10 0
4 1984 0 1 2015 33 2000 15 1 1 15 0
5 1985 0 1 2009 32 2001  8 1 1  8 0
6 1972 1 1 2000 45 1988 12 1 1 12 0
7 1972 0 1 1995 45 1988  7 1 1  7 0
8 1975 1 1 2008 42 1991 17 1 1 17 0
9 1988 1 0 2017 29 2004 13 1 0 13 0
end
[/CODE]
---------- copy up to and including the previous line ---

Listed 9 out of 9 observations

1.2 面板数据处理

假设现在有 2010 年以及 2014 年两个截面数据，首先需要分别将两个截面数据中时变变量（如：受教育程度 edu1，edu2； 是否初婚 firstmarr1，firstmarr2；调研年份 surveyyear1，surveyyear2 ）保存一个数据集，再用 merge 命令将两个数据集合并成类似下图的数据集。（这一步本推文省略）

. use date2,clear
. dataex gender birthy firstmarrmy1 edu1 edu2   ///
         firstmarr1 firstmarr2 surveyyear1 surveyyear2

---------- copy starting from the next line -------------
[CODE]
* Example generated by -dataex-.
* To install: ssc install dataex
clear
input byte gender int(birthy firstmarrmy1)   ///
      float edu1 byte(edu2 firstmarr1 firstmarr2) ///
	  int(surveyyear1 surveyyear2)
1 1971 1996   12 12 1 . 2010 2014
0 1975 1996    9  9 1 . 2010 2014
1 1962 1990 10.5 12 1 . 2010 2014
1 1960 1989 10.5 11 1 . 2010 2014
0 1962 1989  7.5  9 1 . 2010 2014
1 1957 1983   12 12 1 . 2010 2014
1 1994    .   11 14 0 0 2010 2014
end
[/CODE]
---------- copy up to and including the previous line ---

Listed 7 out of 7 observations

首先给每个样本产生一个 id，然后用 reshape 命令将宽数据变为长数据。如下图：

. use date2,clear
. gen id=_n //给每个样本产生一个id
. order id birthy

. reshape long firstmarr edu surveyyear, i(id) //将宽数据变为长数据
(note: j = 1 2)

Data                       wide   ->   long
-------------------------------------------------
Number of obs.                7   ->      14
Number of variables          10   ->       8
j variable (2 values)             ->   _j
xij variables:
          firstmarr1 firstmarr2   ->   firstmarr
                      edu1 edu2   ->   edu
        surveyyear1 surveyyear2   ->   surveyyear
-------------------------------------------------

使用 reshape 命令将宽数据变成长数据之后，数据结构变成下图的样子：

. dataex id _j birthy gender firstmarrmy1 edu firstmarr surveyyear

---------- copy starting from the next line -------------
[CODE]
* Example generated by -dataex-.
* To install: ssc install dataex
clear
input float id byte _j int birthy byte gender ///
      int firstmarrmy1 float edu ///
	  byte firstmarr int surveyyear
1 1 1971 1 1996   12 1 2010
1 2 1971 1 1996   12 . 2014
2 1 1975 0 1996    9 1 2010
2 2 1975 0 1996    9 . 2014
3 1 1962 1 1990 10.5 1 2010
3 2 1962 1 1990   12 . 2014
4 1 1960 1 1989 10.5 1 2010
4 2 1960 1 1989   11 . 2014
5 1 1962 0 1989  7.5 1 2010
5 2 1962 0 1989    9 . 2014
6 1 1957 1 1983   12 1 2010
6 2 1957 1 1983   12 . 2014
7 1 1994 1    .   11 0 2010
7 2 1994 1    .   14 0 2014
end
[/CODE]
---------- copy up to and including the previous line ---

Listed 14 out of 14 observations

由于 2010 年调查时，已经有一些人初婚了，初婚状态变量和初婚时间变量均已经有了，所以她们的第二条记录可以删除。 再删除年龄大于 50 且未婚的样本。

*-如果2010年调查的时候已经初婚了，可以删除2014年初婚状态firstmarr为缺失值的样本
. drop if missing(firstmarr) 
. gen age=surveyyear-birthy 
. drop if age>50 & firstmarr==0 //删除年龄大于50且未婚的样本

删除缺失值和年龄大于 50 且未婚的样本后，数据结构如下图：

. dataex id _j birthy gender firstmarrmy1 edu firstmarr surveyyear age

---------- copy starting from the next line -------------
[CODE]
* Example generated by -dataex-.
* To install: ssc install dataex
clear
input float id byte _j int birthy byte gender ///
      int firstmarrmy1 float edu byte firstmarr ///
	  int surveyyear float age
1 1 1971 1 1996   12 1 2010 39
2 1 1975 0 1996    9 1 2010 35
3 1 1962 1 1990 10.5 1 2010 48
4 1 1960 1 1989 10.5 1 2010 50
5 1 1962 0 1989  7.5 1 2010 48
6 1 1957 1 1983   12 1 2010 53
7 1 1994 1    .   11 0 2010 16
7 2 1994 1    .   14 0 2014 20
end
[/CODE]
---------- copy up to and including the previous line ---

Listed 8 out of 8 observations

开始定义事件，这步可省略 开始设定进入初婚风险年份，也即出生年份 +16。 设定风险结束年份，如果调查年份还未婚，将其风险结束年份替换为调查年份。 再用 stset 命令设定生存分析数据结构。

. rename firstmarr event //是否初婚,定义事件
. gen date0=birthy+16    //进入风险年份
. rename firstmarrmy1 date1 //设置风险结束年份
. replace date1=surveyyear if date1==. 

. stset date1, origin(time date0) id(id) ///
        failure(event==1) time0(date0)

                id:  id
     failure event:  event == 1
obs. time interval:  (date0, date1]
 exit on or before:  failure
    t for analysis:  (time-origin)
            origin:  time date0

--------------------------------------------------------
  8  total observation
  1  entry on or after exit (date0>date1) PROBABLE ERROR
--------------------------------------------------------
  7  observations remaining, representing
  7  subjects
  6  failures in single-failure-per-subject data
 64  total analysis time at risk and under observation
                                 at risk from t   =  0
                      earliest observed entry t   =  0
                           last observed exit t   = 13

如果出现错误，按照错误提示，需要检查数据问题。这里报错的是 date0>date1，因此可以删除 date0 大于或等于 date1 的样本，再次执行 stset命令可以看到没有再报错的信息。此外，一个小 tips， 如果一个样本有多条观测值的话，第二条观测值开始观测的时间 date0 要与上一条的观测结束时间 date1 要保持接续，后面每条观测时间均与前述情况一样，若不是连续可以用 snapspan 命令进行调整。

. drop if date0==date1 //如果报错，可以按照提示，删除date0大于或者等于date1的样本

. stset date1, origin(time date0) id(id) ///
        failure(event==1) time0(date0) //再次执行一次检验

                id:  id
     failure event:  event == 1
obs. time interval:  (date0, date1]
 exit on or before:  failure
    t for analysis:  (time-origin)
            origin:  time date0

-----------------------------------------------------
  7  total observation
  0  exclusion
-----------------------------------------------------
  7  observations remaining, representing
  7  subjects
  6  failures in single-failure-per-subject data
 64  total analysis time at risk and under observation
                               at risk from t =     0
                    earliest observed entry t =     0
                         last observed exit t =    13

 

2. stset 命令解析

在做生存分析之前，首先要告诉 Stata 我们的数据是生存分析的数据，因此需要用到 stset 命令。详情参见 help stset ( [ST] stset )。 在运行 stset 命令之后，Stata 将自动生成四个新变量，分别是 _t0 、 _st 、_d 与 _t 。其中 _t0 表示起始时间点，默认为 0，但如果有左截断，则 _t0 不为 0 ；_st 表示数据的状态是否完好可用；_d 表示事件的结果，是否发生；_t 表示每个样本事件持续的时间。

2.1 stset 最简形式与 failure 选项使用

第一种： 最简命令的格式中只需指定失效时间也即事件持续时间变量 (failtime) 即可：

. stset  failtime

第二种： 加入 failure 选项，命令格式为：

. stset  failtime,  failure(failed)

此处，新增的 failure() 选项用于指定「事件变量」。以初婚风险为例，duration 为事件持续时间，event 为初婚发生事件，则命令为：

. stset duration, failure (event)

2.2 id 选项使用

id 是每一个样本的标识，在生存分析的数据中，有时候一个样本有几行观测值，这时候就需要 id 来帮助识别样本。如果未指定 id，则假定每个观测值代表一个不同的样本，因此构成了一个观测对象一条记录的生存数据集。 当指定 id 时，可以区分一个样本多条观测记录。 当在 stset 命令中使用了 id 选项，它会自动识别出样本中的一些逻辑错误并报错，比如观测时间的重复等。这样我们就知道错在哪里，可以根据报错信息再对样本进行检查。示例见第二节面板数据处理中最后一个图，再次不再赘述。

加入 id 选项，其命令格式为：

. stset  failtime,  failure(failed)  id(id)

2.3 time0 选项使用

time0 是每段观测时期的起始时点。如果没有用 time0 选项，就默认所有 time0 都是从样本第一行的第一个起始点起算。但是如果一个样本有多行观测值，而后面各行的观测值中起始时间点有缺失值，那么就会默认缺失值是第一行的第一个起始时间点，这样就会出现问题。如果用了 time0 选项的话，就不会每一行都以第一行的起始时间点作为起算点，而是以每一行各自开始的时间点为起算。 time0 选项也主要用一个样本多条记录的数据格式中。

2.4 origin 选项使用

origin 是说明一个样本从什么时候开始事件可能发生，也即样本可能面临发生风险的时间。比如说：假设一个女性从 20 岁开始可以登记结婚。 以 Stata 中自有数据集 diet2 为例，加入 origin 选项，其命令格式为：

stset dox,  id(id) failure(fail)  origin(time dob)

其中，dox 是事件退出日期的变量；failure() 选项用于定义事件，括号里面的 fail 变量为事件变量origin() 选项括号内的 time 为固定格式标记，常规的用法是 origin(time exp)。本例中，origin(time dob) 表示事件开始的时间变量为 dob。

2.5 enter 选项与 exit 选项的使用

enter 主要说明一个样本从什么时候进入研究；exit 主要说明一个样本从什么时候离开研究。 以 Stata 中自有数据集 diet2.dta 为例，加入 enter 选项，其命令格式为：

stset dox,  id(id) failure(fail)  origin(time dob)  enter(time doe)

其中，dox 是事件退出日期的变量；failure() 中填入用以标记事件是否发生的变量，本例中为 fail。origin 的括号内 time 为命令，而 dob 是事件可能开始面临风险的时间变量；enter 括号内 time 为命令，而 doe 是样本进入研究的时间。

. webuse diet2
(Diet data with dates split by age)

. stset dox, id(id) failure(fail)  ///
        origin(time dob)  enter(time doe)

                id:  id
     failure event:  fail != 0 & fail < .
obs. time interval:  (dox[_n-1], dox]
 enter on or after:  time doe
 exit on or before:  failure
    t for analysis:  (time-origin)
            origin:  time dob

--------------------------------------------------------------
        755  total observations
          0  exclusions
--------------------------------------------------------------
        755  observations remaining, representing
        337  subjects
         80  failures in single-failure-per-subject data
  1,681,490  total analysis time at risk and under observation
                                    at risk from t =         0
                         earliest observed entry t =    10,985
                              last observed exit t =    25,567

以 Stata 中自有数据集 diet.dta 为例，加入 exit 选项，其命令格式为：

. stset dox,  failure(fail) origin(time dob)  ///
        enter(time doe)  exit(time mdy(12,1,1970))

其中，新增选项 exit() 的含义是：time 为固定格式标记，而 mdy(12,1,1970) 是样本退出研究的时间，说明了样本的 1970 年 12 月 1 日之后会退出分析的风险集。

. webuse diet
(Diet data with dates)

. stset dox, failure(fail) origin(time dob)  ///
        enter(time doe)  exit(time mdy(12,1,1970))

     failure event:  fail != 0 & fail < .
obs. time interval:  (origin, dox]
 enter on or after:  time doe
 exit on or before:  time mdy(12,1,1970)
    t for analysis:  (time-origin)
            origin:  time dob

--------------------------------------------------------------
        337  total observations
          0  exclusions
--------------------------------------------------------------
        337  observations remaining, representing
         54  failures in single-record/single-failure data
  1,166,590  total analysis time at risk and under observation
                                    at risk from t =         0
                         earliest observed entry t =    10,985
                              last observed exit t =    25,567

.

 

3. 参考文献

 

相关课程

 

连享会-直播课 上线了！
http://lianxh.duanshu.com

免费公开课：

 

课程一览

支持回看

专题 嘉宾 直播/回看视频
⭐ 最新专题 因果推断, 空间计量，寒暑假班等
⭕ 数据清洗系列 游万海 直播, 88 元，已上线
研究设计 连玉君 我的特斯拉-实证研究设计，-幻灯片-
面板模型 连玉君 动态面板模型，-幻灯片-
面板模型 连玉君 直击面板数据模型 [免费公开课，2小时]

Note: 部分课程的资料，PPT 等可以前往 连享会-直播课 主页查看，下载。

 

 

关于我们

  • Stata连享会 由中山大学连玉君老师团队创办，定期分享实证分析经验。直播间 有很多视频课程，可以随时观看。
  • 连享会-主页 和 知乎专栏，300+ 推文，实证分析不再抓狂。
  • 公众号关键词搜索/回复 功能已经上线。大家可以在公众号左下角点击键盘图标，输入简要关键词，以便快速呈现历史推文，获取工具软件和数据下载。常见关键词：课程, 直播, 视频, 客服, 模型设定, 研究设计, stata, plus, 绘图, 编程, 面板, 论文重现, 可视化, RDD, DID, PSM, 合成控制法 等
连享会主页 lianxh.cn
连享会主页 lianxh.cn

 

连享会小程序：扫一扫，看推文，看视频……

 

扫码加入连享会微信群，提问交流更方便

✏ 连享会学习群-常见问题解答汇总：
✨ https://gitee.com/arlionn/WD

 

New！ lianxh 命令发布了：
随时搜索连享会推文、Stata 资源，安装命令如下：
  . ssc install lianxh
使用详情参见帮助文件 (有惊喜)：
  . help lianxh