#METHODS OF MATRIX DECOMPOSITION:

#READ DATA:
M=read.table("c:/DATA/Multivariate/Isetosastand.txt",header=F)
M
attach(M)

#CALCULATION COVARIANCE MATRIX S:
S=cov(M)
S

#EIGENVALUES & EIGENVECTORS OF S:
lambda=eigen(S)$values
lambda
E=eigen(S)$vectors
E

#SPECTRAL DECOMPOSITION:
S
E%*%diag(lambda)%*%t(E)

#DIAGONALIZING S:
diag(lambda)
t(E)%*%S%*%E

#EXPRESSION AS PARTIAL SUMS:
P1=(E[,1]*lambda[1])%*%t(E[,1])
P1
P2=(E[,2]*lambda[2])%*%t(E[,2])
P2
P3=(E[,3]*lambda[3])%*%t(E[,3])
P3
P4=(E[,4]*lambda[4])%*%t(E[,4])
P4

P1+P2+P3+P4

#SINGULAR VALUE DECOMPOSITION:
#SINGULAR VALUES:
gamma=svd(M)$d
gamma
#LEFT MATRIX U:
U=svd(M)$u
U
#RIGHT MATRIX V:
V=svd(M)$v
V
#IDENTITY DEMONSTRATED:
U%*%diag(gamma)%*%t(V)
M

#MATRIX SQUARE ROOT:

#FUNCTION FOR SQUARE ROOT OF A SYMMETRIC POSITIVE DEFINITE MATRIX:
matrix.square.root <- function(M) {
M.eig <- eigen(M)
M.square.root <- M.eig$vectors %*% diag(sqrt(M.eig$values)) %*% solve(M.eig$vectors)
return(M.square.root)
}

Ssqrt=matrix.square.root(S)
Ssqrt

#IDENTITY DEMONSTRATED:
Ssqrt%*%Ssqrt

#MATRIX INVERSE SQUARE ROOT:

#FUNCTION FOR INVERSE SQUARE ROOT OF A SYMMETRIC POSITIVE DEFINITE MATRIX:
matrix.inverse.square.root <- function(M) {
M.eig <- eigen(M)
M.square.root <- M.eig$vectors %*% diag(1/sqrt(M.eig$values)) %*% solve(M.eig$vectors)
return(M.square.root)
}

Sinvsqrt=matrix.inverse.square.root(S)
Sinvsqrt

#MATRIX INVERSE:
Sinv=solve(S)
Sinv