r - 如何从超大数据集(尤其是 rxGlm 输出)绘制交互效果

Question

我正在计算 glm 模型关闭庞大的数据数据集。 glm 甚至 speedglm 都需要几天的时间来计算。

我目前有大约 300 万个观察值和总共 400 个变量，其中只有一些用于回归。在我的回归中，我使用 4 个整数自变量( iv1 、 iv2 、 iv3 、 iv4 )，1 个二元自变量作为因子( iv5 )，交互项是 x * y 和 x104016161651616100000000000000000000000000000000000000000000 .最后，我在 y 年和公司 ID ff1 有固定的影响。我有 15 年和 3000 家公司。我通过将它们添加为因素来介绍固定效应。我观察到，特别是 3000 家公司的固定效应使 ff2 stats 和 glm 的计算速度非常慢。

因此，我决定尝试 Microsoft R 的 speedglm (RevoScaleR)，因为它可以处理更多线程和处理器内核。事实上，分析的速度要快得多。此外，我将子样本的结果与标准 rxGlm 之一进行了比较，并且它们匹配。

我使用了以下功能:

mod1 <- rxGlm(formula = dv ~ 
                      iv1 + iv2 + iv3+ 
                      iv4 + iv5 +
                      x * y +
                      ff1  + ff2,
                    family = binomial(link = "probit"), data = dat,
                    dropFirst = TRUE, dropMain = FALSE, covCoef = TRUE, cube = FALSE)

但是，在尝试使用 glm 包绘制交互项时，我遇到了问题。调用以下函数后，我收到以下错误:

> plot(effect("x*y", mod1))
Error in terms.default(model) : no terms component nor attribute

我认为问题在于 effects 不存储绘制交互所需的数据。我相信是因为 rxGlm 对象比 rxGlm 对象小很多，因此可能包含更少的数据(80 MB 对几 GB)。

我现在尝试通过 glm 将 rxGlm 对象转换为 glm 。尽管如此，as.glm() 调用不会产生结果并导致以下错误消息:

Error in dnorm(eta) : 
  Non-numerical argument for mathematical function
In addition: Warning messages:
1: In model.matrix.default(mod, data = list(dv = c(1L, 2L,  :
  variable 'x for y' is absent, its contrast will be ignored

如果我现在将原始 glm 与“转换后的 glm”进行比较，我发现转换后的 glm 包含的项目要少得多。例如，它不包含 effects() 并且为了对比，它仅声明每个变量的 effects。

我现在主要寻找一种以某种格式转置 contr.treatment 输出对象的方法，以便我可以将 if 与 rxGlm 函数一起使用。如果没有办法这样做，我如何使用 effect() 包中的函数获得交互图，例如 RevoScaleR ？ rxLinePlot() 的绘图速度也相当快，但是，我还没有找到一种方法来从中获得典型的交互效果图。我想避免先计算完整的 rxLinePlot() 模型然后再绘图，因为这需要很长时间。

Answer 1

如果你能得到系数，你不能自己滚动吗？
这不会是数据集大小问题

# ex. data
n = 2000
dat <- data.frame( dv = sample(0:1, size = n, rep = TRUE), 
                   iv1 = sample(1:10, size = n, rep = TRUE),
                   iv2 = sample(1:10, size = n, rep = TRUE),
                   iv3 = sample(1:10, size = n, rep = TRUE),
                   iv4 = sample(0:10, size = n, rep = TRUE),
                   iv5 = as.factor(sample(0:1, size = n, rep = TRUE)),
                   x = sample(1:100, size = n, rep = TRUE),
                   y = as.factor(sample(0:1, size = n, rep = TRUE)),
                   ff1  = as.factor(sample(1:15, size = n, rep = TRUE)),
                   ff2  = as.factor(sample(1:100, size = n, rep = TRUE))
                   )

mod1 <- glm(formula = dv ~ 
                      iv1 + iv2 + iv3+ 
                      iv4 + iv5 +
                      x * y +
                      ff1  + ff2,
                    family = binomial(link = "probit"), data = dat)

# coefficients for x, y and their interaction
x1 <- coef(mod1)['x']
y1 <- coef(mod1)['y1']
xy <- coef(mod1)['x:y1']

x <- 1:100
a <- x1*x
b <- x1*x + y1 + xy*x

plot(a~x, type= 'line', col = 'red', xlim = c(0,max(x)), ylim = range(c(a, b)))
lines(b~x, col = 'blue')
legend('topright', c('y = 0', 'y = 1'), col = c('red', 'blue'))