R-面向对象编程（一）

内置类和方法

S3泛型函数和方法在统一各个模型的使用方式上是最有用的。比如我们可以创建一个线性模型，以不同角度查看模型信息：

lm1 = lm(mpg ~ cyl + vs, data = mtcars)

线性模型本质上是由模型拟合产生的数据字段构成的列表，所以 lm1 的类型是 list ，但是它的类是 lm ，因此泛型函数根据 lm 选择方法：

typeof(lm1)
#> [1] "list"

class(lm1)
#> [1] "lm"

甚至没有明确调用S3泛型函数时，S3方法分派也会自动进行。如果我们输入 lm1 ：

lm1
#> 
#> Call:
#> lm(formula = mpg ~ cyl + vs, data = mtcars)
#> 
#> Coefficients:
#> (Intercept)          cyl           vs  
#>      39.625       -3.091       -0.939

实际上， print() 函数被默默地调用了：

print(lm1)
#> 
#> Call:
#> lm(formula = mpg ~ cyl + vs, data = mtcars)
#> 
#> Coefficients:
#> (Intercept)          cyl           vs  
#>      39.625       -3.091       -0.939

为什么打印出来的不像列表呢？因为 print() 是一个泛型函数，它为 lm 选择了一个方法来打印线性模型最重要的信息。我们可以调用 getS3method("print", "lm") 获取实际使用的方法与想象的进行验证：

identical(getS3method("print", "lm"), stats:::print.lm)
#> [1] TRUE

print() 展示模型的一个简要版本， summary() 展示更详细的信息。 summary() 也是一个泛型函数，它为模型的所有类提供了许多方法：

summary(lm1)
#> 
#> Call:
#> lm(formula = mpg ~ cyl + vs, data = mtcars)
#> 
#> Residuals:
#>    Min     1Q Median     3Q    Max 
#> -4.923 -1.953 -0.081  1.319  7.577 
#> 
#> Coefficients:
#>             Estimate Std. Error t value Pr(>|t|)    
#> (Intercept)   39.625      4.225    9.38  2.8e-10 ***
#> cyl           -3.091      0.558   -5.54  5.7e-06 ***
#> vs            -0.939      1.978   -0.47     0.64    
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> Residual standard error: 3.25 on 29 degrees of freedom
#> Multiple R-squared:  0.728,  Adjusted R-squared:  0.71 
#> F-statistic: 38.9 on 2 and 29 DF,  p-value: 6.23e-09

实际上， summary() 的输出结果也是一个对象，包含的数据都可以被访问。在这个例子里，这个对象是一个列表，是 summary.lm 类，它有可供 print() 选择的自己的方法：

lm1summary = summary(lm1)
typeof(lm1summary)
#> [1] "list"

class(lm1summary)
#> [1] "summary.lm"

查看列表成分：

names(lm1summary)
#>  [1] "call"          "terms"         "residuals"     "coefficients" 
#>  [5] "aliased"       "sigma"         "df"            "r.squared"    
#>  [9] "adj.r.squared" "fstatistic"    "cov.unscaled"

还有一些其他有用的且与模型相关的泛型函数，例如 plot() , predict() 。不同的内置模型和第三方扩展包提供的模型都能实现这些泛型函数。

举例，我们可以对线性模型调用 plot() 函数：

oldpar = par(mfrow = c(2, 2))
plot(lm1)

par(oldpar)

为避免依次生成这4个图，我们用 par() 将绘图区域划分为2x2的子区域。

利用 predict() 我们可以使用模型对新数据进行预测，泛型函数 predict() 自动选择正确的方法用新数据进行预测：

predict(lm1, data.frame(cyl = c(6, 8), vs = c(1, 1)))
#>    1    2 
#> 20.1 14.0

这个函数既可以用在样本内，又可以用在样本外。如果我们为模型提供新数据，它就进行样本外预测。

下面我们创建一幅真实值和拟合值的散点图，看一看线性模型的预测效果：

plot(mtcars$mpg, fitted(lm1))

这里的 fitted() 也是泛型函数，等价于 lm1$fitted.values ，拟合值等于用原始数据得到的预测值，即用原始数据构建的模型预测原始数据， predict(lm1, mtcars) 。

真实值与拟合值的差称为残差，可以通过另一个泛型函数 residuals() 获得。

plot(density(residuals(lm1)),
     main = "Density of lm1 residuals")

这些泛型函数不仅适用于 lm 、 glm 和其他内置模型，也适用于其他扩展包提供的模型。

例如我们使用 rpart 包，使用前面的数据和公式拟合一个回归树模型。

if(!require("rpart")) install.packages("rpart")
#> 载入需要的程辑包：rpart
library(rpart)

tree_model = rpart(mpg ~ cyl + vs, data = mtcars)

我们之所以能够使用相同的方法， 是因为这个包的作者希望函数调用的方式与调用R内置函数保持一致 。

typeof(tree_model)
#> [1] "list"
class(tree_model)
#> [1] "rpart"

打印模型：

print(tree_model)
#> n= 32 
#> 
#> node), split, n, deviance, yval
#>       * denotes terminal node
#> 
#> 1) root 32 1130.0 20.1  
#>   2) cyl>=5 21  198.0 16.6  
#>     4) cyl>=7 14   85.2 15.1 *
#>     5) cyl< 7 7   12.7 19.7 *
#>   3) cyl< 5 11  203.0 26.7 *

更详细信息：

summary(tree_model)
#> Call:
#> rpart(formula = mpg ~ cyl + vs, data = mtcars)
#>   n= 32 
#> 
#>       CP nsplit rel error xerror   xstd
#> 1 0.6431      0     1.000  1.089 0.2579
#> 2 0.0893      1     0.357  0.432 0.0811
#> 3 0.0100      2     0.268  0.427 0.0818
#> 
#> Variable importance
#> cyl  vs 
#>  65  35 
#> 
#> Node number 1: 32 observations,    complexity param=0.643
#>   mean=20.1, MSE=35.2 
#>   left son=2 (21 obs) right son=3 (11 obs)
#>   Primary splits:
#>       cyl < 5   to the right, improve=0.643, (0 missing)
#>       vs  < 0.5 to the left,  improve=0.441, (0 missing)
#>   Surrogate splits:
#>       vs < 0.5 to the left,  agree=0.844, adj=0.545, (0 split)
#> 
#> Node number 2: 21 observations,    complexity param=0.0893
#>   mean=16.6, MSE=9.45 
#>   left son=4 (14 obs) right son=5 (7 obs)
#>   Primary splits:
#>       cyl < 7   to the right, improve=0.507, (0 missing)
#>   Surrogate splits:
#>       vs < 0.5 to the left,  agree=0.857, adj=0.571, (0 split)
#> 
#> Node number 3: 11 observations
#>   mean=26.7, MSE=18.5 
#> 
#> Node number 4: 14 observations
#>   mean=15.1, MSE=6.09 
#> 
#> Node number 5: 7 observations
#>   mean=19.7, MSE=1.81

下面对结果进行可视化，得到树图：

oldpar = par(xpd = NA)
plot(tree_model)
text(tree_model, use.n = TRUE)

par(oldpar)