建模过程： 定义模型族 → 从中找到最接近数据的一个模型
1 加载包
library(tidyverse)
library(modelr)
2 查看数据关系：ggplot查看图形
> 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
ggplot(sim1,aes(x,y))+
geom_point()
image.png
3 随机生成模型覆盖到数据上
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/10
) +
geom_point()
4 评价模型：量化数据与模型之间的距离
计算每个数据点与模型之间的 垂直距离
4.1 FUNCTION：根据模型计算预测值
model1 <- function(a,data){
a[1] + data$x * a[2]
}
4.2 FUNCTION：计算均方根误差
实际值-预测值→ 平方 → 均值 → 平方根
measure_distance <- function(a,data){
diff <- data$y - model1(a,data)
sqrt(mean(diff^2))
}
4.3 (purr包map批量计算点到模型的距离和）
models <- models %>%
mutate(dist=map2_dbl(a1,a2,sim1_dist))
4.4 结果可视化
4.4.1 可视化最优线性模型
ggplot(sim1,aes(x,y))+
geom_point(color="grey10")+
geom_abline(
aes(intercept=a1, slope=a2,color=-dist),
data= filter(models,rank(dist)<= 10))
4.4.2 散点图可视化模型参数
ggplot(models,aes(a1,a2))+
geom_point(data= filter(models,rank(dist)<= 10),
size=4,
color="red")+
geom_point(aes(color=-dist))
4.4.3 网格搜索法
grid %>%
ggplot(aes(a1,a2))+
geom_point(
data= filter(grid,rank(dist)<= 10),
size=4,
color="red"
)+
geom_point(aes(color= -dist))
#可视化网格搜索法所得最佳模型
ggplot(sim1,aes(x,y))+
geom_point(size=2,color="grey20")+
geom_abline(
aes(intercept = a1,slope = a2,color= -dist),
data = filter(grid,rank(dist)<= 10))
optim():牛顿-拉夫逊搜索
只要有定义模型与数据集间距离的函数及可以通过修改模型参数使距离最小化的算法就可以找出最小模型
best <- optim(c(0,0),measure_distance,data=sim1)
# 最优模型的参数
best$par
# lm()直接建模
sim1_mod <- lm(y~x,data=sim1)
5. 模型可视化
重点通过预测来理解模型，预测后找出模型未捕获的信息——残差，可以去除数据中显著的模式，以便对剩余的微妙趋势进行研究
5.1 预测
5.1.1预测——首先生成分布均匀的数值网格
modelr::data_grid(dataframe,x,y,...)
第一个参数是数据框，对于随后每个参数，找出其中的唯一值，生成所有组合
grid <- sim1 %>%
data_grid(x) %>%
grid
5.1.2 预测——向数值网格添加预测值
grid <- sim1 %>%
add_predictions(sim1_mod)
grid
5.1.3 预测——绘制预测值
ggplot(sim1,aes(x,y))+
geom_point(aes(y=y))+
geom_line(
aes(y=pred),
data = grid,
colour= "red",
size=1
)
5.2 残差
残差表示模型漏掉的部分，就是预测值与观测值之间的距离
- 生成残差值
# 添加预测是在data_grid生成的数值网格上添加；而残差是直接在原始数据集添加
sim1 <- sim1 %>%
add_residuals(sim1_mod)
sim1
- 可视化残差——绘制频率多边形图
ggplot(sim1,aes(resid))+
geom_freqpoly(binwidth=0.3)
image.png
- 可视化残差——用残差代替原来的变量绘制散点图
> ggplot(sim1,aes(x,resid))+
+ geom_ref_line(h=0)+
+ geom_point()
image.png
6 公式和模型族
6.1 分类变量
> ggplot(sim2,aes(x,y))+
+ geom_point(aes(y=y))+
+ geom_point(
+ data= grid,
+ aes(y=pred),
+ colour= "red",
+ size=4
+ )
image.png
带有分类变量x恶模型会为每个分类预测出均值（均值会使均方根距离最小化）
不能对为观测到的水平i女性预测
6.2.1 交互项（连续变量与分类变量）
ggplot(sim3,aes(x1,y))+
geom_point(aes(color=x2))
image.png
# 添加预测值: 两个变量都要传递到data_grid()，gater_predictions()/spread(predictions()函数
mod1 <- lm(y~x1+x2,data = sim3)
mod2 <- lm(y~x1*x2,data = sim3)
grid <- sim3 %>%
data_grid(x1,x2) %>%
gather_predictions(mod1,mod2)
grid
# 可视化：分层facet_wrap()函数
ggplot(sim3,aes(x1,y,color=x2))+
geom_point()+
geom_line(data = grid,aes(y=pred))+
facet_wrap(~model)
image.png
6.2.2 交互项（两个连续变量）
# A tibble: 300 x 4
x1 x2 rep y
<dbl> <dbl> <int> <dbl>
1 -1 -1 1 4.25
2 -1 -1 2 1.21
3 -1 -1 3 0.353
4 -1 -0.778 1 -0.0467
5 -1 -0.778 2 4.64
6 -1 -0.778 3 1.38
7 -1 -0.556 1 0.975
8 -1 -0.556 2 2.50
9 -1 -0.556 3 2.70
10 -1 -0.333 1 0.558
# ... with 290 more rows
> mod1 <- lm(y~x1+x2,data = sim4)
> mod2 <- lm(y~x1*x2,data = sim4)
> mod2 <- lm(y~x1*x2,data = sim4)
> grid <- sim4 %>%
+ data_grid(
+ x1=seq_range(x1,5),
+ x2=seq_range(x2,5)) %>%
+ gather_predictions(mod1,mod2)
> ggplot(data=grid,aes(x1,x2))+
+ geom_tile(aes(fill=pred))+
+ facet_wrap(~model)
image.png
此时肉眼不能辨别两个模型之间的区别，此时应当分别从x1 x2的角度来看，并表示出多个切面
> ggplot(data=grid,aes(x1,pred,color=x2,group=x2))+
+ geom_line()+
+ facet_wrap(~model)
image.png
image.png
> ggplot(data=grid,aes(x2,pred,color=x1,group=x1))+
+ geom_line()+
+ facet_wrap(~model)
image.png
网友评论