# Plot Leaf Area Index (lai.out) from LPJ-GUESS ( >= v.4) per taxon/year into individual folders # # @author Bikem Ekberzade email{bikemwork@@gmail.com} # # library(tidyverse) library(lubridate) library(patchwork) library(RColorBrewer) #color schemes library(rio) library(rgdal) library(fields) library(maptools) library(sf) library(sp) library(raster) library(maps) library(metR) library(ggmap) ################################################################################ setwd("C:/Users/Documents") LAI <- read.table("lai.out", sep="", header=T) # write.csv(gall13,"LAI_forArcGIS.csv",row.names = F, col.names = T) #to transfer as .csv for further spatial analysis LAI$Year <- as.POSIXct(paste(LAI$Year), format = "%Y") LAI <- as.data.frame(LAI) All_yr <- LAI %>% group_by(year=floor_date(Year, "year")) Per_yr <- All_yr %>% mutate(year=year(Year)) %>% group_by(Year, year) ################################################################################ #Annual loop Years <- c(1961:2100) borders <- sf::st_as_sf(maps::map("world", plot = FALSE, fill = TRUE)) # Alnus glutinosa #to write output dir.create("LAI_loose_scale") setwd("C:/Users/Documents/LAI_loose_scale/") dir.create("Aln_glu") setwd("C:/Users/Documents/LAI_loose_scale/Aln_glu/") for(i in Years){ All_sp <- Per_yr %>% filter(year == i) xmin1 = range(All_sp$Lon)[1] xmax1 = range(All_sp$Lon)[2] ymin1 = range(All_sp$Lat)[1] ymax1 = range(All_sp$Lat)[2] #minimum and maximum value of the trends in grids, #alternately you can find the max/min values per taxa in lai.out and fix the scales for all years for each taxa: #min1 = min(Per_yr$Aln_glu) --> (or 0); max1 = max(Per_yr$Aln_glu) min1<-min(All_sp$Aln_glu) max1<-max(All_sp$Aln_glu) #create breaks and a colour pallete mybreaks<-seq(from=min1, to=max1, length.out = 10) pal11<-rev(colorRampPalette(brewer.pal(9, "Greens"))(length(mybreaks))) #plot the trends Alngl <- ggplot() + geom_sf(data = borders, fill= "white")+ labs(title =paste("Alnus glutinosa",i))+ geom_tile(data = All_sp, aes(fill=All_sp$Aln_glu, x = Lon, y = Lat), colour=NA) + #borders data geom_sf(data = borders,fill = NA,colour = "dark green", size=0.2) + #limiting the borders data with the coordinates of the trend data coord_sf(xlim = c(xmin1, xmax1), ylim = c(ymin1, ymax1)) + #changing color scheme, modifying legend& scale bar scale_fill_gradientn(breaks=mybreaks, limits= c(min1, max1), label=as.character(round(mybreaks, digits=3)), colours=rev(pal11), guide = guide_colourbar(frame.colour = "black", barwidth = 1, barheight=25, ticks = F, title ="Alnus glutinosa", title.position = "top"))+ #and giving uniform ticks to axes scale_x_longitude(ticks = 2) + scale_y_latitude(ticks = 2)+ # #modifying background of the plot theme(legend.position="right", legend.title= element_text(size = 25, color = "dark green"), legend.text = element_text(size = 22), panel.background = element_rect(fill = "#edf0f5"), plot.title = element_text(color="dark green", size=25, face="bold"), axis.line = element_line(colour = "black"), axis.text.x=element_text(size=24), axis.text.y = element_text(size=24), panel.border = element_rect(colour = "black",fill = NA, size = 0.5)) # print(Alngl) png(file=paste(i,"Al_g.png"), width=1980, height=932) print(Alngl) dev.off() } # Abies nordmanniana #to write output dir.create("LAI_loose_scale") setwd("C:/Users/Documents/LAI_loose_scale/") dir.create("Abies_nor") setwd("C:/Users/Documents/LAI_loose_scale/Abies_nor/") for(i in Years){ All_sp <- Per_yr %>% filter(year == i) xmin1 = range(All_sp$Lon)[1] xmax1 = range(All_sp$Lon)[2] ymin1 = range(All_sp$Lat)[1] ymax1 = range(All_sp$Lat)[2] #minimum and maximum value of the trends in grids, #alternately you can find the max/min values per taxa in lai.out and fix the scales for all years for each taxa: #min1 = min(Per_yr$Abi_nor) --> (or 0); max1 = max(Per_yr$Abi_nor) min1<-min(All_sp$Abi_nor) max1<-max(All_sp$Abi_nor) #create breaks and a colour pallete mybreaks<-seq(from=min1, to=max1, length.out = 10) pal11<-rev(colorRampPalette(brewer.pal(9, "Greens"))(length(mybreaks))) #plot the trends Abiesn<- ggplot() + geom_sf(data = borders, fill= "white")+ labs(title =paste("Abies nordmanniana",i))+ geom_tile(data = All_sp, aes(fill=All_sp$Abi_nor, x = Lon, y = Lat), colour=NA) + #borders data geom_sf(data = borders,fill = NA,colour = "dark green", size=0.2) + #limiting the borders data with the coordinates of the trend data coord_sf(xlim = c(xmin1, xmax1), ylim = c(ymin1, ymax1)) + #changing color scheme, modifying legend& scale bar scale_fill_gradientn(breaks=mybreaks, limits= c(min1, max1), label=as.character(round(mybreaks, digits=3)), colours=rev(pal11), guide = guide_colourbar(frame.colour = "black", barwidth = 1, barheight=25, ticks = F, title ="Abies nordmanniana", title.position = "top"))+ #and giving uniform ticks to axes scale_x_longitude(ticks = 2) + scale_y_latitude(ticks = 2)+ # #modifying background of the plot theme(legend.position="right", legend.title= element_text(size = 25, color = "dark green"), legend.text = element_text(size = 22), panel.background = element_rect(fill = "#edf0f5"), plot.title = element_text(color="dark green", size=25, face="bold"), axis.line = element_line(colour = "black"), axis.text.x=element_text(size=24), axis.text.y = element_text(size=24), panel.border = element_rect(colour = "black",fill = NA, size = 0.5)) # print(Abiesn) png(file=paste(i,"Abies_n.png"), width=1980, height=932) print(Abiesn) dev.off() } ################################################################################## # # Add additional taxa as you see fit by copying and pasting the code and # making the necessary adjustment in scale limits, title, etc. including (if necessary) # C3 and Total LAI. # ################################################################################### # BONUS: Total LAI without C3, giving us the LAI for all woody species Total_C3 <- Per_yr %>% mutate(T_C = Total - C3_gr) #to write output dir.create("LAI_loose_scale") setwd("C:/Users/Documents/LAI_loose_scale/") dir.create("Total_C3") setwd("C:/Users/Documents/LAI_loose_scale/Total_C3/") for(i in Years){ All_sp <- Total_C3 %>% filter(year == i) xmin1 = range(All_sp$Lon)[1] xmax1 = range(All_sp$Lon)[2] ymin1 = range(All_sp$Lat)[1] ymax1 = range(All_sp$Lat)[2] #minimum and maximum value of the trends in grids, #alternately you can find the max/min values per taxa in lai.out and fix the scales for all years for each taxa: #min1 = min(Total_C3$T_C) --> (or 0); max1 = max(Total_C3$T_C) min1<-min(All_sp$T_C) max1<-max(All_sp$T_C) #create breaks and a colour pallete mybreaks<-seq(from=min1, to=max1, length.out = 10) pal11<-rev(colorRampPalette(brewer.pal(9, "Greens"))(length(mybreaks))) #plot the trends Total <- ggplot() + geom_sf(data = borders, fill= "white")+ labs(title =paste("Total LAIs, minus C3",i))+ geom_tile(data = All_sp, aes(fill=All_sp$T_C, x = Lon, y = Lat), colour=NA) + #borders data geom_sf(data = borders,fill = NA,colour = "dark green", size=0.2) + #limiting the borders data with the coordinates of the trend data coord_sf(xlim = c(xmin1, xmax1), ylim = c(ymin1, ymax1)) + #changing color scheme, modifying legend& scale bar scale_fill_gradientn(breaks=mybreaks, limits= c(min1, max1), label=as.character(round(mybreaks, digits=3)), colours=rev(pal11), guide = guide_colourbar(frame.colour = "black", barwidth = 1, barheight=25, ticks = F, title ="Total LAIs minus C3", title.position = "top"))+ #and giving uniform ticks to axes scale_x_longitude(ticks = 2) + scale_y_latitude(ticks = 2)+ # #modifying background of the plot theme(legend.position="right", legend.title= element_text(size = 25, color = "dark green"), legend.text = element_text(size = 22), panel.background = element_rect(fill = "#edf0f5"), plot.title = element_text(color="dark green", size=25, face="bold"), axis.line = element_line(colour = "black"), axis.text.x=element_text(size=24), axis.text.y = element_text(size=24), panel.border = element_rect(colour = "black",fill = NA, size = 0.5)) # print(Total) png(file=paste(i,"T_C.png"), width=1980, height=932) print(Total) dev.off() }