set.seed(229)
library(ggplot2)
library(gcookbook)
library(corrplot)
library(igraph)
library(rgl)
library(grid)
library(vcd)
library(MASS)
library(maps)
library(plyr)
library(maptools)

Miscellaneous Graphs

Note: Some problems cannot be completed as written because the code has changed. You will need to find the updated code.

Section 13.1. Making a Correlation Matrix

# Figure 13-1
mcor <- cor(mtcars)
corrplot(mcor)

# Figure 13-2
corrplot(mcor, method="shade", shade.col=NA, tl.col="black", tl.srt=45)

# Figure 13-3
col <- colorRampPalette(c("#BB4444", "#EE9988", "#FFFFFF", "#77AADD", "#4477AA"))
corrplot(mcor, method = "shade", shade.col = NA, tl.col = "black", tl.srt = 45,
 col = col(200), addCoef.col = "black", addcolorlabel = "no", order = "AOE")
## Warning in text.default(pos.xlabel[, 1], pos.xlabel[, 2], newcolnames, srt
## = tl.srt, : "addcolorlabel" is not a graphical parameter
## Warning in text.default(pos.ylabel[, 1], pos.ylabel[, 2], newrownames, col
## = tl.col, : "addcolorlabel" is not a graphical parameter
## Warning in title(title, ...): "addcolorlabel" is not a graphical parameter

Section 13.2. Plotting a Function

# Figure 13-4
p <- ggplot(data.frame(x = c(-3,3)), aes(x))
p + stat_function(fun = dnorm)

# Figure 13-5
myfun <- function(xvar) {
 1/(1 + exp(-xvar + 10))
}
ggplot(data.frame(x = c(0, 20)), aes(x)) + stat_function(fun = myfun)

Section 13.3. Shading a Subregion Under a Function Curve

# Figure 13-6
dnorm_limit <- function(x) {
 y <- dnorm(x)
 y[x < 0 | x > 2] <- NA
 return(y)
}

p <- ggplot(data.frame(x = c(-3, 3)), aes(x))
p + stat_function(fun = dnorm_limit, geom = "area", fill = "blue", alpha = 0.2) +
 stat_function(fun = dnorm)

Section 13.4. Creating a Network Graph

# Figure 13-7
gd <- graph(c(1,2, 2,3, 2,4, 1,4, 5,5, 3,6))
plot(gd)

gu <- graph(c(1,2, 2,3, 2,4, 1,4, 5,5, 3,6), directed = FALSE)
plot(gu, vertex.label = NA)

# Figure 13-8
g <- graph.data.frame(madmen2, directed = TRUE)
par(mar = c(0, 0, 0, 0))
plot(g, layout = layout.fruchterman.reingold, vertex.size = 8, edge.arrow.size = 0.5,
 vertex.label = NA)

# Figure 13-9
g <- graph.data.frame(madmen, directed = FALSE)
par(mar = c(0, 0, 0, 0))
plot(g, layout = layout.circle, vertex.size = 8, vertex.label = NA)

Section 13.5. Using Text Labels in a Network Graph

# Figure 13-10
m <- madmen[1:nrow(madmen) %% 2 == 1, ]
g <- graph.data.frame(m, directed = FALSE)
plot(g, layout = layout.fruchterman.reingold,
 vertex.size = 4, 
 vertex.label = V(g)$name, 
 vertex.label.cex = 0.8,
 vertex.label.dist = 0.4,
 vertex.label.color = "black")

# Figure 13-11
E(g)[c(2, 11, 19)]$label <- "M"
E(g)$color <- "grey70"
E(g)[c(2, 11, 19)]$color <- "red"
plot(g)

Section 13.6. Creating a Heat Map

pres_rating <- data.frame(
 rating = as.numeric(presidents),
 year = as.numeric(floor(time(presidents))),
 quarter = as.numeric(cycle(presidents))
)

# Figure 13-12
p <- ggplot(pres_rating, aes(year, quarter, fill = rating))
p + geom_tile()

# Figure 13-13
p + geom_tile() +
 scale_x_continuous(breaks = seq(1940, 1976, by = 4)) +
 scale_y_reverse() +
 scale_fill_gradient2(midpoint = 50, mid = "grey70", limits = c(0, 100))

Section 13.7. Creating a Three-Dimensional Scatter Plot

# Figure 13-14
plot3d(mtcars$wt, mtcars$disp, mtcars$mpg, type = "s", size = 0.75, lit = FALSE)

# Figure 13-15
interleave <- function(v1, v2) as.vector(rbind(v1,v2))
plot3d(mtcars$wt, mtcars$disp, mtcars$mpg,
 xlab = "Weight", ylab = "Displacement", zlab = "MPG",
 size = .75, type = "s", lit = FALSE)
segments3d(interleave(mtcars$wt, mtcars$wt),
 interleave(mtcars$disp, mtcars$disp),
 interleave(mtcars$mpg, min(mtcars$mpg)),
 alpha = 0.4, col = "blue")

# Figure 13-16
plot3d(mtcars$wt, mtcars$disp, mtcars$mpg,
 xlab = "", ylab = "", zlab = "",
 axes = FALSE,
 size = .75, type = "s", lit = FALSE)
segments3d(interleave(mtcars$wt, mtcars$wt),
 interleave(mtcars$disp, mtcars$disp),
 interleave(mtcars$mpg, min(mtcars$mpg)),
 alpha = 0.4, col = "blue")
rgl.bbox(color = "grey50",
 emission = "grey50", 
 xlen = 0, ylen = 0, zlen = 0)
rgl.material(color = "black")
axes3d(edges = c("x--", "y+-", "z--"),
 ntick = 6 ,
 cex = .75)
mtext3d("Weight", edge = "x--", line = 2)
mtext3d("Displacement", edge = "y+-", line = 3)
mtext3d("MPG", edge = "z--", line = 3)

Section 13.8. Adding a Prediction Surface to a Three-Dimensional Plot

predictgrid <- function(model, xvar, yvar, zvar, res = 16, type = NULL) {
 xrange <- range(model$model[[xvar]])
 yrange <- range(model$model[[yvar]])
 newdata <- expand.grid(x = seq(xrange[1], xrange[2], length.out = res),
 y = seq(yrange[1], yrange[2], length.out = res))
 names(newdata) <- c(xvar, yvar)
 newdata[[zvar]] <- predict(model, newdata = newdata, type = type)
 newdata
}

df2mat <- function(p, xvar = NULL, yvar = NULL, zvar = NULL) {
 if (is.null(xvar)) xvar <- names(p)[1]
 if (is.null(yvar)) yvar <- names(p)[2]
 if (is.null(zvar)) zvar <- names(p)[3]
 x <- unique(p[[xvar]])
 y <- unique(p[[yvar]])
 z <- matrix(p[[zvar]], nrow = length(y), ncol = length(x))
 m <- list(x, y, z)
 names(m) <- c(xvar, yvar, zvar)
 m
}

# Figure 13-17
m <- mtcars
mod <- lm(mpg ~ wt + disp + wt:disp, data = m)
m$pred_mpg <- predict(mod)
mpgrid_df <- predictgrid(mod, "wt", "disp", "mpg")
mpgrid_list <- df2mat(mpgrid_df)

plot3d(m$wt, m$disp, m$mpg, type = "s", size = 0.5, lit = FALSE)
spheres3d(m$wt, m$disp, m$pred_mpg, alpha = 0.4, type = "s", size = 0.5, lit = FALSE)
segments3d(interleave(m$wt, m$wt),
 interleave(m$disp, m$disp),
 interleave(m$mpg, m$pred_mpg),
 alpha = 0.4, col = "red")
surface3d(mpgrid_list$wt, mpgrid_list$disp, mpgrid_list$mpg,
 alpha = 0.4, front = "lines", back = "lines")

# Figure 13-18
plot3d(mtcars$wt, mtcars$disp, mtcars$mpg,
 xlab = "", ylab = "", zlab = "",
 axes = FALSE,
 size=.5, type="s", lit=FALSE)
spheres3d(m$wt, m$disp, m$pred_mpg, alpha = 0.4, type = "s", size = 0.5, lit = FALSE)
segments3d(interleave(m$wt, m$wt),
 interleave(m$disp, m$disp),
 interleave(m$mpg, m$pred_mpg),
 alpha = 0.4, col = "red")
surface3d(mpgrid_list$wt, mpgrid_list$disp, mpgrid_list$mpg,
 alpha = 0.4, front = "lines", back = "lines")
rgl.bbox(color = "grey50", 
 emission = "grey50", 
 xlen = 0, ylen = 0, zlen = 0) 
rgl.material(color = "black")

axes3d(edges = c("x--", "y+-", "z--"),
 ntick = 6,
 cex = .75) 

mtext3d("Weight", edge = "x--", line = 2)
mtext3d("Displacement", edge = "y+-", line = 3)
mtext3d("MPG", edge = "z--", line = 3)

Section 13.9. Saving a Three-Dimensional Plot

plot3d(mtcars$wt, mtcars$disp, mtcars$mpg, type = "s", size = 0.75, lit = FALSE)
rgl.snapshot('3dplot.png', fmt = 'png')

# Save the current viewpoint
view <- par3d("userMatrix")

# Restore the saved viewpoint
par3d(userMatrix = view)

Section 13.10. Animating a Three-Dimensional Plot

Section 13.11. Creating a Dendrogram

# Figure 13-19
c2 <- subset(countries, Year == 2009)
c2 <- c2[complete.cases(c2),]
c2 <- c2[sample(1:nrow(c2), 25), ]
rownames(c2) <- c2$Name
c2 <- c2[,4:7]
c3 <- scale(c2)
hc <- hclust(dist(c3))
plot(hc)

plot(hc, hang = -1)

# Figure 13-20
hc1 <- hclust(dist(c2))
plot(hc1)

Section 13.12. Creating a Vector Field

# Figure 13-21
islice <- subset(isabel, z == min(z))
ggplot(islice, aes(x, y)) +
 geom_segment(aes(xend = x + vx/50, yend = y + vy/50),
 size = 0.25) 
## Warning: Removed 3745 rows containing missing values (geom_segment).

# Figure 13-22
islice <- subset(isabel, z == min(z))
every_n <- function(x, by = 2) {
 x <- sort(x)
 x[seq(1, length(x), by = by)]
}
keepx <- every_n(unique(isabel$x), by = 4)
keepy <- every_n(unique(isabel$y), by = 4)

islicesub <- subset(islice, x %in% keepx & y %in% keepy)

ggplot(islicesub, aes(x, y)) +
 geom_segment(aes(xend = x+vx/50, yend = y+vy/50),
 arrow = arrow(length = unit(0.1, "cm")), size = 0.25)
## Warning: Removed 248 rows containing missing values (geom_segment).

# Figure 13-23
islicesub$speedxy <- sqrt(islicesub$vx^2 + islicesub$vy^2)
geom_segment(aes(xend = x + vx/50, yend = y + vy/50, alpha = speed),
 arrow = arrow(length = unit(0.1, "cm")), size = 0.6)
## mapping: xend = x + vx/50, yend = y + vy/50, alpha = speed 
## geom_segment: arrow = list(angle = 30, length = 0.1, ends = 2, type = 1), lineend = butt, na.rm = FALSE
## stat_identity: na.rm = FALSE
## position_identity
usa <- map_data("usa")
ggplot(islicesub, aes(x, y)) +
 geom_segment(aes(xend = x + vx/50, yend = y + vy/50, colour = speed),
 arrow = arrow(length = unit(0.1, "cm")), size = 0.6) +
 scale_colour_continuous(low = "grey80", high = "darkred") +
 geom_path(aes(long, lat, group = group), data = usa) +
 coord_cartesian(xlim = range(islicesub$x), ylim = range(islicesub$y))
## Warning: Removed 248 rows containing missing values (geom_segment).

Section 13.13. Creating a QQ Plot

# Figure 13-25
qqnorm(heightweight$heightIn)
qqline(heightweight$heightIn)

qqnorm(heightweight$ageYear)
qqline(heightweight$ageYear)

Section 13.14. Creating a Plot of an Empirical Cumulative Distribution Function

# Figure 13-26
ggplot(heightweight, aes(x = heightIn)) + stat_ecdf()

ggplot(heightweight, aes(x = ageYear)) + stat_ecdf()

Section 13.15. Creating a Mosaic Plot

# Figure 13-27
mosaic( ~ Admit + Gender + Dept, data = UCBAdmissions)

# Figure 13-28
mosaic( ~ Dept + Gender + Admit, data = UCBAdmissions,
 highlighting = "Admit", highlighting_fill = c("lightblue", "pink"),
 direction = c("v","h","v"))

# Figure 13-29
mosaic( ~ Dept + Gender + Admit, data = UCBAdmissions,
 highlighting="Admit", highlighting_fill = c("lightblue", "pink"),
 direction=c("v", "v", "h"))

# Figure 13-30
mosaic( ~ Dept + Gender + Admit, data = UCBAdmissions,
 highlighting = "Admit", highlighting_fill = c("lightblue", "pink"),
 direction = c("v", "h", "h"))

Section 13.16. Creating a Pie Chart

# Figure 13-31
fold <- table(survey$Fold)
pie(fold)

Section 13.17. Creating a Map

# Figure 13-32
states_map <- map_data("state")
ggplot(states_map, aes(long, lat, group = group)) +
 geom_polygon(fill = "white", colour = "black")

ggplot(states_map, aes(long, lat, group = group)) +
 geom_path() + coord_map("mercator")

# Figure 13-33
east_asia <- map_data("world", region = c("Japan", "China", "North Korea",
 "South Korea"))
ggplot(east_asia, aes(long, lat, group = group, fill = region)) +
 geom_polygon(colour = "black") +
 scale_fill_brewer(palette = "Set2")

# Figure 13-34
nz1 <- map_data("world", region = "New Zealand")
nz1 <- subset(nz1, long > 0 & lat > -48)
ggplot(nz1, aes(long, lat, group = group)) + geom_path()

nz2 <- map_data("nz")
ggplot(nz2, aes(long, lat, group = group)) + geom_path()

Section 13.18. Creating a Chloropleth Map

# Figure 13-35
crimes <- data.frame(state = tolower(rownames(USArrests)), USArrests)
states_map <- map_data("state")
crime_map <- merge(states_map, crimes, by.x = "region", by.y = "state")
crime_map <- arrange(crime_map, group, order)

ggplot(crime_map, aes(long, lat, group = group, fill = Assault)) +
 geom_polygon(colour = "black") +
 coord_map("polyconic")

# Figure 13-36
ggplot(crimes, aes(map_id = state, fill = Assault)) +
 geom_map(map = states_map, colour = "black") +
 scale_fill_gradient2(low = "#559999", mid = "grey90", high = "#BB650B",
 midpoint = median(crimes$Assault)) +
 expand_limits(x = states_map$long, y = states_map$lat) +
 coord_map("polyconic")

# Figure 13-37
qa <- quantile(crimes$Assault, c(0, 0.2, 0.4, 0.6, 0.8, 1.0))
crimes$Assault_q <- cut(crimes$Assault, qa,
 labels = c("0-20%", "20-40%", "40-60%", "60-80%", "80-100%"),
 include.lowest = TRUE)
pal <- colorRampPalette(c("#559999", "grey80", "#BB650B"))(5)

ggplot(crimes, aes(map_id = state, fill = Assault_q)) +
 geom_map(map = states_map, colour = "black") +
 scale_fill_manual(values = pal) +
 expand_limits(x = states_map$long, y = states_map$lat) +
 coord_map("polyconic") +
 labs(fill = "Assault Rate\nPercentile")

Section 13.19. Making a Map with a Clean Background

theme_clean <- function(base_size = 12) {
require(grid)
 theme_grey(base_size) %+replace%
    theme(
      axis.title = element_blank(),
      axis.text = element_blank(),
      panel.background = element_blank(),
      panel.grid = element_blank(),
      axis.ticks.length = unit(0, "cm"),
      axis.ticks.margin = unit(0, "cm"),
      panel.margin = unit(0, "lines"),
      plot.margin = unit(c(0, 0, 0, 0), "lines"),
      complete = TRUE
 )
}

# Figure 13-38
ggplot(crimes, aes(map_id = state, fill = Assault_q)) +
 geom_map(map = states_map, colour = "black") +
 scale_fill_manual(values = pal) +
 expand_limits(x = states_map$long, y = states_map$lat) +
 coord_map("polyconic") +
 labs(fill = "Assault Rate\nPercentile") +
 theme_clean()
## Warning: `axis.ticks.margin` is deprecated. Please set `margin` property of
## `axis.text` instead
## Warning: `panel.margin` is deprecated. Please use `panel.spacing` property
## instead

Section 13.20. Recreating a Map from a Shapefile

taiwan_shp <- readShapePoly("TWN_adm2.shp")
## Warning: use rgdal::readOGR or sf::st_read
taiwan_map <- fortify(taiwan_shp)
## Regions defined for each Polygons
# Figure 13-39
ggplot(taiwan_map, aes(long, lat, group = group)) + geom_path()


END!