Set directory, read in libraries and source functions:
rm(list = ls())
# Set working directory
setwd("~/Documents/University of Colorado/Fall 2013/CVEN6833/HW2")
# Source functions and libraries
library(fields)
library(maps)
source("SSTxygrid.r")
source("FullSVD.r")
# Define n's
Nrows <- 72 # 72 longitude points on gridded SST Kaplan data
Ncols <- 36 # 36 latitude points on gridded SST Kaplan data
Ntime <- 113 # Monthly SST from Jan 1900 to Dec 2012
Read in the Kaplan SST data. Note, I created my own function (“SSTxygrid.r”) to create the grid, determine undefined points, and set up SST data for analysis.
# Read in data
index1 <- scan("http://civil.colorado.edu/~emgi6568/index1.txt")
KapSST <- readBin("http://iridl.ldeo.columbia.edu/expert/SOURCES/.KAPLAN/.EXTENDED/.v2/.ssta/T/%28Jan%201900%29%28Dec%202012%29RANGE/T/12/boxAverage/T/12/STEP/data.r4",
what = "numeric", n = (Nrows * Ncols * Ntime), size = 4, endian = "swap")
results <- SSTxygrid(KapSST)
SSTdata <- results$SSTdata
xygrid <- results$grid
Perform PCA on the global annual SST anomalies. Note, I created my own function (“FullSVD.r”) which takes a data set, computes PCA, plots EVS, plots first 4 temporal PCs (optional) and plots first 4 spatial PCs (optional).
# Perform PCA on global annual average SSTs
zPCA <- FullSVD(SSTdata, Nlamb = 40, temporal = TRUE, spatial = "SST", xygrid = xygrid)
Figure 1.
The EVS shows the first 40 of 1207 PCs. We see that the first 10 explain the majority of the variance. We will look at the first four.
Figure 2.
Plots of the temporal PCs show (1) a trend, (2) potential ENSO signal, and (3) potential PDO signal. We can better see whether these oscillations are prevalent by looking at the spatial PC plots
Figure 3.
Plots of the spatial PCs clearly show the ENSO signal in the 2nd PC. Notice that the 3rd PC still has some ENSO signal in it.
Since the 2nd mode is ENSO, we will remove it and repeat PCA.
# Remove the 2nd mode because it is ENSO and repeat analysis
Nmodes <- length(zPCA$u[1, ])
Nkeep <- c(2)
E <- matrix(0, nrow = Nmodes, ncol = Nmodes)
E[, Nkeep] <- zPCA$u[, Nkeep]
SSTannavgkeep <- zPCA$pcs %*% t(E)
SSTannrem <- SSTdata - SSTannavgkeep
zPCA <- FullSVD(SSTannrem, Nlamb = 40, temporal = TRUE, spatial = "SST", xygrid = xygrid)
Figure 4.
Now that we have removed the ENSO signal, the variance explained by each PC is shifted.
Figure 5.
Plots of the temporal PCs still show the first PC with the trend, but the other PCs have now moved up, unchanged, with the ENSO signal removed.
Figure 6.
Plots of the spatial PCs clearly show the ENSO signal removed and all other PCs moved up in order.
Repeat PCA analysis for three sub-periods to see trends in PCs.
SSTdata.p1 <- SSTdata[1:36, ]
SSTdata.p2 <- SSTdata[37:76, ]
SSTdata.p3 <- SSTdata[77:113, ]
zPCA <- FullSVD(SSTdata.p1, Nlamb = 40, temporal = TRUE, spatial = "SST", xygrid = xygrid)
zPCA <- FullSVD(SSTdata.p2, Nlamb = 40, temporal = TRUE, spatial = "SST", xygrid = xygrid)
zPCA <- FullSVD(SSTdata.p3, Nlamb = 40, temporal = TRUE, spatial = "SST", xygrid = xygrid)
Comments on sub-period PCA analysis:
The significance of performing PCA on sub periods is to see the shift in principle components over time. In this exercise, we can clearly see that there is a shift in the magnitudes, temporal trend, and spatial pattern of PCs over time. While the sign (negative or positive) is meaningless, a change in the magnitude, periodicity, or trend (pos or neg) is meaningful and we see that through the three time periods. Specifically, ENSO shows up more in the first PC when the subperiods are looked at. Perhaps this is because the long-term trend cannot be detected in shorter time periods.