【问题标题】:Calculating a linear trend line for every row of a table in R为R中表格的每一行计算线性趋势线
【发布时间】:2014-02-14 16:00:38
【问题描述】:

是否有可能在不使用循环的情况下对数据帧的每一行进行线性回归?趋势线的输出(截距+斜率)应作为新列添加到原始数据框中。

为了让我的意图更清楚,我准备了一个非常小的数据示例:

day1 <- c(1,3,1)
day2 <- c(2,2,1)
day3 <- c(3,1,5)
output.intercept <- c(0,4,-1.66667)
output.slope <- c(1,-1,2)
data <- data.frame(day1,day2,day3,output.intercept,output.slope)

输入变量为day1-3;假设这些是连续 3 天不同商店的销售额。我想要做的是计算 3 行的线性趋势线,并将输出参数添加到原始表(参见 output.intercept + output.slope)作为新列。

该解决方案在计算时间方面应该非常有效,因为实际数据帧有许多 100k 行。

最好的,克里斯托夫

【问题讨论】:

  • 什么是响应变量?
  • @SvenHohenstein 响应显示,协变量隐含为1:3(在这种情况下),seq_len(nrow(dat)) 在更一般的情况下。

标签: r regression trendline


【解决方案1】:
design.mat <- cbind(1,1:3)
response.mat <- t(data[,1:3])

reg <- lm.fit(design.mat, response.mat)$coefficients
data <- cbind(data, t(reg))
#  day1 day2 day3 output.intercept output.slope        x1 x2
#1    1    2    3          0.00000            1  0.000000  1
#2    3    2    1          4.00000           -1  4.000000 -1
#3    1    1    5         -1.66667            2 -1.666667  2

但是,如果您有大量数据,由于内存限制,可能需要循环。如果是这种情况,我会使用长格式 data.table 并使用包的 by 语法进行循环。

【讨论】:

  • 哇,完美。多谢!稍后我将使用大数据集进行尝试。什么是精确的“design.mat”?模拟 x 值?
  • 如果你不知道什么是设计矩阵,你应该学习一本关于回归的教科书。
  • 再次感谢,即使使用大数据,该解决方案也能完美运行。然而,当“数据”包含 NA 形式的缺失数据点时,就会出现一个问题。 ("Error in lm.fit(design.mat, response.mat) : NA/NaN/Inf in 'y'") 有什么办法可以解决一些缺失数据点的问题吗?我已经尝试将 na.exclude 函数包含到 lm.fit 语句中,但在这种情况下它不起作用。
  • 在使用lm.fit 之前,您需要从设计和响应矩阵中删除包含NA 值的行。这是lm 为您方便地做的许多事情之一,但总的来说这会花费很多性能。
  • 不幸的是,删除包含 NA 的行不是一种选择,至少在一种特殊情况下,因为在我的一个数据表中,几乎每一列都包含 NA。难道没有另一种可能只对可用数据进行回归并且“简单地”不考虑NA的值吗?否则我将不得不事先处理我的原始数据文件。
【解决方案2】:

使用您的数据,

day1 <- c(1,3,1)
day2 <- c(2,2,1)
day3 <- c(3,1,5)
output.intercept <- c(0,4,-1.66667)
output.slope <- c(1,-1,2)
dat <- data.frame(day1,day2,day3)

我想你想要这样的东西:

fits <- lm.fit(cbind(1, seq_len(nrow(dat))), t(dat))
t(coef(fits))

这给了

R> t(coef(fits))
         x1 x2
[1,]  0.000  1
[2,]  4.000 -1
[3,] -1.667  2

这些可以像这样添加到dat

dat <- cbind(dat, t(coef(fits)))
names(dat)[-(1:3)] <- c("Intercept","Slope")

R> dat
  day1 day2 day3 Intercept Slope
1    1    2    3     0.000     1
2    3    2    1     4.000    -1
3    1    1    5    -1.667     2

如果您可以控制数据最初的排列方式,那么以另一种方式存储数据可能会更容易,将列作为时间序列而不是行,因为这样可以避免在以下情况下转置大矩阵适合通过lm.fit()。理想情况下,您最初希望数据排列如下:

     [,1] [,2] [,3]
day1    1    3    1
day2    2    2    1
day3    3    1    5

即行作为时间点,而不是您现在拥有的单个系列。这是因为 R 期望数据的排列方式。请注意,我们必须在lm.fit() 调用中转置您的dat,这将需要一个大对象的副本。因此,如果您可以在这些数据进入 R 之前控制这些数据的排列/提供方式,那将有助于解决大问题。

lm.fit()lm() 使用的底层精简代码,但我们避免了解析公式和创建模型矩阵的复杂性。如果您想要更高效,您可能需要自己进行 QR 分解(代码在 lm.fit() 中执行此操作),因为 lm.fit() 有一些事情可以作为您可以做的健全性检查如果您确定您的数据不会导致奇异矩阵等。

【讨论】:

  • 非常感谢。我意识到我在 R 中还有很多东西要学,即使是基本的东西。感谢您提供有关数据结构的提示。由于我事先在 R 中进行了一些数据准备,因此我可以控制数据排列。我认为这样会更有效率,因为我的真实数据文件包含 600k 行和只有 100 列。
  • 请注意:我假设语句“fits
【解决方案3】:

我遇到了和 OP 一样的问题。此解决方案适用于具有 NA 的数据。在这种情况下,所有以前的答案都会对我产生错误:

slp = function(x) {
  y = t(x)
  y = y[!is.na(y)]
  len = length(y):1
  b = cov(y,len)/var(len)
  return(b)}

reg_slp <- apply(data,1,slp)

仅获取斜率,但可以轻松添加截距。我怀疑这是否特别有效,但在我的情况下它是有效的。

【讨论】:

    【解决方案4】:

    还是这样?

    day1 <- c(1,3,1)
    day2 <- c(2,2,1)
    day3 <- c(3,1,5)
    data <- data.frame(day1,day2,day3)
    y<-1:3
    
    reg<-apply(data,1,function(x) lm(as.numeric(x)~y))
    data[,c("intercept","slope")]<-rbind(reg[[1]]$coef,reg[[2]]$coef,reg[[3]]$coef)
    

    【讨论】:

    • 说的没错,但效率不高。请注意,lm() 必须解析公式nrow(dat) 次,如果您执行 3 次则很快,如果您执行 100K 次则很慢。此外,这遗漏了lm() 的一个特性,即它接受矩阵响应。所以你根本不需要apply() 或循环;您可以将所有系列放在一个 lm() 调用中:lm(t(data[, 1:3]) ~ I(1:3))。但是,如果您想提高效率,您不想解析公式并生成 model.frame 和 model.matrix 以及所有额外的废话 lm()。使用lm.fit() 进行改进。
    猜你喜欢
    • 1970-01-01
    • 2018-07-18
    • 2021-07-25
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2020-01-16
    • 1970-01-01
    • 2011-02-15
    相关资源
    最近更新 更多