个性化阅读
专注于IT技术分析

R教程:R中的公式

本文概述

当你想到它时, R中的许多函数都使用公式:ggplot2, stats, lattice和dplyr等软件包都使用它们!使用这些R对象的函数的常见示例是glm(), lm(), facet_wrap()等。但是这些公式到底是什么, 为什么要使用它们呢?

这些只是本教程希望回答的一些问题:

  • R中的数据结构
  • R中的公式是什么?
  • 为什么在R中使用公式?
  • 在R中使用公式
    • 如何在R中创建公式
    • 如何连接公式
    • 公式运算符
  • 何时使用公式
    • R建模功能
    • R中的图形功能
    • dplyr中的非标准评估
  • R Formula套装

提示:你是否有兴趣在统计建模的背景下进一步了解公式?看看srcmini的多元和逻辑回归课程。

R中的数据结构

由于公式是R编程语言中的特殊类, 因此, 简短地修改此编程语言中可用的数据类型和数据结构是一个好主意。

请记住, R是一种面向对象的编程语言:该语言围绕对象进行组织。 R中的所有内容都是一个对象。

让我们从头开始:在编程中, 你将使用数据结构来存储数据, 并通过函数对其进行处理。数据结构是计算机内存中组织的数据的接口。正如R语言定义所指出的那样, R不提供对计算机内存的直接访问, 而是提供了许多可以称为”对象”的专用数据结构。每个数据结构都旨在优化存储, 访问或处理的某些方面。

R中的五个主要数据结构是:

  • 原子向量
  • 清单,
  • 矩阵,
  • 数据框, 以及
  • Array
# Create variables
a <- c(1, 2, 3, 4, 5, 6, 7, 8, 9)

b <- list(x = LifeCycleSavings[, 1], y = LifeCycleSavings[, 2])

提示:可以使用typeof()函数返回R对象的类型。对象的类型告诉你有关任何对象的(R内部)类型或存储模式的更多信息:

# Retrieve the types of `a` and `b`
typeof(a)

typeof(b)

‘双’

“清单”

在上面定义变量a和b的示例中, 你可以看到数据结构包含数据元素序列。这些元素可以是相同或不同的数据类型。你可以在R中找到以下6种原子数据类型:

  • 数字(例如100、5、4)包括整数。
  • 字符, 例如” Hello”, ” True”或” 23.4″, 由键盘字符字符串组成。
  • 逻辑, 例如TRUE或FALSE, 由”真值”组成;
  • 原始数据, 例如48 65 6c 6c 6f, 由位组成;
  • 复数(例如2 + 5i)包含复数;最后,
  • double(例如3.14)包含小数。

在R中, 几乎所有对象都具有附加的属性。例如, 你可能已经知道矩阵和数组就是简单的向量, 属性具有dim和附加在向量上的dimnames(可选)。属性用于实现R中使用的类结构。作为一种面向对象的编程语言, 类的概念以及方法是其核心。类是对象的定义。它定义了对象包含哪些信息以及如何使用该对象。

查看以下示例:

# Retrieve the classes of `a` and `b`
class(a)

class(b)

‘数字’

“清单”

请注意, 如果对象不具有class属性, 则它具有隐式类” matrix”, ” array”或mode()函数的结果。

日期和公式是你可能会遇到的一些特殊类。这是本教程的主题!

R中的公式是什么?

在阅读本教程的简介时, 你可能已经看到使用ggplot2之类的软件包或lm()之类的函数时会出现公式。由于通常在这些函数调用中使用公式来表达统计模型的概念, 因此在建模函数以及某些图形函数中经常使用这些R对象是合乎逻辑的。

对?

但是, 公式不仅限于模型。它们是功能强大的通用工具, 可让你捕获两件事:

  • 一个未经评估的表达式, 以及
  • 在其中创建表达式的上下文或环境。

这解释了为什么在函数调用中使用公式来生成”特殊行为”:它们使你可以捕获变量的值而不对其求值, 以便函数可以解释它们。

牢记新的数据结构, 然后可以将这些R对象描述为”语言”对象或未经评估的表达式, 这些表达式具有”公式”类和存储环境的属性。

在上一节中, 你看到对象具有某些(R内部)类型, 这些类型指示如何存储对象。在这种情况下, 公式是”语言”类型的对象。

但这到底是什么意思?

好吧, 在处理R语言本身时, 通常会遇到这种类型的对象。请看以下示例, 以更好地理解这一点:

# Retrieve the object type
typeof(quote(x * 10))

# Retrieve the class
class(quote(x * 10))

‘语言’

‘呼叫’

在上面的示例中, 你要求R返回quote(x * 10)的类型和类。结果, 你看到quote(x * 10)的类型是’language’, 而类是’call’。

这绝对不是公式, 因为你将需要class()返回”公式”!

但是那是什么呢?

代表R中公式的特征是波浪号〜。使用此运算符, 你实际上可以说:”捕获该代码的含义, 而无需立即对其进行评估。”这也解释了为什么你可以将R中的公式视为”引用”运算符。

但是公式到底是什么样子?仔细看看下面的代码行:

# A formula
c <- y ~ x
d <- y ~ x + b

# Double check the class of `c`
class(c)

‘式’

代字号(〜)左侧的变量称为”因变量”, 而右侧的变量称为”独立变量”, 并用加号+联接。

很高兴知道这些变量的名称会根据上下文而变化。你可能已经看到自变量显示为”预测变量(变量)”, “受控变量”, “功能”等。类似地, 你可能会遇到因变量为”响应变量”, “结果变量”或”标签”。

请注意, 即使你在上面的代码块中定义的公式d包含几个变量, 公式的基本结构实际上只是波浪符号〜和至少一个自变量或右侧变量。

请记住, 公式实际上是具有存储环境属性的语言对象:

# Return the type of `d`
typeof(d)

# Retrieve the attributes of `d`
attributes(d)

‘语言’

$class
[1] "formula"

$.Environment
<environment: R_GlobalEnv>

正如你在上面的示例中看到的那样, 公式中包含的变量可以是例如向量。但是, 你经常会看到公式中包含的变量来自数据框, 如以下示例所示:

Sepal.Width ~ Petal.Width + log(Petal.Length) + Species

请注意, 在创建公式本身时, 不会访问已分配给公式中的符号的任何数据值。

既然你知道了公式的样子以及它们在R中的含义, 现在就可以说基础公式对象会有所不同, 具体取决于你使用的是单面还是双面公式。你可以通过查看左侧变量来识别前者。如果没有, 就像〜x中一样, 你就有一个单边公式。

这也意味着一侧公式的长度为2, 而两侧公式的长度为3。

不完全相信吗?看一下下面的代码块。你可以在方括号[[和]]的帮助下访问公式的元素。

e <- ~ x + y + z
f <- y ~ x + b 

# Return the length of `g`
length(e)
length(f)

# Retrieve the elements at index 1 and 2
e[[1]]
e[[2]]
f[[3]]

2

3

`~`



x + y + z



x + b

为什么在R中使用公式?

如你所见, 公式强大的通用工具可让你捕获变量的值而无需求值, 以便函数可以解释它们。那已经是为什么要在R中使用公式的答案的一部分。

同样, 你可以使用这些R对象来表示变量之间的关系。

例如, 在下面的代码块的第一行代码中, 你在第一行代码中说” y是x, a和b的函数”。当然, 你也可以遇到更复杂的公式, 例如在第二行代码中, 你的意思是说”萼片宽度是花瓣宽度的函数, 取决于物种”。

y ~ x + a + b

Sepal.Width ~ Petal.Width | Species
y ~ x + a + b



Sepal.Width ~ Petal.Width | Species

在R中使用公式

既然你对这些特殊的R对象的”什么”和”为什么”有了更多的了解, 是时候学习如何使用基本以及更复杂的公式了!在本节中, 你不仅会看到如何创建和连接基本公式, 还将发现如何在运算符的帮助下构建更复杂的公式。

如何在R中创建公式

你已经知道该怎么做!你已经在本教程中看到了一些示例, 但让我们来概括一下:

y ~ x
~ x + y + z
g <- y ~ x + b

没错-你只需输入公式即可!

但是, 你可能会遇到需要或想要从R对象(例如字符串)创建公式的情况。在这种情况下, 可以使用公式或as.formula()函数:

"y ~ x1 + x2"

h <- as.formula("y ~ x1 + x2")

h <- formula("y ~ x1 + x2")

简单!

如何连接公式

要将多个公式粘合或组合在一起, 你有两个选择。首先, 你可以为每个公式创建单独的变量, 然后使用list():

# Create variables
i <- y ~ x
j <- y ~ x + x1
k <- y ~ x + x1 + x2

# Concatentate
formulae <- list(as.formula(i), as.formula(j), as.formula(k))

# Double check the class of the list elements
class(formulae[[1]])

‘式’

另外, 你也可以使用lapply()函数, 在其中传入一个向量, 其中所有公式都作为第一个参数, 而as.formula作为要应用于该向量的每个元素的函数:

# Join all with `c()`
l <- c(i, j, k)

# Apply `as.formula` to all elements of `f`
lapply(l, as.formula)
[[1]]
y ~ x

[[2]]
y ~ x + x1

[[3]]
y ~ x + x1 + x2

公式运算符

有了这些基本知识, 你就可以深入研究一些更复杂的公式了!在上面, 你已经看到, 代表公式的特征是波浪号〜。除了该符号之外, 你还看到还需要依赖和独立变量, 并且这些变量可以与加号+联接。

但是还有更多!

除+符号外, 还有其他可以为你的公式添加特殊含义的符号:

  • -删除条款;
  • :用于互动;
  • *用于穿越;
  • %in%用于嵌套;和
  • ^用于将极限穿越到指定程度。

在本节的其余部分中, 你将看到所有这些运算符的示例!首先让我们从+和-运算符开始:

# Use multiple independent variables
y ~ x1 + x2

# Ignore objects in an analysis
y ~ x1 - x2

请注意, 在回归建模上下文中, 你通常需要使用:和*符号, 你需要在其中指定交互项。如前所述, 前一个符号用于交互, 这意味着你只需要变量的交互, 而不是变量本身。这与用于交叉的后一个符号形成对比:你使用它来包括两个变量及其相互作用。

小心!通过使用这些运算符, 某些公式可能看起来不同, 但实际上是等效的。考虑以下示例, 它们将产生相同的回归:

y ~ x1 * x2

y ~ x1 + x2 + x1:x2

不确定这两个如何相同?考虑以下R代码块:

# Set seed
set.seed(123)

# Data
x = rnorm(5)
x2 = rnorm(5)
y = rnorm(5)

# Model frame
model.frame(y ~ x * x2, data = data.frame(x = x, y = y, x2=x2))
X x2
1.7150650 -0.56047565 1.2240818
0.4609162 -0.23017749 0.3598138
-1.2650612 1.55870831 0.4007715
-0.6868529 0.07050839 0.1106827
-0.4456620 0.12928774 -0.5558411
model.frame(y ~ x + x2 + x:x2, data = data.frame(x = x, y = y, x2))
X x2
1.7150650 -0.56047565 1.2240818
0.4609162 -0.23017749 0.3598138
-1.2650612 1.55870831 0.4007715
-0.6868529 0.07050839 0.1106827
-0.4456620 0.12928774 -0.5558411

如果你还不知道model.frame()函数, 请不要担心。你将在本教程的后面部分中看到更多有关此的信息!

此外, 这是嵌套的示例, 你可以将其扩展为y〜a + a:b:

y ~ a + b %in% a

所有这些运算符都很酷, 但是如果你想实际执行算术运算呢?假设你要在模型中包含x和x ^ 2。在这种情况下, 你可能会想编写以下公式:y〜x + x ^ 2。

但是结果是否正是你想要的?看一看!

model.frame( y ~ x + x^2, data = data.frame(x = rnorm(5), y = rnorm(5)))
y ~ x + x^2
X
-0.2053091 1.18231565
-0.3030972 0.04779636
-0.7621604 0.86382418
-0.1377784 -1.18333097
-0.3813125 -1.25247842

那不是你期望的结果!

在上面的示例中, 你不保护算术表达式, 因此, R将丢弃x ^ 2项, 因为它被视为x的重复项。

为什么是这样?

好吧, x会给你x的主要作用, 而x ^ 2会给你x的主要作用和二阶相互作用。最后, 你将最终在模型框架中包括x, 因为x的主要作用已经包含在公式中的x项中, 并且没有任何东西可以与x交叉来获得x ^中的二阶相互作用。 2学期。

为避免这种情况, 你有两种解决方案:

  • 你可以预先计算和存储所有变量
  • 你使用I()或”按原样”运算符:y〜x + I(x ^ 2)

请看以下示例, 以了解将I()函数添加到代码中的后果:

model.frame( y ~ x + I(x^2), data = data.frame(x = rnorm(5), y = rnorm(5)))
y ~ x + I(x^2)
X 我(x ^ 2)
1.414090 -0.1996230 0.039849….
1.777646 -1.0675904 1.139749….
1.710137 -1.4071841 1.980167….
1.259111 -1.3747289 1.889879….
-1.490866 0.8323668 0.692834….

最后一行代码实际上告诉R在使用公式之前计算x ^ 2的值。还要注意, 你可以使用”按原样”运算符来缩放模型的变量。你只需要在I()中包装相关的变量名:

y ~ I(2 * x)

当你看到上面的示例时, 这似乎都非常抽象, 因此让我们介绍其他一些情况。例如, 采用多项式回归。在这种情况下, 你将需要使用I()函数。另一方面, 对于仅限于depth = 2交互作用的阶乘方差分析模型, 你不需要此函数, 因为你想扩展为一个包含a, b和c的主要效果以及它们的第二个效果的公式, 订单互动:

# Polynomial Regression
y ~ x + I(x^2) + I(x^3)

# Factorial ANOVA
y ~ (a*b*c)^2

最后, 还有一个其他功能在使用多个变量时会有所帮助, 那就是。操作员。在公式中使用此运算符时, 将引用矩阵中尚未包含在模型中的所有其他变量。例如, 当你要对矩阵或数据框运行回归并且不想键入所有变量时, 这很方便:

y ~ .

如何检查R中的公式

创建公式后, 你可能还需要检查它。在本节中, 将向你介绍一些工具, 这些工具可用于进一步探索你创建的特殊R对象!

请注意, 你已经在上一节中看到了一些检查公式的方法:你看到了诸如attribute(), typeof(), class()等函数。

terms()函数

要检查和比较不同的公式, 可以使用terms()函数:

m <- formula("y ~ x1 + x2")
terms(m)
y ~ x1 + x2
attr(, "variables")
list(y, x1, x2)
attr(, "factors")
   x1 x2
y   0  0
x1  1  0
x2  0  1
attr(, "term.labels")
[1] "x1" "x2"
attr(, "order")
[1] 1 1
attr(, "intercept")
[1] 1
attr(, "response")
[1] 1
attr(, ".Environment")
<environment: R_GlobalEnv>

all.vars

如果你想知道模型中变量的名称, 可以使用all.vars。使用此函数, 你将返回一个字符向量, 其中包含公式中出现的所有名称:

print(all.vars(m))
[1] "y"  "x1" "x2"

update()函数

要修改公式而不将其转换为字符, 可以使用update()函数:

update(y ~ x1 + x2, ~. + x3)
y ~ x1 + x2 + x3

请注意, 你还可以通过使用as.character()将其转换为字符来更新公式。然后, 你可以使用paste()快速构建公式。例如, 如果你想在公式中添加另一个右侧变量, 则只需将其粘贴即可:

as.formula(paste("y ~ x1 + x2", "x3", sep = "+"))

factors <- c("x2", "x3")
as.formula(paste("y~", paste(factors, collapse="+")))
y ~ x1 + x2 + x3



y ~ x2 + x3

在上面的代码块中, 你可以看到可以使用sep或折叠参数来指示一个字符串, 用于在公式中分隔各个术语。

但是, 使用paste()绝对不是调整公式的唯一方法:你还可以使用Reformulate():

reformulate(termlabels = factors, response = 'y')
y ~ x2 + x3

is.formula()

通过将变量传递给is.formula()函数, 仔细检查变量是否为公式。考虑到该函数是plyr库的一部分, 因此在调用is.formula()之前, 需要使该函数在工作区中可用。

# Load `plyr`
library(plyr)

# Check `m`
is.formula(m)

true

何时使用公式

到目前为止, 你已经了解到R公式是通用工具, 不仅限于建模, 还看到了一些可以使用公式的示例。在本节中, 你将更深入地讨论最后一个主题:你将看到在某些情况下可以使用这些工具以发挥自己的优势。当然, 你将介绍数据包的建模和图形功能, 例如网格和统计数据, 但还将介绍dplyr中的非标准评估。

建模功能

R非常适合你需要进行统计建模时使用。如你所知, 统计建模是一种简化的, 数学形式化的方式, 可以逼近现实并可以选择根据这种逼近进行预测。统计模型通常以理想化的形式表示数据生成过程。要进行统计建模, 你需要建模功能。

R中的建模函数是一个典型的示例, 你需要一个公式对象作为参数。在这些函数中可能找到的其他参数是数据, 它允许你指定要在模型持续时间内附加的数据框, 选择要使用的数据的子集, …通常, 如果你想知道要传递给任何特定函数的参数, 不要犹豫使用help()函数还是?在你的R控制台中。

建模函数返回一个模型对象, 其中包含有关拟合的所有信息。诸如print(), summary(), plot(), anova()等的通用R函数将具有为特定对象类定义的方法, 以返回适用于该类型对象的信息。

lm()可能是众所周知的建模函数之一, 它使用了上述所有参数。你可以使用lm()拟合线性模型。你可以使用它执行回归, 方差单层分析和协方差分析。让我们看一个使用lm()并借助print()检查模型的示例:

lm.m <- lm(Sepal.Width ~ Petal.Width + log(Petal.Length) + Species, data = iris, subset = Sepal.Length > 4.6)

print(lm.m)
Call:
lm(formula = Sepal.Width ~ Petal.Width + log(Petal.Length) + 
    Species, data = iris, subset = Sepal.Length > 4.6)

Coefficients:
      (Intercept)        Petal.Width  log(Petal.Length)  Speciesversicolor  
           3.1531             0.6620             0.4612            -1.9265  
 Speciesvirginica  
          -2.3088  

lm()最初使用公式和适当的环境来转换变量之间的关系, 以创建包含数据的数据框。

还有model.frame()方法, 在本教程中你已经看到了一个示例, 该方法最常用于从拟合对象中检索或重新创建模型框架, 而没有其他参数。这使你可以检索数据框的与公式, 子集和权重之外的原始调用参数相对应的列:例如, glm()方法处理offset, etastart和mustart。

在下面的代码块中, 你看到你使用model.frame()函数取回了拟合对象的数据框。请注意, 该代码与你在上面看到的代码块稍有不同:subset参数已被略微修改。

stats::model.frame(formula = Sepal.Width ~ Petal.Width + log(Petal.Length) + Species, data = iris, subset = Sepal.Length > 6.9, drop.unused.levels = TRUE)
分隔宽度 花瓣宽度 log(Petal.Length) 种类
51 3.2 1.4 1.547563 杂色
103 3.0 2.1 1.774952 维吉尼亚
106 3.0 2.1 1.887070 维吉尼亚
108 2.9 1.8 1.840550 维吉尼亚
110 3.6 2.5 1.808289 维吉尼亚
118 3.8 2.2 1.902108 维吉尼亚
119 2.6 2.3 1.931521 维吉尼亚
123 2.8 2.0 1.902108 维吉尼亚
126 3.2 1.8 1.791759 维吉尼亚
130 3.0 1.6 1.757858 维吉尼亚
131 2.8 1.9 1.808289 维吉尼亚
132 3.8 2.0 1.856298 维吉尼亚
136 3.0 2.3 1.808289 维吉尼亚

提示:stats包中还有更多函数可让你使用公式, 例如aggregate()。

对于线性混合效果模型, 该模型允许你对随机效果进行建模以解决由观察者差异等因素引起的变化, 可以将nlme包与lme()函数一起使用。同样在这里, 你会看到公式是你需要提供给建模函数的第一个参数, 并且还有一个数据参数!

# Load packages
library(MASS)
library(nlme)

# Get some data 
data(oats)

# Adjust the data names and columns
names(oats) = c('block', 'variety', 'nitrogen', 'yield')
oats$mainplot = oats$variety
oats$subplot = oats$nitrogen

# Fit a non-linear mixed-effects model 
nlme.m = lme(yield ~ variety*nitrogen, random = ~ 1|block/mainplot, data = oats)

# Retrieve a summary
summary(nlme.m)
Linear mixed-effects model fit by REML
 Data: oats 
       AIC      BIC    logLik
  559.0285 590.4437 -264.5143

Random effects:
 Formula: ~1 | block
        (Intercept)
StdDev:    14.64496

 Formula: ~1 | mainplot %in% block
        (Intercept) Residual
StdDev:    10.29863 13.30727

Fixed effects: yield ~ variety * nitrogen 
                                    Value Std.Error DF   t-value p-value
(Intercept)                      80.00000  9.106958 45  8.784492  0.0000
varietyMarvellous                 6.66667  9.715028 10  0.686222  0.5082
varietyVictory                   -8.50000  9.715028 10 -0.874933  0.4021
nitrogen0.2cwt                   18.50000  7.682957 45  2.407927  0.0202
nitrogen0.4cwt                   34.66667  7.682957 45  4.512152  0.0000
nitrogen0.6cwt                   44.83333  7.682957 45  5.835427  0.0000
varietyMarvellous:nitrogen0.2cwt  3.33333 10.865342 45  0.306786  0.7604
varietyVictory:nitrogen0.2cwt    -0.33333 10.865342 45 -0.030679  0.9757
varietyMarvellous:nitrogen0.4cwt -4.16667 10.865342 45 -0.383482  0.7032
varietyVictory:nitrogen0.4cwt     4.66667 10.865342 45  0.429500  0.6696
varietyMarvellous:nitrogen0.6cwt -4.66667 10.865342 45 -0.429500  0.6696
varietyVictory:nitrogen0.6cwt     2.16667 10.865342 45  0.199411  0.8428
 Correlation: 
                                 (Intr) vrtyMr vrtyVc ntr0.2 ntr0.4 ntr0.6
varietyMarvellous                -0.533                                   
varietyVictory                   -0.533  0.500                            
nitrogen0.2cwt                   -0.422  0.395  0.395                     
nitrogen0.4cwt                   -0.422  0.395  0.395  0.500              
nitrogen0.6cwt                   -0.422  0.395  0.395  0.500  0.500       
varietyMarvellous:nitrogen0.2cwt  0.298 -0.559 -0.280 -0.707 -0.354 -0.354
varietyVictory:nitrogen0.2cwt     0.298 -0.280 -0.559 -0.707 -0.354 -0.354
varietyMarvellous:nitrogen0.4cwt  0.298 -0.559 -0.280 -0.354 -0.707 -0.354
varietyVictory:nitrogen0.4cwt     0.298 -0.280 -0.559 -0.354 -0.707 -0.354
varietyMarvellous:nitrogen0.6cwt  0.298 -0.559 -0.280 -0.354 -0.354 -0.707
varietyVictory:nitrogen0.6cwt     0.298 -0.280 -0.559 -0.354 -0.354 -0.707
                                 vM:0.2 vV:0.2 vM:0.4 vV:0.4 vM:0.6
varietyMarvellous                                                  
varietyVictory                                                     
nitrogen0.2cwt                                                     
nitrogen0.4cwt                                                     
nitrogen0.6cwt                                                     
varietyMarvellous:nitrogen0.2cwt                                   
varietyVictory:nitrogen0.2cwt     0.500                            
varietyMarvellous:nitrogen0.4cwt  0.500  0.250                     
varietyVictory:nitrogen0.4cwt     0.250  0.500  0.500              
varietyMarvellous:nitrogen0.6cwt  0.500  0.250  0.500  0.250       
varietyVictory:nitrogen0.6cwt     0.250  0.500  0.250  0.500  0.500

Standardized Within-Group Residuals:
        Min          Q1         Med          Q3         Max 
-1.81300898 -0.56144838  0.01758044  0.63864476  1.57034166 

Number of Observations: 72
Number of Groups: 
              block mainplot %in% block 
                  6                  18 

请注意, 除了nlme之外, 还有其他软件包, 例如lme4软件包, 该软件包也专用于拟合线性和广义线性混合效果模型。

建模函数和公式存在的另一个示例是nls(), 可用于创建非线性模型:

# Set seed
set.seed(20160227)

# Data
x <- seq(0, 50, 1)
y <- ((runif(1, 10, 20)*x)/(runif(1, 0, 10)+x))+rnorm(51, 0, 1)

# Non-linear model
nls.m <- nls(y ~ a*x/(b+x), start=c(a=4, b=1))

最后一个示例是可用于构建广义线性模型(GLM)的函数。在R中, 你可以使用glm()函数执行此操作。它可能有点老了, 但在这里你还可以使用公式和数据参数:

# Load package
library(MPDiR)

# Get the data
data(Chromatic)

# Model
glm.m <- glm(Thresh ~ Axis:(I(Age^-1) + Age), family = Gamma(link = "identity"), data = Chromatic)

# Get back a summary
summary(glm.m)
Call:
glm(formula = Thresh ~ Axis:(I(Age^-1) + Age), family = Gamma(link = "identity"), data = Chromatic)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-1.2160  -0.3728  -0.0805   0.2311   1.2932  

Coefficients:
                      Estimate Std. Error t value Pr(>|t|)    
(Intercept)          3.282e-04  9.965e-05   3.294  0.00106 ** 
AxisDeutan:I(Age^-1) 7.803e-03  3.686e-04  21.172  < 2e-16 ***
AxisProtan:I(Age^-1) 8.271e-03  3.863e-04  21.410  < 2e-16 ***
AxisTritan:I(Age^-1) 1.166e-02  5.284e-04  22.065  < 2e-16 ***
AxisDeutan:Age       1.521e-05  3.418e-06   4.450 1.06e-05 ***
AxisProtan:Age       1.540e-05  3.434e-06   4.484 9.10e-06 ***
AxisTritan:Age       4.812e-05  5.838e-06   8.241 1.48e-15 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for Gamma family taken to be 0.2054848)

    Null deviance: 543.35  on 510  degrees of freedom
Residual deviance: 100.40  on 504  degrees of freedom
AIC: -4777.6

Number of Fisher Scoring iterations: 6

在继续使用图形函数之前, 还需要了解一件事:在对lm()等函数建模时使用公式时, 会发生从公式到函数的标准转换。为了解释到底发生了什么类型的转换, 你应该返回在本教程前面看到的示例:

lm.m <- lm(Sepal.Width ~ Petal.Width + log(Petal.Length) + Species, data = iris, subset = Sepal.Length > 4.6)

尽管此代码块的目的是拟合线性回归模型, 但该公式用于指定符号模型以及生成预期的设计矩阵。设计矩阵是预测变量或自变量集的二维表示, 其中数据实例在行中, 变量属性在列中。设计矩阵也称为X矩阵。

话虽这么说, 公式方法还定义了设计矩阵中应包括的列。

但是, 这到底意味着什么?让我们看一下下面的代码行:

# A formula
y ~ x

# A converted formula
y = a_1 + a_2 * x

这是一个简单转换的示例:y〜x被转换为y = a_1 + a_2 * x。

要查看和了解R实际发生的情况, 可以使用model_matrix()函数。该函数通过根据对比将因子扩展到一组虚拟变量来创建设计或模型矩阵, 并类似地扩展交互。

不要忘记在数据框dfand中传递公式, 以返回定义模型方程式的标题:

# Load packages
library(tidyverse)
library(modelr)

# A data frame
df <- tribble(
  ~y, ~x1, ~x2, 4, 2, 5, 5, 1, 6
)

# Model matrix
model_matrix(df, y ~ x1)
(截距) x1
1 2
1 1

出现一个额外的(拦截)列。这是R端的默认行为。但是, R将截距添加到模型的方式仅是通过拥有一个充满1的列。如果你不希望这样做, 则需要通过在公式中添加-1来明确删除它, 如下所示:

model_matrix(df, y ~ x1-1)
x1
2
1

还要注意, 当你向模型添加更多变量时, 模型矩阵的增长就不足为奇了。

model_matrix(df, y ~ x1 + x2)
(截距) x1 x2
1 2 5
1 1 6

R中的图形功能

在R中可以找到公式的另一个重要位置是图形函数。那里有很多软件包, 因此本教程仅关注图形, 点阵, ggplot2和ggformula。

图形

基本的R包图形允许你指定散点图或使用公式添加点, 线或文本。看下面的例子:

# Get data
data(airquality)

# Plot
plot(Ozone ~ Wind, data = airquality, pch = as.character(Month))
公式R教程

如果你想了解更多信息, 请随时查看此页面。

格子

点阵是基于网格图形构建的程序包。这是一个通用图形包, 提供基本图形中可用的许多绘图功能的替代实现。具体示例包括带有xyplot()函数的散点图, 带有barchart()函数的条形图和带有bwplot()函数的箱形图。

提示:你想了解更多有关点阵的信息吗?考虑使用R进行带有网格过程的srcmini的数据可视化。

该软件包的特别之处在于, 它使用统计模型的公式符号来描述所需的图, 更具体地说, 是要绘制的变量。它还添加了竖线字符或竖线。指定条件变量:

# Load package
library(lattice)

# Plot histogram
histogram(~ Ozone | factor(Month), data = airquality, layout = c(2, 3), xlab = "Ozone (ppb)")
R公式教程

请注意, 就像建模函数一样, 晶格库中的函数除了公式自变量外, 还具有数据自变量, 就像人们期望的那样!

请记住, 你也可以在上面的代码块中省略data参数, 但随后需要附加数据:

# Load package
library(lattice)

# Attach data
attach(airquality)

# Plot
histogram(~ Ozone | factor(Month), layout = c(2, 3), xlab = "Ozone (ppb)")

ggplot2

你可以在各种ggplot2函数中使用公式:

  • geom_smooth()或stats_smooth(), 用于指定要在平滑函数中使用的公式;这将影响拟合的形式。
  • facet_wrap(), 用于指定绘图面板。
  • facet_grid(), 用于指定需要绘制的行和列(带或不带小平面)。
# Load package
library(ggplot2)

# Plot
ggplot(mpg, aes(displ, hwy)) +
  geom_point() +
  geom_smooth(method = "lm", formula = y ~ splines::bs(x, 3), se = FALSE)
R教程中的公式3

请注意, 要创建此绘图, 公式将使用字母x和y, 而不是变量的名称。使用facet_wrap()时并非如此:

ggplot(mpg, aes(displ, hwy)) +
  geom_point() +
  geom_smooth(span = 0.8) +
  facet_wrap(~drv)
R中的公式

想更多地了解ggplot2?看看这门课程或查阅文档。

ggformula

ggformula软件包当前基于ggplot2构建, 但是提供了基于公式的接口。这类似于你在晶格界面中找到的内容。你还将发现在此程序包中使用了管道运算符, 可从simpeler组件构建更复杂的图形。

提示:如果你想进一步了解R中的管道运算符, 请考虑本教程。

用ggformula创建图的基本方法是

gf_plottype(formula, data = mydata)

请注意, 函数gf_plottype()实际上以gf开头, 这提醒你使用的是基于公式的ggplot2接口。 g代表ggplot2, f代表”公式”。

请看以下示例, 以了解R代码的外观:

# Load package
library(ggformula)

# Plot
gf_point(mpg ~ hp, data = mtcars)
R教程中的公式5

当然, 这只是一个基本图形。你可以使用此可视化包做更多的事情!你可以为这些字形选择其他字形类型和特定属性, 可以绘制一个或两个可变图, 调整字形的位置, 等等。本教程不会对此包进行更详细的介绍, 但是这里的主要要点是, 该包已使公式成为制作图形的主要成分!

如果你想了解的内容不只是本教程中介绍的内容, 请在此处阅读有关ggformula软件包的信息, 或查阅该软件包的RDocumentation页面。

dplyr

dplyr是与非标准评估一起使用的软件包的示例。其他示例是library(magrittr)与library(” magrittr”), 即使一行中有引号, 而另一行中没有引号, 两者也都可以正常工作!比较这些library()与install.packages()的示例, 例如:

# This will work
install.packages("magrittr")

# This won't work
install.packages(magrittr)

在R中, 通常在命名对象的一部分时都需要使用引号, 但是在某些函数(例如library())中却不需要。这些功能被设计为以非标准方式工作。而且, 其中一些功能甚至错过了标准方式!

对于dplyr而言, 你会发现该软件包的大多数功能与其他功能一样工作:所有功能均使用标准评估(SE)。但是, 对于交互使用, 函数还具有非标准评估(NSE), 可以节省一些键入时间。

(让我们面对现实吧, 输入所有这些引号可能是一项艰巨的任务!)

这就是大多数dplyr函数使用非标准评估的原因。他们没有遵循通常的R评估规则。相反, 它们捕获你键入的表达式并以自定义方式对其进行求值。

但是, 这并不意味着这些功能没有标准的评估变体。使用非标准评估的每个函数都具有(并且应该具有)进行实际计算的标准评估转义哈希。标准评估函数应以_结尾。这意味着dplyr包中有多个动词:select(), select _(), mutate(), mutate_()等。

当以交互方式使用时, 将首先使用lazyeval软件包对这些函数进行评估, 然后再将其发送到该函数的标准评估版本。这意味着, 在高级情况下, 使用lazyeval包对select()进行了评估, 然后将其发送到select_()。

综上所述, dplyr和lazyeval可以通过3种方式在标准评估函数中引用变量:

  • 公式
  • quote()和
  • 弦乐
# Load `dplyr`
library(dplyr)

# NSE evaluation
select(iris, Sepal.Length, Petal.Length)

# standard evaluation 
select_(iris, ~Sepal.Length)
select_(iris, ~Sepal.Length, ~Petal.Length) #works
select_(iris, quote(Sepal.Length), quote(Petal.Length)) # yes!
select_(iris, "Sepal.Length", "Petal.Length", "Species")

提示:如果你想了解有关非标准评估的更多信息, 请考虑阅读Hadley Wickham的Advanced R书中有关该主题的章节。

R Formula套装

以前, 你已经看到可以使用诸如.formula, update(), all.vars等功能来创建和检查公式, 这些都是简单的操作和操纵, 但是如何对公式进行高级操纵呢?以下软件包可能会让你感兴趣!

配方套餐

最近, 该软件包在CRAN上发布。对于那些希望将公式推向新高度的人, 此软件包非常理想。该软件包扩展了基类公式。

更具体地说, “公式”对象扩展了基本的”公式”对象:使用此程序包, 你实际上可以定义接受其他公式运算符的公式。分隔多个部分, 或者可能将所有公式运算符(包括竖线字符)保留在左侧以支持多个响应。

你将能够创建的公式示例如下:

  • 多部分公式, 例如y〜x1 + x2 | u1 + u2 + u3 | v1 + v2
  • 多响应公式, 例如y1 + y2〜x1 + x2 + x3
  • 多部分响应, 例如y1 | y2 + y3〜x, 并且
  • 以上三个的组合。
# Load package
library(Formula)

# Create formulas
f1 <- y ~ x1 + x2 | z1 + z2 + z3
F1 <- Formula(f1)

# Retrieve the class of `F1`
class(F1)
  1. ‘式’
  2. ‘式’

请注意, 函数as.formula()和is.formula()在此程序包中也已更新:你将改用is.Formula()和as.Formula()!

在这里阅读更多。

formula.tools

该软件包是最近发布的, 为你提供”用于操纵公式, 表达式, 调用, 赋值和其他R对象的编程实用程序”。这很多, 但从本质上讲, 这全都可以归结为:你可以使用此程序包访问和修改公式结构, 以及提取和替换这些R对象的名称和符号。该软件包由克里斯托弗·布朗(Christopher Brown)撰写。

使用此软件包时, 你可能会发现有用的一些事项如下:

  • get.vars():而不是all.vars(), 此函数将从各种R对象中提取变量名, 但是所有符号等将内插到变量名中。
  • invert():你可以使用此函数反转对象(例如公式)中的运算符。
  • is.one.side():此函数非常便于确定函数是单面还是双面。

请记住, 公式看起来像是单面的:〜x;公式化为x〜y时将是双面的。

还有更多发现!

欢呼!你已经完成了有关R公式的本教程。如果你想了解更多有关它们的信息, 请一定要阅读Hadley Wickham的R for Data Science书籍, 该书的一章完全致力于R中的公式和模型族。

你能想到更多可以在其中找到公式的实例, 还是可以用于操纵公式的更多程序包?随时在Twitter上让我知道:@willems_karlijn。

赞(0)
未经允许不得转载:srcmini » R教程:R中的公式

评论 抢沙发

评论前必须登录!