R语言中的机器学习

Posted 2021-04-18 R语言交流中心

大家平时都会用到一些回归模型，今天我们来看一个集合多个模型建模和可视化的包mixomics。首先看下此包的所包含的方法列表：

然后看下其安装方法：

BiocManager::install("mixOmics")

 

接下来我们就来看下具体如何实现这些功能：

1. PCA分析，主要以线性方式结合多个独立变量进行降维分析。

其中主要的参数

然后我们看实例：

 

library(mixOmics)data(liver.toxicity)X <- liver.toxicity$geneMyResult.pca <- pca(X) plotIndiv(MyResult.pca, group =liver.toxicity$treatment$Dose.Group, legend = TRUE)

 

plotIndiv(MyResult.pca, ind.names = FALSE, group = liver.toxicity$treatment$Dose.Group, pch = as.factor(liver.toxicity$treatment$Time.Group), legend = TRUE, title = 'Liver toxicity: genes, PCA comp 1 - 2', legend.title = 'Dose', legend.title.pch = 'Exposure')

plot(MyResult.pca2)

 

plotLoadings(MyResult.pca)

从上图我们可以看出每个变量对每个组件的贡献(取决于对象)，将其展示在一个柱状图中表示，其中每个柱状图的长度对应于样本对组件的装载重量(重要性)。负载重量可以是正的，也可以是负的。

当然也可以对展示的数据进行设置：

plotLoadings(MyResult.pca, ndisplay = 100, comp=2, name.var = liver.toxicity$gene.ID[, "geneBank"], size.name = rel(0.3))

当然我们也可以从3D角度来看PCA 的分析结果，那样更能展示具体的分布，这里需要用到3主成分的降维：

 

MyResult.pca2 <- pca(X, ncomp = 3)plotIndiv(MyResult.pca2, group = liver.toxicity$treatment$Dose.Group, style="3d", legend = TRUE, title = 'Liver toxicity: genes, PCA comp 1 - 2 - 3')

 

2. 稀疏矩阵PCA分析，此函数相对与PCA多了keepX参数可以设置在每个组件起作用的前几个基因或者样本。

我们直接看下实例：

MyResult.spca <- spca(X, ncomp = 3,keepX = c(15,10,5)) # 1Run the methodplotIndiv(MyResult.spca, group =liver.toxicity$treatment$Dose.Group, #2 Plot the samples pch = as.factor(liver.toxicity$treatment$Time.Group), legend = TRUE, title = 'Liver toxicity: genes, sPCA comp 1 - 2', legend.title = 'Dose', legend.title.pch = 'Exposure')

plotVar(MyResult.spca, cex = 1)

上面这个图为圆形图，其实就是极坐标图，主要体现各变量之间的相关性分布。

在PCA还可以直接看下explainedvariance score（主要体现自变量对因变量的贡献度，介于【0，1】，值越大代表贡献越大），需要用到函数tune.pca(X)：

3. PLS-DA分析，虽然偏最小二乘法最初没有应用于分类和辨别问题。后来经过改造还是被用来进行分类研究。

我们直接看下实例：

data(srbct)X <- srbct$geneY <- srbct$classMyResult.splsda <- splsda(X, Y, keepX =c(50,50))plotIndiv(MyResult.splsda)

selectVar(MyResult.splsda, comp=1)$value #查看组件1中涉及的变量名称值。。

4. PLS分析，偏最小二乘回归(PLS)Wold, Sjostrom和Eriksson 2001)是一种多元方法，它涉及两个数据矩阵X(例如转录组)和Y(例如脂质)。通过对两个矩阵的结构建模，PLS超越了传统的多重回归。与传统的多元回归模型不同，它不局限于不相关的变量。PLS的优点之一是它可以处理许多有噪声的、共线性(相关)和缺失变量，还可以同时在Y中建模几个响应变量。

我们直接看下实例：

data(nutrimouse)X <- nutrimouse$gene Y <- nutrimouse$lipidMyResult.spls <- spls(X,Y, keepX = c(25,25), keepY = c(5,5)) plotIndiv(MyResult.spls)

plotIndiv(MyResult.spls, group =nutrimouse$genotype,
 rep.space = "XY-variate", legend = TRUE, legend.title = 'Genotype', ind.names = nutrimouse$diet, title = 'Nutrimouse: sPLS')

plotArrow(MyResult.spls,group=nutrimouse$diet,legend = TRUE,
 X.label = 'PLS comp 1', Y.label = 'PLS comp 2')

上面的图主要是指通过X，Y的组件将样本进行短箭头的链接。通过大量的短箭头可以看出样本在两个数据集之间良好的一致性。

下面看下各变量之间的关系的可视化：

X11()cim(MyResult.spls, comp = 1)

network(MyResult.spls, comp = 1)# network(MyResult.spls,comp = 1, cutoff = 0.6, save = 'jpeg', name.save = 'PLSnetwork')就可以进行对图像的保存。

另外我们还可以对模型的组件进行判别，最终确定需要多少各组件合适。判别用到了交叉验证的方法，具体实例如下：

set.seed(30) # for reproducbility in thisvignette, otherwise increase nrepeatperf.pls <- perf(MyResult.spls,validation = "Mfold", folds = 5, progressBar = FALSE, nrepeat= 10)plot(perf.pls$Q2.total)abline(h = 0.0975)

以上实例表明一个组件就可以。

最后是最优变量数的获取，需要用到和PCA类似的函数tune.pls。我们看下实例：

list.keepX <- c(2:10, 15, 20)# tuning based on MAEset.seed(30) # for reproducbility in thisvignette, otherwise increase nrepeattune.spls.MAE <- tune.spls(X, Y, ncomp =3, test.keepX =list.keepX, validation ="Mfold", folds = 5, nrepeat = 10,progressBar = FALSE, measure = 'MAE')plot(tune.spls.MAE, legend.position = 'topright')

图中Y轴MAE指的是平均绝对误差。平均绝对误差是绝对误差的平均值，能更好地反映预测值误差的实际情况.值越低代表此时模型越优化。

5. DIABLO分析，相当于是pls的扩展，可以在X中引入多个矩阵。

我们直接看下实例：

data(breast.TCGA)# extract training data and name each dataframeX <- list(mRNA =breast.TCGA$data.train$mrna, miRNA = breast.TCGA$data.train$mirna, protein = breast.TCGA$data.train$protein)Y <- breast.TCGA$data.train$subtypelist.keepX <- list(mRNA = c(16, 17),miRNA = c(18,5), protein = c(5, 5))MyResult.diablo <- block.splsda(X, Y,keepX=list.keepX)plotIndiv(MyResult.diablo, ind.names = FALSE, legend=TRUE, cex=c(1,2,3), title = 'BRCA with DIABLO')

 

plotVar(MyResult.diablo, var.names =c(FALSE, FALSE, TRUE), legend=TRUE, pch=c(16,16,1))

 

plotDiablo(MyResult.diablo, ncomp = 1)

从上图我们可以看出DIABLO在mRNA和蛋白质数据集之间提取了很强的相关性。

接下来看下变量之间的相关性：

 

circosPlot(MyResult.diablo, cutoff=0.7)

circos图基于相似矩阵构建的，表示不同类型变量之间的相关性图。

 

cimDiablo(MyResult.diablo, color.blocks =c('darkorchid', 'brown1', 'lightgreen'), comp = 1, margin=c(8,20),legend.position = "right")

 

最后还提供了ROC曲线的绘制：

Myauc.diablo <- auroc(MyResult.diablo,roc.block = "miRNA", roc.comp = 2)

 

至此整个包的功能介绍完了，具体的应用希望大家各显神通.

欢迎大家学习交流！