肿瘤单基因研究的若干问题之三：目的基因的筛选1.0

肿瘤单基因研究的若干问题之三：目的基因的筛选1.0 by 单细胞天地

前言

Hello小伙伴们大家好！临床菜鸟刚刚上路，本着菜就多练的原则，几乎住在了科室。

医疗同行应该都知道了三明医改五年内全国推行的消息，不知是喜还是忧？

无论怎样，我们还是要做好自己，要么努力提升职业技能，好好学临床+科研；要么提前规划改行；要么准备雅思托福考试。

总之，我们能做的就是无条件相信国家，然后做好自己。

好了，不制造焦虑了。

言归正传，开始我们今天的话题：如何筛选目的基因？

如果你身处大佬课题组，大概率在硕博入学前课题就固定了（有目的基因）那么恭喜你，可以暂时划过本期推文了（后期将单出一期有目的基因的推文，大佬课题组的也可以忽略）。

如果是一般的课题组，或是放养课题组，抑或是没有课题组，那么本期推文或许适合你和我一起来探讨。目的基因的筛选，从测序数据入手，或许你可以这样做。。。

1.从基因表达与预后入手

如果没有目的基因，最简单粗暴的方式就是通过差异基因分析和生存分析筛选目的基因。表达与预后的问题前两期已经和大家聊过。虽然，但是，，，还是推荐大家首选在肿瘤中高表达且是预后危险因素的基因进行研究。

在进行差异基因分析时，有几个常用的R包，它们分别是DESeq2、edgeR和limma。以下是关于这些R包的一些信息：

DESeq2：这是一个基于负二项式分布的差异表达分析工具，适用于RNA-seq数据。它通过估计样本之间的大小差异来标准化计数数据，并使用局部回归来估计分散度。DESeq2可以处理批次效应和样本间的差异性。安装DESeq2的代码如下：
```
if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")
BiocManager::install("DESeq2")
```
edgeR：edgeR也是一个用于RNA-seq数据差异表达分析的R包，它使用经验贝叶斯方法来估计差异表达。edgeR适合于比较小的样本集，并且可以处理技术重复。安装edgeR的代码如下：
```
BiocManager::install("edgeR")
```
limma：limma包提供了用于微阵列和RNA-seq数据的线性模型和统计测试。它使用Bayesian方法来估计差异表达，并且可以用于处理大规模的基因表达数据。limma还支持voom变换，这是一种用于改善计数数据的统计方法。安装limma的代码如下：
```
BiocManager::install("limma")
```

这些R包都有各自的优势和特点，选择哪个包取决于数据类型、样本大小和研究目的。在使用这些工具之前，需要准备好表达矩阵和分组信息，然后按照各自的流程进行分析。分析结果通常包括差异表达基因列表、统计显著性和差异倍数等。

在使用这些R包时，可以参考相应的用户指南和文献，以确保分析的正确性和结果的可解释性。此外，还可以使用其他工具和R包，如pheatmap、ggplot2等，来进行结果的可视化和进一步分析。

具体的分析代码，微信搜一搜一大堆，在这里就不罗列了。

那么生存分析（KM/COX）的批量筛选，可以在差异表达基因的基础上进行，下面的代码还是很好用的：

library(limma)
library(survival)

pFilter=0.05                         
expFile="diffGeneExp.txt"            
cliFile="time.txt"     
setwd("")    
rt=read.table(expFile, header=T, sep="\t", check.names=F)
rt=as.matrix(rt)
rownames(rt)=rt[,1]
exp=rt[,2:ncol(rt)]
dimnames=list(rownames(exp), colnames(exp))
data=matrix(as.numeric(as.matrix(exp)), nrow=nrow(exp), dimnames=dimnames)
data=avereps(data)
data=data[rowMeans(data)>0,]
data=t(data)
cli=read.table(cliFile, header=T, sep="\t", check.names=F, row.names=1)
cli$futime=cli$futime/365
sameSample=intersect(row.names(data), row.names(cli))
data=data[sameSample,,drop=F]
cli=cli[sameSample,,drop=F]
rt=cbind(cli, data)

outTab=data.frame()
for(gene in colnames(rt[,3:ncol(rt)])){
 if(sd(rt[,gene])<0.1){next}
 group=ifelse(rt[,gene]<=median(rt[,gene]), "low", "high")
 
    #km
 diff=survdiff(Surv(futime, fustat) ~group, data=rt)
 pValue=1-pchisq(diff$chisq,df=1)
 
 #cox
 cox=coxph(Surv(futime, fustat) ~ rt[,gene], data = rt)
 coxSummary=summary(cox)
 coxP=coxSummary$coefficients[,"Pr(>|z|)"]

 if((pValue<pFilter) & (coxP<pFilter)){
     outTab=rbind(outTab,
                  cbind(gene=gene,
                        KM=pValue,
                        HR=coxSummary$conf.int[,"exp(coef)"],
                        HR.95L=coxSummary$conf.int[,"lower .95"],
                        HR.95H=coxSummary$conf.int[,"upper .95"],
                  coxPvalue=coxP) )
 }
}

outTab=outTab[order(as.numeric(as.vector(outTab[,"KM"]))),]
write.table(outTab, file="survival.txt", sep="\t", row.names=F, quote=F)

由于上述方法简单粗暴，我们在具体执行时，可能会面临以下问题：一是筛到的基因太多；二是筛到的基因太少。对于这两个问题，应该如何解决，最终顺利筛选到你心仪的基因？后期我们会专门出一期推文。同时也欢迎大家给我们私信或留言讨论。

2. 从基因集入手

如果课题组有固定的研究方向，或者你有感兴趣的研究方向，或者想追热门，那么可以从某一大方向入手。如表观遗传？代谢？程序性细胞死亡？免疫？等等等等

逐渐缩小范围：以程序性细胞死亡（https://mp.weixin.qq.com/s/SlESiT9PSbD_eHKI7AngEQ）为例，其下包括：细胞凋亡（Apoptosis）、坏死性凋亡 (Necroptosis)、细胞焦亡（Pyroptosis）、铁死亡（Ferroptosis）、PARP1依赖的死亡（Parthanatos）、内亡或者侵入性细胞死亡（Entotic cell death）、依赖中性粒细胞外诱捕网的死亡（Netotic cell death）、溶酶体依赖性细胞死亡（LCD），、自噬依赖性细胞死亡（ADCD）、免疫原性细胞死亡（immunogenic cell death，ICD）、铜死亡（Cuprotosis）、双硫死亡（Disulfidptosis）......

选定了方向之后，就要筛选该方向的相关基因，那么就要先获取下载基因集。基因集的下载方式有很多，主要有以下几种:

1.专用数据库搜索：比如你对铁死亡感兴趣，可以去铁死亡数据库下载基因集，FerrDb数据库（http://www.zhounan.org/ferrdb/current/）于2020年1月份发布，目前更新至2.0版本，数据库发布时整合了Pubmed中784篇铁死亡相关研究文献数据，提供253个铁死亡调节因子、111个铁死亡标志基因和95个铁死亡相关疾病信息，其中铁死亡调节因子包含108个铁死亡驱动基因和69个铁死亡抑制基因，以及35个铁死亡诱导剂和41个铁死亡抑制剂。如果对线粒体感兴趣，可以去专门收录线粒体相关基因标签的数据库——MitoCarta3.0（https://www.broadinstitute.org/mitocarta/mitocarta30-inventory-mammalian-mitochondrial-proteins-and-pathways）。

2.如果没有特定的基因集数据库，可以使用多个通用数据库和工具，以下是一些常用的资源：

MSigDB (Molecular Signatures Database): 提供了大量注释的基因集，用于基因集富集分析（GSEA）。这些基因集可以根据不同的生物学主题进行分类，如Hallmark基因集、KEGG通路、GO术语等。MSigDB 官网
GeneCards: 一个全面的基因信息数据库，提供基因的基本信息、功能、表达数据以及与疾病的关系。GeneCards 官网
KEGG (Kyoto Encyclopedia of Genes and Genomes): 一个集成的数据库资源，用于理解生物系统的高级功能和用途，特别是大规模分子数据集。KEGG 官网
Reactome: 提供人类生物过程的详细路径信息，包括代谢通路、信号传导等。Reactome 官网
Gene Ontology (GO): 提供基因和基因产品的标准化描述，包括生物过程、分子功能和细胞组分。GO 官网
NCBI (National Center for Biotechnology Information): 提供多种基因相关信息，包括基因序列、基因表达数据等。NCBI GEO 官网
CNGBdb (China National GeneBank DataBase): 提供生物大数据共享和应用服务的平台，包括多个科学数据库。CNGBdb 官网
Bio-db Search 2.0: 提供多种生物数据库的搜索，包括蛋白质、基因、代谢通路等。Bio-db Search 官网
OmicShare Tools: 提供多种生物信息学分析工具，包括KEGG富集分析、GO富集分析等。OmicShare Tools 官网