smooth PRC diagram of treatment and plot lines or points without species loading

Source: R/plot_smoothPRC_cdt.R

plot_smoothPRC_cdt creates a PRC diagram of the treatments without species loadings.

plot_smoothPRC_cdt(
  object,
  axis = 1,
  mod_prc = NULL,
  with_classical_lines = TRUE,
  with_unconstrained_scores = TRUE,
  flip = FALSE,
  sample_times_only = FALSE
)

Arguments

object

a result of smoothPRC.

axis

the number of the axis to plot

mod_prc

a results of doPRC, i.e. a classical PRC with the same data and design. Default NULL, in which the classical PRC is computed using doPRC and the data in object.

with_classical_lines

logical to add or suppress classical PRC lines. Default TRUE, i.e. with classical PRC lines.

with_unconstrained_scores

logical to add or suppress the unconstrained scores in the smooth PRC showing the scatter of the individual samples (cosms). Default TRUE, i.e. with individual points.

flip

logical. Should the axis be reversed? Default FALSE. Can be numeric with -1 meaning: reverse, and 1 meaning do not reverse.

sample_times_only

logical to PRClines based on the sampled times only, i.e. without interpolation. Default FALSE.

Examples


## ---- Packages ----

suppressPackageStartupMessages({
  library(LMMsolver)
  library(dplyr)
  library(ggplot2)
  library(permute)
  library(PRC)
})

#the chlorpyrifos experiment from van den Brink & ter Braak 1999
data(pyrifos, package = "vegan") #log-transformed species data from package vegan
Y <- pyrifos
Design <- data.frame(Time=gl(11, 12, labels=c(-4, -1, 0.1, 1, 2, 4, 8, 12, 15, 19, 24)),
                     Treatment=factor(rep(c(0.1, 0, 0, 0.9, 0, 44, 6, 0.1, 44, 0.9, 0, 6), 11)),
                     cosm = gl(12, 1, length=132))
mod_prc <- doPRC(pyrifos ~ Treatment:Time + Condition(Time),  data = Design)
#> 
#> Some constraints or conditions were aliased because they were redundant. This
#> can happen if terms are constant or linearly dependent (collinear):
#> 'Treatment44:Time-4', 'Treatment44:Time-1', 'Treatment44:Time0.1',
#> 'Treatment44:Time1', 'Treatment44:Time2', 'Treatment44:Time4',
#> 'Treatment44:Time8', 'Treatment44:Time12', 'Treatment44:Time15',
#> 'Treatment44:Time19', 'Treatment44:Time24'
cntr <- how(plots = Plots(strata =Design$cosm,type = "free"),
            within = Within(type = "none"), nperm = 99)

anova(mod_prc, permutations = cntr)
#> Permutation test for rda under reduced model
#> Plots: Design$cosm, plot permutation: free
#> Permutation: none
#> Number of permutations: 99
#> 
#> Model: rda(formula = pyrifos ~ Treatment:Time + Condition(Time), data = data, scale = scale)
#>          Df Variance     F Pr(>F)  
#> Model    44   96.684 1.312   0.03 *
#> Residual 77  128.959               
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


b_mod_prc <- vegan::scores(mod_prc,choices= 1, display = "sp")

smoothPRC_model  <- set_smoothPRC_model(data = Design)
smoothPRC_model$spline
#> ~spl2D(time, dose, nseg = c(6, 6), degree = 3, pord = 2)
#> <environment: 0x0000028814ec1d10>
smooth_PRC <- smoothPRC(Y, lmm_model = smoothPRC_model)
summary(smooth_PRC$obj)
#>       Length Class    Mode
#> axis1 32     LMMsolve list
#> axis2 32     LMMsolve list
cor(smooth_PRC$B[,"B1"], b_mod_prc)
#>           RDA1
#> [1,] 0.9887614

cntr$nperm <- 19
an <- anova(smooth_PRC, permutations = cntr, verbose = TRUE)
#> iteration 10 out of  19 
#> [1] 0.1
an$pval
#> [1] 0.1

#graph
plot_species_scores_bk(smooth_PRC$B, threshold = 20, scoresname = "B1")
#> 20  species with names out of 178 species

plot_smoothPRC_cdt(smooth_PRC,flip = FALSE,sample_times_only = TRUE)


plotPRC(mod_prc)
#> The analysis is based on model
#>  pyrifos~Treatment:Time + Condition(Time) 
#> The current score, condition and treatment(s) are:
#>  score    (y-axis): PRC1 
#> condition (x-axis): Time 
#> treatments (color ): Treatment 
#> Set these arguments yourself if you wish otherwise.
#> 34  species with names out of 178 species

plotsmoothPRC(smooth_PRC, flip = FALSE, threshold = 1)
#> 122  species with names out of 178 species


# illustration that smoothPRC extends doPRC
classical_PRC_model <- set_smoothPRC_model(
  fixed= Yb ~ Treatment:Time,
  fixedH0 = Yb ~Time,
  spline = NULL,
  splineH0 = NULL,
  data = Design)
classical_PRC <- smoothPRC(Y, lmm_model = classical_PRC_model)
plotsmoothPRC(classical_PRC, flip = FALSE, threshold = 1, title ="Classical PRC via smoothPRC")
#> 122  species with names out of 178 species

# the dotted and solid lines overlap