RcppArmadillo中的fastLm和fastLmPure函数的区别

Question

这里是一个例子：

require(Rcpp)
require(RcppArmadillo)
require(zoo)
require(repmis)

myData <- source_DropboxData(file = "example.csv", 
                              key = "cbrmkkbssu5bn96", sep = ",", header = TRUE)

dolm = function(x) coef(fastLmPure(as.matrix(x[,2]), x[,1]))

myCoef = rollapply(myData, 260, dolm, by.column = FALSE)

summary(myCoef) # 80923 NA's

dolm2 = function(x) coef(fastLm(x[,1] ~ x[,2] + 0, data = as.data.frame(x)))

myCoef2 = rollapply(myData, 260, dolm2, by.column = FALSE)

summary(myCoef2) #  0 NA's 

在上面的示例中，使用fastLmPure的第一种方法会在输出中产生NAs，而使用fastLm的第二种方法则不会。

这里是用R编写的fastLm和fastLmPure函数的链接：

https://github.com/RcppCore/RcppArmadillo/blob/master/R/fastLm.R

这里是链接到底层的C ++编写的fastLm函数：

https://github.com/RcppCore/RcppArmadillo/blob/master/src/fastLm.cpp

从这些链接和RcppArmadillo文档来看，我不清楚是什么原因导致输出结果不同？为什么第二个输出没有NAs？最重要的问题是，哪个例程/代码部分防止了第二种方法中出现NAs，它是如何实现的？

- Bobby Digital

---

我邀请你阅读源代码。它真的并不复杂... - Dirk Eddelbuettel

@DirkEddelbuettel 谢谢，我已经找到问题所在了。顺便说一句，提供已经在问题中提供的链接非常有帮助 ;) - Bobby Digital

Answer 1

你正在使用两个不同的接口调用两个不同的函数。

特别地，当通过公式 y ~ X 使用 fastLm() 时，将依赖 R 内部（并且非常缓慢!!）函数为您创建向量和矩阵，这些向量和矩阵对应于 fastLm(X, y)。

以下是一个设置示例：

R> data(mtcars)
R> lm(mpg ~ cyl + disp + hp + wt - 1, data=mtcars)

Call:
lm(formula = mpg ~ cyl + disp + hp + wt - 1, data = mtcars)

Coefficients:
    cyl     disp       hp       wt  
 5.3560  -0.1206  -0.0313   5.6913  

R> fastLm(mpg ~ cyl + disp + hp + wt - 1, data=mtcars)

Call:
fastLm.formula(formula = mpg ~ cyl + disp + hp + wt - 1, data = mtcars)

Coefficients:
      cyl      disp        hp        wt 
 5.356014 -0.120609 -0.031306  5.691273 
R> fastLm(mtcars[, c("cyl","disp","hp","wt")], mtcars[,"mpg"])

Call:
fastLm.default(X = mtcars[, c("cyl", "disp", "hp", "wt")], y = mtcars[, 
    "mpg"])

Coefficients:
      cyl      disp        hp        wt 
 5.356014 -0.120609 -0.031306  5.691273 
R> 

现在我们需要在左右两侧都添加一个NA。为了方便索引，我们将使用整个行：

R> mtcars[7, ] <- NA
R> lm(mpg ~ cyl + disp + hp + wt - 1, data=mtcars)

Call:
lm(formula = mpg ~ cyl + disp + hp + wt - 1, data = mtcars)

Coefficients:
    cyl     disp       hp       wt  
 5.3501  -0.1215  -0.0332   5.8281  

R> fastLm(mpg ~ cyl + disp + hp + wt - 1, data=mtcars)

Call:
fastLm.formula(formula = mpg ~ cyl + disp + hp + wt - 1, data = mtcars)

Coefficients:
      cyl      disp        hp        wt 
 5.350102 -0.121478 -0.033184  5.828065 
R> fastLm(na.omit(mtcars[, c("cyl","disp","hp","wt")]), na.omit(mtcars[,"mpg"]))

Call:
fastLm.default(X = na.omit(mtcars[, c("cyl", "disp", "hp", "wt")]), 
    y = na.omit(mtcars[, "mpg"]))

Coefficients:
      cyl      disp        hp        wt 
 5.350102 -0.121478 -0.033184  5.828065 
R> 

这里的关键是：只要我们在缺失值上保持一致，不管采用哪种方法，结果仍然相同。

- Dirk Eddelbuettel