r回归模型中的元组变量答案

【问题标题】：tuple variable in r regression modelr回归模型中的元组变量
【发布时间】：2018-08-16 09:56:13
【问题描述】：

我试图在我的模型中使用一天中的时间作为自变量。由于时间是一个循环变量，我将其转换为 (sin(pi * hour / 12), cos(pi * hour / 12))。

我四处搜索，但我仍然不知道如何在 R 中使用 (sin, cos) 格式的向量/元组值创建列。不知道lm、glm、glm.nb（MASS）和glmer（lme4）的模型能不能支持这种数据。

对不起，我是新手。如果向量类型的变量不应该包含在回归模型中，我会去 Cross Validated (stats) 寻求有关处理循环变量的建议。请帮助并分享您的经验，谢谢！

【问题讨论】：

你能用数学术语描述你追求什么吗？这听起来像是一个编程问题，所以简历可能跑题了？
我认为您可以将这两个组件作为单独的预测变量包含在模型中，进行拟合、预测和反向转换到时间。
@RomanLuštrik 感谢您的帮助！没错，这主要是一个编程/软件问题，所以我决定在这里问。例如，我想知道 R 是否支持数据框单元格中的元组值，以及是否可以从其他 2 列创建一列。此外，如果回归包接受这种数据。研究设计是我的问题的背景。
@MikkoMarttila 感谢您的帮助！您在此处对 SO 的建议“交叉验证”了我在 CV 上阅读的答案：sin 和 cos 作为单独的变量进入回归模型。不幸的是，我看到的答案没有提到如何用 2 个变量一起解释结果。我也没有机会阅读提到的 R 中的循环统计一书......

标签： r regression

【解决方案1】：

这是一个很好的问题。

在内部，R 使用矩阵将模型拟合到数据。您需要生成一个块矩阵，而不是将您的模型视为元组的函数；例如，这就是样条线的实现方式。

在您的示例中，一列包含sin(x) 和另一列cox(x)，并且块用类属性标记；在内部管理这个对象的函数是makepredictcall 和predict。

任何使用标准 model.frame / model.matrix 处理的模型都应该与此兼容。

sincos <- function(x, period=168/2/pi) {
  structure(cbind(`_sin`=sin(x/period), 
                  `_cos`=cos(x/period)),
            class="sincos", 
            period=period)
}

这里我们将类和周期设置为属性。

makepredictcall.sincos <- function(var, call){
  if (as.character(call)[1L] != "sincos")
    return(call)
  call["period"] <- attr(var, "period")
  call
}

根据需要设置通话时间。

predict.sincos <- function(object, newx, ...)
{
  if(missing(newx))
    return(object)
  sincos(newx, period=attr(object, "period"))
}

使用拟合模型中的周期调用我们的函数。

这是一个使用lm的简短示例：

#FAKE DATA EXAMPLE
N <- 1000
hr <- sample(168, N, replace = TRUE)
Y = 5 + sinpi(hr * 2/168) + cospi(hr * 2/168) + rnorm(N)

lm(Y~sinpi(hr*2/168)+cospi(hr*2/168))
#> 
#> Call:
#> lm(formula = Y ~ sinpi(hr * 2/168) + cospi(hr * 2/168))
#> 
#> Coefficients:
#>       (Intercept)  sinpi(hr * 2/168)  cospi(hr * 2/168)  
#>            5.0078             1.0243             0.9637

我们的自定义函数完全匹配：

lm(Y~sincos(hr))
#> 
#> Call:
#> lm(formula = Y ~ sincos(hr))
#> 
#> Coefficients:
#>    (Intercept)  sincos(hr)_sin  sincos(hr)_cos  
#>         5.0078          1.0243          0.9637

其他函数也将能够分辨出这两列是模型中的一个术语：

anova(lm(Y~sincos(hr)))
#> Analysis of Variance Table
#> 
#> Response: Y
#>             Df  Sum Sq Mean Sq F value    Pr(>F)    
#> sincos(hr)   2 1051.05  525.53  553.24 < 2.2e-16 ***
#> Residuals  997  947.06    0.95                      
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

我会在接下来的几天内将此添加到 stackoverflow 包中，如果其他人觉得这有帮助的话。

【讨论】：