|
建模过程：选定模型族- 拟合模型 1.准备工作[AppleScript] 纯文本查看 复制代码 library(tidyverse)
## -- Attaching packages -------------------
## √ ggplot2 3.1.0 √ purrr 0.2.5
## √ tibble 1.4.2 √ dplyr 0.7.7
## √ tidyr 0.8.2 √ stringr 1.3.1
## √ readr 1.1.1 √ forcats 0.3.0
## -- Conflicts --- tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
[AppleScript] 纯文本查看 复制代码 library(modelr)
options(na.action = na.warn)
2.简单模型认识一下sim [AppleScript] 纯文本查看 复制代码 sim1
## # A tibble: 30 x 2
## x y
## <int> <dbl>
## 1 1 4.20
## 2 1 7.51
## 3 1 2.13
## 4 2 8.99
## 5 2 10.2
## 6 2 11.3
## 7 3 7.36
## 8 3 10.5
## 9 3 10.5
## 10 4 12.4
## # ... with 20 more rows
class(sim1)
## [1] "tbl_df" "tbl" "data.frame"
str(sim1)[mw_shl_code=applescript,true]ggplot(sim1, aes(x, y)) +
geom_point()
## Classes 'tbl_df', 'tbl' and 'data.frame': 30 obs. of 2 variables:
## $ x: int 1 1 1 2 2 2 3 3 3 4 ...
## $ y: num 4.2 7.51 2.13 8.99 10.24 ...[/mw_shl_code] 散点图
确定模型的基本形式 —线性—y= a_0 + a_1 * x [AppleScript] 纯文本查看 复制代码 models <- tibble(
a1 = runif(250, -20, 40),
a2 = runif(250, -5, 5)
)
ggplot(sim1, aes(x, y)) +
geom_abline(
aes(intercept = a1, slope = a2),
data = models, alpha = 1/4 )+
geom_point()计算出的y值：预测值 y实际值：响应变量 计算预测值与实际值之间的差异 [AppleScript] 纯文本查看 复制代码 model1 <- function(a, data) {
a[1] + data$x * a[2]
}
model1(c(7, 1.5), sim1)
## [1] 8.5 8.5 8.5 10.0 10.0 10.0 11.5 11.5 11.5 13.0 13.0 13.0 14.5 14.5
## [15] 14.5 16.0 16.0 16.0 17.5 17.5 17.5 19.0 19.0 19.0 20.5 20.5 20.5 22.0
## [29] 22.0 22.0这里的a有两个元素，分别是模型的两个参数。 给出了30个距离，实际上应该用均方根误差来表示，量化为一个数值。 [AppleScript] 纯文本查看 复制代码 measure_distance <- function(mod, data) {
diff <- data$y - model1(mod, data)
sqrt(mean(diff ^ 2))
}
measure_distance(c(7, 1.5), sim1)
## [1] 2.665212使用purrr计算所有模型的均方根误差 [AppleScript] 纯文本查看 复制代码 sim1_dist <- function(a1, a2) {
measure_distance(c(a1, a2), sim1)
}
models <- models %>%
mutate(dist = purrr::map2_dbl(a1, a2, sim1_dist))
models
## # A tibble: 250 x 3
## a1 a2 dist
## <dbl> <dbl> <dbl>
## 1 -19.1 0.176 34.1
## 2 -14.1 -3.37 50.6
## 3 6.48 3.08 8.73
## 4 9.80 2.08 6.12
## 5 5.49 -3.74 34.9
## 6 36.9 -4.82 20.5
## 7 34.0 1.90 29.0
## 8 22.3 1.65 16.1
## 9 2.46 -4.19 40.3
## 10 37.0 3.49 40.9
## # ... with 240 more rows将最好的（dist最小的）10个模型覆盖到数据上。 [AppleScript] 纯文本查看 复制代码 ggplot(sim1, aes(x, y)) +
geom_point(size = 2, color = "grey30") +
geom_abline(
aes(intercept = a1, slope = a2, color = -dist),
data = filter(models, rank(dist) <= 10)
)intercept是截距，slope是斜率。 a1和a2组成散点图，红色边框显示前10个最佳模型 [AppleScript] 纯文本查看 复制代码 ggplot(models, aes(a1, a2)) +
geom_point(
data = filter(models, rank(dist) <= 10),
size = 4, color = "red" )+
geom_point(aes(colour = -dist))更加系统化的搜索模型参数的方法： 网格搜索法 这里的“系统”是相对于之前生成随机数的操作说的。 生成网格，显示最佳10个模型 [AppleScript] 纯文本查看 复制代码 grid <- expand.grid(
a1 = seq(-5, 20, length = 25),
a2 = seq(1, 3, length = 25)
) %>%
mutate(dist = purrr::map2_dbl(a1, a2, sim1_dist))
grid %>%
ggplot(aes(a1, a2)) +
geom_point(
data = filter(grid, rank(dist) <= 10),
size = 4, colour = "red" )+
geom_point(aes(color = -dist))第三种方法：牛顿—拉夫逊搜索 [AppleScript] 纯文本查看 复制代码 best <- optim(c(0, 0), measure_distance, data = sim1)
best$par
## [1] 4.222248 2.051204
ggplot(sim1, aes(x, y)) +
geom_point(size = 2, color = "grey30") +
geom_abline(intercept = best$par[1], slope = best$par[2])线性模型拟合函数—lm() [AppleScript] 纯文本查看 复制代码 sim1_mod <- lm(y ~ x, data = sim1)
coef(sim1_mod)
## (Intercept) x
## 4.220822 2.051533 emmm 所以lm这么厉害，我们为啥要学之前那三种奇怪的办法。有一种被大神耍了的感觉，是不是只是为了强调一下lm很厉害， |
发表于 2018-11-26 14:24:39
