TCGA-12.支持向量机生存模型

花花花写于2020-02-07 闭门不出的第n天，我们将原定于下周一的学习小组提前到今天开始了，我也收到公司通知，我们可能五一之前都没办法出去做线下培训，要转型线上培训了，如果不是这场灾难，此刻我应该是在广州的酒店，准备明天的课程了。

1.准备输入数据

输入数据是TCGA的表达矩阵expr、临床信息meta和group_list。保存为forest.Rdata了，在生信星球公众号后台聊天窗口回复“森林”即可获得。

load("forest.Rdata")
exprSet = expr[,group_list=="tumor"]

## 必须保证生存资料和表达矩阵，两者一致
all(substring(colnames(exprSet),1,12)==meta$ID)
#> [1] TRUE

library(ROCR)
library(genefilter)
library(Hmisc)
library(e1071)

2.构建支持向量机模型

2.1.切割数据

用R包caret切割数据，生成的结果是一组代表列数的数字，用这些数字来给表达矩阵和meta取子集即可。

#library(caret)
#set.seed(12345679)
#sam<- createDataPartition(meta$event, p = .5,list = FALSE)
load("sam.Rdata")

train <- exprSet[,sam]
test <- exprSet[,-sam]
train_meta <- meta[sam,]
test_meta <- meta[-sam,]

2.2 train数据集建模

x=t(log2(train+1))
y=as.factor(train_meta$event)
model = svm(x,y,kernel = "linear")
summary(model) 
#> 
#> Call:
#> svm.default(x = x, y = y, kernel = "linear")
#> 
#> 
#> Parameters:
#>    SVM-Type:  C-classification 
#>  SVM-Kernel:  linear 
#>        cost:  1 
#> 
#> Number of Support Vectors:  179
#> 
#>  ( 108 71 )
#> 
#> 
#> Number of Classes:  2 
#> 
#> Levels: 
#>  0 1

2.3.模型预测

用训练集构建模型，预测测试集的生死。不同于其他模型，这个预测结果是分类变量，直接预测生死(0,1)，而不是prob(百分数)。

x=t(log2(test+1))
y=as.factor(test_meta$event)
pred = predict(model, x)
table(pred,y)
#>     y
#> pred   0   1
#>    0 151  42
#>    1  36  32

目前留有一个疑问，端详了一下模型内部结构，还不知道如何找到构建模型用到的基因。由于公众号文章一经发出无法修改，此问题后期如需更新，我将在简书（原文链接跳转）更新。

微信公众号生信星球同步更新我的文章，欢迎大家扫码关注！

我们有为生信初学者准备的学习小组,点击查看◀️
想要参加我的线上线下课程，也可加好友咨询🔼
如果需要提问，请先看生信星球答疑公告