花花花写于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
目前留有一个疑问，端详了一下模型内部结构，还不知道如何找到构建模型用到的基因。由于公众号文章一经发出无法修改，此问题后期如需更新，我将在简书（原文链接跳转）更新。
微信公众号生信星球同步更新我的文章，欢迎大家扫码关注！
我们有为生信初学者准备的学习小组,点击查看◀️
想要参加我的线上线下课程，也可加好友咨询🔼
如果需要提问，请先看生信星球答疑公告
网友评论