0x00 前情提要

同学问了一个关于使用 ezcox 进行批量Cox 模型处理的问题,用的是R语言(什么鬼?)。火线在菜鸟教程上扫了一遍概念和基本语法,又结合一个案例练习,算是稍稍了解一下这个语言,记录一下解决的思路和过程,以备复盘和参考。

发来了代码和报错两部分信息,分两步同时进行:

  1. 直接通过报错信息,试试能不能google出相关问题和解决方式;
    1. 报的错误是一个通用性的错误,不太好定位解决具体方案,但是起到了理解R语言部分语法的作用
  2. 找了一个案例,直接从搭环境到导出csv数据,实操了一遍;

    大佬如是说:cox函数,判断一个变量和人的生存状况(用生存时间time和生存状态status 0指生存1指死亡)的相关性,导出一个HR值(距离1越远就越影响越大)和P值

0x01 解决问题流程

使用 ezcox 进行批量Cox 模型处理的问题 - 图1

0x02 R语言初见

  • R 语言是为数学研究工作者设计的一种数学编程语言,主要用于统计分析、绘图、数据挖掘。
  • R 语言与 C 语言都是贝尔实验室的研究成果,但两者有不同的侧重领域,R 语言是一种解释型的面向数学理论研究工作者的语言,而 C 语言是为计算机软件工程师设计的。R 语言是解释运行的语言(与 C 语言的编译运行不同),它的执行速度比 C 语言慢得多,不利于优化。但它在语法层面提供了更加丰富的数据结构操作并且能够十分方便地输出文字和图形信息,所以它广泛应用于数学尤其是统计学领域。
  • R 语言官方网站:https://cran.r-project.org/
  • 官方镜像站列表:https://cran.r-project.org/mirrors.html

    0x03 部署环境

  • 安装R语言官方的IDE - cran or 安装更加美观的 R-Studio😳

  • 安装包 ( 详见例题 )
  • 如果使用xlsx函数,需要安装java 环境(java 9+jdk)

    0x04 基本语法

    略:基本上是用到了再去查

  • 菜鸟教程-R语言教程

  • Google

    0x05 例题实操

    1. 安装ezcox包并调用

    ```r

    library(ezcox)

Welcome to ‘ezcox’ package!

You are using ezcox version 1.0.2

Project home : https://github.com/ShixiangWang/ezcox Documentation: https://shixiangwang.github.io/ezcox

Cite as : arXiv:2110.14232

  1. ```r
  2. > library(survival) #安装survival数据包

2. 查看所用数据集

  1. > lung #查看lung数据集
  2. --------------------------------------------------------------------------
  3. inst time status age sex ph.ecog ph.karno pat.karno meal.cal wt.loss
  4. 1 3 306 2 74 1 1 90 100 1175 NA
  5. 2 3 455 2 68 1 0 90 90 1225 15
  6. 3 3 1010 1 56 1 0 90 90 NA 15
  7. 4 5 210 2 57 1 1 90 60 1150 11
  8. 5 1 883 2 60 1 0 100 90 NA 0
  9. 6 12 1022 1 74 1 1 50 80 513 0
  10. 7 7 310 2 68 2 2 70 60 384 10
  11. 8 11 361 2 71 2 2 60 80 538 1
  12. 9 1 218 2 53 1 1 70 80 825 16
  13. 10 7 166 2 61 1 2 70 70 271 34
  1. > str(lung)
  2. ----------------------------------------------------------------
  3. 'data.frame': 228 obs. of 10 variables:
  4. $ inst : num 3 3 3 5 1 12 7 11 1 7 ...
  5. $ time : num 306 455 1010 210 883 ...
  6. $ status : num 2 2 1 2 2 1 2 2 2 2 ...
  7. $ age : num 74 68 56 57 60 74 68 71 53 61 ...
  8. $ sex : num 1 1 1 1 1 1 2 2 1 1 ...
  9. $ ph.ecog : num 1 0 0 1 0 1 2 2 1 2 ...
  10. $ ph.karno : num 90 90 90 90 100 50 70 60 70 70 ...
  11. $ pat.karno: num 100 90 90 60 90 80 60 80 80 70 ...
  12. $ meal.cal : num 1175 1225 NA 1150 NA ...
  13. $ wt.loss : num NA 15 15 11 0 0 10 1 16 34 ...

3. 分类变量因子化

  1. lung$sex <- factor(lung$sex)
  2. lung$ph.ecog <- factor(lung$ph.ecog)
  1. str(lung)
  2. ----------------------------------------------------
  3. 'data.frame': 228 obs. of 10 variables:
  4. $ inst : num 3 3 3 5 1 12 7 11 1 7 ...
  5. $ time : num 306 455 1010 210 883 ...
  6. $ status : num 2 2 1 2 2 1 2 2 2 2 ...
  7. $ age : num 74 68 56 57 60 74 68 71 53 61 ...
  8. $ sex : Factor w/ 2 levels "1","2": 1 1 1 1 1 1 2 2 1 1 ...
  9. $ ph.ecog : Factor w/ 4 levels "0","1","2","3": 2 1 1 2 1 2 3 3 2 3 ...
  10. $ ph.karno : num 90 90 90 90 100 50 70 60 70 70 ...
  11. $ pat.karno: num 100 90 90 60 90 80 60 80 80 70 ...
  12. $ meal.cal : num 1175 1225 NA 1150 NA ...
  13. $ wt.loss : num NA 15 15 11 0 0 10 1 16 34 ...

4. 批量完成单因素Cox回归分析

  1. > results <- ezcox(lung, time = "time",status = "status",covariates = c("age", "sex", "ph.ecog","ph.karno","pat.karno","meal.cal","wt.loss"))
  2. -----------------------------------------------------------------------------------------------------------------------------------------------
  3. => Processing variable age
  4. ==> Building Surv object...
  5. ==> Building Cox model...
  6. ==> Done.
  7. => Processing variable sex
  8. ==> Building Surv object...
  9. ==> Building Cox model...
  10. ==> Done.
  11. => Processing variable ph.ecog
  12. ==> Building Surv object...
  13. ==> Building Cox model...
  14. ==> Done.
  15. => Processing variable ph.karno
  16. ==> Building Surv object...
  17. ==> Building Cox model...
  18. ==> Done.
  19. => Processing variable pat.karno
  20. ==> Building Surv object...
  21. ==> Building Cox model...
  22. ==> Done.
  23. => Processing variable meal.cal
  24. ==> Building Surv object...
  25. ==> Building Cox model...
  26. ==> Done.
  27. => Processing variable wt.loss
  28. ==> Building Surv object...
  29. ==> Building Cox model...
  30. ==> Done.

5. 查看结果

  1. > results
  2. -------------------------------------------------------------------------------------------------------------
  3. # A tibble: 9 x 12
  4. Variable is_control contrast_level ref_level n_contrast n_ref beta HR lower_95 upper_95 p.value
  5. <chr> <lgl> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
  6. 1 age FALSE age age 228 228 0.0187 1.02 1 1.04 4.19e-2
  7. 2 sex FALSE 2 1 90 138 -0.531 0.588 0.424 0.816 1.49e-3
  8. 3 ph.ecog FALSE 1 0 113 63 0.369 1.45 0.98 2.13 6.34e-2
  9. 4 ph.ecog FALSE 2 0 50 63 0.916 2.5 1.61 3.88 4.48e-5
  10. 5 ph.ecog FALSE 3 0 1 63 2.21 9.1 1.22 67.9 3.14e-2
  11. 6 ph.karno FALSE ph.karno ph.karno 227 227 -0.0164 0.984 0.972 0.995 4.96e-3
  12. 7 pat.karno FALSE pat.karno pat.karno 225 225 -0.0199 0.98 0.97 0.991 2.82e-4
  13. 8 meal.cal FALSE meal.cal meal.cal 181 181 -0.000124 1 0.999 1 5.93e-1
  14. 9 wt.loss FALSE wt.loss wt.loss 214 214 0.00132 1 0.989 1.01 8.28e-1
  15. # ... with 1 more variable: global.pval <dbl>

6. 表格形式查看(富文本格式)

  1. > library(knitr)
  2. > knitr::kable(results)
  3. ----------------------------------------------------------------------------------------------------------------------------------------
  4. |Variable |is_control |contrast_level |ref_level | n_contrast| n_ref| beta| HR| lower_95| upper_95| p.value| global.pval|
  5. |:---------|:----------|:--------------|:---------|----------:|-----:|---------:|-----:|--------:|--------:|--------:|-----------:|
  6. |age |FALSE |age |age | 228| 228| 0.018700| 1.020| 1.000| 1.040| 4.19e-02| 0.039500|
  7. |sex |FALSE |2 |1 | 90| 138| -0.531000| 0.588| 0.424| 0.816| 1.49e-03| 0.001110|
  8. |ph.ecog |FALSE |1 |0 | 113| 63| 0.369000| 1.450| 0.980| 2.130| 6.34e-02| 0.000356|
  9. |ph.ecog |FALSE |2 |0 | 50| 63| 0.916000| 2.500| 1.610| 3.880| 4.48e-05| 0.000356|
  10. |ph.ecog |FALSE |3 |0 | 1| 63| 2.210000| 9.100| 1.220| 67.900| 3.14e-02| 0.000356|
  11. |ph.karno |FALSE |ph.karno |ph.karno | 227| 227| -0.016400| 0.984| 0.972| 0.995| 4.96e-03| 0.005970|
  12. |pat.karno |FALSE |pat.karno |pat.karno | 225| 225| -0.019900| 0.980| 0.970| 0.991| 2.82e-04| 0.000412|
  13. |meal.cal |FALSE |meal.cal |meal.cal | 181| 181| -0.000124| 1.000| 0.999| 1.000| 5.93e-01| 0.590000|
  14. |wt.loss |FALSE |wt.loss |wt.loss | 214| 214| 0.001320| 1.000| 0.989| 1.010| 8.28e-01| 0.829000|

7. 结果输出为CSV格式

  1. write.csv(results,"D:\\results.csv")
  2. ⬇️⬇️结果⬇️⬇️ D盘导出一个 results.csv 文件

image.png

0x06 解决问题

当我稍微理解R语言并开始着手调试同学的命令时,我其实才真正理解他的问题

  1. library(ezcox) #安装程序包 ezcox
  2. library(dplyr) #安装程序包 dplyr - 用做数据转换和筛选
  3. data<-read.xlsx("C:\\Users\\Robot\\Desktop\\survival_group1.xlsx","Sheet1", header = TRUE, sep = ",") #以表格的形式读取"survival_group1.xlsx"中的数据
  4. data<-data[,2:100] #将data数据中的 2到100列筛选保存到 data 中(去掉的time和status 是自变量)
  5. #变量因子化
  6. data[,2:99]<-lapply(data[,2:99],as.factor) #将所有data中的所有列的数据 变量因子化(个人认为这步有问题,因子化分为10级,
  7. #可是上百个随机数据10级是不够用的,建议去掉)
  8. results <- ezcox(data, time = "time",status = "status",covariates = c(3:99),return_models = TRUE )
  9. write.csv(results,"C:\\Users\\Robot\\Desktop\\results.csv")

0x07 参考链接