本文为笔者从官网学习的代码实录，不对的地方请多指教！
官网地址：https://www.tidymodels.org/
官网介绍：The tidymodels framework is a collection of packages for modeling and machine learning using tidyverse principles.
本次内容为get started部分的pipeline
内容目录如下图
内容目录
载入需要的包
library(tidymodels) # for the parsnip package, along with the rest of tidymodels
# Helper packages
library(readr) # for importing data
library(broom.mixed) # for converting bayesian models to tidy tibbles
library(dotwhisker) # for visualizing regression results
#读入测试数据
urchins <-read_csv("https://tidymodels.org/start/models/urchins.csv") %>% #读入数据
setNames(c("food_regime", "initial_volume", "width")) %>% #定义列名称
mutate(food_regime = factor(food_regime, levels = c("Initial", "Low", "High")))#设定因子变量
glimpse(urchins)
测试数据变量情况
food_regime列为三种不同的喂养策略，initial_volume为金枪鱼初始的体积，width为金枪鱼最终喂养后的宽度数据。研究的主要目的是看看不同的喂养策略对于金枪鱼最终宽度的影响。从尝试可以知道，金枪鱼的初始体积initial_volume也会影响最终width的结果。
#对数据进行可视化
ggplot(data = urchins,#数据集
aes(x = initial_volume, #全局映射
y = width,
group = food_regime,
color = food_regime)) +
geom_point() + #绘制点图
geom_smooth(method = lm, se = FALSE) +#绘制平滑曲线
scale_color_viridis_d(option = "plasma", end = .7) #色盲友好颜色
数据可视化情况，不同分组对于最终连续变量的影响
从图中可以看出，三组共同的趋势为：金枪鱼的初始体积越大，最终喂养后的宽度越大。不同的喂养策略产生的直线的斜率有点不同
可以看到，因为本研究的结局变量为数值型变量，所以应该用线性回归模型进行拟合分析
选定模型后，我们还需要对模型内部的engine进行选择，其定义如下：The engine value is often a mash-up of the software that can be used to fit or train the model as well as the estimation method. 个人认为engine的作用主要是确定损失函数。
linear_reg()#查看线性回归默认的engine
线性回归默认的engine
linear_reg()的engine可选项目
lm_mod <- linear_reg()#定义需要的模型，默认参数
lm_fit <- lm_mod %>% #对模型进行拟合
fit(width ~ initial_volume * food_regime, data = urchins)
tidy(lm_fit)#查看模型
模型拟合结果
#以下对于估计值及标准误进行可视化
tidy(lm_fit) %>%
dwplot(dot_args = list(size = 2, color = "black"),
whisker_args = list(color = "black"),
vline = geom_vline(xintercept = 0, colour = "grey50", linetype = 2))+
可视化结果
#构建测试数据
new_points <- expand.grid(initial_volume = 20,
food_regime = c("Initial", "Low", "High"))
new_points
测试数据
#进行点数据预测
m%ean_pred <- predict(lm_fit, new_data = new_points)
mean_pred
均值估计结果
#95%CI估计
conf_int_pred <- predict(lm_fit,
new_data = new_points,
type = "conf_int")
conf_int_pred
95%CI估计
#构建可视化需要的数据集
plot_data <-
new_points %>%
bind_cols(mean_pred) %>%
bind_cols(conf_int_pred)
plot_data
可视化需要的数据集
#画图
ggplot(plot_data, aes(x = food_regime)) +
geom_point(aes(y = .pred)) +
geom_errorbar(aes(ymin = .pred_lower,
ymax = .pred_upper),
width = .2) +
labs(y = "urchin size")
可视化预测结果
利用其他engine进行数据拟合及分析
#以下利用贝叶斯模型进行数据分析拟合
# 设定数据的先验分布，这是后面贝叶斯engine的参数
prior_dist <- rstanarm::student_t(df = 1)
#设定种子数
set.seed(123)
# 定义模型
bayes_mod <-
linear_reg() %>%
set_engine("stan",
prior_intercept = prior_dist,
prior = prior_dist)
# 训练模型
bayes_fit <-
bayes_mod %>%
fit(width ~ initial_volume * food_regime, data = urchins)
print(bayes_fit, digits = 5)#展示模型
贝叶斯模型拟合结果
tidy(bayes_fit, conf.int = TRUE)
整洁展示
#贝叶斯模型进行可视化
bayes_plot_data <-
new_points %>%
bind_cols(predict(bayes_fit, new_data = new_points)) %>%
bind_cols(predict(bayes_fit, new_data = new_points, type = "conf_int"))
ggplot(bayes_plot_data, aes(x = food_regime)) +
geom_point(aes(y = .pred)) +
geom_errorbar(aes(ymin = .pred_lower, ymax = .pred_upper), width = .2) +
labs(y = "urchin size") +
ggtitle("Bayesian model with t(1) prior distribution")
贝叶斯模型结果
以下网址为tidymodels包提供的可以拟合的模型和engine
https://www.tidymodels.org/find/parsnip/
