R语言基础| 回归分析

写在前面

R语言拥有丰富的数据处理、统计分析和机器学习工具包，涵盖了从简单的描述统计到复杂的模型建立的各个方面。再加上数据的处理可以完美的衔接后续的可视化，这使得它成为处理各种类型和规模的数据集的理想选择。回归分析是统计学中一种用于探究自变量与因变量之间关系的方法。在 R 语言中，回归分析可以通过内置的lm()函数来实现

更多R语言教程和测试数据可见：R语言基础学习手册

8.1 定义与分类

回归分析是一种用于建立和探索变量之间关系的统计方法。它用于预测一个或多个自变量（也称为解释变量）对因变量（也称为响应变量）的影响程度。

不同类型的回归分析：

本章重点讨论最小二乘（OLS)回归法，是一种常用的统计学方法，用于建立因变量与一个或多个自变量之间的线性关系模型。包括简单线性回归，多项式回归和多元线性回归。

简单线性回归：模型包含一个自变量和一个因变量。

多项式回归：只有一个自变量但同时包含变量的幂(比如，x,x2,x3)。

多元线性回归：当不止一个自变量时。

8.2 OLS回归

8.2.1 函数lm()拟合回归模型

8.2.1.1语法

myfit <- lm(formula,data)

myfit:是你为线性回归模型分配的变量名，可以根据需要自定义这个变量名。

lm():是R中用于拟合线性模型的函数。它接受两个主要参数：formula 和 data。

formula:是描述回归模型的公式，它指定了因变量和自变量之间的关系。公式通常采用类似于Y ~ X1 + X2+….Xk 的格式，其中 Y 是因变量，X1 和 X2 是自变量。

data:是包含你的因变量和自变量的数据框或数据集。你需要指定数据的名称，以便模型可以从中获取变量。

8.2.1.2 formula

其中formula形式为

Y ~ X1 + X2+....Xk

～左边为因变量，右边为各个自变量，自变量之间用符号+来分隔，下表中的符号也用不同方式修改这一表达式：

8.2.1.3 拟合lm()函数后用下面函数获取信息

8.2.2 简单线性回归实例展示

先回忆一下公式：

myfit <- lm(y~x,data)

数据集：系统已安装的women提供了年龄在30-39间女性的身高体重信息。现想通过身高来预测体重，通过获得一个方程可以分辨出那些过重或过轻的个体。

women

##    height weight
## 1      58    115
## 2      59    117
## 3      60    120
## 4      61    123
## 5      62    126
## 6      63    129
## 7      64    132
## 8      65    135
## 9      66    139
## 10     67    142
## 11     68    146
## 12     69    150
## 13     70    154
## 14     71    159
## 15     72    164

这里用身高预测体重，因此身高为自变量，体重为因变量。

fit <- lm(weight~height,women)
summary(fit)

## 
## Call:
## lm(formula = weight ~ height, data = women)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -1.7333 -1.1333 -0.3833  0.7417  3.1167 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -87.51667    5.93694  -14.74 1.71e-09 ***
## height        3.45000    0.09114   37.85 1.09e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.525 on 13 degrees of freedom
## Multiple R-squared:  0.991,  Adjusted R-squared:  0.9903 
## F-statistic:  1433 on 1 and 13 DF,  p-value: 1.091e-14

通过summary()函数显示，斜率为3.45，截距为-87.52，因此该方程为：

weight=3.45*height-87.52

Multiple R-squared:0.991表示这个模型可以解释99.1%的数据，这个值在0-1之间，越大越接近于1说明模型拟合的越好。

可以利用plot()函数对模型作图，生成四幅图：

1）残差拟合图

2）正态的QQ图

3）大小位列图

4）残差影响图

都是用来评价拟合模型的。

plot(fit)

abline()绘制拟合曲线，必须在高级绘图的基础上使用：

plot(women$height,women$weight)
abline(fit)

8.2.3 多项式回归实例展示

一般来说，n次多项式生成一个带有n-1个弯曲的曲线，但几乎高于三次的没有必要。

二次和三次项语法：

fit2 <- lm(y~x+I(x^2),data)
fit3 <- lm(y~x+I(x^2)+I(x^3),data)

从上面的拟合曲线可以看到，只有2个数据点在拟合曲线上，因此我们可以通过添加一个二项式来提高回归的精度。

fit2 <- lm(weight~height+I(height^2),women)
summary(fit2)

## 
## Call:
## lm(formula = weight ~ height + I(height^2), data = women)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.50941 -0.29611 -0.00941  0.28615  0.59706 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 261.87818   25.19677  10.393 2.36e-07 ***
## height       -7.34832    0.77769  -9.449 6.58e-07 ***
## I(height^2)   0.08306    0.00598  13.891 9.32e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3841 on 12 degrees of freedom
## Multiple R-squared:  0.9995, Adjusted R-squared:  0.9994 
## F-statistic: 1.139e+04 on 2 and 12 DF,  p-value: < 2.2e-16

plot(women$height,women$weight)
lines(women$height,fitted(fit2))

可以绘制好几条线：

plot(women$height,women$weight)
lines(women$height,women$weight,col="red")
lines(women$height,fitted(fit),col="black")
lines(women$height,fitted(fit2),col="blue")

可以看到用多项式回归后的蓝色预测值[fitted(fit2]和红色真实值[women$weight]几乎重合了，说明拟合的效果非常好。

8.2.4多元线性回归

当自变量大于1个时，需要用到多元线性回归。当做该回归时注意几点：

1.lm()函数需要的是数据框，as.data.frame()函数

2.首先检查自变量之间的相关性

3.当存在明显的相关性时可以考虑引用交互项

举例一

要研究一个州的犯罪率和其他因素的关系，包括人口、文盲率、收入和结霜天数，数据来自R中自带的state.x77数据集。

states <- as.data.frame(state.x77[,c("Murder","Population","Illiteracy","Income","Frost")])
library(psych)
corr.test(states)

## Call:corr.test(x = states)
## Correlation matrix 
##            Murder Population Illiteracy Income Frost
## Murder       1.00       0.34       0.70  -0.23 -0.54
## Population   0.34       1.00       0.11   0.21 -0.33
## Illiteracy   0.70       0.11       1.00  -0.44 -0.67
## Income      -0.23       0.21      -0.44   1.00  0.23
## Frost       -0.54      -0.33      -0.67   0.23  1.00
## Sample Size 
## [1] 50
## Probability values (Entries above the diagonal are adjusted for multiple tests.) 
##            Murder Population Illiteracy Income Frost
## Murder       0.00       0.09       0.00   0.43  0.00
## Population   0.01       0.00       0.46   0.43  0.09
## Illiteracy   0.00       0.46       0.00   0.01  0.00
## Income       0.11       0.15       0.00   0.00  0.43
## Frost        0.00       0.02       0.00   0.11  0.00
## 
##  To see confidence intervals of the correlations, print with the short=FALSE option

fit <- lm(Murder~Population+Illiteracy+Income+Frost,states)
summary(fit)

## 
## Call:
## lm(formula = Murder ~ Population + Illiteracy + Income + Frost, 
##     data = states)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.7960 -1.6495 -0.0811  1.4815  7.6210 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 1.235e+00  3.866e+00   0.319   0.7510    
## Population  2.237e-04  9.052e-05   2.471   0.0173 *  
## Illiteracy  4.143e+00  8.744e-01   4.738 2.19e-05 ***
## Income      6.442e-05  6.837e-04   0.094   0.9253    
## Frost       5.813e-04  1.005e-02   0.058   0.9541    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.535 on 45 degrees of freedom
## Multiple R-squared:  0.567,  Adjusted R-squared:  0.5285 
## F-statistic: 14.73 on 4 and 45 DF,  p-value: 9.133e-08

这个例子中不难看出自变量之间存在着显著相关性，但没有考虑交互项的引用，或者改用其他方法拟合，导致Multiple R-squared: 0.567，即只能符合56.7%的真实值。

举例二，引入交互项

以mtcar数据框中的汽车数据为例，研究汽车重量和发动机功率（自变量）对每加仑汽车行驶英里数（因变量）之间的关系。同时，考虑到汽车重量和发动机功率之间也存在关系，对这两个自变量引入交互项：

library(psych)
data <- as.data.frame(mtcars[,c("mpg","hp","wt")])
corr.test(data)

## Call:corr.test(x = data)
## Correlation matrix 
##       mpg    hp    wt
## mpg  1.00 -0.78 -0.87
## hp  -0.78  1.00  0.66
## wt  -0.87  0.66  1.00
## Sample Size 
## [1] 32
## Probability values (Entries above the diagonal are adjusted for multiple tests.) 
##     mpg hp wt
## mpg   0  0  0
## hp    0  0  0
## wt    0  0  0
## 
##  To see confidence intervals of the correlations, print with the short=FALSE option

fit <- lm(mpg~hp+wt+hp:wt,data)
summary(fit)

## 
## Call:
## lm(formula = mpg ~ hp + wt + hp:wt, data = data)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.0632 -1.6491 -0.7362  1.4211  4.5513 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 49.80842    3.60516  13.816 5.01e-14 ***
## hp          -0.12010    0.02470  -4.863 4.04e-05 ***
## wt          -8.21662    1.26971  -6.471 5.20e-07 ***
## hp:wt        0.02785    0.00742   3.753 0.000811 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.153 on 28 degrees of freedom
## Multiple R-squared:  0.8848, Adjusted R-squared:  0.8724 
## F-statistic: 71.66 on 3 and 28 DF,  p-value: 2.981e-13

结果显示：hp:wt的Pr(>|t|) 具有显著性，说明因变量与其中一个自变量的关系依赖于另外一个自变量的水平。

注：我这里还是先看了相关性，也可以不做这一步，直接看交互项的Pr(>|t|)，若其小于0.05即说明交互显著。

8.2.5 图形展示交互的结果（effects包中的effect()函数）

语法：

plot(effect(term，mod,,xlevels),line=c()),multiline=TRUE)

term:模型要绘制的项；

mod:通过lm()拟合的模型,注意mod后面有2个，，

xlevels:是一个列表，指定变量要设定的常量值；

line=c():1实线，2虚线，3点线；

multiline=TRUE：添加直线。

举例：

library(effects)

## Warning: 程辑包'effects'是用R版本4.3.2 来建造的

## lattice theme set by effectsTheme()
## See ?effectsTheme for details.

plot(effect("hp:wt",fit,,list(wt=c(2.2,3.2,4.2))),line=c(1,2,3),multiline=TRUE)

从图上可以看到，随着汽车重量（wt）的增加，发动机功率与每加仑汽油行驶英里数（mpg）的关系减弱了。当wt=4.2时，直线几乎是水平的，表明随着发动机功率（hp）的增加，mpg不再发生改变。

8.3 回归模型的诊断

主要为了评价回归模型适用性的工具。虽然summary()函数对模型进行了整体描述，但没有提供关于模型在多大程度上满足统计假设的任何信息。例如OLS回归需要满足：

1）正态性：自变量值固定时，因变量成正态分布。

2）独立性：要求因变量值相互独立。例如前面用的身高体重的例子中，A的体重通常不会影响B的体重，但如果在家族中做这样的统计则可能不那么独立了。

3）线性：因变量与自变量呈线性关系（直线/曲线）。

4）同方差性（不变方差）：因变量的方差不随自变量的水平不同而变化。

8.3.1 评价OLS拟合的标准方法——plot()函数生成四幅图

8.3.1.1语法展示

fit <- lm(weight~height,women)
par(mfrow=c(2,2))#将4张图放在一张大图中，成2x2排列
plot(fit)

8.3.1.2 如何根据这四幅图来评价OLS回归

1)Residuals vs Fitted:残差vs预测(拟合)图，评价线性：

点为残差值的分布，残差值随着拟合值的增加而随机分布在0附近，无明显模式或趋势则认为线性合理。若点呈曲线分布（如上述例子），可能存在二次项。

2)Nomal Q-Q:正态Q- Q图，评价正态性：

若满足正态假设，图上的点应该落在呈45度角的直线上；若非如此则违反正态性。

3)Scale-Location:位置-尺度关系图，评价同方差性：

水平线周围的点随机分布则满足该假设。

4)Residuals vs Leverage:残差-杠杆关系图，评价离群点、高杆杠值、强影响点：

独立性无法从图中获得，只能从数据的收集方式验证。

8.4 异常观测值

8.4.1 离群点

指那些模型预测效果不佳的观测点。残差值（观测值-预测值）大于2或小于-2则被认为是离群点。

正残差值说明模型对因变量的预测过小，负残差值说明模型对因变量预测过大。

8.4.1.2 检验方法-car包的outlierTest()函数

语法：

library(car)
outlierTest(fit)

举例：

states <- as.data.frame(state.x77[,c("Murder","Population","Illiteracy","Income","Frost")])
fit <- lm(Murder~Population+Illiteracy+Income+Frost,states)
library(car)

## Warning: 程辑包'car'是用R版本4.3.2 来建造的

## 
## 载入程辑包：'car'

## The following object is masked from 'package:psych':
## 
##     logit

## The following object is masked from 'package:dplyr':
## 
##     recode

outlierTest(fit)

##        rstudent unadjusted p-value Bonferroni p
## Nevada 3.542929         0.00095088     0.047544

当P<0.05时为离群点，这里Nevada的P=0.048，因此Nevada为离群点。

8.4.2 高杆杠值点

与其他自变量有关的离群点，是由许多异常的自变量值组合起来的，与因变量没有关系。

8.4.2.1 检验方法-帽子统计量(hat statistic)

对于一个给定的数据集，帽子均值为p/n，

p:模型估计的参数数目(包含截距项)，也可以理解为自变量的个数或特征的数量，length(cofficients(fit)))；

n:观测值数量，代表了样本中的数据点数量，也可以理解为样本量(length(fitted(fit)))；

若观测点的帽子值大于帽子均值的2或3倍，即可判定为高杠杆值点。

语法(hat.plot()函数也是car包中的)：

hat.plot <- function(fit){
  p <- length(coefficients(fit))
  n <- length(fitted(fit))
  plot(hatvalues(fit),main="Index Plot of Hat Values")
  abline(h=c(2,3)*p/n,col="red",lty=2)
  identify(1:n,hatvalues(fit),labels=names(hatvalues(fit)))#需配合Plots中Finish键选中显示
}
hat.plot(fit)

注：h表示”Hat Values”（帽值）的阈值。帽值是通过杠杆统计量（Leverage Statistics）计算得出的，用于检测在线性回归模型中数据点对回归方程的影响程度。

举例：

states <- as.data.frame(state.x77[,c("Murder","Population","Illiteracy","Income","Frost")])
print(states)

##                Murder Population Illiteracy Income Frost
## Alabama          15.1       3615        2.1   3624    20
## Alaska           11.3        365        1.5   6315   152
## Arizona           7.8       2212        1.8   4530    15
## Arkansas         10.1       2110        1.9   3378    65
## California       10.3      21198        1.1   5114    20
## Colorado          6.8       2541        0.7   4884   166
## Connecticut       3.1       3100        1.1   5348   139
## Delaware          6.2        579        0.9   4809   103
## Florida          10.7       8277        1.3   4815    11
## Georgia          13.9       4931        2.0   4091    60
## Hawaii            6.2        868        1.9   4963     0
## Idaho             5.3        813        0.6   4119   126
## Illinois         10.3      11197        0.9   5107   127
## Indiana           7.1       5313        0.7   4458   122
## Iowa              2.3       2861        0.5   4628   140
## Kansas            4.5       2280        0.6   4669   114
## Kentucky         10.6       3387        1.6   3712    95
## Louisiana        13.2       3806        2.8   3545    12
## Maine             2.7       1058        0.7   3694   161
## Maryland          8.5       4122        0.9   5299   101
## Massachusetts     3.3       5814        1.1   4755   103
## Michigan         11.1       9111        0.9   4751   125
## Minnesota         2.3       3921        0.6   4675   160
## Mississippi      12.5       2341        2.4   3098    50
## Missouri          9.3       4767        0.8   4254   108
## Montana           5.0        746        0.6   4347   155
## Nebraska          2.9       1544        0.6   4508   139
## Nevada           11.5        590        0.5   5149   188
## New Hampshire     3.3        812        0.7   4281   174
## New Jersey        5.2       7333        1.1   5237   115
## New Mexico        9.7       1144        2.2   3601   120
## New York         10.9      18076        1.4   4903    82
## North Carolina   11.1       5441        1.8   3875    80
## North Dakota      1.4        637        0.8   5087   186
## Ohio              7.4      10735        0.8   4561   124
## Oklahoma          6.4       2715        1.1   3983    82
## Oregon            4.2       2284        0.6   4660    44
## Pennsylvania      6.1      11860        1.0   4449   126
## Rhode Island      2.4        931        1.3   4558   127
## South Carolina   11.6       2816        2.3   3635    65
## South Dakota      1.7        681        0.5   4167   172
## Tennessee        11.0       4173        1.7   3821    70
## Texas            12.2      12237        2.2   4188    35
## Utah              4.5       1203        0.6   4022   137
## Vermont           5.5        472        0.6   3907   168
## Virginia          9.5       4981        1.4   4701    85
## Washington        4.3       3559        0.6   4864    32
## West Virginia     6.7       1799        1.4   3617   100
## Wisconsin         3.0       4589        0.7   4468   149
## Wyoming           6.9        376        0.6   4566   173

fit <- lm(Murder~Population+Illiteracy+Income+Frost,states)
hat.plot <- function(fit){
  p <- length(coefficients(fit))
  n <- length(fitted(fit))
  plot(hatvalues(fit),main="Index Plot of Hat Values")
  abline(h=c(2,3)*p/n,col="red",lty=2)
  identify(1:n,hatvalues(fit),labels=names(hatvalues(fit)))
}
hat.plot(fit)

## integer(0)

注：1:n表示一个从1到n的整数序列。这个序列用于表示x轴上的位置，其中n是回归模型中的观测值数量(样本量，从n <- length(fitted(fit))获得)。

通过使用1:n作为x轴的位置，identify()函数将在图形上绘制由帽子值和位置确定的数据点，并允许用户交互式地选择这些数据点来显示其标签。

这里identify()函数运行需要到脚本里进行，在右下模块中Plots的右上角有一个Finish的按钮，先选中需要展示标签的圆圈或点后再点Finish即可完成点的标签展示。

8.4.3 强影响点

即对模型参数估计值的影响有些比例失衡的点，例如若移除模型的一个观测点时模型会发生巨大的改变，那么就需要检测一下数据中是否存在强影响点。

8.4.3.1 检验方法1-Cook距离(库克距离)或称D统计量

通常，Cooks’D值大于4/(n-k-1),则表明它是强影响点，n:样本量，k:自变量数目。

语法+实例展示：

cutoff <- 4/(nrow(states)-length(coefficients(fit))-1)#计算Cooks’D值
plot(fit,which = 4,cook.levels = cutoff)#绘制回归诊断图，注释见下方
abline(h=cutoff,lty=2,col="red")

注释：1. which参数表示要绘制的回归诊断图的类型。它可以取不同的值来绘制不同的图形，其中包括：

which = 1：绘制标准化残差（Standardized Residuals）图

which = 2：绘制QQ图（Quantile-Quantile Plot）

which = 3：绘制斜方差图（Scale-Location Plot）

which = 4：绘制Cook距离（Cook’s Distance）图

which = 5：绘制杠杆值（Leverage）图

which = 6：绘制标准化残差与杠杆值的组合图

在上述代码中，which = 4表示要绘制Cook距离图。

2.cook.levels参数表示Cook距离的阈值，用于判断哪些观测值具有较大的影响。

在上述代码中，cutoff是通过公式计算得到的阈值，它的计算方式是 4 / (nrow(states) - length(coefficients(fit)) - 1)。cook.levels = cutoff将阈值传递给plot()函数，用于在Cook距离图中绘制阈值线。

8.4.3.2 检验方法2-变量添加图(added variable plot)

Cooks’D图有助于识别强影响点，但并不提供关于这些点如何影响模型的信息。变量添加图可以做到。

变量添加图：一个因变量和k个自变量，对于每个自变量Xk，绘制Xk在其他k-1个自变量上回归的残差值相对于因变量在其他k-1个自变量上回归的残差值的关系图。

语法(car中的avPlots()函数）：

library(car)
avPlots(fit)#显示出全部的子图
avPlots(fit,ask=FALSE,id=list(method="identify"))#交互式显示子图，两行代码选一种即可

备注：

1.ask参数是一个逻辑值，用于确定是否在每个子图绘制后等待用户输入。如果将ask设置为FALSE，则函数会自动绘制所有的子图而无需用户交互。默认情况下，ask的取值为TRUE，即绘制每个子图后会等待用户输入，用户可以通过按下回车键或点击图形窗口的任意位置来继续绘制下一个子图。在交互模式下，这可以方便用户观察和分析每个子图。

2.id参数是一个列表，用于指定标识数据点的方法。通过method参数可以选择不同的标识方法。在上述代码中，id=list(method=“identify”)指定了使用”identify”方法来标识数据点。这意味着在每个子图中，用户可以通过鼠标点击或光标移动来标识具体的数据点，并显示其观测值标签。这有助于用户识别特定的观测值或异常值，并进一步分析其对回归模型的影响。

举例:

library(car)
avPlots(fit,ask=FALSE)

变量添加图就可以看出每个自变量与因变量之间的强影响点，并可以预测出如果删除某一强影响点后的模型变化趋势。

8.4.4 综合反应异常观测值-car包的influencePlot()函数

该函数可以将离群点、杠杆值和强影响点信息整合到一张图形中。

语法：

library(car)
influencePlot(fit)

举例：

library(car)
influencePlot(fit,main="Influence Plot",
              sub="Circle size is proportional to Cook's distance")

##                 StudRes        Hat       CookD
## Alaska        1.7536917 0.43247319 0.448050997
## California   -0.2761492 0.34087628 0.008052956
## Nevada        3.5429286 0.09508977 0.209915743
## Rhode Island -2.0001631 0.04562377 0.035858963

#main 主标题，sub 下标题，同样可以添加id=list(method="identify可以交互式显示

注：

纵坐标超过+2或小于-2的为离群点(Nevada和Rhode Island);

水平轴超过0.2或0.3的为高杠杆点(California和Alaska)；

圆圈大小与影响成比例，圆圈很大的点可能是强影响点(Nevada和Alaska)

8.5 改进措施

1）删除观测点；例如删除强影响点或离群点

2）变量变换；例如用Ya或log来变换，见P191

3）添加或删除变量；

4）使用其他回归方法。

8.6 补充lines()函数的用法:

lines()函数用于在一个已有的图形上添加线条。lines()函数通常与绘图函数（如plot()）一起使用，以在同一图形中添加多个线条，绘制曲线、折线、平滑线等。

lines(x,y=NULL,type="l")

type：线条的类型。默认值为”l”，表示绘制折线。其他常用选项包括”b”（绘制实心点和折线）、“o”（绘制实心点和折线，但不绘制线条的起始和结束点）、“s”（绘制平滑线）等。

8.7 补充car包中的scatterplotMatrix()函数

R语言中car包（Caroline的R包）提供的函数，用于创建散点图矩阵（scatterplot matrix）。散点图矩阵是一种可视化工具，用于同时展示多个变量之间的散点图。

语法：

scatterplotMatrix(data,diagonal = "histogram", smooth = FALSE)

该函数会生成散点图矩阵(需安装car包)，例如：

states <- as.data.frame(state.x77[,c("Murder","Population","Illiteracy","Income","Frost")])
cor(states)

##                Murder Population Illiteracy     Income      Frost
## Murder      1.0000000  0.3436428  0.7029752 -0.2300776 -0.5388834
## Population  0.3436428  1.0000000  0.1076224  0.2082276 -0.3321525
## Illiteracy  0.7029752  0.1076224  1.0000000 -0.4370752 -0.6719470
## Income     -0.2300776  0.2082276 -0.4370752  1.0000000  0.2262822
## Frost      -0.5388834 -0.3321525 -0.6719470  0.2262822  1.0000000

library(car)
scatterplotMatrix(states,smooth=FALSE,main="Scatter Plot Matrix")

8.8 补充par(mfrow=c())函数

用于设置绘图设备中的绘图布局。通过指定参数mfrow的值，可以将绘图设备分割成一个指定行数和列数的矩阵，并按照这个布局进行图形的排列。

语法：

par(mfrow=c(rows,cols))

例如，par(mfrow=c(2, 2))将绘图设备分割成一个2行2列的矩阵，共有4个绘图区域。par(mfrow=c(2, 2))函数被调用后，后续的4个plot()函数调用将会自动按照2行2列的布局排列在绘图设备中，每个图形占据一个绘图区域。

需要注意的是，在使用完分割布局后，为了恢复默认的绘图布局，可以使用par(mfrow=c(1, 1))函数将绘图设备恢复为单一的绘图区域。

8.9 补充identify()函数

用于在图形设备上标识数据点，并通过交互式地选择数据点来显示其标签或索引。它通常用于在散点图中标识和标记特定的数据点。

语法：

identify(x,y,labels,pos=1,cex=1,col="black",...)

x和y:是数据点的x和y坐标。

labels是一个可选参数，用于指定要显示的标签。如果未提供labels参数，identify()函数将默认显示数据点的索引。

pos是标签的位置参数，用于指定标签的位置。默认值为1，表示标签位于数据点的上方。

cex是标签的字体缩放参数，默认为1，表示正常大小。

col是标签的颜色，默认为黑色。

…表示其他图形参数，如颜色、点形状等，可以根据需要进行设置

完整教程请查看

R语言基础学习手册