This function produces figures that illustrate statistical correlations between temperature variation during certain phases and the timing of phenological event, based on a PLS analysis conducted with the PLS_pheno or the PLS_chill_force function.
plot_PLS(
PLS_output,
PLS_results_path,
VIP_threshold = 0.8,
colorscheme = "color",
plot_bloom = TRUE,
fonttype = "serif",
add_chill = c(NA, NA),
add_heat = c(NA, NA),
plot_titles_Temp = "Mean temperature",
plot_titles_chill_force = c("Chill Accumulation", "Heat Accumulation"),
axis_labels_Temp = expression("Mean temperature ("^"o" * "C)"),
axis_labels_chill_force = c("Chill Portions per day", "GDH per day"),
chill_force_same_scale = TRUE
)a PLS_output object - the output of the PLS_pheno function. This object is a list with a list element called PLS_summary (and an optional object called PLS_output). This element is a data.frame with the following columns: Date, JDay, Coef, VIP, Tmean, Tstdev. Date is the day of the year in MDD format. JDay is the Julian day (day of the year) of the year in which the biological event is observed; since the analysis will often start in the year before the event, this column often starts with negative numbers. Coef is the coefficient of the PLS regression output. VIP is the Variable Importance in the Projection, another output of the PLS regression. Tmean is the mean observed temperature of the respective day of the year, for the duration of the phenology record. Tstdev is the standard deviation of temperature on a given day of the year over the length of the phenology record.
the path where analysis outputs should be saved. Should include the file name, but without suffix.
the VIP threshold, above which a variable is considered important. Defaults to 0.8.
color scheme used for plotting. For grayscale image, this should be set to "bw". Otherwise a color plot is produced.
boolean variable specifying whether the range of bloom dates should be shown in the plots. If set to TRUE, this range is shown by a semi-transparent gray rectangle. The median bloom date is shown as a dashed line. This only works if the full range of bloom dates is visible in the plot, and it should be set to FALSE if anything other than Julian dates are used as dependent variables.
font style to be used for the figure. Can be 'serif' (default) or 'sans'.
option for indicating the chilling period in the plot. This should be a numeric vector: c(start_chill,end_chill).
option for indicating the forcing period in the plot. This should be a numeric vector: c(start_heat,end_heat).
title for the bottom plot, which relates PLS outputs to values of the input variable (temperature in the original version). Only affects the output for PLS_Temp_pheno objects.
titles for the bottom plots, which relate PLS outputs to values of the input variables (chill and heat accumulation). Only affects the output for PLS_chillforce_pheno objects.
y-axis label for the bottom plot, which relates PLS outputs to values of the input variable (temperature in the original version). Only affects the output for PLS_Temp_pheno objects.
y-axis labels for the bottom plots, which relate PLS outputs to values of the input variables (chill and heat accumulation). Only affects the output for PLS_chillforce_pheno objects.
Boolean parameter indicating whether the two sets of VIP scores and model coefficients resulting from a PLS_chillforce_pheno analysis should be shown on the same scale in the separate output diagrams. Since this is generally advisable for comparison, this defaults to TRUE.
Eike Luedeling
Ths figure illustrates results from the PLS_pheno function, which uses Partial Least Squares (or Projection to Latent Structures) regression to examine the relationship between mean daily temperatures and the timing of an annual biological event. It produces a plot (as a bmp image) with three panels: the top panel shows the value of the VIP score for each day of the year; the middle panel shows the model coefficients and the bottom panel shows the mean temperature and its standard deviation. In the top plot, all days with VIP scores above VIP_threshold are shown in blue. In the other two panels, values for the same days are shown in red, which high VIP scores coincide with negative model coefficients, and in green for positive coefficients. This function does not produce an output, but as side effects it produces a bmp image and a table that summarizes all data used for making the figure in the specified folder.
The method is described here:
Luedeling E and Gassner A, 2012. Partial Least Squares Regression for analyzing walnut phenology in California. Agricultural and Forest Meteorology 158, 43-52.
Wold S, 1995. PLS for multivariate linear modeling. In: van der Waterbeemd H (ed) Chemometric methods in molecular design: methods and principles in medicinal chemistry, vol 2. Chemie, Weinheim, pp 195-218.
Wold S, Sjostrom M, Eriksson L, 2001. PLS-regression: a basic tool of chemometrics. Chemometr Intell Lab 58(2), 109-130.
Mevik B-H, Wehrens R, Liland KH, 2011. PLS: Partial Least Squares and Principal Component Regression. R package version 2.3-0. http://CRAN.R-project.org/package0pls.
Some applications:
Guo L, Dai J, Wang M, Xu J, Luedeling E, 2015. Responses of spring phenology in temperate zone trees to climate warming: a case study of apricot flowering in China. Agricultural and Forest Meteorology 201, 1-7.
Luedeling E, Kunz A and Blanke M, 2013. Identification of chilling and heat requirements of cherry trees - a statistical approach. International Journal of Biometeorology 57,679-689.
Yu H, Luedeling E and Xu J, 2010. Stronger winter than spring warming delays spring phenology on the Tibetan Plateau. Proceedings of the National Academy of Sciences (PNAS) 107 (51), 22151-22156.
Yu H, Xu J, Okuto E and Luedeling E, 2012. Seasonal Response of Grasslands to Climate Change on the Tibetan Plateau. PLoS ONE 7(11), e49230.
weather<-fix_weather(KA_weather[which(KA_weather$Year>2004),])
#Plots look much better with weather<-fix_weather(KA_weather)
#but that takes to long to run for passing CRAN checks
PLS_results<-PLS_pheno(
weather_data=weather$weather,
split_month=6, #last month in same year
bio_data=KA_bloom)
PLS_results_path<-paste(getwd(),"/PLS_output",sep="")
#plot_PLS(PLS_results,PLS_results_path)
#plot_PLS(PLS_results,PLS_results_path,add_chill=c(307,19),add_heat=c(54,109))
dc<-daily_chill(stack_hourly_temps(weather,50.4), 11)
plscf<-PLS_chill_force(daily_chill_obj=dc, bio_data_frame=KA_bloom, split_month=6)
#plot_PLS(plscf,PLS_results_path)
#plot_PLS(plscf,PLS_results_path,add_chill=c(307,19),add_heat=c(54,109))
Run the code above in your browser using DataLab