生存分析凭什么不需要矫正P值

生存分析是大数据时代筛选目标基因的超级有效策略。各种数据挖掘文章本质上都是要把目标基因集缩小，比如表达量矩阵通常是2万多个蛋白编码基因，不管是表达芯片还是RNA-seq测序的，采用何种程度的差异分析，最后都还有成百上千个目标基因。如果是临床队列，通常是会跟生存分析进行交集，或者多个数据集差异结果的交集，比如：多个数据集整合神器-RobustRankAggreg包 ，这样的基因集就是100个以内的数量了，但是仍然有缩小的空间，比如lasso等统计学算法，最后搞成10个左右的基因组成signature即可顺利发表。

虽然生存分析如此重要而且如此常见，但是仍然有一些未解之谜，不同数据库来源，病人的不同时期的记录信息，以及不同的阈值分组，拿到的结果居然是可以不一样的！虽然大家都倾向于做各种花式分析，然后挑选具有统计学显著意义的生存分析结果。

生存分析最重要的是病人分组

我在生信技能树多次分享过生存分析的细节；

• 基因表达量高低分组的cox和连续变量cox回归计算的HR值差异太大?
• 学徒作业-两个基因突变联合看生存效应
• TCGA数据库里面你的基因生存分析不显著那就TMA吧
• 对“不同数据来源的生存分析比较”的补充说明
• 批量cox生存分析结果也可以火山图可视化
• 既然可以看感兴趣基因的生存情况，当然就可以批量做完全部基因的生存分析
• 多测试几个数据集生存效应应该是可以找到统计学显著的！
• 我不相信kmplot这个网页工具的结果(生存分析免费做)
• 为什么不用TCGA数据库来看感兴趣基因的生存情况
• 200块的代码我的学徒免费送给你，GSVA和生存分析
• 集思广益-生存分析可以随心所欲根据表达量分组吗
• 生存分析时间点问题
• 寻找生存分析的最佳基因表达分组阈值
• apply家族函数和for循环还是有区别的（批量生存分析出图bug）
• TCGA数据库生存分析的网页工具哪家强
• KM生存曲线经logRNA检验后也可以计算HR值

可以看到，有基因表达量高低分组，基因突变与否分组，多个基因表达量和突变联合分组，甲基化高低分组，gsea和gsva等基因集得分进行分组，五花八门，其中200块的代码我的学徒免费送给你，GSVA和生存分析 视频值得看！

TCGA数据库的RNA-seq表达矩阵全部基因高低表达分组后批量生存分析

虽然说可以各种花式分析，然后挑选具有统计学显著意义的生存分析结果，但是最开始基本上都是对全部基因分组后批量生存分析，可以是表达量高低，包括mRNA，lncRNA，miRNA的表达量，以及甲基化信号值高低等等，一个基因可以把病人分组，只要是有分组，就可以进行一次生存分析。

比如我们可以下载TCGA数据库的RNA-seq表达矩阵，读入到R里面构建成为 expr 这个数据变量，然后整理好临床表型，构建成为phe这个变量，接下来就可以使用下面的代码对RNA-seq表达矩阵全部基因高低表达分组后批量生存分析：

## 批量生存分析 使用  logrank test 方法
mySurv=with(phe,Surv(time, event))
log_rank_p <- apply(expr , 1 , function(gene){
  # gene=as.numeric(expr[1,])
  phe$group=ifelse(gene>median(gene),'high','low')  
  data.survdiff=survdiff(mySurv~group,data=phe)
  p.val = 1 - pchisq(data.survdiff$chisq, length(data.survdiff$n) - 1)
  #cat(p.val) 
  return(p.val)
})
log_rank_p=sort(log_rank_p)
head(log_rank_p)
boxplot(log_rank_p)  

但是最近挑选具有统计学显著意义的生存分析结果的基因时候，发现很多基因都是表达量相关的，也就是说，它们尽管说是不同的基因在不同病人表达量不同，但是它们对病人的分组效果其实是类似的。

那么，我们一下子对几万个基因进行批量生存分析，每一次每一个基因的生存分析都是独立的P值，为什么我们没有对这样的P值进行矫正呢？

大家耳熟能详的矫正P值有，adjust.p ， q值，以及FDR，他们的作用都是把P值的放大，这样之前那些小于0.05或者0.01的具有统计学显著的基因就不再显著啦，就是把筛选标准严格一点而已。

生存分析凭什么不需要矫正P值？

难道就是因为我们希望统计学显著的生存结果，就选择性展示它吗？