6

I want to plot a bunch of rasters and I created a code to adjust breaks for each one and plot them trough a for loop. But i'm getting a problematic color scale bar, and my efforts haven't being effective to solve that. Example:

I have precipitation ranging from 0 to 11.000...but most part of the data is between 0 and 5.000... and very few up to 11.000. So I need to change the breaks to capture this variation... more breaks where I have more data.

Then I created a breaks object for that.
But when I plot the raster, the scale color bar gets awful, very messy...

#get predictors (These are a way lighter version of mine)
predictors_full<-getData('worldclim', var='bio', res=10)

predic_legends<-c(
"Annual Mean Temperature [°C*10]",
"Mean Diurnal Range [°C]",
"Isothermality",
"Temperature Seasonality [standard deviation]",
"Max Temperature of Warmest Month [°C*10]",
"Min Temperature of Coldest Month [°C*10]",
"Temperature Annual Range [°C*10]",
"Mean Temperature of Wettest Quarter [°C*10]",
"Mean Temperature of Driest Quarter [°C*10]",
"Mean Temperature of Warmest Quarter [°C*10]",
"Mean Temperature of Coldest Quarter [°C*10]",
"Annual Precipitation [mm/year]",
"Precipitation of Wettest Month [mm/month]",
"Precipitation of Driest Month [mm/month]",
"Precipitation Seasonality [coefficient of variation]",
"Precipitation of Wettest Quarter [mm/quarter]",
"Precipitation of Driest Quarter [mm/quarter]",
"Precipitation of Warmest Quarter [mm/quarter]",
"Precipitation of Coldest Quarter [mm/quarter]",
)

# Crop rasters and rename
xmin=-120; xmax=-35; ymin=-60; ymax=35
limits <- c(xmin, xmax, ymin, ymax)
predictors <- crop(predictors_full,limits)

predictor_names<-c("mT_annual","mT_dayn_rg","Isotherm","T_season",
"maxT_warm_M","minT_cold_M","rT_annual","mT_wet_Q","mT_dry_Q",
"mT_warm_Q","mT_cold_Q","P_annual","P_wet_M","P_dry_M","P_season",
"P_wet_Q","P_dry_Q","P_warm_Q","P_cold_Q")

names(predictors)<-predictor_names

#Set a palette
Blues_up<-c('#fff7fb','#ece7f2','#d0d1e6','#a6bddb','#74a9cf','#3690c0','#0570b0','#045a8d','#023858','#233159')
colfunc_blues<-colorRampPalette(Blues_up)

#Create a loop to plot all my Predictor rasters
for (i in 1:19) {
#save a figure
png(file=paste0(predictor_names[[i]],".png"),units="in", width=12, height=8.5, res=300)

#Define a plot area
par(mar = c(2,2, 3, 3), mfrow = c(1,1))

#extract values from rasters
vmax<- maxValue(predictors[[i]])
vmin<-minValue(predictors[[i]])
vmedn=(maxValue(predictors[[i]])-minValue(predictors[[i]]))/2

#breaks 
break1<-c((seq(from=vmin,to= vmedn, length.out = 40)),(seq(from=(vmedn+(vmedn/5)),to=vmax,length.out = 5)))

#plot without the legend because the legend would come out with really messy, with too many marks and uneven spaces
plot(predictors[[i]], col =colfunc_blues(45) , breaks=break1,  margin=FALSE, 
            main =predic_legends[i],legend.shrink=1)
dev.off()
}

Precipitation raster according to the above code This figure is the i=12 from all rasters in the loop

Then I wrote a different code to set different breaks to the color bar

#Plot the raster with no color scale bar    
plot(predictors[[i]], col =colfunc_blues(45) , breaks=break1,  margin=FALSE, 
        main =predic_legends[i],legend=FALSE)

#breaks for the color scale
def_breaks = seq(vmax,vmin,length.out=(10))

#plot only the legend
image.plot(predictors_full[[i]], zlim = c(vmin,vmax), 
             legend.only = TRUE, col = colfunc_greys(30),
             axis.args = list(at = def_breaks, labels =def_breaks,cex.axis=0.5)) 

But that doesn't work, because the colors don't really match the numbers in the map... Look at the color for 6.000 in each map... It's different.

Different breaks for plot and scale

Any tips on how to proceed on that? I'm new to R so I struggle a lot to reach my goals... Also, I'm getting a lot of decimal places in the numbers... how to change that for 2 decimal places?

EDIT: @jbaums taught me to use log... I liked but it's not yet what I seek

levelplot(predictors[[12]]+1, col.regions=colorRampPalette(brewer.pal(9, 'Blues')), zscaleLog=TRUE, at=seq(1, 4, len=100), margin=FALSE)

Using log

Thai
  • 493
  • 5
  • 15
  • What do you actually want it to look like? – jbaums Aug 07 '17 at 05:13
  • 1
    I dont want the colors equally distributed along numbers because i have very few data on high and low numbers and too much in the middle... If they are equally distributed, I get a map with almost no light and dark blue, and everything would be "middle blue" with few blue variation. I want one shade of light blue for one extreme, one dark blue for the other extreme, and the rest of the blue along the data where there are indeed data. The codes I posted were my attempt to do that... but I was introduced to programming and R last month... I lack basis, although I have been reading a lot – Thai Aug 07 '17 at 11:17
  • @jbaums... Did you understand this explanation? Please, let me know if I was not clear enought, and I'll try to do better! Thank you in advance for your attentio! – Thai Aug 07 '17 at 11:18
  • 2
    You could try a log scale. E.g.: `library(rasterVis); library(RColorBrewer); levelplot(predictors[[12]]+1, col.regions=colorRampPalette(brewer.pal(9, 'Blues')), zscaleLog=TRUE, at=seq(1, 4, len=100), margin=FALSE)` – jbaums Aug 07 '17 at 11:46
  • Thank you @jbaums... I really liked to learn this possinility! But I'll keep trying other options because using log I loose the easy information of precipitation amount in the scale, although I can see where rains more.. – Thai Aug 07 '17 at 12:59

1 Answers1

5

You can avoid log scale (as some users said you) using classIntervals() function from classInt package.

Using levelplot() (in my opinion the result is better than raster::plot() function):

# Normal breaks
break1 <- classIntervals(predictors[[12]][!is.na(predictors[[12]])], n = 50, style = "equal")

levelplot(predictors[[12]], col.regions=colorRampPalette(brewer.pal(9, 'Blues')), at=break1$brks, margin=FALSE,main =predic_legends[12])

enter image description here

# Using quantiles
break1 <- classIntervals(predictors[[12]][!is.na(predictors[[12]])], n = 50, style = "quantile")

levelplot(predictors[[12]], col.regions=colorRampPalette(brewer.pal(9, 'Blues')), at=break1$brks, margin=FALSE,main =predic_legends[12])

enter image description here

Also, you have more options to choose, such like sd, pretty, kmeans, hclust and others.


Adding polygons and point to the plot

First, I'll save the plot above to p, the line is too long for this example:

p <- levelplot(predictors[[12]], col.regions=colorRampPalette(brewer.pal(9, 'Blues')), at=break1$brks, margin=FALSE,main =predic_legends[12])

I'll use the same data than your, wrld_simpl data, as polygons to add into the plot and I'll create points to be added to the plot also.

library(maptools)
library(rgeos)

data(wrld_simpl)
pts <- gCentroid(wrld_simpl, byid = T)

To add lines, polygons, points or even text, you can use layer() function and a panel.spplot object:

p + layer(sp.polygons(wrld_simpl)) + layer(sp.points(pts))

enter image description here

Finally, you can also change color, fill, symbology, and so on:

p + layer(sp.polygons(wrld_simpl,col='firebrick')) + layer(sp.points(pts,pch = 12,col='red'))

enter image description here

Check ?panel.spplot for more information.

aldo_tapia
  • 1,153
  • 16
  • 27
  • oh my! Good darwin! Exactly what I needed! @aldo_tapia, muito obrigada!!! – Thai Aug 08 '17 at 23:34
  • levelplot indeed produces better results then plot... but I was introduced to R last month, thus I'm still learning and I have some troubles adding one plot over another with levelplot, because I always plot in the same map a raster, lines for countries from wrld_simpl, spatialPoints and a polygon... with reguar plot it is easy, I just do `plot(x, add=TRUE)` for everything I want to plot... with levelplot is different. But I'll dedicate some time to learn levelplot and ggplot as soon as possible! Thank you very much, @aldo_tapia – Thai Aug 08 '17 at 23:45
  • 1
    @Thai Let me help you. I'm far of my machine right now, but I will improve my answer as soon as posible – aldo_tapia Aug 09 '17 at 00:45
  • I'm already using and exploring everything you taught me! Thanks, @aldo_tapia !!! – Thai Aug 09 '17 at 23:08