From 513dfbdcdda94a901b5901b486ff5500c7d158b1 Mon Sep 17 00:00:00 2001 From: ranke Date: Wed, 10 May 2006 15:44:14 +0000 Subject: The inverse prediction works in a variety of cases and is tested with Examples 7 and 8 from Massart! I need to compare with the DIN and draper examples, and finish the package vignette. git-svn-id: http://kriemhild.uft.uni-bremen.de/svn/chemCal@6 5fad18fb-23f0-0310-ab10-e59a3bee62b4 --- R/calplot.R | 52 ++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 52 insertions(+) create mode 100644 R/calplot.R (limited to 'R/calplot.R') diff --git a/R/calplot.R b/R/calplot.R new file mode 100644 index 0000000..cea1149 --- /dev/null +++ b/R/calplot.R @@ -0,0 +1,52 @@ +calplot <- function(object, xlim = "auto", ylim = "auto", + xlab = "Concentration", ylab = "Response", alpha=0.05) +{ + UseMethod("calplot") +} + +calplot.default <- function(object, xlim = "auto", ylim = "auto", + xlab = "Concentration", ylab = "Response", alpha=0.05) +{ + stop("Calibration plots only implemented for univariate lm objects.") +} + +calplot.lm <- function(object, xlim = "auto", ylim = "auto", + xlab = "Concentration", ylab = "Response", alpha=0.05) +{ + if (length(object$coef) > 2) + stop("More than one independent variable in your model - not implemented") + + if (alpha <= 0 | alpha >= 1) + stop("Alpha should be between 0 and 1 (exclusive)") + + m <- object + level <- 1 - alpha + x <- m$model$x + y <- m$model$y + newdata <- data.frame(x = seq(0,max(x),length=250)) + pred.lim <- predict(m, newdata, interval = "prediction",level=level) + conf.lim <- predict(m, newdata, interval = "confidence",level=level) + if (xlim == "auto") xlim = c(0,max(x)) + if (ylim == "auto") ylim = range(c(pred.lim,y)) + plot(1, + type = "n", + xlab = xlab, + ylab = ylab, + xlim = xlim, + ylim = ylim + ) + points(x,y, pch = 21, bg = "yellow") + matlines(newdata$x, pred.lim, lty = c(1, 4, 4), + col = c("black", "red", "red")) + matlines(newdata$x, conf.lim, lty = c(1, 3, 3), + col = c("black", "green4", "green4")) + + legend(min(x), + max(pred.lim, na.rm = TRUE), + legend = c("Fitted Line", "Confidence Bands", + "Prediction Bands"), + lty = c(1, 3, 4), + lwd = 2, + col = c("black", "green4", "red"), + horiz = FALSE, cex = 0.9, bg = "gray95") +} -- cgit v1.2.1