全文下载链接：http://tecdat.cn/?p=9024

最近咱们被客户要求撰写对于狭义相加模型（GAM）的钻研报告，包含一些图形和统计输入。

这篇文章探讨了为什么应用狭义相加模型 是一个不错的抉择。为此，咱们首先须要看一下线性回归，看看为什么在某些状况下它可能不是最佳抉择。

回归模型

假如咱们有一些带有两个属性Y和X的数据。如果它们是线性相关的，则它们可能看起来像这样：

为了查看这种关系，咱们能够应用回归模型。线性回归是一种应用X来预测变量Y的办法。将其利用于咱们的数据将预测成红线的一组值：

这就是“直线方程式”。依据此等式，咱们能够从直线在y轴上开始的地位（“截距”或）开始形容，并且每个单位的x都减少了多少y（“斜率”），咱们将它称为x的系数，或称为）。还有一点天然的稳定，如果没有的话，所有的点都将是完满的。咱们将此称为“残差”（）。

数学上是：

或者，如果咱们用理论数字代替，则会失去以下后果：

这篇文章通过思考每个数据点和线之间的差别（“残差）而后最小化这种差别来估算模型。

咱们在线的上方和下方都有正误差和负误差，因而，通过对它们进行平方并最小化“平方和”，使它们对于预计都为正。这称为“一般最小二乘法”或OLS。

3非线性关系如何？

因而，如果咱们的数据看起来像这样，咱们该怎么办：

咱们刚刚看到的模型的要害假如之一是y和x线性相关。如果咱们的y不是正态分布的，则应用狭义线性模型 _（Nelder＆Wedderburn，1972）_，其中y通过链接函数进行变换，但再次假如f（y）和x线性相关。如果不是这种状况，并且关系在x的范畴内变动，则可能不是最合适的。咱们在这里有一些抉择：

  • 咱们能够应用线性拟合，然而如果这样做的话，咱们会在数据的某些局部下面或者上面。
  • 咱们能够分为几类。我在上面的图中应用了三个，这是一个正当的抉择。同样，咱们可能处于数据某些局部之下或之上，而在类别之间的边界左近仿佛是精确的。例如，如果x = 49时，与x = 50相比，y是否有很大不同？
  • 咱们能够应用多项式之类的变换。上面，我应用三次多项式，因而模型适宜：。这些的组合使函数能够润滑地近似变动。这是一个很好的抉择，但可能会极其稳定，并可能在数据中引起相关性，从而升高拟合度。

4样条曲线

多项式的进一步细化是拟合“分段”多项式，咱们在数据范畴内将多项式链在一起以形容形态。“样条线”是分段多项式，以绘图员用来绘制曲线的工具命名。物理样条曲线是一种柔性条，能够蜿蜒成形，并由砝码固定。在结构数学样条曲线时，咱们有多项式函数，二阶导数间断，固定在“结”点上。

上面是一个ggplot2 对象，该 对象的 geom_smooth 的公式蕴含ns 函数中的“天然三次样条”  。这种样条曲线为“三次”，并且应用10个结

5润滑函数

样条曲线能够是润滑的或“摇晃的”，这能够通过扭转节点数（k）或应用润滑惩办来管制。如果咱们减少结的数目，它将更“摇晃”。这可能会更靠近数据，而且误差也会更小，但咱们开始“适度拟合”关系，并拟合咱们数据中的噪声。当咱们联合润滑惩办时，咱们会惩办模型中的复杂度，这有助于缩小适度拟合。

6狭义相加模型（GAM）

狭义加性模型（GAM）（Hastie，1984）应用润滑函数（如样条曲线）作为回归模型中的预测因子。

这些模型是严格可加的，这意味着咱们不能像失常回归那样应用交互项，然而咱们能够通过从新参数化作为一个更润滑的模型来实现同样的成果。事实并非如此，但实质上，咱们正转向一种模型，如：

摘自Wood _（2017）_的GAM的更正式示例 是：

其中：

  • i≡E（Yi），Y的冀望
  • Yi〜EF（i，i），Yi是一个响应变量，依据均值i和形态参数的指数族散布。
  • Ai是任何严格参数化模型重量的模型矩阵的一行，其中为对应的参数向量。
  • fi是协变量xk的润滑函数，其中k是每个函数的根底。

如果您要建设回归模型，但狐疑润滑拟合会做得更好，那么GAM是一个不错的抉择。它们适宜于非线性或有噪声的数据。

7 gam拟合

那么，如何 为上述S型数据建设 GAM模型？

在这里，我将应用三次样条回归 ：

gam(Y ~ s(X, bs="cr")

下面的设置意味着：

  • s()指定光滑器。还有其余选项，然而s是一个很好的默认选项
  • bs=“cr”通知它应用三次回归样条（'basis'）。
  • s函数计算出要应用的默认结数，然而您能够将其更改为k=10，例如10个结。

8模型输入：

查看模型摘要：

## ## Family: gaussian ## Link function: identity ## Parametric coefficients:##             Estimate Std. Error t value Pr(>|t|)    ## (Intercept)  43.9659     0.8305   52.94   <2e-16 ***## ---## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## ## Approximate significance of smooth terms:##        edf Ref.df     F p-value    ## s(X) 6.087  7.143 296.3  <2e-16 ***## ---## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## ## R-sq.(adj) =  0.876   Deviance explained = 87.9%## GCV = 211.94  Scale est. = 206.93    n = 300
  • 显示了咱们截距的模型系数，所有非润滑参数将在此处显示
  • 每个润滑项的总体含意如下。
  • 这是基于“无效自由度”（edf）的，因为咱们应用的样条函数能够扩大为许多参数，但咱们也在惩办它们并缩小它们的影响。

9查看模型：

该 gam.check() 函数可用于查看残差图，但它也能够测试润滑器以查看是否有足够的结来形容数据。然而如果p值很低，则须要更多的结。

## ## Method: GCV   Optimizer: magic## Smoothing parameter selection converged after 4 iterations.## The RMS GCV score gradient at convergence was 1.107369e-05 .## The Hessian was positive definite.## Model rank =  10 / 10 ## ## Basis dimension (k) checking results. Low p-value (k-index<1) may## indicate that k is too low, especially if edf is close to k'.## ##        k'  edf k-index p-value## s(X) 9.00 6.09     1.1    0.97

10它比线性模型好吗？

让咱们比照具备雷同数据的一般线性回归模型：

anova(my_lm, my_gam)
## Analysis of Variance Table## ## Model 1: Y ~ X## Model 2: Y ~ s(X, bs = "cr")##   Res.Df   RSS     Df Sum of Sq      F    Pr(>F)    ## 1 298.00 88154                                      ## 2 292.91 60613 5.0873     27540 26.161 < 2.2e-16 ***## ---## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

咱们的方差分析函数在这里执行了f测验，咱们的GAM模型显著优于线性回归。

11小结

所以，咱们看了什么是回归模型，咱们是如何解释一个变量y和另一个变量x的。其中一个根本假如是线性关系，但状况并非总是这样。当关系在x的范畴内变动时，咱们能够应用函数来扭转这个形态。一个很好的办法是在“结”点处将润滑曲线链接在一起，咱们称之为“样条曲线”

咱们能够在惯例回归中应用这些样条曲线，然而如果咱们在GAM的背景中应用它们，咱们同时预计了回归模型以及如何使咱们的模型更润滑。

下面的示例显示了基于样条的GAM，其拟合度比线性回归模型好得多。

12用GAM进行建模用电负荷工夫序列

我曾经筹备了一个文件，其中蕴含四个用电工夫序列来进行剖析。数据操作将由data.table程序包实现。

将提及的智能电表数据读到data.table。

DT <- as.data.table(read_feather("ind"))

应用GAM回归模型。将工作日的字符转换为整数，并应用recode包中的函数从新编码工作日：1.星期一，…，7星期日。

DT[, week_num := as.integer(car::recode(week,    "'Monday'='1';'Tuesday'='2';'Wednesday'='3';'Thursday'='4';    'Friday'='5';'Saturday'='6';'Sunday'='7'"))]

将信息存储在日期变量中，以简化工作。

n_type <- unique(DT[, type])n_date <- unique(DT[, date])n_weekdays <- unique(DT[, week])period <- 48

让咱们看一下用电量的一些数据并对其进行剖析。

data_r <- DT[(type == n_type[1] & date %in% n_date[57:70])] ggplot(data_r, aes(date_time, value)) +  geom_line() +  theme(panel.border = element_blank(),        panel.background = element_blank(),        panel.grid.minor = element_line(colour = "grey90"),        panel.grid.major = element_line(colour = "grey90"),        panel.grid.major.x = element_line(colour = "grey90"),        axis.text = element_text(size = 10),        axis.title = element_text(size = 12, face = "bold")) +  labs(x = "Date", y = "Load (kW)")

在绘制的工夫序列中能够看到两个次要的季节性：每日和每周。咱们在一天中有48个测量值，在一周中有7天，因而这将是咱们用来对因变量–电力负荷进行建模的自变量。

训练咱们的第一个GAM。通过平滑函数s对自变量建模，对于每日季节性，应用三次样条回归，对于每周季节性，应用P样条。

gam_1 <- gam(Load ~ s(Daily, bs = "cr", k = period) +               s(Weekly, bs = "ps", k = 7),             data = matrix_gam,             family = gaussian)

首先是可视化。

layout(matrix(1:2, nrow = 1))plot(gam_1, shade = TRUE)

咱们在这里能够看到变量对电力负荷的影响。在左图中，白天的负载峰值约为下午3点。在左边的图中，咱们能够看到在周末负载量缩小了。

点击题目查阅往期内容

R语言对用电负荷工夫序列数据进行K-medoids聚类建模和GAM回归

左右滑动查看更多

01

02

03

04

让咱们应用summary函数对第一个模型进行诊断。

## ## Family: gaussian ## Link function: identity ## ## Formula:## Load ~ s(Daily, bs = "cr", k = period) + s(Weekly, bs = "ps", ##     k = 7)## ## Parametric coefficients:##             Estimate Std. Error t value Pr(>|t|)    ## (Intercept)  2731.67      18.88   144.7   <2e-16 ***## ---## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## ## Approximate significance of smooth terms:##              edf Ref.df     F p-value    ## s(Daily)  10.159 12.688 119.8  <2e-16 ***## s(Weekly)  5.311  5.758 130.3  <2e-16 ***## ---## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## ## R-sq.(adj) =  0.772   Deviance explained = 77.7%## GCV = 2.4554e+05  Scale est. = 2.3953e+05  n = 672

EDF：预计的自由度–能够像对给定变量进行平滑解决那样来解释（较高的EDF值示意更简单的样条曲线）。P值：给定变量对因变量的统计显着性，通过F测验进行测验（越低越好）。调整后的R平方（越高越好）。咱们能够看到R-sq.（adj）值有点低。

让咱们绘制拟合值：

咱们须要将两个自变量的交互作用包含到模型中。

第一种交互类型对两个变量都应用了一个平滑函数。

gam_2 <- gam(Load ~ s(Daily, Weekly),           summary(gam_2)$r.sq
## [1] 0.9352108

R方值表明后果要好得多。

summary(gam_2)$s.table
##                     edf   Ref.df        F p-value## s(Daily,Weekly) 28.7008 28.99423 334.4754       0

仿佛也很好，p值为0，这意味着自变量很重要。拟合值图：

当初，让咱们尝试上述变量交互。这能够通过function实现te，也能够定义根本函数。

## [1] 0.9268452

与以前的模型类似gam_2。

summary(gam_3)$s.table
##                       edf   Ref.df        F p-value## te(Daily,Weekly) 23.65709 23.98741 354.5856       0

十分类似的后果。让咱们看一下拟合值：

与gam_2模型相比，只有一点点差别，看起来te拟合更好。

## [1] 0.9727604
summary(gam_4)$sp.criterion
##   GCV.Cp ## 34839.46
summary(gam_4)$s.table
##                       edf   Ref.df        F p-value## te(Daily,Weekly) 119.4117 149.6528 160.2065       0

咱们能够在这里看到R方略有回升。
让咱们绘制拟合值：

这仿佛比gam_3模型好得多。

## [1] 0.965618
summary(gam_4_fx)$s.table
##                  edf Ref.df        F       p-value## te(Daily,Weekly) 335    335 57.25389 5.289648e-199

咱们能够看到R平方比模型gam_4低，这是因为咱们适度拟合了模型。证实lambda和EDF的预计工作失常。

因而，让咱们在案例（模型）中尝试ti办法。

## [1] 0.9717469
summary(gam_5)$sp.criterion
##   GCV.Cp ## 35772.35
summary(gam_5)$s.table
##                        edf     Ref.df          F p-value## s(Daily)         22.583649  27.964970  444.19962       0## s(Weekly)         5.914531   5.995934 1014.72482       0## ti(Daily,Weekly) 85.310314 110.828814   41.22288       0

而后应用t2。

## [1] 0.9738273
summary(gam_6)$sp.criterion
##   GCV.Cp ## 32230.68
summary(gam_6)$s.table
##                       edf   Ref.df        F p-value## t2(Daily,Weekly) 98.12005 120.2345 86.70754       0

我还输入了最初三个模型的GCV得分值，这也是在一组拟合模型中抉择最佳模型的良好规范。咱们能够看到，对于t2相应模型gam_6，GCV值最低。

在统计中宽泛应用的其余模型抉择规范是AIC（Akaike信息准则）。让咱们看看三个模型：

AIC(gam_4, gam_5, gam_6)
##             df      AIC## gam_4 121.4117 8912.611## gam_5 115.8085 8932.746## gam_6 100.1200 8868.628

最低值在gam_6模型中。让咱们再次查看拟合值。

咱们能够看到的模型的拟合值gam_4和gam_6十分类似。能够应用软件包的更多可视化和模型诊断性能来比拟这两个模型。

第一个是function gam.check，它绘制了四个图：残差的QQ图，线性预测变量与残差，残差的直方图以及拟合值与因变量的关系图。让咱们诊断模型gam_4和gam_6。

gam.check(gam_4)

## ## Method: GCV   Optimizer: magic## Smoothing parameter selection converged after 7 iterations.## The RMS GCV score gradiant at convergence was 0.2833304 .## The Hessian was positive definite.## The estimated model rank was 336 (maximum possible: 336)## Model rank =  336 / 336 ## ## Basis dimension (k) checking results. Low p-value (k-index<1) may## indicate that k is too low, especially if edf is close to k'.## ##                      k'    edf k-index p-value## te(Daily,Weekly) 335.00 119.41    1.22       1
gam.check(gam_6)

## ## Method: GCV   Optimizer: magic## Smoothing parameter selection converged after 9 iterations.## The RMS GCV score gradiant at convergence was 0.05208856 .## The Hessian was positive definite.## The estimated model rank was 336 (maximum possible: 336)## Model rank =  336 / 336 ## ## Basis dimension (k) checking results. Low p-value (k-index<1) may## indicate that k is too low, especially if edf is close to k'.## ##                      k'    edf k-index p-value## t2(Daily,Weekly) 335.00  98.12    1.18       1

咱们能够再次看到模型十分类似，只是在直方图中能够看到一些差别。

layout(matrix(1:2, nrow = 1))plot(gam_4, rug = FALSE, se = FALSE, n2 = 80, main = "gam n.4 with te()")plot(gam_6, rug = FALSE, se = FALSE, n2 = 80, main = "gam n.6 with t2()")

该模型gam_6 有更多的“波浪形”的轮廓。因而，这意味着它对因变量的拟合度更高，而平滑因子更低。

vis.gam(gam_6, n.grid = 50, theta = 35, phi = 32, zlab = "",        ticktype = "detailed", color = "topo", main = "t2(D, W)")

咱们能够看到最高峰值是Daily变量的值靠近30（下午3点），而Weekly变量的值是1（星期一）。

vis.gam(gam_6, main = "t2(D, W)", plot.type = "contour",        color = "terrain", contour.col = "black", lwd = 2)

再次能够看到，电力负荷的最高值是星期一的下午3:00，直到星期四都十分类似，而后负荷在周末缩小。

本文摘选 《 R语言狭义相加模型（GAM）在电力负荷预测中的利用 》 ，点击“浏览原文”获取全文残缺材料。

点击题目查阅往期内容

分位数回归、GAM样条曲线、指数平滑和SARIMA对电力负荷工夫序列预测
实现狭义相加模型GAM和一般最小二乘(OLS)回归
R语言非参数模型厘定保险费率：部分回归、狭义相加模型GAM、样条回归
R语言狭义加性模型GAMs剖析温度、臭氧环境数据绘制偏回归图与偏残差图
R语言狭义相加（加性）模型（GAMs）与润滑函数可视化
R语言里的非线性模型：多项式回归、部分样条、平滑样条、 狭义相加模型GAM剖析
R语言用规范最小二乘OLS，狭义相加模型GAM ，样条函数进行逻辑回归LOGISTIC分类
R语言ISLR工资数据进行多项式回归和样条回归剖析
R语言中的多项式回归、部分回归、核平滑和平滑样条回归模型
R语言用泊松Poisson回归、GAM样条曲线模型预测骑自行车者的数量
R语言分位数回归、GAM样条曲线、指数平滑和SARIMA对电力负荷工夫序列预测
R语言中的多项式回归、B样条曲线(B-spline Curves)回归
R语言狭义相加模型 (GAMs)剖析预测CO2工夫序列数据
R语言中实现狭义相加模型GAM和一般最小二乘(OLS)回归
在r语言中应用GAM（狭义相加模型）进行电力负荷工夫序列剖析
R语言用泊松Poisson回归、GAM样条曲线模型预测骑自行车者的数量
Python用狭义加性模型GAM进行工夫序列剖析
R语言狭义线性模型GLM、多项式回归和狭义可加模型GAM预测泰坦尼克号幸存者
R语言中的狭义线性模型（GLM）和狭义相加模型（GAM）：多元（平滑）回归剖析保险资金投资组合信用风险敞口
R语言对用电负荷工夫序列数据进行K-medoids聚类建模和GAM回归
对用电负荷工夫序列数据进行K-medoids聚类建模和GAM回归