Homework4

To build and visualize eigenimagery that accounts for 80% of the variability

#loading required libraries
library(dplyr)

## 
## Attaching package: 'dplyr'

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

library(jpeg)
library(RSpectra)
library("EBImage")

img_files <- list.files(pattern = "\\.jpg$")
height=100; 
width=240;
matrx = matrix(1, nrow = 3, ncol = 3)
matrx[2, 2] = -8
image_plot = function(path, add=FALSE)
{ jpg_img = readJPEG(path, native=T) # read the file
  resoln = dim(jpg_img)[2:1] # get the resolution, [x, y]
  if (!add) # initialize an empty plot area
    plot(1,1,xlim=c(1,resoln[1]),ylim=c(1,resoln[2]),asp=1,type='n',xaxs='i',yaxs='i',xaxt='n',yaxt='n',xlab='',ylab='',bty='n')
  rasterImage(jpg_img,1,1,resoln[1],resoln[2])
}
img=array(rep(0,length(img_files)*height*width*3), dim=c(length(img_files), height, width,3))

for (i in 1:length(img_files)){

  tmp= EBImage::resize(readJPEG(img_files[i]),100, 240)
  img[i,,,]=array(tmp,dim=c(1, 100, 240,3))
}
flat=matrix(0, length(img_files), prod(dim(img)))
for (i in 1:length(img_files)){
  r=as.vector(img[i,,,1]); g=as.vector(img[i,,,2]);b=as.vector(img[i,,,3])
  flat[i,] <- t(c(r, g, b))
}
shoes=as.data.frame(t(flat))

par(mfrow=c(3,3))
par(mai=c(.3,.3,.3,.3))
for (i in 1:length(img_files)){ 
image_plot(writeJPEG(img[i,,,]))
}

### Compute EigenValues

scaled=scale(shoes, center = TRUE, scale = TRUE)
Sigma_=cor(scaled)
myeigen=eigs(Sigma_,5,which="LM")
cumsum(myeigen$values) / sum(eigen(Sigma_)$values)

## [1] 0.6908865 0.7901734 0.8431089 0.8704645 0.8894482

From the above EigenValues the 80% of the variability has already been surpassed.

scaling=diag(myeigen$values[1:5]^(-1/2)) / (sqrt(nrow(scaled)-1))
eigenshoes=scaled%*%myeigen$vectors[,1:5]%*%scaling

newdata=img
dim(newdata)=c(length(img_files),height*width*3)
mypca=princomp(t(as.matrix(newdata)), scores=TRUE, cor=TRUE)

pcaScores=t(mypca$scores)
dim(pcaScores)=c(length(img_files),height,width,3)
par(mfrow=c(5,5))
par(mai=c(.001,.001,.001,.001))
for (i in 1:length(img_files)){
image_plot(writeJPEG(pcaScores[i,,,], bg="white"))  
}


for (i in 1:length(img_files)){

  tmp= EBImage::gblur(EBImage::resize(readJPEG(img_files[i]),100, 240),sigma = 8, radius=5)
  img[i,,,]=array(tmp,dim=c(1, 100, 240,3))
}
flat=matrix(0, length(img_files), prod(dim(img)))
for (i in 1:length(img_files)){
  r=as.vector(img[i,,,1]); g=as.vector(img[i,,,2]);b=as.vector(img[i,,,3])
  flat[i,] <- t(c(r, g, b))
}
shoes=as.data.frame(t(flat))
par(mfrow=c(3,3))

par(mai=c(.3,.3,.3,.3))
for (i in 1:length(img_files)){ 
image_plot(writeJPEG(img[i,,,]))
}

### Compute new EigenValues

scaled=scale(shoes, center = TRUE, scale = TRUE)
Sigma_=cor(scaled)
myeigen=eigs(Sigma_,5,which="LM")
cumsum(myeigen$values) / sum(eigen(Sigma_)$values)

## [1] 0.7782187 0.8817021 0.9256703 0.9432350 0.9557032

Variability has increased. Blurring images might achieves greater refinement

pcaScores=t(mypca$scores)
dim(pcaScores)=c(length(img_files),height,width,3)
par(mfrow=c(5,5))
par(mai=c(.001,.001,.001,.001))
for (i in 1:length(img_files)){
image_plot(writeJPEG(pcaScores[i,,,], bg="white"))  
}