对那些基于降维（特征、成分）的回归方法作个简介
对基于降维（特征、成分）的回归方法作个简介mp.weixin.qq.com
很多生物学数据集具有高维的特征，体现在变量数量非常多且存在高度共线性，难以通过常规的多元回归分析变量间关系。在前文基于相似或相异度矩阵的多元回归（MRM）中，提到了一种替代方法是可以根据多变量数据计算所有样本之间的相似度或距离，然后基于相似或相异度矩阵执行回归分析，解释一组变量对另一组变量的整体效应。本篇列举另一种方法，基于降维的回归。备注：这个方法的名称是白鱼同学瞎取的，因为实在不知道该怎么称呼，有些片面......基于降维的回归、基于特征的回归、基于成分的回归等等，怎么叫都行吧......算了名称无所谓了，看了下文的方法过程描述就知道是怎么一回事了。对基于降维（特征、成分）的回归方法作个简介很多生物学数据集具有高维的特征，体现在变量数量非常多且存在高度共线性，难以通过常规的多元回归分析变量间关系。在前文基于相似或相异度矩阵的多元回归（MRM）中，提到了一种替代方法是可以根据多变量数据计算所有样本之间的相似度或距离，然后基于相似或相异度矩阵执行回归分析，解释一组变量对另一组变量的整体效应。本篇列举另一种方法，基于降维的回归。备注：这个方法的名称是白鱼同学瞎取的，因为实在不知道该怎么称呼，有些片面......基于降维的回归、基于特征的回归、基于成分的回归等等，怎么叫都行吧......算了名称无所谓了，看了下文的方法过程描述就知道是怎么一回事了。
基于降维（特征、成分）的回归的计算原理
这类方法是个统称。简单来说，就是首先自变量进行降维处理，将原始的大量自变量减少为较少的一组不相关成分，这些成分代表了原始自变量的特征组合，并对这些成分（而非原始的自变量）执行回归的过程。
基于降维的回归方法中，最具代表性的一个子类别就是主成分回归（Principal Component Regression，PCR），它是主成分分析（PCA）和普通最小二乘（OLS）回归的结合，最简单也最容易理解。就以主成分回归为例，包含三个步骤：
（1）使用所有的自变量运行主成分分析（PCA），对自变量进行降维，获得一组互不相关（正交）的主成分，这些主成分承载了原始变量的方差；
（2）使用主成分代替原始的自变量，运行普通最小二乘（OLS）回归，也就是通过线性回归建立这些主成分与响应变量的线性关系；
（3）计算与输入变量相对应的模型参数。
其它类别的基于降维的回归方法在主成分回归的基础上延伸得来。实际的应用中，根据具体情况选择合适的方法，总之非常灵活多样：
（1）降维方法不局限于PCA，也可以是主坐标分析（PCoA）、非度量多维标度（NMDS）等；
（2）回归方法不局限于线性回归，也可以是其它类型的非线性回归等。
基于降维（特征、成分）的回归的特点
下文简单概括这类方法的特点。
在一定程度上降低模型的复杂性
降维后的成分有效代表了原始变量的整体特征。通过回归拟合重要成分与响应变量的关系，该过程减少了模型中变量的数量，因此也相应地降低了模型的复杂性，这总的来说是一个优点。
但是很多情况下，成分的选择需要谨慎。并非每个成分都与响应变量有关，对于不重要的成分最好不予考虑，并且成分的选择不宜过多，否则难以体现降维的优势。一般而言，在前几个降维后的成分中进行选择就可以了，它们通常承载了变量集的绝大部分方差，并且与响应变量最相关的成分也通常包含在内。
这种“降低模型的复杂性”是如何体现的，可以参考下文中的文献实例一。
利于更好的趋势特征建模
很多分析中，为了更好的描述趋势响应关系，通常会排除一些不相干变量。而这类基于降维的回归方法，在这种情况下会显得很有用。
降维的一个明显的优势是，可以通过其中2-3个成分绘制散点图，利于观察对象的分布趋势。这样就可以将趋势明显的一组成分提取出来，并用它们建模以解释生物学问题。特别是当某一成分即有效代表了数据集中对象的整体响应趋势时，使用该成分的建模结果往往表现非常良好。
这种“利于更好的趋势特征建模”是如何体现的，可以参考下文中的文献实例二。
避免多重共线性以及缓解过拟合问题
此外，由于降维后的各成分之间相互正交（不相关），直接避开了原始变量之间的多重共线性问题，相比常规的基于变量的回归，可以无需考虑变量选择。
同时使用少数成分代替原始大量变量，也在另一层面上减少了应用于回归中的变量数量，缓解了当变量数量远大于对象数量时出现的过拟合问题。
可以实现高维数据结构的可视化
由于结合了降维过程，就可以像PCA等那样，允许用图表来描述多维数据结构，例如成分载荷图、对象和变量的相关图以及得分图等，更有助于描述变量之间的关系，并提供了关于对象相似度和数据集整体结构的信息。而对于常规的基于变量的回归而言，大量的变量间关系则难以通过图表直接描述。
参考下文中的文献实例三，展示了某种类型的可视化方式。
缺点是可能不利于解释变量效应
在这类方法中，成分实际上视为一组虚拟变量或者潜变量。尽管基于成分的回归可以降低模型的复杂性，但由于成分并非具体的实际变量，可能会难以解释原始自变量的效应。
因此，使用这类方法时，应尽可能避免从原始变量角度解释问题，尽可能基于成分所代表的特征去解释。但若仍期望寻找在某成分上具有代表性的变量，通常可以使用相关性分析建立成分和变量的关联，尽可能选择一些相关系数绝对值高且自身方差更大的变量。
列举一些文献中的实例
就以群落分析为例吧，白鱼同学比较了解这种。很多研究期望了解物种多样性和环境或功能的关系，然而环境中的物种数不胜数，如果仅从生物统计的角度入手，单独分析各物种丰度和环境或功能的关联肯定难以实现。这时候，一种可选方法是通过对物种多度数据进行降维，特别是在基于β多样性的有关分析中非常有用。
很多同学咨询过，看到很多微生物组的相关文献中，使用微生物群落的β多样性进行建模，然后描述β多样性和环境、功能等的关系。然后就很好奇，说一般β多样性不是计算的群落物种组成相似或相异度矩阵吗？结果是一个对称矩阵，不是一组变量，该怎样用于建立群落β多样性与环境、功能等的关系？类似这种：
然后白鱼同学通常回复两种可选思路：
（1）有些回归方法是可以基于相似或相异度矩阵进行建模的，就可以对β多样性求个相异矩阵，环境变量也求个相异矩阵，然后执行矩阵式回归；
（2）通过本篇提到的这类基于降维的回归去解决。
关于（1）就不多提了。关于（2），简单找了几篇文献中的例子，在这里帮助大家了解这类基于降维的回归是如何实际应用的。
实例一
β多样性平时都是怎样表示的，最常见的就是降维方法，如PCA、PCoA、NMDS等，对不对？所以，一种常见的分析β多样性和环境因子的思路就建立在该方法之上。可以首先基于群落物种组成多度数据，选择一种合适的方法进行降维，随后就可以将样本在各成分轴上的坐标提取出来，这些成分就可视为代表“群落结构β多样性”变量。由于各成分之间正交，因此各成分也反映了潜在的不同方式的物种-环境关联模式。最后用于建立这些成分和环境因子回归。
Jiao等（2018）研究了退耕还林过程中土壤微生物群落的结构组成及其对土壤营养物质循环的驱动效应。这里节选了文中章节“Potential drivers of soil multi-nutrient cycling in reforested ecosystems”的一小部分内容，该部分阐述了土壤养分循环的潜在微生物驱动力。
作者建立了微生物因子与土壤养分因子的回归模型，以评估不同类型微生物的多样性对土壤养分循环的相对重要性。微生物因子以微生物群落的α和β多样性指数衡量，α多样性指数直接计算Richness和Shannon指数；对于β多样性指数的获得，首先基于群落物种组成丰度执行NMDS，随后提取两个轴上的坐标以量化各群落的β多样性指数，两组坐标即代表了两类“群落结构β多样性”成分。土壤养分因子使用土壤多养分循环指数衡量，并通过随机森林（RF）建立了包含各类型微生物群落α和β多样性与土壤多养分循环指数的回归。
结果显示，在全范围的垂直土壤剖面中，贡献土壤养分循环能力的不同类型微生物多样性指数的重要程度上，细菌>古细菌>真菌。当考虑土壤分层时，不同类型微生物多样性指数对土壤多种养分指数的贡献在表层（Superficial）和深层（Deep）土壤中存在区别，古细菌β多样性更好地预测了表层土壤的动态，但细菌β多样性在深层土壤中的作用更高。
节选展示了Jiao等（2018）的图5 A部分。
图A，展示了细菌、古细菌和真菌群落的α和β多样性在贡献土壤多养分循环指数中的重要性；Whole不考虑土壤分层，Superficial代表表层土壤，Deep代表深层土壤；微生物的重要性以随机森林中的“percentage of increase of mean square error”指示值衡量。
实例二
有时为了更好地描述趋势响应关系，可能就需要只选择1-2个少数重要成分。如果观察到样本在这些成分上的分布规律与响应变量的变化趋势存在比较一致的特征，那么使用这些成分的建模就会更加理想。还是以β多样性为例，首先观察各样本点在PCA、PCoA、NMDS等图中的分布，如果发现这些样本沿某一特征轴存在明显的与环境梯度一致的趋势，并且该特征轴的方差贡献也不是很低时，就可以将这些样本在该特征轴上的坐标提取出来。这组坐标就可视为代表“有规律的β多样性梯度”的成分变量，随后用于建立它们和环境因子回归。
Edwards等（2018）研究了水稻根系微生物组与水稻发育时间的关系。文中章节“Developmental stage correlates with microbiota succession”的部分内容提到，在影响水稻根系微生物群落组建的因素中，观察到植物生长时间与基因型之间存在显著的统计学交互效应，并且植物发育阶段（营养生长、生殖生长等特定阶段）比植物年龄（单纯的生长天数）更重要。
作者首先对水稻根系（包括根际、根表和根内）微生物群落进行PCoA分析后，发现样本沿第二PCoA轴上的分布特征与水稻的生长时间的趋势相吻合，因此第二PCoA轴就代表了一组与水稻生长时间密切相关的“根系微生物成分”。随后，为了探讨水稻基因型和生长时间与根系微生物群落β多样性的关系，将样本在该轴上的坐标提取出作为响应变量，并考虑水稻基因型和生长时间的交互效应作为自变量，建立二者的回归。
而样本在第一PCoA轴上的分布与植物生长时间无规律，就未考虑。其它PCoA轴的方差贡献很低，或者趋势也不明显，故也不考虑。
回归结果显示水稻基因型和生长时间的交互效应是显著的。作者假设，如果植物发育速率对根系相关微生物群有影响，那么发育较快的水稻品种比发育较慢的水稻品种具有更大的回归系数（斜率），结果观察到在水稻生长的前84天中（水稻主要处于营养生长期），不同基因型水稻的斜率显著不同。但在随后的生长期（水稻主要处于生殖生长期）中则没有，这在某种程度上体现了发育阶段比植物年龄对根系微生物群落的组建更重要。
节选展示了Edwards等（2018）的图6 B、C部分。
图B，显示了水稻根际、根表和根内微生物群落组成的PCoA。
图C，使用第二PCoA轴的坐标代表“时间响应型根系微生物成分”，建立其与不同品种水稻生长时间的关系，图中展示了回归的斜率。
实例三
下面的这种图在文献中很常见，使用非约束排序方法（常见如PCA、PCoA、NMDS等，下图使用的NMDS）描述群落物种组成结构，是β多样性可视化的一种方式。但是与常规的NMDS图不同的是，除了样方点外，有时还在图中添加了环境变量的梯度信息，以类似等高线的方式表现出来，直观地呈现了不同样方之间在环境组成上的差异。
通常，要描述群落结构与环境梯度的关系，单独分析各物种肯定不现实，因为群落中每种物种对环境的响应都不同。这时候，就可以在β多样性角度上，分析整体的趋势。首先使用NMDS实现对物种多度数据的降维，将样方投影在降维后的二维空间中，以此表示各群落之间β多样性的差异。简单来说，点越近表示群落物种组成越相似，点越远表示群落物种组成差异越大，想必这点大家都没疑问。随后为了表示群落组成结构与环境梯度的关系，将两个NMDS轴上的样方坐标点提取出来，这样就获得了两组“群落结构β多样性”成分，然后建立这两个成分（作为自变量）与各环境变量（作为响应变量）的回归。最后将回归结果绘制在图中，这样就将各样方所在的环境梯度信息表现出来了。通常而言，回归若是线性的，环境梯度的展示方式就以带箭头的向量从(0,0)点引出；回归若是平滑的，就是这种等高线类型。最后，NMDS图既呈现了各样方间群落结构的β多样性水平，又反映出了群落所处的环境特征。
节选展示了Ho等（2019）的图4部分，展示了鱼类群落结构组成与环境梯度的关系。
关于R语言的计算方法
如上文所述，这类基于成分的回归方法是一个统称，大致分为两步，第一步降维，第二步回归。无论降维方法或回归方法的选择都丰富多样，实际使用时最好参考相关文献确定分析思路。
对于执行，通常也是两步分开执行。白鱼同学懒得写R代码了，直接摆一些相关方法的链接。
先前整理了一些常见的降维分析方法，可见专辑“多变量的降维分析”（注意RDA、CCA、db-RDA等约束排序方法除外，原因参考下文）。
以及常见的回归方法，可见专辑“回归与建模”。
其它补充
RDA、CCA、db-RDA（CAP）等方法其实也是降维+回归的组合，尽管很相似，但和本篇探讨的方法过程存在不同，实际使用时注意区分。特别是db-RDA，确实很容易搞混，因此下面作了简单的比较。
关于上文介绍的方法
本篇探讨的基于降维的回归过程中，降维后的成分作为自变量使用。就以物种和环境的关系为例，对物种多度数据降维后，以物种成分为自变量，环境作为响应变量看待。这种情况下，解释的是物种对环境功能的贡献（参考上文示例文献一），或者物种对环境特征的预测或指示作用（参考上文示例文献二）。
关于RDA、CCA、db-RDA等方法
RDA、CCA、db-RDA（CAP）等方法和上述过程相比，某些操作过程恰好相反，它们以物种作为响应变量，环境作为自变量，主要用于解释物种变化是如何随环境变化作出响应的。
不妨先来看下RDA的过程，对于给定的两个矩阵，物种多度变量矩阵和环境变量矩阵，RDA过程如下。（比较和上文方法的区别，它是先回归，后降维，物种作为响应变量，环境作为自变量）
（1）执行物种spe1与环境变量env1的线性回归，将回归模型拟合的物种丰度值存储在拟合值矩阵，并将物种丰度的残差存储在残差矩阵。如此对物种组成矩阵中的所有物种重复相同的操作，最终获得包含所有物种丰度拟合值及残差的两个矩阵。
（2）回归过程执行完毕后，使用PCA对拟合值矩阵进行降维获得约束轴。所以，RDA本身就是线性模型的延伸，RDA的整体方差解释率也就等于线性回归的R2。残差通常不考虑，但如果感兴趣，同样使用PCA对残差值进行降维获得残差轴。
然后再来看下db-RDA（CAP）的过程。（比较和上文方法的区别，它是先降维，再回归，最后还有一步降维；降维后的物种成分作为响应变量，环境作为自变量）
（1）基于物种多度变量矩阵，计算相异度矩阵，并执行PCoA。
（2）将PCoA获得的特征轴（成分）提取出，作为响应变量，环境作为自变量，执行普通RDA的过程。
参考文献
https://www.xlstat.com/en/solutions/features/principal-component-regressionhttps://iq.opengenus.org/principal-component-regression/Ho J K, Ramchunder S J, Memory A, et al. Native and introduced fish community structure in a freshwater swamp forest: Implications for conservation and management. Aquatic Conservation-marine and Freshwater Ecosystems, 2019, 29(1): 47-58.Jiao S, Chen W, Wang J, et al. Soil microbiomes with distinct assemblies through vertical soil profiles drive the cycling of multiple nutrients in reforested ecosystems. Microbiome, 2018, 6(1): 1-13.Edwards J, Santosmedellin C, Liechty Z, et al. Compositional shifts in root-associated bacterial and archaeal microbiota track the plant life cycle in field-grown rice. PLOS Biology, 2018, 16(2).
友情链接
R包randomForest的随机森林回归模型以及对重要变量的选择mp.weixin.qq.com
559
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
