vignettes/plotting.Rmd
plotting.Rmd
#> tweaking rgl
Reproducing the examples in this vignette.
All spectra used in this manual are installed automatically with package hyperSpec.
Terminology Throughout the documentation of the package, the following terms are used:
- wavelength indicates any type of spectral abscissa.
- intensity indicates any type of spectral ordinate.
- extra data indicates non-spectroscopic data.
For some plots of the faux_cell
dataset, the pre-processed spectra and their cluster averages \(\pm\) one standard deviation are more suitable:
set.seed(1)
faux_cell <- generate_faux_cell()
faux_cell_preproc <- faux_cell - spc_fit_poly_below(faux_cell)
faux_cell_preproc <- faux_cell_preproc / rowMeans(faux_cell)
faux_cell_preproc <- faux_cell_preproc - quantile(faux_cell_preproc, 0.05)
cluster_cols <- c("dark blue", "orange", "#C02020")
cluster_meansd <- aggregate(faux_cell_preproc, faux_cell$region, mean_pm_sd)
cluster_means <- aggregate(faux_cell_preproc, faux_cell$region, mean)
Package hyperSpec comes with 6 major predefined plotting functions:
plot()
plot_spc()
plot_c()
plot_c()
is a package lattice function)
levelplot()
levelplot()
.
More details in sections 1.4, and 5.
plot_matrix()
plot_map()
levelplot()
for map or image plots.
More details in sections 1.5, 2.9, and 7.
(plot_map()
is a package lattice function)
plot_voronoi()
plot_map()
that produces Voronoi tesselations.
More details in sections 1.6, and 2.10.
(plot_voronoi()
is a package lattice function)
NOTE. Functions plot_map()
, plot_voronoi()
, and
levelplot()
are package lattice functions. Therefore, in loops,
functions, R Markdown chunks, etc. lattice objects need to be printed
explicitly by, e.g., print(plot_map(object))
(R
FAQ: Why do lattice/trellis graphics not work?).
plot_spc()
Function plot_spc()
plots the spectra (Fig. 1.1), i.e., the intensities ($spc
) over the wavelengths (@wavelength
).
plot_spc(flu)
plot_matrix()
Function plot_matrix()
plots the spectra, i.e., the colour coded intensities ($spc
) over the wavelengths (@wavelength
) and the row number (Fig. 1.2).
plot_matrix(flu)
plot_c()
Function plot_c()
plots an intensity over a single other data column, e.g., concentration (calibration plot, e.g., Fig. 1.3), depth (depth profile plot), or time (time series plot).
plot_c(flu)
#> Warning in plot_c(flu): Intensity at first wavelengh only is used.
The following warning is expected:
Intensity at first wavelengh only is used.
levelplot()
Function levelplot()
plots a false colour map, defined by a formula (Fig. 5.1).
levelplot(spc ~ x * y, faux_cell[, , 1200], aspect = "iso")
The following warning is expected:
Only first wavelength is used for plotting
.
plot_map()
Function plot_map()
is a specialized version of levelplot()
.
It uses a single value (e.g., average intensity or cluster membership) over two data columns (default are $x
and $y
) to plot a false colour map (Fig. @ref(fig:plot_map)).
plot_map(faux_cell[, , 1200])
plot_voronoi()
Function plot_voronoi()
is a special version of plot_map()
that produces Voronoi diagram of the hyperSpec
object (Fig. 1.5).
plot_voronoi(sample(faux_cell, 300), region ~ x * y)
plot()
The hyperSpec’s plot()
method uses its second argument to determine which of the specialized plots to produce.
This allows some handy abbreviations.
All further arguments are handed over to the function actually producing the plot.
"spc"
Command plot(x, "spc")
is equivalent to plot_spc(flu)
(Fig. 2.1).
plot(flu, "spc")
"spcmeansd"
Command plot(x, "spcmeansd")
plots mean spectrum \(\pm\) 1 standard deviation (Fig. 2.2).
plot(faux_cell_preproc, "spcmeansd")
"spcprctile"
Code plot(x, "spcprctile")
plots median, 16th and 84th percentile for each wavelength (Fig. 2.3).
For Gaussian distributed data, 16th, 50th and 84th percentile are equal to mean \(\pm\) standard deviation.
Spectroscopic data frequently are not Gaussian distributed.
The percentiles give a better idea of the true distribution.
They are also less sensitive to outliers.
plot(faux_cell_preproc, "spcprctile")
"spcprctl5"
The code plot(x, "spcprctl5")
is like "spcprctile"
plus 5th and 95th percentile (Fig. 2.4).
plot(faux_cell_preproc, "spcprctl5")
"c"
Command plot(x, "c")
is equivalent to plot_c(flu)
(Fig. 2.5).
plot(flu, "c")
#> Warning in plot_c(x, ...): Intensity at first wavelengh only is used.
"ts"
Function plot(x, "ts")
plots a time series plot and is equivalent to plot_c(laser, spc ~ t)
(Fig. 2.6).
plot(laser[, , 405], "ts")
"depth"
Code plot(x, "depth")
plots a depth profile plot and is the same as plot_c(laser, spc ~ z)
(Fig. 2.7).
depth.profile <- new("hyperSpec",
spc = as.matrix(rnorm(20) + 1:20),
data = data.frame(z = 1:20),
labels = list(
spc = "I / a.u.",
z = expression(`/`(z, mu * m)),
.wavelength = expression(lambda)
)
)
plot(depth.profile, "depth")
"mat"
Code plot(x, "mat")
plots the spectra matrix (Fig. ??).
plot(laser, "mat")
The code is equivalent to:
plot_matrix(laser)
A lattice alternative is:
levelplot(spc ~ .wavelength * .row, data = laser)
"map"
Code plot(x, "map")
is equivalent to plot_map(faux_cell)
(Fig. 2.8).
plot(faux_cell[, , 1200], "map")
"voronoi"
Use plot(x, "voronoi")
for a Voronoi plot (Fig. @ref(fig:plot_voronoi)).
See ?plot_voronoi
and ?latticeExtra::panel.voronoi
.
plot_spc()
Function plot_spc()
offers a variety of parameters for customized plots.
If only one wavelength range is needed (Fig. 3.1), the extract command is handiest.
plot_spc(paracetamol[, , 700 ~ 1200])
Numbers connected with the tilde (~
) are interpreted as having the same units as the wavelengths.
If wl.range
already contains indices use wl.index = TRUE
.
For more details refer to vignette("hyperspec", package = "hyperSpec")
.
To plot severel wavelength ranges (Fig. 3.2), use wl.range = list (600 ~ 1800, 2800 ~ 3100)
.
Cut the wavelength axis appropriately with xoffset = 750
.
If available, the package package plotrix[3,4] is used to produce the cut mark.
To create a plot with reversed abscissa (Fig. 3.3), use wl.reverse = TRUE
plot_spc(paracetamol, wl.reverse = TRUE)
To have different colours of spectra (Fig. 3.4), use col = vector_of_colours
.
plot_spc(flu, col = palette_matlab_dark(6))
To plot dots instead of lines (Fig. 3.5), use, e.g., lines.args = list (pch = 20, type = "p")
To plot additional spectra onto an existing plot (Fig. 3.7), use add = TRUE
Argument func
may be used to calculate summary characteristics before plotting.
To plot, e.g., the standard deviation of the spectra, use plot_spc(..., func = sd)
(Fig. 3.8).
plot_spc(faux_cell_preproc, func = sd)
To plot a different line at zero intensity \((I = 0)\), argument zeroline
may be used (Fig. 3.9).
Argument zeroline
takes a list with parameters that are passed to function abline()
, NA
suppresses the line.
Function plot_spc()
uses base graphics.
After plotting the spectra, more content may be added to the graphic by abline()
, lines()
, points()
, etc. (Fig. 3.10).
plot(laser, "spcmeansd")
abline(
v = c(405.0063, 405.1121, 405.2885, 405.3591),
col = c("black", "blue", "red", "darkgreen")
)
To stack spectra, use stacked = TRUE
(Fig. 3.11).
plot_spc(
cluster_means,
col = cluster_cols,
stacked = TRUE
)
The spectra to be stacked can be grouped: stacked = "factor"
.
Alternatively, the name of an extra data column can be used for grouping (Fig. 3.12).
op <- par(las = 1, mgp = c(3.1, .7, 0))
plot(
cluster_meansd,
stacked = ".aggregate",
fill = ".aggregate",
col = cluster_cols
)
yoffset
Stacking values can also be given manually as numeric values in yoffset
(Fig. 3.13).
It is possible to obtain a denser stacking (Fig. 3.14).
yoffsets <- apply(cluster_means[[]], 2, diff)
yoffsets <- -apply(yoffsets, 1, min)
plot(cluster_means,
yoffset = c(0, cumsum(yoffsets)),
col = cluster_cols
)
Function plot_spc()
allows fine grained customization of almost all aspects of the plot (see example in Fig. 3.15).
This is possible by giving arguments to the functions that actually perform the plotting plot()
for setting up the plot area, lines()
for the plotting of the lines, axis()
for the axes, etc.
The arguments for these functions should be given in lists as plot.args
, lines.args
, axis.args
, etc.
yoffset <- apply(faux_cell_preproc, 2, quantile, c(0.05, 0.95))
yoffset <- range(yoffset)
plot(faux_cell_preproc[1],
plot.args = list(ylim = c(0, 2) * yoffset),
lines.args = list(type = "n")
)
yoffset <- (0:1) * diff(yoffset)
for (i in 1:3) {
plot(faux_cell_preproc, "spcprctl5",
yoffset = yoffset[i],
col = "gray", add = TRUE
)
plot(faux_cell_preproc[i],
yoffset = yoffset[i],
col = palette_matlab_dark(3)[i], add = TRUE,
lines.args = list(lwd = 2)
)
}
plot_c()
Spectra intensities of one wavelength can be plotted over the concentration for univariate calibration (Fig. 4.1).
plot_c(flu[, , 450])
The default is to use the first intensity only.
A function to compute a summary of the intensities before the drawing can be used via argument func
(Fig. 4.2).
plot_c(flu, func = range, groups = .wavelength)
If func()
returns more than one value, the different results are accessible by .wavelength
.
Lattice conditioning (operator |
) can be used to plot more traces separately (Fig. 4.3).
plot_c(flu[, , c(405, 445)], spc ~ c | .wavelength,
cex = .3, scales = list(alternating = c(1, 1))
)
Argument groups
may be used as a grouping parameter to plot more traces in one panel (Fig. 4.4).
Arguments of lattice function xyplot()
can be given to plot_c()
(Fig. 4.5).
plot_c(flu[, , 450],
ylab = expression(I["450 nm"] / a.u.),
xlim = range(0, flu$c + .01),
ylim = range(0, flu$spc + 10),
pch = 4
)
As plot_c()
uses the package lattice function xyplot()
, additions to the plot must be made via the panel function (Fig. 4.6).
panelcalibration <- function(x, y, ..., clim = range(x), level = .95) {
panel.xyplot(x, y, ...)
lm <- lm(y ~ x)
panel.abline(coef(lm), ...)
cx <- seq(clim[1], clim[2], length.out = 50)
cy <- predict(lm, data.frame(x = cx),
interval = "confidence",
level = level
)
panel.lines(cx, cy[, 2], col = "gray")
panel.lines(cx, cy[, 3], col = "gray")
}
Abscissae other than c
may be specified by explicitly giving the model formula (Fig. 4.7).
plot_c(laser[, , c(405.0063, 405.1121, 405.2885, 405.3591)],
spc ~ t,
groups = .wavelength,
type = "b",
col = c("black", "blue", "red", "darkgreen")
)
levelplot()
Package hyperSpec’s function levelplot()
can use two special column names:
Besides that, it behaves exactly like lattice levelplot()
.
Particularly, the data is given as the second argument (Fig. 5.1).
levelplot(spc ~ x * y, data = faux_cell[, , 800])
If the colour-coded value is a factor, the display is adjusted to this fact (Fig. 5.2).
levelplot(region ~ x * y, data = faux_cell)
plot_matrix()
It is often useful to plot the spectra against an additional coordinate, e.g., the time for time
series, the depth for depth profiles, etc.
This can be done by plot(object, "mat")
.
The actual plotting is done by image()
, but levelplot()
can produce spectra matrix plots as well and these plots can be grouped or conditioned.
Argument col
can be used to provide a different colour palette (Fig. 6.1).
plot(laser, "mat", col = heat.colors(20))
This is the same as:
plot_matrix(laser, col = heat.colors(20))
Different extra data column can be used as y-axis (Fig. 6.2).
plot_matrix(laser, y = "t")
Alternatively, y values and axis label can be given separately.
plot_matrix(laser, y = laser$t, ylab = labels(laser, "t"))
Contour lines may also be added (Fig. 6.3).
plot_matrix(flu, col = palette_matlab_dark(20))
plot_matrix(flu, col = "white", contour = TRUE, add = TRUE)
In levelplot()
, colour-coded points may be set via special panel function (Fig. 6.4).
library("latticeExtra")
barb <- collapse(barbiturates)
barb <- wl_sort(barb)
levelplot(spc ~ .wavelength * z, barb,
panel = panel.levelplot.points,
cex = .33, col.symbol = NA,
col.regions = palette_matlab
)
plot_map()
Function plot_map()
is a specialized version of levelplot()
(Fig. 7.3).
The spectral intensities may be summarized by a function before plotting (default: mean()
).
The same scale is used for x and y axes (aspect = "iso"
).
plot_map(faux_cell[, , 1200])
Specify the colour-coded variable, abscissa and ordinate as formula: colour.coded ~ abscissa * ordinate
(Fig. 7.2).
plot_map(faux_cell[, , 1200], spc ~ y * x)
The plotting of colour maps is done via R package lattice (aka Trellis graphic approach), which is highly customizable.
Use function trellis.par.get()
and trellis.par.set()
to get/set the settings for the current graphics device.
my_theme <- trellis.par.get()
names(my_theme) # note how many parameters are tunable
#> [1] "grid.pars" "fontsize" "background" "panel.background"
#> [5] "clip" "add.line" "add.text" "plot.polygon"
#> [9] "box.dot" "box.rectangle" "box.umbrella" "dot.line"
#> [13] "dot.symbol" "plot.line" "plot.symbol" "reference.line"
#> [17] "strip.background" "strip.shingle" "strip.border" "superpose.line"
#> [21] "superpose.symbol" "superpose.polygon" "regions" "shade.colors"
#> [25] "axis.line" "axis.text" "axis.components" "layout.heights"
#> [29] "layout.widths" "box.3d" "par.title.text" "par.xlab.text"
#> [33] "par.ylab.text" "par.zlab.text" "par.main.text" "par.sub.text"
Any of these parameters can be fine-tuned to produce the desired output.
For example, parameter my_theme$region
is responsible for the appearance of color maps, and it contains elements $alpha
and $col
.
By changing these parameters you can create your own theme for plotting and pass it to the plotting function via par.settings
.
Fig. 8.1 uses a customized lattice theme.
my_theme$regions$col <- grDevices::terrain.colors
plot_map(faux_cell[, , 1200], par.settings = my_theme)
It is possible to persistently (i.e. inside of the current R session) set lattice parameters, so they would apply to all further plots.
This is done via a call to trellis.par.set()
, for example trellis.par.set(my_theme)
.
The current settings can be visualized via a call to show.settings()
Graphical parameters for trellis may be displayed via show.settings()
(Fig. 8.2).
# Display current trellis parameters
show.settings()
show.settings(my_theme)
An overview of different colour palettes, and ways to create your own, can be found in the R colour cheatsheet.
A map of the average intensity at particular wavelengths can be plotted if the wavelengths of interested are explicitly extracted (Fig. 8.3).
Logical conditions may be used to create subplots (Fig. 8.4).
plot_map(
faux_cell[, , 1500],
spc ~ y * x | x > 5,
col.regions = palette_matlab(20)
)
.wavelength
Function plot_map()
automatically applies the function in func
before plotting.
Argument func
defaults to function mean()
.
In order to suppress this, use func = NULL
.
This allows conditioning on the wavelengths.
Fig. 8.5 demonstrates an example to plot the maps of principal components scores.
pca <- prcomp(~spc, data = faux_cell_preproc$.)
scores <- decomposition(faux_cell, pca$x,
label.wavelength = "PC",
label.spc = "score / a.u."
)
plot_map(
scores[, , 1:3],
spc ~ y * x | as.factor(.wavelength),
func = NULL,
col.regions = palette_matlab(20)
)
Alternatively, use levelplot()
directly (Fig. 8.6).
levelplot(
spc ~ y * x | as.factor(.wavelength),
scores[, , 1:3],
aspect = "iso",
col.regions = palette_matlab(20)
)
Fig. 8.7 shows an example of a Voronoi plot.
Voronoi uses panel.voronoi()
from package latticeExtra[5].
plot_voronoi(
sample(faux_cell, 300), region ~ x * y,
col.regions = palette_matlab(20)
)
If the spectra come from a rectangular grid, missing positions can be marked with the following panel function:
mark.missing <- function(x, y, z, ...) {
panel.levelplot(x, y, z, ...)
miss <- expand.grid(x = unique(x), y = unique(y))
miss <- merge(miss, data.frame(x, y, TRUE), all.x = TRUE)
miss <- miss[is.na(miss[, 3]), ]
panel.xyplot(miss[, 1], miss[, 2], pch = 4, ...)
}
Fig. 8.8 shows the result.
plot_map(sample(faux_cell[, , 1200], length(faux_cell) - 20),
col.regions = palette_matlab(20),
col = "black",
panel = mark.missing
)
The panel function used by plot_map()
defaults to panel.levelplot.raster()
which assumes an evenly spaced measurement grid.
Even if the spectra are measured on a nominally evenly spaced grid, the actual stage position may slightly vary due to positioning inaccuracy and some manufacturers (e.g., Kaiser) record the position reported by the stage rather than the position requested by the stage control.
This leads to weird-looking output with holes, and possibly wrong columns (Fig. 8.9).
uneven <- faux_cell[, , 1200]
uneven$x <- uneven$x + round(rnorm(nrow(uneven), sd = 0.05), digits = 1)
uneven$y <- uneven$y + round(rnorm(nrow(uneven), sd = 0.05), digits = 1)
plot_map(uneven)
#> Warning in (function (x, y, z, subscripts, at = pretty(z), ..., col.regions = regions$col, : 'x'
#> values are not equispaced; output may be wrong
#> Warning in (function (x, y, z, subscripts, at = pretty(z), ..., col.regions = regions$col, : 'y'
#> values are not equispaced; output may be wrong
The symptom of this situation are warnings about values in x and/or y not being equispaced; and that the output, therefore, may be wrong.
One possibility to obtain a correct map is using plot_voronoi()
instead which will construct a mosaic-like image with the respective “pixel” areas being centred around the actually recorded $x
and $y
position (Fig. 8.10).
plot_voronoi(uneven, backend = "deldir")
Another possibility that underlines a point shape of the measurements is switching to latticeExtra::panel.levelplot.points()
(Fig. 8.11).
plot_map(
uneven,
panel = panel.levelplot.points,
cex = 0.75,
col.symbol = NA
)
Alternatively, the measurement raster positions can be rounded to their nominal raster (e.g. Fig. 8.12).
rx <- raster_make(uneven$x, startx = -11.55, d = 1, tol = 0.3)
uneven$x <- rx$x
ry <- raster_make(uneven$y, startx = -4.77, d = 1, tol = 0.3)
uneven$y <- ry$x
plot_map(uneven)
Package rgl[6] offers fast 3d plotting in R. As package rgl’s axis annotations are sometimes awkward, they may better be set manually:
library(rgl)
laser <- laser[, , 404.8 ~ 405.6] / 10000
laser$t <- laser$t / 3600
cols <- rep(palette_matlab(nrow(laser)), nwl(laser))
surface3d(
y = wl(laser), x = laser$t,
z = laser$spc, col = cols
)
aspect3d(c(1, 1, 0.25))
axes3d(c("x+-", "y--", "z--"))
axes3d("y--", nticks = 25, labels = FALSE)
mtext3d("t / h", "x+-", line = 2.5)
mtext3d("lambda / nm", "y--", line = 2.5)
mtext3d("I / a.u.", edge = "z--", line = 2.5)
Package hyperSpec offers basic interaction, identify_spc()
for spectra plots, and map.identify()
and map.sel.poly()
for maps.
The first two identify points in spectra plots and map plots, respectively.
Function map.sel.poly()
selects the part of a hyperSpec
object that lies inside the user defined polygon.
identify_spc()
Finding Out Wavelength, Intensity and Spectrum
Function identify_spc()
allows to measure points in graphics produced by plot_spc()
.
It works correctly with reversed and cut wavelength axes.
identify_spc(plot_spc(paracetamol, wl.range = c(600 ~ 1800, 2800 ~ 3200), xoffset = 800))
The result is a data.frame with the indices of the spectra, the wavelength, and its intensity.
map.identify()
finding a spectrum in a map plot
Function map.identify()
returns the spectra indices of the clicked points.
map.identify(faux_cell[, , 1200])
map.sel.poly()
selecting spectra inside a polygon in a map plot
Function map.sel.poly()
returns a logical indicating which spectra are inside the polygon drawn by the user:
map.sel.poly(faux_cell[, , 1200])
For base graphics (as produced by plot_spc()
), locator()
may be useful as well.
It returns the clicked coordinates. Note that these are not transformed according to xoffset
& Co.
For lattice graphics, grid.locator()
may be used instead.
If it is not called in the panel function, a preceding call to trellis.focus()
is needed:
plot(laser, "mat")
trellis.focus()
grid.locator()
Function identify()
(or panel.identify()
for lattice graphics) allows to identify points of the plot directly.
Note that the returned indices correspond to the plotted object.
Methods plot_map
, plot_voronoi
, levelplot
, and plot_c
use package lattice functions.
Therefore, in loops, functions, Sweave, R Markdown chunks, etc. the lattice object needs to be printed explicitly by print(plot_map(object))
(R FAQ: Why do lattice/trellis graphics not work?).
The same holds for package ggplot2 graphics.
sessioninfo::session_info("hyperSpec")
#> ─ Session info ───────────────────────────────────────────────────────────────────────────────────
#> setting value
#> version R version 4.4.0 (2024-04-24)
#> os Ubuntu 22.04.4 LTS
#> system x86_64, linux-gnu
#> ui X11
#> language en
#> collate C.UTF-8
#> ctype C.UTF-8
#> tz UTC
#> date 2024-05-27
#> pandoc 3.1.11 @ /opt/hostedtoolcache/pandoc/3.1.11/x64/ (via rmarkdown)
#>
#> ─ Packages ───────────────────────────────────────────────────────────────────────────────────────
#> ! package * version date (UTC) lib source
#> brio 1.1.5 2024-04-24 [2] RSPM
#> callr 3.7.6 2024-03-25 [2] RSPM
#> cli 3.6.2 2023-12-11 [2] RSPM
#> colorspace 2.1-0 2023-01-23 [2] RSPM
#> crayon 1.5.2 2022-09-29 [2] RSPM
#> deldir 2.0-4 2024-02-28 [2] RSPM
#> desc 1.4.3 2023-12-10 [2] RSPM
#> diffobj 0.3.5 2021-10-05 [2] RSPM
#> digest 0.6.35 2024-03-11 [2] RSPM
#> dplyr 1.1.4 2023-11-17 [2] RSPM
#> evaluate 0.23 2023-11-01 [2] RSPM
#> fansi 1.0.6 2023-12-08 [2] RSPM
#> farver 2.1.2 2024-05-13 [2] RSPM
#> fs 1.6.4 2024-04-25 [2] RSPM
#> generics 0.1.3 2022-07-05 [2] RSPM
#> ggplot2 * 3.5.1 2024-04-23 [2] RSPM
#> glue 1.7.0 2024-01-09 [2] RSPM
#> gtable 0.3.5 2024-04-22 [2] RSPM
#> hyperSpec * 0.200.0.9000 2024-05-27 [1] local
#> hySpc.testthat 0.2.1 2020-06-24 [2] RSPM
#> interp 1.1-6 2024-01-26 [2] RSPM
#> isoband 0.2.7 2022-12-20 [2] RSPM
#> jpeg 0.1-10 2022-11-29 [2] RSPM
#> jsonlite 1.8.8 2023-12-04 [2] RSPM
#> labeling 0.4.3 2023-08-29 [2] RSPM
#> lattice * 0.22-6 2024-03-20 [4] CRAN (R 4.4.0)
#> latticeExtra * 0.6-30 2022-07-04 [2] RSPM
#> lazyeval 0.2.2 2019-03-15 [2] RSPM
#> lifecycle 1.0.4 2023-11-07 [2] RSPM
#> magrittr 2.0.3 2022-03-30 [2] RSPM
#> MASS 7.3-60.2 2024-05-25 [4] local
#> Matrix 1.7-0 2024-03-22 [4] CRAN (R 4.4.0)
#> mgcv 1.9-1 2023-12-21 [4] CRAN (R 4.4.0)
#> munsell 0.5.1 2024-04-01 [2] RSPM
#> nlme 3.1-164 2023-11-27 [4] CRAN (R 4.4.0)
#> pillar 1.9.0 2023-03-22 [2] RSPM
#> pkgbuild 1.4.4 2024-03-17 [2] RSPM
#> pkgconfig 2.0.3 2019-09-22 [2] RSPM
#> pkgload 1.3.4 2024-01-16 [2] RSPM
#> png 0.1-8 2022-11-29 [2] RSPM
#> praise 1.0.0 2015-08-11 [2] RSPM
#> processx 3.8.4 2024-03-16 [2] RSPM
#> ps 1.7.6 2024-01-18 [2] RSPM
#> R6 2.5.1 2021-08-19 [2] RSPM
#> RColorBrewer 1.1-3 2022-04-03 [2] RSPM
#> Rcpp 1.0.12 2024-01-09 [2] RSPM
#> R RcppEigen <NA> <NA> [?] <NA>
#> rematch2 2.1.2 2020-05-01 [2] RSPM
#> rlang 1.1.3 2024-01-10 [2] RSPM
#> rprojroot 2.0.4 2023-11-05 [2] RSPM
#> scales 1.3.0 2023-11-28 [2] RSPM
#> testthat 3.2.1.1 2024-04-14 [2] RSPM
#> tibble 3.2.1 2023-03-20 [2] RSPM
#> tidyselect 1.2.1 2024-03-11 [2] RSPM
#> utf8 1.2.4 2023-10-22 [2] RSPM
#> vctrs 0.6.5 2023-12-01 [2] RSPM
#> viridisLite 0.4.2 2023-05-02 [2] RSPM
#> waldo 0.5.2 2023-11-02 [2] RSPM
#> withr 3.0.0 2024-01-16 [2] RSPM
#> xml2 1.3.6 2023-12-04 [2] RSPM
#>
#> [1] /tmp/RtmpxUEHmP/temp_libpath17ce4b3a89e0
#> [2] /home/runner/work/_temp/Library
#> [3] /opt/R/4.4.0/lib/R/site-library
#> [4] /opt/R/4.4.0/lib/R/library
#>
#> R ── Package was removed from disk.
#>
#> ──────────────────────────────────────────────────────────────────────────────────────────────────