1

I am trying to recreate the following plsr biplot:

plsr loading plots

plsr code

df.metric <- plsr(y ~ LMA + LDMC + Thick + Carbon + Nitrogen + Tough, scale 
= TRUE, validation = "LOO", method = "oscorespls", data = df)

extract fungal taxa loadings

df2<-df.metric$Yloadings
comp1a <- df2[,1]
comp2a <- df2[,2]
namesa <- df2[,0]
df2<-as.data.frame(cbind(namesa,comp1a, comp2a))

extract leaf traits loadings

df1<-df.metric$loadings
comp1 <- df1[,1]
comp2 <- df1[,2]
names <- df1[,0]
df1<-as.data.frame(cbind(names, comp1, comp2))

Generate two plots, one for fungal taxa and one for leaf traits

#generate fungal taxa plot
plot.fungal.taxa<-ggplot(data=df2, aes(comp1a,comp2a))+
ylab("")+
xlab("")+
theme_bw()+ 
theme(panel.border = element_rect(colour = "black", fill=NA, 
size=1),panel.grid.major = element_blank(), 
                 panel.grid.minor = element_blank(), 
                 axis.line = element_line(colour = "black"))+
geom_text(aes(label=rownames(df2)), color="red")+
scale_x_continuous(breaks = c(0.10,0.05,0,-0.05,-0.10,-0.15))+ 
scale_y_continuous(breaks = c(0.10,0.05,0,-0.05,-0.10,-0.15))+
coord_fixed(ylim=c(0.10, -0.15),xlim=c(0.10, -0.15))+
theme(axis.ticks = element_line(colour = "red")) +
theme(axis.text.y=element_text(angle = 90, hjust = 0.65)) +
theme(axis.text.y = element_text(margin=margin(10,10,10,5,"pt")))

#generate leaf traits plot
plot.leaf.traits<-ggplot(data=df1, aes(comp1,comp2))+
ylab("Comp 2")+
xlab("Comp 1")+
theme_bw() + 
theme(panel.border = element_rect(colour = "black", fill=NA, size=1),
                 panel.grid.major = element_blank(), 
                 panel.grid.minor = element_blank(), 
                 axis.line = element_line(colour = "black"))+
geom_text(aes(label=rownames(df1)), color="black")+
scale_x_continuous(breaks = c(-0.8,-0.6,-0.4,-0.2,0,0.2,0.4,0.6))+ 
scale_y_continuous(breaks = c(-0.8,-0.6,-0.4,-0.2,0,0.2,0.4,0.6))+
coord_fixed(ylim=c(0.6, -0.8),xlim=c(0.6, -0.8))+
theme(axis.ticks = element_line(colour = "black")) +
theme(axis.text.y=element_text(angle = 90, hjust = 0.65)) +
theme(axis.text.y = element_text(margin=margin(10,10,10,5,"pt")))

function to overlay plots

ggplot_dual_axis = function(plot.leaf.traits, plot.fungal.taxa, which.axis = 
"x")
{
# Update plot with transparent panel
plot.fungal.taxa = plot.fungal.taxa + theme(panel.background = 
element_rect(fill = NA))
grid.newpage()
# Increase right margin if which.axis == "y"
if(which.axis == "y") plot.leaf.traits = plot.leaf.traits + 
theme(plot.margin = unit(c(0.7, 1.5, 0.4, 0.4), "cm"))
# Extract gtable
g1 = ggplot_gtable(ggplot_build(plot.leaf.traits))
g2 = ggplot_gtable(ggplot_build(plot.fungal.taxa))
# Overlap the panel of the second plot on that of the first
pp = c(subset(g1$layout, name == "panel", se = t:r))
g = gtable_add_grob(g1, g2$grobs[[which(g2$layout$name=="panel")]], pp$t, 
pp$l, pp$b, pp$l)
# Steal axis from second plot and modify
axis.lab = ifelse(which.axis == "x", "axis-b", "axis-l")
ia = which(g2$layout$name == axis.lab)
ga = g2$grobs[[ia]]
ax = ga$children[[2]]
# Switch position of ticks and labels
if(which.axis == "x") ax$heights = rev(ax$heights) else ax$widths = 
rev(ax$widths)
ax$grobs = rev(ax$grobs)
if(which.axis == "x") 
ax$grobs[[2]]$y = ax$grobs[[2]]$y - unit(1, "npc") + unit(0.15, "cm") else
ax$grobs[[1]]$x = ax$grobs[[1]]$x - unit(1, "npc") + unit(0.15, "cm")
# Modify existing row to be tall enough for axis
if(which.axis == "x") g$heights[[2]] = g$heights[g2$layout[ia,]$t]
# Add new row or column for axis label
if(which.axis == "x") {
g = gtable_add_grob(g, ax, 2, 4, 2, 4) 
g = gtable_add_rows(g, g2$heights[1], 1)
g = gtable_add_grob(g, g2$grob[[6]], 2, 4, 2, 4)
} else {
g = gtable_add_cols(g, g2$widths[g2$layout[ia, ]$l], length(g$widths) - 1)
g = gtable_add_grob(g, ax, pp$t, length(g$widths) - 1, pp$b) 
g = gtable_add_grob(g, g2$grob[[7]], pp$t, length(g$widths), pp$b - 1)
}
# Draw it
grid.draw(g)
}

Run function on individual plots

ggplot_dual_axis(plot.leaf.traits, plot.fungal.taxa, "y")

this is what I end up getting:

plsr loading plot using ggplot2

My question is how to I get the top x axis to match on top of the plot? Currently it sits on top and adjacent to the plot. I used a previous code I found here (Plotting Partial Least Squares Regression (plsr) biplot with ggplot2). Any help would be amazing!

Tellez
  • 11
  • 1

1 Answers1

0

In the future when you're posting a question, you should try to include a minimal reproducible example. You almost did that, except you didn't include any data to work with.

If you have a question related to code from a particular package, you should be able to grab data that comes with that package or R, just look at the help section for any of the functions. For example, below I just copied and pasted the line in the PLSR package to make an example PLSR model. Alternatively you could have grabbed the table from the post your reference at the end of your question.

The main problem is that the code you're building off of broke when ggplot updated. You can follow that conversation here and here. Below is some code that should work, with the package version numbers noted at the top.

#Make a ggplot object for a Partial Least Squares Regression (PLSR) plot

#####################################################
## Note that this code may break as ggplot updates,##
## as is noted on some of the below posts.         ##
#####################################################

#Mostly taken from the posts below

#Links to posts-------------

#https://stackoverflow.com/questions/48664746/how-to-set-two-x-axis-and-two-y-axis-using-ggplot2

#https://stackoverflow.com/questions/39137287/plotting-partial-least-squares-regression-plsr-biplot-with-ggplot2

#https://stackoverflow.com/questions/21026598/ggplot2-adding-secondary-transformed-x-axis-on-top-of-plot

#https://stackoverflow.com/questions/36754891/ggplot2-adding-secondary-y-axis-on-top-of-a-plot/36759348#36759348


#load libraries------
library(pls) #version 2.7.1
library(ggplot2) #version 3.1.0
library(grid) #version 3.5.1
library(gtable) #version 0.2.0
library(cowplot) #version 0.9.3
library(ggplotify) #version 0.0.3

#Read data into PLSR model-----
dens1 <- plsr(density ~ NIR, ncomp = 5, data = yarn)

#Extract information from plsr (AKA mvr) model----
df1<-as.data.frame(dens1$loadings[,1:2])
names(df1) <- c("comp1", "comp2")

df2<-as.data.frame(dens1$scores[,1:2])
names(df2) <- c("comp1a", "comp2a")

#make ggplot objects------

#Plot Loadings - colored red

p1 <- ggplot(data=df1, 
           aes(x = comp1, y = comp2)) +
  geom_text(aes(label = rownames(df1)), 
            color="red") +
  theme_bw() + 
  theme(panel.border = element_rect(colour = "black", 
                                    fill=NA, 
                                    size=1),
        panel.grid.major = element_blank(), 
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black"),
        axis.ticks = element_line(colour = "red"),
        axis.text.y = element_text(margin = margin(10,10,10,5,"pt"), 
                                   angle = 90, 
                                   hjust = 0.65, 
                                   colour = "red"),
        axis.text.x = element_text(colour = "red")) +
  scale_y_continuous(limits = c(min(df1), max(df1))) +
  scale_x_continuous(limits = c(min(df1), max(df1)))



#Plot 2 - scores in black
p2 <- ggplot(data=df2, 
             aes(x = comp1a, y = comp2a)) +
  geom_text(aes(label = rownames(df2)), 
            color="black") +
  theme_bw() + 
  theme(panel.border = element_rect(colour = "black", 
                                    fill=NA, 
                                    size=1),
        panel.grid.major = element_blank(), 
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black"),
        axis.ticks = element_line(colour = "black"),
        axis.text.y = element_text(margin = margin(10,10,10,5,"pt"), 
                                   angle = 90, 
                                   hjust = 0.65, 
                                   colour = "black"),
        axis.text.x = element_text(colour = "black")) +
  scale_y_continuous(limits = c(min(df2), max(df2))) +
  scale_x_continuous(limits = c(min(df2), max(df2)))

#Final plot----

#Overlay plots in order to get two graphs with different axes on same plot
#rename plots in case you want to make adjustments without regenerating plots

plot1 <- p1
plot2 <- p2

# Update plot with transparent panel
plot2 = plot2 + 
    theme(panel.background = element_rect(fill = "transparent")) 

#clean plot space
grid.newpage()

# Extract gtables from ggplot objects
g1 = ggplot_gtable(ggplot_build(plot1))
g2 = ggplot_gtable(ggplot_build(plot2))

# Get the location of the plot panel in g1.
# These are used later when transformed elements of g2 are put back into g1

pp <- c(subset(g1$layout, name == "panel", se = t:r))

# Overlap panel for second plot on that of the first plot

g1 <- gtable_add_grob(g1, g2$grobs[[which(g2$layout$name == "panel")]], pp$t, pp$l, pp$b, pp$l)

#Note from stack overflow post:  
# Get the location of the plot panel in g1.
# These are used later when transformed elements of g2 are put back into g1
pp <- c(subset(g1$layout, name == "panel", se = t:r))

# Overlap panel for second plot on that of the first plot
g1 <- gtable_add_grob(g1, g2$grobs[[which(g2$layout$name == "panel")]], pp$t, pp$l, pp$b, pp$l)

# Then proceed as before:

# ggplot contains many labels that are themselves complex grob; 
# usually a text grob surrounded by margins.
# When moving the grobs from, say, the left to the right of a plot,
# Make sure the margins and the justifications are swapped around.
# The function below does the swapping.
# Taken from the cowplot package:
# https://github.com/wilkelab/cowplot/blob/master/R/switch_axis.R 

hinvert_title_grob <- function(grob){

  # Swap the widths
  widths <- grob$widths
  grob$widths[1] <- widths[3]
  grob$widths[3] <- widths[1]
  grob$vp[[1]]$layout$widths[1] <- widths[3]
  grob$vp[[1]]$layout$widths[3] <- widths[1]

  # Fix the justification
  grob$children[[1]]$hjust <- 1 - grob$children[[1]]$hjust 
  grob$children[[1]]$vjust <- 1 - grob$children[[1]]$vjust 
  grob$children[[1]]$x <- unit(1, "npc") - grob$children[[1]]$x
  grob
}

# Get the y axis title from g2

# Which grob contains the y axis title?
index <- which(g2$layout$name == "ylab-l") 

# Extract that grob
ylab <- g2$grobs[[index]]                

# Swap margins and fix justifications
ylab <- hinvert_title_grob(ylab)         

# Put the transformed label on the right side of g1
g1 <- gtable_add_cols(g1, g2$widths[g2$layout[index, ]$l], pp$r)
g1 <- gtable_add_grob(g1, ylab, pp$t, pp$r + 1, pp$b, pp$r + 1, clip = "off", name = "ylab-r")

# Get the y axis from g2 (axis line, tick marks, and tick mark labels)

# Which grob
index <- which(g2$layout$name == "axis-l")  

# Extract the grob
yaxis <- g2$grobs[[index]]                  

# yaxis is a complex of grobs containing the axis line, the tick marks, and the tick mark labels.
# The relevant grobs are contained in axis$children:
#   axis$children[[1]] contains the axis line;
#   axis$children[[2]] contains the tick marks and tick mark labels.

# First, move the axis line to the left
yaxis$children[[1]]$x <- unit.c(unit(0, "npc"), unit(0, "npc"))

# Second, swap tick marks and tick mark labels
ticks <- yaxis$children[[2]]
ticks$widths <- rev(ticks$widths)
ticks$grobs <- rev(ticks$grobs)

# Third, move the tick marks
ticks$grobs[[1]]$x <- ticks$grobs[[1]]$x - unit(1, "npc") + unit(1, "mm")

# Fourth, swap margins and fix justifications for the tick mark labels
ticks$grobs[[2]] <- hinvert_title_grob(ticks$grobs[[2]])

# Fifth, put ticks back into yaxis
yaxis$children[[2]] <- ticks

# Put the transformed yaxis on the right side of g1
g1 <- gtable_add_cols(g1, g2$widths[g2$layout[index, ]$l], pp$r)
g1 <- gtable_add_grob(g1, yaxis, pp$t, pp$r + 1, pp$b, pp$r + 1, clip = "off", name = "axis-r")

#Draw it for a dummy check

grid.newpage()
grid.draw(g1)

# function that can vertically invert a title grob, with margins treated properly
# title grobs are used a lot in the new ggplot2 version (>1.0.1)
vinvert_title_grob <- function(grob) {
  heights <- grob$heights
  grob$heights[1] <- heights[3]
  grob$heights[3] <- heights[1]
  grob$vp[[1]]$layout$heights[1] <- heights[3]
  grob$vp[[1]]$layout$heights[3] <- heights[1]

  grob$children[[1]]$hjust <- 1 - grob$children[[1]]$hjust
  grob$children[[1]]$vjust <- 1 - grob$children[[1]]$vjust
  grob$children[[1]]$y <- unit(1, "npc") - grob$children[[1]]$y
  grob
  }

# Copy title xlab from g2 and swap margins
index <- which(g2$layout$name == "xlab-b")
xlab <- g2$grobs[[index]]
xlab <- vinvert_title_grob(xlab)

# Put xlab at the top of g1
g1 <- gtable_add_rows(g1, g2$heights[g2$layout[index, ]$t], pp$t-1)
g1 <- gtable_add_grob(g1, xlab, pp$t, pp$l, pp$t, pp$r, clip = "off", name="xlab-t")

# Get "feet" axis (axis line, tick marks and tick mark labels) from g2
index <- which(g2$layout$name == "axis-b")
xaxis <- g2$grobs[[index]]

# Move the axis line to the bottom (Not needed in your example)
xaxis$children[[1]]$y <- unit.c(unit(0, "npc"), unit(0, "npc"))

# Swap axis ticks and tick mark labels
ticks <- xaxis$children[[2]]
ticks$heights <- rev(ticks$heights)
ticks$grobs <- rev(ticks$grobs)

# Move tick marks
ticks$grobs[[2]]$y <- ticks$grobs[[2]]$y - unit(1, "npc") + unit(3, "pt")

# Swap tick mark labels' margins
ticks$grobs[[1]] <- vinvert_title_grob(ticks$grobs[[1]])

# Put ticks and tick mark labels back into xaxis
xaxis$children[[2]] <- ticks

# Add axis to top of g1
g1 <- gtable_add_rows(g1, g2$heights[g2$layout[index, ]$t], pp$t)
g1 <- gtable_add_grob(g1, xaxis, pp$t+1, pp$l, pp$t+1, pp$r, clip = "off", name = "axis-t")

#remove title and axes titles if you want
g1 <- gtable_remove_grobs(g1, c("title", "xlab-t", "xlab-b", "ylab-r", "ylab-l"))

# Draw it
grid.newpage()
my_PLS = ggplotify::as.ggplot(g1)

#save plot in square format----
ggsave(paste0("my_PLS_",Sys.Date(),".png"), width = 6, height = 6, units = "in", plot = my_PLS)

Should look like this:

PLSR_ggplot_example

sformel
  • 1
  • 2