#' @name rsp.pls.plot
#' @title Plots for use with respeciate profile Positive Least Squares models
#' @aliases pls.plot pls_plot pls_plot_species pls_plot_profile
#' @description
#' The \code{pls_plot} functions are intended for use with PLS models built
#' using \code{rsp_pls_profile} (documented separately). They generate some
#' plots commonly used with source apportionment model outputs.
#' @param pls A \code{rsp_pls_profile} output, intended for use with
#' \code{pls_} functions.
#' @param id numeric or character
#' identifying the species or profile to plot. If numeric, these are treated
#' as indices of the species or profile, respectively, in the PLS model; if
#' character, species is treated as the name of species and profile is treated
#' as the profile code. Both can be concatenated to produce multiple plots and
#' the special case \code{id = -1} is a short cut to all species or profiles,
#' respectively.
#' @param plot.type numeric, the plot type if
#' multiple options are available.
#' @param ... other arguments, typically passed on to the associated
#' \code{lattice} plot.
#' @param output character, output method, one of: 'plot' to return just the
#' requested plot; 'data' to return just the data; and, c('plot', 'data') to
#' plot then return the data invisibly (default).
#' @param log (for \code{pls_plot_profile} only) logical, if \code{TRUE} this
#' applies 'log' scaling to the primary Y axes of the plot.
#' @return \code{pls_plot}s produce various plots commonly used in source
#' apportionment studies.
######################################################
# notes on recent changes
######################################################
# shifted my type to type to plot.type because it conflicts with type in lattice::xyplot....
# all plots use .rsp_plot_output
# so need to update them all if .rsp function formals change
# also pls_plot_profile uses rsp_plot_profile with output forced...
# for first plot layer
# CHECK: are we doing same for rsp_plot???
# all plots have old versions
# pls_plot... old
# can hopefully loose these at some point.
# to think about
#########################
# log axes - possible improvement
# NOTE: maybe log note here and do in rsp_plot or xxx.r ???
# need to check barplot handling of origin for this...
# zero handling is currently a based on a hard-coded offset in plot(..., log="y", off.set)
# but like it automatically estimated...
# should check log axes handling when only 1 or 2 species in profile
# I think it is suspect
# also like a tidier handling of spacing when plotting log axes
# GENERAL NOTES
# TO DO
# these all need code tidying
# to think about
###################################
# common approaches for ordering subsetting and renaming
# for both species and profiles when plotting.
# also like to make this and rsp_plot consistent
# check individual function notes
# plot.type=2 for pls_plot pls_plot_species or pls_plot_profile
# the standard per species summary:
# rsp_x species amount versus sample n#/index axes,
# black line of actual .actual spcies
# stacked bar plot of prediction per [profile]
# (like what CRC plots but automated plot)
# plot.type=3 summary as pie chart ???
####################################
###################################
## pls_plots
###################################
###################################
#these are all draft
###################################
###################################
## pls_plot
###################################
###################################
#' @rdname rsp.pls.plot
#' @export
## this replaces previous pls_plot (now pls_plot.old)
## now imports via data.table::
## need this to kill as.data.table load warning message
## #' @import data.table
#test
#devtools::load_all()
#d1 <- readRDS("C:\\Users\\trakradmin\\OneDrive - University of Leeds\\Documents\\pkg\\respeciate\\test\\my.working.rds")
#ref <- rsp(c("4868", "4914", "8948", "91155", "91163", "95441", "95529"))
#mod <- rsp_pls_profile(d1, ref, power=2)
#pls_plot(mod)
pls_plot <- function (pls, plot.type = 1, ...,
output="default"){
#current using lattice/latticeExtra for the panelling/layers...
########################
# to do
########################
# id
#basic plots in development...
# plot element ordering
# currently as it comes...
# because *I think* pls_report kills all pre-model handling...
# no id handling
# that maybe need to be in plot type...
# maybe also want to do it at end
# so missing case locations are retained (if needed for plot)...
############################
# nags
############################
# pls_plot(..., horizontal=FALSE) errors
# should flip x and y...
# type = 1
############################
#in development
#pls fit summary
# simple proportional fit plot
#think like
# https://latticeextra.r-forge.r-project.org/man/postdoc.html
# (but without the 100 percent (proportion=1) limit...)
# type = 2
############################
#to do?
#################
#setup
#################
.x.args <- list(...)
dat <- pls_report(pls)
.ord.pro.c <- .rsp_profile_code_order(dat)
######################################
#option to not do name simplification?
#
dat$.species <- .rsp_tidy_species_name(dat$.species)
.sp.ref <- unique(dat$.species)
#type
if(!plot.type %in% c(1)){
stop("pls_plot: plot.type unknown, check ?pls_plot",
call. = FALSE)
}
############################
#type 1
############################
if(plot.type==1){
##################################
# note
# maybe use .rsp_get_prop_from_pls
# but check naming change...???
##################################
.tmp <- names(dat)
.tmp <- .tmp[grep("^.x_", .tmp)]
.refs <- c(.tmp, "pred")
#make summary pls. prop.table
.ans <- lapply(.sp.ref, function(x){
.tmp <- subset(dat, .species==x)
.d2 <- .tmp[1, c(".species", .refs)]
for(.ref in .refs){
#use only paired cases to calculate skew...
.tmp2 <- .tmp[c(.ref, ".value")]
.tmp2[.tmp2==0] <- NA
.tmp2 <- na.omit(.tmp2)
.d2[, .ref] <- sum(.tmp2[,.ref], na.rm=TRUE) / sum(.tmp2[,".value"], na.rm=TRUE)
}
.d2
})
.ans <- do.call(rbind, .ans)
#barchart formula
.for <- paste(paste("`", .tmp, "`", sep=""), collapse="+")
.for <- as.formula(paste(".species~", .for, sep=""))
.tmp <- gsub("^.x_", "", .tmp)
#plot lists
gr.ls <- list(h=-1, v=-1, col="grey", lty=3)
if("grid" %in% names(.x.args) && is.list(.x.args$grid)){
gr.ls <- modifyList(gr.ls, .x.args$grid)
.x.args$grid <- NULL
}
bar.ls <- list(v=c(0.5,1,2), lty=c(3,2,3), col="red")
if("bars" %in% names(.x.args) && is.list(.x.args$bars)){
bar.ls <- modifyList(bar.ls, .x.args$bars)
.x.args$bars <- NULL
}
if("col" %in% names(.x.args)){
#could allow function as col input???
.cols <- rep(.x.args$col, length.out=length(.tmp))
} else {
.cols <- rainbow(length(.tmp))
.x.args$col <- .cols
}
ky.ls <- list(space="right", text=list(text=.tmp),
rect=list(col=.cols))
if("key" %in% names(.x.args) && is.list(.x.args$key)){
if(any(c("x", "y") %in% names(.x.args$key))){
ky.ls$space <- NULL
}
.x.args$key <- modifyList(ky.ls, .x.args$key)
} else {
.x.args$key <- ky.ls
}
#####################################
#note
#maybe use .rsp_panelpal for grid and bars
# see pls_plot_profile ???
#####################################
pl.ls <- list(x=.for, data=.ans, origin=0, stack=TRUE,
grid=TRUE, bars=TRUE, xlim=c(-0.025, NA),
xlab="mean(model) / mean(measurements)",
#prepanel=function(...){
# ans <- lattice::prepanel.default.bwplot(...)
# print(ans)
# ans
#},
panel=function(...){
.temp <- list(...)
if(.temp$grid){
do.call(panel.grid, gr.ls)
}
panel.barchart(...)
if(.temp$bars){
do.call(panel.abline, bar.ls)
}
}
)
pl.ls <- modifyList(pl.ls, .x.args)
p <- do.call(barchart, pl.ls)
}
#output
############################
.rsp_plot_output(as.data.frame(.ans), pl.ls, p, output)
}
####################################
####################################
## pls_plot_profile
####################################
####################################
#' @rdname rsp.pls.plot
#' @export
## now imports from xxx.r
## #' @import data.table
#############################
#this needs a lot of work
#############################
#test
#devtools::load_all()
#d1 <- readRDS("C:\\Users\\trakradmin\\OneDrive - University of Leeds\\Documents\\pkg\\respeciate\\test\\my.working.rds")
#ref <- rsp(c("4868", "4914", "8948", "91155", "91163", "95441", "95529"))
#mod <- rsp_pls_profile(d1, ref, power=2)
#pls_plot_profile(mod)
# log scale may need work
# but that is in rsp_plot_profile/unexported functions...
pls_plot_profile <- function (pls, plot.type=1, log = FALSE, ...,
output="default")
{
#new version of pls_plot
#to do
##############################
# log (needs better axes control) but that is in rsp_profile_plot
# id needs to be enabled...
#setup
#############################
.x.args <- list(...)
dat <- pls_report(pls)
.ord.pro.c <- .rsp_profile_code_order(dat)
######################################
#option to not do name simplification?
#
dat$.species <- .rsp_tidy_species_name(dat$.species)
.sp.ref <- unique(dat$.species)
#plot.type control
if(!plot.type %in% c(1)){
stop("pls_plot_profile: plot.type unknown, check ?pls_plot_profile",
call. = FALSE)
}
############################
#type 1
############################
if(plot.type==1){
#make first plot and output .ans
##############################
#get profiles .m_ columns
.ans <- .rsp_get_m_from_pls(dat)
.p1.prof <- unique(.ans$.profile.id)
#send to rsp_plot_profile with any user arguments
# to make first plot
#set cols
#set cols
p1.ls <- list(rsp=.ans, layout =c(1, length(.p1.prof)), log=log,
multi.profile = "panel", id=1:length(.p1.prof),
order=FALSE, silent=TRUE)
if(!"col" %in% names(p1.ls)){
#maybe need better handling
p1.ls$col <- trellis.par.get("superpose.symbol")$col[1]
}
#issue with species_code not being known made this tricky...
p1.ls <- modifyList(p1.ls, .x.args)
p1.ls$output <- "plot"
p1 <- do.call(rsp_plot_profile, p1.ls)
#make second plot and .ans2
######################################
.ans2 <- .rsp_get_prop_from_pls(dat)
.ans2$.pc <- .ans2$.prop * 100
#could do this in the panel so any missing is greyed out ???
.ans2$.pc[is.na(.ans2$.pc)] <- 0
p2.ls <- .rsp_panelPal("tc", list(x =.pc~factor(.species)|factor(.profile.id),
data=.ans2,
type=c("h", "p"), pch=18, layout=c(1,7),
ylab="Total Contribution (%)",
scales=list(x=list(rot=90))),
#note: function is cheat to use .rsp... outside lattice
# could make it the default if no panel set in call???
function(...){list(...)}, ...)
#if tc layer not turned off..
if(!is.null(p2.ls)){
if(!"col" %in% names(p2.ls)){
#maybe need better handling
p2.ls$col <- trellis.par.get("superpose.symbol")$col[2]
}
p2 <- do.call(xyplot, p2.ls)
p1 <- update(doubleYScale(p1, p2, add.ylab2 = TRUE),
par.settings = simpleTheme(col=c(p1.ls$col[1], p2.ls$col[1]))
)
}
}
#output
############################
#list is void here...
.rsp_plot_output(list(profile = .ans, tc = .ans2), list(), p1, output)
}
####################################
####################################
## pls_plot_species
####################################
####################################
#' @rdname rsp.pls.plot
#' @export
## now imports from xxx.r
## #' @import data.table
#############################
#this needs a lot of work
#############################
#test
#devtools::load_all()
#d1 <- readRDS("C:\\Users\\trakradmin\\OneDrive - University of Leeds\\Documents\\pkg\\respeciate\\test\\my.working.rds")
#ref <- rsp(c("4868", "4914", "8948", "91155", "91163", "95441", "95529"))
#mod <- rsp_pls_profile(d1, ref, power=2)
#pls_plot_species(mod)
# id enabled but
# species order is always as supplied...
# probably actually alphabetic
# look like order(character(unique(.profile.id)))
# to do
# limit default output to < 7 plots?
# .rsp_panelpal handling like other plots
# layer .mod ???
# log ??? (not sure it is needed/useful)
# decide how to reorder or rename species, profiles and x data
# (do this in plots and data ???)
# decide how to modify .index
pls_plot_species <- function (pls, id, plot.type=1, ...,
output = "default")
{
#new version of pls_plot
#to do
##############################
# most stuff
# log not sure about doing them...
# id
#setup
#############################
.x.args <- list(...)
dat <- pls_report(pls)
.ord.pro.c <- .rsp_profile_code_order(dat)
.sp.ref <- unique(dat$.species)
my.species <- if (missing(id)) {
.sp.ref
#default option (print the lot...)
############################
#possibly a warning if lots of species to plot
##################
} else {
id
}
if (is.numeric(my.species)) {
if (all(my.species == -1)) {
my.species <- .sp.ref
}
else {
my.species <- .sp.ref[my.species]
}
}
if(length(my.species)>6 & missing(id)){
#to think about
# option to turn off warning???
# (using in older versions of code)
#if(!silent){
warning("RSP/PLS> ", length(my.species), " species... ",
"just showing first 6 to reduce plot clutter",
"\n\t (maybe use id to force larger range if sure)",
sep="", call.=FALSE)
#}
my.species <- my.species[1:6]
}
######################################
#option to not do name simplification?
#
dat$.species <- .rsp_tidy_species_name(dat$.species)
.sp.ref <- unique(dat$.species)
my.species <- .rsp_tidy_species_name(my.species)
if (!any(my.species %in% .sp.ref)) {
stop("pls_plot_species> unknown species, please check", call. = FALSE)
}
.tmp <- dat[c(".species", ".profile.id", ".value")]
.tmp <- data.table::as.data.table(.tmp)
.tmp <- data.table::dcast(.tmp, .profile.id ~ .species,
mean,
na.rm=TRUE,
value.var = ".value")
.tmp2 <- data.table::melt(.tmp, id.vars=".profile.id", variable.name=".species",
value.name=".value")
.tmp <- dat[c(".species", ".profile.id", "pred")]
.tmp <- data.table::as.data.table(.tmp)
.tmp <- data.table::dcast(.tmp, .profile.id ~ .species,
mean,
na.rm=TRUE,
value.var = "pred")
.tmp <- data.table::melt(.tmp, id.vars=".profile.id", variable.name=".species",
value.name="pred")
.tmp <- data.table::merge.data.table(.tmp2, .tmp)
.tmp <- as.data.frame(.tmp)
.tmp$.index <- as.numeric(factor(.tmp$.profile.id, levels=.ord.pro.c,
ordered = TRUE))
.tmp<- .tmp[order(.tmp$.index),]
#plot.type control
if(!plot.type %in% c(1,2)){
stop("pls_plot_species: plot.type unknown, check ?pls_plot_profile",
call. = FALSE)
}
.tmp <- subset(.tmp, .species %in% my.species)
############################
#type 1
############################
if(plot.type==1){
.mc <- if ("mod.col" %in% names(.x.args)) {
.x.args$mod.col
} else {
"red"
}
plt <- list(x=pred~.value | .species, data=.tmp,
#prepanel forces x and y lims to same range
prepanel=function(...){
.tmp <- prepanel.default.xyplot(...)
.tmp$xlim <- range(c(.tmp$xlim, .tmp$ylim))
.tmp$ylim <- .tmp$xlim
.tmp
},
panel= function(x, y, xlim, ylim, ...){
#user control of grid - like loa...
.rsp_panelPal("grid", list(h=-1,v=-1, col="grey", lty=3),
lattice::panel.grid, ...)
#TO DO
#user control of y=x
panel.ablineq(a = 0, b = 1, adj = c(0,1),
col.line ="grey", lty=2, label="")
#user control of main plotted data via
# standard lattice
panel.xyplot(x=x, y=y, ...)
#CURRENTLY JUST col via mod.col
#user control of model
panel.ablineq(lm(y ~ x + 0), cex = 0.8,
x = min(c(x, y), na.rm=TRUE),
y = max(c(x, y), na.rm=TRUE),
r.squared = TRUE, adj = c(0,1),
sep = " (", sep.end = ")",
offset=0, varStyle = NULL,
col.line = .mc, col.text = .mc, digits = 2)
},
xlab="Measurement", ylab="model",
scales=list(y=list(relation="free",
rot=90),
x=list(relation="free")))
plt <- modifyList(plt, .x.args)
p <- do.call(xyplot, plt)
}
if(plot.type==2){
#xlab
if(!"xlab" %in% names(.x.args)){
.x.args$xlab <- "Sample [index]"
}
if(!"ylab" %in% names(.x.args)){
.x.args$ylab <- "Measurement, Model"
} else {
if(length(.x.args$ylab)>1){
if(!"key.text" %in% names(.x.args)){
.x.args$key.text <- .x.args$ylab[1:2]
}
.x.args$ylab <- paste(.x.args$ylab[1], .x.args$ylab[2], sep =", ")
}
}
if(!"key.text" %in% names(.x.args)){
.x.args$key.text <- c("Measurement", "Model")
}
if(!"col" %in% names(.x.args)){
.x.args$col <- trellis.par.get("superpose.line")$col[1:2]
}
if("mod.col" %in% names(.x.args)){
.x.args$col <- c(.x.args$col[1], .x.args$mod.col)
}
p2.ls <- list(x= .value + pred ~ .index | .species, data=.tmp,
auto.key = list(text=.x.args$key.text,
space="top", columns=2),
type="l",
panel= function(...){
.rsp_panelPal("grid", list(h=-1,v=-1, col="grey", lty=3),
panel.grid, ...)
lattice::panel.xyplot(...)
},
scale=list(relation="free"),
par.settings = simpleTheme(col=.x.args$col))
p2.ls <- modifyList(p2.ls, .x.args)
p <- do.call(xyplot, p2.ls)
##plot(p)
}
#output
############################
#this needs working up based on input from Dennis...
## plot(p)
## return(invisible(.tmp))
.rsp_plot_output(as.data.frame(.tmp), list(p2.ls=plt, p2.ls=p2.ls), p, output)
}
######################################
######################################
## unexported
######################################
######################################
#old
# holding until new versions are fully up and running...
## pls_plot.old
############################
# currently not exporting...
############################
# #' @rdname rsp.pls.plot
# #' @export
## now imports via data.table::
## need this to kill the as.data.table load message
## #' @import data.table
##
#############################
#this needs a lot of work
#############################
# this uses unexported rsp_profile_pie function below...
# both pls_plot and rsp_profile_pie need work...
pls_plot.old <- function (pls, n, plot.type = 1, ...){
#current using lattice/latticeExtra for the panelling/layers...
#basic plots finished but...
# currently not passing arguments generally
# the par setting seem to be dropped when using plot(p)
# ahead of end of function
############################
# nags
############################
# type = 1
############################
# note sure about the layer naming
# zero is not bottom of barchart...
# type = 2
############################
# the label positioning is messy (see not about nudge)
# cex setting too small if only one panel...
#wondering about
# https://latticeextra.r-forge.r-project.org/man/postdoc.html
# as an alternative to type=2
# (but 100 percent measured rather than proportion...)
#################
#setup
#################
.x.args <- list(...)
dat <- pls_report(pls)
.ord.pro.c <- .rsp_profile_code_order(dat)
#dat$.species <- .rsp_tidy_species_name(dat$.species)
.sp.ref <- unique(dat$.species)
#species
# now defaulting to all plots
species <- if (missing(n)) {
species <- .sp.ref
}
else {
n
}
if (is.numeric(species)) {
if (all(species == -1)) {
species <- .sp.ref
}
else {
species <- .sp.ref[species]
}
}
if (!any(species %in% .sp.ref)) {
stop("RSP_PLS> unknown species, please check", call. = FALSE)
}
################################
#note:
# could condition here BUT currently
# holding on to everything until just before plot
# might not need to do this....
#################################
.sp.ord <- unique(dat$.species)
#####################################
#messy at moment...
.sp.m.pro <- names(dat)[grep("^.n_", names(dat))]
.sp.pro <- gsub("^.n_", "", .sp.m.pro)
#line col....
.col <- lattice::trellis.par.get("superpose.line")$col[1]
#bar cols
.cols <- if ("col" %in% names(.x.args)) {
#could include if you supply a function..?
#could use col.regions?
.cols <- .x.args$col
}
else {
.cols <- heat.colors(n = length(.sp.m.pro))
}
if (length(.cols) != length(.sp.m.pro)) {
stop("pls_plot> halted; expecting ", length(.sp.m.pro),
"colours; given ", length(.cols), sep = "", call. = FALSE)
}
######################
# build x_[profile]
######################
#now done in pls_report
for (i in .sp.pro) {
dat[, paste(".x_", i, sep = "")] <- dat[, paste(".m_", i, sep = "")] *
dat[, paste(".n_", i, sep = "")]
}
.sp.x.pro <- names(dat)[grep("^.x_", names(dat))]
.rep <- dat[c(".species", ".species.id", ".profile.id",
.sp.x.pro)]
.rep <- data.table::melt(data.table::as.data.table(.rep),
id = c(".species.id", ".species", ".profile.id"))
.tot <- data.table::as.data.table(dat)
.cs <- c(".value", "pred", .sp.x.pro)
.tot <- .tot[, lapply(.SD, function(x) sum(x, na.rm = TRUE)),
.SDcols = .cs, by = c(".species.id", ".species")]
###########################
# now plotting as panels
# using
###########################
######################################################
# now using rsp_ function to track all pls model cases
# previous method only tracked valid cases for the plotted data
# so no gaps where models dropped/not built...
#########################################################
.rep$.index <- as.numeric(factor(.rep$.profile.id, levels = .ord.pro.c,
ordered = TRUE))
dat$.index <- as.numeric(factor(dat$.profile.id, levels = .ord.pro.c,
ordered = TRUE))
.tmp <- dat[c(".species.id", ".profile.id", ".index", ".value", "pred")]
.rep <- data.table::merge.data.table(.rep, .tmp)
.rep$variable <- gsub("^x_", "", .rep$variable)
#print(names(.rep))
#return(dat)
.rep <- subset(as.data.frame(.rep), .species %in% species)
if (1 %in% plot.type) {
#lattice sets panel order based
.sp <- if(is.factor(.rep$.species)){
levels(.rep$.species)
} else {
sort(unique(.rep$.species))
}
.sp <- .sp[.sp %in% .rep$.species]
#.y.scale <- lapply(unique(.rep$.species), function(x){
.y.scale <- lapply(.sp, function(x){
.tmp <- .rep[.rep$.species==x,]
c(0, max(c(.tmp$.value, .tmp$pred), na.rm=TRUE))
})
###############################################
#use loa method to generalise this?
###############################################
p2 <- lattice::xyplot(.value ~ .index | .species, .rep,
panel=lattice::panel.xyplot,
type="s", xlab="Sample [index]",
ylab="Measurement",
scales=list(relation="free"),
ylim=.y.scale)
p <- lattice::barchart(value ~ factor(.index) | .species, .rep,
groups=.rep$variable, stack=TRUE,
panel=function(x, y, col, groups, ..., subscripts){
#grid control like loa
.rsp_panelPal("grid", list(h=-1,v=-1, col="grey", lty=3),
lattice::panel.grid, ...)
lattice::panel.barchart(x=x, y=y, col=col,
groups=groups,
subscripts=subscripts, ...)
.y <- .rep$.value[subscripts]
#col needs to be from option[1]
lattice::panel.xyplot(x=x, y=.y,
col=.col,
type="l",
subscripts=subscripts,...)
},
scales=list(relation="free"),
#auto.key=list(space="top", columns=2,
# col.line=.cols,
# points=FALSE, rectangles=TRUE),
ylim=.y.scale,
col=.cols,
border=NA,
#par.settings = list(superpose.polygon = list(col = .cols,
# pch =c (15, 15)),
#superpose.symbol = list(fill = .cols)),
auto.key=list(space="top", columns = 3,
cex=0.8,
points=FALSE,
rectangles=TRUE)) #,
#xscale.components = function(lim,...){
# lim <- as.numeric(as.character(lim))
# ans <- lattice::xscale.components.default(lim=lim,...)
# print(ans)
# ans
#})
plot(update(latticeExtra::doubleYScale(p2, p, add.axis = FALSE),
par.settings = list(superpose.polygon = list(col = .cols),
superpose.symbol = list(fill = .cols))))
#p2 <- lattice::xyplot(.value ~ factor(.index) | .species, dat,
# type="l", scales=list(relation="free"))
#plot(cheat(p, latticeExtra::as.layer(p2)))
#plot(latticeExtra::doubleYScale(p, p2, add.axis=FALSE, add.ylab2=FALSE))
}
if (2 %in% plot.type) {
p <- lattice::xyplot(value ~ .index | .species, .rep,
groups=.rep$variable,
totals=.rep$.value,
scales=list(relation="free",
draw=FALSE),
ylab="", xlab="",
col = .cols,
auto.key=list(space="top", columns = 3,
cex=0.8,
points=FALSE,
rectangles=TRUE),
ylim=c(-2,2), xlim=c(-2,2),
between = list(x = 0.2, y = 0.2),
panel=.rsp_panel.pie,
par.settings = list(superpose.polygon = list(col = .cols),
axis.line = list(col = 'transparent'),
superpose.symbol = list(fill = .cols))
)
plot(p)
}
invisible(.rep)
}
########################
#currently not exporting
########################
pls_plot_profile.old <- function (pls, n, log = FALSE, ...)
{
#########################
#previous plot used base r graphics
#this moved to lattice/latticeExtra
#so we can panel outputs
#########################
#setup
.x.args <- list(...)
.plt.args <- .x.args[names(.x.args %in% c())]
dat <- pls_report(pls)
.sp.ord <- unique(dat$.species.id)
.sp.m.pro <- names(dat)[grep("^.m_", names(dat))]
.sp.pro <- gsub("^.m_", "", .sp.m.pro)
#defaulting n to all profiles as one plot
profile <- if (missing(n)) {
profile <- .sp.pro
} else {
n
}
if (is.numeric(profile)) {
if (all(profile == -1)) {
profile <- .sp.pro
}
else {
profile <- .sp.pro[profile]
}
}
if (!any(profile %in% .sp.pro)) {
stop("RSP_PLS> unknown profile(s), please check", call. = FALSE)
}
#########################
#build x_[profile]
#########################
n_profile <- paste(".n_", profile, sep = "")
m_profile <- paste(".m_", profile, sep = "")
dat <- dat[c(".species.id", ".species", ".profile.id",
n_profile, m_profile, "pred", ".value")]
for (i in profile) {
dat[, paste(".x_", i, sep = "")] <- dat[, paste(".m_", i, sep = "")] *
dat[, paste(".n_", i, sep = "")]
}
.rep <- data.table::as.data.table(dat)
.cols <- c(".value", "pred", paste(".x_", profile, sep = ""))
.rep <- .rep[, lapply(.SD, function(x) sum(x, na.rm = TRUE)),
.SDcols = .cols, by = c(".species.id", ".species")]
.rep <- as.data.frame(.rep)
#########################
# y2 setup
#########################
# by default this is .value
# but might want mod prediction
if ("y2" %in% names(.x.args) && .x.args$y2 == "pred") {
for (i in profile) {
.rep[, paste("pc_", i, sep = "")] <-
(.rep[, paste(".x_", i, sep = "")]/.rep$pred) * 100
}
}
else {
for (i in profile) {
.rep[, paste("pc_", i, sep = "")] <-
(.rep[, paste(".x_", i, sep = "")]/.rep$.value) * 100
}
}
#might not need all of following now we
#we are not pulling apart to plot one at time...
dat <- dat[!duplicated(dat$.species), ]
dat$.profile <- dat$.profile[1]
dat$.profile.id <- dat$.profile.id[1]
dat <- merge(.rep, dat[c(".species.id", ".species", ".profile.id",
m_profile)], )
dat <- dat[order(ordered(dat$.species.id, levels = .sp.ord)), ]
################################
# build pc_[profile]
################################
rownames(dat) <- 1:nrow(dat)
.ref <- names(dat)[grep("pc_", names(dat))]
.oth <- c(".species.id", ".species", ".profile.id", ".value", "pred")
.temp <- data.table::as.data.table(dat[c(.oth, gsub("^pc_", ".x_", .ref))])
.d1 <- data.table::melt(.temp, measure.vars = gsub("^pc_", ".x_", .ref),
variable.name = "pls_profile", value.name = "loading")
.temp <- data.table::as.data.table(dat[c(.oth, .ref)])
.d2 <- data.table::melt(.temp, measure.vars = .ref,
variable.name = "pls_profile", value.name = "percent_contr")
.d2$pls_profile <- gsub("^pc_", ".x_", .d2$pls_profile)
dat <- as.data.frame(merge(.d1, .d2, all=T))
#############################
############################
#now using lattice to handle logs
###############
#.dat <- dat
#don't need local version of dat because not changing data ahead of plot
#if(log){
# .dat$loading <- log10(.dat$loading)
# .ylim <- lapply(profile, function(x){
# .temp <- subset(.dat, pls_profile==x)
# .temp <- range(.temp$loading, na.rm=TRUE, finite=TRUE)
# if(.temp[1] == .temp[2]){
# .temp <- c(.temp[1]-1, .temp[1]+1)
# }
# range(c(floor(.temp), ceiling(.temp)))
# })
#} else {
# .ylim <- lapply(profile, function(x){
# .temp <- subset(.dat, pls_profile==x)
# .temp <- range(.temp$loading, na.rm=TRUE, finite=TRUE)
# range(pretty(.temp))
# })
#}
######################
#plot
######################
#now using lattice/latticeExtra
##
#think there is more here that can be generalized...
p1.ls <- list(x = loading~.species | pls_profile,
data=dat, ylab="Source Loading",
panel = function(...){
.rsp_panelPal("grid", list(h=-1,v=-1, col="grey", lty=3),
panel.grid, ...)
panel.barchart(...)
},
between=list(y=.2),
scales=list(x=list(rot=90),
y=list(rot=c(0,90),
relation="free")),
layout=c(1,length(profile)))
if(log){
p1.ls$scales$y$log=10
p1.ls$yscale.components <- .rsp_yscale.component.log10
}
p1.ls <- modifyList(p1.ls, .x.args)
if(!"col" %in% names(p1.ls)){
p1.ls$col <- trellis.par.get("superpose.line")$col[1]
}
p1 <- do.call(barchart, p1.ls)
if("mod" %in% names(.x.args) && !.x.args$mod){
#if mod FALSE just plot 1
plot(p1)
} else {
#add mod layer (total contributions) as y2
.col2 <- if("mod.col" %in% names(.x.args)){
.x.args$mod.col
} else {
trellis.par.get("superpose.line")$col[2]
}
p2.ls <- list(x = percent_contr ~ factor(.species) | pls_profile,
pch=16, type=c("h", "p"), col= c(.col2, .col2),
ylab="Total Contribution (%)",
data=dat)
.tmp <- .x.args[grepl("^mod[.]", names(.x.args))]
if(length(.tmp)>0){
names(.tmp) <- gsub("^mod[.]", "", names(.tmp))
p2.ls <- modifyList(p2.ls, .tmp)
}
p2 <- do.call(xyplot, p2.ls)
plot(update(doubleYScale(p1, p2, add.ylab2 = TRUE),
par.settings = simpleTheme(col=c(p1.ls$col[1], .col2))))
}
############
#output
############
#could pass plot and data as list???
return(invisible(dat))
}
## #' @rdname rsp.pls.plot
## #' @export
## now imports from xxx.r
## #' @import data.table
#############################
#this needs a lot of work
#############################
#currently not exporting
pls_plot_species.old <- function (pls, id, plot.type = 1, ...)
{
#not including NAs....
###########################
# setup
###########################
.x.args <- list(...)
dat <- pls_report(pls)
.ord.pro.c <- .rsp_profile_code_order(dat)
.sp.ref <- unique(dat$.species)
species <- if (missing(id)) {
.sp.ref
#default option (print the lot...)
############################
#possibly a warning if lots of species to plot
##################
} else {
id
}
if (is.numeric(species)) {
if (all(species == -1)) {
species <- .sp.ref
}
else {
species <- .sp.ref[species]
}
}
if (!any(species %in% .sp.ref)) {
stop("RSP_PLS> unknown species, please check", call. = FALSE)
}
############################
#if not earlier, then here?
#possibly a warning if lots of species to plot
##################
#########################
#could drop a lot of this??
#########################
.xlb <- if ("xlab" %in% names(.x.args)) {
.x.args$xlab
} else {
"Measurement"
}
.ylb <- if ("ylab" %in% names(.x.args)) {
.x.args$ylab
} else {
"Model"
}
.bc <- if ("col" %in% names(.x.args)) {
.x.args$col
} else {
par("col")
}
.mc <- if ("mod.col" %in% names(.x.args)) {
.x.args$mod.col
} else {
"red"
}
dat <- subset(dat, .species %in% species)
# lims <- range(c(d2$.value, d2$pred), na.rm = TRUE, finite = TRUE)
# mod <- lm(pred ~ 0 + .value, d2)
# .sum <- paste("y = ", signif(summary(mod)$coefficients[1,
# 1], 3), "x (adj.R2 = ", signif(summary(mod)$adj.r.squared,
# 3), ")", sep = "")
.lims <- lapply(species, function(x){
.d <- subset(dat, .species==x)
range(c(.d$pred, .d$.value), finite=TRUE, na.rm=TRUE)
})
if (1 %in% plot.type) {
.mc <- if ("mod.col" %in% names(.x.args)) {
.x.args$mod.col
} else {
"red"
}
p1.ls <- list(x=pred~.value | .species, data=dat,
#prepanel forces x and y lims to same range
prepanel=function(...){
.tmp <- prepanel.default.xyplot(...)
.tmp$xlim <- range(c(.tmp$xlim, .tmp$ylim))
.tmp$ylim <- .tmp$xlim
.tmp
},
panel= function(x, y, xlim, ylim, ...){
#user control of grid - like loa...
.rsp_panelPal("grid", list(h=-1,v=-1, col="grey", lty=3),
lattice::panel.grid, ...)
#TO DO
#user control of y=x
panel.ablineq(a = 0, b = 1, adj = c(0,1),
col.line ="grey", lty=2, label="")
#user control of main plotted data via
# standard lattice
panel.xyplot(x=x, y=y, ...)
#CURRENTLY JUST col via mod.col
#user control of model
panel.ablineq(lm(y ~ x + 0), cex = 0.8,
x = min(c(x, y), na.rm=TRUE),
y = max(c(x, y), na.rm=TRUE),
r.squared = TRUE, adj = c(0,1),
sep = " (", sep.end = ")",
offset=0, varStyle = NULL,
col.line = .mc, col.text = .mc, digits = 2)
},
xlab="Measurement", ylab="model",
scales=list(y=list(relation="free",
rot=90),
x=list(relation="free")))
p1.ls <- modifyList(p1.ls, .x.args)
p <- do.call(xyplot, p1.ls)
plot(p)
# plot(d2$.value, d2$pred, type = "n", main = i, col = .bc,
# xlab = .xlb, ylab = .ylb, xlim = lims, ylim = lims)
# grid()
# abline(a = 0, b = 1, col = "grey")
# points(d2$.value, d2$pred)
# abline(mod, col = .mc, lty = 2)
# text(lims[1], lims[2], .sum, adj = c(0, 1), cex = 0.75)
}
if (2 %in% plot.type) {
#xlab
if(!"xlab" %in% names(.x.args)){
.x.args$xlab <- "Sample [index]"
}
if(!"ylab" %in% names(.x.args)){
.x.args$ylab <- "Measurement, Model"
} else {
if(length(.x.args$ylab)>1){
if(!"key.text" %in% names(.x.args)){
.x.args$key.text <- .x.args$ylab[1:2]
}
.x.args$ylab <- paste(.x.args$ylab[1], .x.args$ylab[2], sep =", ")
}
}
if(!"key.text" %in% names(.x.args)){
.x.args$key.text <- c("Measurement", "Model")
}
if(!"col" %in% names(.x.args)){
.x.args$col <- trellis.par.get("superpose.line")$col[1:2]
}
if("mod.col" %in% names(.x.args)){
.x.args$col <- c(.x.args$col[1], .x.args$mod.col)
}
#ylab
#can to two terms for
#if("ylab" %in% names(.x.args)){
# if(length(.x.args$ylab)>1){
# if(!"key.text" %in% names(.x.args)){
# .x.args$key.text <- .x.args$ylab[1:2]
# }
# .x.args$ylab <- paste(.x.args$ylab[1], .x.args$ylab[2], sep =", ")
# } else {
# if(!"key.text" %in% names(.x.args)){
# .x.args$key.text <- c("Measurement", "Model")
# }
# }
#} else {
# if(!"key.text" %in% names(.x.args)){
# .x.args$key.text <- c("Measurement", "Model")
# }
# .x.args$ylab <- "Measurement, Model"
#}
#########################
#previous code
#########################
#plot(d2$.value, type = "n", main = i, col = .bc,
# ylab = .ylb, xlab = .xlb, ylim = lims)
#lines(d2$.value)
#lines(d2$pred, col = .mc)
########################
#using standardised index
#make 'ordered profile codes' at top
# before any subsetting...
# .ord.pro.c <- rsp_profile_code_order(dat)
dat$.index <- as.numeric(factor(dat$.profile.id, levels=.ord.pro.c,
ordered = TRUE))
dat<- dat[order(dat$.index),]
p2.ls <- list(x= .value + pred ~ .index | .species, data=dat,
auto.key = list(text=.x.args$key.text,
space="top", columns=2),
type="l",
panel= function(...){
.rsp_panelPal("grid", list(h=-1,v=-1, col="grey", lty=3),
panel.grid, ...)
lattice::panel.xyplot(...)
},
scale=list(relation="free"),
par.settings = simpleTheme(col=.x.args$col))
p2.ls <- modifyList(p2.ls, .x.args)
p <- do.call(xyplot, p2.ls)
plot(p)
######################
# or any with any missing are plot on different x scale
# maybe find longest, take range for that
#xyplot(.value + .pred ~ )
}
invisible(dat)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.