# svd2.r                 September 2007, revised 2008, 2009
 #
 # some examples of use of svd to do Singular Value Decomposition
 #
 # second example -- another singular design matrix
 #
  a <- 3            # number of levels of factor A (alpha's)
  b <- 4            # number of levels of factor B (beta's)
  N <- a*b          # number of observations
  One <- matrix(1,N,1)        # intercept column
  XA  <- matrix(c(rep(1,b),rep(0,N)),N,a)        # factor A cols
  XB  <- matrix(rep(diag(1,b),a),N,b,byrow=TRUE) # factor B cols
  X <- cbind(One,XA,XB)
 #                    write out design matrix
  X
 #                    next svd
  svd(X)
 #                    inner product matrix
  xpx <- t(X) %*% X
  print("X'X")
  xpx
  expx <- eigen(xpx,symmetric=TRUE) # eigenvectors and values of X'X
  sqrt(expx$values)
 # done
  rm(list=ls())
 q()