R/PCA.R

In sq77/PCA7:

#' A PCA Function
#'
#' This function allows you to perfrom singular value decomposition.
#' @keywords PCA
#' @examples
#' PCA()

PCA <- function(C,variants){
  # Need call from previous function instead of re-extract from matrix C.
  NR = dim(C)[1]
  NC = dim(C)[2]
  STEST = 1e-8
  # Initialized matrix PR, PC, and vector EV.
  # Determine the limit for the number of components.
  NXL = min(NR-1, NC-1,7)
  PR = matrix(0, nrow = NR, ncol = NXL)
  PC = matrix(0, nrow =NXL, ncol = NC)
  EV = rep(0, NXL)
  # Store temporary information NX ITER STEST SFT
  temp = matrix(0, ncol = 4, nrow = NXL)
  colnames(temp) = c("NX","ITER","STEST","SFT")

  if(variants == 7){
    CA.P = ca(C)
    PR = CA.P$rowcoord[,1:NXL]
    PC = t(CA.P$colcoord[,1:NXL])
    EV = CA.P$sv[1:NXL]
  }
  else{
    for(NX in 1:NXL){
      SFT = 1
      ITER = 0
      # Initialize VR with uniform random function
      VR = runif(NR, min = -0.5, max = 0.5)
      SR = VR
      while(ITER<200 & SFT>STEST){
        ITER = ITER + 1
        # Compute VC from VR
        VC = t(VR)%*%C
        EVP = sum(VC^2)
        SVP = sqrt(EVP)
        VC = VC/SVP
        # Recalculate VR from VC by reverse sum
        VR = C%*%t(VC)
        EVP = sum(VR^2)
        SVP = sqrt(EVP)
        VR = VR/SVP
        SFT = max(abs(VR-SR))
        # Retain VR as SR
        SR = VR
      }
      # Store eigenvalue first
      EV[NX] = EVP

      SQVP = sqrt(SVP)
      VR = VR*SQVP
      VC = VC*SQVP

      PR[,NX] = VR
      PC[NX,] = VC

      # Write NX, ITER, STEST, SFT into a matrix
      temp[NX,] = c(NX, ITER, STEST, SFT)
      # Remove NX from data copy
      C = C - VR%*%VC
    }
  }
  return(list(PC = PC,PR = PR, EV = EV))
}