if (FALSE) {
if(require("mvoutlier") && require("RColorBrewer"))
{
data(bsstop)
Data.1 <- bsstop[, 1:14]
colnames(Data.1)
Data.1.scaled <- scale(as.matrix(Data.1[5:14])) # standardised data...
rownames(Data.1.scaled) <- Data.1[, 1]
#compute principal components:
pca <- princomp(Data.1.scaled, cor = FALSE, scores = TRUE)
# use covariance matrix to match the following...
pca$loadings
data(bss.background)
backdrop <- function()
plot(bss.background, asp = 1, type = "l", xaxt = "n", yaxt = "n",
xlab = "", ylab = "", bty = "n", col = "grey")
pc1 <- pca$scores[, 1]
backdrop()
points(Data.1$XCOO[pc1 > 0], Data.1$YCOO[pc1 > 0], pch = 16, col = "blue")
points(Data.1$XCOO[pc1 < 0], Data.1$YCOO[pc1 < 0], pch = 16, col = "red")
#Geographically Weighted PCA and mapping the local loadings
# Coordinates of the sites
Coords1 <- as.matrix(cbind(Data.1$XCOO,Data.1$YCOO))
d1s <- SpatialPointsDataFrame(Coords1,as.data.frame(Data.1.scaled))
pca.gw <- gwpca(d1s,vars=colnames(d1s@data),bw=1000000,k=10)
local.loadings <- pca.gw$loadings[, , 1]
# Mapping the winning variable with the highest absolute loading
# note first component only - would need to explore all components..
lead.item <- colnames(local.loadings)[max.col(abs(local.loadings))]
df1p = SpatialPointsDataFrame(Coords1, data.frame(lead = lead.item))
backdrop()
colour <- brewer.pal(8, "Dark2")[match(df1p$lead, unique(df1p$lead))]
plot(df1p, pch = 18, col = colour, add = TRUE)
legend("topleft", as.character(unique(df1p$lead)), pch = 18, col =
brewer.pal(8, "Dark2"))
backdrop()
#Glyph plots give a view of all the local loadings together
glyph.plot(local.loadings, Coords1, add = TRUE)
#it is not immediately clear how to interpret the glyphs fully,
#so inter-actively identify the full loading information using:
check.components(local.loadings, Coords1)
# GWPCA with an optimal bandwidth
bw.choice <- bw.gwpca(d1s,vars=colnames(d1s@data),k=2)
pca.gw.auto <- gwpca(d1s,vars=colnames(d1s@data),bw=bw.choice,k=2)
# note first component only - would need to explore all components..
local.loadings <- pca.gw.auto$loadings[, , 1]
lead.item <- colnames(local.loadings)[max.col(abs(local.loadings))]
df1p = SpatialPointsDataFrame(Coords1, data.frame(lead = lead.item))
backdrop()
colour <- brewer.pal(8, "Dark2")[match(df1p$lead, unique(df1p$lead))]
plot(df1p, pch = 18, col = colour, add = TRUE)
legend("topleft", as.character(unique(df1p$lead)), pch = 18,
col = brewer.pal(8, "Dark2"))
# GWPCPLOT for investigating the raw multivariate data
gw.pcplot(d1s, vars=colnames(d1s@data),focus=359, bw = bw.choice)
}
}
Run the code above in your browser using DataLab