swp(x, n1 [, n2, ...] [, quiet:T, diag:d, tolerance:tol, keepswept:T]), x a REAL matrix, n1, n2, ... positive integers, or vectors of positive integers, tol > 0 a REAL scalar, d a REAL vector with positive elements |

swp(x,Cols) uses a form of the Beaton SWP operator on the REAL matrix x to compute a REAL result of the same size. If x is m by n, Cols should be a vector whose elements are positive integers <= min(m,n). A single SWP of x on row and column k produces a matrix y of the same size as x with y[i,j] = x[i,j] - x[i,k]x[k,j]/x[k,k], for i and j not equal to k y[i,k] = x[i,k]/x[k,k], for i not equal to k y[k,j] = -x[k,j]/x[k,k], for j not equal to k y[k,k] = 1/x[k,k] The element x[k,k] is sometimes called a "pivot". swp(x,Cols) does successive SWPs of x on Cols[1], Cols[2], ... rows and columns. If, attempting to SWP row and column M, abs(pivot) is found to be too small in comparison with the abs(x[M,M]), the message WARNING: tolerance failure pivoting column M; not swept is printed and that SWP is skipped. The threshold for finding tolerance failure may be modified by keyword 'tolerance'; see below. If x is n by n and non-singular, swp(x,run(n)) computes the inverse of x (it can fail for certain non-positive definite but invertible matrices). swp(x,Cols1, Cols2, ...), where Cols1, Cols2, ..., are vectors of positive integers is equivalent to swp(x,vector(Cols1, Cols2, ...)). For example, swp(cp,1,2,3,4) is equivalent to swp(cp,run(4)). swp() can be very useful in regression and analysis of variance. swp(x,Cols,...,tolerance:Tol), where Tol is a small REAL scalar does the same, but Tol is used in the test for what constitutes a small pivot. Column M will be skipped if abs(pivot) < Tol*abs(x[M,M]), where pivot is the value of x[M,M] just before a SWP of row and column M is attempted. The default value of Tol is 10^(-12). swp(x,Cols,...,quiet:T) does the same, but no warning message is printed when a too small pivot is found. swp(x,Cols,...,diag:d [,quiet:T]), where d is a REAL vector, does the same, but the pivot for row and clumn M is compared with abs(d[M]) to test whether a row and column will be swept. The length of d must always be min(nrows(x), ncols(x)). This feature allows the sequence, say, Cmd> d <- diag(x); x1 <- swp(x,1); x1 <- swp(x1,2,diag:d) to be completely equivalent to Cmd> x1 <- swp(x,1,2) Without diag:d, the comparison element for swp(x1,2) would be x1[2,2] instead of x[2,2]. swp(x,Cols,..., keepswept:T [,tolerance:Tol, ,diag:d, quiet:F]) does the same, except that the result is a structure with components 'matrix' and 'sweptcols'. Component 'matrix' is the swept version of x. Component 'sweptcols' is a vector of integers, with abs(sweptcols[i]) the number of the i-th row and column swept. sweptcols[i] < 0 if and only if row and column abs(sweptcols[i]) failed the tolerance check. No warning message is printed unless 'quiet:F' is an argument. The use of 'keepswept:T' may allow you to compute a generalized inverse of a singular square matrix. Cmd> tmp <- swp(a, run(nrows(a)), keepswept:T) Cmd> b <- tmp$matrix; J <- (-tmp$sweptcols)[tmp$sweptcols < 0] Cmd> if(!isnull(J)){b[J,] <- b[,J] <- 0;;} Now, even if a is singular, a %*% b %*% a should be a within rounding error and b %*% a %*% b should be b within rounding error. See also solve(), qr().

Gary Oehlert 2003-01-15