搜索
查看: 8684|回复: 11

TCGA表达批量做生存分析--菜鸟视野

[复制链接]

5

主题

13

帖子

284

积分

中级会员

Rank: 3Rank: 3

积分
284
发表于 2017-2-26 18:54:03 | 显示全部楼层 |阅读模式
本帖最后由 ybb7 于 2017-2-26 20:22 编辑

感谢jimmy的原贴,http://www.biotrainee.com/thread-929-1-1.html
这个是修改版,从菜鸟的视角来解决问题。
##上述数据第一行,第一列列都是如下形式


##生存时间,生成行列。得到的数据类型,X为基因,举例子

[AppleScript] 纯文本查看 复制代码
rm(list=ls())
setwd('C:\\Users\\Administrator\\Desktop\\')
###读取自己的TCGA数据
dat=read.csv('3expforRsuoxiao.csv',header =TRUE,row.names= 1,sep = ",")##基因不能太多
rt=read.csv('3osforR.csv',header =TRUE,row.names= 1,sep = ",")##,header =TRUE第一行,row.names= 1第一列  
os=rt[,2]  ##(时间)
STATUS=rt[,1] ###(生存的event)
log_rank_p <- apply(dat,  1, function(values1){  ## for a matrix 1 indicates rows, 2 indicates columns, c(1, 2) indicates rows and columns.
  group=ifelse(values1>median(values1),'high','low')  ##分高低组  ifelse(test, yes, no)test, an object which can be coerced to logical mode.
  survival_dat <- data.frame(group=group,os=os,STATUS=STATUS,stringsAsFactors = F)  ##行,列的生存时间, 即 基因表达值分成了高低组,和相应
  library(survival)###表格可以任意添加参数,如性别,年龄等。
  my.surv <- Surv(survival_dat$os,survival_dat$STATUS)  ###my.surv <- Surv(OS_MONTHS,OS_STATUS=='LIVING') 
  ##这个Surv函数第一个参数必须是数值型的时间,第二个参数是逻辑向量,1,0表示死亡与否
  kmfit2 <- survfit(my.surv~survival_dat$group)  ##### 生存数据对象~表达高低 创建KM生存曲线 $绝对引用?
  data.survdiff=survdiff(my.surv~group,data = survival_dat) ##data = survival_dat 单个的表达高低与生存数据对象
  ##用于不同组的统计检验 Two-sample test  survdiff(Surv(futime, fustat) ~ rx,data=ovarian)#生存数据对象~表达高低
  p.val = 1 - pchisq(data.survdiff$chisq, length(data.survdiff$n) - 1)
})
names(log_rank_p[log_rank_p<0.05])

经验分享:
1:模仿别人的代码,先把原代码 每行跑一遍,理解每行代码的意思,和注意一些小细节。
2:help()  不会的函数,百度,google,看每行的参数
3:学会引申别人的代码,被自己所用
简洁版
[AppleScript] 纯文本查看 复制代码
rm(list=ls())
setwd('C:\\Users\\Administrator\\Desktop')
dat=read.csv('3expforRsuoxiao.csv',header =TRUE,row.names= 1,sep = ",")
rt=read.csv('3osforR.csv',header =TRUE,row.names= 1,sep = ",")
os=rt[,2]   
STATUS=rt[,1]
log_rank_p <- apply(dat,  1, function(values1){ 
  group=ifelse(values1>median(values1),'high','low')
  survival_dat <- data.frame(group=group,os=os,STATUS=STATUS,stringsAsFactors = F) 
  library(survival)
  my.surv <- Surv(survival_dat$os,survival_dat$STATUS) 
  kmfit2 <- survfit(my.surv~survival_dat$group)
  data.survdiff=survdiff(my.surv~group,data = survival_dat)
  p.val = 1 - pchisq(data.survdiff$chisq, length(data.survdiff$n) - 1)
})
names(log_rank_p[log_rank_p<0.05])







本帖子中包含更多资源

您需要 登录 才可以下载或查看,没有帐号?立即注册

x
回复

使用道具 举报

4

主题

24

帖子

184

积分

注册会员

Rank: 2

积分
184
发表于 2017-3-29 04:07:02 | 显示全部楼层
pdf(paste(gene_name,".pdf",sep=""))
plot(kimfit2)
dev.off()
回复 支持 0 反对 1

使用道具 举报

6

主题

36

帖子

457

积分

中级会员

Rank: 3Rank: 3

积分
457
发表于 2017-2-26 21:28:34 | 显示全部楼层
你这个好像有个严重的问题,临床样本matrix和基因表达的matrix TCGA id没有对应,除非整理表格前已经整理好了
回复 支持 反对

使用道具 举报

5

主题

13

帖子

284

积分

中级会员

Rank: 3Rank: 3

积分
284
 楼主| 发表于 2017-2-26 21:42:45 | 显示全部楼层
017中大普外 发表于 2017-2-26 21:28
你这个好像有个严重的问题,临床样本matrix和基因表达的matrix TCGA id没有对应,除非整理表格前已经整理好 ...

我的是自己整理好的。
回复 支持 反对

使用道具 举报

0

主题

1

帖子

37

积分

新手上路

Rank: 1

积分
37
发表于 2017-3-4 18:43:39 | 显示全部楼层
ybb7 发表于 2017-2-26 21:42
我的是自己整理好的。

好像下载下来的数据是,临床资料一个表,基因信息一个表,是把它们都整理到同一个表格上吗?如何整理呀?
回复 支持 反对

使用道具 举报

0

主题

2

帖子

31

积分

新手上路

Rank: 1

积分
31
发表于 2017-3-22 15:25:09 | 显示全部楼层
楼主 我有一个问题  我按照你这个处理出来了我自己的数据 但是我不会批量输出km曲线图  我的导师要求我把200多个miRNA的生存曲线图都做出来。。。能不能帮帮我  我自己也正在摸索方法。。。但是感觉实在是太复杂了  
回复 支持 反对

使用道具 举报

0

主题

2

帖子

31

积分

新手上路

Rank: 1

积分
31
发表于 2017-3-23 10:34:10 | 显示全部楼层
我按照楼主的方法做出来了  但是我plot(kimfit2)只能看到一张km图 前面的都被覆盖了
请问有人知道该怎么输出所有的km图么?
回复 支持 反对

使用道具 举报

2

主题

17

帖子

145

积分

注册会员

Rank: 2

积分
145
发表于 2017-5-8 20:51:20 | 显示全部楼层
请问“##生存时间,生成行列。得到的数据类型,X为基因,举例子”  X为基因是什么意思?是指X那一列都是基因吗?这里的生存事件和生存时间排列的顺序是按sample1,2,3.。。。。的顺序吗?
回复 支持 反对

使用道具 举报

5

主题

13

帖子

284

积分

中级会员

Rank: 3Rank: 3

积分
284
 楼主| 发表于 2017-5-8 22:31:14 | 显示全部楼层
miniaowai 发表于 2017-5-8 20:51
请问“##生存时间,生成行列。得到的数据类型,X为基因,举例子”  X为基因是什么意思?是指X那一列都是基 ...

1:x基因是图1的gene_1中的某一个
2:Sample是图1中的样本, 这里的生存事件和生存时间排列的顺序是按sample1,2,3.。。。。的顺序吗
是的
回复 支持 反对

使用道具 举报

2

主题

17

帖子

145

积分

注册会员

Rank: 2

积分
145
发表于 2017-5-9 21:08:47 | 显示全部楼层
再请教楼主一个问题,小弟在阅读文献中看到
"   The association between the lncRNA expression
and patient’s disease free survival (DFS) and disease
specific survival (DSS) was assessed by univariable Cox
regression analysis along with a permutation test using
Biometric Research Branch-Array (BRB-Array) Tools"   请问这个通过BRB-Array Tools进行univariable Cox
regression analysis along with a permutation test  得到的结果应该与我们用上述代码在R语言中实现的一样吧,  所用表达矩阵、目的都是一样的,都是为了筛选与肿瘤预后相关的 lncRNA  这几天一直困惑着我 还请楼主帮忙解答下 感谢!
回复 支持 反对

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

QQ|手机版|小黑屋|生信技能树    

GMT+8, 2018-9-23 01:45 , Processed in 0.116431 second(s), 25 queries .

Powered by Discuz! X3.2

© 2001-2013 Comsenz Inc.