The hardware and bandwidth for this mirror is donated by dogado GmbH, the Webhosting and Full Service-Cloud Provider. Check out our Wordpress Tutorial.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]dogado.de.

Generalized inverse

Michael Friendly

2022-12-08

In matrix algebra, the inverse of a matrix is defined only for square matrices, and if a matrix is singular, it does not have an inverse.

The generalized inverse (or pseudoinverse) is an extension of the idea of a matrix inverse, which has some but not all the properties of an ordinary inverse.

A common use of the pseudoinverse is to compute a ‘best fit’ (least squares) solution to a system of linear equations that lacks a unique solution.

library(matlib)

Construct a square, singular matrix [See: Timm, EX. 1.7.3]

A <-matrix(c(4, 4, -2,
             4, 4, -2,
            -2, -2, 10), nrow=3, ncol=3, byrow=TRUE)
det(A)
## [1] 0

The rank is 2, so inv(A) won’t work

R(A)
## [1] 2

In the echelon form, this rank deficiency appears as the final row of zeros

echelon(A)
##      [,1] [,2] [,3]
## [1,]    1    1    0
## [2,]    0    0    1
## [3,]    0    0    0

inv() will throw an error

try(inv(A))
## Error in Inverse(X, tol = sqrt(.Machine$double.eps), ...) : 
##   X is numerically singular

A generalized inverse does exist for any matrix, but unlike the ordinary inverse, the generalized inverse is not unique, in the sense that there are various ways of defining a generalized inverse with various inverse-like properties. The function matlib::Ginv() calculates a Moore-Penrose generalized inverse.

(AI <- Ginv(A))
##         [,1] [,2]    [,3]
## [1,] 0.27778    0 0.05556
## [2,] 0.00000    0 0.00000
## [3,] 0.05556    0 0.11111

We can also view this as fractions:

Ginv(A, fractions=TRUE)
##      [,1] [,2] [,3]
## [1,] 5/18    0 1/18
## [2,]    0    0    0
## [3,] 1/18    0  1/9

Properties of generalized inverse (Moore-Penrose inverse)

The generalized inverse is defined as the matrix \(A^-\) such that

A %*% AI %*% A
##      [,1] [,2] [,3]
## [1,]    4    4   -2
## [2,]    4    4   -2
## [3,]   -2   -2   10
AI %*% A %*% AI
##         [,1] [,2]    [,3]
## [1,] 0.27778    0 0.05556
## [2,] 0.00000    0 0.00000
## [3,] 0.05556    0 0.11111

In addition, both \(A * A^-\) and \(A^- * A\) are symmetric, but neither product gives an identity matrix, A %*% AI != AI %*% A != I

zapsmall(A %*% AI)
##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    1    0    0
## [3,]    0    0    1
zapsmall(AI %*% A)
##      [,1] [,2] [,3]
## [1,]    1    1    0
## [2,]    0    0    0
## [3,]    0    0    1

Rectangular matrices

For a rectangular matrix, \(A^- = (A^{T} A)^{-1} A^{T}\) is the generalized inverse of \(A\) if \((A^{T} A)^-\) is the ginv of \((A^{T} A)\) [See: Timm: EX 1.6.11]

A <- cbind( 1, matrix(c(1, 0, 1, 0, 0, 1, 0, 1), nrow=4, byrow=TRUE))
A
##      [,1] [,2] [,3]
## [1,]    1    1    0
## [2,]    1    1    0
## [3,]    1    0    1
## [4,]    1    0    1

This \(4 \times 3\) matrix is not of full rank, because columns 2 and 3 sum to column 1.

R(A)
## [1] 2
(AA <- t(A) %*% A)
##      [,1] [,2] [,3]
## [1,]    4    2    2
## [2,]    2    2    0
## [3,]    2    0    2
(AAI <- Ginv(AA))
##      [,1] [,2] [,3]
## [1,]  0.5 -0.5    0
## [2,] -0.5  1.0    0
## [3,]  0.0  0.0    0

The generalized inverse of \(A\) is \((A^{T} A)^- A^{T}\), AAI * t(A)

AI <- AAI  %*%  t(A)

Show that it is a generalized inverse:

A %*% AI %*% A
##      [,1] [,2] [,3]
## [1,]    1    1    0
## [2,]    1    1    0
## [3,]    1    0    1
## [4,]    1    0    1
AI %*% A %*% AI
##      [,1] [,2] [,3] [,4]
## [1,]  0.0  0.0  0.5  0.5
## [2,]  0.5  0.5 -0.5 -0.5
## [3,]  0.0  0.0  0.0  0.0

These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.
Health stats visible at Monitor.