Language Reference

INVUPDT Function

INVUPDT (matrix, vector<, scalar>);

The INVUPDT function updates a matrix inverse.

The arguments to the INVUPDT function are as follows:

matrix: is an $n \times n$ nonsingular matrix. In most applications matrix is symmetric positive definite.
vector: is an $n \times 1$ or $1 \times n$ vector.
scalar: is a numeric scalar.

The Sherman-Morrison-Woodbury formula is

$(\mb{A} + \mb{U}\mb{V}^{\prime })^{-1} = \mb{A}^{-1} - \mb{A}^{-1} \mb{U} (\mb{I} + \mb{V}^{\prime }\mb{A}^{-1}\mb{U})^{-1} \mb{V}^{\prime }\mb{A}^{-1}$

where $\mb{A}$ is an $n \times n$ nonsingular matrix and $\mb{U}$ and $\mb{V}$ are $n\times k$ . The formula shows that a rank k update to $\mb{A}$ corresponds to a rank k update of $\mb{A}^{-1}$ .

The INVUPDT function implements the Sherman-Morrison-Woodbury formula for rank-one updates with $\mb{U} = w \mb{X}$ and $\mb{V} = \mb{X}$ , where $\mb{X}$ is an $n \times 1$ vector and w is a scalar.

If $\mb{M} = \mb{A}^{-1}$ , then you can call the INVUPDT function as follows:

R = invupdt(M, X, w);

This statement computes the following matrix:

$\mb{R} = \mb{M} - w\mb{MX} (\mb{I} + w\mb{X}^{\prime }\mb{MX})^{-1} \mb{X}^{\prime }\mb{M}$

The matrix $\mb{R}$ is equivalent to $(\mb{A} + w\mb{XX}^{\prime })^{-1}$ . If $\mb{A}$ is symmetric positive definite, then so is $\mb{R}$ .

If w is not specified, then it is given a default value of 1.

A common use of the INVUPDT function is in linear regression. If $\mb{Z}$ is a design matrix, $\mb{M} = (\mb{Z}^{\prime }\mb{Z})^{-1}$ is the associated inverse crossproduct matrix, and $\mb{v}$ is a new observation to be used in estimating the parameters of a linear model, then the inverse crossproducts matrix that includes the new observation can be updated from $\mb{M}$ by using the following statement:

M2 = invupdt(M, v);

If w is 1, the function adds an observation to the inverse; if w is $-1$ , the function removes an observation from the inverse. If weighting is used, w is the weight.

To perform the computation, the INVUPDT function uses about $2n^2$ multiplications and additions, where n is the row dimension of the positive definite argument matrix.

The following program demonstrates adding or removing observations from a linear fit and updating the inverse crossproduct matrix:

X = {0, 1, 1, 1, 2, 2, 3, 4, 4};
Y = {1, 1, 2, 6, 2, 3, 3, 3, 4};

/* find linear fit */
Z = j(nrow(X), 1, 1) || X;           /* design matrix */
M = inv(Z`*Z);

b = M*Z`*Y;                          /* LS estimate */
resid = Y - Z*b;                     /* residuals */
print "Original Fit", b resid;

/* residual for observation (1,6) seems too large.
   Take obs number 4 out of data set and refit. */
v = z[4,];
M = invupdt(M, v, -1);               /* update inverse crossprod */

keepObs = (1:3) || (5:nrow(X));
Z = Z[keepObs, ];
Y = Y[keepObs, ];
b = M*Z`*Y;                          /* new LS estimate */
print "After deleting observation 4", b;

/* Add a new obs (x,y) = (0,2) and refit. */
obs = {0 2};
v = 1 || obs[1];                     /* new row in design matrix */
M = invupdt(M, v);

Z = Z // v;
Y = Y // obs[2];
b = M*Z`*Y;                          /* new LS estimate */
print "After adding observation (0,2)", b;

Figure 24.171: Refitting Linear Regression Models

Original Fit

b	resid
2.0277778	-1.027778
0.375	-1.402778
	-0.402778
	3.5972222
	-0.777778
	0.2222222
	-0.152778
	-0.527778
	0.4722222

After deleting observation 4

b
1
0.6470588

After adding observation (0,2)

b
1.3
0.5470588