diff options
| -rw-r--r-- | R/inverse.predict.lm.R | 25 | ||||
| -rw-r--r-- | man/inverse.predict.Rd | 2 | ||||
| -rw-r--r-- | vignettes/chemCal.Rmd | 6 |
3 files changed, 15 insertions, 18 deletions
diff --git a/R/inverse.predict.lm.R b/R/inverse.predict.lm.R index 927e672..77d548f 100644 --- a/R/inverse.predict.lm.R +++ b/R/inverse.predict.lm.R @@ -23,10 +23,9 @@ inverse.predict.lm <- function(object, newdata, ..., ws <- ifelse(length(object$weights) > 0, mean(object$weights), 1) } if (length(object$weights) > 0) { - wx <- split(object$weights,object$model[[xname]]) - w <- sapply(wx,mean) + w <- object$weights } else { - w <- rep(1,length(split(object$model[[yname]],object$model[[xname]]))) + w <- rep(1, nrow(object$model)) } .inverse.predict(object = object, newdata = newdata, ws = ws, alpha = alpha, var.s = var.s, w = w, xname = xname, yname = yname) @@ -40,8 +39,7 @@ inverse.predict.rlm <- function(object, newdata, ..., if (ws == "auto") { ws <- mean(object$w) } - wx <- split(object$weights,object$model[[xname]]) - w <- sapply(wx,mean) + w <- object$w .inverse.predict(object = object, newdata = newdata, ws = ws, alpha = alpha, var.s = var.s, w = w, xname = xname, yname = yname) } @@ -57,14 +55,13 @@ inverse.predict.rlm <- function(object, newdata, ..., ybars <- mean(newdata) m <- length(newdata) - yx <- split(object$model[[yname]], object$model[[xname]]) - n <- length(yx) + n <- nrow(object$model) df <- n - length(object$coef) - x <- as.numeric(names(yx)) - ybar <- sapply(yx,mean) - yhatx <- split(object$fitted.values, object$model[[xname]]) - yhat <- sapply(yhatx, mean) - se <- sqrt(sum(w * (ybar - yhat)^2)/df) + xi <- object$model[[2]] + yi <- object$model[[1]] + yihat <- object$fitted.values + + se <- sqrt(sum(w * (yi - yihat)^2)/df) if (var.s == "auto") { var.s <- se^2/ws @@ -72,14 +69,14 @@ inverse.predict.rlm <- function(object, newdata, ..., b1 <- object$coef[[xname]] - ybarw <- sum(w * ybar)/sum(w) + ybarw <- sum(w * yi)/sum(w) # This is the adapted form of equation 8.28 (see package vignette) sxhats <- 1/b1 * sqrt( (var.s / m) + se^2 * (1/sum(w) + (ybars - ybarw)^2 * sum(w) / - (b1^2 * (sum(w) * sum(w * x^2) - sum(w * x)^2))) + (b1^2 * (sum(w) * sum(w * xi^2) - sum(w * xi)^2))) ) if (names(object$coef)[1] == "(Intercept)") { diff --git a/man/inverse.predict.Rd b/man/inverse.predict.Rd index 347d670..26ba6b8 100644 --- a/man/inverse.predict.Rd +++ b/man/inverse.predict.Rd @@ -51,7 +51,7 @@ This is elaborated in the package vignette. } \note{ - The function was validated with examples 7 and 8 from Massart et al. (1997). + The function was validated with examples 7 and 8 from Massart et al. (1997). } \references{ Massart, L.M, Vandenginste, B.G.M., Buydens, L.M.C., De Jong, S., Lewi, P.J., diff --git a/vignettes/chemCal.Rmd b/vignettes/chemCal.Rmd index cea1678..2515abb 100644 --- a/vignettes/chemCal.Rmd +++ b/vignettes/chemCal.Rmd @@ -107,9 +107,9 @@ with s_e = \sqrt{ \frac{\sum w_i (y_i - \hat{y_i})^2}{n - 2}} \end{equation} -where $w_i$ is the weight for calibration standard $i$, $y_i$ is the mean $y$ -value (!) observed for standard $i$, $\hat{y_i}$ is the estimated value for -standard $i$, $n$ is the number calibration standards, $w_s$ is the weight +where $w_i$ is the weight for calibration standard $i$, $y_i$ is the $y$ +value observed for standard $i$, $\hat{y_i}$ is the estimated value for +standard $i$, $n$ is the number of calibration samples, $w_s$ is the weight attributed to the sample $s$, $m$ is the number of replicate measurements of sample $s$, $\bar{y_s}$ is the mean response for the sample, $\bar{y_w} = \frac{\sum{w_i y_i}}{\sum{w_i}}$ is the weighted mean of responses |