[Rd] solve.default(a, b) should not ignore names(dimnames(.))

[Mikael Jagan]

On 2022-11-10 1:10 pm, Mikael Jagan wrote:
> Hello,
> 
> Currently, solve.default() is defined such that
> 
>   > identical(dimnames(solve.default(a, b)),
>               list(dimnames(a)[[2L]], dimnames(b)[[2L]]))
> 
> is always TRUE, i.e., ignoring names(dimnames(a)) and names(dimnames(b)).
> Hence we see:
> 
>   > a <- b <- diag(2L)
>   > dimnames(a) <- list(A1 = c("a11", "a12"), A2 = c("a21", "a22"))
>   > dimnames(b) <- list(B1 = c("b11", "b12"), B2 = c("b21", "b22"))
>   > a
>        A2
> A1    a21 a22
>     a11   1   0
>     a12   0   1
>   > b
>        B2
> B1    b21 b22
>     b11   1   0
>     b12   0   1
>   > solve.default(a, b)
>       b21 b22
> a21   1   0
> a22   0   1
> 
> I claim that solve.default() should be changed to instead give:
> 
>   > identical(dimnames(solve.default(a, b)),
>               c(dimnames(a)[2L], dimnames(b)[2L]))
> 
> This would make solve.default() consistent with `%*%`, which
> _does_ respect names(dimnames(.)) :
> 
>   > a %*% b
>        B2
> A1    b21 b22
>     a11   1   0
>     a12   0   1
> 
> If others agree, then I would submit a minimal patch to the R-level
> solve.default() in src/library/base/R/solve.R and to the C-level
> La_solve() and La_solve_cmplx() in src/modules/lapack/Lapack.c ...
> 

Well, taking a closer look, I see more related issues, which might
be considered simultaneously.  (Sorry if this is too much for one thread.)

* qr.solve(a, b) and solve.qr(qr(a), b) also ignore names(dimnames(.))

   > qr.solve(a, b)
       b21 b22
   a21   1   0
   a22   0   1
   > solve.qr(qr(a), b)
       b21 b22
   a21   1   0
   a22   0   1


* dimnames(qr.solve(a)) and dimnames(solve.qr(qr(a))) do not agree
   with dimnames(solve.default(a)), i.e., when 'b' is missing

   > solve.default(a)
       a11 a12
   a21   1   0
   a22   0   1
   > qr.solve(a)
       [,1] [,2]
   a21    1    0
   a22    0    1
   > solve.qr(qr(a))
       [,1] [,2]
   a21    1    0
   a22    0    1


* More controversially: qr.qy(qr, y) and qr.qty(qr, y) are currently
   documented as retaining the 'dimnames' of 'y', and they _do_.  But
   _conventionally_ the rownames of matrix products are taken from the
   first factor, in this case 'Q'.  It is not currently documented what
   the 'dimnames' of the implicitly stored 'Q' and 'R' factors should be,
   leading to inconsistencies like this:

   > a
        A2
   A1    a21 a22
     a11   1   0
     a12   0   1
   > qr.Q(qr(a)) %*% qr.R(qr(a))
         A2
          a21 a22
     [1,]   1   0
     [2,]   0   1

   Hence I propose to enforce {and partly document} the following:
   1. that Q=qr.Q(qr(X)) gets the rownames of X and NULL colnames
   2. that R=qr.R(qr(X)) gets the (permuted) colnames of X and NULL rownames
   3. that dimnames(qr.qy (qr(X), y)) <=> dimnames(  Q  %*% y)
   4. that dimnames(qr.qty(qr(X), y)) <=> dimnames(t(Q) %*% y)

   I'm not sure what the level of disruption on CRAN would be,
   but internal consistency here would certainly be nice ...

> Mikael