生存曲线求助

数据见附件，想 作图20210812.xlsx (11.09 KB)

2021-8-13 00:20:26 上传

绘制成类似如下曲线

edragon1983 发表于 2021-8-14 17:52
谢谢大侠！ 还差最后一个关键环节。number censored如何绘制，恳请出手相助！
其他环节通过r与stata也基 ...

已经改进。

1. /*0. 环境设置*/
   
2. option missing = ' ';
   
3. filename RawData "D:\作图20210812.xlsx ";
   
4. 
   
5. /*1. 导入数据*/
   
6. proc import datafile = RawData dbms = excel out = RawData replace;
   
7. run;
   
8. 
   
9. /*2. 执行生存分析*/
   
10. /*2.1 计算各时点的生存率,并执行log-rank检验*/
    
11. ods output ProductLimitEstimates = ProductLimitEstimates SurvivalPlot = SurvivalPlot HomTests = HomTests;
    
12. proc lifetest data = RawData plots=(survival(atrisk=0 to 135 by 15));
    
13.   strata group;
    
14.   time time*_COL1(1);
    
15. run;
    
16. ods output close;
    
17. 
    
18. data SurvivalPlot;
    
19.   set SurvivalPlot;
    
20. 
    
21.   /*获取删失数量*/
    
22.   do i = 1 to rec;
    
23.     set ProductLimitEstimates(keep=STRATUM time Censor rename=(STRATUM=TmpSTRATUM time=Tmptime))nobs=rec point=i;
    
24.     if StratumNum=TmpSTRATUM and time>=Tmptime then CountCensor = sum(CountCensor,Censor);
    
25.   end;
    
26.   if mod(time,15) = 0 then NumRiskAndCensor = cats(AtRisk,'(',CountCensor,')');
    
27.   if Survival then Survival = Survival * 100;
    
28.   if Censored then Censored = Censored * 100;
    
29. run;
    
30. 
    
31. data _null_;
    
32.   set HomTests;
    
33.   if Test = '对数秩' then call symputx('LogRank',cats(put(ProbChiSq,PValue6.3)));
    
34. run;
    
35. 
    
36. /*2.2 计算HR*/
    
37. ods output HazardRatios = HazardRatios;
    
38. proc phreg data = RawData;
    
39.   model time*_COL1(1) = group;
    
40.   hazardratio group;
    
41. run;
    
42. ods output close;
    
43. 
    
44. data _null_;
    
45.   set HazardRatios;
    
46.   call symputx('HR',cats(put(HazardRatio,8.2)));
    
47.   call symputx('CI',catx('-',put(WaldLower,8.2),put(WaldUpper,8.2)));
    
48. run;
    
49. 
    
50. /*3. 对计算结果进行绘图呈现*/
    
51. proc sgplot data = SurvivalPlot noborder nowall;
    
52.   styleattrs datacontrastcolors = (cx5b86c3 cxc7533b);
    
53.   step x = time y = Survival /group = StratumNum name = 'step';
    
54.   scatter x = time y = Censored /markerattrs = (symbol = plus);
    
55.   xaxistable NumRiskAndCensor /class = StratumNum colorgroup = StratumNum
    
56.     title = "Number at risk(Number censored)" titleattrs = (weight = bold);
    
57. 
    
58.   xaxis values=(0 to 135 by 15) display = (nolabel);
    
59.   yaxis values=(0 to 100 by 10);
    
60. 
    
61.   keylegend 'step' 'stat' /position = topright location = inside across = 1 noborder;
    
62.   inset "(*ESC*){unicode '0020'x}" "(*ESC*){unicode '0020'x}"
    
63.     "HR &HR. (95% CI &CI.), log rank P=&LogRank" /position = topright;
    
64.   label Survival = 'Disease free survival(%)';
    
65. run;

复制代码

Kaplan Meier plot

proc lifetest data=Graph20210812 plots=survival(strata=individual);
time time * censor(0);
strata Group;
run;

完成！

1. /*0. 环境设置*/
   
2. option missing = ' ';
   
3. filename RawData "D:\作图20210812.xlsx ";
   
4. 
   
5. /*1. 导入数据*/
   
6. proc import datafile = RawData dbms = excel out = RawData replace;
   
7. run;
   
8. 
   
9. /*2. 执行生存分析*/
   
10. /*2.1 计算各时点的生存率,并执行log-rank检验*/
    
11. ods output SurvivalPlot = SurvivalPlot HomTests = HomTests;
    
12. proc lifetest data = RawData plots=(survival(atrisk=0 to 135 by 15));
    
13.   strata group;
    
14.   time time*_COL1(1);
    
15. run;
    
16. ods output close;
    
17. 
    
18. data SurvivalPlot;
    
19.   set SurvivalPlot;
    
20. 
    
21.   if mod(time,15) then AtRisk = .;
    
22.   if Survival then Survival = Survival * 100;
    
23.   if Censored then Censored = Censored * 100;
    
24. run;
    
25. 
    
26. data _null_;
    
27.   set HomTests;
    
28.   if Test = '对数秩' then call symputx('LogRank',cats(put(ProbChiSq,PValue6.3)));
    
29. run;
    
30. 
    
31. /*2.2 计算HR*/
    
32. ods output HazardRatios = HazardRatios;
    
33. proc phreg data = RawData;
    
34.   model time*_COL1(1) = group;
    
35.   hazardratio group;
    
36. run;
    
37. ods output close;
    
38. 
    
39. data _null_;
    
40.   set HazardRatios;
    
41.   call symputx('HR',cats(put(HazardRatio,8.2)));
    
42.   call symputx('CI',catx('-',put(WaldLower,8.2),put(WaldUpper,8.2)));
    
43. run;
    
44. 
    
45. /*3. 对计算结果进行绘图呈现*/
    
46. proc sgplot data = SurvivalPlot noborder nowall;
    
47.   step x = time y = Survival /group = StratumNum name = 'step';
    
48.   scatter x = time y = Censored /markerattrs = (symbol = plus);
    
49.   xaxistable AtRisk /class = StratumNum colorgroup = StratumNum
    
50.     title = 'Number at risk' titleattrs = (weight = bold);
    
51. 
    
52.   xaxis values=(0 to 135 by 15) display = (nolabel);
    
53.   yaxis values=(0 to 100 by 10);
    
54. 
    
55.   keylegend 'step' 'stat' /position = topright location = inside across = 1 noborder;
    
56.   inset "(*ESC*){unicode '0020'x}" "(*ESC*){unicode '0020'x}"
    
57.     "HR &HR. (95% CI &CI.), log rank P=&LogRank" /position = topright;
    
58.   label Survival = 'Disease free survival(%)';
    
59. run;

复制代码

whymath 发表于 2021-8-13 23:53
完成！

谢谢大侠！ 还差最后一个关键环节。number censored如何绘制，恳请出手相助！
其他环节通过r与stata也基本能实现，就卡在number censored的绘制了。多谢多谢！

whymath 发表于 2021-8-13 23:53
完成！

再请教一个问题：
ods output SurvivalPlot = SurvivalPlot HomTests = HomTests;
当中 SurvivalPlot = SurvivalPlot无法实现，proc lifetest结果里没有看到SurvivalPlot。需要在什么方面注意呢

ods output SurvivalPlot无法输出数据集

edragon1983 发表于 2021-8-14 21:49
再请教一个问题：
ods output SurvivalPlot = SurvivalPlot HomTests = HomTests;
当中 SurvivalPlot = ...

这是本应有的结果，却未给出。应该观察程序日志，不要放过任何提示信息。

谢谢大侠出手相助！多谢多谢！