灰色预测模型代码_生信审稿人最常问的验证!临床预测模型中的PCA主成分分析!这点你注意到了没!(附代码)...
在高维度数据不断出现的现在,输入特征的数量呈指数形式迅速增长,形成了所谓的维数灾难。今天,我们将重点讲讲主成分分析,从原始数据集中找出一个更小的,但能最大程度保留原来大部分信息的变量集合。
所谓主成分分析,即Principle Component Analysis (PCA),其实就是寻找主成分的过程,且通常被认为是一种特殊的非监督学习算法,可以对复杂或多变的变量进行处理分析。一般而言,成分可以认为是特征的规范化线性组合。其中,第一主成分是能够最大程度解释数据中方差的特征线性组合,而第二主成分是在方向上与第一主成分垂直的条件下,最大程度解释数据中方差的另一种线性组合。其后的每一个主成分都遵循相同的原则。此外,主成分得分是对于每个主成分和每个观测变量计算而得到的。
下面,我们一起来看下如何进行主成分分析并建立模型。
1.R包的安装与读取
1.R包的安装与读取
rm(list = ls()) #清空环境变量
options(stringsAsFactors = F)###1. R包的安装与读取
if(!require(psych))install.packages("psych")
if(!require(ggplot2))install.packages("ggplot2")
if(!require(ggpubr))install.packages("ggpubr")
if(!require(ROCR))install.packages("ROCR")
library(psych) #PCA packages
library(ggplot2)
library(ggpubr)
library(ROCR)
psych包常被用于进行主成分分析,包括:判断主成分的个数;提取主成分;获取主成分得分;以及列出主成分方程,解释主成分意义。
2.数据的读取与预处理
2.1 读取表达数据
rt <- read.table("exp.txt", header=T, sep=" ", check.names=F, row.names = 1)
exp <- t(log2(rt+1))
str(exp)
与之前的模型构建一样,首先使用read.table()函数读取表达数据,并通过结构函数str()检查数据。
2.2 数据标准化
exp.scale <- scale(exp)
exp.cor <- cor(exp.scale)
cor.plot(exp.cor)
在分析之前,我们需要对数据进行标准化,使数据的均值为0,标准差为1。同时,完成标准化后,使用psych包提供的cor.plot()函数,将输入基因之间的相关性进行可视化展示。
2.3 读取临床数据
cli <- read.table("cliData.txt",header=T,sep=" ",check.names=F, row.names = 1)
head(cli)
str(cli)
随后,将提前清洗准备好的临床数据进行读取。
3.主成分分析
3.1 主成分抽取
pca <- principal(exp.scale, rotate = "none")
通过psych包的principal()函数抽取主成分。首先,我们将rotate的参数设置为none,得到其中的各个成分。
plot(pca$values, type = "b", ylab = "Eigenvalues", xlab = "Component")
接着,为了确定要保留的成分的数量,使用碎石图来评估能解释大部分数据方差的主成分,其中x轴表示主成分的数量,y轴表示相应的特征。在碎石图中,需要找出使变化率降低的那个点,也就是常说的“肘点”或者弯曲点。从图中可以看出,4个主成分是比较令人信服的。
3.2 正交旋转与解释
正交旋转又称为“方差最大法”,而旋转的意义是使变量在某个主成分上的载荷最大化,减少主成分之间的相关性,从而有助于对主成分的解释。
pca.rotate <- principal(exp.scale, nfactors = 4, rotate = "varimax")
pca.rotate
在输出的结果中,4个主成分的变量载荷分别标注为RC1到RC4,而SS loading中的值是每个主成分的特征值。同时,Cumulative Var结果表示这4个旋转后的主成分可以解释86%的全部方差。一般而言,我们选择的主成分至少应该解释大于70%的全部方差。可见,4个主成分已经可以用于后续的建模分析。
3.3 根据主成分建立因子得分
pca.scores <- data.frame(pca.rotate$scores)
head(pca.scores)
检查旋转和的主成分载荷,并将其作为每个患者的因子得分。而且,这些得分说明了每个患者与旋转后的主成分之间的相关程度。
exp.pca <-cbind(cli, pca.scores)
最后,将得分与临床信息进行合并。
4.模型构建与评价
4.1 回归分析
#通过lm()函数建立线性模型
lm_all <- lm(fustat ~ RC1+RC2+RC3+RC4, data = exp.pca)
summary(lm_all)
首先,将所有的因子作为输入,然后查看结果摘要。结果显示,整体模型在统计学上是具有显著性的,P值为0.0348。而且,其中RC1和RC4具有显著性。
#使用RC1和RC4变量建立线性模型
lm_2 <- lm(fustat ~ RC1+RC4, data = exp.pca)
summary(lm_2)
随后,将不显著的主成分因子排出模型,只保留RC1和RC4,再一次构建模型。可以发现,整体模型的统计学显著性得到了提升,P值为0.0057。
4.2 模型评价
exp.pca$predict <- lm_2$fitted.values
ggboxplot(exp.pca, x = "fustat", y = "predict", color = "fustat", palette = "jco", add = "jitter") + stat_compare_means()
提取模型得到的预测值,且箱型图结果显示两组之间具有显著性差异,P值为0.0035。
4.3 ROC分析
pred <- prediction(exp.pca$predict, exp.pca$fustat)
ROC <- performance(pred,"tpr","fpr")
auc <- performance(pred,"auc")@y.values[[1]]
auc
#[1] 0.7101449
plot(ROC, col="red", #曲线的颜色 xlab="False positive rate", ylab="True positive rate", #x轴和y轴的名称 lty=1,lwd=3, main=paste("AUC=", auc))
abline(0, 1, lty=2, lwd=3) #绘制对角线
同时,ROC曲线分析结果显示,AUC值为0.71,表明其模型的预测能力较为一般。
好了,PCA的相关内容就介绍到这里了,简单的讲呢,就是使用一些相关的变量构造一个新的变量。在今后的文章分析中,遇到其他模型不理想时,或者其他高维度数据时候,不妨试试主成分分析来进行分析~~~
欢迎大家关注解螺旋生信频道-挑圈联靠公号~
本文首发于“挑圈联考”微信公众号
转载请注明:解螺旋·临床医生科研成长平台
灰色预测模型代码_生信审稿人最常问的验证!临床预测模型中的PCA主成分分析!这点你注意到了没!(附代码)...相关推荐
- c语言自学门槛,初学C语言的人最常问的几个问题
初学C语言的人最常问的几个问题 C语言是一门通用计算机编程语言,应用广泛.对于新手来说学习C语言并不是那么容易,下面是C语言初学者最常问的几个问题,欢迎阅读! 1.多久能学会编程? 这是一个没有答案的 ...
- 引用另一模板的宏_生信人值得拥有的编程模板Shell
前言 "工欲善其事必先利其器",生信工程师每天写代码.搭流程,而且要使用至少三门编程语言,没有个好集成开发环境(IDE,Integrated Development Environ ...
- cpan安装_生信平台搭建(五):安装perl模块
perl模块也是生物信息分析中经常需要配置的东西,尽管很多人觉得python很流行,但是依然有大量的生物软件依赖于perl模块,如果配置不正确就无法运行,典型的就是circos,里面调用大量perl的 ...
- 编程心得体会_生信编程语言的经验之谈
对于生信狗来说,日常的数据分析任务大概可以分成三大类: 纯文本处理 数值计算 图表可视化 其中,对这几个任务的可选的解决方案,以及各自的最优解决方案分别为: 纯文本处理: Linux下的文本处理命令, ...
- .md是什么文件_生信中常见的数据文件格式
TCGA | GEO | 文献阅读 | 数据库 | 理论知识 R语言 | Bioconductor | 服务器与Linux 前面我们介绍了各种测序技术的原理:illumina.Sanger.第三代和第 ...
- 生物信息学软件_生信软件操作视频教程大赛
楔子 朋友圈偶然看到由信息中心生命科学图书馆联合营养与健康院团委.研究生会.中科院创新创业俱乐部举办的2019年度生物软件操作视频征集大赛,虽然是生物软件操作大赛,但是里面列出来的几乎都是生物信息学软 ...
- 属于服务器端运行的程序_生信分析云平台产品开发 - 5 生信分析pipeline服务器端运行...
在上文 [生信分析云平台产品开发 - 4 生信分析pipeline的图形化] 讨论了生信分析pipeline的图形化,如何用图形的方式显示生信pipeline,但是pipeline脚本按照变量的形式保 ...
- java肝癌晚期_生信分析43.肿瘤浸润免疫与肝癌(HCCDB+oncomine)
生信论文的套路 ONCOMINE从全景.亚型两个维度做表达差异分析: 临床标本从蛋白水平确认(或HPA数据库),很重要: Kaplan-Meier Plotter从临床意义的角度阐明其重要性: cBi ...
- cluego使用说明_生信分析绘图神器,你值得拥有!
GO和KEGG分析是最常用的生信分析方法,在SCI论文中也经常见到,那么你能想到的GO和KEGG分析结果的展示方法有哪些? 条形图? 条形图2? 饼状图? 表格? 相比于上面这些,这样的网络图展示起来 ...
最新文章
- 机器翻译:谷歌翻译是如何对几乎所有语言进行翻译的?
- 【正一专栏】内马尔要走快走、走好不送!
- UC伯克利博士尤洋回国创业,曾破ImageNet纪录!已获超千万融资
- Python实现 灰色关联分析 与结果可视化
- [译] 曝光!UX 行话大全
- 【收藏】gitee:使用k8s部署nacos
- Spark _12_每个网址的每个地区访问量 ,由大到小排序
- EasyDarwin开源手机直播方案:EasyPusher手机直播推送,EasyDarwin流媒体服务器,EasyPlayer手机播放器...
- 常州win8如何禁用应用商店_Win8系统当中Windows defnedder安全软件应该如何禁用?...
- BAT面试进阶:最全Memcached面试30题含答案
- tomcat安装-tomcat8.5
- [转]asp 无法连接 access,出现 -2147467259 未指定的错误
- python车牌识别系统开源代码_TensorFlow车牌识别完整版代码(含车牌数据集)
- 迅雷mac版精简安装教程
- electron选mysql的优缺点_大型Electron应用本地数据库技术选型
- 神经网络背后的数学原理:反向传播过程及公式推导
- C++ 侯捷视频学习(草稿)
- jhu研究生录取 计算机,背景一般获约翰霍普金斯大学JHU信息安全硕士录取
- “百度杯”CTF比赛2017年2月场WP--web
- 为什么顶尖高手,都是深度思考者?