Python 散点图线性拟合_基于R软件实现多元线性回归

一个多元线性回归在R中的实现示例 在一元回归 中，只包含一个预测变量和响应变量间的关系。与此相比，当存在两个或以上的预测变量时，称为多元回归(Multiple Regression)。 如果只考虑变量间的线性关系时，就是多元线性回归( Multiple Linear Regression )。其基本思想是确定一个方程，给出一组预测变量( X )和一个响应变量( Y )之间的线性关系：

式中β₁-β_n代表了预测变量X₁-X_n的回归系数(斜率)，β₀是线性方程的截距。回归的目的是找到一组最优的模型参数(斜率和截距值的组合)，使响应变量的观测值和预测值之间的残差平方和最小化。

从技术上来说，简单多项式回归也可以视为多元线性回归的特例，例如二次回归有两个预测变量(X₁=X和X₂=X²)，三次回归有三个预测变量(X₁=X、X₂=X²和X₃=X³)。

本篇继续通过一个示例，展示在R语言中拟合多元线性回归的方法过程。

下文的示例数据、R代码等，可见网盘附件(提取码，54ua)：

https://pan.baidu.com/s/1ixB_7VQ5dXHu7ckyv1r-GQ

示例数据概要

节选自马里兰州生物资源调查(https://dnr.maryland.gov/streams/Pages/mbss.aspx)。

就节选部分数据为例，记录了所调查的马里兰州河流中每75米长的区段水域内，鱼类物种Rhinichthys cataractae的丰度，并测量了每段水域中相应的环境特征。

其中第一列代表了调查河流区段的位置信息，其余各列依次为：

fish，水域中鱼类物种Rhinichthys cataractae的个体数量，代表了物种丰度；

acre，水域流域面积(英亩，acre)；

do2，水域溶解氧含量(毫克/升，mg/L)；

depth，水域最大深度(厘米，cm)；

no3，水域硝酸盐浓度(毫克/升，mg/L)；

so4，水域硫酸盐浓度(毫克/升，mg/L)；

temp，水域温度(摄氏度，℃)。

一个生物学目的是探索可能影响鱼类物种Rhinichthys cataractae丰度的环境因素，并对物种丰度变化的原因作出假设。

线性关系无疑是最简单的类型，尽管实际中变量间关系并非线性这样简单，但由于线性关系比较直观且易于理解，大多数实际分析中我们更期望首选将问题简化为线性模型去解释。

简单线性回归初步探索变量间关系

因此，不妨首先站在单个预测变量与响应变量角度上，通过一元线性回归初步观察每个环境因素和物种丰度的线性关系。

#读取鱼类物种丰度和水体环境数据
dat  
#绘制二维散点图观测各环境变量与鱼类物种丰度的关系
library(ggplot2)
dat_plot  
p geom_point() +
facet_wrap(~variable, ncol = 3, scale = 'free') +
geom_smooth(method = 'lm')
p
#拟合各环境变量与鱼类物种丰度的一元回归，并提取各自 R2 和 p 值
env R2_adj p_value  
for (i in env) {
    fit_stat     R2_adj     p_value }
env_stat env_stat  #数据框中存储了各环境变量与鱼类物种丰度的一元回归的 R2 和 p 值
#在散点图中添加各环境变量与鱼类物种丰度的一元回归的 R2 和 p 值作为标识
#注：该 R2 和 p 值仅为单个变量一元回归的 R2 和 p 值，和下文即将提到的多元回归的 R2 和 p 值存在区别
env_stat$fish for (i in env) env_stat[i,'value'] env_stat$variable env_stat$label  
env_stat dat_plot$label dat_plot  
p + geom_text(data = dat_plot, aes(label = label), size = 3)

通过分别查看单个环境变量和物种丰度的关系，可以看到acre(流域面积)、no3(硝酸盐浓度)以及depth(水域深度)与物种丰度间的线性关系较明显，其余三个环境因素则很弱。

         

多元线性回归描述变量间线性响应

接下来拟合环境因素与物种丰度的多元线性回归，以综合看待环境因素的效应，为可能导致物种丰度变化的原因提供见解。

类似一元线性回归，多元线性回归在R中也是通过函数lm()执行。

#首先不妨使用全部环境变量拟合与鱼类物种丰度的多元线性回归
#自变量之间用“+”相连
fit1 summary(fit1)  #展示拟合方程的简单统计
#如果数据集中所有的自变量都使用到，且不考虑交互作用，则上式可以简化如下
#使用“.”代替所有自变量
fit1 summary(fit1)  #展示拟合方程的简单统计

整体而言，全模型p<0.001是显著的，表明多元线性回归中存在部分环境因素可以解释物种丰度变化。回归方程校正后的R²=0.25，算不上很高，这是可预期的，正如上所提到的，物种丰度和大部分环境因素间关系实际上并非线性这样简单。尽管如此，仍可以对结果作出一些解释。

关注单个预测变量的p值，p值越小代表该预测变量(环境因素)和响应变量(物种丰度)间线性关系越明显。例如，acre(流域面积)和no3(硝酸盐浓度)的p值显著，且回归系数大于0，表明在其它环境因素保持不变时，它们的增加有助于物种丰度的增加；depth(水域深度)接近临界显著性，暗示水深的增加也有助于提升物种的丰度；temp(水域温度)的p值不显著，表明物种的丰度并非随着温度的增加而呈线性增长。

通过比较，可以看到环境因素与物种丰度各自的一元回归和整体的多元线性回归之间趋势基本一致，尽管存在微小差异。

如果期望从中提取某些关键的信息，或者用作其它统计，这是一些参考用法。

#提取或统计重要数值，例如
coefficients(fit1)  #各变量的斜率和截距项
#也可以 names(summary(fit)) 后查看主要的内容项，然后从中提取，例如
summary(fit1)$adj.r.squared  #校正后 R2

此外，如期望作图观测变量间关系，由于多元线性回归是多维结构，很难直接在三维空间中表现出来。一种可选替代方案是通过多个一元回归，单独观测每个环境因素和物种丰度的关系，这个在上文已经展示了。但需注意的是，两种方法还是存在一些差异的，这个在上文比较二者就看出来了。

#作图，三维空间中最多展示 1 个响应变量和 2-3 个预测变量间的线性关系，例如
library(car)
scatter3d(fish~acre+no3, data = dat)
#要么就考虑使用多个一元回归散点图替代，但是它和多元线性回归是存在一些区别的

模型比较和选择的一个参考

由于未观察到物种丰度与do2(水域溶解氧含量)、so4(水域硫酸盐浓度)和temp(水域温度)呈线性响应，因此不妨：

(1)考虑转化数据，并继续观察转化后的变量间线性关系是否明显；

(2)直接把它们从回归方程中剔除。

转化数据需要思考很多问题，懒得思考，就直接选择把它们剔除了，来看新拟合的多元线性回归方程。

#只考虑使用 3 个和鱼类物种丰度线性关系较为明显的环境变量拟合多元线性回归
#这里将处于线性关系临界值的 depth 也算在内
fit2 summary(fit2)  #展示拟合方程的简单统计

整体而言，去除了3个不显著的预测变量(环境因素)后，方程的R²并未发生明显变化，也就是精度未降低，且全模型p<0.001仍然是显著的。

如果和上文使用全部环境因素拟合的最初的回归方程相比，肯定更期望选择后者，因为在二者精度相差不大时，后者更简约、更容易理解且解释轻松。

如果期望通过一些统计检验方法比较二者，以选择更优的一个，如下提供两个示例参考。

前后两个多元线性回归可视为嵌套模型，所谓嵌套模型，即其中之一的一些项完全包含在另一模型中。此时可以通过anova()比较两个嵌套模型的拟合优度，以检验不含do2(水域溶解氧含量)、so4(水域硫酸盐浓度)和temp(水域温度)的回归方程与包含这3个变量的回归方程的预测效果是否一样好。

还可以通过如AIC(Akaike Information Criterion，赤池信息准则)等值进行比较。AIC会在回归方程复杂性与其拟合优度之间进行权衡，AIC值较小的回归方程优先选择，表明较少的预测变量已经获得了足够的拟合度。

#anova() 比较两个嵌套模型的拟合优度
anova(fit2, fit1)
#AIC 比较两个回归方程
AIC(fit2, fit1)

anova()对是否应该添加变量do2、so4和temp到多元线性回归方程中进行了检验，由于p>0.05，表明添加这3个变量前后两回归方程的差别不大，表明可以将这3个变量直接删除。类似地，AIC值也显示剔除这3个变量更佳。

其它

本篇未涉及预测变量(环境因素)间的交互作用，带交互项的多元回归在放下节展示。需注意，考虑交互项时，就不再是线性方程了。

本篇未探讨预测变量(环境因素)间可能存在的共线性问题，它比较复杂，关于共线性以及变量选择等过程，将在下下节细讲。

参考资料

http://www.biostathandbook.com/multi p leregression.html Robert I. Kabacoff. R 语言实战(第二版)(王小宁   刘撷芯   黄俊文   等   译) .  人民邮电出版社 , 2016. 友情链接 这张统计图，强烈推荐，包括了直方图、密度曲线、散点图、趋势线、相关系数及显著性 2020-07-27 工作科研学习，效率太低，是因为浏览器没这么用 2020-07-24 听说你的Pubmed插件Scholarscope不能用了 2020-07-22 【修订版】对中文期刊友好的，临床基线资料table快速生成，不学就OUT了 2020-07-09 看过许多统计教程，这篇是我最推荐的 2020-05-25 End 如你有以下问题： 数据清洗、数据整理、统计分析、样本量计算、诊断试验、ROC曲线、临床预测模型、Graphpad作图、R语言、研究设计、投稿选刊、文章修回等问题，均可向刘老师咨询(微信号：tbwhere2，咨询形式可微信沟通)。