plotBeta with ggplot2 - Githubissues

Hi all,
I wrote a some code converting the plotBeta function from the Hmsc package to make an equivalent plot with ggplot2. I know that previously someone had posted about the difficulty of adjusting the plot parameters for plotBeta (see issue #64). Since I managed to get this working in ggplot2, I wanted to make the code available for the HMSC community to use. I hope it is alright posting it here, even though, strictly speaking, this is not an issue.
#Converting HMSC::plotBeta to ggplot figure

library(Hmsc)
library(tidyverse)

#'[Note: m = Hmsc fitted model]

m<-your.fitted.model

#Get posterior estimates of Beta as usual
postBeta = getPostEstimate(m, parName="Beta")

#Get species and covariate names from model
spNames<- m$spNames
covNames<-m$covNames

#pull the data from postBeta for plotting
mbeta = postBeta$mean
betaP = postBeta$support

#Plot param == "Sign" ----

#'[Note: This will plot the sign of the beta parameters for all parameters
#'[      with a support level greater than the one specified

#Choose the support level for Beta plot
supportLevel <- 0.95

#format the data as sign and narrow to the covariates that have the desired support level 
toPlot = sign(mbeta)
toPlot = toPlot * ((betaP > supportLevel) + (betaP < (1 - supportLevel)) > 0)

#format the data as matrix and add column and row names
betaMat = matrix(toPlot, nrow = m$nc, ncol = ncol(m$Y))
colnames(betaMat)<- spNames
rownames(betaMat) <- covNames

#remove intercept
betaMat<- as.data.frame(betaMat) %>%
  dplyr::slice(-1)

#update covariate names
rownames(betaMat) <-c('Depth (log)', 'Rugosity', 'Sediment', 'SST','Chla', 'FSLE', 
  'MLD', 'Salinity', 'SST Fprob', 'Chl Fprob')

#reorder species alphabetically, so displays from top to bottom A-Z
betaMat<-betaMat[, rev(order(colnames(betaMat)))]

#reformat for ggplot
betaMatmelt<-as.data.frame(melt(as.matrix(betaMat)))

colors = colorRampPalette(c('#370f71', '#fdfdf8', '#fdfdf8', '#ff6700'))
colorLevels<-3
cols_to_use= colors(colorLevels)

fname<-paste0(draft_fig_dir, "/fall_fprob_plotBeta_sign.tif")

tiff(filename = fname, width = 7.5, height = 8, units = "in", compression = "lzw", 
     res = 300)

ggplot(betaMatmelt, aes(x = Var1, y = Var2, fill = factor(value))) +
  labs(x = "Environmental Niche", y = "Species", fill = "Sign") +
  #geom_raster() +
  geom_tile(color = 'gray60')+
  scale_fill_manual(breaks = levels(factor(betaMatmelt$value)),
                    values = cols_to_use,
                    labels = c('\u2013', '0', '+')) +
  theme(plot.background = element_rect(fill = "white"),
        legend.background = element_rect(fill = "white"),
        panel.border = element_rect(fill = NA, color = NA),
        legend.margin = margin(l = 1, unit = 'cm'),
        legend.title = element_text(hjust = 0.1),
        legend.key.width = unit(1, 'cm'),
        legend.key.height = unit(4, 'cm'),
        legend.text = element_text(size = 12),
        axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5))+
  scale_x_discrete(expand = c(0, 0)) + #because his x and y axis are numbers
  scale_y_discrete(expand = c(0, 0))

dev.off()

#plot param = "Mean" ------

#'[Note: This will plot the magnitude of the beta parameters for all parameters
#'[      with a support level greater than the one specified

#Choose the support level for Beta plot
supportLevel <- 0.95

#format data as mean for all values with .95 support level
toPlot = mbeta
toPlot = toPlot * ((betaP > supportLevel) + (betaP < (1 - supportLevel)) > 0)

#replace zeros with NA
toPlot[toPlot == 0] <- NA

#format the data as matrix and add column and row names
betaMat = matrix(toPlot, nrow = m$nc, ncol = ncol(m$Y))
colnames(betaMat)<- spNames
rownames(betaMat) <- covNames

#remove intercept
betaMat<- as.data.frame(betaMat) %>%
  dplyr::slice(-1)

#update covariate names
rownames(betaMat) <-c('Depth (log)', 'Rugosity', 'Sediment', 'SST','Chla', 'FSLE', 
                      'MLD', 'Salinity', 'SST Fprob', 'Chl Fprob')

#reorder species alphabetically, so displays from top to bottom A-Z
betaMat<-betaMat[, rev(order(colnames(betaMat)))]

#reformat for ggplot
betaMatmelt<-as.data.frame(melt(as.matrix(betaMat)))

colors = colorRampPalette(c('#370f71', '#fdfdf8', '#fdfdf8', '#ff6700'))
colorLevels<-3
cols_to_use= colors(colorLevels)

#figure out the scale - use this to determine what to set n.breaks to - may
#                       need to play with the n.breaks value a bit and see what works best
min(betaMatmelt$value)
max(betaMatmelt$value)
length(seq(-1,6, 1))

fname<-paste0(draft_fig_dir, "/fall_fprob_plotBeta_mean.tif")

tiff(filename = fname, width = 7.5, height = 8, units = "in", compression = "lzw", 
     res = 300)

ggplot(betaMatmelt, aes(x = Var1, y = Var2, fill = value)) +
  labs(x = "Environmental Niche", y = "Species", fill = "Mean") +
  geom_tile(color = 'gray60')+
  scale_fill_steps2(n.breaks = 7,
                     high = cols_to_use[3],
                     mid = 'white',
                     low = cols_to_use[1],
                    nice.breaks = TRUE,
                    na.value = "white",
                     ) +
  theme(plot.background = element_rect(fill = "white"),
        legend.background = element_rect(fill = "white"),
        panel.border = element_rect(fill = NA, color = NA),
        legend.margin = margin(l = 1, unit = 'cm'),
        legend.title = element_text(hjust = 0.08),
        legend.key.width = unit(1, 'cm'),
        legend.key.height = unit(2.25, 'cm'),
        legend.text = element_text(size = 12),
        axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5))+
  scale_x_discrete(expand = c(0, 0)) + 
  scale_y_discrete(expand = c(0, 0))

dev.off()
Enjoy! Chandra Goetsch
hmsc-r / HMSC

plotBeta with ggplot2 #139
