diff --git a/.Rbuildignore b/.Rbuildignore index 5af5293d..f1d67fcb 100644 --- a/.Rbuildignore +++ b/.Rbuildignore @@ -8,9 +8,11 @@ ^Meta$ ^cran-comments\.md$ ^CRAN-RELEASE$ +^CRAN-SUBMISSION$ ^_pkgdown\.yml$ ^docs$ ^pkgdown$ ^codecov\.yml$ -^vignettes/spectragryph.* ^vignettes/advanced.* +^vignettes/app.* +^vignettes/spectragryph.* diff --git a/CRAN-SUBMISSION b/CRAN-SUBMISSION new file mode 100644 index 00000000..e7ae3d08 --- /dev/null +++ b/CRAN-SUBMISSION @@ -0,0 +1,3 @@ +Version: 1.0.0 +Date: 2023-09-01 14:53:15 UTC +SHA: 1924ab057e7bd36723eccdd4649204fef0b2366c diff --git a/DESCRIPTION b/DESCRIPTION index d13ebc54..43b56724 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -2,7 +2,7 @@ Package: OpenSpecy Type: Package Title: Analyze, Process, Identify, and Share Raman and (FT)IR Spectra Version: 1.0.0 -Date: 2023-08-29 +Date: 2023-09-01 Authors@R: c(person("Win", "Cowger", role = c("cre", "aut", "dtc"), email = "wincowger@gmail.com", comment = c(ORCID = "0000-0001-9226-3104")), @@ -31,13 +31,16 @@ Authors@R: c(person("Win", "Cowger", role = c("cre", "aut", "dtc"), person("Mary C", "Norris", role = c("ctb")), person("Christine M", "Knauss", role = c("ctb"), comment = c(ORCID = "0000-0003-4404-8922")), - person("Aleksandra","Karapetrova", role = c("ctb", "dtc", "rev"), + person("Aleksandra", "Karapetrova", role = c("ctb", "dtc", "rev"), comment = c(ORCID = "0000-0002-9856-1644")), - person("Vesna","Teofilovic", role = c("ctb"), + person("Vesna", "Teofilovic", role = c("ctb"), comment = c(ORCID = "0000-0002-3557-1482")), - person("Laura A. T.","Markley", role = c("ctb"), + person("Laura A. T.", "Markley", role = c("ctb"), comment = c(ORCID = "0000-0003-0620-8366")), - person("Shreyas","Patankar", role = c("ctb", "dtc"))) + person("Shreyas", "Patankar", role = c("ctb", "dtc")), + person("Katherine", "Lasdin", role = c("ctb")), + person("National Renewable Energy Laboratory", role = c("fnd")), + person("Possibility Lab", role = c("fnd"))) Description: Raman and (FT)IR spectral analysis tool for plastic particles and other environmental samples (Cowger et al. 2021, ). With read_any(), Open Specy provides a diff --git a/NAMESPACE b/NAMESPACE index e62d5632..90eaf4d6 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -4,6 +4,8 @@ S3method(adj_intens,OpenSpecy) S3method(adj_intens,default) S3method(ai_classify,OpenSpecy) S3method(ai_classify,default) +S3method(as.data.frame,OpenSpecy) +S3method(as.data.table,OpenSpecy) S3method(as_OpenSpecy,OpenSpecy) S3method(as_OpenSpecy,data.frame) S3method(as_OpenSpecy,default) @@ -32,6 +34,8 @@ S3method(heatmap_spec,default) S3method(interactive_plot,OpenSpecy) S3method(interactive_plot,default) S3method(lines,OpenSpecy) +S3method(make_rel,OpenSpecy) +S3method(make_rel,default) S3method(match_spec,OpenSpecy) S3method(match_spec,default) S3method(plot,OpenSpecy) diff --git a/NEWS.md b/NEWS.md index 12f37b71..7b89b11e 100644 --- a/NEWS.md +++ b/NEWS.md @@ -6,7 +6,8 @@ - The Shiny app has been outsourced to an own GitHub repository: https://github.com/wincowgerDEV/OpenSpecy-shiny - Spectra are now stored in dedicated `OpenSpecy` objects, which can be managed - with a set of new functions including `c_spec()` for concatenating spectra + with a set of new functions including `c_spec()` for concatenating spectra or + converting them back to tables - Various functions have been renamed and improved, for instance, to facilitate reading (and writing) spectral files - New functions include `def_features()` to identify microplastics in spectral diff --git a/R/adj_intens.R b/R/adj_intens.R index a1ded912..1a887602 100644 --- a/R/adj_intens.R +++ b/R/adj_intens.R @@ -9,7 +9,7 @@ #' absorbance units. For example, see \code{\link{subtr_baseline}()}. #' To run those functions properly, you will need to first convert any spectra #' from transmittance or reflectance to absorbance using this function. -#' The transmittance adjustment uses the \eqn{log10(1 / T)} +#' The transmittance adjustment uses the \eqn{log(1 / T)} #' calculation which does not correct for system and particle characteristics. #' The reflectance adjustment uses the Kubelka-Munk equation #' \eqn{(1 - R)^2 / 2R}. We assume that the reflectance intensity @@ -62,7 +62,7 @@ adj_intens.OpenSpecy <- function(x, type = "none", make_rel = TRUE, ...) { adj <- switch(type, "reflectance" = (1 - spec/100)^2 / (2 * spec/100), - "transmittance" = log10(1/adj_neg(spec, ...)), + "transmittance" = log(1/adj_neg(spec, ...)), "none" = adj_neg(spec, ...) ) diff --git a/R/as_OpenSpecy.R b/R/as_OpenSpecy.R index 542648c3..c337699f 100644 --- a/R/as_OpenSpecy.R +++ b/R/as_OpenSpecy.R @@ -318,7 +318,6 @@ check_OpenSpecy <- function(x) { stop("object 'x' is not of class 'OpenSpecy'", call. = F) if(!identical(names(x), c("wavenumber", "spectra", "metadata"))) stop("names of the object components are incorrect", call. = F) - if(!(cw <- is.vector(x$wavenumber))) warning("Wavenumber is not a vector", call. = F) if(!(cs <- is.data.table(x$spectra))) diff --git a/R/conform_spec.R b/R/conform_spec.R index 04762c48..f1e3ea7c 100644 --- a/R/conform_spec.R +++ b/R/conform_spec.R @@ -60,7 +60,7 @@ conform_spec.OpenSpecy <- function(x, range = NULL, res = 5, type = "interp", wn <- conform_res(range, res = res) } else { - wn <- range + wn <- range[range >= min(x$wavenumber) & range <= max(x$wavenumber)] } if(type == "interp") @@ -70,7 +70,8 @@ conform_spec.OpenSpecy <- function(x, range = NULL, res = 5, type = "interp", if(type == "roll") { join <- data.table("wavenumber" = wn) # Rolling join option - spec <- x$spectra[,"wavenumber" := x$wavenumber] + spec <- x$spectra + spec$wavenumber <- x$wavenumber spec <- spec[join, roll = "nearest", on = "wavenumber"] spec <- spec[,-"wavenumber"] } diff --git a/R/data_norm.R b/R/data_norm.R index a6ea0cf1..eacfc482 100644 --- a/R/data_norm.R +++ b/R/data_norm.R @@ -75,15 +75,6 @@ adj_neg <- function(y, na.rm = FALSE) { } } -#' @rdname data_norm -#' -#' @export -make_rel <- function(y, na.rm = FALSE) { - r <- range(y, na.rm = na.rm) - - (y - r[1]) / (r[2] - r[1]) -} - #' @rdname data_norm #' @importFrom data.table fifelse #' diff --git a/R/def_features.R b/R/def_features.R index a587f5e3..68f63568 100644 --- a/R/def_features.R +++ b/R/def_features.R @@ -154,6 +154,13 @@ def_features.OpenSpecy <- function(x, features, ...) { features_dt <- rbindlist(lapply(seq_along(convex_hulls), function(i) { hull <- convex_hulls[[i]] id <- names(convex_hulls)[i] + if(nrow(hull) == 1) + return(data.table(feature_id = id, + area = 1, + perimeter = 4, + feret_min = 1, + feret_max = 1) + ) # Calculate Feret dimensions dist_matrix <- as.matrix(dist(hull)) diff --git a/R/gen_OpenSpecy.R b/R/gen_OpenSpecy.R index 25b9f267..8369ec07 100644 --- a/R/gen_OpenSpecy.R +++ b/R/gen_OpenSpecy.R @@ -2,7 +2,7 @@ #' @rdname gen_OpenSpecy #' #' @description -#' Methods to visualize \code{OpenSpecy} objects. +#' Methods to visualize and convert \code{OpenSpecy} objects. #' #' @details #' \code{head()} shows the first few lines of an \code{OpenSpecy} object. @@ -19,6 +19,8 @@ #' \code{head()}, \code{print()}, and \code{summary()} return a textual #' representation of an \code{OpenSpecy} object. #' \code{plot()} and \code{lines()} return a plot. +#' \code{as.data.frame()} and \code{as.data.table()} convert \code{OpenSpecy} +#' objects into tabular data. #' #' @examples #' data("raman_hdpe") @@ -38,7 +40,9 @@ #' @seealso #' \code{\link[utils]{head}()}, \code{\link[base]{print}()}, #' \code{\link[base]{summary}()}, \code{\link[graphics]{matplot}()}, and -#' \code{\link[graphics]{matlines}()} +#' \code{\link[graphics]{matlines}()}, +#' \code{\link[base]{as.data.frame}()}, +#' \code{\link[data.table]{as.data.table}()} #' #' @importFrom utils head #' @importFrom graphics matplot matlines @@ -93,7 +97,7 @@ summary.OpenSpecy <- function(object, ...) { cat("\n$spectra\n") sl <- length(object$spectra) - sr <- range(object$spectra) + sr <- range(object$spectra, na.rm = T) array(c(sl, sr), c(1,3), list("", c("Number", "Min. Intensity", "Max. Intensity"))) |> print() @@ -104,3 +108,21 @@ summary.OpenSpecy <- function(object, ...) { t(array(c(xr, yr), c(2,2), list(c("Min.", "Max."), c("x", "y")))) |> print() names(object$metadata) |> print() } + +#' @rdname gen_OpenSpecy +#' +#' @method as.data.frame OpenSpecy +#' @export +as.data.frame.OpenSpecy <- function(x, ...) { + data.frame(wavenumber = x$wavenumber, x$spectra) +} + +#' @rdname gen_OpenSpecy +#' +#' @method as.data.table OpenSpecy +#' +#' @importFrom data.table data.table +#' @export +as.data.table.OpenSpecy <- function(x, ...) { + data.table(wavenumber = x$wavenumber, x$spectra) +} diff --git a/R/interactive_plots.R b/R/interactive_plots.R index 5f437bf6..5d1519de 100644 --- a/R/interactive_plots.R +++ b/R/interactive_plots.R @@ -89,7 +89,7 @@ plotly_spec.OpenSpecy <- function(x, x2 = NULL, p <- plot_ly(dt, type = "scatter", mode = "lines", ...) |> add_trace(x = ~wavenumber, y = ~make_rel(intensity, na.rm = T), color = ~id, line = line, - name = "Your Spectra", showlegend = F) |> + name = "x1", showlegend = F) |> layout(xaxis = list(title = "wavenumber [cm-1]", autorange = "reversed"), yaxis = list(title = "intensity [-]"), @@ -100,7 +100,7 @@ plotly_spec.OpenSpecy <- function(x, x2 = NULL, # Add dummy trace for Your Spectra p <- p |> add_trace(x = NULL, y = NULL, - line = line, name = "Your Spectra", showlegend = T) + line = line, name = "x1", showlegend = T) if (!is.null(x2)) { dt2 <- cbind(wavenumber = x2$wavenumber, x2$spectra) |> @@ -109,14 +109,14 @@ plotly_spec.OpenSpecy <- function(x, x2 = NULL, p <- p |> add_trace(data = dt2, x = ~wavenumber, y = ~make_rel(intensity, na.rm = T), color = ~id, type = "scatter", mode = "lines", - name = "Library Spectra", + name = "x2", line = line2, showlegend = F) # Add dummy trace for Library Spectra p <- p |> add_trace(x = NULL, y = NULL, line = line2, - name = "Library Spectra", showlegend = T) + name = "x2", showlegend = T) } return(p) diff --git a/R/make_rel.R b/R/make_rel.R new file mode 100644 index 00000000..74cd4b40 --- /dev/null +++ b/R/make_rel.R @@ -0,0 +1,57 @@ +#' @rdname make_rel +#' @title Make spectral intensities relative +#' +#' @description +#' \code{make_rel()} converts intensities \code{x} into relative values between +#' 0 and 1 using the standard normalization equation. +#' If \code{na.rm} is \code{TRUE}, missing values are removed before the +#' computation proceeds. +#' +#' @details +#' \code{make_rel()} is used to retain the relative height proportions between +#' spectra while avoiding the large numbers that can result from some spectral +#' instruments. +#' +#' @param x a numeric vector or an \R OpenSpecy object +#' @param na.rm logical. Should missing values be removed? +#' @param \ldots further arguments passed to \code{make_rel()}. +#' +#' @return +#' \code{make_rel()} return numeric vectors (if vector provided) or an +#' \code{OpenSpecy} object with the normalized intensity data. +#' +#' @examples +#' make_rel(c(-1000, -1, 0, 1, 10)) +#' +#' @author +#' Win Cowger, Zacharias Steinmetz +#' +#' @seealso +#' \code{\link[base]{min}()} and \code{\link[base]{round}()}; +#' \code{\link{adj_intens}()} for log transformation functions; +#' \code{\link{conform_spec}()} for conforming wavenumbers of an +#' \code{OpenSpecy} object to be matched with a reference library +#' +#' +#' @export +make_rel <- function(x, ...) { + UseMethod("make_rel") +} + +#' @rdname make_rel +#' +#' @export +make_rel.default <- function(x, na.rm = FALSE, ...) { + r <- range(x, na.rm = na.rm) + + return((x - r[1]) / (r[2] - r[1])) +} + +#' @rdname make_rel +#' +#' @export +make_rel.OpenSpecy <- function(x, na.rm = FALSE, ...) { + x$spectra <- x$spectra[, lapply(.SD, make_rel, na.rm = na.rm, ...)] + + return(x) +} diff --git a/R/manage_spec.R b/R/manage_spec.R index b5b645e1..355c776c 100644 --- a/R/manage_spec.R +++ b/R/manage_spec.R @@ -12,6 +12,8 @@ #' spectra having all the same wavenumber range. #' @param res defaults to \code{NULL}, the resolution you want the output #' wavenumbers to be. +#' @param size the number of spectra to sample. +#' @param prob probabilities to use for the sampling. #' @param \ldots further arguments passed to submethods. #' #' @return @@ -115,14 +117,10 @@ sample_spec.default <- function(x, ...) { #' @rdname manage_spec #' #' @export -sample_spec.OpenSpecy <- function(x, ...) { +sample_spec.OpenSpecy <- function(x, size = 1, prob = NULL, ...) { # replace = false is mandatory currently because we don't have a way to # rename and recoordinate duplicates. - cols <- sample(1:ncol(x$spectra), ...) + cols <- sample(1:ncol(x$spectra), size = size, replace = FALSE, prob = prob, ...) - as_OpenSpecy( - x = x$wavenumber, - spectra = x$spectra[, cols, with = F], - metadata = x$metadata[cols, ] - ) + filter_spec(x, cols) } diff --git a/R/match_spec.R b/R/match_spec.R index 86eefcd0..aacdf1ed 100644 --- a/R/match_spec.R +++ b/R/match_spec.R @@ -22,6 +22,8 @@ #' If \code{NULL} (default), all matches will be returned. #' @param cor_matrix a correlation matrix for object and library, #' can be returned by \code{cor_spec()} +#' @param order an \code{OpenSpecy} used for sorting, ideally the unprocessed +#' one; \code{NULL} skips sorting. #' @param add_library_metadata name of a column in the library metadata to be #' joined; \code{NULL} if you don't want to join. #' @param add_object_metadata name of a column in the object metadata to be @@ -53,19 +55,23 @@ #' #' @examples #' data("test_lib") +#' #' unknown <- read_extdata("ftir_ldpe_soil.asp") |> #' read_any() |> #' conform_spec(range = test_lib$wavenumber, #' res = spec_res(test_lib)) |> #' process_spec() -#' matches <- cor_spec(unknown, test_lib) +#' cor_spec(unknown, test_lib) #' #' test_lib_extract <- filter_spec(test_lib, #' logic = grepl("polycarbonate", test_lib$metadata$polymer_class, #' ignore.case = TRUE) #' ) #' -#' matches2 <- cor_spec(unknown, library = test_lib_extract) +#' cor_spec(unknown, library = test_lib_extract) +#' +#' match_spec(unknown, test_lib, add_library_metadata = "sample_name", +#' top_n = 1) #' #' @author #' Win Cowger, Zacharias Steinmetz @@ -136,16 +142,24 @@ match_spec.default <- function(x, ...) { #' #' @export match_spec.OpenSpecy <- function(x, library, na.rm = T, top_n = NULL, - add_library_metadata = NULL, + order = NULL, add_library_metadata = NULL, add_object_metadata = NULL, fill = NULL, ...) { if(is_OpenSpecy(library)) { - cor_spec(x, library = library) |> + res <- cor_spec(x, library = library) |> ident_spec(x, library = library, top_n = top_n, add_library_metadata = add_library_metadata, add_object_metadata = add_object_metadata) } else { - ai_classify(x, library, fill) + res <- ai_classify(x, library, fill) + } + + if(!is.null(order)) { + .r <- NULL + match <- match(colnames(order$spectra), res$object_id) + setorder(res[, .r := order(match)], .r)[, .r := NULL] } + + return(res) } #' @rdname match_spec diff --git a/R/process_spec.R b/R/process_spec.R index 2d31a6ff..0c51bc9e 100644 --- a/R/process_spec.R +++ b/R/process_spec.R @@ -34,6 +34,9 @@ #' \code{\link{subtr_baseline}()}. #' @param make_rel logical; if \code{TRUE} spectra are automatically normalized #' with \code{\link{make_rel}()}. +#' @param make_rel_args named list of arguments passed to +#' \code{\link{make_rel}()}. +#' @param na.rm Whether to allow NA or set all NA values to #' @param \ldots further arguments passed to subfunctions. #' #' @return @@ -105,6 +108,8 @@ process_spec.OpenSpecy <- function(x, active = TRUE, polynomial = 3, window = 11, derivative = 1, abs = TRUE), make_rel = TRUE, + make_rel_args = list( + na.rm = TRUE), ...) { if(active) { if(adj_intens) @@ -124,7 +129,8 @@ process_spec.OpenSpecy <- function(x, active = TRUE, x <- do.call("smooth_intens", c(list(x, make_rel = F), smooth_intens_args)) if(make_rel) - x$spectra <- x$spectra[, lapply(.SD, make_rel)] + x <- do.call("make_rel", c(list(x), + make_rel_args)) } return(x) diff --git a/R/read_ext.R b/R/read_ext.R index b4e41353..4602cf39 100644 --- a/R/read_ext.R +++ b/R/read_ext.R @@ -303,7 +303,7 @@ read_jdx <- function(file, share = NULL, other_info = NULL, license = "CC BY-NC"), ...) { - jdx <- read.jdx(file, ...) + jdx <- read.jdx(file, encoding = 'latin1') x <- jdx@wavelength y <- as.numeric(unname(jdx@data$spc[1,])) diff --git a/R/sig_noise.R b/R/sig_noise.R index aae90245..ecef02c8 100644 --- a/R/sig_noise.R +++ b/R/sig_noise.R @@ -10,8 +10,11 @@ #' Options include \code{"sig"} (mean intensity), \code{"noise"} (standard #' deviation of intensity), \code{"sig_times_noise"} (absolute value of #' signal times noise), \code{"sig_over_noise"} (absolute value of signal / -#' noise), or \code{"tot_sig"} (total signal = signal * number of data -#' points). +#' noise), \code{"run_sig_over_noise"} (absolute value of signal / +#' noise where signal is estimated as the max intensity and noise is +#' estimated as the height of a low intensity region.), +#' \code{"log_tot_sig"} (sum of the inverse log intensities, useful for spectra in log units), +#' or \code{"tot_sig"} (sum of intensities). #' @param na.rm logical; indicating whether missing values should be removed #' when calculating signal and noise. Default is \code{TRUE}. #' @param \ldots further arguments passed to subfunctions; currently not used. @@ -45,7 +48,7 @@ sig_noise.default <- function(x, ...) { #' @rdname sig_noise #' #' @export -sig_noise.OpenSpecy <- function(x, metric = "sig_over_noise", +sig_noise.OpenSpecy <- function(x, metric = "run_sig_over_noise", na.rm = TRUE, ...) { vapply(x$spectra, function(y) { if(length(y[!is.na(y)]) < 20) { diff --git a/R/smooth_intens.R b/R/smooth_intens.R index c52326c3..a2c4a58c 100644 --- a/R/smooth_intens.R +++ b/R/smooth_intens.R @@ -60,7 +60,7 @@ smooth_intens.default <- function(x, ...) { #' #' @export smooth_intens.OpenSpecy <- function(x, polynomial = 3, window = 11, - derivative = 0, abs = FALSE, + derivative = 1, abs = TRUE, make_rel = TRUE, ...) { filt <- x$spectra[, lapply(.SD, .sgfilt, p = polynomial, n = window, m = derivative, abs = abs, ...)] diff --git a/_pkgdown.yml b/_pkgdown.yml index c27a0f6c..b92617f1 100644 --- a/_pkgdown.yml +++ b/_pkgdown.yml @@ -3,3 +3,11 @@ template: bootstrap: 5 development: mode: auto +articles: +- title: Articles + navbar: ~ + contents: + - app + - sop + - advanced + - spectragryph diff --git a/cran-comments.md b/cran-comments.md index 72c6bdb7..299d0275 100644 --- a/cran-comments.md +++ b/cran-comments.md @@ -16,16 +16,15 @@ ## Comments Possibly misspelled words in DESCRIPTION: - preprocessing (46:5) - wavenumber (47:62) + wavenumber (50:62) -Both words are spelled correctly. +The word is spelled correctly. - installed size is 6.5Mb - sub-directories of 1Mb or more: - doc 4.9Mb + installed size is 6.1Mb + sub-directories of 1Mb or more: + doc 4.5Mb We tried to reduce the vignette size as much as possible but would like to keep -the current number of screenshots for better comprehensibility. +the current number of figures for better comprehensibility. For further details, see NEWS.md diff --git a/inst/extdata/fitr_nitrocellulose.jdx b/inst/extdata/fitr_nitrocellulose.jdx index 3a09a662..357d5c98 100644 --- a/inst/extdata/fitr_nitrocellulose.jdx +++ b/inst/extdata/fitr_nitrocellulose.jdx @@ -1023,1027 +1023,3 @@ 607.64+66946033+67829741+67681322+66557188+64851073+63203464+62214539+62111525 599.92+62840942 ##END= -##TITLE=Nitrocellulose -##JCAMP-DX=4.24 -##DATA TYPE=INFRARED SPECTRUM -##ORIGIN=Pacific Northwest National Laboratory Under IARPA Contract -##OWNER=Public domain -##DATE=March 2017 -##CAS REGISTRY NO=9004-70-0 -##MOLFORM=C 6 H 9 ( N O 2 ) O 5 ; C 6 H 8 ( N O 2 ) 2 O 5 ; C 6 H 7 ( N O 2 ) 3 O 5 -##$NIST SOURCE=IARPA-IR-S -##$NIST DOC FILE=META.MIR.PUR.Nitrocellulose.SLD.UNKNOWN.PNNL161235.pdf -##SPECTROMETER/DATA SYSTEM=Bruker IR Cube FTIR -##APERTURE=6 mm -##EXTERNAL DIFFUSE REFLECTANCE ACCESSORY=A 562-G integrating sphere -##BEAMSPLITTER=Ge on KBr -##DETECTOR (DIA. DET. PORT IN SPHERE)=2×2 mm, 60° field of view MCT -##SPHERE DIAMETER=75 mm -##ACQUISITION MODE=double-sided, forward-backward -##SCANNER SPEED=40 kHz -##COADDED SCANS=2048 -##PHASE RESOLUTION=32.00 -##PHASE CORRECTION=Mertz -##ZEROFILLING=4× -##SPECTRAL RANGE=7,500 to 600 cm-1 saved; 7500 to 600 cm-1 reported -##STATE=solid -##RESOLUTION=0.96450084 -##SPECTRAL RESOLUTION=4 cm-1 -##WAVENUMBER ACCURACY=± 0.4 cm-1 -##APODIZATION FUNCTION=Blackman-Harris 3-term -##LOW PASS FILTER=Open -##SWITCH GAIN ON=512 points -##XUNITS=1/cm -##YUNITS=Reflectance -##XLABEL=Wavenumber -##YLABEL=Diffuse-only reflectance -##XFACTOR=1 -##YFACTOR=5.6371038e-010 -##DELTAX=-0.96450084 -##FIRSTX=7498.994 -##LASTX=599.91952 -##FIRSTY=0.60527939 -##MAXX=7498.99 -##MINX=599.916 -##MAXY=0.60527941 -##MINY=0.009793684 -##NPOINTS=7154 -##XYDATA=(X++(Y..Y)) -7498.99+1073741786+1073642077+1073318312+1072754315+1071941097+1070878342 -7493.21+1069574930+1068046089+1066309582+1064383067+1062284730+1060036354 -7487.42+1057662257+1055283508+1052931086+1050637878+1048434968+1046350588 -7481.63+1044407895+1042622221+1041000015+1039540220+1038236174+1037075930 -7475.85+1036040666+1035133765+1034357027+1033714996+1033214017+1032857369 -7470.06+1032639340+1032543755+1032549147+1032669370+1032880842+1033160303 -7464.27+1033477089+1033789857+1034048805+1034202123+1034203497+1034019199 -7458.48+1033632310+1033042090+1032261016+1031313196+1030232254+1029056360 -7452.70+1027819352+1026537828+1025204916+1023791962+1022261852+1020584134 -7446.91+1018624629+1016371283+1013806332+1010982856+1007986924+1004931780 -7441.12+1001939125+999113957+996527753+994213292+992170044+990372527+988778235 -7434.37+987331657+985968927+984627240+983260069+981849018+980408996+978982402 -7427.62+977626228+976399476+975355541+974656942+974346394+974429503+974863762 -7420.87+975555806+976377060+977195459+977903257+978428026+978724723+978763951 -7414.12+978529110+978030352+977320333+976494427+975664714+974917899+974279780 -7407.37+973707958+973116258+972413746+971534760+970444196+969124607+967559922 -7400.61+965733962+963643555+961315454+958809292+956204373+953584545+951034925 -7393.86+948650149+946541027+944817103+943546893+942718766+942230264+941921726 -7387.11+941640044+941299785+940902745+940509829+940186911+939963490+939828253 -7380.36+939757304+939748739+939835972+940067111+940454952+940973271+941560002 -7373.61+942134044+942610386+942904122+942970630+942734944+942143984+941189396 -7366.86+939916543+938411387+936768350+935057749+933315426+931556820+929799800 -7360.11+928076087+926419199+924841614+923318694+921793871+920205078+918530743 -7353.35+916786200+915015117+913269200+911588839+909996767+908464014+906999143 -7346.60+905622245+904364829+903263903+902356157+901670880+901221501+900997446 -7339.85+900960755+901054120+901220021+901423034+901662739+901970008+902389676 -7333.10+902958537+903689915+904569641+905616747+906779846+908016644+909302080 -7326.35+910625370+911978477+913343109+914681942+915936186+917029288+917879725 -7319.60+918420249+918616178+918474597+918038118+917369019+916529684+915568224 -7312.85+914502085+913356221+912132112+910817387+909396502+907809929+906079872 -7306.09+904243761+902342623+900407649+898453219+896478382+894475103+892438199 -7299.34+890373698+888303064+886261772+884290953+882425712+880682914+879052196 -7292.59+877495017+875952853+874373999+872692951+870876666+868923769+866855091 -7285.84+864699023+862478191+860203010+857875067+855498062+853088596+850680241 -7279.09+848321317+846067759+843974710+842086207+840421917+838964924+837657495 -7272.34+836410282+835126327+833729973+832186488+830503800+828718549+826876940 -7265.58+825019682+823175641+821362898+819595780+817891734+816274497+814769870 -7258.83+813394505+812144068+810986520+809863917+808700818+807416492+805940149 -7252.08+804229336+802287331+800170807+797979158+795828218+793819281+792017429 -7245.33+790444972+789088692+787913856+786876478+785928816+785018585+784086255 -7238.58+783066217+781894924+780524317+778935365+777139011+775207420+773179927 -7231.83+771108870+769042625+767019413+765059221+763158241+761288031+759401801 -7225.08+757445045+755367697+753134018+750726720+748145380+745400994+742511539 -7218.32+739499376+736389460+733207115+729974175+726705179+723405995+720080324 -7211.57+716697503+713277569+709830301+706368971+702906266+699453870+696022146 -7204.82+692620187+689255922+685936490+682669027+679461464+676322576+673262039 -7198.07+670289633+667414611+664644637+661985528+659440455+657010001+654691786 -7191.32+652480628+650369339+648325721+646355959+644456037+642619345+640838587 -7184.57+639106150+637414105+635754574+634119574+632501227+630891498+629282192 -7177.82+627664426+626029373+624368256+622673092+620937325+619156250+617326853 -7171.06+615447708+613518338+611538902+609509981+607455895+605363004+603233316 -7164.31+601069687+598875871+596656201+594415385+592158550+589890771+587617123 -7157.56+585342630+583072260+580810772+578562766+576332682+574124696+571942722 -7150.81+569790249+567670500+565586121+563539489+561532456+559566342+557642100 -7144.06+555759994+553919919+552121133+550362474+548642514+546959404+545311398 -7137.31+543696435+542112664+540558235+539031456+537530740+536054715+534602057 -7130.55+533171709+531762720+530374138+529005223+527655288+526323698+525009713 -7123.80+523712646+522431704+521166040+519914863+518677325+517452635+516239892 -7117.05+515038464+513847556+512666481+511494553+510331031+509175069+508025610 -7110.30+506881702+505741970+504605199+503470120+502335411+501199909+500062663 -7103.55+498922772+497779816+496633582+495484281+494332496+493179018+492024960 -7096.80+490871535+489720067+488571877+487428127+486289876+485158128+484033834 -7090.05+482917893+481864865+480774988+479697218+478632665+477582492+476547810 -7083.29+475529570+474528775+473546062+472582065+471636994+470711114+469804214 -7076.54+468916029+468046136+467193849+466358479+465539182+464734951+463944837 -7069.79+463167781+462402620+461648245+460903491+460167144+459437986+458714750 -7063.04+457996325+457281495+456569202+455858654+455149005+454439832+453730817 -7056.29+453021908+452313263+451605306+450898670+450194202+449492854+448795840 -7049.54+448104589+447420581+446745455+446080902+445428826+444791025+444169348 -7042.79+443565700+442981903+442419598+441880316+441365433+440876007+440412908 -7036.03+439976614+439567441+439185336+438830036+438501064+438197653+437918880 -7029.28+437663632+437430642+437218483+437025646+436850520+436691334+436546370 -7022.53+436413777+436291704+436178328+436071799+435970371+435872380+435776266 -7015.78+435680548+435583958+435485333+435383641+435278037+435167727+435052131 -7009.03+434930666+434802963+434668652+434527494+434379331+434224084+434061673 -7002.28+433892177+433715677+433532383+433342507+433146314+432944225+432736638 -6995.52+432524056+432306979+432086043+431861908+431635263+431406899+431177557 -6988.77+430948083+430719322+430492095+430267273+430045703+429828177+429615488 -6982.02+429408377+429207557+429013690+428827356+428649111+428479404+428318685 -6975.27+428167218+428025214+427892780+427769888+427656406+427552018+427456300 -6968.52+427368777+427288708+427215354+427147841+427085219+427026456+426970495 -6961.77+426916252+426862591+426808427+426752652+426694126+426631795+426564626 -6955.02+426491615+426411811+426324340+426228437+426123362+426008400+425883023 -6948.26+425746835+425599544+425441125+425271656+425091402+424900786+424700389 -6941.51+424490926+424273241+424048419+423817571+423581832+423342498+423100784 -6934.76+422857908+422615164+422373742+422134830+421899567+421669036+421444399 -6928.01+421226794+421017330+420817039+420626714+420447094+420278656+420121823 -6921.26+419976700+419843340+419721558+419610932+419510773+419420263+419338449 -6914.51+419264328+419196895+419135065+419077809+419024122+418973025+418923593 -6907.75+418874954+418826263+418776672+418725364+418671676+418615001+418554864 -6901.00+418490867+418422746+418350449+418273976+418193537+418109529+418022508 -6894.25+417933293+417842730+417751823+417661683+417573446+417488276+417407387 -6887.50+417331944+417263057+417201810+417149259+417106277+417073789+417052589 -6880.75+417043337+417046721+417063295+417093509+417137760+417196364+417269507 -6874.00+417357268+417459568+417576275+417707070+417851506+418009000+418178866 -6867.25+418360362+418552617+418754785+418965967+419185264+419411777+419644661 -6860.49+419883070+420126237+420373396+420623885+420877018+421132186+421388808 -6853.74+421646355+421904325+422162243+422419737+422676491+422932241+423186669 -6846.99+423439669+423691269+423941546+424190555+424438481+424685507+424931846 -6840.24+425177736+425423414+425669119+425915220+426161929+426409326+426657542 -6833.49+426906683+427156802+427407926+427660054+427913134+428167059+428421646 -6826.74+428676708+428931982+429187203+429441975+429695901+429948584+430199522 -6819.99+430448267+430694263+430937007+431175971+431410652+431640602+431865318 -6813.23+432084430+432297568+432504415+432704706+432898203+433084775+433264289 -6806.48+433436692+433602063+433760430+433911977+434056915+434195482+434328022 -6799.73+434454853+434576423+434693103+434805395+434913775+435018639+435120489 -6792.98+435219776+435316948+435412454+435506718+435600241+435693527+435787104 -6786.23+435881659+435977826+436076319+436177879+436283271+436393290+436508728 -6779.48+436630404+436759138+436895696+437040819+437195168+437359455+437534290 -6772.72+437720201+437917796+438127497+438349702+438584754+438832838+439094086 -6765.97+439368445+439655757+439955652+440267706+440591312+440925835+441270429 -6759.22+441624329+441986634+442356447+442732815+443114867+443501782+443892795 -6752.47+444287192+444684284+445083544+445484496+445886770+446290048+446694067 -6745.72+447098878+447504324+447910616+448317912+448726371+449136417+449548312 -6738.97+449973848+450413605+450868218+451338375+451824603+452327485+452847443 -6732.22+453384848+453939963+454512948+455103856+455701107+456304227+456912581 -6725.46+457525322+458141448+458759952+459379620+459999287+460617739+461233653 -6718.71+461845813+462453057+463054326+463648617+464232704+464806165+465368841 -6711.96+465920731+466461837+466992474+467513014+468023879+468525809+469019545 -6705.21+469505932+469988406+470467603+470944368+471419441+471893721+472368160 -6698.46+472843709+473321214+473801574+474285688+474774242+475267978+475767318 -6691.71+476272737+476784501+477302661+477827113+478357698+478893940+479435257 -6684.96+479980856+480529839+481081148+481633567+482185881+482736661+483284428 -6678.20+483827807+484365211+484895267+485416547+485927624+486427333+486914513 -6671.45+487388159+487847530+488291940+488721018+489134500+489532385+489914833 -6664.70+490282320+490635320+490974628+491301089+491615707+491919646+492214016 -6657.95+492500033+492779071+493052241+493320917+493586262+493849439+494111454 -6651.20+494373363+494636065+494900300+495166755+495435854+495708125+495983832 -6644.45+496263187+496546296+496833370+497124304+497419308+497718277+498021476 -6637.69+498328905+498640827+498957454+499279209+499606358+499939480+500279052 -6630.94+500625708+500980136+501343023+501715162+502097399+502490526+502895391 -6624.19+503312890+503743819+504188916+504648921+505124312+505615562+506122991 -6617.44+506646650+507186434+507742025+508313001+508898569+509497988+510110095 -6610.69+510733728+511367564+512010176+512659979+513315385+513974809+514636560 -6603.94+515298998+515960484+516619538+517274575+517924377+518567624+519203258 -6597.19+519830379+520448355+521056656+521655071+522243440+522822082+523391261 -6590.43+523951610+524503818+525048730+525587350+526120790+526650053+527176250 -6583.68+527700596+528224044+528747544+529272049+529798298+530326874+530858252 -6576.93+531392749+531930523+532471470+533015431+533561982+534110541+534660370 -6570.18+535210675+535760398+536308588+536854028+537395557+537932115+538462436 -6563.43+538985566+539500449+540006239+540502247+540987841+541462544+541926198 -6556.68+542378591+542819828+543250122+543669895+544079571+544479888+544871482 -6549.93+545255305+545632149+546003072+546369078+546731067+547090041+547447060 -6543.17+547802863+548158295+548514151+548871011+549229404+549589806+549952481 -6536.42+550317800+550685868+551056897+551430939+551808100+552188434+552571992 -6529.67+552958829+553349048+553742757+554140061+554541172+554946248+555355606 -6522.92+555769458+556188121+556612017+557041412+557476835+557918654+558367240 -6516.17+558823069+559286511+559757884+560237398+560725318+561221750+561726693 -6509.42+562240043+562761534+563290850+563827568+564371105+564920828+565476102 -6502.66+566036029+566599973+567167037+567736427+568307403+568879119+569450835 -6495.91+570022023+570591994+571160168+571726227+572289749+572850468+573408227 -6489.16+573962867+574514493+575063159+575609075+576152454+576693665+577240586 -6482.41+577776035+578310479+578844553+579378574+579913071+580448467+580984867 -6475.66+581522536+582061580+582601998+583143473+583685795+584228698+584771548 -6468.91+585313869+585855028+586394177+586930683+587463753+587992698+588516622 -6462.16+589034888+589546863+590051860+590549402+591039120+591520431+591993178 -6455.40+592457308+592912767+593359555+593797779+594228020+594650437+595065610 -6448.65+595474228+595877031+596274758+596668255+597058687+597446739+597776213 -6441.90+598105106+598434844+598766804+599102146+599442036+599787317+600138679 -6435.15+600496755+600861915+601234371+601614176+602001277+602395514+602796572 -6428.40+603203974+603617244+604035749+604458694+604885340+605314841+605746245 -6421.65+606178653+606611220+607043153+607473659+607902208+608328325+608751640 -6414.90+609172153+609589706+610004562+610416934+610827138+611235650+611642999 -6408.14+612049714+612456376+612863619+613272184+613682758+614096028+614512629 -6401.39+614933354+615358678+615789131+616225134+616666795+617114165+617567033 -6394.64+618024977+618487521+618953765+619422864+619893813+620365397+620836558 -6387.89+621306133+621772853+622235661+622693393+623144993+623589297+624025406 -6381.14+624452211+624868865+625274522+625668548+626050468+626419910+626776876 -6374.39+627121523+627454275+627775767+628086737+628388191+628681134+628966780 -6367.63+629246241+629520680+629791312+630059195+630325333+630590678+630856129 -6360.88+631122426+631390362+631660624+631933688+632210135+632490284+632774397 -6354.13+633062581+633354836+633650951+633950819+634253965+634560072+634868716 -6347.38+635179475+635491873+635805592+636120210+636435463+636751086+637066973 -6340.63+637383072+637699488+638016432+638334170+638653176+638973874+639297057 -6333.88+639623254+639953310+640287965+640628119+640974511+641328040+641689341 -6327.13+642059154+642438165+642826958+643226059+643635893+644056829+644489185 -6320.37+644933172+645388737+645855827+646334125+646823049+647321913+647829818 -6313.62+648345546+648867884+649395402+649926462+650459532+650992972+651525089 -6306.87+652054300+652578963+653097600+653608729+654110976+654603020+655083644 -6300.12+655551739+656006405+656446902+656872491+657282906+657677990+658057795 -6293.37+658422638+658773259+659110505+659435486+659749523+660054202+660351003 -6286.62+660641725+660928059+661211855+661494806+661778602+662064883+662355129 -6279.86+662650662+662952751+663262452+663580613+663907919+664244795+664591505 -6273.11+664948153+665314476+665690105+666074403+666466368+666864835+667268483 -6266.36+667675726+668084820+668493967+668901157+669304541+669702110+670092118 -6259.61+670472663+670842211+671199336+671542661+671871237+672184110+672480701 -6252.86+672760532+673023287+673269018+673497778+673709938+673906132+674087152 -6246.11+674254004+674408009+674550700+674683716+674809067+674928602+675044541 -6239.36+675158948+675274042+675391991+675514804+675644489+675782898+675931722 -6232.60+676092441+676266483+676454905+676658765+676878696+677115334+677368995 -6225.85+677639786+677927600+678232121+678552819+678888849+679239259+679602674 -6219.10+679977774+680362971+680756469+681156416+681560858+681967837+682375186 -6212.35+682781002+683183170+683579892+683969319+684349811+684719994+685078493 -6205.60+685424197+685756262+686073894+686376511+686663797+686935486+687191580 -6198.85+687432288+687657929+687868873+688065754+688249313+688420341+688579791 -6192.10+688728668+688868187+688999512+689123699+689242071+689355684+689465544 -6185.34+689572708+689677968+689782119+689885687+689989309+690093301+690197927 -6178.59+690303451+690409981+690517567+690626159+690735860+690846566+690958276 -6171.84+691070938+691184393+691298589+691413365+691528459+691643712+691758912 -6165.09+691873688+691987831+692101127+692213208+692323966+692433139+692540673 -6158.34+692646409+692750295+692852331+692952516+693050850+693147441+693242233 -6151.59+693335281+693426637+693516249+693604221+693690608+693775567+693859416 -6144.83+693942419+694025052+694107896+694191692+694277286+694365681+694457936 -6138.08+694555214+694594759+694669885+694747389+694828648+694915140+695029335 -6131.33+695157858+695303192+695467559+695653179+695861797+696095104+696354317 -6124.58+696640227+696953260+697293572+697660688+698053975+698472479+698914774 -6117.83+699379167+699863704+700366004+700883160+701412053+701949088+702490299 -6111.08+703031563+703568545+704096804+704612004+705110074+705587157+706039919 -6104.33+706465455+706861438+707226228+707558610+707858002+708124193+708357289 -6097.57+708557818+708726573+708864559+708973256+709054461+709110343+709143438 -6090.82+709156761+709153589+709137306+709111559+709079838+709045474+709011426 -6084.07+708980234+708953853+708933552+708920017+708913250+708912986+708918484 -6077.32+708928794+708942909+708959827+708978542+708998051+709017242+709034900 -6070.57+709049809+709060329+709064823+709061281+709047747+709022053+708982137 -6063.82+708926097+708852346+708759721+708647588+708516105+708365801+708198103 -6057.07+708014809+707818192+707610896+707395670+707175369+706952952+706731223 -6050.31+706512931+706300771+706097070+705904101+705723979+705558555+705409573 -6043.56+705278565+705166908+705075710+705005977+704958290+704933019+704930111 -6036.81+704948932+704988531+705047267+705123080+705213432+705315309+705425274 -6030.06+705539787+705654986+705767067+705872380+705967437+706049171+706115045 -6023.31+706162996+706191545+706199687+706187104+706153639+706099713+706025804 -6016.56+705932756+705821468+705692999+705548457+705389218+705216762+705032781 -6009.80+704839072+704637909+704431565+704222841+704014541+703809835+703611791 -6003.05+703423581+703248217+703088502+702946657+702824796+702724294+702646207 -5996.30+702591172+702559504+702551150+702565848+702603278+702662808+702743961 -5989.55+702846208+702968862+703111447+703273171+703453134+703650015+703862281 -5982.80+704087763+704323872+704567700+704815969+705065242+705312137+705553110 -5976.05+705784989+706004604+706209256+706396515+706564319+706711133+706835585 -5969.30+706936669+707013751+707066196+707093634+707095855+707073069+707025646 -5962.54+706954327+706860433+706745551+706611636+706460856+706295696+706118587 -5955.79+705931804+705737672+705538095+705334764+705128896+704921335+704712347 -5949.04+704501932+704289455+704074070+703854720+703629819+703397781+703156649 -5942.29+702904363+702638647+702357230+702057838+701738408+701396880+701031720 -5935.54+700641606+700225639+699783503+699315250+698821461+698303406+697762618 -5928.79+697201211+696621776+696027063+695420189+694804381+694183022+693559495 -5922.04+692937343+692319843+691710114+691110801+690524176+689951878+689394807 -5915.28+688853331+688327240+687815688+687317511+686831019+686354465+685885736 -5908.53+685422717+684963399+684505666+684047722+683587822+683124327+682655810 -5901.78+682252214+681830908+681388983+680924907+680438308+679930087+679401723 -5895.03+678784117+678160854+677536111+676913008+676179515+675426461+674663574 -5888.28+673901797+673151968+672424079+671726484+671178876+670694656+670265102 -5881.53+669877685+669517864+669170150+668819370+668451091+668052148+667610381 -5874.77+667114848+666555398+665922619+665207683+664403082+663503848+662508923 -5868.02+661422271+660253197+659015554+657726682+656406300+655075291+653754856 -5861.27+652465561+651226332+650053240+648958077+647947343+647020988+646173036 -5854.52+645392437+644665289+643977369+643316623+642675332+642050800+641445513 -5847.77+640866184+640322329+639824999+639384977+639011094+638709692+638483945 -5841.02+638334276+638258463+638251431+638305727+638412732+638564834+638757961 -5834.27+638993383+639278289+639624946+640048208+640561822+641175040+641889923 -5827.51+642700814+643595026+644554636+645558232+646582394+647431879+648232091 -5820.76+648970712+649636163+650218295+650709968+651110127+651595721+652034685 -5814.01+652433047+652798577+653140529+653469315+653795829+654130855+654478991 -5807.26+654847588+655242355+655667150+656123666+656611587+657128478+657670324 -5800.51+658231625+658805773+659385790+659964590+660535355+661092003+661629249 -5793.76+662143022+662630255+663089203+663518968+663919603+664291689+664636178 -5787.01+664954021+665245906+665512203+665752595+665966341+666152014+666307922 -5780.25+666432215+666523201+666579505+666600388+666585902+666537158+666456164 -5773.50+666345987+666210380+666053890+665881276+665697295+665506282+665311674 -5766.75+665115851+664919869+664723041+664523463+664317648+664101100+663868480 -5760.00+663614131+663332503+663018413+662667844+662277942+661847436+661376857 -5753.25+660868583+660326737+659756871+659165805+658561205+657950948+657342912 -5746.50+656744444+656161996+655601012+655065669+654558769+654081687+653634634 -5739.74+653216553+652825487+652458688+652113142+651785571+651472750+651171772 -5732.99+650879834+650594504+650313616+650035054+649756861+649477136+649193921 -5726.24+648905526+648609940+648305419+647990061+647661961+647319376+646960665 -5719.49+646584561+646190112+645777000+645345385+644895847+644429497+643947657 -5712.74+643451859+642943691+642424737+641896743+641361241+640819924+640274272 -5705.99+639725501+639174826+638623094+638071098+637519630+636969379+636421189 -5699.24+635875960+635334537+634797661+634265913+633739716+633219124+632703924 -5692.48+632193799+631688115+631186238+630687268+630190308+629694352+629198448 -5685.73+628701858+628204263+627705664+627206641+626708094+626211451+625718350 -5678.98+625230694+624750387+624279490+623819802+623372960+622940340+622522841 -5672.23+622120778+621733942+621361750+621002987+620656119+620319349+619990562 -5665.48+619667538+619347897+619028996+618708245+618382948+618050301+617707821 -5658.73+617353023+616983528+616597379+616192726+615768301+615323204+614857277 -5651.97+614370995+613865735+613343714+612807790+612261556+611708926+611154180 -5645.22+610601708+610055792+609520449+608999275+608495071+608010059+607545348 -5638.47+607101308+606677253+606271913+605883279+605509131+605147090+604794724 -5631.72+604449918+604110927+603776325+603445317+603117693+602793612+602473759 -5624.97+602159035+601850285+601548302+601253457+600965590+600684120+600407937 -5618.22+600135560+599865139+599594824+599322818+599047798+598769024+598486655 -5611.47+598202066+597917635+597637064+597365005+597106955+596869049+596657259 -5604.71+596477454+596334816+596233521+596176423+596164845+596198469+596275181 -5597.96+596391279+596541478+596719537+596918269+597130270+597348193+597565058 -5591.21+597774468+597970979+598149885+598307538+598441136+598548881+598629452 -5584.46+598682320+598707380+598704736+598674602+598617398+598533602+598424007 -5577.71+598289880+598133179+597956546+597763736+597559348+597348880+597138570 -5570.96+596935134+596745654+596577164+596436376+596329212+596260483+596233521 -5564.21+596249857+596309016+596408408+596543751+596708805+596896117+597097280 -5557.45+597303308+597505211+597694215+597862335+598002700+598109864+598179967 -5550.70+598210895+598202436+598155965+598074389+597961727+597823107+597664080 -5543.95+597490725+597308964+597124772+596943487+596769921+596608197+596461329 -5537.20+596331538+596219934+596126727+596051284+595992336+595947980+595916047 -5530.45+595894160+595879727+595870316+595863655+595857787+595851125+595842561 -5523.70+595831617+595818347+595803385+595787895+595773568+595762254+595756280 -5516.94+595757919+595769550+595793288+595830930+595883692+595952209+596036217 -5510.19+596134499+596244940+596364634+596489773+596615916+596738411+596852025 -5503.44+596951787+597032781+597090460+597120912+597120912+597088293+597021785 -5496.69+596921282+596787685+596623000+596430190+596212849+595975365+595722603 -5489.94+595459531+595191278+594922708+594658367+594402380+594158076+593928258 -5483.19+593714882+593519323+593342320+593184033+593044091+592921966+592816652 -5476.44+592726988+592651545+592588791+592536768+592493469+592456567+592423419 -5469.68+592391170+592356594+592316308+592266771+592204175+592124926+592025375 -5462.93+591902404+591753210+591575890+591369334+591133436+590869413+590579643 -5456.18+590267720+589938299+589597194+589250802+588905890+588569278+588247206 -5449.43+587945064+587666978+587415378+587190900+586992010+586815378+586655557 -5442.68+586505623+586357540+586202371+586031237+585835572+585607922+585342154 -5435.93+585034144+584681672+584285002+583846461+583370225+582862268+582329674 -5429.18+581780163+581221505+580660785+580104084+579555630+579017855+578490918 -5422.42+577972969+577459989+576946269+576424725+575887532+575326548+574815683 -5415.67+574196967+573530300+572811399+572038572+571212507+570336218+569414938 -5408.92+568455804+567467275+566458815+565440152+564420643+563408852+562412023 -5402.17+561435548+560483023+559555716+558652305+557769301+556900942+556039508 -5395.42+555175960+554300305+553401864+552469904+551494064+550464879+549374157 -5388.67+548215498+546984094+545677246+544294004+542835266+541303676+539703304 -5381.91+538039491+536318421+534546862+532731793+530880298+528999143+527094621 -5375.16+525172124+523236198+521290175+519336009+517374760+515406109+513428576 -5368.41+511439359+509434704+507409749+505359152+503277204+501158249+498996841 -5361.66+496787851+494526786+492209945+489834632+487399473+484904572+482351252 -5354.91+479742421+477082148+474375458+471628218+468846877+466038098+463208437 -5348.16+460364448+457512000+454656485+451802293+448952911+446110508+443276247 -5341.41+440450259+437631859+434819643+432014407+429213347+426413370+423610962 -5334.65+420802447+418021292+415222003+412399426+409549674+406669682+403757017 -5327.90+400809882+398654448+395579452+392434669+389236596+385995197+382717636 -5321.15+379412716+376093654+372776072+369472342+366186877+362914762+359645026 -5314.40+356366275+353073251+349761035+346430368+343085295+339736177+336404665 -5307.65+333124196+329930458+326841980+323842136+320878638+317884477+314818257 -5300.90+311678603+308498055+305313093+302134475+298941318+295703303+292416305 -5294.15+289127114+285927692+282919495+280166359+277663502+275339815+273091306 -5287.39+270824399+268485592+266063491+263570493+261020266+258416536+255756343 -5280.64+252963901+250379891+247843092+245376793+242999604+240725586+238564706 -5273.89+236523361+234604406+232806625+231123700+229544396+228054466+226640375 -5267.14+225293136+224009815+222792367+221644600+220569460+219568283+218641874 -5260.39+217792177+217022444+216335237+215729223+215196867+214725455+214301016 -5253.64+213912726+213554703+213224330+212918937+212633198+212359169+212088524 -5246.88+211815698+211539871+211264574+210995409+210737188+210492131+210259603 -5240.13+210036896+209820203+209605373+209388574+209166713+208937688+208700178 -5233.38+208453442+208197732+207935533+207672725+207418403+207182783+206974020 -5226.63+206795933+206647572+206524878+206422618+206335016+206254445+206168983 -5219.88+206061727+205913485+205708660+205441200+205116868+204750386+204359769 -5213.13+203961500+203569311+203196366+202858235+202573764+202362424+202239360 -5206.38+202210428+202269944+202401626+202581443+202780783+202969522+203119086 -5199.62+203205935+203214539+203138462+202979210+202743801+202442955+202091461 -5192.87+201710216+201327649+200977820+200694474+200503025+200414365+200423326 -5186.12+200512250+200657241+200833609+201018410+201190549+201330900+201423882 -5179.37+201460189+201438698+201366243+201255510+201122864+200987191+200869665 -5172.62+200793007+200778864+200843192+200991355+201216018+201499272+201818728 -5165.87+202154110+202491211+202821782+203140564+203442124+203719351+203964170 -5159.11+204169735+204332582+204453597+204537063+204588331+204611210+204606452 -5152.36+204572101+204505593+204406584+204279793+204136375+203993485+203872232 -5145.61+203794358+203779184+203841238+203988767+204223356+204540235+204929133 -5138.86+205375353+205860643+206363855+206861966+207331937+207753627+208113064 -5132.11+208404883+208632758+208807474+208943662+209056575+209160448+209268537 -5125.36+209393848+209549016+209744469+209985588+210270494+210590386+210932720 -5118.61+211285985+211643599+212005045+212373919+212753975+213145675+213545015 -5111.85+213945028+214338619+214720697+215088553+215440430+215774015+216086016 -5105.10+216373143+216633559+216868161+217080717+217277704+217467620+217660443 -5098.35+217866576+218095469+218353743+218643394+218960947+219298404+219645933 -5091.60+219995616+220344427+220694837+221052979+221425157+221814544+222219752 -5084.85+222635666+223055915+223475899+223894589+224313807+224735959+225161177 -5078.10+225586237+226006010+226416426+226817166+227212356+227609052+228014656 -5071.35+228434191+228868450+229313309+229760600+230200067+230622034+231019920 -5064.59+231391662+231739772+232069881+232388200+232699593+233007392+233315243 -5057.84+233629280+233958701+234313631+234701419+235123757+235576176+236050535 -5051.09+236538244+237033037+237531610+238032721+238535391+239037982+239538405 -5044.34+240035419+240529789+241024635+241525111+242037192+242566375+243116178 -5037.59+243686546+244273594+244870343+245468863+246061727+246642562+247205317 -5030.84+247743251+248248697+248714280+249134899+249507831+249843041+250146690 -5024.08+250426547+250690782+250947589+251205586+251473997+251762419+252079655 -5017.33+252432682+252827792+253260756+253732551+254242147+254786107+255358114 -5010.58+255949365+256549286+257146961+257732370+258297213+258835094+259341227 -5003.83+259812123+260245563+260641360+261002027+261333087+261642683+261940701 -4997.08+262237661+262543847+262868959+263222066+263611545+264044932+264527935 -4990.33+265078716+265692251+266358469+267063756+267791539+268523366+269240523 -4983.58+269926460+270569205+271162200+271703808+272164713+272556730+272895615 -4976.82+273199739+273489932+273788399+274115098+274483034+274894903+275344150 -4970.07+275819646+276326705+276855122+277396439+277941483+278478227+278993110 -4963.32+279475373+279921633+280337176+280732683+281118118+281498002+281870828 -4956.57+282232737+282582539+282924861+283268821+283623963+283995969+284384709 -4949.82+284785793+285194279+285607866+286027904+286457299+286896554+287340805 -4943.07+287780165+288203216+288602344+288977259+289334859+289685375+290037318 -4936.32+290394125+290754448+291115987+291479905+291852255+292240968+292650194 -4929.56+293076391+293509249+293937560+294355720+294767008+295181627+295611233 -4922.81+296063229+296538303+297032250+297540524+298062359+298600927+299159743 -4916.06+299737751+300326755+300913353+301483695+302027999+302542723+303029903 -4909.31+303494720+303943464+304382032+304816053+305251528+305695356+306154913 -4902.56+306637149+307147433+307689437+308264537+308871516+309505616+310158934 -4895.81+310822007+311486057+312144450+312792429+313425789+314039879+314629571 -4889.05+315190687+315721430+316223624+316702478+317165127+317619079+318070705 -4882.30+318524843+318984901+319452520+319927435+320407689+320890164+321371317 -4875.55+321847183+322313083+322764049+323196484+323610335+324010230+324404521 -4868.80+324802697+325212267+325637512+326080864+326545998+327040474+327575764 -4862.05+328163393+328809785+329512747+330262073+331044337+331848461+332668711 -4855.30+333503234+334349257+335197632+336030464+336823275+337551613+338198719 -4848.55+338760654+339245799+339669748+340048205+340392535+340709453+341004801 -4841.79+341287567+341571495+341872923+342205992+342578660+342991665+343440675 -4835.04+343920162+344425582+344953708+345501131+346063543+346636766+347218236 -4828.29+347807557+348404438+349005126+349600527+350178059+350726830+351243458 -4821.54+351734444+352213298+352694979+353190301+353704285+354238385+354794638 -4814.79+355378487+355998234+356661227+357368946+358114017+358881134+359651608 -4808.04+360409182+361143996+361853513+362540165+363207625+363857375+364487510 -4801.29+365094463+365676964+366239507+366792508+367348205+367914608+368491082 -4794.53+369069116+369637422+370188282+370721960+371246174+371771604+372305784 -4787.78+372849215+373396638+373942502+374487176+375038591+375608113+376202695 -4781.03+376819137+377443986+378059133+378649432+379207984+379736454+380240922 -4774.28+380726067+381191835+381634764+382052634+382449383+382835849+383226413 -4767.53+383633868+384065721+384523718+385005743+385508387+386028980+386566516 -4760.78+387120892+387691154+388273840+388862871+389451531+390034958+390611908 -4754.02+391184179+391754521+392324492+392893776+393461157+394026741+394593328 -4747.27+395165362+395745352+396330338+396911121+397476176+398017651+398535574 -4740.52+399037293+399532218+400025768+400516517+400998172+401464813+401915540 -4733.77+402354980+402789398+403220696+403643060+404044911+404414670+404747185 -4727.02+405046524+405324293+405595057+405871796+406163311+406474281+406806848 -4720.27+407162016+407540209+407940262+408358687+408789880+409227126+409664054 -4713.52+410095458+410517531+410928052+411326625+411714995+412096491+412475265 -4706.76+412855070+413237809+413677592+414054727+414420549+414784176+415147935 -4700.01+415501755+415826233+416105853+416343469+416563559+416800408+417077252 -4693.26+417396417+417739716+418080214+418395044+418671570+418905829+419095494 -4686.51+419235356+419321108+419360654+419381775+419424598+419521347+419678867 -4679.76+419879898+420103081+420342415+420608316+420910774+421240486+421563537 -4673.01+421832663+422006705+422066313+422019684+421896924+421735729+421566921 -4666.26+421408396+421269485+421159122+421089548+421072234+421112704+421208766 -4659.50+421352435+421531023+421726345+421916009+422080641+422214979+422334328 -4652.75+422470675+422658648+422921218+423262587+423672183+424134383+424637001 -4646.00+425173797+425742473+426341073+426965419+427604780+428236528+428824554 -4639.25+429326775+429712607+429979565+430156620+430289557+430416229+430549060 -4632.50+430676922+430785962+430885460+431019825+431254586+431642479+432190854 -4625.75+432851997+433544939+434193605+434757365+435236139+435653295+436030826 -4618.99+436373755+436669209+436896991+437045313+437123241+437164213+437216844 -4612.24+437325752+437511901+437763791+438044970+438313329+438540979+438724564 -4605.49+438882587+439042936+439228397+439446875+439689381+439935720+440164534 -4598.74+440363001+440531756+440683223+440833606+440991497+441149996+441288061 -4591.99+441383144+441428003+441439396+441449970+441487057+441553644+441626787 -4585.24+441674501+441677250+441638841+441579814+441521421+441472386+441425835 -4578.49+441366121+441279972+441165169+441032443+440899401+440781558+440685813 -4571.73+440609974+440547828+440494616+440450048+440417587+440402176+440407595 -4564.98+440434478+440481214+440546717+440634056+440751978+440912116+441123007 -4558.23+441382351+441670985+441953010+442186687+442344631+442433740+442499931 -4551.48+444038077+443576802+442690997+441830409+441534717+442105693+443364986 -4544.73+444684284+445300991+444736306+443049548+440767442+438541429+436782505 -4537.98+435543435+434716260+434328710+434631221+435867992+437948486+440368657 -4531.22+442492159+442859223+443932182+444641408+445001863+445101044+445043946 -4524.47+444910771+444754335+444623063+444575376+444665305+444913468+445290153 -4517.72+445723778+446129646+446440828+446624545+446679687+446621902+446469483 -4510.97+446239189+445947146+445605935+445217195+444770724+444258167+443691209 -4504.22+443105483+442541618+442020311+441534030+441061442+440588166+440111454 -4497.47+439628318+439122872+438568576+437945261+437252741+436509309+435737144 -4490.72+434949831+434153082+433354191+432564896+431789294+431008590+430180120 -4483.96+429257491+428213848+427046625+425761506+424356826+422825130+421168084 -4477.21+419402843+417552325+415632128+413651528+411627630+409595088+407599236 -4470.46+405679911+403858973+402136503+400490241+398878926+397255292+395589444 -4463.71+393888755+392194991+390558669+389010796+387552137+386162445+384814440 -4456.96+383480788+382138439+380778406+379419774+378113561+376927570+375919268 -4450.21+375117761+374526536+374141788+373964680+373995820+374217178+374581096 -4443.46+375018739+375463624+375870471+376218106+376500659+376720723+376891937 -4436.70+377042611+377207374+377407930+377638488+377871768+378084695+378281417 -4429.95+378494053+378758711+379089005+379466853+379856729+380226145+380556756 -4423.20+380841689+381079464+381273913+381436773+381585465+381733205+381879306 -4416.45+382007141+382091360+382107908+382042087+381891281+381662996+381369366 -4409.70+381023582+380642271+380252315+379891913+379597755+379383374+379226594 -4402.95+379076554+378875655+378577558+378151520+377579064+376854929+375993998 -4396.19+375033859+374024818+373011282+372016251+371038508+370061663+369067002 -4389.44+368042972+366985424+365893327+364766574+363610956+362444605+361294221 -4382.69+360178994+359093003+358004210+356872197+355669183+354386548+353026885 -4375.94+351594978+350100078+348563572+347021143+345515326+344083102+342746040 -4369.19+341503665+340331235+339186031+338023514+336815556+335558721+334269558 -4362.44+332975504+331710053+330511664+329417822+328453877+327625354+326922683 -4355.69+326333389+325851284+325473488+325190934+324979488+324798388+324597384 -4348.93+324331986+323982289+323564922+323128284+322731984+322423445+322225930 -4342.18+322139781+322148769+322224318+322329234+322425877+322487733+322506660 -4335.43+322491302+322458550+322424899+322401056+322390165+322389266+322393707 -4328.68+322401558+322414326+322432222+322448928+322450223+322416943+322327728 -4321.93+322160506+321895108+321521039+321046706+320503089+319936554+319395026 -4315.18+318914349+318510463+318177420+317891297+317620084+317337318+317032719 -4308.43+316712576+316389604+316070836+315752411+315423571+315073874+314697480 -4301.67+314293990+313869247+313437130+313019022+312638186+312309558+312031101 -4294.92+311782172+311529568+311237286+310875086+310424042+309879050+309249549 -4288.17+308558669+307839820+307129933+306456340+305827632+305230619+304637042 -4281.42+304015735+303342750+302602622+301787977+300897334+299938676+298933705 -4274.67+297916496+296925932+295996721+295153739+294409408+293761271+293190163 -4267.92+292662169+292136897+291578821+290965762+290290266+289555769+288771681 -4261.16+287950692+287110405+286277098+285487829+284785264+284204878+283759860 -4254.41+283434853+283194039+282998850+282823354+282664115+282537602+282474424 -4247.66+282512410+282685765+283011274+283477835+284043181+284648309+285238054 -4240.91+285777467+286257272+286690685+287106757+287543369+288037237+288610248 -4234.16+289257513+289945750+290624313+291242157+291761401+292163675+292452071 -4227.41+292650247+292794921+292920192+293040995+293149428+293228492+293270998 -4220.66+293287123+293294842+293306446+293324157+293345146+293364813+293373615 -4213.90+293353816+293286356+293164575+293001212+292821064+292644749+292477632 -4207.15+292310331+292125848+291906128+291635311+291302532+290904435+290447231 -4200.40+289946226+289423915+288906204+288417094+287972738+287578024+287228116 -4193.65+286912414+286616987+286325261+286019498+285685530+285317832+284920475 -4186.90+284504799+284087458+283689995+283336704+283047488+282829302+282673367 -4180.15+282559621+282463322+282357771+282214947+282009925+281728745+281373789 -4173.40+280960360+280507015+280029193+279540956+279060411+278611005+278216503 -4166.64+277893611+277645501+277458162+277298950+277122027+276882930+276558320 -4159.89+276158637+275723902+275304737+274940581+274646793+274413724+274212878 -4153.14+274005106+273749753+273415705+272991306+272488001+271934498+271364791 -4146.39+270807561+270281047+269792810+269340126+268909224+268476657+268014352 -4139.64+267497486+266911390+266254556+265538801+264789236+264043505+263346200 -4132.89+262736049+262231343+261824496+261490925+261205437+260951590+260716433 -4126.13+260479081+260208634+259876675+259472842+259007708+258500253+257964963 -4119.38+257411143+256853780+256319865+255840800+255436174+255104321+254824992 -4112.63+254568951+254304161+253998716+253625176+253169558+252638233+252056842 -4105.88+251461600+250890280+250376322+249944020+249601460+249335190+249112774 -4099.13+248896173+248657341+248385229+248080523+247745921+247382108+246991254 -4092.38+246581235+246163709+245746765+245329582+244903941+244459875+243990220 -4085.63+243492308+242969627+242434364+241907718+241413242+240965026+240560056 -4078.87+240181388+239808826+239428228+239032774+238618843+238183897+237728940 -4072.12+237262670+236800126+236355294+235931925+235520638+235104539+234670042 -4065.37+234213393+233738558+233250057+232750083+232243791+231746883+231283282 -4058.62+230872021+230514527+230195256+229894674+229601018+229308789+229007917 -4051.87+228677173+228289623+227824542+227276220+226655998+225990044+225315790 -4045.12+224673971+224097761+223603629+223191205+222850443+222568048+222327180 -4038.37+222103945+221868496+221592934+221260605+220869777+220430826+219961000 -4031.61+219479729+219004933+218548826+218114118+217692521+217268188+216826646 -4024.86+216366350+215905168+215474491+215101770+214792016+214522851+214256237 -4018.11+213956937+213603236+213186212+212706460+212176034+211623628+211092898 -4011.36+210629773+210265088+210004541+209830340+209710237+209606312+209481225 -4004.61+209304210+209058241+208742459+208386326+208036524+207736484+207509270 -3997.86+207350163+207233074+207130972+207017742+206879690+206717305+206538888 -3991.10+206354762+206174601+206006823+205858489+205735108+205640447+205575763 -3984.35+205538874+205522789+205514846+205497968+205454682+205373846+205256254 -3977.60+205115823+204974203+204852276+204763589+204712955+204698694+204715783 -3970.85+204758250+204821176+204903056+205006069+205132345+205278353+205432834 -3964.10+205581050+205711053+205816062+205891135+205930284+205929054+205890249 -3957.35+205825340+205748615+205668732+205585980+205496064+205394980+205280851 -3950.60+205154206+205020396+204893909+204797134+204752381+204771321+204850002 -3943.84+204971520+205111250+205240592+205331221+205363127+205333270+205258039 -3937.09+205167330+205094703+205070132+205115268+205239270+205435372+205681327 -3930.34+205944148+206187064+206376530+206489285+206519023+206478936+206398563 -3923.59+206315851+206267979+206283205+206374046+206533734+206737739+206952502 -3916.84+207146106+207297533+207399978+207459376+207489563+207506508+207522276 -3910.09+207541586+207559389+207561702+207529968+207451274+207329003+207186034 -3903.33+207056481+206971998+206951723+207000044+207108675+207259217+207426888 -3896.58+207585704+207713354+207794917+207826994+207822408+207809217+207819276 -3889.83+207870531+207957777+208060487+208159773+208245314+208309112+208338784 -3883.08+208323704+208269474+208201552+208149424+208124127+208111306+208085216 -3876.33+208030656+207949226+207845789+207711570+207525540+207274232+206968257 -3869.58+206638902+206322169+206049660+205849620+205746051+205744677+205818732 -3862.83+205911145+205957775+205915375+205775750+205560193+205302289+205032662 -3856.07+204772934+204538887+204340882+204179211+204040009+203901085+203743115 -3849.32+203555433+203330889+203061341+202743365+202390193+202036214+201726301 -3842.57+201495162+201351096+201273869+201225627+201166057+201062119+200892716 -3835.82+200651095+200346984+200005588+199660743+199342860+199069122+198840467 -3829.07+198646031+198469663+198291167+198081968+197804251+197422953+196926971 -3822.32+196343148+195727762+195142314+194632467+194223624+193924879+193729399 -3815.57+193607922+193511503+193387567+193200110+192935981+192595272+192180125 -3808.81+191696751+191163351+190609016+190059147+189522496+188993709+188469243 -3802.06+187956608+187470036+187021040+186617550+186269559+185989582+185781387 -3795.31+185628017+185491934+185328373+185100221+184785219+184375121+183872662 -3788.56+183290306+182649227+181976480+181299781+180641057+180010935+179407326 -3781.81+178818877+178231486+177633124+177015690+176375047+175712411+175034681 -3775.06+174351902+173672071+172998425+172331493+171671500+171016887+170360199 -3768.30+169687954+168987530+168255941+167501988+166722976+165848854+165069577 -3761.55+164412321+163825986+163219946+162525020+161730610+160881918+160043350 -3754.80+159255681+158508310+157745066+156895528+155917559+154829295+153719447 -3748.05+153039893+152145417+151261026+150411223+149611870+148855010+148109635 -3741.30+147339901+146525983+145673656+144804926+143940162+143086712+142241575 -3734.55+141401566+140571391+139762495+138984566+138235504+137498086+136745838 -3727.80+135953636+135105512+134195241+133222229+132188366+131097379+129955797 -3721.04+128772119+127554434+126307644+125032716+123728432+122394278+121031271 -3714.29+119640403+118219756+116763622+115265141+113721338+112135835+110516583 -3707.54+108870230+107198202+105498425+103770740+102020732+100258326+98492490 -3700.79+96726753+94958849+93183682+91396250+89592673+87770334+85929003+84072941 -3693.07+82212432+80363745+78547235+76784882+75097807+73504434+72018806+70649488 -3685.36+69399222+68265888+67244171+66326935+65505179+64767080+64097571+63479999 -3677.64+62898555+62339084+61788059+61232719+60664347+60083670+59502946+58941077 -3669.93+58413175+57922479+57461317+57020199+56597115+56199004+55834723+55506233 -3662.21+55204799+54914675+54621366+54319433+54014800+53720030+53446216+53197002 -3654.49+52967898+52750237+52537039+52326512+52120729+51921337+51726313+51530426 -3646.78+51328583+51119882+50909857+50708205+50522717+50353264+50190880+50023149 -3639.06+49842426+49651869+49465724+49304017+49184211+49114415+49091708+49106062 -3631.35+49145257+49196096+49241503+49259194+49226419+49130150+48975831+48786910 -3623.63+48593363+48417804+48269922+48151977+48066634+48018280+48006794+48021102 -3615.91+48040035+48041046+48009352+47939962+47832154+47683112+47487474+47243596 -3608.20+46960957+46662473+46379096+46139974+45963133+45850765+45790076+45758543 -3600.48+45731448+45689137+45621651+45529707+45422570+45314352+45219335+45147807 -3592.77+45104227+45089444+45105053+45155053+45241698+45358431+45486157+45597510 -3585.05+45668089+45687026+45657823+45588893+45483391+45338896+45157214+44954967 -3577.33+44764490+44622347+44553298+44561410+44632829+44744311+44868188+44972262 -3569.62+45021125+44984990+44852985+44640323+44381956+44115851+43867351+43644858 -3561.90+43446804+43271272+43120053+42996500+42902623+42839323+42809188+42817822 -3554.19+42870935+42968040+43097042+43233759+43347449+43409410+43400535+43315424 -3546.47+43163200+42966078+42756082+42569094+42436345+42375531+42385830+42449288 -3538.75+42537393+42618777+42664723+42652193+42566186+42402748+42171935+41898283 -3531.04+41616628+41363489+41167127+41040633+40981576+40978040+41016482+41085320 -3523.32+41172446+41260941+41329207+41357799+41338918+41281656+41208942+41147926 -3515.61+41119598+41132581+41182726+41257316+41340947+41420487+41487551+41539153 -3507.89+41577853+41611824+41653695+41717114+41811725+41939127+42091847+42254889 -3500.17+42408289+42530566+42604093+42621384+42587406+42513890+42408540+42269183 -3492.46+42090535+41878933+41660139+41470418+41336817+41263172+41234157+41233486 -3484.74+41258869+41319790+41422126+41556405+41700712+41834138+41948356+42048871 -3477.03+42147143+42251525+42362888+42474853+42575855+42652318+42693453+42696341 -3469.31+42667868+42621212+42568744+42515767+42458266+42384654+42280232+42133550 -3461.59+41943495+41723871+41501696+41308661+41170309+41099855+41100030+41168855 -3453.88+41302221+41490399+41713915+41946007+42162716+42353422+42523247+42685272 -3446.16+42849447+43016653+43180075+43328786+43448721+43522386+43532203+43468556 -3438.45+43337149+43159238+42964115+42779787+42628187+42525584+42483431+42505824 -3430.73+42584611+42696840+42807985+42880841+42886954+42816048+42678855+42501241 -3423.01+42312026+42131118+41965257+41813985+41680803+41580109+41533913+41560376 -3415.30+41661900+41820408+42003818+42180929+42335354+42468195+42587700+42695155 -3407.58+42779615+42826450+42832903+42816728+42809211+42837701+42908098+43002920 -3399.87+43093229+43153122+43168107+43136373+43066788+42975495+42879737+42792085 -3392.15+42713470+42631250+42527070+42390135+42225824+42053180+41894335+41765521 -3384.43+41676610+41636903+41658239+41749486+41906673+42108834+42324579+42524083 -3376.72+42687552+42806637+42882645+42926410+42958266+43002246+43072805+43162724 -3369.00+43244726+43289621+43287787+43255815+43222597+43207037+43206281+43203819 -3361.29+43188927+43168725+43161508+43178707+43212846+43241924+43245783+43221338 -3353.57+43185765+43167089+43190543+43269620+43404097+43582286+43784731+43988947 -3345.85+44174990+44330030+44449416+44534233+44588095+44616215+44626544+44630479 -3338.14+44640686+44666763+44711480+44769744+44830080+44877585+44897880+44881679 -3330.42+44828497+44747470+44654306+44566046+44496280+44452700+44436813+44445487 -3322.71+44473970+44519449+44582743+44666175+44768627+44881484+44988862+45071545 -3314.99+45111444+45094483+45012904+44868611+44675923+44459890+44248765+44064271 -3307.27+43916584+43807471+43737900+43712510+43737222+43813065+43931880+44077217 -3299.56+44229656+44373256+44499505+44606817+44696244+44766010+44808443+44811896 -3291.84+44767477+44676835+44555165+44425965+44310475+44219089+44150208+44096547 -3284.13+44053704+44025188+44021068+44052105+44123986+44235557+44381315+44555767 -3276.41+44755674+44978015+45214828+45449185+45656756+45813232+45903191+45924463 -3268.69+45885896+45801307+45684502+45548366+45407420+45280146+45187495+45146489 -3260.98+45162429+45225570+45316463+45417673+45523706+45640713+45775884+45926083 -3253.26+46076552+46210953+46322970+46418123+46504090+46579021+46628931+46636690 -3245.55+46594521+46508808+46393744+46260140+46109522+45938067+45746228+45545055 -3237.83+45353197+45186926+45052149+44945504+44863311+44810105+44799257+44844992 -3230.11+44952644+45114884+45316205+45541183+45778987+46020512+46251099+46446252 -3222.40+46576698+46621959+46583452+46485731+46362645+46237634+46113384+45979317 -3214.68+45830467+45681141+45560912+45496298+45493542+45535965+45597084+45658272 -3206.97+45717217+45781676+45855533+45929588+45984862+46005788+45993066+45965330 -3199.25+45947742+45956231+45989908+46036699+46085619+46134271+46186025+46241728 -3191.53+46296103+46342951+46383884+46431102+46500538+46600654+46726242+46861869 -3183.82+46991594+47106711+47205842+47288574+47349161+47376722+47362536+47308793 -3176.10+47231162+47152785+47094393+47068381+47080092+47132623+47228740+47367318 -3168.39+47538075+47719899+47884689+48003946+48054815+48024499+47913700+47738713 -3160.67+47530224+47327482+47169343+47085993+47094852+47200965+47399333+47675731 -3152.95+48005274+48351419+48669903+48919503+49075454+49137003+49123502+49061963 -3145.24+48975702+48881798+48796254+48738503+48728451+48777694+48882720+49026026 -3137.52+49183890+49335060+49465925+49571672+49654678+49721113+49776820+49824117 -3129.81+49861812+49887982+49902454+49906172+49898139+49874038+49828872+49762043 -3122.09+49680177+49595241+49519712+49463067+49431624+49430342+49463153+49531056 -3114.37+49630392+49754506+49898863+50065087+50258967+50482101+50723156+50956910 -3106.66+51153361+51290749+51362742+51375962+51343210+51280245+51208040+51154871 -3098.94+51151610+51220497+51365177+51569532+51805387+52042742+52256564+52429486 -3091.23+52552893+52627906+52665337+52683573+52704734+52750650+52839240+52980629 -3083.51+53173420+53403383+53647307+53881093+54087784+54260458+54398788+54503278 -3075.79+54572334+54604494+54603222+54579356+54548061+54522020+54505591+54494182 -3068.08+54478586+54450665+54406665+54347704+54279672+54214698+54172766+54179732 -3060.36+54259569+54423672+54664288+54957601+55275140+55595131+55905454+56197894 -3052.65+56461283+56681882+56851443+56974639+57066419+57139066+57189608+57199831 -3044.93+57150168+57036284+56874620+56693857+56521302+56377058+56279899+56256003 -3037.21+56337281+56545529+56873047+57276113+57688941+58049931+58322427+58497156 -3029.50+58579452+58575639+58491063+58338895+58148385+57961562+57818726+57743785 -3021.78+57740216+57798140+57903949+58043904+58201253+58353652+58477113+58555318 -3014.07+58586576+58582003+58556203+58517590+58465534+58395385+58305932+58202410 -3006.35+58092424+57978804+57856070+57713948+57545933+57358105+57168771+56997657 -2998.63+56853472+56730137+56614502+56498013+56381836+56271698+56166827+56055877 -2990.92+55924374+55766272+55588463+55404006+55220276+55032984+54831593+54612438 -2983.20+54389153+54191564+54052604+53992199+54010005+54092156+54224957+54402112 -2975.49+54618438+54857502+55088542+55280090+55420455+55526032+55628160+55748766 -2967.77+55885774+56020250+56137201+56240783+56351073+56486647+56648773+56819868 -2960.05+56976582+57106783+57218104+57332636+57472089+57643521+57834936+58023563 -2952.34+58191023+58334560+58465627+58596330+58723789+58825540+58870980+58841215 -2944.62+58741572+58596641+58432439+58261582+58083859+57899244+57718124+57557240 -2936.91+57425228+57313128+57200314+57070503+56922809+56768480+56619068+56479291 -2929.19+56351053+56242482+56170931+56156260+56210919+56335510+56521394+56754721 -2921.47+57017000+57283660+57526477+57721970+57860617+57949145+58003467+58038287 -2913.76+58062025+58080661+58105297+58155053+58251035+58405741+58616519+58868595 -2906.04+59145690+59440304+59756726+60105874+60496781+60930888+61401686+61898296 -2898.33+62409928+62929152+63452818+63980382+64510174+65036668+65551604+66049272 -2890.61+66531502+67007414+67488111+67981424+68491820+69023739+69582113+70167290 -2882.89+70770053+71374462+71969136+72557459+73156449+73784377+74448256+75142289 -2875.18+75857562+76592759+77355832+78155714+78992147+79852114+80714942+81561420 -2867.46+82380519+83170429+83935140+84679616+85406196+86113189+86795485+87446662 -2859.74+88061823+88639705+89182482+89693182+90172901+90621210+91040322+91439549 -2852.03+91834963+92242947+92671839+93118047+93569058+94010890+94434093+94835362 -2844.31+95215399+95575847+95917739+96241886+96550914+96851271+97152507+97462334 -2836.60+97779510+98090567+98375242+98618845+98822381+99001326+99173405+99345967 -2828.88+99514828+99675633+99836649+100019937+100249199+100535500+100873493 -2822.13+101249914+101654197+102080572+102521591+102962307+103383104+103768500 -2815.38+104113464+104421897+104701596+104962494+105219486+105492689+105799971 -2808.63+106145656+106516447+106890879+107255683+107614618+107982045+108368094 -2801.87+108769992+109177096+109582542+109988199+110397650+110805448+111195152 -2795.12+111549409+111862282+112142338+112403996+112656877+112903124+113143172 -2788.37+113381264+113623201+113867888+114102768+114309826+114479480+114621167 -2781.62+114761955+114933776+115158598+115440861+115769794+116126429+116489646 -2774.87+116838959+117156485+117430607+117659922+117854239+118031109+118209976 -2768.12+118407914+118638115+118908615+119219321+119558523+119903541+120227662 -2761.37+120510414+120744594+120935434+121094131+121231536+121356489+121477557 -2754.61+121603846+121742453+121894185+122052512+122207904+122354825+122495745 -2747.86+122638000+122786560+122939137+123089150+123233507+123377255+123528630 -2741.11+123688067+123842244+123970846+124062202+124124705+124184367+124270846 -2734.36+124404576+124594188+124843236+125156665+125539840+125991242+126495419 -2727.61+127022184+127533313+127993715+128382270+128696967+128950880+129161044 -2720.86+129338020+129484623+129604436+129711692+129832337+129994325+130215684 -2714.11+130500181+130841603+131229034+131646640+132070219+132468607+132813109 -2707.35+133088671+133296853+133448901+133556753+133630913+133685315+133740416 -2700.60+133817657+133928204+134065476+134208154+134330636+134413719+134449973 -2693.85+134444951+134415728+134386320+134378799+134403211+134453264+134512001 -2687.10+134563468+134600568+134623222+134630306+134614882+134568041+134487853 -2680.35+134384721+134277993+134186729+134120737+134077041+134042016+133997065 -2673.60+133925904+133821001+133686121+133533187+133375058+133219916+133071410 -2666.84+132933582+132814854+132725864+132672493+132649747+132641935+132628784 -2660.09+132593429+132527370+132430502+132307967+132166267+132011205+131848992 -2653.34+131689331+131547512+131442344+131390057+131397974+131462500+131572611 -2646.59+131715738+131881639+132062461+132250883+132440296+132628269+132821317 -2639.84+133035896+133292042+133601571+133958273+134338025+134709846+135049524 -2633.09+135345771+135595956+135798996+135955737+136076964+136189665+136330718 -2626.34+136530308+136797742+137121044+137478260+137848443+138213075+138551127 -2619.58+138838175+139054128+139193422+139267979+139300599+139315706+139334805 -2612.83+139377218+139460248+139595445+139783497+140012495+140261068+140503363 -2606.08+140713566+140870215+140960725+140984079+140949054+140868232+140752570 -2599.33+140609694+140444679+140259019+140047282+139795683+139486867+139109005 -2592.58+138662216+138157484+137610193+137034274+136356240+135761830+135199405 -2585.83+134646894+134086597+133509171+132918066+132332445+131773893+131236753 -2579.08+130668975+129990611+129142632+128123824+126986101+125797586+124609374 -2572.32+123446235+122316522+121225072+120183821+119219348+118369717+117664337 -2565.57+117095462+116605731+116109670+115539962+114887332+114208928+113601314 -2558.82+113153931+112906230+112828487+112837713+112845233+112814411+112795815 -2552.07+112915191+113314914+114079810+115190834+116537188+117977184+119405470 -2545.32+120780874+122106305+123392138+124641187+125859177+127061108+128251593 -2538.57+128990623+129874499+130663173+131364297+131985154+132534468+133025295 -2531.81+133473498+133891210+134281059+134636228+134946035+135201970+135399406 -2525.06+135535343+135606014+135607508+135539387+135408406+135229659+135024161 -2518.31+134813983+134617182+134444805+134300885+134183980+134088077+134001823 -2511.56+133907956+133785830+133618080+133397871+133131785+132836913+132534798 -2504.81+132247090+131993759+131790679+131644353+131545886+131471316+131391405 -2498.06+131286422+131155508+131014231+130883290+130778267+130708494+130682946 -2491.31+130715195+130820046+131003499+131357728+131557967+131775334+131999045 -2484.55+132217893+132414338+132596059+132759197+132901386+133021555+133119625 -2477.80+133196205+133252258+133288962+133307717+133310321+133298915+133275824 -2471.05+133243271+133203012+133156118+133102919+133043099+132975890+132900209 -2464.30+132814999+132719387+132613003+132496072+132369572+132235234+132095424 -2457.55+131952971+131811046+131672994+131542251+131422214+131315791+131224924 -2450.80+131150195+131090652+131043916+131006367+130973390+130939700+130899785 -2444.05+130848621+130782219+130698040+130594934+130472914+130332774+130175887 -2437.29+130003973+129819014+129623045+129418168+129206510+128990438+128772542 -2430.54+128555492+128341812+128133551+127932018+127737728+127550416+127369184 -2423.79+127192816+127020069+126849820+126681263+126513751+126346806+126179862 -2417.04+126012204+125842973+125671006+125495114+125313974+125126464+124931804 -2410.29+124729808+124520965+124306294+124087050+123864647+123640592+123416602 -2403.54+123194702+122977136+122766205+122564183+122373170+122195018+122030915 -2396.78+121880928+121743590+121616019+121494237+121373724+121249894+121118675 -2390.03+120976830+120822521+120655577+120477755+120292505+120104373+119918476 -2383.28+119739583+119571463+119416255+119274145+119143442+119020894+118902443 -2376.53+118783925+118661972+118534850+118402944+118268924+118137163+118012659 -2369.78+117899507+117799626+117711746+117631175+117550036+117458270+117345159 -2363.03+117201278+117020429+116801053+116546757+116265789+115969728+115671736 -2356.28+115384702+115119595+114884107+114681965+114512814+114373057+114256786 -2349.52+114157037+114066474+113977907+113884529+113780088+113659363+113518483 -2342.77+113355345+113170029+112964900+112744678+112515878+112286060+112062534 -2336.02+111851127+111655092+111474508+111306400+111145456+110985319+110820146 -2329.27+110645978+110461733+110269425+110073555+109880110+109695111+109523363 -2322.52+109367421+109227162+109099902+108981107+108865372+108747476+108623282 -2315.77+108490352+108348104+108197575+108040867+107880439+107718404+107556099 -2309.02+107393926+107231284+107066930+106899206+106726618+106548188+106363790 -2302.26+106174179+105980899+105786107+105592246+105401888+105217252+105039721 -2295.51+104869380+104704854+104543507+104381797+104215923+104042330+103858184 -2288.76+103661838+103453339+103234624+103009485+102782985+102560839+102348494 -2282.01+102150563+101970137+101808446+101664585+101535640+101417116+101303780 -2275.26+101190457+101072786+100947759+100814010+100671821+100523089+100370981 -2268.51+100219461+100072699+99934469+99807605+99693668+99592881+99504129 -2261.75+99425124+99352523+99282155+99209336+99129148+99036794+98927978+98799482 -2254.04+98649482+98477931+98286450+98078248+97857709+97629953+97400413+97174130 -2246.32+96955078+96745680+96546764+96357965+96178048+96005308+95837736+95673389 -2238.61+95510786+95349261+95188945+95030512+94874683+94722092+94573335+94429130 -2230.89+94290279+94157586+94031760+93913243+93802226+93698565+93601460+93509416 -2223.17+93420195+93331092+93239313+93142260+93037786+92924463+92801789+92670300 -2215.46+92531408+92387052+92239398+92090528+91942371+91796673+91654901+91518190 -2207.74+91387441+91263412+91146792+91038247+90938181+90846362+90761898+90683151 -2200.03+90608039+90534294+90459618+90381862+90299229+90210464+90115136+90013655 -2192.31+89907106+89796837+89684036+89569437+89453233+89335086+89214130+89089011 -2184.59+88958004+88819331+88671650+88514639+88349227+88177657+88003106+87829355 -2176.88+87660237+87499240+87348975+87210738+87084093+86966831+86855273+86744745 -2169.16+86630332+86507684+86373703+86227060+86068449+85900533+85727734+85555681 -2161.45+85390528+85238122+85103242+84989113+84897182+84827132+84776966+84743157 -2153.73+84721006+84705251+84690871+84673748+84651021+84621157+84583806+84539674 -2146.01+84490395+84438438+84386839+84338842+84297459+84265190+84243897+84234923 -2138.30+84239245+84257603+84290447+84337797+84399131+84473371+84558938+84653810 -2130.58+84755462+84860822+84966386+85068508+85163750+85249257+85322843+85383060 -2122.87+85429221+85461404+85480509+85488102+85486113+85476584+85461437+85442464 -2115.15+85421271+85399324+85377938+85358298+85341565+85328930+85321667+85321039 -2107.43+85328328+85344750+85371350+85408794+85457122+85515653+85583087+85657737 -2099.72+85737772+85821344+85906726+85992412+86077311+86160810+86242742+86323241 -2092.00+86402450+86480398+86556825+86631244+86702920+86770928+86834073+86890886 -2084.29+86939756+86979130+87007936+87025885+87033604+87032626+87025297+87014690 -2076.57+87004520+86998790+87001275+87014968+87041488+87080927+87131839+87191455 -2068.85+87256046+87321325+87382989+87437272+87481403+87513950+87534959+87546015 -2061.14+87550218+87551982+87556589+87569515+87595651+87638587+87700013+87779388 -2053.42+87873778+87978067+88085448+88188350+88279574+88353338+88406094+88437041 -2045.71+88448183+88444000+88430763+88415458+88404574+88403219+88414539+88439546 -2037.99+88477380+88525596+88580724+88638971+88696855+88751706+88801858+88846558 -2030.27+88885714+88919741+88949367+88975457+88998950+89020798+89042091+89064249 -2022.56+89089091+89118763+89155559+89201706+89259161+89329436+89413430+89511276 -2014.84+89622279+89745019+89877447+90017138+90161462+90307768+90453565+90596672 -2007.13+90735160+90867377+90991895+91107696+91214271+91311654+91400526+91482068 -1999.41+91558007+91630555+91702310+91776015+91854286+91939278+92032591+92135122 -1991.69+92247024+92367755+92496205+92630979+92770729+92914299+93060731+93209224 -1983.98+93358993+93509172+93658716+93806297+93950178+94088395+94219151+94341296 -1976.26+94454785+94560904+94662358+94762993+94867407+94980254+95105353+95244879 -1968.55+95398785+95564593+95737782+95912379+96081927+96240452+96383461+96508785 -1960.83+96616946+96711052+96796163+96878274+96963186+97055468+97157668+97269861 -1953.11+97389641+97512467+97632576+97743983+97841478+97921428+97982108+98023860 -1945.40+98048986+98061357+98065837+98067542+98071144+98080343+98097565+98123860 -1937.68+98158714+98200050+98244248+98286436+98320999+98342054+98344116+98322816 -1929.97+98275717+98202944+98107617+97995589+97874878+97754583+97643778+97550248 -1922.25+97479464+97433753+97412129+97410510+97422630+97441055+97458297+97467892 -1914.53+97465090+97447485+97414978+97369524+97314528+97254179+97192859+97134704 -1906.82+97083276+97041240+97010226+96990711+96982047+96982569+96989786+97000743 -1899.10+97012387+97022075+97027970+97029384+97026707+97021295+97015110+97010292 -1891.39+97008673+97011250+97017674+97025948+97032530+97032913+97022286+96996374 -1883.67+96951998+96887664+96803782+96702685+96588232+96465162+96338140+96211177 -1875.95+96086970+95966741+95850134+95735396+95619893+95500728+95375510+95242864 -1868.24+95102657+94956041+94805234+94653219+94503284+94358551+94221411+94093186 -1860.52+93973955+93862542+93756628+93652907+93547389+93435725+93313612+93177094 -1852.81+93022791+92848214+92651835+92433292+92193250+91933265+91655654+91363439 -1845.09+91060187+90749950+90436951+90125379+89819120+89521526+89235264+88961930 -1837.37+88701831+88453893+88215682+87983796+87754296+87523176+87286882+87042697 -1829.66+86789082+86525865+86254261+85976677+85696238+85416156+85139093+84866770 -1821.94+84599792+84337645+84078677+83820198+83558686+83290202+83010979+82717765 -1814.23+82408097+82080281+81733327+81366964+80981542+80577980+80157612+79721939 -1806.51+79272613+78811424+78340402+77861813+77378083+76891497+76404080+75917382 -1798.79+75432403+74949446+74467995+73986644+73503270+73015331+72520306+72016076 -1791.08+71500995+70974032+70434671+69882847+69318862+68743188+68156325+67558723 -1783.36+66950832+66333193+65706646+65072434+64432372+63788874+63144934+62503894 -1775.65+61869219+61244099+60631139+60032037+59447448+58876948+58319182+57772109 -1767.93+57233178+56699540+56168268+55636580+55104918+54568111+54024501+53473562 -1760.21+52959215+52387459+51804432+51216578+50628291+50039571+49448248+48853904 -1752.50+48260171+47673702+46990906+46481416+45991646+45500173+44995071+44477439 -1744.78+43954762+43400518+42844785+42276703+41705384+41141575+40593435+40063924 -1737.07+39551612+39053917+38568974+38094209+37624730+37154685+36680845+36204067 -1729.35+35726281+35246562+34761222+34267979+33769485+33272184+32780689+32292108 -1721.63+31795009+31274476+30720269+30132679+29522835+28908824+28310449+27744232 -1713.92+27218535+26729817+26262739+25795677+25309505+24795413+24258446+23715380 -1706.20+23187435+22690546+22228068+21790262+21360728+20924665+20473817+20007867 -1698.49+19534886+19071291+18638992+18258815+17943992+17698554+17521227+17411404 -1690.77+17373609+17418860+17562583+17819501+18197226+18692276+19291249+19975936 -1683.05+20728771+21536657+22393139+23298594+24256989+25269624+26330098+27424597 -1675.34+28536623+29652190+30762699+31865693+32964287+34065398+35176257+36300313 -1667.62+37434569+38569926+39694048+40794720+41861299+42884142+43854180+44654637 -1659.91+45519256+46337450+47110468+47839321+48524750+49166662+49763309+50311429 -1652.19+50807686+51250205+51638757+51973828+52255896+52485889+52666437+52802758 -1644.47+52902447+52974332+53027157+53068556+53104182+53137096+53167284+53191603 -1636.76+53204463+53199206+53169950+53113186+53028439+52917749+52784528+52632129 -1629.04+52462667+52276076+52069926+51840138+51582697+51295784+50981883+50649068 -1621.33+50310980+49985219+49690506+49443384+49255345+49130834+49066322+49050610 -1613.61+49066391+49093221+49111580+49106954+49072910+49012307+48936613+48863473 -1605.89+48813176+48804651+48851826+48961411+49132552+49358269+49627960+49929962 -1598.18+50253532+50589998+50933089+51278678+51624145+51967626+52307631+52642987 -1590.46+52973076+53297891+53617651+53931972+54239066+54535429+54816235+55076362 -1582.75+55311774+55521056+55706756+55876132+56040565+56213371+56406558+56627705 -1575.03+56878096+57152369+57439311+57723504+57987752+58216373+58398689+58531342 -1567.31+58618621+58670742+58701049+58723260+58749463+58788565+58844877+58917504 -1559.60+59001379+59089669+59176459+59257929+59332294+59399741+59462819+59526724 -1551.88+59598353+59684257+59788573+59911822+60050627+60198149+60344918+60480703 -1544.17+60596986+60689413+60758988+60811168+60853516+60893028+60934193+60978067 -1536.45+61022264+61061618+61089830+61101474+61093874+61067843+61027234+60977644 -1528.73+60925033+60874478+60829322+60790867+60758518+60730135+60702709+60673387 -1521.02+60640708+60605128+60568173+60530624+60491065+60445579+60411247+60355577 -1513.30+60265549+60133941+59962251+59760434+59544375+59309409+59098267+58929360 -1505.59+58681653+58653871+58671112+58695233+58688552+58630668+58525599+58390382 -1497.87+58231890+58032901+57765653+57421593+57029497+56645522+56328333+56117951 -1490.15+56029793+56057252+56177731+56357219+56557298+56743460+56890361+56981909 -1482.44+57007682+56962751+56852554+56695225+56515942+56335682+56162994+55992039 -1474.72+55804978+55576118+55276921+54882221+54372502+53732051+52946014+52003319 -1467.01+50907656+49687436+48396581+47107134+45893770+44816473+43907464+43168229 -1459.29+42580130+42124860+41804875+41644599+41670408+41881531+42240297+42687317 -1451.57+43170142+43660752+44148708+44617725+45030618+45336034+45488724+45466173 -1443.86+45272764+44935694+44496587+43998365+43466970+42896843+42257222+41522120 -1436.14+40706080+39872088+39104882+38471281+37999466+37680935+37485498+37375314 -1428.43+37314688+37274267+37233575+37184754+37140940+37140751+37235121+37454917 -1420.71+37779673+38139771+38460743+38719999+38960083+39243479+39590222+39953932 -1412.99+40251939+40423615+40473754+40472677+40518117+40686069+40997228+41410558 -1405.28+41840300+42188420+42379017+42380487+42207096+41905483+41534147+41139232 -1397.56+40728433+40258006+39651427+38849562+37857610+36749937+35627797+34566463 -1389.85+33594657+32713926+31931732+31277469+30788911+30485812+30352782+30344373 -1382.13+30406569+30502544+30628499+30810713+31084339+31464259+31923862+32395624 -1374.41+32796543+33064277+33180637+33171500+33090024+32997475+32945125+32958637 -1366.70+33031446+33136958+33256367+33401295+33608115+33904606+34279395+34679372 -1358.98+35035902+35293079+35417646+35390581+35197794+34833896+34316195+33693715 -1351.27+33039075+32425065+31899146+31470365+31114949+30794178+30471329+30118246 -1343.55+29713920+29245627+28717369+28158265+27616730+27139825+26751893+26530116 -1335.83+26335748+26168985+26035979+25942882+25887040+25852413+25813606+25746632 -1328.12+25639789+25496646+25328711+25143780+24939538+24708501+24450470+24181705 -1320.40+23932111+23731280+23593329+23510037+23455260+23396674+23307780+23174253 -1312.69+22993376+22769058+22507139+22214332+21900106+21578279+21265070+20974560 -1304.97+20660415+20410795+20250000+20268021+20520136+21056091+21909042+23099147 -1297.25+24645841+26572449+28890806+31570151+34510740+37546432+40485088+43169643 -1289.54+45525993+47571987+49389191+51079901+52732192+54402085+56110351+57848517 -1281.82+59589359+61296524+62925769+64415520+65675679+66590417+67041613+66942757 -1274.11+66266097+65054716+63425690+61564520+59697712+58034864+56700366+55697239 -1266.39+54933757+54299792+53739036+53263160+52904556+52657800+52464170+52249232 -1258.67+51975057+51659917+51351445+51080809+50837566+50585379+50305244+50029771 -1250.96+49836168+49799556+49943473+50228776+50592275+51000668+51470177+52035562 -1243.24+52700739+53418371+54114638+54733360+55261156+55719332+56139230+56544478 -1235.53+56946560+57345291+57726491+58058866+58296832+58391605+58306500+58028890 -1227.81+57573107+56973371+56271566+55507370+54714592+53920020+53138048+52361666 -1220.09+51559088+50685399+49706789+48623201+47476018+46339184+45300692+44440937 -1212.38+43810686+43414082+43207034+43118308+43086158+43088900+43149124+43309113 -1204.66+43594349+43987639+44428510+44838730+45159081+45368148+45479096+45549506 -1196.95+45585466+45531395+45299284+44825424+44111463+43219075+42228422+41197698 -1189.23+40151764+39097696+38042180+36993794+35951684+34903319+33841806+32791940 -1181.51+32009683+31251990+30649558+30190854+29851923+29585348+29340935+29082562 -1173.80+28800160+28512399+28258867+28085036+28026818+28101491+28307948+28633441 -1166.08+29060410+29568477+30131489+30714295+31273950+31766585+32156204+32420640 -1158.37+32552494+32557176+32451641+32264835+32037052+31814302+31636943+31526838 -1150.65+31480073+31469612+31456296+31401917+31278417+31071948+30784679+30436525 -1142.93+30064558+29714730+29424603+29204433+29030067+28855781+28641946+28379849 -1135.22+28096581+27836500+27633093+27491720+27394773+27323185+27276196+27273659 -1127.50+27339402+27478865+27668382+27864834+28029789+28152565+28257380+28389765 -1119.79+28589296+28863245+29175243+29455372+29627355+29639647+29485080+29199853 -1112.07+28843760+28473680+28124631+27806712+27515007+27242280+26984374+26736712 -1104.35+26489089+26228346+25951045+25677526+25436491+25310091+25070670+25359923 -1096.64+25825200+26333362+26804806+27249253+27721049+28247120+28795422+29307169 -1088.92+29749543+30139502+30527184+30963247+31474396+32054145+32667789+33271576 -1081.21+33839288+34378734+34924486+35517376+36188824+36960184+37842318+38823852 -1073.49+39859572+40880300+41823395+42658418+43388226+44035171+44633893+45229551 -1065.77+45857546+46498159+47044231+47332795+47247115+46820489+46252497+45806964 -1058.06+45653931+45766401+45944408+45947458+45641738+45064708+44375328+43739030 -1050.34+43233901+42833930+42460077+42050355+41592001+41112517+40653681+40239403 -1042.63+39867750+39518632+39169193+38806670+38430956+38047034+37657174+37263911 -1034.91+36884258+36559615+36341554+36261019+36301017+36413259+36566048+36780448 -1027.19+37116736+37615808+38242372+38879545+39380507+39636138+39612793+39350111 -1019.48+38939706+38496862+38134673+37925549+37872367+37909636+37950460+37935261 -1011.76+37854177+37737180+37636519+37621799+37780919+38196462+38909191+39872029 -1004.05+40959781+42022005+42958818+43759110+44465987+45107459+45650306+46025366 -996.33+46195922+46207652+46164594+46160675+46229496+46359327+46535880+46763765 -988.61+47050079+47384758+47749098+48138505+48560753+49005999+49427579+49769689 -980.90+50020436+50231796+50477082+50785544+51123874+51448121+51777089+52213116 -973.18+52881227+53830247+54970382+56094041+56969564+57454768+57567609+57474296 -965.47+57397783+57493250+57777878+58131322+58387204+58433946+58256355+57898894 -957.75+57397446+56746983+55930890+54975454+53969766+53024071+52197537+51454144 -950.03+50685032+49786511+48742620+47656751+46709331+46066864+45798442+45847567 -942.32+46069808+46307933+46459608+46503921+46481115+46447987+46441352+46468999 -934.60+46517892+46552649+46496037+46222435+45603650+44596743+43300921+41916272 -926.89+40624237+39625848+38612624+37755165+37078886+36606697+36344557+36263200 -919.17+36294246+36348757+36351734+36274408+36145105+36030431+35994844+36059393 -911.45+36182031+36271070+36224930+35976691+35520267+34906970+34219363+33540940 -903.74+32937182+32450491+32100788+31882833+31760215+31665518+31517044+31252000 -896.02+30860597+30397994+29960292+29633144+29442812+29342300+29243836+29077170 -888.31+28835417+28578441+28391628+28326881+28363939+28417623+28387604+28220266 -880.59+27943508+27652687+27457222+27404230+27504152+27718438+27983265+28248372 -872.87+28519384+28797725+29082265+29368305+29648697+29916405+30167284+30401880 -865.16+30625625+30847281+31076071+31388116+31793102+32269530+32789911+33344971 -857.44+33885927+34494645+35229591+35905666+37112523+38529730+40141756+41901627 -849.73+43728575+45518489+47172941+48633437+49897514+51007749+52024658+53003832 -842.01+53984857+54985964+55998165+56759730+57383931+57806976+58054901+58268495 -834.29+58654954+59375124+60429797+61628318+62676310+63337764+63564423+63504424 -826.58+63385636+63360741+63422432+63441174+63283686+62916603+62421182+61921545 -818.86+61496221+61147298+60877716+60606826+60355703+60079751+59725217+59255649 -811.15+58683252+58079187+57544723+57154490+56909691+56738088+56545053+56279694 -803.43+55967673+55686243+55502605+55422741+55388615+55323878+55193049+55030050 -795.71+54917077+54930123+55088734+55341318+55591219+55742792+55738199+55567402 -788.00+55253847+54832584+54337963+53804778+53272154+52775617+52325484+51887974 -780.28+51393306+50778869+50044435+49277236+48618006+48179458+47967057+47857548 -772.57+47655357+47197429+46440747+45479932+44485323+43605538+42898287+42331263 -764.85+41836738+41383304+40823761+40387397+40115913+39909179+39674669+39494448 -757.13+39591292+40091293+40831559+41443782+41654521+41496674+41222417+41031834 -749.42+40909451+40712031+40377075+40026034+39888121+40152025+40874306+41975163 -741.70+43281529+44585330+45701548+46506674+46933587+46957920+46621451+46088606 -733.99+45653405+45628498+46157830+47090318+48022222+48492084+48213307+47210564 -726.27+45806686+44481626+43718048+43574861+43902111+44417427+44788036+44759970 -718.55+44244155+43311677+42123326+40858892+39689974+38780371+38683110+38799410 -710.84+38991338+39235635+39497752+39749315+39982342+40208105+40441753+40685931 -703.12+40925377+41134818+41293723+41397420+41458073+41494760+41521717+41543105 -695.41+41557895+41570240+41415230+41356556+41425513+41777698+42356349+43119412 -687.69+43946937+44649865+45037471+45022136+44681213+44211988+43815338+43605442 -679.97+43601050+43777957+44120361+44621495+45237088+45846341+46276426+46395938 -672.26+46212836+45895277+45699846+45858947+46485034+47514744+48731617+49873939 -664.54+50749795+51275364+51455337+51353397+51066581+50699081+50327861+49980455 -656.83+49641015+49288847+48957096+48774896+48947818+48838113+48967838+49261765 -649.11+49548714+49629364+49387826+48857595+48766163+48677777+48609124+48564999 -641.39+48393713+48224773+48056963+47888413+47716132+47537064+47350763+47315355 -633.68+47281576+47199455+47137649+47085289+47016501+46912714+46788934+46672171 -625.96+46570770+46471150+46350102+46183775+46046235+45968271+45952457+45986276 -618.25+46038999+46092426+46145888+46199322+46252715+46306076+46359482+46414016 -610.53+46527107+46762146+47367860+48110571+48813661+49185602+49373092+49482463 -602.81+49571242+49675667+49768566+49077520 -##END= diff --git a/man/OpenSpecy-package.Rd b/man/OpenSpecy-package.Rd index 03bcc458..2064357d 100644 --- a/man/OpenSpecy-package.Rd +++ b/man/OpenSpecy-package.Rd @@ -78,6 +78,9 @@ Other contributors: \item Vesna Teofilovic (\href{https://orcid.org/0000-0002-3557-1482}{ORCID}) [contributor] \item Laura A. T. Markley (\href{https://orcid.org/0000-0003-0620-8366}{ORCID}) [contributor] \item Shreyas Patankar [contributor, data contributor] + \item Katherine Lasdin [contributor] + \item National Renewable Energy Laboratory [funder] + \item Possibility Lab [funder] } } diff --git a/man/adj_intens.Rd b/man/adj_intens.Rd index 5007a61c..5db1a9db 100644 --- a/man/adj_intens.Rd +++ b/man/adj_intens.Rd @@ -38,7 +38,7 @@ Many of the Open Specy functions will assume that the spectrum is in absorbance units. For example, see \code{\link{subtr_baseline}()}. To run those functions properly, you will need to first convert any spectra from transmittance or reflectance to absorbance using this function. -The transmittance adjustment uses the \eqn{log10(1 / T)} +The transmittance adjustment uses the \eqn{log(1 / T)} calculation which does not correct for system and particle characteristics. The reflectance adjustment uses the Kubelka-Munk equation \eqn{(1 - R)^2 / 2R}. We assume that the reflectance intensity diff --git a/man/data_norm.Rd b/man/data_norm.Rd index 117a1dde..97bcc171 100644 --- a/man/data_norm.Rd +++ b/man/data_norm.Rd @@ -4,7 +4,6 @@ \alias{adj_res} \alias{conform_res} \alias{adj_neg} -\alias{make_rel} \alias{mean_replace} \alias{is_empty_vector} \title{Normalization and conversion of spectral data} @@ -15,8 +14,6 @@ conform_res(x, res = 5) adj_neg(y, na.rm = FALSE) -make_rel(y, na.rm = FALSE) - mean_replace(y, na.rm = TRUE) is_empty_vector(x) diff --git a/man/gen_OpenSpecy.Rd b/man/gen_OpenSpecy.Rd index bca31257..f43b4f3f 100644 --- a/man/gen_OpenSpecy.Rd +++ b/man/gen_OpenSpecy.Rd @@ -6,6 +6,8 @@ \alias{plot.OpenSpecy} \alias{lines.OpenSpecy} \alias{summary.OpenSpecy} +\alias{as.data.frame.OpenSpecy} +\alias{as.data.table.OpenSpecy} \title{Generic Open Specy Methods} \usage{ \method{head}{OpenSpecy}(x, ...) @@ -17,6 +19,10 @@ \method{lines}{OpenSpecy}(x, ...) \method{summary}{OpenSpecy}(object, ...) + +\method{as.data.frame}{OpenSpecy}(x, ...) + +\method{as.data.table}{OpenSpecy}(x, ...) } \arguments{ \item{x}{an \code{OpenSpecy} object.} @@ -29,9 +35,11 @@ \code{head()}, \code{print()}, and \code{summary()} return a textual representation of an \code{OpenSpecy} object. \code{plot()} and \code{lines()} return a plot. +\code{as.data.frame()} and \code{as.data.table()} convert \code{OpenSpecy} +objects into tabular data. } \description{ -Methods to visualize \code{OpenSpecy} objects. +Methods to visualize and convert \code{OpenSpecy} objects. } \details{ \code{head()} shows the first few lines of an \code{OpenSpecy} object. @@ -56,7 +64,9 @@ plot(raman_hdpe) \seealso{ \code{\link[utils]{head}()}, \code{\link[base]{print}()}, \code{\link[base]{summary}()}, \code{\link[graphics]{matplot}()}, and -\code{\link[graphics]{matlines}()} +\code{\link[graphics]{matlines}()}, +\code{\link[base]{as.data.frame}()}, +\code{\link[data.table]{as.data.table}()} } \author{ Zacharias Steinmetz, Win Cowger diff --git a/man/make_rel.Rd b/man/make_rel.Rd new file mode 100644 index 00000000..7be18bfb --- /dev/null +++ b/man/make_rel.Rd @@ -0,0 +1,49 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/make_rel.R +\name{make_rel} +\alias{make_rel} +\alias{make_rel.default} +\alias{make_rel.OpenSpecy} +\title{Make spectral intensities relative} +\usage{ +make_rel(x, ...) + +\method{make_rel}{default}(x, na.rm = FALSE, ...) + +\method{make_rel}{OpenSpecy}(x, na.rm = FALSE, ...) +} +\arguments{ +\item{x}{a numeric vector or an \R OpenSpecy object} + +\item{na.rm}{logical. Should missing values be removed?} + +\item{\ldots}{further arguments passed to \code{make_rel()}.} +} +\value{ +\code{make_rel()} return numeric vectors (if vector provided) or an +\code{OpenSpecy} object with the normalized intensity data. +} +\description{ +\code{make_rel()} converts intensities \code{x} into relative values between +0 and 1 using the standard normalization equation. +If \code{na.rm} is \code{TRUE}, missing values are removed before the +computation proceeds. +} +\details{ +\code{make_rel()} is used to retain the relative height proportions between +spectra while avoiding the large numbers that can result from some spectral +instruments. +} +\examples{ +make_rel(c(-1000, -1, 0, 1, 10)) + +} +\seealso{ +\code{\link[base]{min}()} and \code{\link[base]{round}()}; +\code{\link{adj_intens}()} for log transformation functions; +\code{\link{conform_spec}()} for conforming wavenumbers of an +\code{OpenSpecy} object to be matched with a reference library +} +\author{ +Win Cowger, Zacharias Steinmetz +} diff --git a/man/manage_spec.Rd b/man/manage_spec.Rd index 50bd2612..063fe063 100644 --- a/man/manage_spec.Rd +++ b/man/manage_spec.Rd @@ -22,7 +22,7 @@ sample_spec(x, ...) \method{sample_spec}{default}(x, ...) -\method{sample_spec}{OpenSpecy}(x, ...) +\method{sample_spec}{OpenSpecy}(x, size = 1, prob = NULL, ...) } \arguments{ \item{x}{a list of \code{OpenSpecy} objects.} @@ -35,6 +35,10 @@ spectra having all the same wavenumber range.} \item{res}{defaults to \code{NULL}, the resolution you want the output wavenumbers to be.} +\item{size}{the number of spectra to sample.} + +\item{prob}{probabilities to use for the sampling.} + \item{\ldots}{further arguments passed to submethods.} } \value{ diff --git a/man/match_spec.Rd b/man/match_spec.Rd index f90eb855..81572e90 100644 --- a/man/match_spec.Rd +++ b/man/match_spec.Rd @@ -35,6 +35,7 @@ match_spec(x, ...) library, na.rm = T, top_n = NULL, + order = NULL, add_library_metadata = NULL, add_object_metadata = NULL, fill = NULL, @@ -83,6 +84,9 @@ when calculating correlations. Default is \code{TRUE}.} \item{top_n}{integer; specifying the number of top matches to return. If \code{NULL} (default), all matches will be returned.} +\item{order}{an \code{OpenSpecy} used for sorting, ideally the unprocessed +one; \code{NULL} skips sorting.} + \item{add_library_metadata}{name of a column in the library metadata to be joined; \code{NULL} if you don't want to join.} @@ -134,19 +138,23 @@ matrix as a named vector. } \examples{ data("test_lib") + unknown <- read_extdata("ftir_ldpe_soil.asp") |> read_any() |> conform_spec(range = test_lib$wavenumber, res = spec_res(test_lib)) |> process_spec() -matches <- cor_spec(unknown, test_lib) +cor_spec(unknown, test_lib) test_lib_extract <- filter_spec(test_lib, logic = grepl("polycarbonate", test_lib$metadata$polymer_class, ignore.case = TRUE) ) -matches2 <- cor_spec(unknown, library = test_lib_extract) +cor_spec(unknown, library = test_lib_extract) + +match_spec(unknown, test_lib, add_library_metadata = "sample_name", + top_n = 1) } \seealso{ diff --git a/man/process_spec.Rd b/man/process_spec.Rd index 383f8c0d..c0d4317b 100644 --- a/man/process_spec.Rd +++ b/man/process_spec.Rd @@ -27,6 +27,7 @@ process_spec(x, ...) smooth_intens = TRUE, smooth_intens_args = list(polynomial = 3, window = 11, derivative = 1, abs = TRUE), make_rel = TRUE, + make_rel_args = list(na.rm = TRUE), ... ) } @@ -75,6 +76,11 @@ to the spectra.} \item{make_rel}{logical; if \code{TRUE} spectra are automatically normalized with \code{\link{make_rel}()}.} +\item{make_rel_args}{named list of arguments passed to +\code{\link{make_rel}()}.} + +\item{na.rm}{Whether to allow NA or set all NA values to} + \item{\ldots}{further arguments passed to subfunctions.} } \value{ diff --git a/man/sig_noise.Rd b/man/sig_noise.Rd index 3688cba2..99d0e588 100644 --- a/man/sig_noise.Rd +++ b/man/sig_noise.Rd @@ -10,7 +10,7 @@ sig_noise(x, ...) \method{sig_noise}{default}(x, ...) -\method{sig_noise}{OpenSpecy}(x, metric = "sig_over_noise", na.rm = TRUE, ...) +\method{sig_noise}{OpenSpecy}(x, metric = "run_sig_over_noise", na.rm = TRUE, ...) } \arguments{ \item{x}{an \code{OpenSpecy} object.} @@ -19,8 +19,11 @@ sig_noise(x, ...) Options include \code{"sig"} (mean intensity), \code{"noise"} (standard deviation of intensity), \code{"sig_times_noise"} (absolute value of signal times noise), \code{"sig_over_noise"} (absolute value of signal / -noise), or \code{"tot_sig"} (total signal = signal * number of data -points).} +noise), \code{"run_sig_over_noise"} (absolute value of signal / +noise where signal is estimated as the max intensity and noise is +estimated as the height of a low intensity region.), +\code{"log_tot_sig"} (sum of the inverse log intensities, useful for spectra in log units), +or \code{"tot_sig"} (sum of intensities).} \item{na.rm}{logical; indicating whether missing values should be removed when calculating signal and noise. Default is \code{TRUE}.} diff --git a/man/smooth_intens.Rd b/man/smooth_intens.Rd index afed0c95..2d2fc55f 100644 --- a/man/smooth_intens.Rd +++ b/man/smooth_intens.Rd @@ -14,8 +14,8 @@ smooth_intens(x, ...) x, polynomial = 3, window = 11, - derivative = 0, - abs = FALSE, + derivative = 1, + abs = TRUE, make_rel = TRUE, ... ) diff --git a/tests/testthat/test-adj_intens.R b/tests/testthat/test-adj_intens.R index 39ed63e4..f290aa83 100644 --- a/tests/testthat/test-adj_intens.R +++ b/tests/testthat/test-adj_intens.R @@ -1,10 +1,16 @@ data("raman_hdpe") + raman_hdpe$spectra$intensity2 <- raman_hdpe$spectra$intensity * 2 +raman_hdpe$metadata <- rbind(raman_hdpe$metadata, raman_hdpe$metadata) + +test_that("adj_intens() handles input errors correctly", { + adj_intens(1:1000) |> expect_error() +}) test_that("adj_intens() works as expected", { - expect_silent(adj <- adj_intens(raman_hdpe)) - expect_silent(adj_intens(raman_hdpe, type = "reflectance")) - expect_silent(adj_intens(raman_hdpe, type = "transmittance")) + adj <- adj_intens(raman_hdpe) |> expect_silent() + adj_intens(raman_hdpe, type = "reflectance") |> expect_silent() + adj_intens(raman_hdpe, type = "transmittance") |> expect_silent() expect_equal( cor(raman_hdpe$spectra$intensity, adj$spectra$intensity), 1, ignore_attr = F) @@ -13,6 +19,7 @@ test_that("adj_intens() works as expected", { ignore_attr = F) expect_s3_class(adj, "OpenSpecy") + expect_true(check_OpenSpecy(adj)) expect_equal(nrow(adj$spectra), nrow(raman_hdpe$spectra)) expect_equal(adj$wavenumber, raman_hdpe$wavenumber) diff --git a/tests/testthat/test-adj_range.R b/tests/testthat/test-adj_range.R index 6e09a892..96c00bdf 100644 --- a/tests/testthat/test-adj_range.R +++ b/tests/testthat/test-adj_range.R @@ -1,5 +1,17 @@ library(data.table) +test_that("flatten_range() error handling", { + test <- as_OpenSpecy(x = 1:10, spectra = data.table(V1 = 1:10)) + + expect_s3_class(test, "OpenSpecy") + expect_true(check_OpenSpecy(test)) + + expect_error(flatten_range(test)) + expect_error(flatten_range(test, min = c(1000), + max = c(2000, 3000))) + expect_error(flatten_range(test, min = c(2000), max = c(1000))) +}) + test_that("restrict_range() provides correct range", { test_noise <- as_OpenSpecy(x = seq(400,4000, by = 10), spectra = data.table(intensity = rnorm(361))) @@ -11,7 +23,9 @@ test_that("restrict_range() provides correct range", { max = c(1500, 2500)) |> expect_silent() - is_OpenSpecy(single_range) |> expect_true() + check_OpenSpecy(single_range) |> expect_true() + check_OpenSpecy(double_range) |> expect_true() + expect_identical(single_range$wavenumber, seq(1000,2000, by = 10)) expect_identical(double_range$wavenumber, c(seq(1000,1500, by = 10), seq(2000,2500, by = 10))) @@ -23,6 +37,8 @@ test_that("flatten_range() function test", { make_rel = F) |> expect_silent() + expect_true(check_OpenSpecy(flat_sam)) + expect_equal(flat_sam$spectra$V1[4:5], c(4.5, 4.5)) expect_equal(flat_sam$spectra$V1[7:10], c(8.5, 8.5, 8.5, 8.5)) @@ -30,6 +46,8 @@ test_that("flatten_range() function test", { flat_hdpe <- flatten_range(raman_hdpe, min = c(500, 1000), max = c(700, 1500)) |> expect_silent() + expect_true(check_OpenSpecy(flat_hdpe)) + expect_equal(flat_hdpe$spectra$intensity[1:50], make_rel(raman_hdpe$spectra$intensity)[1:50]) expect_equal(flat_hdpe$spectra$intensity[60:100] |> unique() |> round(6), @@ -39,6 +57,7 @@ test_that("flatten_range() function test", { flat_map <- flatten_range(tiny_map, min = c(1000, 2000), max = c(1200, 2400), make_rel = F) |> expect_silent() + expect_true(check_OpenSpecy(flat_map)) expect_false(all.equal(flat_map$spectra, tiny_map$spectra) |> isTRUE()) expect_equal(flat_map$spectra[1:20], tiny_map$spectra[1:20]) @@ -47,12 +66,3 @@ test_that("flatten_range() function test", { |> as.numeric(), c(-0.8694, -1.246, -0.8304, -1.1909, -0.7857)) }) - -test_that("flatten_range() error handling", { - test <- as_OpenSpecy(x = 1:10, spectra = data.table(V1 = 1:10)) - - expect_error(flatten_range(test)) - expect_error(flatten_range(test, min = c(1000), - max = c(2000, 3000))) - expect_error(flatten_range(test, min = c(2000), max = c(1000))) -}) diff --git a/tests/testthat/test-as_OpenSpecy.R b/tests/testthat/test-as_OpenSpecy.R index 47f19d15..9383e82d 100644 --- a/tests/testthat/test-as_OpenSpecy.R +++ b/tests/testthat/test-as_OpenSpecy.R @@ -42,6 +42,15 @@ test_that("as_OpenSpecy() generates OpenSpecy objects", { expect_s3_class(ost, "OpenSpecy") expect_s3_class(osl, "OpenSpecy") + expect_true(check_OpenSpecy(osf)) + expect_true(check_OpenSpecy(ost)) + expect_true(check_OpenSpecy(osl)) + + expect_true(is_OpenSpecy(osf)) + expect_true(is_OpenSpecy(ost)) + expect_true(is_OpenSpecy(osl)) + expect_false(is_OpenSpecy(df)) + expect_equal(names(osf), c("wavenumber", "spectra", "metadata")) expect_equal(names(ost), c("wavenumber", "spectra", "metadata")) expect_equal(names(osl), c("wavenumber", "spectra", "metadata")) @@ -50,11 +59,6 @@ test_that("as_OpenSpecy() generates OpenSpecy objects", { expect_equal(ost$spectra, osf$spectra) expect_equal(ost$wavenumber, osf$wavenumber) expect_equal(ost$metadata, osf$metadata) - - expect_true(is_OpenSpecy(osf)) - expect_true(is_OpenSpecy(ost)) - expect_true(is_OpenSpecy(osl)) - expect_false(is_OpenSpecy(df)) }) test_that("check_OpenSpecy() work as expected", { @@ -85,6 +89,8 @@ test_that("check_OpenSpecy() work as expected", { test_that("'OpenSpecy' objects are read correctly", { os <- as_OpenSpecy(df) + expect_true(check_OpenSpecy(os)) + expect_equal(range(os$wavenumber) |> round(2), c(301.04, 3198.12)) expect_equal(range(os$spectra) |> round(2), c(26, 816)) expect_length(os$wavenumber, 964) @@ -101,6 +107,8 @@ test_that("'OpenSpecy' objects are transcribed to and from 'hyperSpec' objects", openhyper <- as_OpenSpecy(hyper) expect_true(is_OpenSpecy(openhyper)) + expect_true(check_OpenSpecy(openhyper)) + expect_equal(openhyper$wavenumber, hyper@wavelength) expect_equal(unlist(openhyper$spectra$V1),unname(t(hyper$spc)[,1])) }) diff --git a/tests/testthat/test-conform_spec.R b/tests/testthat/test-conform_spec.R index 0cc8b172..3526dc1b 100644 --- a/tests/testthat/test-conform_spec.R +++ b/tests/testthat/test-conform_spec.R @@ -27,13 +27,24 @@ test_that("conform_spec() handles errors correctly", { }) test_that("conform_spec() conforms wavenumbers correctly", { - sam <- as_OpenSpecy(wn <- seq(1000, 2000, 5), + wn <- seq(1000, 2000, 5) + sam <- as_OpenSpecy(x = wn, data.table(intensity = rnorm(length(wn)))) + expect_true(check_OpenSpecy(sam)) + new_wavenumbers <- c(1100, 1200, 1300, 1400, 1500, 1600, 1700, 1800, 1900) + wider_wavenumbers <- seq(800, 2500, by = 100) conf_new <- conform_spec(sam, new_wavenumbers) |> expect_silent() + expect_true(check_OpenSpecy(conf_new)) + conf_roll <- conform_spec(sam, new_wavenumbers,res = NULL, type = "roll") |> expect_silent() + expect_true(check_OpenSpecy(conf_roll)) + + conf_wider <- conform_spec(x = sam,range = wider_wavenumbers, res = NULL) |> + expect_silent() + expect_true(check_OpenSpecy(conf_wider)) expect_equal(length(conf_new$wavenumber), length(conf_new$spectra[[1]])) expect_equal(range(conf_new$wavenumber), range(new_wavenumbers)) @@ -46,4 +57,8 @@ test_that("conform_spec() conforms wavenumbers correctly", { conform_spec(raman_hdpe)$spectra$intensity[c(63, 143, 283, 325, 402)] |> round(2) |> expect_equal(c(78.84, 65.00, 105.73, 109.41, 116.00)) + + expect_true(all(conf_wider$wavenumber <= max(sam$wavenumber)) & all(conf_wider$wavenumber >= min(sam$wavenumber))) + + }) diff --git a/tests/testthat/test-data_norm.R b/tests/testthat/test-data_norm.R index fd88d10b..c3a16071 100644 --- a/tests/testthat/test-data_norm.R +++ b/tests/testthat/test-data_norm.R @@ -4,13 +4,6 @@ test_that("adj_neg() give correct output", { c(1425.14, 542.03, 1.00, 458.30, 396.34)) }) -test_that("make_rel() give correct output", { - rel <- make_rel(c(-1000, -1, 0, 1, 10)) - expect_equal(range(rel), c(0, 1)) - expect_equal(round(make_rel(c(965.86, 82.75, -458.28, -0.98, -62.94)), 4), - c(1.0000, 0.3799, 0.0000, 0.3211, 0.2776)) -}) - test_that("adj_res() function test", { expect_equal(adj_res(seq(500, 4000, 4), 5), round(seq(500, 4000, 4)/5)*5) }) @@ -24,8 +17,3 @@ test_that("adj_neg() function test", { expect_equal(adj_neg(c(-1000, -1, 0, 1, 10)), c(1, 1000, 1001, 1002, 1011)) expect_equal(adj_neg(c(1, 2, 3)), c(1, 2, 3)) }) - -test_that("make_rel() function test", { - expect_equal(make_rel(c(1, 2, 3)), c(0, 0.5, 1)) - expect_equal(make_rel(c(10, 20, 30)), c(0, 0.5, 1)) -}) diff --git a/tests/testthat/test-def_features.R b/tests/testthat/test-def_features.R index 646c43b8..240bece1 100644 --- a/tests/testthat/test-def_features.R +++ b/tests/testthat/test-def_features.R @@ -1,8 +1,17 @@ map <- read_extdata("CA_tiny_map.zip") |> read_any() +test_that("def_features() handles input errors correctly", { + def_features(1:1000) |> expect_error() +}) + +test_that("collapse_spec() handles input errors correctly", { + collapse_spec(1:1000) |> expect_error() +}) + test_that("features are identified when given logical", { map$metadata$particles <- map$metadata$y == 0 id_map <- def_features(map, map$metadata$particles) + expect_true(check_OpenSpecy(id_map)) expect_length(unique(id_map$metadata$feature_id), 2) expect_equal(max(id_map$metadata$area, na.rm = T), 13) expect_equal(max(id_map$metadata$feret_max, na.rm = T), 13) @@ -13,6 +22,7 @@ test_that("features are identified when given logical", { test_that("particles are identified when given character", { map$metadata$particles <- ifelse(map$metadata$y == 1, "particle", "not_particle") id_map <- def_features(map, map$metadata$particles) + expect_true(check_OpenSpecy(id_map)) expect_length(unique(id_map$metadata$feature_id), 3) expect_equal(max(id_map$metadata$area, na.rm = T), 182) expect_equal(round(max(id_map$metadata$feret_max, na.rm = T)), 19) @@ -38,6 +48,8 @@ test_that("the original spectrum remains unmodified and metadata is amended", { id_map <- def_features(map,ifelse(map$metadata$x == 1, "particle", "not_particle")) + expect_true(check_OpenSpecy(id_map)) + expect_equal(id_map$wavenumber, map$wavenumber) expect_equal(id_map$spectra, map$spectra) expect_contains(id_map$metadata, map$metadata) @@ -50,15 +62,22 @@ test_that("the original spectrum remains unmodified and metadata is amended", { test_that("collapse particles returns expected values", { particles <- ifelse(map$metadata$y == 1, "particleA", "particleB") id_map <- def_features(map, particles) + expect_true(check_OpenSpecy(id_map)) + test_collapsed <- collapse_spec(id_map) + expect_true(check_OpenSpecy(test_collapsed)) + expect_equal(test_collapsed$metadata |> nrow(), 3) expect_equal(test_collapsed$metadata$feret_max |> round(2), c(13, 13, 18.69)) expect_equal(test_collapsed$metadata$centroid_x |> unique(), 6) particles <- map$metadata$y == 1 id_map <- def_features(map, particles) + expect_true(check_OpenSpecy(id_map)) + test_collapsed <- collapse_spec(id_map) + expect_true(check_OpenSpecy(test_collapsed)) expect_equal(test_collapsed$metadata |> nrow(), 2) expect_equal(test_collapsed$metadata$feret_max |> round(2), c(NA, 13)) diff --git a/tests/testthat/test-gen_OpenSpecy.R b/tests/testthat/test-gen_OpenSpecy.R index 766dceaa..b40f999c 100644 --- a/tests/testthat/test-gen_OpenSpecy.R +++ b/tests/testthat/test-gen_OpenSpecy.R @@ -1,5 +1,6 @@ library(data.table) data("raman_hdpe") +csv <- read_extdata("raman_hdpe.csv") |> read.csv() test_that("print() and summary() work as expected", { print(raman_hdpe) |> expect_output( @@ -8,7 +9,7 @@ test_that("print() and summary() work as expected", { summary(raman_hdpe) |> expect_output("Length.*964.*spectra.*26") }) -test_that("head returns the first few lines of the OpenSpecy object", { +test_that("head() returns the first few lines of the OpenSpecy object", { head <- head(raman_hdpe) nrow(head) |> expect_equal(6) @@ -16,7 +17,17 @@ test_that("head returns the first few lines of the OpenSpecy object", { raman_hdpe$spectra))) }) -test_that("ploting works without errors", { +test_that("plotting works without errors", { plot(raman_hdpe) |> expect_silent() }) +test_that("conversion to tables works as expected", { + df <- as.data.frame(raman_hdpe) |> expect_silent() + dt <- as.data.table(raman_hdpe) |> expect_silent() + + expect_s3_class(df, "data.frame") + expect_s3_class(dt, "data.table") + + expect_equal(df, csv) + expect_equal(dt, as.data.table(csv)) +}) diff --git a/tests/testthat/test-interactive_plots.R b/tests/testthat/test-interactive_plots.R index d1c1df78..dcb22960 100644 --- a/tests/testthat/test-interactive_plots.R +++ b/tests/testthat/test-interactive_plots.R @@ -1,10 +1,8 @@ map <- read_zip(read_extdata("CA_tiny_map.zip")) data("raman_hdpe") -test_that("heatmap_spec() generates 'plotly' object", { - heatmap_spec(map, z = map$metadata$y) |> - expect_silent() |> - expect_s3_class("plotly") +test_that("plotly_spec() handles input errors correctly", { + plotly_spec(1:1000) |> expect_error() }) test_that("plotly_spec() generates 'plotly' object", { @@ -13,6 +11,16 @@ test_that("plotly_spec() generates 'plotly' object", { expect_s3_class("plotly") }) +test_that("heatmap_spec() handles input errors correctly", { + heatmap_spec(1:1000) |> expect_error() +}) + +test_that("heatmap_spec() generates 'plotly' object", { + heatmap_spec(map, z = map$metadata$y) |> + expect_silent() |> + expect_s3_class("plotly") +}) + test_that("interactive_plot() generates 'plotly' object", { interactive_plot(map, x2 = raman_hdpe, select = 2) |> suppressWarnings() |> diff --git a/tests/testthat/test-io_spec.R b/tests/testthat/test-io_spec.R index dbc972fb..e0602e3b 100644 --- a/tests/testthat/test-io_spec.R +++ b/tests/testthat/test-io_spec.R @@ -29,7 +29,7 @@ test_that("read_spec() gives expected output", { jsn <- read_extdata("raman_hdpe.json") |> read_spec() |> expect_silent() rds <- read_extdata("raman_hdpe.rds") |> read_spec() |> expect_silent() - read_spec(read_extdata("raman_hdpe.csv")) |> expect_error() + read_extdata("raman_hdpe.csv") |> read_spec() |> expect_error() read_extdata("raman_hdpe.yml") |> read_spec(share = tmp) |> expect_message() read_extdata("raman_hdpe.json") |> read_spec(share = tmp) |> expect_message() @@ -38,6 +38,10 @@ test_that("read_spec() gives expected output", { expect_s3_class(yml, "OpenSpecy") expect_s3_class(jsn, "OpenSpecy") expect_s3_class(rds, "OpenSpecy") + + expect_true(check_OpenSpecy(yml)) + expect_true(check_OpenSpecy(jsn)) + expect_true(check_OpenSpecy(rds)) jsn$metadata$file_name <- yml$metadata$file_name <- rds$metadata$file_name <- NULL @@ -62,7 +66,7 @@ test_that("read_spec() and write_spec() work nicely together", { expect_equal(yml[1:2], raman_hdpe[1:2]) }) -test_that("as_hyperspec function", { +test_that("as_hyperspec() function works", { hyperspec_object <- as_hyperSpec(raman_hdpe) # Verify the class of the output diff --git a/tests/testthat/test-make_rel.R b/tests/testthat/test-make_rel.R new file mode 100644 index 00000000..e43ea4b1 --- /dev/null +++ b/tests/testthat/test-make_rel.R @@ -0,0 +1,21 @@ +test_that("make_rel() function test", { + expect_equal(make_rel(c(1, 2, 3)), c(0, 0.5, 1)) + expect_equal(make_rel(c(10, 20, 30)), c(0, 0.5, 1)) +}) + +test_that("make_rel() gives correct output", { + rel <- make_rel(c(-1000, -1, 0, 1, 10)) + expect_equal(range(rel), c(0, 1)) + expect_equal(round(make_rel(c(965.86, 82.75, -458.28, -0.98, -62.94)), 4), + c(1.0000, 0.3799, 0.0000, 0.3211, 0.2776)) +}) + +test_that("make_rel() gives correct output with OpenSpecy objects", { + data("raman_hdpe") + rel <- make_rel(raman_hdpe) |> expect_silent() + + expect_s3_class(rel, "OpenSpecy") + expect_true(check_OpenSpecy(rel)) + + expect_equal(range(rel$spectra), c(0, 1)) +}) diff --git a/tests/testthat/test-manage_lib.R b/tests/testthat/test-manage_lib.R index fc57f016..e2170609 100644 --- a/tests/testthat/test-manage_lib.R +++ b/tests/testthat/test-manage_lib.R @@ -30,6 +30,18 @@ test_that("check_lib() finds test library", { expect_warning() }) +test_that("load_lib() works with test library", { + skip_on_cran() + skip_if_offline(host = "api.osf.io") + + tl <- load_lib(type = "test", path = tmp) |> + expect_silent() + expect_true(check_OpenSpecy(tl)) + expect_type(tl, "list") + expect_s3_class(tl, "OpenSpecy") + expect_identical(tl, test_lib, ignore_attr = T) +}) + test_that("get_lib() downloads complete library", { skip_on_cran() skip_if_offline(host = "api.osf.io") @@ -50,15 +62,24 @@ test_that("check_lib() finds complete library", { expect_warning() }) -test_that("load_lib() works as expected", { +test_that("load_lib() works with complete library", { skip_on_cran() skip_if_offline(host = "api.osf.io") + skip_if_not(testthat:::on_ci(), "Not on CI") - tl <- load_lib(type = "test", path = tmp) |> + tl <- load_lib(type = "derivative", path = tmp) |> expect_silent() expect_type(tl, "list") expect_s3_class(tl, "OpenSpecy") - expect_identical(tl, test_lib, ignore_attr = T) + expect_true(check_OpenSpecy(tl)) +}) + +test_that("rm_lib() works as expected", { + skip_on_cran() + skip_if_offline(host = "api.osf.io") + + rm_lib(type = "test", path = tmp) |> + expect_silent() }) # Tidy up diff --git a/tests/testthat/test-manage_spec.R b/tests/testthat/test-manage_spec.R index 936f9e3b..46125d98 100644 --- a/tests/testthat/test-manage_spec.R +++ b/tests/testthat/test-manage_spec.R @@ -18,7 +18,8 @@ test_that("c_spec() merges identical files without range specification", { specs <- lapply(c(read_extdata("raman_hdpe.yml"), read_extdata("raman_hdpe.yml")), read_spec) same <- c_spec(specs) |> expect_silent() - + expect_true(check_OpenSpecy(same)) + expect_equal(same$wavenumber, raman_hdpe$wavenumber) expect_equal(same$spectra$intensity, raman_hdpe$spectra$intensity) }) @@ -26,6 +27,7 @@ test_that("c_spec() merges identical files without range specification", { test_that("c_spec() merges different files with common range", { diff <- c_spec(specs, range = "common", res = 5) |> expect_silent() + expect_true(check_OpenSpecy(diff)) diff$wavenumber[1:2] |> expect_equal(c(655, 660)) diff$spectra$intensity[1:2] |> round(2) |> expect_equal(c(53.87, 59.00)) @@ -36,6 +38,8 @@ test_that("c_spec() merges different files with specified range", { spec <- c_spec(specs, range = c(1000, 2000), res = 5) |> expect_silent() + expect_true(check_OpenSpecy(spec)) + spec$wavenumber |> expect_equal(seq(1000, 2000, 5)) spec$spectra$intensity |> expect_equal( @@ -47,6 +51,7 @@ test_that("sample_spec() returns a subset of the spectra", { tiny_map <- read_any(read_extdata("CA_tiny_map.zip")) sampled <- sample_spec(tiny_map, size = 5) expect_s3_class(sampled, "OpenSpecy") + expect_true(check_OpenSpecy(sampled)) expect_equal(ncol(sampled$spectra), 5) }) diff --git a/tests/testthat/test-match_spec.R b/tests/testthat/test-match_spec.R index c90c3a9f..baebfcc6 100644 --- a/tests/testthat/test-match_spec.R +++ b/tests/testthat/test-match_spec.R @@ -2,11 +2,11 @@ tmp <- file.path(tempdir(), "OpenSpecy-testthat") dir.create(tmp, showWarnings = F) -# Create test data for cor_spec function data("test_lib") -unknown <- read_extdata("ftir_ldpe_soil.asp") |> read_any() |> - conform_spec(range = test_lib$wavenumber, res = spec_res(test_lib)) |> +unknown <- read_extdata("ftir_ldpe_soil.asp") |> read_any() +preproc <- conform_spec(unknown, range = test_lib$wavenumber, + res = spec_res(test_lib)) |> process_spec(smooth_intens = T, make_rel = T) # Create a subset of test_lib for filtering @@ -14,19 +14,24 @@ test_lib_extract <- filter_spec( test_lib, logic = test_lib$metadata$polymer_class == "polycarbonates" ) -# Match_spec function with AI -test_that("match_spec returns correct structure with AI", { +test_that("ai_classify() handles input errors correctly", { + ai_classify(1:1000) |> expect_error() +}) + +test_that("match_spec() returns correct structure with AI", { skip_on_cran() skip_if_offline(host = "api.osf.io") get_lib("model", path = tmp) lib <- load_lib(type = "model", path = tmp) + expect_error(check_OpenSpecy(lib)) + set.seed(47) rn <- runif(n = length(unique(lib$variables_in))) fill <- as_OpenSpecy(as.numeric(unique(lib$variables_in)), spectra = data.frame(rn)) - matches <- match_spec(x = unknown, library = lib, na.rm = T, fill = fill) |> + matches <- match_spec(x = preproc, library = lib, na.rm = T, fill = fill) |> expect_silent() nrow(matches) |> expect_equal(1) names(matches) |> expect_contains(c("x", "y", "z", "value", "name")) @@ -34,13 +39,23 @@ test_that("match_spec returns correct structure with AI", { grepl("polyamide", matches$name) |> expect_true() }) +test_that("match_spec() handles input errors correctly", { + match_spec(1:1000) |> expect_error() +}) + # Match_spec function -test_that("match_spec returns correct structure", { - matches <- match_spec(x = unknown, library = test_lib, na.rm = T, top_n = 5, +test_that("match_spec() returns correct structure", { + matches <- match_spec(x = preproc, library = test_lib, na.rm = T, top_n = 5, add_library_metadata = "sample_name", add_object_metadata = "col_id") |> expect_silent() + order <- match_spec(x = preproc, library = test_lib, na.rm = T, top_n = 5, + order = unknown, + add_library_metadata = "sample_name", + add_object_metadata = "col_id") |> + expect_silent() + nrow(matches) |> expect_equal(5) names(matches) |> expect_contains(c("object_id", "library_id", "match_val")) round(matches$match_val, 2) |> expect_equal(c(0.67, 0.57, 0.55, 0.50, 0.43)) @@ -48,21 +63,28 @@ test_that("match_spec returns correct structure", { c("poly(ethylene)", "polystyrene", "poly(vinyl chloride)", "poly(dimethylsiloxane) (pdms)", NA) ) + expect_equal(matches, order) +}) + +test_that("cor_spec() handles input errors correctly", { + cor_spec(1:1000) |> expect_error() + restrict_range(raman_hdpe, 300, 310) |> cor_spec(test_lib) |> + expect_error() }) # Write the tests for cor_spec function -test_that("cor_spec returns a data.table with correct columns", { - matches <- cor_spec(unknown, library = test_lib) |> +test_that("cor_spec() returns a data.table with correct columns", { + matches <- cor_spec(preproc, library = test_lib) |> expect_silent() - unknown2 <- unknown - unknown2$wavenumber[1:3] <- unknown2$wavenumber[1:3] + 1 + preproc2 <- preproc + preproc2$wavenumber[1:3] <- preproc2$wavenumber[1:3] + 1 - matches2 <- cor_spec(unknown2, library = test_lib) |> + matches2 <- cor_spec(preproc2, library = test_lib) |> expect_warning() inherits(matches, "matrix") |> expect_true() expect_identical(dim(matches), c(ncol(test_lib$spectra), - ncol(unknown$spectra))) + ncol(preproc$spectra))) top_matches <- max_cor_named(cor_matrix = matches, na.rm = T) |> expect_silent() @@ -79,7 +101,7 @@ test_that("cor_spec returns a data.table with correct columns", { expect_silent() expect_equal(nrow(test_metadata), 1) - full_test <- ident_spec(matches, unknown, library = test_lib, top_n = 5, + full_test <- ident_spec(matches, preproc, library = test_lib, top_n = 5, add_library_metadata = "sample_name") |> expect_silent() @@ -92,24 +114,35 @@ test_that("cor_spec returns a data.table with correct columns", { ) }) +test_that("filter_spec() handles input errors correctly", { + filter_spec(1:1000) |> expect_error() +}) + # Write the tests for filter_spec function -test_that("filter_spec returns OpenSpecy object with filtered spectra", { +test_that("filter_spec() returns erroneous OpenSpecy object when removing all spectra", { os_filtered <- filter_spec(test_lib, logic = rep(F, ncol(test_lib$spectra))) |> expect_silent() + expect_warning(check_OpenSpecy(os_filtered)) expect_equal(ncol(os_filtered$spectra), 0) expect_equal(nrow(os_filtered$metadata), 0) }) # Write the tests for filter_spec function -test_that("filter_spec returns OpenSpecy object with filtered spectra", { +test_that("filter_spec() returns OpenSpecy object with filtered spectra", { logic <- rep(F,ncol(test_lib$spectra)) logic[1] <- TRUE os_filtered <- filter_spec(test_lib, logic = logic) |> expect_silent() + expect_true(check_OpenSpecy(os_filtered)) + expect_equal(ncol(os_filtered$spectra), 1) expect_equal(nrow(os_filtered$metadata), 1) }) +test_that("get_metadata() handles input errors correctly", { + get_metadata(1:1000) |> expect_error() +}) + # Tidy up unlink(tmp, recursive = T) diff --git a/tests/testthat/test-process_spec.R b/tests/testthat/test-process_spec.R index acc934be..707bdcc9 100644 --- a/tests/testthat/test-process_spec.R +++ b/tests/testthat/test-process_spec.R @@ -18,16 +18,18 @@ test_that("process_spec() returns expected values", { derivative = 1, abs = T, make_rel = T)) - + expect_true(check_OpenSpecy(conf)) + proc <- process_spec(raman_hdpe, smooth_intens = TRUE, smooth_intens_args = list( polynomial = 3, window = 11, derivative = 1 - )) |> - expect_silent() - + )) |> expect_silent() + + expect_true(check_OpenSpecy(proc)) + expect_equal(proc, conform_spec(raman_hdpe) |> smooth_intens(derivative = 1, make_rel = T)) }) diff --git a/tests/testthat/test-read_envi.R b/tests/testthat/test-read_envi.R index 9f54bf92..669a2441 100644 --- a/tests/testthat/test-read_envi.R +++ b/tests/testthat/test-read_envi.R @@ -9,6 +9,7 @@ test_that("ENVI files are read", { expect_message() |> expect_warning() expect_s3_class(tiny_map, "OpenSpecy") + expect_true(check_OpenSpecy(tiny_map)) expect_equal(ncol(tiny_map$spectra), 208) expect_length(tiny_map$wavenumber, 427) diff --git a/tests/testthat/test-read_ext.R b/tests/testthat/test-read_ext.R index 6547a3f6..67b6004d 100644 --- a/tests/testthat/test-read_ext.R +++ b/tests/testthat/test-read_ext.R @@ -32,6 +32,8 @@ test_that("read_text() gives expected output", { expect_warning() |> expect_warning() expect_s3_class(csv, "OpenSpecy") + expect_true(check_OpenSpecy(csv)) + expect_equal(names(csv), c("wavenumber", "spectra", "metadata")) expect_equal(csv$wavenumber, raman_hdpe$wavenumber) expect_equal(csv$spectra, raman_hdpe$spectra) @@ -47,6 +49,8 @@ test_that("read_asp() gives expected output", { expect_message() |> expect_warning() expect_s3_class(asp, "OpenSpecy") + expect_true(check_OpenSpecy(asp)) + expect_equal(names(asp), c("wavenumber", "spectra", "metadata")) expect_length(asp$wavenumber, 1798) range(asp$wavenumber) |> round(1) |> @@ -65,6 +69,8 @@ test_that("read_spa() gives expected output", { expect_message() |> expect_warning() expect_s3_class(spa, "OpenSpecy") + expect_true(check_OpenSpecy(spa)) + expect_equal(names(spa), c("wavenumber", "spectra", "metadata")) expect_length(spa$wavenumber, 1738) range(spa$wavenumber) |> round(1) |> @@ -74,7 +80,7 @@ test_that("read_spa() gives expected output", { }) test_that("read_jdx() gives expected output", { - suppressWarnings(jdx <- read_jdx(read_extdata("fitr_nitrocellulose.jdx"), + suppressWarnings(jdx <- read_jdx(file = read_extdata("fitr_nitrocellulose.jdx"), encoding = "latin1")) |> capture_messages() |> expect_match("JDX file inconsistency.*") @@ -82,6 +88,8 @@ test_that("read_jdx() gives expected output", { expect_error() expect_s3_class(jdx, "OpenSpecy") + expect_true(check_OpenSpecy(jdx)) + expect_equal(names(jdx), c("wavenumber", "spectra", "metadata")) expect_length(jdx$wavenumber, 7154) range(jdx$wavenumber) |> round(1) |> @@ -99,6 +107,8 @@ test_that("read_spc() gives expected output", { expect_message() |> expect_warning() expect_s3_class(spc, "OpenSpecy") + expect_true(check_OpenSpecy(spc)) + expect_equal(names(spc), c("wavenumber", "spectra", "metadata")) expect_length(spc$wavenumber, 559) range(spc$wavenumber) |> round(1) |> diff --git a/tests/testthat/test-read_opus.R b/tests/testthat/test-read_opus.R index 23dba994..d901ac27 100644 --- a/tests/testthat/test-read_opus.R +++ b/tests/testthat/test-read_opus.R @@ -5,8 +5,10 @@ dir.create(tmp, showWarnings = F) test_that("opus files are read correctly", { single <- read_extdata("ftir_ps.0") |> read_opus() |> expect_silent() + multi <- c(read_extdata("ftir_ps.0"), read_extdata("ftir_ps.0")) |> read_opus() |> expect_silent() + read_extdata("raman_hdpe.csv") |> read_opus() |> expect_error() read_extdata("ftir_ps.0") |> read_opus(share = tmp) |> expect_message() |> @@ -14,6 +16,8 @@ test_that("opus files are read correctly", { expect_s3_class(single, "OpenSpecy") expect_s3_class(multi, "OpenSpecy") + expect_true(check_OpenSpecy(single)) + expect_true(check_OpenSpecy(multi)) expect_equal(names(single), c("wavenumber", "spectra", "metadata")) expect_equal(names(multi), c("wavenumber", "spectra", "metadata")) diff --git a/tests/testthat/test-smooth_intens.R b/tests/testthat/test-smooth_intens.R index 954f2167..4f88636b 100644 --- a/tests/testthat/test-smooth_intens.R +++ b/tests/testthat/test-smooth_intens.R @@ -1,11 +1,16 @@ data("raman_hdpe") +test_that("smooth_intens() handles input errors correctly", { + smooth_intens(1:1000) |> expect_error() +}) + test_that("smooth_intens() works as expected", { smt <- smooth_intens(raman_hdpe, polynomial = 3) |> expect_silent() cor(smt$spectra$intensity, - smooth_intens(raman_hdpe, polynomial = 1)$spectra$intensity) |> round(4) |> - expect_equal(0.9756, ignore_attr = F) + smooth_intens(raman_hdpe, polynomial = 1)$spectra$intensity) |> + round(4) |> + expect_equal(0.8043, ignore_attr = F) expect_s3_class(smt, "OpenSpecy") expect_equal(nrow(smt$spectra), nrow(raman_hdpe$spectra)) expect_equal(smt$wavenumber, raman_hdpe$wavenumber) diff --git a/tests/testthat/test-subtr_baseline.R b/tests/testthat/test-subtr_baseline.R index 52c6b542..9c36d5c5 100644 --- a/tests/testthat/test-subtr_baseline.R +++ b/tests/testthat/test-subtr_baseline.R @@ -7,6 +7,8 @@ test_that("polynomial subtr_baseline() works as expected", { subtr_baseline(raman_hdpe, degree = 1)$spectra$intensity) |> round(4) |> expect_equal(0.9763, ignore_attr = F) expect_s3_class(poly, "OpenSpecy") + expect_true(check_OpenSpecy(poly)) + expect_equal(nrow(poly$spectra), nrow(raman_hdpe$spectra)) expect_equal(poly$wavenumber, raman_hdpe$wavenumber) expect_equal(range(poly$spectra), c(0, 1)) @@ -20,6 +22,7 @@ test_that("manual subtr_baseline() works as expected", { man <- subtr_baseline(raman_hdpe, type = "manual", baseline = bl) |> expect_silent() + expect_true(check_OpenSpecy(man)) cor(raman_hdpe$spectra$intensity, man$spectra$intensity) |> expect_equal(1, ignore_attr = F) diff --git a/tests/testthat/test-workflows.R b/tests/testthat/test-workflows.R index d4456292..59d8e8e2 100644 --- a/tests/testthat/test-workflows.R +++ b/tests/testthat/test-workflows.R @@ -11,6 +11,7 @@ test_that("Raman batch analysis with test library", { is_OpenSpecy(batch) |> expect_true() plot(batch) |> expect_silent() plotly_spec(batch) |> expect_silent() + expect_true(check_OpenSpecy(batch)) get_lib(type = "test", path = tmp) |> expect_no_error() check_lib(type = "test", path = tmp) |> expect_silent() @@ -19,6 +20,8 @@ test_that("Raman batch analysis with test library", { filter_spec(lib, lib$metadata$SpectrumType == "Raman") |> expect_silent() batch2 <- conform_spec(batch, lib$wavenumber, res = spec_res(lib$wavenumber)) |> expect_silent() + + expect_true(check_OpenSpecy(batch2)) plotly_spec(batch2) |> expect_silent() @@ -27,6 +30,8 @@ test_that("Raman batch analysis with test library", { batch3 <- process_spec(batch2, subtr_baseline = T) |> expect_silent() plotly_spec(x = batch3, x2 = batch) |> expect_silent() + expect_true(check_OpenSpecy(batch3)) + matches <- cor_spec(batch3, library = lib) |> expect_silent() test_max_cor <- max_cor_named(matches) |> expect_silent() sig_noise(batch3, metric = "run_sig_over_noise") |> @@ -52,6 +57,7 @@ test_that("Raman batch analysis with complete library", { batch2 <- conform_spec(batch, range = lib$wavenumber, res = spec_res(lib$wavenumber)) |> expect_silent() + expect_true(check_OpenSpecy(batch2)) plotly_spec(batch2) |> expect_silent() @@ -60,6 +66,8 @@ test_that("Raman batch analysis with complete library", { batch3 <- process_spec(batch2, subtr_baseline = T) |> expect_silent() plotly_spec(x = batch3, x2 = batch) |> expect_silent() + expect_true(check_OpenSpecy(batch3)) + matches <- cor_spec(batch3, library = lib) |> expect_silent() test_max_cor <- max_cor_named(matches) |> expect_silent() test_sn <- sig_noise(batch3, metric = "run_sig_over_noise") |> diff --git a/vignettes/advanced.Rmd b/vignettes/advanced.Rmd index 6509573f..c33b2491 100644 --- a/vignettes/advanced.Rmd +++ b/vignettes/advanced.Rmd @@ -39,7 +39,7 @@ you publish your results. The authors were kind enough to make their data open access so that it could be used in Open Specy and we should return the favor by citing them. -# Pearson's R +# Pearson's **r** Correlation values are used to identify the closest matches available in the current Open Specy spectral libraries to improve material identification and diff --git a/vignettes/app.Rmd b/vignettes/app.Rmd new file mode 100644 index 00000000..6d59da22 --- /dev/null +++ b/vignettes/app.Rmd @@ -0,0 +1,437 @@ +--- +title: "Open Specy App Tutorial" +author: > + Win Cowger, Zacharias Steinmetz +date: "`r Sys.Date()`" +output: rmarkdown::html_vignette +vignette: > + %\VignetteIndexEntry{Open Specy App Tutorial} + %\VignetteEncoding{UTF-8} + %\VignetteEngine{knitr::rmarkdown} +--- + +```{r, include = FALSE} +knitr::opts_chunk$set( + collapse = TRUE, + comment = "#>", + warning = FALSE +) +``` + +# Document Overview + +This document outlines a common workflow for using the Open Specy +app and highlights some topics that users are often requesting a tutorial on. If +the document is followed sequentially from beginning to end, the user will have +a better understanding of every procedure involved in using the Open Specy R + app as a tool for interpreting spectra. It takes approximately 45 +minutes to read through and follow along with this standard operating procedure +the first time. Afterward, knowledgeable users should be able to thoroughly +analyze spectra at an average speed of 1 min^-1^ or faster with the new batch +and automated procedures. If you are looking for documentation about the R package please see the [package vignette](http://wincowger.com/OpenSpecy-package/articles/sop.html) + +# Running the App +To get started with the Open Specy user interface, access +[https://openanalysis.org/openspecy/](https://openanalysis.org/openspecy/) +or start the Shiny GUI directly from your own computer in R. +If using the package, you just need to read in the library and run the command `run_app()`. + +```{r setup} +library(OpenSpecy) +``` + +```{r eval = F} +run_app() +``` + +# Reading Data + +Once the app is open, if you have your own data, click **Browse** at the top left hand corner of the +Analyze Spectra tab. + +```{r, fig.align="center", out.width="98%", echo=FALSE} +knitr::include_graphics("app/mainpage.jpg") +``` + +Before uploading, indicate if you would like to share the uploaded data or not +using the slider. If selected, any data uploaded to the tool will automatically +be shared under [CC-BY 4.0 +license](https://creativecommons.org/licenses/by/4.0/) and will be available for +researchers and other ventures to use to improve spectral analysis, build +machine learning tools, etc. Some users may choose not to share if they need to +keep their data private. If switched off, none of the uploaded data will be +stored or shared in Open Specy. + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/uploadfile.jpg") +``` + +Open Specy allows for upload of native Open Specy .y(a)ml, .json, or .rds files. +In addition, .csv, .asp, .jdx, .0, .spa, .spc, and .zip files can be imported. +.zip files can either contain multiple files with individual spectra in them of +the non-zip formats or it can contain a .hdr and .dat file that form an ENVI +file for a spectral map. Open Specy and .csv files should always load correctly +but the other file types are still in development, though most of the time these +files work perfectly. If uploading a .csv file, label the column with the +wavenumbers `wavenumber` and name the column with the intensities `intensity`. +Wavenumber units must be cm^-1^. Any other columns are not used by the software. +Always keep a copy of the original file before alteration to preserve metadata +and raw data for your records. + +```{r sample_data, echo=FALSE, out.width="50%"} +knitr::kable(head(raman_hdpe), caption = "Sample data `raman_hdpe`") +``` + +It is best practice to cross check files in the proprietary software they came +from and Open Specy before use in Open Specy. Due to the complexity of some +proprietary file types, we haven't been able to make them fully compatible yet. +If your file is not working, please contact the administrator and share the file +so that we can work on integrating it. + +The specific steps to converting your instrument's native files to .csv can be +found in its software manual or you can check out +[Spectragryph](https://www.effemm2.de/spectragryph/), which supports many +spectral file conversions see [Spectragryph Tutorial](http://wincowger.com/OpenSpecy-package/articles/spectragryph.html). + +If you don't have your own data, you can use a test dataset. + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/samplefile.jpg") +``` + +A .csv file of an HDPE Raman spectrum will download on your computer. This file +can also be used as a template for formatting .csv data into an Open Specy +accepted format. + +# Visualization + +## Spectra + +In the app after spectral data are uploaded, it will appear in the main window. +This plot is selectable, zoomable, and provides information on hover. You can +also save a .png file of the plot view using the camera icon at the top right +when you hover over the plot. + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/spectra_vis.jpg") +``` + + +## Maps + +Spectral maps can also be visualized as overlaid spectra but in addition the +spatial information can be plotted as a heatmap. A similar plot should appear in +the app if you upload multiple spectra or a spectral map. It is important to +note that when multiple spectra are uploaded in batch they are prescribed `x` +and `y` coordinates, this can be helpful for visualizing summary statistics and +navigating vast amounts of data but shouldn't be confused with data which +actually has spatial coordinates. Hovering over the map will reveal information +about the signal and noise and correlation values. Clicking the map will provide +the selected spectrum underneath. In the app, the colors of the heatmap are +either the signal and noise or the spectral identifications depending on whether +the identification is turned on or not. Pixels are black if the spectra does not +pass the signal-noise threshold and/or the correlation threshold. + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/map_vis.jpg") +``` + +# Processing + +The goal of this processing is to increase the signal to noise ratio (S/N) without +distorting the shape or relative size of the peaks. After uploading data, you can preprocess the data using intensity adjustment, +baseline subtraction, smoothing, flattening, and range selection and save your +preprocessed data. Once the process button is selected the default processing +will initiate. This is an absolute derivative transformation, it is kind of +magic, it does something similar to smoothing, baseline subtraction, and +intensity correction simultaneously and really quickly. By clicking +preprocessing you should see the spectrum update with the processed spectra. + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/processing.jpg") +``` +To view the raw data again, just deselect preprocessing. Toggling preprocessing +on and off can help you to make sure that the spectra are processed correctly. + + +## Threshold Signal-Noise + +Considering whether you have enough signal to analyze spectra is important. +Classical spectroscopy would recommend your highest peak to be at least 10 times +the baseline of your processed spectra before you begin analysis. If your +spectra is below that threshold, you may want to consider recollecting it. In +practice, we are rarely able to collect spectra of that good quality and more +often use 5. The Signal Over Noise technique searches your spectra for high and +low regions and conducts division on them to derive the signal to noise ratio. +Signal Times Noise multiplies the mean signal by the standard deviation of the +signal and Total Signal sums the intensities. The latter can be really useful +for thresholding spectral maps to identify particles which we will discuss +later. + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/snr_threshold.jpg") +``` + +If analyzing spectra in batch, we recommend looking at the heatmap and +optimizing the percent of spectra that are above your signal to noise threshold +to determine the correct settings instead of looking through spectra +individually. Good Signal tells you what percent of your data are above your +signal threshold. + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/signal_settings.jpg") +``` + +## Intensity Adjustment + +```{r, fig.align='center', echo=FALSE} +knitr::include_graphics("app/intensityadjustment.jpg") +``` + +Open Specy assumes that intensity units are in absorbance units but Open Specy +can adjust reflectance or transmittance spectra to absorbance units using this +selection. The transmittance adjustment uses the $\log_{10} 1/T$ calculation +which does not correct for system or particle characteristics. The reflectance +adjustment use the Kubelka-Munk equation $\frac{(1-R)^2}{2R}$. If none is +selected, Open Specy assumes that the uploaded data is an absorbance spectrum +and does not apply an adjustment. + +## Conforming + +Conforming spectra is essential before comparing to a reference library and can +be useful for summarizing data when you don't need it to be highly resolved +spectrally. In the app, conforming happens behind the scenes without any user +input to make sure that the spectra the user uploads and the library spectra +will be compatible. We set the spectral resolution to 5 because +this tends to be pretty good for a lot of applications and is in between 4 and 8 +which are commonly used wavenumber resolutions. + +## Smoothing + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/smoothingpoly.jpg") +``` + +The value on the slider is +the polynomial order of the [Savitzky-Golay (SG) +filter](https://en.wikipedia.org/wiki/Savitzky%E2%80%93Golay_filter). Higher +numbers lead to more wiggly fits and thus less smooth, lower numbers lead to +more smooth fits. The SG filter is fit to a moving window of 11 data points by +default where the center point in the window is replaced with the polynomial +estimate. Larger windows will produce smoother fits. The derivative order is set +to 1 by default which transforms the spectra to their first derivative. A zero +order derivative will have no derivative transformation. When smoothing is done +well, peak shapes and relative heights should not change. The absolute value is +primarily useful for first derivative spectra where the absolute value results +in an absorbance-like spectrum which is why we set it as the default. + +## Baseline Correction + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/baselinecorrectionpoly.jpg") +``` + +The goal of baseline correction is to get all non-peak regions of the spectra to +zero absorbance. The higher the polynomial order, the more wiggly the fit to the +baseline. If the baseline is not very wiggly, a more wiggly fit could remove +peaks which is not desired. The baseline correction algorithm used in Open Specy +is called "iModPolyfit" (Zhao et al. 2007). This algorithm iteratively fits +polynomial equations of the specified order to the whole spectrum. During the +first fit iteration, peak regions will often be above the baseline fit. The data +in the peak region is removed from the fit to make sure that the baseline is +less likely to fit to the peaks. The iterative fitting terminates once the +difference between the new and previous fit is small. An example of a good +baseline fit below. For those who have been with OpenSpecy a while, you will +notice the app no longer supports manual baseline correction, it was a hard +choice but it just didn't make sense to keep it now that we are moving toward +high throughput automated methods. It does still exist in the R function though. + +## Range Selection + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/rangeselection.jpg") +``` + +Sometimes an instrument operates with high noise at the ends of the spectrum +and, a baseline fit produces distortions, or there are regions of interest for +analysis. Range selection accomplishes those goals. You should look into the +signal to noise ratio of your specific instrument by wavelength to determine +what wavelength ranges to use. Distortions due to baseline fit can be assessed +from looking at the preprocess spectra. Additionally, you can restrict the range +to examine a single peak or a subset of peaks of interests. The maximum and +minimum wavenumbers will specify the range to keep in the app. + +## Flattening Ranges + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/flattten_example.jpg") +``` + +Sometimes there are peaks that really shouldn't be in your spectra and can +distort your interpretation of the spectra but you don't necessarily want to +remove the regions from the analysis because you believe those regions should be +flat instead of having a peak. One way to deal with this is to replace the peak +values with the mean of the values around the peak. This is the purpose of the +Flatten Range. By default it is set to flatten the CO2 region for +FTIR spectra because that region often needs to be flattened when atmospheric +artifacts occur in spectra. + +## Min-Max Normalization + +Min-Max normalization is used throughout the app but is not one of the options +the user can specify. Often we regard spectral intensities as arbitrary and +min-max normalization allows us to view spectra on the same scale without +drastically distorting their shapes or relative peak intensities. In the app, it +is primarily used for plotting and as a processing step before correlation +analysis. + +## Downloading Processed Data + +You can download a csv file of your processed data using the Your Spectra option +with download. This will download a version of your spectra with whatever +processing options you have selected at the time. + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/download_yourspec.jpg") +``` + +# Identifying Spectra + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/identification.jpg") +``` + +After uploading data and preprocessing it (if desired) you can now identify the +spectrum. To identify the spectrum go to the **Identification** box. Pro tip: if +you select **Identification** without uploading data to the app, you'll be able +to explore the library by itself. + +## Reading Libraries + +These options define the strategy for identification. The **ID Library** will +inform which library is used. Both (default) will search both FTIR and Raman +libraries. Deriv will search against a derivative transformed library. No +Baseline will search against a baseline corrected library. This should be in +line with how you choose to process your spectra. Cor options use a simple +Pearson correlation search algorithm. AI is currently experimental and uses +either a multinomial model or correlation on mediod spectra from the library. +Correlation thresholding will set the minimum value from matching to use as a +'positive identification' this will black out pixels in a spectral map view if +they are lower than the threshold. + +## Matches + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/matches.jpg") +``` + +Top matches in the app can be assessed by clicking the cog in the right hand +corner of the Spectra box. This will open a side window with the matches sorted +from most to least similar. Clicking on rows in the table will add the selected +match to the spectra viewer. Using the table's filter options, you can restrict +the range of Pearson\'s r values or search for specific material types. + +## Selection Metadata + +```{r, fig.align="center", out.width="98%", echo=FALSE} +knitr::include_graphics("app/selectionmetadata.jpg") +``` + +Whatever match is selected from the match table may have additional metadata +about it. That metadata will be displayed below the plot. Some of this metadata +may assist you in interpreting the spectra. For example, if the spectra has +metadata which says it is a liquid and you are analyzing a solid particle, that +spectrum may not be the best match. + +## Match Plot + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/matchplot.jpg") +``` + +This plot is dynamically updated by selecting matches from the match table. The +red spectrum is the spectrum that you selected from the reference library and +the white spectrum is the spectrum that you are trying to identify. Whenever a +new dataset is uploaded, the plot and data table in this tab will be updated. +These plots can be saved as a .png by clicking the camera button at the top of +the plot. + +## Sharing Reference Data + +If you have reference data or AI models that you think would be useful for +sharing with the spectroscopy community through OpenSpecy please contact the +website administrator to discuss options for collaborating. + +## Download Top Matches + +To download the top matches and associated metadata, use the top matches download option. +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/download.jpg") +``` + +# Characterizing Particles + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/download_threshold.jpg") +``` + +Thresholding only works and shouldn't be used for single particles. You can +threshold spectral maps using the signal-noise values and/or the correlation +threshold. Whatever selection you have active will inform how your spectra is +thresholded. It should get thresholded in coincidence with what you see on the +spectral map. Pro tip: The file that gets downloaded is an OpenSpecy with median +spectra for each of your thresholded particles. You can upload that back to the +app and analyze that independently. + +## Additional App Features + +Accessibility is extremely important to us and we are making strives to improve +the accessibility of Open Specy for all spectroscopists. Please reach out if you +have ideas for improvement. + +We added a Google translate plugin to all pages in the app so that you can +easily translate the app. We know that not all languages will be fully supported +but we will continue to try and improve the translations. + +```{r, fig.align="center", echo=FALSE} +knitr::include_graphics("app/translate.jpg") +``` + +# References + +Chabuka BK, Kalivas JH (2020). “Application of a Hybrid Fusion Classification +Process for Identification of Microplastics Based on Fourier Transform Infrared +Spectroscopy.†*Applied Spectroscopy*, **74**(9), 1167–1183. doi: +[10.1177/0003702820923993](https://doi.org/10.1177/0003702820923993). + +Cowger W, Gray A, Christiansen SH, Christiansen SH, Christiansen SH, De Frond H, +Deshpande AD, Hemabessiere L, Lee E, Mill L, et al. (2020). “Critical Review of +Processing and Classification Techniques for Images and Spectra in Microplastic +Research.†*Applied Spectroscopy*, **74**(9), 989–1010. doi: +[10.1177/0003702820929064](https://doi.org/10.1177/0003702820929064). + +Cowger W, Steinmetz Z, Gray A, Munno K, Lynch J, Hapich H, Primpke S, +De Frond H, Rochman C, Herodotou O (2021). “Microplastic Spectral Classification +Needs an Open Source Community: Open Specy to the Rescue!†+*Analytical Chemistry*, **93**(21), 7543–7548. doi: +[10.1021/acs.analchem.1c00123](https://doi.org/10.1021/acs.analchem.1c00123). + +Primpke S, Wirth M, Lorenz C, Gerdts G (2018). “Reference Database Design for +the Automated Analysis of Microplastic Samples Based on Fourier Transform +Infrared (FTIR) Spectroscopy.†*Analytical and Bioanalytical Chemistry*, +**410**(21), 5131–5141. doi: +[10.1007/s00216-018-1156-x](https://doi.org/10.1007/s00216-018-1156-x). + +Renner G, Schmidt TC, Schram J (2017). “A New Chemometric Approach for Automatic +Identification of Microplastics from Environmental Compartments Based on FT-IR +Spectroscopy.†*Analytical Chemistry*, **89**(22), 12045–12053. doi: +[10.1021/acs.analchem.7b02472](https://doi.org/10.1021/acs.analchem.7b02472). + +Savitzky A, Golay MJ (1964). “Smoothing and Differentiation of Data by +Simplified Least Squares Procedures.†*Analytical Chemistry*, **36**(8), +1627–1639. + +Zhao J, Lui H, McLean DI, Zeng H (2007). “Automated Autofluorescence Background +Subtraction Algorithm for Biomedical Raman Spectroscopy.†+*Applied Spectroscopy*, **61**(11), 1225–1232. doi: +[10.1366/000370207782597003](https://doi.org/10.1366/000370207782597003). diff --git a/vignettes/sop/baselinecorrectionpoly.jpg b/vignettes/app/baselinecorrectionpoly.jpg similarity index 100% rename from vignettes/sop/baselinecorrectionpoly.jpg rename to vignettes/app/baselinecorrectionpoly.jpg diff --git a/vignettes/app/download.jpg b/vignettes/app/download.jpg new file mode 100644 index 00000000..54c96ddf Binary files /dev/null and b/vignettes/app/download.jpg differ diff --git a/vignettes/app/download_threshold.jpg b/vignettes/app/download_threshold.jpg new file mode 100644 index 00000000..b22ee18e Binary files /dev/null and b/vignettes/app/download_threshold.jpg differ diff --git a/vignettes/app/download_yourspec.jpg b/vignettes/app/download_yourspec.jpg new file mode 100644 index 00000000..a86364c9 Binary files /dev/null and b/vignettes/app/download_yourspec.jpg differ diff --git a/vignettes/sop/flattten_example.jpg b/vignettes/app/flattten_example.jpg similarity index 100% rename from vignettes/sop/flattten_example.jpg rename to vignettes/app/flattten_example.jpg diff --git a/vignettes/sop/identification.jpg b/vignettes/app/identification.jpg similarity index 100% rename from vignettes/sop/identification.jpg rename to vignettes/app/identification.jpg diff --git a/vignettes/sop/intensityadjustment.jpg b/vignettes/app/intensityadjustment.jpg similarity index 100% rename from vignettes/sop/intensityadjustment.jpg rename to vignettes/app/intensityadjustment.jpg diff --git a/vignettes/sop/mainpage.jpg b/vignettes/app/mainpage.jpg similarity index 100% rename from vignettes/sop/mainpage.jpg rename to vignettes/app/mainpage.jpg diff --git a/vignettes/sop/map_vis.jpg b/vignettes/app/map_vis.jpg similarity index 100% rename from vignettes/sop/map_vis.jpg rename to vignettes/app/map_vis.jpg diff --git a/vignettes/sop/matches.jpg b/vignettes/app/matches.jpg similarity index 100% rename from vignettes/sop/matches.jpg rename to vignettes/app/matches.jpg diff --git a/vignettes/app/matchplot.jpg b/vignettes/app/matchplot.jpg new file mode 100644 index 00000000..64a8c832 Binary files /dev/null and b/vignettes/app/matchplot.jpg differ diff --git a/vignettes/sop/processing.jpg b/vignettes/app/processing.jpg similarity index 100% rename from vignettes/sop/processing.jpg rename to vignettes/app/processing.jpg diff --git a/vignettes/sop/rangeselection.jpg b/vignettes/app/rangeselection.jpg similarity index 100% rename from vignettes/sop/rangeselection.jpg rename to vignettes/app/rangeselection.jpg diff --git a/vignettes/sop/samplefile.jpg b/vignettes/app/samplefile.jpg similarity index 100% rename from vignettes/sop/samplefile.jpg rename to vignettes/app/samplefile.jpg diff --git a/vignettes/sop/selectionmetadata.jpg b/vignettes/app/selectionmetadata.jpg similarity index 100% rename from vignettes/sop/selectionmetadata.jpg rename to vignettes/app/selectionmetadata.jpg diff --git a/vignettes/sop/sig_over_noise.jpg b/vignettes/app/sig_over_noise.jpg similarity index 100% rename from vignettes/sop/sig_over_noise.jpg rename to vignettes/app/sig_over_noise.jpg diff --git a/vignettes/sop/signal_settings.jpg b/vignettes/app/signal_settings.jpg similarity index 100% rename from vignettes/sop/signal_settings.jpg rename to vignettes/app/signal_settings.jpg diff --git a/vignettes/sop/smoothingpoly.jpg b/vignettes/app/smoothingpoly.jpg similarity index 100% rename from vignettes/sop/smoothingpoly.jpg rename to vignettes/app/smoothingpoly.jpg diff --git a/vignettes/sop/snr_threshold.jpg b/vignettes/app/snr_threshold.jpg similarity index 100% rename from vignettes/sop/snr_threshold.jpg rename to vignettes/app/snr_threshold.jpg diff --git a/vignettes/sop/spectra_vis.jpg b/vignettes/app/spectra_vis.jpg similarity index 100% rename from vignettes/sop/spectra_vis.jpg rename to vignettes/app/spectra_vis.jpg diff --git a/vignettes/app/translate.jpg b/vignettes/app/translate.jpg new file mode 100644 index 00000000..26a79dde Binary files /dev/null and b/vignettes/app/translate.jpg differ diff --git a/vignettes/sop/uploadfile.jpg b/vignettes/app/uploadfile.jpg similarity index 100% rename from vignettes/sop/uploadfile.jpg rename to vignettes/app/uploadfile.jpg diff --git a/vignettes/sop.Rmd b/vignettes/sop.Rmd index 04d4c4da..55136b0a 100644 --- a/vignettes/sop.Rmd +++ b/vignettes/sop.Rmd @@ -1,13 +1,11 @@ --- -title: "Open Specy Tutorial" +title: "Open Specy Package Tutorial" author: > - Jessica Meyers, Jeremy Conkle, Win Cowger, Zacharias Steinmetz, - Andrew Gray, Chelsea Rochman, Sebastian Primpke, Jennifer Lynch, - Hannah Hapich, Hannah De Frond, Keenan Munno, Bridget O’Donnell + Win Cowger, Zacharias Steinmetz date: "`r Sys.Date()`" output: rmarkdown::html_vignette vignette: > - %\VignetteIndexEntry{Open Specy Tutorial} + %\VignetteIndexEntry{Open Specy Package Tutorial} %\VignetteEncoding{UTF-8} %\VignetteEngine{knitr::rmarkdown} --- @@ -23,22 +21,22 @@ knitr::opts_chunk$set( # Document Overview This document outlines a common workflow for using the Open Specy package and -app and highlights some topics that users are often requesting a tutorial on. If -the document is followed sequentially from beginning to end, the user will have -a better understanding of every procedure involved in using the Open Specy R -package and app as a tool for interpreting spectra. It takes approximately 45 -minutes to read through and follow along with this standard operating procedure -the first time. Afterward, knowledgeable users should be able to thoroughly -analyze spectra at an average speed of 1 min^-1^ or faster with the new batch -and automated procedures. - -# R Package and App +highlights some topics that users are often requesting a tutorial on. If the +document is followed sequentially from beginning to end, the user will have a +better understanding of every procedure involved in using the Open Specy R +package as a tool for interpreting spectra. It takes approximately 45 minutes to +read through and follow along with this standard operating procedure the first +time. Afterward, knowledgeable users should be able to thoroughly analyze +spectra at an average speed of 1 min^-1^ or faster with the new batch and +automated procedures. The Open Specy R package is the backbone of the Shiny app. The choice is yours as to which you start with, we use both on a regular basis. The tutorial will -talk through the R functions and app features step by step together. +talk through the R functions you can use to programatically analyze spectra. If +you are looking for a tutorial about how to use the app see [this +tutorial](http://wincowger.com/OpenSpecy-package/articles/app.html). -## Installation +# Installation After installing the R package, you just need to read in the library. @@ -46,36 +44,34 @@ After installing the R package, you just need to read in the library. library(OpenSpecy) ``` -## Running the App +# Running the App To get started with the Open Specy user interface, access -[https://openanalysis.org/openspecy/](https://openanalysis.org/openspecy/) -or start the Shiny GUI directly from your own computer in R. +[https://openanalysis.org/openspecy/](https://openanalysis.org/openspecy/) or +start the Shiny GUI directly from your own computer in R. If you are looking for +a tutorial about how to use the app see [this +tutorial](http://wincowger.com/OpenSpecy-package/articles/app.html). ```{r eval=FALSE} run_app() ``` -## Reading Data +# Read Data -If you have your own data, click **Browse** at the top left hand corner of the -Analyze Spectra tab. +The following line of code will read in your data when using the package and +interprets which reading function to use based on the file extension. -```{r, fig.align="center", out.width="98%", echo=FALSE} -knitr::include_graphics("sop/mainpage.jpg") +```{r, eval=FALSE} +read_any("path/to/your/data") ``` -Before uploading, indicate if you would like to share the uploaded data or not -using the slider. If selected, any data uploaded to the tool will automatically -be shared under [CC-BY 4.0 -license](https://creativecommons.org/licenses/by/4.0/) and will be available for -researchers and other ventures to use to improve spectral analysis, build -machine learning tools, etc. Some users may choose not to share if they need to -keep their data private. If switched off, none of the uploaded data will be -stored or shared in Open Specy. +These file type specific functions will also read in spectral +data accordingly if you have a particular format in mind. -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/uploadfile.jpg") +```{r, eval=FALSE} +read_text(".csv") +read_asp(".asp") +read_opus(".0") ``` Open Specy allows for upload of native Open Specy .y(a)ml, .json, or .rds files. @@ -90,10 +86,6 @@ Wavenumber units must be cm^-1^. Any other columns are not used by the software. Always keep a copy of the original file before alteration to preserve metadata and raw data for your records. -```{r sample_data, echo=FALSE, out.width="50%"} -knitr::kable(head(raman_hdpe), caption = "Sample data `raman_hdpe`") -``` - It is best practice to cross check files in the proprietary software they came from and Open Specy before use in Open Specy. Due to the complexity of some proprietary file types, we haven't been able to make them fully compatible yet. @@ -103,37 +95,11 @@ so that we can work on integrating it. The specific steps to converting your instrument's native files to .csv can be found in its software manual or you can check out [Spectragryph](https://www.effemm2.de/spectragryph/), which supports many -spectral file conversions (see Mini Tutorial section: File conversion in -Spectragryph to Open Specy accepted format). - -The following line of code will read in your data when using the package and -interprets which reading function to use based on the file extension. - -```{r, eval=FALSE} -read_any("path/to/your/data") -``` - -These file type specific functions will also read in spectral -data accordingly if you have a particular format in mind. - -```{r, eval=FALSE} -read_text(".csv") -read_asp(".asp") -read_opus(".0") -``` +spectral file conversions. For instructions, see [Spectragryph +Tutorial](http://wincowger.com/OpenSpecy-package/articles/spectragryph.html). If you don't have your own data, you can use a test dataset. -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/samplefile.jpg") -``` - -A .csv file of an HDPE Raman spectrum will download on your computer. This file -can also be used as a template for formatting .csv data into an Open Specy -accepted format. - -The following line of code does the same to load the example file into R: - ```{r} data("raman_hdpe") ``` @@ -164,7 +130,7 @@ batch_example <- read_extdata("testdata_zipped.zip") |> You will notice now that the R package reads in files into an object with class `OpenSpecy`. This is a class we created for high throughput spectral analysis which now also preserves spectral metadata. You can even create these from -scratch if you'd like: +scratch if you'd like. ```{r} scratch_OpenSpecy <- as_OpenSpecy(x = seq(1000,2000, by = 5), @@ -172,6 +138,36 @@ scratch_OpenSpecy <- as_OpenSpecy(x = seq(1000,2000, by = 5), metadata = list(file_name = "fake_noise")) ``` +Open Specy objects are lists with three components, `wavenumber` is a vector of +the wavenumber values for the spectra and corresponds to the rows in `spectra` +which is a `data.table` where each column is a set of spectral intensities. +`metadata` is a data.table which holds additional information about the spectra. +Each row corresponds to a column in `spectra`. + +```{r} +# Access the wavenumbers +scratch_OpenSpecy$wavenumber +``` + +```{r} +# Access the spectra +scratch_OpenSpecy$spectra +``` + +```{r} +# Access the metadata +scratch_OpenSpecy$metadata +``` + +```{r} +# Performs checks to ensure that OpenSpecy objects are adhering to our standards; +# returns TRUE if it passes. +check_OpenSpecy(scratch_OpenSpecy) + +# Checks only the object type to make sure it has OpenSpecy type +is_OpenSpecy(scratch_OpenSpecy) +``` + We have some generic functions built for inspecting the spectra: ```{r} @@ -186,19 +182,46 @@ summary(scratch_OpenSpecy) # summarizes the contents of the spectra head(scratch_OpenSpecy) # shows the top wavenumbers and intensities ``` -## Visualization +# Save Data + +Open Specy objects can be saved most accurately as .rds, .yml, or .json files. +.rds will be the most reproducible as it is a native R file format and floating +point errors can happen with .json or .yml. -### Spectra +```{r, eval=F} +write_spec(scratch_OpenSpecy, "test_scratch_OpenSpecy.rds") +write_spec(scratch_OpenSpecy, "test_scratch_OpenSpecy.yml", digits = 5) +write_spec(scratch_OpenSpecy, "test_scratch_OpenSpecy.json", digits = 5) +``` -In the app after spectral data are uploaded, it will appear in the main window. -This plot is selectable, zoomable, and provides information on hover. You can -also save a .png file of the plot view using the camera icon at the top right -when you hover over the plot. +# Format Conversions -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/spectra_vis.jpg") +Another great spectroscopy R package is +[hyperSpec](https://CRAN.R-project.org/package=hyperSpec). We +actually depend on their functions for several of ours. They are currently +making some awesome new features and we want to integrate well with them so that +both packages can be easily used together. That is why we created the +`as_hyperSpec` function and made sure that `as_OpenSpecy` can convert from +`hyperSpec` objects. + +```{r eval=F} +hyperspecy <- as_hyperSpec(scratch_OpenSpecy) + +identical(as_OpenSpecy(hyperspecy)$wavenumber, scratch_OpenSpecy$wavenumber) + +# There are some challenges with reproducing the names and floating point +# numbers but the underlying data is preserved. +identical(unname(unlist(signif(as_OpenSpecy(hyperspecy)$spectra, 2))), + unname(unlist(signif(scratch_OpenSpecy$spectra, 2)))) + +# Metadata is not perfectly preserved yet. +identical(as_OpenSpecy(hyperspecy)$metadata, scratch_OpenSpecy$metadata) ``` +# Visualization + +## Spectra + In R, we have two ways to visualize your spectra, one is quick and efficient and the other is interactive. Here is an example of quick and efficient plotting. @@ -209,39 +232,30 @@ plot(scratch_OpenSpecy) # quick and efficient This is an example of an interactive plot. ```{r, fig.align="center", out.width="98%"} -plotly_spec(scratch_OpenSpecy) # this will min-max normalize your data even - # if it isn't already but are not changing your - # underlying data +# This will min-max normalize your data even if it isn't already but are not +# changing your underlying data +plotly_spec(scratch_OpenSpecy) ``` With interactive plots you can also plot two different datasets simultaneously to compare. ```{r, fig.align="center", out.width="98%"} -plotly_spec(scratch_OpenSpecy, rds_example) # interactive comparison +plotly_spec(scratch_OpenSpecy, rds_example) ``` -### Maps +## Maps Spectral maps can also be visualized as overlaid spectra but in addition the -spatial information can be plotted as a heatmap. A similar plot should appear in -the app if you upload multiple spectra or a spectral map. It is important to +spatial information can be plotted as a heatmap. It is important to note that when multiple spectra are uploaded in batch they are prescribed `x` and `y` coordinates, this can be helpful for visualizing summary statistics and navigating vast amounts of data but shouldn't be confused with data which -actually has spatial coordinates. Hovering over the map will reveal information -about the signal and noise and correlation values. Clicking the map will provide -the selected spectrum underneath. In the app, the colors of the heatmap are -either the signal and noise or the spectral identifications depending on whether -the identification is turned on or not. Pixels are black if the spectra does not -pass the signal-noise threshold and/or the correlation threshold. - -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/map_vis.jpg") -``` +actually has spatial coordinates. -The same plot can be created in R but the user can control what values form the -colors of the heatmap by specifying `z`. +In R the user can control what values form the colors of the heatmap by +specifying `z`, it is handy to put the information you want in the metadata for +this reason. This example just shows the x values of the spectra. ```{r, fig.align="center", out.width="98%"} heatmap_spec(spectral_map, @@ -258,107 +272,216 @@ that spectrum. interactive_plot(spectral_map, select = 100, z = spectral_map$metadata$x) ``` -## Processing -After uploading data, you can preprocess the data using intensity adjustment, -baseline subtraction, smoothing, flattening, and range selection and save your -preprocessed data. Once the process button is selected the default processing -will initiate. This is an absolute derivative transformation, it is kind of -magic, it does something similar to smoothing, baseline subtraction, and -intensity correction simultaneously and really quickly. By clicking -preprocessing you should see the spectrum update with the processed spectra. +# Combining OpenSpecy Objects + +Sometimes you have several OpenSpecy objects that you want to combine into one. +The easiest way to do that is by having spectra which all are in the same exact +format with the same series of wavenumbers. The default settings of `c_spec` +assume that is the case. + +```{r} +#Same range example. +data("raman_hdpe") +raman_hdpe1 <- raman_hdpe +raman_hdpe2 <- raman_hdpe + +raman_hdpe1$spectra <- raman_hdpe1$spectra - 100 + +raman_hdpe2$spectra <- raman_hdpe2$spectra - 200 + +combined_spectra <- c_spec(list(raman_hdpe, raman_hdpe1, raman_hdpe2)) + +plot(combined_spectra) -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/processing.jpg") ``` -To view the raw data again, just deselect preprocessing. Toggling preprocessing -on and off can help you to make sure that the spectra are processed correctly. +If you have different wavenumber ranges for the spectra you want to combine, +you can set `range = "common"` and `res` equal to the wavenumber resolution you +want and the function will collapse all the spectra to whatever their common +range is using linear interpolation. + +```{r} +asp_example <- read_extdata("ftir_ldpe_soil.asp") |> + read_any() +ps_example <- read_extdata("ftir_ps.0") |> + read_any() # preserves some metadata + +combined <- c_spec(list(asp_example, ps_example), range = "common", res = 5) + +plot(combined) + +``` + + +# Filtering OpenSpecy Objects + +OpenSpecy objects can have any number of spectra in them. To create an OpenSpecy +with a subset of the spectra that is in an Open Specy object you can use the +`filter_spec` function. Filtering is allowed by index number, name, or using a +logical vector. Filtering will update the `spectra` and `metadata` items of the +`OpenSpecy` but not the `wavenumber`. + +```{r} +spectral_map <- read_extdata("CA_tiny_map.zip") |> + read_any() # preserves some metadata + +# Extract the 150th spectrum. +filter_spec(spectral_map, 150) |> + plot() +``` + +```{r} +# Extract the spectrum with column name "8_5". +filter_spec(spectral_map, "8_5") |> + plot() +``` + +```{r} +# Extract the spectra with column names "8_5", "10_3", or "12_10" using a +# logical argument based on the spectra column names. +filter_spec(spectral_map, names(spectral_map$spectra) %in% c("8_5", "10_3", "12_10")) |> + plot() +``` + +```{r} +filter_spec(spectral_map, spectral_map$metadata$y == 1) |> + plot() #extracts the spectra with y values of 1 for their metadata. +``` + +# Sampling OpenSpecy Objects + +Sometimes you have a large suite of examples of spectra of the same material and +you want to reduce the number of spectra you run through the analysis for +simplicity or you are running simulations or other proceedures that require you +to first sample from the spectra contained in your OpenSpecy objects before +doing analysis. the `sample_spec` function serves this purpose. + +```{r} +spectral_map <- read_any(read_extdata("CA_tiny_map.zip")) +sample_spec(spectral_map,size = 10) |> + plot() +``` + +# Processing + +The goal of this all processing is to increase the signal to noise ratio (S/N) +of the spectra without distorting the shape, position, or relative size of the +peaks. After loading data, you can preprocess the data using intensity +adjustment, baseline subtraction, smoothing, flattening, and range selection. +The default settings is an absolute derivative transformation, it is kind of +magic, it does something similar to smoothing, baseline subtraction, and +intensity correction simultaneously and really quickly. The `process_spec()` function is a monolithic function for all processing procedures which is optimized by default to result in high signal to noise in most cases, same as the app. ```{r} -processed <- process_spec(raman_hdpe) +processed <- process_spec(raman_hdpe, + active = TRUE, + adj_intens = FALSE, + adj_intens_args = list(type = "none"), + conform_spec = TRUE, + conform_spec_args = list(range = NULL, res = 5, + type = "interp"), + restrict_range = FALSE, + restrict_range_args = list(min = 0, max = 6000), + flatten_range = FALSE, + flatten_range_args = list(min = 2200, max = 2420), + subtr_baseline = FALSE, + subtr_baseline_args = list(type = "polynomial", + degree = 8, raw = FALSE, + baseline = NULL), + smooth_intens = TRUE, + smooth_intens_args = list(polynomial = 3, window = 11, + derivative = 1, abs = TRUE), + make_rel = TRUE) summary(processed) summary(raman_hdpe) ``` -You can compare the processed and unprocessed data in an overlay plot: +You can compare the processed and unprocessed data in an overlay plot. ```{r eval=FALSE} plotly_spec(raman_hdpe, processed) ``` We want people to use the `process_spec()` function for most processing -operations, all processing functions can be tuned using its parameters. However, -we recognize that nesting of functions and order of operations can be useful for -users to control so you can also use individual functions for each operation if -you'd like. +operations. All other processing functions can be tuned using its parameters in +the single function and we have set defaults and nested the processing functions +in a way that provides typically quality results. However, we recognize that +nesting of functions and order of operations can be useful for users to control +so you can also use individual functions for each operation if you'd like. See +explanations of each processing sub-function below. -#### Threshold Signal-Noise +## Threshold Signal-Noise Considering whether you have enough signal to analyze spectra is important. Classical spectroscopy would recommend your highest peak to be at least 10 times the baseline of your processed spectra before you begin analysis. If your spectra is below that threshold, you may want to consider recollecting it. In practice, we are rarely able to collect spectra of that good quality and more -often use 5. The Signal Over Noise technique searches your spectra for high and -low regions and conducts division on them to derive the signal to noise ratio. -Signal Times Noise multiplies the mean signal by the standard deviation of the -signal and Total Signal sums the intensities. The latter can be really useful -for thresholding spectral maps to identify particles which we will discuss -later. +often use 5. The Run Signal Over Noise technique searches your spectra for high +and low regions and conducts division on them to derive the signal to noise +ratio. This is a good way to automatically calculate the signal to noise ratio. +In the example below you can see that our signal to noise ratio is increased by +the processing, the goal of processing is generally to maximize this. Signal +Times Noise multiplies the mean signal by the standard deviation of the signal +and Total Signal sums the intensities. The latter can be really useful for +thresholding spectral maps to identify particles which we will discuss later. -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/snr_threshold.jpg") +```{r} +sig_noise(processed, metric = "run_sig_over_noise") > + sig_noise(raman_hdpe, metric = "run_sig_over_noise") ``` If analyzing spectra in batch, we recommend looking at the heatmap and optimizing the percent of spectra that are above your signal to noise threshold to determine the correct settings instead of looking through spectra -individually. Good Signal tells you what percent of your data are above your -signal threshold. +individually. Setting the `min_sn` will threshold the heatmap image to only +color spectra which have a `sn` value over the threshold. -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/signal_settings.jpg") -``` +```{r} +spectral_map_p <- spectral_map |> + process_spec(flatten_range = T) + +spectral_map_p$metadata$sig_noise <- sig_noise(spectral_map_p) -#### Intensity Adjustment +heatmap_spec(spectral_map_p, sn = spectral_map_p$metadata$sig_noise, min_sn = 5) -```{r, fig.align='center', echo=FALSE} -knitr::include_graphics("sop/intensityadjustment.jpg") ``` +## Intensity Adjustment + Open Specy assumes that intensity units are in absorbance units but Open Specy -can adjust reflectance or transmittance spectra to absorbance units using this -selection. The transmittance adjustment uses the $\log_{10} 1/T$ calculation -which does not correct for system or particle characteristics. The reflectance -adjustment use the Kubelka-Munk equation $\frac{(1-R)^2}{2R}$. If none is -selected, Open Specy assumes that the uploaded data is an absorbance spectrum -and does not apply an adjustment. +can adjust reflectance or transmittance spectra to absorbance units. The +transmittance adjustment uses the $\log_{10} 1/T$ calculation which does not +correct for system or particle characteristics. The reflectance adjustment use +the Kubelka-Munk equation $\frac{(1-R)^2}{2R}$. This is the respective R code for a scenario where the spectra doesn't need intensity adjustment: ```{r, fig.align="center", fig.width=6} -raman_adj <- raman_hdpe |> - adj_intens(type = "none") +data("raman_hdpe") +trans_raman_hdpe <- raman_hdpe +trans_raman_hdpe$spectra <- 2 - trans_raman_hdpe$spectra^2 + +rev_trans_raman_hdpe <- trans_raman_hdpe |> + adj_intens(type = "transmittance") -plot(raman_adj) +plotly_spec(trans_raman_hdpe, rev_trans_raman_hdpe) ``` -#### Conforming +## Conforming Conforming spectra is essential before comparing to a reference library and can be useful for summarizing data when you don't need it to be highly resolved -spectrally. In the app, conforming happens behind the scenes without any user -input to make sure that the spectra the user uploads and the library spectra -will be compatible. In code, we set the default spectral resolution to 5 because -this tends to be pretty good for a lot of applications and is in between 4 and 8 -which are commonly used wavenumber resolutions. There are several ways that this -is function can be specified. +spectrally. We set the default spectral resolution to 5 because this tends to be +pretty good for a lot of applications and is in between 4 and 8 which are +commonly used wavenumber resolutions. There are several ways that this function +can be specified. ```{r} conform_spec(raman_hdpe) |> # default convert res to 5 wavenumbers. @@ -375,64 +498,51 @@ conform_spec(spectral_map, range = ps_example$wavenumber, res = 10, summary() ``` -#### Smoothing - -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/smoothingpoly.jpg") -``` +## Smoothing -The first step of the Open Specy preprocessing routing is spectral smoothing. -The goal of this function is to increase the signal to noise ratio (S/N) without -distorting the shape or relative size of the peaks. The value on the slider is -the polynomial order of the [Savitzky-Golay (SG) -filter](https://en.wikipedia.org/wiki/Savitzky%E2%80%93Golay_filter). Higher -numbers lead to more wiggly fits and thus less smooth, lower numbers lead to -more smooth fits. The SG filter is fit to a moving window of 11 data points by -default where the center point in the window is replaced with the polynomial -estimate. Larger windows will produce smoother fits. The derivative order is set -to 1 by default which transforms the spectra to their first derivative. A zero -order derivative will have no derivative transformation. When smoothing is done -well, peak shapes and relative heights should not change. The absolute value is -primarily useful for first derivative spectra where the absolute value results -in an absorbance-like spectrum which is why we set it as the default. +The [Savitzky-Golay +filter](https://en.wikipedia.org/wiki/Savitzky%E2%80%93Golay_filter) is used for +smoothing. Higher polynomial numbers lead to more wiggly fits and thus less +smoothing, lower numbers lead to more smooth fits. The SG filter is fit to a +moving window of 11 data points by default where the center point in the window +is replaced with the polynomial estimate. Larger windows will produce smoother +fits. The derivative order is set to 1 by default which transforms the spectra +to their first derivative. A zero order derivative will have no derivative +transformation. When smoothing is done well, peak shapes and relative heights +should not change. The absolute value is primarily useful for first derivative +spectra where the absolute value results in an absorbance-like spectrum which is +why we set it as the default. You'll notice a new function we are using +`c_spec()` which is used to combine spectral objects into one OpenSpecy object. -Examples of smoothing with the R package: +Examples of smoothing: -```{r smooth_intens, fig.cap = "Sample `raman_hdpe` spectrum with different smoothing polynomials from Cowger et al. (2020).", fig.width=6, fig.align="center"} -library(ggplot2) +```{r smooth_intens, fig.cap = "Sample `raman_hdpe` spectrum with different smoothing polynomials.", fig.width=6, fig.align="center"} data("raman_hdpe") -none <- adj_intens(raman_hdpe, make_rel = T) -p1 <- smooth_intens(raman_hdpe, polynomial = 1) -p4 <- smooth_intens(raman_hdpe, polynomial = 4) +none <- make_rel(raman_hdpe) +p1 <- smooth_intens(raman_hdpe, polynomial = 1, derivative = 0, abs = F) +p4 <- smooth_intens(raman_hdpe, polynomial = 4, derivative = 0, abs = F) compare_derivative <- c_spec(list(none, p1, p4)) plot(compare_derivative) ``` -Derivative transformation can happen with the same function. You'll notice a new -function we are using `c_spec()` which is used to combine spectral objects into -one OpenSpecy object. +Derivative transformation can happen with the same function. -```{r compare_derivative, fig.cap = "Sample `raman_hdpe` spectrum with different derivatives from Cowger et al. (2020).", fig.width=6, fig.align="center"} +```{r compare_derivative, fig.cap = "Sample `raman_hdpe` spectrum with different derivatives.", fig.width=6, fig.align="center"} -none <- adj_intens(raman_hdpe, make_rel = T) +none <- make_rel(raman_hdpe) d1 <- smooth_intens(raman_hdpe, derivative = 1, abs = T) -d2 <- smooth_intens(raman_hdpe, derivative = 2) +d2 <- smooth_intens(raman_hdpe, derivative = 2, abs = T) compare_derivative <- c_spec(list(none, d1, d2)) plot(compare_derivative) ``` -#### Baseline Correction - -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/baselinecorrectionpoly.jpg") -``` +## Baseline Correction -The second step of Open Specy's preprocessing routine is baseline correction. The goal of baseline correction is to get all non-peak regions of the spectra to zero absorbance. The higher the polynomial order, the more wiggly the fit to the baseline. If the baseline is not very wiggly, a more wiggly fit could remove @@ -443,16 +553,19 @@ first fit iteration, peak regions will often be above the baseline fit. The data in the peak region is removed from the fit to make sure that the baseline is less likely to fit to the peaks. The iterative fitting terminates once the difference between the new and previous fit is small. An example of a good -baseline fit below. For those who have been with OpenSpecy a while, you will -notice the app no longer supports manual baseline correction, it was a hard -choice but it just didn't make sense to keep it now that we are moving toward -high throughput automated methods. It does still exist in the R function though. +baseline fit below. Manual baseline correction can also be specified by +providing a baseline `OpenSpecy` object. ```{r subtr_baseline, fig.cap = "Sample `raman_hdpe` spectrum with different degrees of background subtraction (Cowger et al., 2020).", fig.width=6, fig.align="center"} data("raman_hdpe") -none <- adj_intens(raman_hdpe) -d2 <- subtr_baseline(raman_hdpe, degree = 2) +alternative_baseline <- smooth_intens(raman_hdpe, polynomial = 1, window = 51, + derivative = 0, abs = F, make_rel = F) |> + flatten_range(min = 2700, max = 3200, make_rel = F) + +none <- make_rel(raman_hdpe) +d2 <- subtr_baseline(raman_hdpe, type = "manual", + baseline = alternative_baseline) d8 <- subtr_baseline(raman_hdpe, degree = 8) compare_subtraction <- c_spec(list(none, d2, d8)) @@ -460,11 +573,7 @@ compare_subtraction <- c_spec(list(none, d2, d8)) plot(compare_subtraction) ``` -#### Range Selection - -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/rangeselection.jpg") -``` +## Range Selection Sometimes an instrument operates with high noise at the ends of the spectrum and, a baseline fit produces distortions, or there are regions of interest for @@ -472,14 +581,13 @@ analysis. Range selection accomplishes those goals. You should look into the signal to noise ratio of your specific instrument by wavelength to determine what wavelength ranges to use. Distortions due to baseline fit can be assessed from looking at the preprocess spectra. Additionally, you can restrict the range -to examine a single peak or a subset of peaks of interests. The maximum and -minimum wavenumbers will specify the range to keep in the app. In the R -function, multiple ranges can be specified simultaneously. +to examine a single peak or a subset of peaks of interests. Multiple ranges can +be specified simultaneously. -```{r restrict_range, fig.cap = "Sample `raman_hdpe` spectrum with different degrees of background subtraction (Cowger et al., 2020).", fig.width=6, fig.align="center"} +```{r restrict_range, fig.cap = "Sample `raman_hdpe` spectrum with different degrees of range restriction.", fig.width=6, fig.align="center"} data("raman_hdpe") -none <- adj_intens(raman_hdpe) +none <- make_rel(raman_hdpe) r1 <- restrict_range(raman_hdpe, min = 1000, max = 2000) r2 <- restrict_range(raman_hdpe, min = c(1000, 1800), max = c(1200, 2000)) @@ -490,11 +598,7 @@ compare_ranges <- c_spec(list(none, r1, r2), range = "common") plot(compare_ranges) ``` -#### Flattening Ranges - -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/flattten_example.jpg") -``` +## Flattening Ranges Sometimes there are peaks that really shouldn't be in your spectra and can distort your interpretation of the spectra but you don't necessarily want to @@ -502,15 +606,15 @@ remove the regions from the analysis because you believe those regions should be flat instead of having a peak. One way to deal with this is to replace the peak values with the mean of the values around the peak. This is the purpose of the `flatten_range` function. By default it is set to flatten the CO2 region for -FTIR spectra because that region often needs to be flattened when atmouspheric +FTIR spectra because that region often needs to be flattened when atmospheric artefacts occur in spectra. Like `restrict_range`, the R function can accept multiple ranges. ```{r flatten_range, fig.cap = "Sample `raman_hdpe` spectrum with different degrees of background subtraction (Cowger et al., 2020).", fig.width=6, fig.align="center"} single <- filter_spec(spectral_map, 120) # Function to filter spectra by index # number or name or a logical vector. -none <- adj_intens(single) -f1 <- flatten_range(single) +none <- make_rel(single) +f1 <- flatten_range(single) #default flattening the CO2 region. f2 <- flatten_range(single, min = c(1000, 2500), max = c(1200, 3000)) compare_flats <- c_spec(list(none, f1, f2)) @@ -518,157 +622,227 @@ compare_flats <- c_spec(list(none, f1, f2)) plot(compare_flats) ``` -#### Min-Max Normalization - -Min-Max normalization is used throughout the app but is not one of the options -the user can specify. Often we regard spectral intensities as arbitrary and -min-max normalization allows us to view spectra on the same scale without -drastically distorting their shapes or relative peak intensities. In the app, it -is primarily used for plotting and as a processing step before correlation -analysis. In the package, most of the processing functions will min-max -transform your spectra relative by default if you do not specify otherwise. This -plot shows the two spectra on such different scales that one looks like a -straight line. - -```{r adj_intens, fig.cap = "Sample `raman_hdpe` spectrum with different degrees of background subtraction (Cowger et al., 2020).", fig.width=6, fig.align="center"} -none <- adj_intens(raman_hdpe, make_rel = F) -relative <- adj_intens(raman_hdpe, make_rel = T) +## Min-Max Normalization -compare_rel <- c_spec(list(none, relative)) - -plot(compare_rel) -``` +Often we regard spectral intensities as arbitrary and min-max normalization +allows us to view spectra on the same scale without drastically distorting their +shapes or relative peak intensities. In the package, most of the processing +functions will min-max transform your spectra by default if you do not specify +otherwise. This plot shows two plots that are nearly identical except for the +min-max transformed spectrum has a y axis that ranges from 0-1. -## Identifying Spectra +```{r make_rel, fig.cap = "Sample `raman_hdpe` spectrum with one being relative and the other untransformed.", fig.width=6, fig.align="center"} +none <- raman_hdpe +relative <- make_rel(raman_hdpe) +plot(none) +plot(relative) -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/identification.jpg") ``` -After uploading data and preprocessing it (if desired) you can now identify the -spectrum. To identify the spectrum go to the **Identification** box. Pro tip: if -you select *Identification** without uploading data to the app, you'll be able -to explore the library by itself. +# Identifying Spectra -### Reading Libraries +## Reading Libraries -These options define the strategy for identification. The ID Library will inform -which library is used. Both (default) will search both FTIR and Raman libraries. -Deriv will search against a derivative transformed library. No Baseline will -search against a baseline corrected library. This should be in line with how you -choose to process your spectra. Cor options use a simple Pearson correlation -search algorithm. AI is currently experimental and uses either a multinomial -model or correlation on mediod spectra from the library. Correlation -thresholding will set the minimum value from matching to use as a 'positive -identification' - -In the R package, you'll download and load the libraries into your working -environment. +Reference libraries are spectra with known identities. The Open Specy library +now has over 10,000 spectra in it an is getting so large that we cannot fit it +within the R package size limit of 5 MB. We host the reference libraries on +[OSF](https://osf.io/x7dpz/) and have a function to pull the libraries down +automatically. Running get_lib by itself will download all libraries to your +package director or you can specify which libraries you want and where you want +them. ```{r eval=F} -get_lib() #will load all libraries unless specified otherwise. -lib <- load_lib(type = "derivative") +get_lib() +get_lib(type = "derivative", path = "your/path") ``` -### Matches +After download the libraries you can load them into your active environment. You +should load them one at a time. Libraries are also `OpenSpecy` objects so you +can use any Open Specy object as a library. -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/matches.jpg") +```{r eval=F} +lib <- load_lib(type = "derivative") ``` -Top matches in the app can be assessed by clicking the cog in the right hand -corner of the Spectra box. This will open a side window with the matches sorted -from most to least similar. Clicking on rows in the table will add the selected -match to the spectra viewer. Using the table's filter options, you can restrict -the range of Pearson\'s r values or search for specific material types. - -The same table can be returned using the OpenSpecy package commands in the R -console. +## Matches + +Before attempting to use a reference library to identify spectra it is really +important to understand what format the reference library is in. All the +OpenSpecy reference libraries are in Absorbance units. `derivative` has been +absolute first derivative transformed, `nobaseline` has been baseline corrected, +`raw` is the rawest form of the reference spectra (not recommended except for +advanced uses). The previously mentioned libraries all have Raman and FTIR +spectra in them. `mediod` is the mediod compressed library version of the +absolute first derivative for FTIR, `model` is an exception because it is a +multinomial regression approach for FTIR to identification of absolute +derivative spectra. All of the libraries have wavenumbers at 5 cm^-1 resolution. +Whichever library you choose, you first need to get your spectra into a similar +enough format to use for comparison. That means conforming the wavenumbers to +the same range and processing the spectra using the same processing procedures. +In this example we use the `data("test_lib")` which is a subsampled version of +the absolute derivative library and `data("raman_hdpe")` which is an unprocessed +Raman spectrum in absorbance units of HDPE plastic. With single spectra it is +easy to look at the spectra but when doing in batch, also refer to the +Thresholding Signal-Noise section for guidance on making sure your batch spectra +are processed to quality specs. ```{r} data("test_lib") +data("raman_hdpe") -processed <- process_spec(ps_example, - conform_spec_args = list(range = test_lib$wavenumber, - res = NULL)) +processed <- process_spec(x = raman_hdpe, + conform_spec_args = list(range = test_lib$wavenumber, + res = NULL, + type = "interp" + )) -test <- match_spec(processed, test_lib, add_object_metadata = "col_id", - add_library_metadata = "sample_name") +# Check to make sure that the signal to noise ratio of the processed spectra is +# greater than 10. +print(sig_noise(processed) > 10) +plotly_spec(raman_hdpe, processed) ``` -### Selection Metadata +After your spectra is processed similarly to the library specifications, you can +identify the spectra using `match_spec()`. The `add_library_metadata` and +`add_object_metadata` options specify the column name in the metadata that you +want to add metadata from and `top_n` specifies how many matches you want. In +this example we just identified a single spectrum with the library but you can +also send an OpenSpecy object with multiple spectra. The output `matches` is a +data.table with at least 3 columns, `object_id` tells you the column names of +the spectra in `x`, `library_id` tells you the column names from the library +that it matched to. `match_val` is the value of the Pearson correlation +coefficient (default) or other correlation if specified in `...` or if using the +model identification option `match_val` will be the model confidence. The output +in this example returned the correct material type, HDPE, as the top match. If +using Pearson correlation, 0.7 is a good threshold to use for a positive ID. In +this example, only our top match is greater than the threshold so we would +disregard the other matches. If no matches were above our threshold, we would +proclaim that the spectrum is of an unknown identity. You'll also notice in this +example that we matched to a library with both Raman and FTIR spectra but the +Raman spectra had the highest hits, this is the rationale for lazily matching to +a library with both. If you want to just match to a library with FTIR or Raman +spectra, you can first filter the library using `filter_spec()` using +`SpectrumType`. -```{r, fig.align="center", out.width="98%", echo=FALSE} -knitr::include_graphics("sop/selectionmetadata.jpg") +```{r} +matches <- match_spec(x = processed, library = test_lib, + add_library_metadata = "sample_name", top_n = 5) +print(matches[,c("object_id", "library_id", "match_val", "SpectrumType", + "SpectrumIdentity")]) ``` -Whatever match is selected from the match table may have additional metadata -about it. That metadata will be displayed below the plot. Some of this metadata -may assist you in interpreting the spectra. For example, if the spectra has -metadata which says it is a liquid and you are analyzing a solid particle, that -spectrum may not be the best match. +## Library Metadata -The R command for manual metadata selection using `sample_name == 5381` as -example is: +The libraries we have created have over 100 variables of metadata in them and +this can be onerous to read through especially given that many of the variables +are `NA` values. We created `get_metadata()` to remedy this by removing columns +from the metadata which are all blank values. The function below will return the +metadata for the top match in `matches`. Remember, similar `filter_spec()`, you +can specify `logic` for more than one thing at a time. -```{r, eval=FALSE} -find_spec(sample_name == 5381, library = spec_lib, which = "raman") +```{r} +get_metadata(x = test_lib, logic = matches[[1,"library_id"]], rm_empty = T) ``` -### Match Plot +## Plot Matches -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/matchplot.png") -``` +Overlaying unknown spectra and the best matches can be extremely useful to +identify peaks that don't fit to the reference library which may need further +investigation. The example below shows great correspondence between the best +match and the unknown spectrum. All major peaks are accounted for and the +correct relative height. There are two small peaks in the unknown spectrum near +500 that are not accounted for which could be investigated further but we would +call this a positive id to HDPE. -This plot is dynamically updated by selecting matches from the match table. The -red spectrum is the spectrum that you selected from the reference library and -the white spectrum is the spectrum that you are trying to identify. Whenever a -new dataset is uploaded, the plot and data table in this tab will be updated. -These plots can be saved as a .png by clicking the camera button at the top of -the plot. +```{r} +plotly_spec(processed, filter_spec(test_lib, logic = matches[[1,"library_id"]])) +``` -### Sharing Reference Data +## Sharing Reference Data If you have reference data or AI models that you think would be useful for sharing with the spectroscopy community through OpenSpecy please contact the -website administrator to discuss options for collaborating. +package administrator to discuss options for collaborating. + +# Characterizing Particles + +Sometimes the spectroscopy task we want to perform is to identify particles in a +spectral map. This is especially common for microplastic analysis where a +spectral map is used to image a sample and spectral information is used to +differentiate microplastic particles from nonplastic particles. In addition to +the material id, one often wants to measure the shape and size of the particles. +In a brute force technique, one could first identify every spectrum in the map, +then use thresholding and image analysis to measure the particles. However, more +often than not, particles are well separated on the image surface and background +spectra is quite different from particle spectra and therefore we can use +thresholding a priori to identify and measure the particles, then pass an +exemplary spectrum for each particle to the identification routine. It is +important to note here that this is at the bleeding edge of theory and technique +so we may be updating these functions in the near future. + +## Brute Force + +```{r eval = F} +data("test_lib") +test_map <- read_any(read_extdata("CA_tiny_map.zip")) -## Characterizing Particles +test_map_processed <- process_spec(test_map, conform_spec_args = list( + range = test_lib$wavenumber, res = NULL) + ) -### Thresholding +identities <- match_spec(test_map_processed, test_lib, order = test_map, + add_library_metadata = "sample_name", top_n = 1) -### Collapse Spectra +features <- ifelse(identities$match_val > 0.7, + tolower(identities$polymer_class), "unknown") -## Additional App Features +id_map <- def_features(x = test_map_processed, features = features) -Accessibility is extremely important to us and we are making strives to improve -the accessibility of Open Specy for all spectroscopists. Please reach out if you -have ideas for improvement. +id_map$metadata$identities <- features +# Also should probably be implemented automatically in the function when a +# character value is provided. -We added a Google translate plugin to all pages in the app so that you can -easily translate the app. We know that not all languages will be fully supported -but we will continue to try and improve the translations. +# Collapses spectra to their median for each particle +test_collapsed <- collapse_spec(id_map) +``` -### Download Data +## A Priori Particle Thresholding -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/download.jpg") -``` +```{r eval = F} +data("test_lib") +test_map <- read_any(read_extdata("CA_tiny_map.zip")) + +# Characterize the total signal as a threshold value. +snr <- sig_noise(test_map,metric = "log_tot_sig") + +# Use this to find your particles and the sig_noise value to use for +# thresholding. +heatmap_spec(test_map, z = snr) + +# Set define the feature regions based on the threshold. 400 appeared to be +# where I suspected my particle to be in the previous heatmap. +id_map <- def_features(x = test_map, features = snr > 400) + +# Check that the thresholding worked as expected. +heatmap_spec(id_map, z = id_map$metadata$feature_id) + +# Collapse the spectra to their medians based on the threshold. Important to +# note here that the particles with id -88 are anything from the FALSE values +# so they should be your background. +collapsed_id_map <- id_map |> + collapse_spec() -After you have the preprocessing parameters set, we recommend that you download -the preprocessed data for your records. The download data button will append the -uploaded data to three columns created by the preprocessing parameters. -"Wavelength" and "Absorbance" are columns from the data uploaded by the user. -"NormalizedIntensity" is the max-min normalized value (Equation 1) of the -"Absorbance". "Smoothed" is the Savitzky-Golay filter specified by the slider -explained above. "BaselineRemoved" is the smoothed and baseline corrected value -that is visible on the center plot. +# Process the collapsed spectra. +id_map_processed <- process_spec(collapsed_id_map, conform_spec_args = list( + range = test_lib$wavenumber, res = NULL) + ) -## Conceptual diagram of Open Specy +# Check the spectra. +plot(id_map_processed) -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/flowchart.png") +# Get the matches of the collapsed spectra for the particles. +matches <- match_spec(id_map_processed, test_lib, + add_library_metadata = "sample_name", top_n = 1) ``` # References diff --git a/vignettes/sop.html b/vignettes/sop.html index 2b89eca7..ce4fb386 100644 --- a/vignettes/sop.html +++ b/vignettes/sop.html @@ -11,9 +11,9 @@ - + -Open Specy Tutorial +Open Specy Package Tutorial -
                               # if it isn't already but are not changing your
-                               # underlying data
+
plotly_spec(scratch_OpenSpecy) # this will min-max normalize your data even
+
+ +
                               # if it isn't already but are not changing your
+                               # underlying data

With interactive plots you can also plot two different datasets simultaneously to compare.

-
plotly_spec(scratch_OpenSpecy, rds_example) # interactive comparison
-
- +
plotly_spec(scratch_OpenSpecy, rds_example) # interactive comparison
+
+ -
-

Maps

+
+

Maps

Spectral maps can also be visualized as overlaid spectra but in -addition the spatial information can be plotted as a heatmap. A similar -plot should appear in the app if you upload multiple spectra or a -spectral map. It is important to note that when multiple spectra are -uploaded in batch they are prescribed x and y -coordinates, this can be helpful for visualizing summary statistics and -navigating vast amounts of data but shouldn’t be confused with data -which actually has spatial coordinates. Hovering over the map will -reveal information about the signal and noise and correlation values. -Clicking the map will provide the selected spectrum underneath. In the -app, the colors of the heatmap are either the signal and noise or the -spectral identifications depending on whether the identification is -turned on or not. Pixels are black if the spectra does not pass the -signal-noise threshold and/or the correlation threshold.

-

-

The same plot can be created in R but the user can control what -values form the colors of the heatmap by specifying z.

-
heatmap_spec(spectral_map,
-             z = spectral_map$metadata$x) # will show more advanced example 
-
- -
                                          # later in tutorial
+addition the spatial information can be plotted as a heatmap. It is +important to note that when multiple spectra are uploaded in batch they +are prescribed x and y coordinates, this can +be helpful for visualizing summary statistics and navigating vast +amounts of data but shouldn’t be confused with data which actually has +spatial coordinates.

+

In R the user can control what values form the colors of the heatmap +by specifying z, it is handy to put the information you +want in the metadata for this reason. This example just shows the x +values of the spectra.

+
heatmap_spec(spectral_map,
+             z = spectral_map$metadata$x) # will show more advanced example 
+
+ +
                                          # later in tutorial

An interactive plot of the heatmap and spectra overlayed can be generated with the interactive_plot() function. A user can hover over the heatmap to identify the next row id they are interested in and update the value of select to see that spectrum.

-
interactive_plot(spectral_map, select = 100, z = spectral_map$metadata$x)
-
- +
interactive_plot(spectral_map, select = 100, z = spectral_map$metadata$x)
+
+ +
+
+
+

Combining OpenSpecy Objects

+

Sometimes you have several OpenSpecy objects that you want to combine +into one. The easiest way to do that is by having spectra which all are +in the same exact format with the same series of wavenumbers. The +default settings of c_spec assume that is the case.

+
#Same range example. 
+data("raman_hdpe")
+raman_hdpe1 <- raman_hdpe
+raman_hdpe2 <- raman_hdpe
+
+raman_hdpe1$spectra <- raman_hdpe1$spectra - 100 
+
+raman_hdpe2$spectra <- raman_hdpe2$spectra - 200 
+
+combined_spectra <- c_spec(list(raman_hdpe, raman_hdpe1, raman_hdpe2))
+
+plot(combined_spectra)
+

+

If you have different wavenumber ranges for the spectra you want to +combibne, you can set range = "common" and res +equal to the wavenumber resolution you want and the function will +collapse all the spectra to whatever their common range is using linear +interpolation.

+
asp_example <- read_extdata("ftir_ldpe_soil.asp") |> 
+  read_any() 
+ps_example <- read_extdata("ftir_ps.0") |> 
+  read_any() # preserves some metadata
+
+combined <- c_spec(list(asp_example, ps_example), range = "common", res = 5)
+
+plot(combined)
+

+
+

Filtering OpenSpecy Objects

+

OpenSpecy objects can have any number of spectra in them. To create +an OpenSpecy with a subset of the spectra that is in an Open Specy +object you can use the filter_spec function. Filtering is +allowed by index number, name, or using a logical vector. Filtering will +update the spectra and metadata items of the +OpenSpecy but not the wavenumber.

+

+spectral_map <- read_extdata("CA_tiny_map.zip") |> 
+  read_any() # preserves some metadata
+
+filter_spec(spectral_map, 150) |>
+    plot() #extracts the 150th spectrum. 
+

+
filter_spec(spectral_map, "8_5") |>
+    plot() #extracts the spectrum with column name "8_5". 
+

+
filter_spec(spectral_map, names(spectral_map$spectra) %in% c("8_5", "10_3", "12_10")) |>
+    plot() #extracts the spectra with column names "8_5", "10_3", or "12_10" using a logical argument based on the spectra column names.
+

+
filter_spec(spectral_map, spectral_map$metadata$y == 1) |>
+    plot() #extracts the spectra with y values of 1 for their metadata. 
+

-
-

Processing

-

After uploading data, you can preprocess the data using intensity -adjustment, baseline subtraction, smoothing, flattening, and range -selection and save your preprocessed data. Once the process button is -selected the default processing will initiate. This is an absolute +

+

Sampling OpenSpecy Objects

+

Sometimes you have a large suite of examples of spectra of the same +material and you want to reduce the number of spectra you run through +the analysis for simplicity or you are running simulations or other +proceedures that require you to first sample from the spectra contained +in your OpenSpecy objects before doing analysis. the +sample_spec function serves this purpose.

+
spectral_map <- read_any(read_extdata("CA_tiny_map.zip"))
+sample_spec(spectral_map,size = 10) |>
+    plot()
+

+
+
+

Processing

+

The goal of this all processing is to increase the signal to noise +ratio (S/N) of the spectra without distorting the shape, position, or +relative size of the peaks. After loading data, you can preprocess the +data using intensity adjustment, baseline subtraction, smoothing, +flattening, and range selection. The default settings is an absolute derivative transformation, it is kind of magic, it does something similar to smoothing, baseline subtraction, and intensity correction -simultaneously and really quickly. By clicking preprocessing you should -see the spectrum update with the processed spectra.

-

-

To view the raw data again, just deselect preprocessing. Toggling -preprocessing on and off can help you to make sure that the spectra are -processed correctly.

+simultaneously and really quickly.

The process_spec() function is a monolithic function for all processing procedures which is optimized by default to result in high signal to noise in most cases, same as the app.

-
processed <- process_spec(raman_hdpe)
-summary(processed)
-#> $wavenumber
-#>  Length Min. Max. Res.
-#>     579  305 3195    5
-#> 
-#> $spectra
-#>  Number Min. Intensity Max. Intensity
-#>       1              0              1
-#> 
-#> $metadata
-#>   Min. Max.
-#> x    1    1
-#> y    1    1
-#> [1] "x"                 "y"                 "user_name"        
-#> [4] "spectrum_type"     "spectrum_identity" "organization"     
-#> [7] "license"           "session_id"        "file_id"
-summary(raman_hdpe)
-#> $wavenumber
-#>  Length   Min.    Max. Res.
-#>     964 301.04 3198.12 2.54
-#> 
-#> $spectra
-#>  Number Min. Intensity Max. Intensity
-#>       1             26            816
-#> 
-#> $metadata
-#>   Min. Max.
-#> x    1    1
-#> y    1    1
-#> [1] "x"                 "y"                 "user_name"        
-#> [4] "spectrum_type"     "spectrum_identity" "organization"     
-#> [7] "license"           "session_id"        "file_id"
+
processed <- process_spec(raman_hdpe,
+                          active = TRUE,
+                          adj_intens = FALSE,
+                          adj_intens_args = list(type = "none"),
+                          conform_spec = TRUE,
+                          conform_spec_args = list(range = NULL, res = 5, type = "interp"),
+                          restrict_range = FALSE,
+                          restrict_range_args = list(min = 0, max = 6000),
+                          flatten_range = FALSE,
+                          flatten_range_args = list(min = 2200, max = 2420),
+                          subtr_baseline = FALSE,
+                          subtr_baseline_args = list(type = "polynomial", degree = 8, raw = FALSE, baseline =
+                            NULL),
+                          smooth_intens = TRUE,
+                          smooth_intens_args = list(polynomial = 3, window = 11, derivative = 1, abs = TRUE),
+                          make_rel = TRUE)
+summary(processed)
+#> $wavenumber
+#>  Length Min. Max. Res.
+#>     579  305 3195    5
+#> 
+#> $spectra
+#>  Number Min. Intensity Max. Intensity
+#>       1              0              1
+#> 
+#> $metadata
+#>   Min. Max.
+#> x    1    1
+#> y    1    1
+#> [1] "x"                 "y"                 "user_name"        
+#> [4] "spectrum_type"     "spectrum_identity" "organization"     
+#> [7] "license"           "session_id"        "file_id"
+summary(raman_hdpe)
+#> $wavenumber
+#>  Length   Min.    Max. Res.
+#>     964 301.04 3198.12 2.54
+#> 
+#> $spectra
+#>  Number Min. Intensity Max. Intensity
+#>       1             26            816
+#> 
+#> $metadata
+#>   Min. Max.
+#> x    1    1
+#> y    1    1
+#> [1] "x"                 "y"                 "user_name"        
+#> [4] "spectrum_type"     "spectrum_identity" "organization"     
+#> [7] "license"           "session_id"        "file_id"

You can compare the processed and unprocessed data in an overlay -plot:

-
plotly_spec(raman_hdpe, processed)
+plot.

+
plotly_spec(raman_hdpe, processed)

We want people to use the process_spec() function for -most processing operations, all processing functions can be tuned using -its parameters. However, we recognize that nesting of functions and -order of operations can be useful for users to control so you can also -use individual functions for each operation if you’d like.

-
-

Threshold Signal-Noise

+most processing operations. All other processing functions can be tuned +using its parameters in the single function and we have set defaults and +nested the processing functions in a way that provides typically quality +results. However, we recognize that nesting of functions and order of +operations can be useful for users to control so you can also use +individual functions for each operation if you’d like. See explanations +of each processing sub-function below.

+
+

Threshold Signal-Noise

Considering whether you have enough signal to analyze spectra is important. Classical spectroscopy would recommend your highest peak to be at least 10 times the baseline of your processed spectra before you begin analysis. If your spectra is below that threshold, you may want to consider recollecting it. In practice, we are rarely able to collect -spectra of that good quality and more often use 5. The Signal Over Noise -technique searches your spectra for high and low regions and conducts -division on them to derive the signal to noise ratio. Signal Times Noise -multiplies the mean signal by the standard deviation of the signal and -Total Signal sums the intensities. The latter can be really useful for -thresholding spectral maps to identify particles which we will discuss -later.

-

+spectra of that good quality and more often use 5. The Run Signal Over +Noise technique searches your spectra for high and low regions and +conducts division on them to derive the signal to noise ratio. This is a +good way to automatically calculate the signal to noise ratio. In the +example below you can see that our signal to noise ratio is increased by +the processing, the goal of processing is generally to maximize this. +Signal Times Noise multiplies the mean signal by the standard deviation +of the signal and Total Signal sums the intensities. The latter can be +really useful for thresholding spectral maps to identify particles which +we will discuss later.

+
sig_noise(processed, metric = "run_sig_over_noise") > sig_noise(raman_hdpe, metric = "run_sig_over_noise")
+#> intensity 
+#>      TRUE

If analyzing spectra in batch, we recommend looking at the heatmap and optimizing the percent of spectra that are above your signal to noise threshold to determine the correct settings instead of looking -through spectra individually. Good Signal tells you what percent of your -data are above your signal threshold.

-

+through spectra individually. Setting the min_sn will +threshold the heatmap image to only color spectra which have a +sn value over the threshold.

+
spectral_map_p <- spectral_map |>
+    process_spec(flatten_range = T)
+
+spectral_map_p$metadata$sig_noise <- sig_noise(spectral_map_p)
+
+heatmap_spec(spectral_map_p, sn = spectral_map_p$metadata$sig_noise, min_sn = 5)
+
+
-
-

Intensity Adjustment

-

+
+

Intensity Adjustment

Open Specy assumes that intensity units are in absorbance units but Open Specy can adjust reflectance or transmittance spectra to absorbance -units using this selection. The transmittance adjustment uses the \(\log_{10} 1/T\) calculation which does not +units. The transmittance adjustment uses the \(\log_{10} 1/T\) calculation which does not correct for system or particle characteristics. The reflectance -adjustment use the Kubelka-Munk equation \(\frac{(1-R)^2}{2R}\). If none is selected, -Open Specy assumes that the uploaded data is an absorbance spectrum and -does not apply an adjustment.

+adjustment use the Kubelka-Munk equation \(\frac{(1-R)^2}{2R}\).

This is the respective R code for a scenario where the spectra doesn’t need intensity adjustment:

-
raman_adj <- raman_hdpe |> 
-  adj_intens(type = "none")
-
-plot(raman_adj)
-

+
data("raman_hdpe")
+trans_raman_hdpe <- raman_hdpe
+trans_raman_hdpe$spectra <- 2 - trans_raman_hdpe$spectra^2
+    
+rev_trans_raman_hdpe <- trans_raman_hdpe |> 
+                            adj_intens(type = "transmittance")
+
+plotly_spec(trans_raman_hdpe, rev_trans_raman_hdpe)
+
+
-
-

Conforming

+
+

Conforming

Conforming spectra is essential before comparing to a reference library and can be useful for summarizing data when you don’t need it to -be highly resolved spectrally. In the app, conforming happens behind the -scenes without any user input to make sure that the spectra the user -uploads and the library spectra will be compatible. In code, we set the -default spectral resolution to 5 because this tends to be pretty good -for a lot of applications and is in between 4 and 8 which are commonly -used wavenumber resolutions. There are several ways that this is -function can be specified.

-
conform_spec(raman_hdpe) |> # default convert res to 5 wavenumbers.
-  summary()
-#> $wavenumber
-#>  Length Min. Max. Res.
-#>     579  305 3195    5
-#> 
-#> $spectra
-#>  Number Min. Intensity Max. Intensity
-#>       1       41.65339       662.0608
-#> 
-#> $metadata
-#>   Min. Max.
-#> x    1    1
-#> y    1    1
-#> [1] "x"                 "y"                 "user_name"        
-#> [4] "spectrum_type"     "spectrum_identity" "organization"     
-#> [7] "license"           "session_id"        "file_id"
-
-# Force one spectrum to have the exact same wavenumbers as another
-conform_spec(asp_example, range = ps_example$wavenumber, res = NULL) |> 
-  summary()
-#> $wavenumber
-#>  Length     Min.     Max.     Res.
-#>    2126 399.2239 4497.537 1.928618
-#> 
-#> $spectra
-#>  Number Min. Intensity Max. Intensity
-#>       1             NA             NA
-#> 
-#> $metadata
-#>   Min. Max.
-#> x    1    1
-#> y    1    1
-#> [1] "x"          "y"          "file_name"  "license"    "session_id"
-#> [6] "file_id"    "col_id"
-
-# Specify the wavenumber resolution and use a rolling join instead of linear
-# approximation (faster for large datasets). 
-conform_spec(spectral_map, range = ps_example$wavenumber, res = 10,
-             type = "roll") |> 
-  summary()
-#> $wavenumber
-#>  Length Min. Max. Res.
-#>     329  720 4000   10
-#> 
-#> $spectra
-#>  Number Min. Intensity Max. Intensity
-#>     208      -1.317072       1.168225
-#> 
-#> $metadata
-#>   Min. Max.
-#> x    0   12
-#> y    0   15
-#>  [1] "x"             "y"             "file_name"     "license"      
-#>  [5] "description"   "samples"       "lines"         "bands"        
-#>  [9] "header offset" "data type"     "interleave"    "z plot titles"
-#> [13] "pixel size"    "session_id"    "file_id"       "col_id"
+be highly resolved spectrally. We set the default spectral resolution to +5 because this tends to be pretty good for a lot of applications and is +in between 4 and 8 which are commonly used wavenumber resolutions. There +are several ways that this function can be specified.

+
conform_spec(raman_hdpe) |> # default convert res to 5 wavenumbers.
+  summary()
+#> $wavenumber
+#>  Length Min. Max. Res.
+#>     579  305 3195    5
+#> 
+#> $spectra
+#>  Number Min. Intensity Max. Intensity
+#>       1       41.65339       662.0608
+#> 
+#> $metadata
+#>   Min. Max.
+#> x    1    1
+#> y    1    1
+#> [1] "x"                 "y"                 "user_name"        
+#> [4] "spectrum_type"     "spectrum_identity" "organization"     
+#> [7] "license"           "session_id"        "file_id"
+
+# Force one spectrum to have the exact same wavenumbers as another
+conform_spec(asp_example, range = ps_example$wavenumber, res = NULL) |> 
+  summary()
+#> $wavenumber
+#>  Length     Min.     Max.     Res.
+#>    1736 651.8728 3998.025 1.928618
+#> 
+#> $spectra
+#>  Number Min. Intensity Max. Intensity
+#>       1   0.0009781419      0.5156938
+#> 
+#> $metadata
+#>   Min. Max.
+#> x    1    1
+#> y    1    1
+#> [1] "x"          "y"          "file_name"  "license"    "session_id"
+#> [6] "file_id"    "col_id"
+
+# Specify the wavenumber resolution and use a rolling join instead of linear
+# approximation (faster for large datasets). 
+conform_spec(spectral_map, range = ps_example$wavenumber, res = 10,
+             type = "roll") |> 
+  summary()
+#> $wavenumber
+#>  Length Min. Max. Res.
+#>     329  720 4000   10
+#> 
+#> $spectra
+#>  Number Min. Intensity Max. Intensity
+#>     208      -1.317072       1.168225
+#> 
+#> $metadata
+#>   Min. Max.
+#> x    0   12
+#> y    0   15
+#>  [1] "x"             "y"             "file_name"     "license"      
+#>  [5] "description"   "samples"       "lines"         "bands"        
+#>  [9] "header offset" "data type"     "interleave"    "z plot titles"
+#> [13] "pixel size"    "session_id"    "file_id"       "col_id"
-
-

Smoothing

-

-

The first step of the Open Specy preprocessing routing is spectral -smoothing. The goal of this function is to increase the signal to noise -ratio (S/N) without distorting the shape or relative size of the peaks. -The value on the slider is the polynomial order of the Savitzky-Golay -(SG) filter. Higher numbers lead to more wiggly fits and thus less -smooth, lower numbers lead to more smooth fits. The SG filter is fit to -a moving window of 11 data points by default where the center point in -the window is replaced with the polynomial estimate. Larger windows will -produce smoother fits. The derivative order is set to 1 by default which -transforms the spectra to their first derivative. A zero order -derivative will have no derivative transformation. When smoothing is -done well, peak shapes and relative heights should not change. The -absolute value is primarily useful for first derivative spectra where -the absolute value results in an absorbance-like spectrum which is why -we set it as the default.

-

Examples of smoothing with the R package:

-
library(ggplot2)
-data("raman_hdpe")
-
-none <- adj_intens(raman_hdpe, make_rel = T)
-p1 <- smooth_intens(raman_hdpe, polynomial = 1)
-p4 <- smooth_intens(raman_hdpe, polynomial = 4)
-
-compare_derivative <- c_spec(list(none, p1, p4))
-    
-plot(compare_derivative)
+
+

Smoothing

+

The Savitzky-Golay +filter is used for smoothing. Higher polynomial numbers lead to more +wiggly fits and thus less smoothing, lower numbers lead to more smooth +fits. The SG filter is fit to a moving window of 11 data points by +default where the center point in the window is replaced with the +polynomial estimate. Larger windows will produce smoother fits. The +derivative order is set to 1 by default which transforms the spectra to +their first derivative. A zero order derivative will have no derivative +transformation. When smoothing is done well, peak shapes and relative +heights should not change. The absolute value is primarily useful for +first derivative spectra where the absolute value results in an +absorbance-like spectrum which is why we set it as the default. You’ll +notice a new function we are using c_spec() which is used +to combine spectral objects into one OpenSpecy object.

+

Examples of smoothing:

+
data("raman_hdpe")
+
+none <- make_rel(raman_hdpe)
+p1 <- smooth_intens(raman_hdpe, polynomial = 1, derivative = 0, abs = F)
+p4 <- smooth_intens(raman_hdpe, polynomial = 4, derivative = 0, abs = F)
+
+compare_derivative <- c_spec(list(none, p1, p4))
+    
+plot(compare_derivative)
-Sample `raman_hdpe` spectrum with different smoothing polynomials from Cowger et al. (2020). +Sample `raman_hdpe` spectrum with different smoothing polynomials.

Sample raman_hdpe spectrum with different smoothing -polynomials from Cowger et al. (2020). +polynomials.

-

Derivative transformation can happen with the same function. You’ll -notice a new function we are using c_spec() which is used -to combine spectral objects into one OpenSpecy object.

-

-none <- adj_intens(raman_hdpe, make_rel = T)
-d1 <- smooth_intens(raman_hdpe, derivative = 1, abs = T)
-d2 <- smooth_intens(raman_hdpe, derivative = 2)
-
-compare_derivative <- c_spec(list(none, d1, d2))
-    
-plot(compare_derivative)
+

Derivative transformation can happen with the same function.

+

+none <- make_rel(raman_hdpe)
+d1 <- smooth_intens(raman_hdpe, derivative = 1, abs = T)
+d2 <- smooth_intens(raman_hdpe, derivative = 2, abs = T)
+
+compare_derivative <- c_spec(list(none, d1, d2))
+    
+plot(compare_derivative)
-Sample `raman_hdpe` spectrum with different derivatives from Cowger et al. (2020). +Sample `raman_hdpe` spectrum with different derivatives.

-Sample raman_hdpe spectrum with different derivatives from -Cowger et al. (2020). +Sample raman_hdpe spectrum with different derivatives.

-
-

Baseline Correction

-

-

The second step of Open Specy’s preprocessing routine is baseline -correction. The goal of baseline correction is to get all non-peak -regions of the spectra to zero absorbance. The higher the polynomial -order, the more wiggly the fit to the baseline. If the baseline is not -very wiggly, a more wiggly fit could remove peaks which is not desired. -The baseline correction algorithm used in Open Specy is called -“iModPolyfit†(Zhao et al. 2007). This algorithm iteratively fits -polynomial equations of the specified order to the whole spectrum. -During the first fit iteration, peak regions will often be above the -baseline fit. The data in the peak region is removed from the fit to -make sure that the baseline is less likely to fit to the peaks. The -iterative fitting terminates once the difference between the new and -previous fit is small. An example of a good baseline fit below. For -those who have been with OpenSpecy a while, you will notice the app no -longer supports manual baseline correction, it was a hard choice but it -just didn’t make sense to keep it now that we are moving toward high -throughput automated methods. It does still exist in the R function -though.

-
data("raman_hdpe")
-
-none <- adj_intens(raman_hdpe)
-d2 <- subtr_baseline(raman_hdpe, degree = 2)
-d8 <- subtr_baseline(raman_hdpe, degree = 8)
-
-compare_subtraction <- c_spec(list(none, d2, d8))
-
-plot(compare_subtraction)
+
+

Baseline Correction

+

The goal of baseline correction is to get all non-peak regions of the +spectra to zero absorbance. The higher the polynomial order, the more +wiggly the fit to the baseline. If the baseline is not very wiggly, a +more wiggly fit could remove peaks which is not desired. The baseline +correction algorithm used in Open Specy is called “iModPolyfit†(Zhao et +al. 2007). This algorithm iteratively fits polynomial equations of the +specified order to the whole spectrum. During the first fit iteration, +peak regions will often be above the baseline fit. The data in the peak +region is removed from the fit to make sure that the baseline is less +likely to fit to the peaks. The iterative fitting terminates once the +difference between the new and previous fit is small. An example of a +good baseline fit below. Manual baseline correction can also be +specified by providing a baseline OpenSpecy object.

+
data("raman_hdpe")
+
+alternative_baseline <- smooth_intens(raman_hdpe, polynomial = 1, window = 51, derivative = 0, abs = F, make_rel = F) |> flatten_range(min = 2700, max = 3200, make_rel = F)
+
+none <- make_rel(raman_hdpe)
+d2 <- subtr_baseline(raman_hdpe, type = "manual", baseline = alternative_baseline)
+d8 <- subtr_baseline(raman_hdpe, degree = 8)
+
+compare_subtraction <- c_spec(list(none, d2, d8))
+
+plot(compare_subtraction)
-Sample `raman_hdpe` spectrum with different degrees of background subtraction (Cowger et al., 2020). +Sample `raman_hdpe` spectrum with different degrees of background subtraction (Cowger et al., 2020).

Sample raman_hdpe spectrum with different degrees of background subtraction (Cowger et al., 2020).

-
-

Range Selection

-

+
+

Range Selection

Sometimes an instrument operates with high noise at the ends of the spectrum and, a baseline fit produces distortions, or there are regions of interest for analysis. Range selection accomplishes those goals. You @@ -2789,31 +2884,29 @@

Range Selection

by wavelength to determine what wavelength ranges to use. Distortions due to baseline fit can be assessed from looking at the preprocess spectra. Additionally, you can restrict the range to examine a single -peak or a subset of peaks of interests. The maximum and minimum -wavenumbers will specify the range to keep in the app. In the R -function, multiple ranges can be specified simultaneously.

-
data("raman_hdpe")
-
-none <- adj_intens(raman_hdpe)
-r1 <- restrict_range(raman_hdpe, min = 1000, max = 2000)
-r2 <- restrict_range(raman_hdpe, min = c(1000, 1800), max = c(1200, 2000))
-
-compare_ranges <- c_spec(list(none, r1, r2), range = "common")
-# Common argument crops the ranges to the most common range between the spectra
-# when joining. 
-
-plot(compare_ranges)
+peak or a subset of peaks of interests. Multiple ranges can be specified +simultaneously.

+
data("raman_hdpe")
+
+none <- make_rel(raman_hdpe)
+r1 <- restrict_range(raman_hdpe, min = 1000, max = 2000)
+r2 <- restrict_range(raman_hdpe, min = c(1000, 1800), max = c(1200, 2000))
+
+compare_ranges <- c_spec(list(none, r1, r2), range = "common")
+# Common argument crops the ranges to the most common range between the spectra
+# when joining. 
+
+plot(compare_ranges)
-Sample `raman_hdpe` spectrum with different degrees of background subtraction (Cowger et al., 2020). +Sample `raman_hdpe` spectrum with different degrees of range restriction.

-Sample raman_hdpe spectrum with different degrees of -background subtraction (Cowger et al., 2020). +Sample raman_hdpe spectrum with different degrees of range +restriction.

-
-

Flattening Ranges

-

+
+

Flattening Ranges

Sometimes there are peaks that really shouldn’t be in your spectra and can distort your interpretation of the spectra but you don’t necessarily want to remove the regions from the analysis because you @@ -2822,18 +2915,18 @@

Flattening Ranges

values around the peak. This is the purpose of the flatten_range function. By default it is set to flatten the CO2 region for FTIR spectra because that region often needs to be -flattened when atmouspheric artefacts occur in spectra. Like +flattened when atmospheric artefacts occur in spectra. Like restrict_range, the R function can accept multiple ranges.

-
single <- filter_spec(spectral_map, 120) # Function to filter spectra by index
-                                         # number or name or a logical vector. 
-none <- adj_intens(single)
-f1 <- flatten_range(single)
-f2 <- flatten_range(single, min = c(1000, 2500), max = c(1200, 3000))
-
-compare_flats <- c_spec(list(none, f1, f2))
-
-plot(compare_flats)
+
single <- filter_spec(spectral_map, 120) # Function to filter spectra by index
+                                         # number or name or a logical vector. 
+none <- make_rel(single)
+f1 <- flatten_range(single) #default flattening the CO2 region. 
+f2 <- flatten_range(single, min = c(1000, 2500), max = c(1200, 3000))
+
+compare_flats <- c_spec(list(none, f1, f2))
+
+plot(compare_flats)
Sample `raman_hdpe` spectrum with different degrees of background subtraction (Cowger et al., 2020).

@@ -2842,143 +2935,239 @@

Flattening Ranges

-
-

Min-Max Normalization

-

Min-Max normalization is used throughout the app but is not one of -the options the user can specify. Often we regard spectral intensities -as arbitrary and min-max normalization allows us to view spectra on the -same scale without drastically distorting their shapes or relative peak -intensities. In the app, it is primarily used for plotting and as a -processing step before correlation analysis. In the package, most of the -processing functions will min-max transform your spectra relative by -default if you do not specify otherwise. This plot shows the two spectra -on such different scales that one looks like a straight line.

-
none <- adj_intens(raman_hdpe, make_rel = F)
-relative <- adj_intens(raman_hdpe, make_rel = T)
-
-compare_rel <- c_spec(list(none, relative))
-
-plot(compare_rel)
+
+

Min-Max Normalization

+

Often we regard spectral intensities as arbitrary and min-max +normalization allows us to view spectra on the same scale without +drastically distorting their shapes or relative peak intensities. In the +package, most of the processing functions will min-max transform your +spectra by default if you do not specify otherwise. This plot shows two +plots that are nearly identical except for the min-max transformed +spectrum has a y axis that ranges from 0-1.

+
none <- raman_hdpe
+relative <- make_rel(raman_hdpe)
+plot(none)
-Sample `raman_hdpe` spectrum with different degrees of background subtraction (Cowger et al., 2020). +Sample `raman_hdpe` spectrum with one being relative and the other untransformed.

-Sample raman_hdpe spectrum with different degrees of -background subtraction (Cowger et al., 2020). +Sample raman_hdpe spectrum with one being relative and the +other untransformed. +

+
+
plot(relative)
+
+Sample `raman_hdpe` spectrum with one being relative and the other untransformed. +

+Sample raman_hdpe spectrum with one being relative and the +other untransformed.

-
-

Identifying Spectra

-

-

After uploading data and preprocessing it (if desired) you can now -identify the spectrum. To identify the spectrum go to the -Identification box. Pro tip: if you select -*Identification** without uploading data to the app, you’ll be able to -explore the library by itself.

-
-

Reading Libraries

-

These options define the strategy for identification. The ID Library -will inform which library is used. Both (default) will search both FTIR -and Raman libraries. Deriv will search against a derivative transformed -library. No Baseline will search against a baseline corrected library. -This should be in line with how you choose to process your spectra. Cor -options use a simple Pearson correlation search algorithm. AI is -currently experimental and uses either a multinomial model or -correlation on mediod spectra from the library. Correlation thresholding -will set the minimum value from matching to use as a ‘positive -identification’

-

In the R package, you’ll download and load the libraries into your -working environment.

-
get_lib() #will load all libraries unless specified otherwise. 
-lib <- load_lib(type = "derivative")
+
+

Identifying Spectra

+
+

Reading Libraries

+

Reference libraries are spectra with known identities. The Open Specy +library now has over 10,000 spectra in it an is getting so large that we +cannot fit it within the R package size limit of 5 MB. We host the +reference libraries on OSF and have +a function to pull the libraries down automatically. Running get_lib by +itself will download all libraries to your package director or you can +specify which libraries you want and where you want them.

+
get_lib() 
+get_lib(type = "derivative", path = "your/path")
+

After download the libraries you can load them into your active +environment. You should load them one at a time. Libraries are also +OpenSpecy objects so you can use any Open Specy object as a +library.

+
lib <- load_lib(type = "derivative")
-
-

Matches

-

-

Top matches in the app can be assessed by clicking the cog in the -right hand corner of the Spectra box. This will open a side window with -the matches sorted from most to least similar. Clicking on rows in the -table will add the selected match to the spectra viewer. Using the -table’s filter options, you can restrict the range of Pearson's r values -or search for specific material types.

-

The same table can be returned using the OpenSpecy package commands -in the R console.

-
data("test_lib")
-
-processed <- process_spec(ps_example, 
-                        conform_spec_args = list(range = test_lib$wavenumber, 
-                                                 res = NULL))
-
-test <- match_spec(processed, test_lib, add_object_metadata = "col_id",
-                   add_library_metadata = "sample_name")
+
+

Matches

+

Before attempting to use a reference library to identify spectra it +is really important to understand what format the reference library is +in. All the OpenSpecy reference libraries are in Absorbance units. +derivative has been absolute first derivative transformed, +nobaseline has been baseline corrected, raw is +the rawest form of the reference spectra (not recommended except for +advanced uses). The previously mentioned libraries all have Raman and +FTIR spectra in them. mediod is the mediod compressed +library version of the absolute first derivative for FTIR, +model is an exception because it is a multinomial +regression approach for FTIR to identification of absolute derivative +spectra. All of the libraries have wavenumbers at 5 cm^-1 resolution. +Whichever library you choose, you first need to get your spectra into a +similar enough format to use for comparison. That means conforming the +wavenumbers to the same range and processing the spectra using the same +processing procedures. In this example we use the +data("test_lib") which is a subsampled version of the +absolute derivative library and data("raman_hdpe") which is +an unprocessed Raman spectrum in absorbance units of HDPE plastic. With +single spectra it is easy to look at the spectra but when doing in +batch, also refer to the Thresholding Signal-Noise section for guidance +on making sure your batch spectra are processed to quality specs.

+
data("test_lib")
+data("raman_hdpe")
+
+processed <- process_spec(x = raman_hdpe, 
+                          conform_spec_args = list(range = test_lib$wavenumber, 
+                                                   res = NULL,
+                                                   type = "interp"
+                                                   ))
+
+print(sig_noise(processed) > 10) #Check to make sure that the signal to noise ratio of the processed spectra is greater than 10.
+#> intensity 
+#>      TRUE
+plotly_spec(raman_hdpe, processed)
+
+ +

After your spectra is processed similarly to the library +specifications, you can identify the spectra using +match_spec(). The add_library_metadata and +add_object_metadata options specify the column name in the +metadata that you want to add metadata from and top_n +specifies how many matches you want. In this example we just identified +a single spectrum with the library but you can also send an OpenSpecy +object with multiple spectra. The output matches is a +data.table with at least 3 columns, object_id tells you the +column names of the spectra in x, library_id +tells you the column names from the library that it matched to. +match_val is the value of the Pearson correlation +coefficient (default) or other correlation if specified in +... or if using the model identification option +match_val will be the model confidence. The output in this +example returned the correct material type, HDPE, as the top match. If +using Pearson correlation, 0.7 is a good threshold to use for a positive +ID. In this example, only our top match is greater than the threshold so +we would disregard the other matches. If no matches were above our +threshold, we would proclaim that the spectrum is of an unknown +identity. You’ll also notice in this example that we matched to a +library with both Raman and FTIR spectra but the Raman spectra had the +highest hits, this is the rationale for lazily matching to a library +with both. If you want to just match to a library with FTIR or Raman +spectra, you can first filter the library using +filter_spec() using SpectrumType.

+
matches <- match_spec(x = processed, library = test_lib,
+                   add_library_metadata = "sample_name", top_n = 5)
+print(matches[,c("object_id", "library_id", "match_val", "SpectrumType", "SpectrumIdentity")])
+#>    object_id                       library_id match_val SpectrumType
+#> 1: intensity 0031bb13faea1e04b52ffbeca009e8ab 0.9741244        Raman
+#> 2: intensity 0b10e6cc01c8e53f2221647a9c828164 0.4946028        Raman
+#> 3: intensity 101b6ae86864958ddd95445d0ab01fe4 0.3833523         FTIR
+#> 4: intensity 00002e4e3fac430aa1fdfea6e26f85e4 0.3225426         FTIR
+#> 5: intensity 0062d71f901713f2ac7a89c72fb186d4 0.3169888        Raman
+#>     SpectrumIdentity
+#> 1:              HDPE
+#> 2: Polyvinylchloride
+#> 3:           Nitrile
+#> 4:                PS
+#> 5:      polyisoprene
-
-

Selection Metadata

-

-

Whatever match is selected from the match table may have additional -metadata about it. That metadata will be displayed below the plot. Some -of this metadata may assist you in interpreting the spectra. For -example, if the spectra has metadata which says it is a liquid and you -are analyzing a solid particle, that spectrum may not be the best -match.

-

The R command for manual metadata selection using -sample_name == 5381 as example is:

-
find_spec(sample_name == 5381, library = spec_lib, which = "raman")
+
+

Library Metadata

+

The libraries we have created have over 100 variables of metadata in +them and this can be onerous to read through especially given that many +of the variables are NA values. We created +get_metadata() to remedy this by removing columns from the +metadata which are all blank values. The function below will return the +metadata for the top match in matches. Remember, similar +filter_spec(), you can specify logic for more +than one thing at a time.

+
get_metadata(x = test_lib, logic = matches[[1,"library_id"]], rm_empty = T)
+#>     x y                      sample_name                      SpectrumID
+#> 1: 37 1 0031bb13faea1e04b52ffbeca009e8ab HDPE_Sarah_#18_1s_20ac_10x_25mW
+#>      Organization SpectrumType SpectrumIdentity Polymer.Category LibraryType
+#> 1: J. Lynch, NIST        Raman             HDPE      polyolefins    Polymers
+#>                polymer_class plastic_or_not
+#> 1: Polyolefins (POLYALKENES)        plastic
+#>                                                   url_polymer_class
+#> 1: https://www.polymerdatabase.com/polymer%20index/polyalkenes.html
+#>           polymer                                                url_polymer
+#> 1: POLY(ETHYLENE) https://www.polymerdatabase.com/polymers/polyethylene.html
+#>                                                               url_more_info
+#> 1: https://www.polymerdatabase.com/polymer%20classes/Polyolefin%20type.html
+#>    snr_deriv                          file_id
+#> 1:  1339.684 77f837c1640910f2879184d61ef1df48
-
-

Match Plot

-

-

This plot is dynamically updated by selecting matches from the match -table. The red spectrum is the spectrum that you selected from the -reference library and the white spectrum is the spectrum that you are -trying to identify. Whenever a new dataset is uploaded, the plot and -data table in this tab will be updated. These plots can be saved as a -.png by clicking the camera button at the top of the plot.

+
+

Plot Matches

+

Overlaying unknown spectra and the best matches can be extremely +useful to identify peaks that don’t fit to the reference library which +may need further investigation. The example below shows great +correspondence between the best match and the unknown spectrum. All +major peaks are accounted for and the correct relative height. There are +two small peaks in the unknown spectrum near 500 that are not accounted +for which could be investigated further but we would call this a +positive id to HDPE.

+
plotly_spec(processed, filter_spec(test_lib, logic = matches[[1,"library_id"]]))
+
+
-
-

Sharing Reference Data

+
+

Sharing Reference Data

If you have reference data or AI models that you think would be useful for sharing with the spectroscopy community through OpenSpecy -please contact the website administrator to discuss options for +please contact the package administrator to discuss options for collaborating.

-
-

Characterizing Particles

-
-

Thresholding

-
-
-

Collapse Spectra

-
-
-
-

Additional App Features

-

Accessibility is extremely important to us and we are making strives -to improve the accessibility of Open Specy for all spectroscopists. -Please reach out if you have ideas for improvement.

-

We added a Google translate plugin to all pages in the app so that -you can easily translate the app. We know that not all languages will be -fully supported but we will continue to try and improve the -translations.

-
-

Download Data

-

-

After you have the preprocessing parameters set, we recommend that -you download the preprocessed data for your records. The download data -button will append the uploaded data to three columns created by the -preprocessing parameters. “Wavelength†and “Absorbance†are columns from -the data uploaded by the user. “NormalizedIntensity†is the max-min -normalized value (Equation 1) of the “Absorbanceâ€. “Smoothed†is the -Savitzky-Golay filter specified by the slider explained above. -“BaselineRemoved†is the smoothed and baseline corrected value that is -visible on the center plot.

-
+
+

Characterizing Particles

+

Sometimes the spectroscopy task we want to perform is to identify +particles in a spectral map. This is especially common for microplastic +analysis where a spectral map is used to image a sample and spectral +information is used to differentiate microplastic particles from +nonplastic particles. In addition to the material id, one often wants to +measure the shape and size of the particles. In a brute force technique, +one could first identify every spectrum in the map, then use +thresholding and image analysis to measure the particles. However, more +often than not, particles are well separated on the image surface and +background spectra is quite different from particle spectra and +therefore we can use thresholding a priori to identify and measure the +particles, then pass an exemplary spectrum for each particle to the +identification routine. It is important to note here that this is at the +bleeding edge of theory and technique so we may be updating these +functions in the near future.

+
+

Brute Force

+
data("test_lib")
+test_map <- read_any(read_extdata("CA_tiny_map.zip"))
+
+test_map_processed <- process_spec(test_map, conform_spec_args = list(range = test_lib$wavenumber, res = NULL))
+
+identities <- match_spec(test_map_processed, test_lib, add_library_metadata = "sample_name", top_n = 1)[order(factor(object_id, levels = colnames(test_map$spectra)))]#Should integrate this sorting into match_spec, tried but couldn't get it. 
+
+features <- ifelse(identities$match_val > 0.7, tolower(identities$polymer_class), "unknown")
+
+id_map <- def_features(x = test_map_processed, features = features)
+
+id_map$metadata$identities <- features #also should probably be implemented automatically in the function when a character value is provided. 
+
+test_collapsed <- collapse_spec(id_map) #Collapses spectra to their median for each particle. 
-
-

Conceptual diagram of Open Specy

-

+
+

A Priori Particle Thresholding

+
data("test_lib")
+test_map <- read_any(read_extdata("CA_tiny_map.zip"))
+
+snr <- sig_noise(test_map,metric = "log_tot_sig") #Characterize the total signal as a threshold value. 
+
+heatmap_spec(test_map, z = snr) #Use this to find your particles and the sig_noise value to use for thresholding. 
+
+id_map <- def_features(x = test_map, features = snr > 400) #Set define the feature regions based on the threshold. 400 appeared to be where I suspected my particle to be in the previous heatmap. 
+
+heatmap_spec(id_map, z = id_map$metadata$feature_id) #Check that the thresholding worked as expected. 
+
+collapsed_id_map <- id_map |>
+    collapse_spec() #Collapse the spectra to their medians based on the threshold. Important to note here that the particles with id -88 are anything from the FALSE values so they should be your background. 
+
+id_map_processed <- process_spec(
+collapsed_id_map, conform_spec_args = list(range = test_lib$wavenumber, res = NULL)) #Process the collapsed spectra. 
+
+plot(id_map_processed) #check the spectra. 
+
+matches <- match_spec(id_map_processed, test_lib, add_library_metadata = "sample_name", top_n = 1) #Get the matches of the collapsed spectra for the particles. 
diff --git a/vignettes/sop/deselection.jpg b/vignettes/sop/deselection.jpg deleted file mode 100644 index 7031c600..00000000 Binary files a/vignettes/sop/deselection.jpg and /dev/null differ diff --git a/vignettes/sop/download.jpg b/vignettes/sop/download.jpg deleted file mode 100644 index 11b72b34..00000000 Binary files a/vignettes/sop/download.jpg and /dev/null differ diff --git a/vignettes/sop/flowchart.png b/vignettes/sop/flowchart.png deleted file mode 100644 index 6b8dda43..00000000 Binary files a/vignettes/sop/flowchart.png and /dev/null differ diff --git a/vignettes/sop/hoverexample.jpg b/vignettes/sop/hoverexample.jpg deleted file mode 100644 index 1fd96b87..00000000 Binary files a/vignettes/sop/hoverexample.jpg and /dev/null differ diff --git a/vignettes/sop/matchplot.png b/vignettes/sop/matchplot.png deleted file mode 100644 index 0c8d13a4..00000000 Binary files a/vignettes/sop/matchplot.png and /dev/null differ diff --git a/vignettes/spectragryph.Rmd b/vignettes/spectragryph.Rmd index 83c50d6d..7c6a8742 100644 --- a/vignettes/spectragryph.Rmd +++ b/vignettes/spectragryph.Rmd @@ -1,5 +1,5 @@ --- -title: "File conversion in SpectraGryph" +title: "File Conversion in SpectraGryph" author: > Jessica Meyers, Jeremy Conkle, Win Cowger, Zacharias Steinmetz, Andrew Gray, Chelsea Rochman, Sebastian Primpke, Jennifer Lynch, @@ -7,7 +7,7 @@ author: > date: "`r Sys.Date()`" output: rmarkdown::html_vignette vignette: > - %\VignetteIndexEntry{File conversion in SpectraGryph} + %\VignetteIndexEntry{File Conversion in SpectraGryph} %\VignetteEncoding{UTF-8} %\VignetteEngine{knitr::rmarkdown} --- @@ -20,6 +20,9 @@ knitr::opts_chunk$set( ) ``` +[Spectragryph](https://www.effemm2.de/spectragryph/) supports many spectral file +conversions which facilitate data import to Open Specy. + 1. Download Spectragryph from [https://www.effemm2.de/spectragryph/down.html](https://www.effemm2.de/spectragryph/down.html) @@ -31,14 +34,10 @@ knitr::include_graphics("spectragryph/spectragryph-1.png") ``` 3. Click File, Save/export data, save data as, and save it as an spc -file. ¸ +file. ```{r, fig.align="center", out.width="98%", echo=FALSE} knitr::include_graphics("spectragryph/spectragryph-2.png") ``` 4. Then upload that .spc file to Open Specy. - -```{r, fig.align="center", echo=FALSE} -knitr::include_graphics("sop/uploadfile.jpg") -```