This page provides a concise reference and a sequence of practice problems for matrix calculations in base R. It is suitable for review or for use alongside lectures on linear algebra, regression, and numerical computation.
Tip: Use Ctrl/Cmd + Enter in RStudio to run selected lines. Use set.seed() to obtain reproducible random matrices.
# Create matrices
A <- matrix(1:9, nrow = 3, ncol = 3) # fills by column by default
B <- matrix(seq(2, 20, by = 2), nrow = 3) # another 3x3
C <- matrix(rnorm(12), nrow = 3, ncol = 4) # 3x4 random
# Bind and diag
X <- cbind(1, A) # add intercept column
Y <- rbind(A, 10:12)
I3 <- diag(3) # 3x3 identity
# Elementwise vs matrix multiplication
A + B # addition
A * B # elementwise product
A %*% B # matrix product
# Transpose, inverse, determinant
At <- t(A)
if (det(A) != 0) Ainv <- solve(A)
# Norms and rank
norm_A <- norm(A, type = "F")
rank_A <- qr(A)$rank
A[1, 2] # element (row 1, col 2)
A[ , 2] # column 2 (as vector)
A[1:2, 2:3] # submatrix
A[diag(3) == 1] <- 5 # set diagonal to 5
rowSums(A); colMeans(A)
apply(A, 1, sd) # row-wise sd
# Linear system Ax = b
set.seed(1)
A <- matrix(rnorm(9), 3, 3)
b <- rnorm(3)
x <- solve(A, b) # solves Ax = b
# Eigen and SVD
E <- eigen(A) # E$values, E$vectors
sv <- svd(A) # U, D, V
# Least squares (normal equations)
X <- cbind(1, matrix(rnorm(20), 10, 2))
y <- rnorm(10)
beta_hat <- solve(t(X) %*% X, t(X) %*% y)
matrix(data, nrow, ncol, byrow = FALSE) — build a matrix.cbind(...), rbind(...) — column/row bind.diag(n) — identity; diag(v) — diagonal with entries v.t(A) — transpose; solve(A) — inverse or solve system; det(A) — determinant.%*% — matrix product; * — elementwise product.rowSums/colSums, rowMeans/colMeans, apply(A, MARGIN, FUN).eigen(A), svd(A), qr(A)$rank, crossprod(X) (= t(X)%*%X), tcrossprod(X).Work through these in order. Each problem includes an optional solution. Unless stated otherwise, assume all matrices are conformable and invertible where needed.
Create two 3×3 matrices A and B. Compute A + B, A * B (elementwise), and A %*% B. Check whether A %*% B equals B %*% A.
set.seed(42)
A <- matrix(sample(1:9), 3, 3)
B <- matrix(sample(1:9), 3, 3)
A + B
A * B
A %*% B
all.equal(A %*% B, B %*% A) # generally FALSE
Compute t(A), det(A), and (if invertible) solve(A). Verify numerically that A %*% solve(A) is the identity.
At <- t(A)
DA <- det(A)
if (abs(DA) > 1e-12) {
Ainv <- solve(A)
I3 <- A %*% Ainv
round(I3, 8)
}
Replace the diagonal entries of A with 5. Extract the submatrix formed by rows 1–2 and columns 2–3.
diag(A) <- 5
subA <- A[1:2, 2:3]
A; subA
Compute rowSums(A), colMeans(A), and the row‑wise standard deviation using apply.
rowSums(A)
colMeans(A)
apply(A, 1, sd)
Generate a random invertible 3×3 matrix A and vector b. Solve Ax = b for x.
set.seed(1)
A <- matrix(rnorm(9), 3, 3); b <- rnorm(3)
x <- solve(A, b)
# check
max(abs(A %*% x - b))
Simulate data with two predictors (plus intercept) and solve for \(\hat{\beta}\) using (X'X)^{-1} X'y.
set.seed(7)
X <- cbind(1, matrix(rnorm(200), 100, 2))
y <- 3 + 2*X[,2] - 1*X[,3] + rnorm(100, sd=0.5)
beta_hat <- solve(t(X) %*% X, t(X) %*% y)
round(beta_hat, 3)
Compute the eigenvalues and eigenvectors of a symmetric positive definite matrix and verify the reconstruction A = V %*% D %*% t(V).
set.seed(3)
M <- matrix(rnorm(16), 4, 4)
A <- crossprod(M) # SPD
E <- eigen(A)
V <- E$vectors; D <- diag(E$values)
max(abs(A - V %*% D %*% t(V)))
Construct a nearly collinear design matrix and compare solve vs qr.solve for least squares.
set.seed(4)
x1 <- rnorm(100)
x2 <- x1 + rnorm(100, sd=1e-4) # nearly collinear
X <- cbind(1, x1, x2)
y <- 1 + 2*x1 - 2*x2 + rnorm(100, sd=0.2)
# Normal equations (may be unstable)
beta_ne <- solve(t(X)%*%X, t(X)%*%y)
# QR-based solution (preferred)
beta_qr <- qr.solve(X, y)
qr(X)$rank
cbind(beta_ne, beta_qr)
Build a 4×4 block matrix \(\begin{bmatrix}I_2 & A\\ 0 & B\end{bmatrix}\) for given 2×2 matrices A, B. Verify its determinant equals det(A) * det(B).
I2 <- diag(2)
A <- matrix(c(2,1,0,3), 2, 2)
B <- matrix(c(4,2,1,2), 2, 2)
Z <- matrix(0, 2, 2)
Blk <- rbind(cbind(I2, A), cbind(Z, B))
all.equal(det(Blk), det(A) * det(B))
Show that crossprod(X) equals t(X) %*% X and use it to compute a sample covariance matrix.
set.seed(6)
X <- matrix(rnorm(300), 100, 3)
all.equal(crossprod(X), t(X) %*% X)
# sample covariance (unbiased)
Xc <- scale(X, center = TRUE, scale = FALSE)
S <- crossprod(Xc) / (nrow(X) - 1)
S
Let u be a non‑zero vector. Form the projection matrix P = u %*% t(u) / as.numeric(t(u) %*% u) and verify P %*% P = P.
u <- c(1,2,3)
P <- (u %*% t(u)) / as.numeric(t(u) %*% u)
max(abs(P %*% P - P))
Use svd to reconstruct a matrix and verify the maximum absolute reconstruction error.
set.seed(9)
A <- matrix(rnorm(25), 5, 5)
SV <- svd(A)
Ahat <- SV$u %*% diag(SV$d) %*% t(SV$v)
max(abs(A - Ahat))
?matrix, ?solve, ?eigen, ?svd, ?qr, ?crossprod