lm和biglm的QR分解不同吗? [英] QR decomposition different in lm and biglm?

问题描述

我正在尝试从biglm中使用的QR分解中恢复R矩阵.为此，我使用了vcov.biglm中的部分代码，并将其放入类似这样的函数中:

I'm trying to recover the R matrix from the QR decomposition used in biglm. For this I am using a portion of the code in vcov.biglm and put it into a function like so:

qr.R.biglm <- function (object, ...) {
    # Return the qr.R matrix from a biglm object
    object$qr <- .Call("singcheckQR", object$qr)
    p <- length(object$qr$D)
    R <- diag(p)
    R[row(R) > col(R)] <- object$qr$rbar
    R <- t(R)
    R <- sqrt(object$qr$D) * R
    dimnames(R) <- list(object$names, object$names)
    return(R)
}

更具体地说，我试图从基本包中获得与使用qr.R相同的结果，该基本包用于"qr"类的QR分解，例如lm类(lm $ qr)中包含的那些分解.基本功能的代码如下:

More specifically, I'm trying to get the same result as using qr.R from the base package, which is used on QR decompositions of class "qr" such as those contained in the lm class (lm$qr). The code for the base function is as follows:

qr.R <- function (qr, complete = FALSE) {
    if (!is.qr(qr)) 
        stop("argument is not a QR decomposition")
    R <- qr$qr
    if (!complete) 
        R <- R[seq.int(min(dim(R))), , drop = FALSE]
    R[row(R) > col(R)] <- 0
    R
}

对于除样本外的样本，我设法获得相同的结果.

I manage to get the same result for a sample regression, except for the signs.

x <- as.data.frame(matrix(rnorm(100 * 10), 100, 10))
y <- seq.int(1, 100)
fit.lm <- lm("y ~ .", data =  cbind(y, x))
R.lm <- qr.R(fit.lm$qr)

library(biglm)
fmla <- as.formula(paste("y ~ ", paste(colnames(x), collapse = "+")))
fit.biglm <- biglm(fmla, data = cbind(y, x))
R.biglm <- qr.R.biglm(fit.biglm)

将两者进行比较，很明显，绝对值匹配，但符号不匹配.

Comparing both, it's clear that the absolute values match, but not the signs.

mean(abs(R.lm) - abs(R.biglm) < 1e-6)
[1] 1
mean(R.lm - R.biglm < 1e-6)
[1] 0.9338843

我不太清楚为什么会这样.我希望能够获得与biglm的lm相同的R矩阵结果.

I can't quite figure out why this is. I would like to be able to get the same result for the R matrix as lm from biglm.

推荐答案

两个R矩阵之间的区别在于 biglm 显然执行其旋转，使得R的对角元素均为正，而 lm (或者实际上是它调用的例程)没有施加此类约束. (一种策略或另一种策略不应具有数值上的优势，因此区别仅是惯例之一，即AFAIKT.)

The difference between the two R matrices is that biglm apparently performs its rotations such that R's diagonal elements are all positive, while lm (or, really, the routines it calls) imposes no such constraint. (There should be no numerical advantage to one strategy or the other, so the difference is just one of convention, AFAIKT.)

您可以通过施加附加约束来使 lm 的结果与 biglm 的结果相同.我会使用一个反射矩阵，该矩阵将列乘以1或-1，这样对角元素都以正数结尾:

You can make lm's results identical to biglm's by imposing that additional constraint yourself. I'd use a reflection matrix that multiplies columns by either 1 or -1, such that the diagonal elements all end up positive:

## Apply the necessary reflections
R.lm2 <- diag(sign(diag(R.lm))) %*% R.lm

## Show that they did the job
mean(R.lm2 - R.biglm < 1e-6)
# [1] 1

这篇关于lm和biglm的QR分解不同吗?的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！