R中的`qr.qy()`函数

Question

我找不到有关此功能的确切操作的文档。我有一个矩阵的 QR 分解X：

X = matrix(c(1,1,1,1,-1.5,-0.5,0.5,1.5,2.25,0.25,0.25, 2.25,-3.275,-0.125,0.125,3.375), 
nrow=4, byrow=F)

     [,1] [,2] [,3]   [,4]
[1,]    1 -1.5 2.25 -3.375
[2,]    1 -0.5 0.25 -0.125
[3,]    1  0.5 0.25  0.125
[4,]    1  1.5 2.25  3.375

该函数qr(X)产生一个列表：

$qr (rounding output)
     [,1]       [,2]          [,3]          [,4]
[1,] -2.0          0          -2.5             0
[2,]  0.5     -2.236             0        -4.583
[3,]  0.5      0.447             2             0 
[4,]  0.5      0.894        -0.929        -1.341

$rank
[1] 4

$qraux
[1] 1.500000 1.000000 1.368524 1.341641

$pivot
[1] 1 2 3 4

attr(,"class")
[1] "qr"

我选择 的对角线元素qr(X)$qr，我命名为z：

z = qr(X)$qr
z = Q * (row(Q) == col(Q))

     [,1]      [,2] [,3]      [,4]
[1,]   -2  0.000000    0  0.000000
[2,]    0 -2.236068    0  0.000000
[3,]    0  0.000000    2  0.000000
[4,]    0  0.000000    0 -1.341641

到现在为止还挺好。现在下一个电话我不明白：

(raw = qr.qy(qr(X), z))

     [,1] [,2] [,3] [,4]
[1,]    1 -1.5    1 -0.3
[2,]    1 -0.5   -1  0.9
[3,]    1  0.5   -1 -0.9
[4,]    1  1.5    1  0.3

---

取得一些进展：

所以，多亏了答案和一些阅读，我认为对象qr(X)$qr完全包含在上三角形中的 R ：

     [,1]       [,2]          [,3]          [,4]
[1,] -2.0          0          -2.5             0
[2,]          -2.236             0        -4.583
[3,]                             2             0 
[4,]                                      -1.341

的下方三角形qr(X)$qr包含有关 Q 的信息：

     [,1]       [,2]          [,3]          [,4]
[1,]                                            
[2,]  0.5                                       
[3,]  0.5      0.447                             
[4,]  0.5      0.894        -0.929      

以某种方式调用qr.Q(qr(X))返回 Q 使用内部函数qr.qy()与 qr() 和 1 的对角矩阵作为输入。

但是这个操作是如何进行的呢？Q 右上角的其余部分如何填充？我认为它利用了$qraux，但它是如何做到的：

     [,1]       [,2] [,3]       [,4]
[1,] -0.5  0.6708204  0.5  0.2236068
[2,] -0.5  0.2236068 -0.5 -0.6708204
[3,] -0.5 -0.2236068 -0.5  0.6708204
[4,] -0.5 -0.6708204  0.5 -0.2236068

简而言之，qr.qy() 具体是如何工作的？

我刚刚发现：“qy.qr()：返回矩阵乘法的结果：Q %*% y，其中 Q 是由 qr 表示的 order-nrow(x) 正交（或酉）变换。”

Answer 1

Q来自 QR 分解的矩阵仅隐含在函数的返回值中qr。列表元素qr是 Q 矩阵的紧凑表示；它包含在列表元素的下三角部分qr和向量中qraux。RQR 分解的上三角矩阵是qr返回值中列表元素的上三角部分。

qr.qy经过一些中间步骤后的 R 函数最终会调用 Lapack 子例程dormqr，它不会Q显式生成矩阵。它使用包含在列表元素qr和 中的信息qraux。请参阅http://www.netlib.org/lapack/explore-html/da/d82/dormqr_8f.html。

所以qr.qy不会将 的紧凑形式Q转换为实际的Q. 它使用紧凑形式来计算Q %*% z。

R 函数qr.Q（您已经完成）使用qr.qy对角线上带有 1 的对角矩阵来生成Q.

为什么要这样做？出于效率原因。

使用以下代码，您可以检查这一点：

library(rbenchmark)

benchmark(qr.qy(XQR,z), {Q <- qr.Q(qr(X)); Q %*% z},  { Q %*% z},
          replications=10000, 
          columns=c("test","replications","elapsed","relative")  )

带输出

                                     test replications elapsed relative
3                           {    Q %*% z}        10000   0.022    1.000
2       {    Q <- qr.Q(qr(X));   Q %*% z}        10000   0.486   22.091
1                           qr.qy(XQR, z)        10000   0.152    6.909

教训：仅Q当您确实需要以显式形式使用它并且需要使用不同的输入矩阵多次生成它时才生成。如果Q是固定的并且没有改变，那么你可以使用 Q %*% z.