r - 从 GAM 平滑对象确定导数

Question

我有一个非常简单的时间序列数据集，由单个变量(“AVERAGE”)的年平均值组成。我希望研究时间序列的“趋势”分量的变化率(一阶导数)和加速度(二阶导数)以及相关的标准误差。我使用MGCV的GAM和PREDICT函数获得了“趋势”，如下所示:

A <- gam(AVERAGE ~ s(YEAR), data=DF, na.action=na.omit)
B <- predict(A, type="response", se.fit=TRUE)

我已经通过 2 个单独的方法确定了导数，应用高自由度三次平滑样条并通过一阶和二阶差异(轻微平滑)和自举来近似误差，两者都产生可比较的结果。

我注意到“gam.fit3”函数有助于确定二阶导数，但不会直接调用。我还注意到，使用类型为“lpmatrix”的“predict.gam”有助于平滑的导数。我想直接使用“GAM”函数来计算一阶和二阶导数，但不够熟练来计算或提取这些导数。我试图在“Predict.gam”帮助页面末尾为一个变量重新配置伍德的示例，但没有成功。任何让我朝着正确方向前进的帮助都会很棒。谢谢菲尔。

Answer 1

来自 predict.gam 的示例用途 finite differences逼近平滑项的导数

这是对单个预测器模型执行此操作的示例。这比帮助中的示例更直接。

A <- gam(AVERAGE ~ s(YEAR), data=DF, na.action=na.omit)
# new data for prediction
newDF <- with(DF, data.frame(YEAR = unique(YEAR)))
# prediction of smoothed estimates at each unique year value
# with standard error    
B <- predict(A,  newDF, type="response", se.fit=TRUE)


# finite difference approach to derivatives following
# example from ?predict.gam

eps <- 1e-7
X0 <- predict(A, newDF, type = 'lpmatrix')


newDFeps_p <- newDF + eps

X1 <- predict(A, newDFeps_p, type = 'lpmatrix')

# finite difference approximation of first derivative
# the design matrix
Xp <- (X0 - X1) / eps

# first derivative
fd_d1 <- Xp %*% coef(A)

# second derivative
newDFeps_m <- newDF - eps

X_1 <- predict(A, newDFeps_m, type = 'lpmatrix')
# design matrix for second derivative
Xpp <- (X1 + X_1 - 2*X0)  / eps^2
# second derivative
fd_d2 <- Xpp %*% coef(A)

如果您使用 bootstrap 来获得置信区间，您应该能够获得这些近似值的置信区间。