fastQR (version 1.0.0)

qrupdate: Fast updating of the QR factorization

Description

qrupdate provides the update of the QR factorization after the addition of \(m>1\) rows or columns to the matrix \(X\in\mathbb{R}^{n\times p}\) with \(n>p\). The QR factorization of the matrix \(X\) returns the matrices \(Q\in\mathbb{R}^{n\times n}\) and \(R\in\mathbb{R}^{n\times p}\) such that \(X=QR\). The \(Q\) and \(R\) matrices are factorized as \(Q=\begin{bmatrix}Q_1&Q_2\end{bmatrix}\) and \(R=\begin{bmatrix}R_1\\R_2\end{bmatrix}\), with \(Q_1\in\mathbb{R}^{n\times p}\), \(Q_2\in\mathbb{R}^{n\times (n-p)}\) such that \(Q_1^{\top}Q_2=Q_2^\top Q_1=0\) and \(R_1\in\mathbb{R}^{p\times p}\) upper triangular matrix and \(R_2\in\mathbb{R}^{(n-p)\times p}\). qrupdate accepts in input the matrices \(Q\) and either the complete matrix \(R\) or the reduced one, \(R_1\). See Golub and Van Loan (2013) for further details on the method.

Usage

qrupdate(Q, R, k, U, type = NULL, fast = NULL, complete = NULL)

Value

A named list containing

Q

the updated Q matrix.

R

the updated R matrix.

Arguments

Q

a \(n\times p\) matrix.

R

a \(p\times p\) upper triangular matrix.

k

position where the columns or the rows are added.

U

either a \(n\times m\) matrix or a \(p\times m\) matrix of columns or rows to be added.

type

either 'row' of 'column', for adding rows or columns. Default is 'column'.

fast

fast mode: disable to check whether the provided matrices are valid inputs. Default is FALSE.

complete

logical expression of length 1. Indicates whether an arbitrary orthogonal completion of the \(Q\) matrix is to be made, or whether the \(R\) matrix is to be completed by binding zero-value rows beneath the square upper triangle.

References

golub_van_loan.2013fastQR

bjorck.2015fastQR

bjorck.2024fastQR

bernardi_etal.2024fastQR

Examples

## Add one column
## generate sample data
set.seed(1234)
n <- 12
p <- 5
X <- matrix(rnorm(n * p, 1), n, p)

## get the initial QR factorization
output <- qr(X, complete = TRUE)
Q      <- output$Q
R      <- output$R

## create column u to be added
k  <- p+1
u  <- matrix(rnorm(n), n, 1)
X1 <- cbind(X, u)

## update the QR decomposition
out <- fastQR::qrupdate(Q = Q, R = R,
                       k = k, U = u,
                       type = "column",
                       fast = FALSE,
                       complete = TRUE)

## check
round(out$Q %*% out$R - X1, 5)
max(abs(out$Q %*% out$R - X1))

## Add m columns
## create data: n > p
set.seed(1234)
n <- 10
p <- 5
X <- matrix(rnorm(n * p, 1), n, p)

## get the initial QR factorization
output <- fastQR::qr(X, complete = TRUE)
Q      <- output$Q
R      <- output$R

## create the matrix of two columns to be added
## in position 2
k  <- 2
m  <- 2
U  <- matrix(rnorm(n*m), n, m)
X1 <- cbind(X[,1:(k-1)], U, X[,k:p])

# update the QR decomposition
out <- fastQR::qrupdate(Q = Q, R = R,
                        k = k, U = U, type = "column",
                       fast = FALSE, complete = TRUE)

## check
round(out$Q %*% out$R - X1, 5)
max(abs(out$Q %*% out$R - X1))

## Add one row
## create data: n > p
set.seed(1234)
n <- 12
p <- 5
X <- matrix(rnorm(n * p, 1), n, p)

## get the initial QR factorization
output <- fastQR::qr(X, complete = TRUE)
Q      <- output$Q
R      <- output$R
R1     <- R[1:p,]

## create the row u to be added
u  <- matrix(data = rnorm(p), p, 1)
k  <- n+1
if (k<=n) {
  X1 <- rbind(rbind(X[1:(k-1), ], t(u)), X[k:n, ])
} else {
  X1 <- rbind(rbind(X, t(u)))
}

## update the QR decomposition
out <- fastQR::qrupdate(Q = Q, R = R,
                        k = k, U = u,
                        type = "row",
                        complete = TRUE)

## check
round(out$Q %*% out$R - X1, 5)
max(abs(out$Q %*% out$R - X1))

## Add m rows
## create data: n > p
set.seed(1234)
n <- 12
p <- 5
X <- matrix(rnorm(n * p, 1), n, p)

## get the initial QR factorization
output <- fastQR::qr(X, complete = TRUE)
Q      <- output$Q
R      <- output$R
R1     <- R[1:p,]

## create the matrix of rows U to be added:
## two rows in position 5
m  <- 2
U  <- matrix(data = rnorm(p*m), p, m)
k  <- 5
if (k<=n) {
  X1 <- rbind(rbind(X[1:(k-1), ], t(U)), X[k:n, ])
} else {
  X1 <- rbind(rbind(X, t(U)))
}

## update the QR decomposition
out <- fastQR::qrupdate(Q = Q, R = R,
                        k = k, U = U,
                        type = "row",
                        complete = FALSE)

## check
round(out$Q %*% out$R - X1, 5)
max(abs(out$Q %*% out$R - X1))