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

$(\bA + \bU \bV ^{\prime })^{-1} = \bA ^{-1} - \bA ^{-1} \bU (\bI + \bV ^{\prime }\bA ^{-1}\bU )^{-1} \bV ^{\prime }\bA ^{-1}$

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

The INVUPDT function implements the Sherman-Morrison-Woodbury formula for rank-one updates with $\bU = w \bX$ and $\bV = \bX$ , where $\bX$ is an $n \times 1$ vector and w is a scalar.

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

R = invupdt(M, X, w);

This statement computes the following matrix:

$\bR = \bM - w\mb{MX} (\bI + w\bX ^{\prime }\mb{MX})^{-1} \bX ^{\prime }\bM$

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

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, $\bM = (\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 $\bM$ 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 25.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