R での QR 分解は qr(x, LAPACK=TRUE), qr(x, LAPACK=FALSE) がある(後者がデフォルト，すなわち qr(x)）
R factor を取り出すには，qr.R() がある。
Julia でも同じ関数名なのだが，特別な引数はない。
qr(x) とした場合には Q factor, R factor の2要素が返される。
ここで混乱が始まる。しかるべく取り出した R の R factor と Julia の R factor が違う。
R の qr(x, LAPACK=TRUE) は，Julia の LAPACK.geqp3!(x) に対応する。 pivoted QR decomposition
R の qr(x, LAPACK=FALSE) は，Julia の LAPACK.geqrf!(x) に対応する。 unpivoted QR decomposition
Julia の qr（x） の結果の R factor は qr.R(qr(x)) と同じである。
=====#
using LinearAlgebra, RCall
function qrR(R, complete = false)
mn = minimum(size(R))
!complete && (R = R[1:mn, :])
[R[i, j] = 0 for i = 1:mn, j = 1:mn if i > j]
R
end
★★ pivotted
a1 = reshape([1.,3,2,1,2,3,4,3,2,5,4,3], 4, :);
LAPACK.geqp3!(a1); # pivoted QR decomposition
qrR(a1) # LAPACK.geqp3!() の引数は破壊される
#=====
3×3 Matrix{Float64}:
-7.34847 -5.98764 -3.81032
0.0 1.46566 -0.555939
0.0 0.0 -0.415227
=====#
R"""
a1 = matrix(c(1.0,3,2,1,2,3,4,3,2,5,4,3), 4)
res = qr(a1, LAPACK = TRUE)
print(qr.R(res));
"""
#=====
[,1] [,2] [,3]
[1,] -7.348469 -5.987642 -3.8103174
[2,] 0.000000 1.465656 -0.5559386
[3,] 0.000000 0.000000 -0.4152274
=====#
★★ unpivotted
a1 = reshape([1.,3,2,1,2,3,4,3,2,5,4,3], 4, :);
LAPACK.geqrf!(a1); # unpivoted QR decomposition
qrR(a1) # LAPACK.geqrf!() の引数は破壊される
#=====
3×3 Matrix{Float64}:
-3.87298 -5.68038 -7.22957
0.0 2.39444 1.22506
0.0 0.0 0.482243
=====#
R"""
a1 = matrix(c(1.0,3,2,1,2,3,4,3,2,5,4,3), 4)
res2 = qr(a1, LAPACK = FALSE)
print(qr.R(res2))
"""
#=====
[,1] [,2] [,3]
[1,] -3.872983 -5.680376 -7.2295689
[2,] 0.000000 2.394438 1.2250613
[3,] 0.000000 0.000000 0.4822428
=====#
qr(a1)
#=====
LinearAlgebra.QRCompactWY{Float64, Matrix{Float64}}
Q factor:
4×4 LinearAlgebra.QRCompactWYQ{Float64, Matrix{Float64}}:
-0.258199 0.222738 -0.289346 -0.894427
-0.774597 -0.584688 0.241121 5.55112e-16
-0.516398 0.445477 -0.578691 0.447214
-0.258199 0.640373 0.723364 2.10942e-15
R factor:
3×3 Matrix{Float64}:
-3.87298 -5.68038 -7.22957
0.0 2.39444 1.22506
0.0 0.0 0.482243
=====#
