Assignment Prompt

With the attached data file, build and visualize eigenimagery that accounts for 80% of the variability. Provide full R code and discussion.

Load Libraries

library(foreach)
library(jpeg)
library(EBImage)
## Warning: package 'EBImage' was built under R version 4.2.3

Load Images from Folder

filepath <- "~/Downloads/jpg"
files <- list.files(filepath)
print(files)
##  [1] "RC_2500x1200_2014_us_53446.jpg" "RC_2500x1200_2014_us_53455.jpg"
##  [3] "RC_2500x1200_2014_us_53469.jpg" "RC_2500x1200_2014_us_53626.jpg"
##  [5] "RC_2500x1200_2014_us_53632.jpg" "RC_2500x1200_2014_us_53649.jpg"
##  [7] "RC_2500x1200_2014_us_53655.jpg" "RC_2500x1200_2014_us_53663.jpg"
##  [9] "RC_2500x1200_2014_us_53697.jpg" "RC_2500x1200_2014_us_54018.jpg"
## [11] "RC_2500x1200_2014_us_54067.jpg" "RC_2500x1200_2014_us_54106.jpg"
## [13] "RC_2500x1200_2014_us_54130.jpg" "RC_2500x1200_2014_us_54148.jpg"
## [15] "RC_2500x1200_2014_us_54157.jpg" "RC_2500x1200_2014_us_54165.jpg"
## [17] "RC_2500x1200_2014_us_54172.jpg"
height=1200
width=2500
scale=20
plot_jpeg = function(path, add=FALSE) #initialize function
{
  require('jpeg')
  jpg = readJPEG(path, native=T) # read the file
  res = dim(jpg)[2:1] # get the resolution, [x is 2, y is 1]
  if (!add) # initialize an empty plot area if add==FALSE
    plot(1,1,xlim=c(1,res[1]),ylim=c(1,res[2]), #set the X Limits by size
         asp=1, #aspect ratio
         type='n', #don't plot
         xaxs='i',yaxs='i',#prevents expanding axis windows +6% as normal
         xaxt='n',yaxt='n',xlab='',ylab='', # no axes or labels
         bty='n') # no box around graph
  rasterImage(jpg,1,1,res[1],res[2]) #image, xleft,ybottom,xright,ytop
}

Load Data into Array

im=array(rep(0,length(files)*height/scale*width/scale*3),
         #set dimension to N, x, y, 3 colors, 4D array)
         dim=c(length(files), height/scale, width/scale,3)) 

for (i in 1:length(files)){
  #define file to be read
  tmp=paste0("~/Downloads/jpg/", files[i])
  #read the file
  temp=EBImage::resize(readJPEG(tmp),height/scale, width/scale)
  #assign to the array
  im[i,,,]=array(temp,dim=c(1, height/scale, width/scale,3))
}

Actual Plot

par(mfrow=c(3,3)) #set graphics to 3 x 3 table
par(mai=c(.3,.3,.3,.3)) #set margins 
for (i in 1:17){
  plot_jpeg(writeJPEG(im[i,,,]))
}

height=1200
width=2500
scale=20
newdata=im
dim(newdata)=c(length(files),height*width*3/scale^2)
mypca=princomp(t(as.matrix(newdata)), scores=TRUE, cor=TRUE)
mypca
## Call:
## princomp(x = t(as.matrix(newdata)), cor = TRUE, scores = TRUE)
## 
## Standard deviations:
##    Comp.1    Comp.2    Comp.3    Comp.4    Comp.5    Comp.6    Comp.7    Comp.8 
## 3.4319009 1.3118000 0.9316975 0.6809679 0.5683219 0.5255886 0.4877556 0.4528049 
##    Comp.9   Comp.10   Comp.11   Comp.12   Comp.13   Comp.14   Comp.15   Comp.16 
## 0.4060420 0.3924175 0.3791956 0.3684830 0.3602187 0.3374253 0.3241916 0.3181866 
##   Comp.17 
## 0.2808942 
## 
##  17  variables and  22500 observations.
sum(mypca$sdev^2/sum(mypca$sdev^2))
## [1] 1
mycomponents=mypca$sdev^2/sum(mypca$sdev^2)
sum(mycomponents[1:17])
## [1] 1
sum(mycomponents[1:17]) 
## [1] 1

Calculate the cumulative proportion of variance explained

cumulative_variance <- cumsum(mypca$sdev^2) / sum(mypca$sdev^2)
cumulative_variance
##    Comp.1    Comp.2    Comp.3    Comp.4    Comp.5    Comp.6    Comp.7    Comp.8 
## 0.6928202 0.7940449 0.8451073 0.8723847 0.8913841 0.9076338 0.9216282 0.9336889 
##    Comp.9   Comp.10   Comp.11   Comp.12   Comp.13   Comp.14   Comp.15   Comp.16 
## 0.9433872 0.9524455 0.9609037 0.9688907 0.9765235 0.9832209 0.9894033 0.9953587 
##   Comp.17 
## 1.0000000

Eigenshoes

mypca2=t(mypca$scores)
dim(mypca2)=c(length(files),height/scale,width/scale,3)
par(mfrow=c(5,5))
par(mai=c(.001,.001,.001,.001))

for (i in 1:17){
  plot_jpeg(writeJPEG(mypca2[i,,,], quality=1,bg="white"))
}