diff --git a/.Rbuildignore b/.Rbuildignore new file mode 100644 index 0000000..91114bf --- /dev/null +++ b/.Rbuildignore @@ -0,0 +1,2 @@ +^.*\.Rproj$ +^\.Rproj\.user$ diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..807ea25 --- /dev/null +++ b/.gitignore @@ -0,0 +1,3 @@ +.Rproj.user +.Rhistory +.RData diff --git a/DESCRIPTION b/DESCRIPTION new file mode 100644 index 0000000..b20220b --- /dev/null +++ b/DESCRIPTION @@ -0,0 +1,26 @@ +Package: basicTrendline +Version: 2.0.3 +Date: 2018-07-26 +Title: Add Trendline and Confidence Interval of Basic Regression Models to Plot +Authors@R: c( + person("Weiping", "Mei", email = "meiweipingg@163.com", role = c("aut", "cre", "cph"), comment = c(ORCID = "0000-0001-6400-9862")), + person("Guangchuang", "Yu", email = "guangchuangyu@gmail.com", role = c("aut"), comment = c(ORCID = "0000-0002-6485-8781")), + person("Jiangshan", "Lai", role = c("ctb")), + person("Qiang", "Rao", role = "ctb"), + person("Yu", "Umezawa", role = "ctb") + ) +Maintainer: Weiping Mei +Description: Plot, draw regression line and confidence interval, and show regression equation, R-square and P-value, as simple as possible, by using different models ("line2P", "line3P", "log2P", "exp2P", "exp3P", "power2P", "power3P") built in the 'trendline()' function. +Depends: R (>= 2.1.0) +Imports: + graphics, + stats, + scales, + investr +BugReports: https://github.com/PhDMeiwp/basicTrendline/issues +License: GPL-3 +URL: https://github.com/PhDMeiwp/basicTrendline +LazyData: true +RoxygenNote: 6.0.1 + + diff --git a/NAMESPACE b/NAMESPACE new file mode 100644 index 0000000..a3f4150 --- /dev/null +++ b/NAMESPACE @@ -0,0 +1,12 @@ +# Generated by roxygen2: do not edit by hand + +export(SSexp2P) +export(SSexp3P) +export(SSpower2P) +export(SSpower3P) +export(trendline) +export(trendline_summary) +import(graphics) +import(investr) +import(scales) +import(stats) diff --git a/NEWS b/NEWS new file mode 100644 index 0000000..98742df --- /dev/null +++ b/NEWS @@ -0,0 +1,25 @@ +Changes in version 2.0.3 +------------------------------------------------------------------------------ +Date: 2018-07-26 +* add several arguments to `trendline()` function, including show.equation, show.Rpvalue, Rname, Pname, xname, yname, yhat, CI.fill, CI.level, CI.alpha, CI.color, CI.lty, CI.lwd, ePos.x, ePos.y, las. +* enable to draw confidence interval for regression models (arguments CI.fill, CI.level, etc.) +* add 'show.equation' and show.Rpvale' arguments to enable to choose which parameter to show +* add 'Rname' and 'Pname' arguments to specify the character of R-square and P-vlaue (i.e. R^2 or r^2; P or p) +* add 'xname' and 'ynameto' arguments to specify the character of 'x' and 'y' in the equation +* add 'yhat' argument to enable to add a hat symbol on the top of 'y' in the equation +* add 'ePos.x' and 'ePos.y' arguments to specify the x and y co-ordinates of equation's position +* deleted the 'ePos' argument +* add "Residual Sum of Squares" to the output of 'trendline_summary()' function + + +Changes in version 1.2.0 +------------------------------------------------------------------------------ +Date: 2018-07-13 +* change the expression for `model` of `exp3P` using a supscript +* add `trendline_summary()` function +* add `SSexp2P()` function +* add `SSpower2P` function +* add `Pvalue.corrected` argument in `trendline()` and `trendline_summary()` , for P-vlaue calculation for non-linear regression. +* add `Details` in `trendline()` and `trendline_summary()` +* add `...` argument in `trendline()` as same as those in `plot()` + diff --git a/R/SSexp2P.R b/R/SSexp2P.R new file mode 100644 index 0000000..ddeefc7 --- /dev/null +++ b/R/SSexp2P.R @@ -0,0 +1,45 @@ +#' Self-Starting Nls 'exp2P' Regression Model +#' +#' This selfStart model evaluates the power regression function (formula as: y=a*exp(b*x)). It has an initial attribute that will evaluate initial estimates of the parameters 'a' and 'b' for a given set of data. +#' +#' @usage SSexp2P(predictor, a, b) +#' @param predictor a numeric vector of values at which to evaluate the model. +#' @param a,b The numeric parameters responsing to the exp2P model. +#' @export +#' @examples +#' library(basicTrendline) +#' x<-1:5 +#' y<-c(2,4,8,20,25) +#' xy<-data.frame(x,y) +#' getInitial(y ~ SSexp2P(x,a,b), data = xy) +#' ## Initial values are in fact the converged values +#' +#' fitexp2P <- nls(y~SSexp2P(x,a,b), data=xy) +#' summary(fitexp2P) +#' +#' @author Weiping Mei \email{meiweipingg@163.com} +#' @seealso \code{\link{trendline}}, \code{\link{SSexp3P}}, \code{\link{SSpower3P}}, \code{\link[stats]{nls}}, \code{\link[stats]{selfStart}} + + +# selfStart method for exp2P model (formula as y = a *exp(b*x)) +SSexp2P<-selfStart( + function(predictor,a,b){a*exp(b*predictor)}, + function(mCall,LHS, data) + { + xy <- sortedXyData(mCall[["predictor"]],LHS, data) + + if (min(y)>0){ + lmFit <- lm(log(xy[,"y"]) ~ xy[,"x"]) + coefs <- coef(lmFit) + a <- exp(coefs[1]) #intercept + b <- coefs[2] #slope + value <- c(a, b) + names(value) <- mCall[c("a","b")] + value + }else{stop(" +>>Try to use other selfStart functions. +Because the 'SSexp2P' function need ALL x values greater than 0.") + } + },c("a","b")) + +# getInitial(y~SSexp2P(x,a,b),data = xy) diff --git a/R/SSexp3P.R b/R/SSexp3P.R new file mode 100644 index 0000000..d8534a8 --- /dev/null +++ b/R/SSexp3P.R @@ -0,0 +1,52 @@ +#' Self-Starting Nls 'exp3P' Regression Model +#' +#' This selfStart model evaluates the exponential regression function (formula as: y=a*exp(b*x)+c). It has an initial attribute that will evaluate initial estimates of the parameters a, b, and c for a given set of data. +#' +#' @usage SSexp3P(predictor, a, b, c) +#' @param predictor a numeric vector of values at which to evaluate the model. +#' @param a,b,c Three numeric parameters responsing to the exp3P model. +#' @export +#' @examples +#' library(basicTrendline) +#' x<-1:5 +#' y<-c(2,4,8,16,28) +#' xy<-data.frame(x,y) +#' getInitial(y ~ SSexp3P(x,a,b,c), data = xy) +#' ## Initial values are in fact the converged values +#' +#' fitexp3P <- nls(y~SSexp3P(x,a,b,c), data=xy) +#' summary(fitexp3P) +#' +#' @author Weiping Mei \email{meiweipingg@163.com} +#' @seealso \code{\link{trendline}}, \code{\link{SSexp3P}}, \code{\link{SSpower3P}}, \code{\link[stats]{nls}}, \code{\link[stats]{selfStart}} + +# selfStart method for exp3P model (formula as y = a *exp(b*x)+ c) +SSexp3P<-selfStart( + function(predictor,a,b,c){a*exp(b*predictor)+c}, + function(mCall,LHS, data) + { + xy <- sortedXyData(mCall[["predictor"]],LHS, data) + y=xy[,"y"] + x=xy[,"x"] + adjy=y-min(y)+1 + xadjy=data.frame(x,adjy) + + lmFit <- lm(log(adjy) ~ x) + coefs <- coef(lmFit) + get.b <- coefs[2] #slope + + nlsFit<-nls(adjy~cbind(1+exp(b*x),exp(b*x)), + start = list(b=get.b),data = xadjy,algorithm = "plinear", + nls.control(maxiter = 5000000,minFactor = 10^(-10))) + + coef<-coef(nlsFit) + b<-coef[1] + c<-coef[2]+min(y)-1 + a<-coef[3]+coef[2] + + value <- c(a,b,c) + names(value) <- mCall[c("a","b","c")] + value + },c("a","b","c")) + + # getInitial(y~SSexp3P(x,a,b,c),data = z) diff --git a/R/SSpower2P.R b/R/SSpower2P.R new file mode 100644 index 0000000..beaefca --- /dev/null +++ b/R/SSpower2P.R @@ -0,0 +1,47 @@ +#' Self-Starting Nls 'power2P' Regression Model +#' +#' This selfStart model evaluates the power regression function (formula as: y=a*x^b). It has an initial attribute that will evaluate initial estimates of the parameters 'a' and 'b' for a given set of data. +#' +#' @usage SSpower2P(predictor, a, b) +#' @param predictor a numeric vector of values at which to evaluate the model. +#' @param a,b The numeric parameters responsing to the exp2P model. +#' @export +#' @examples +#' library(basicTrendline) +#' x<-1:5 +#' y<-c(2,4,8,20,25) +#' xy<-data.frame(x,y) +#' getInitial(y ~ SSpower2P(x,a,b), data = xy) +#' ## Initial values are in fact the converged values +#' +#' fitpower2P <- nls(y~SSpower2P(x,a,b), data=xy) +#' summary(fitpower2P) +#' +#' @author Weiping Mei \email{meiweipingg@163.com} +#' @seealso \code{\link{trendline}}, \code{\link{SSexp3P}}, \code{\link{SSpower3P}}, \code{\link[stats]{nls}}, \code{\link[stats]{selfStart}} + + +# selfStart method for power2P model (formula as y = a *x^b) +SSpower2P<-selfStart( + function(predictor,a,b){a*predictor^b}, + function(mCall,LHS, data) +{ + xy <- sortedXyData(mCall[["predictor"]],LHS, data) + + if (min(x)>0){ + lmFit <- lm(log(xy[,"y"]) ~ log(xy[,"x"])) # both x and adjy values should be greater than 0. + coefs <- coef(lmFit) + a <- exp(coefs[1]) #intercept + b <- coefs[2] #slope + + value <- c(a,b) + names(value) <- mCall[c("a","b")] + value + + }else{stop(" +>>Try to use other selfStart functions. +Because the 'SSpower2P' function need ALL x values greater than 0.") + } + },c("a","b")) + +# getInitial(y~SSpower2P(x,a,b),data = xy) diff --git a/R/SSpower3P.R b/R/SSpower3P.R new file mode 100644 index 0000000..fc72975 --- /dev/null +++ b/R/SSpower3P.R @@ -0,0 +1,60 @@ +#' Self-Starting Nls 'power3P' Regression Model +#' +#' This selfStart model evaluates the power regression function (formula as: y=a*x^b+c). It has an initial attribute that will evaluate initial estimates of the parameters a, b, and c for a given set of data. +#' +#' @usage SSpower3P(predictor, a, b, c) +#' @param predictor a numeric vector of values at which to evaluate the model. +#' @param a,b,c Three numeric parameters responsing to the exp3P model. +#' @export +#' @examples +#' library(basicTrendline) +#' x<-1:5 +#' y<-c(2,4,8,20,25) +#' xy<-data.frame(x,y) +#' getInitial(y ~ SSpower3P(x,a,b,c), data = xy) +#' ## Initial values are in fact the converged values +#' +#' fitpower3P <- nls(y~SSpower3P(x,a,b,c), data=xy) +#' summary(fitpower3P) +#' +#' @author Weiping Mei \email{meiweipingg@163.com} +#' @seealso \code{\link{trendline}}, \code{\link{SSexp3P}}, \code{\link{SSpower3P}}, \code{\link[stats]{nls}}, \code{\link[stats]{selfStart}} + +# selfStart method for power3P model (formula as y = a *x^b+ c) +SSpower3P<-selfStart( + function(predictor,a,b,c){a*predictor^b+c}, + function(mCall,LHS, data) + { + xy <- sortedXyData(mCall[["predictor"]],LHS, data) + y=xy[,"y"] + x=xy[,"x"] + + if (min(x)>0){ + + adjy=y-min(y)+1 + xadjy=data.frame(x,adjy) + + lmFit <- lm(log(adjy) ~ log(x)) # both x and adjy values should be greater than 0. + coefs <- coef(lmFit) + get.b <- coefs[2] #slope + + nlsFit<-nls(adjy~cbind(1+x^b,x^b), + start = list(b=get.b),data = xadjy,algorithm = "plinear", + nls.control(maxiter = 5000000,minFactor = 10^(-10))) + + coef<-coef(nlsFit) + b<-coef[1] + c<-coef[2]+min(y)-1 + a<-coef[3]+coef[2] + + value <- c(a,b,c) + names(value) <- mCall[c("a","b","c")] + value + + }else{stop(" +>>Try to use other selfStart functions. +Because the 'SSpower3P' function need ALL x values greater than 0.") + } + },c("a","b","c")) + + # getInitial(y~SSpower3P(x,a,b,c),data = xy) diff --git a/R/trendline.R b/R/trendline.R new file mode 100644 index 0000000..69eb6d3 --- /dev/null +++ b/R/trendline.R @@ -0,0 +1,478 @@ +#' Add Trendline and Show Equation to Plot +#' +#' Plot, draw regression line and confidence interval, and show regression equation, R-square and P-value, as simple as possible, +#' by using different models built in the 'trendline()' function. The function includes the following models in the latest version: +#' "line2P" (formula as: y=a*x+b), "line3P" (y=a*x^2+b*x+c), "log2P" (y=a*ln(x)+b), "exp2P" (y=a*exp(b*x)),"exp3P" (y=a*exp(b*x)+c), "power2P" (y=a*x^b), and "power3P" (y=a*x^b+c). +#' Besides, the summarized result of each fitted model is also output by default. +#' +#' @param x,y the x and y arguments provide the x and y coordinates for the plot. Any reasonable way of defining the coordinates is acceptable. +#' @param model select which model to fit. Default is "line2P". The "model" should be one of c("line2P", "line3P", "log2P", "exp2P", "exp3P", "power2P", "power3P"), their formulas are as follows:\cr "line2P": y=a*x+b \cr "line3P": y=a*x^2+b*x+c \cr "log2P": y=a*ln(x)+b \cr "exp2P": y=a*exp(b*x) \cr "exp3P": y=a*exp(b*x)+c \cr "power2P": y=a*x^b \cr "power3P": y=a*x^b+c +#' @param Pvalue.corrected if P-value corrected or not, the value is one of c("TRUE", "FALSE"). +#' @param linecolor color of regression line. +#' @param lty line type. lty can be specified using either text c("blank","solid","dashed","dotted","dotdash","longdash","twodash") or number c(0, 1, 2, 3, 4, 5, 6). Note that lty = "solid" is identical to lty=1. +#' @param lwd line width. Default is 1. +#' @param show.equation whether to show the regression equation, the value is one of c("TRUE", "FALSE"). +#' @param show.Rpvalue whether to show the R-square and P-value, the value is one of c("TRUE", "FALSE"). +#' @param Rname to specify the character of R-square, the value is one of c(0, 1), corresponding to c(r^2, R^2). +#' @param Pname to specify the character of P-value, the value is one of c(0, 1), corresponding to c(p, P). +#' @param xname to specify the character of "x" in equation, see Examples [case 5]. +#' @param yname to specify the character of "y" in equation, see Examples [case 5]. +#' @param yhat whether to add a hat symbol (^) on the top of "y" in equation. Default is FALSE. +#' @param CI.fill fill the confidance interval? (TRUE by default, see 'CI.level' to control) +#' @param CI.level level of confidence interval to use (0.95 by default) +#' @param CI.alpha alpha value of fill color of confidence interval. +#' @param CI.color line or fill color of confidence interval. +#' @param CI.lty line type of confidence interval. +#' @param CI.lwd line width of confidence interval. +#' @param summary summarizing the model fits. Default is TRUE. +#' @param ePos.x,ePos.y equation position. Default as ePos.x = "topleft". If no need to show equation, set ePos.x = NA. It's same as those in \code{\link[graphics]{legend}}. +#' @param text.col the color used for the equation text. +#' @param eDigit the numbers of digits for equation parameters. Default is 5. +#' @param eSize font size in percentage of equation. Default is 1. +#' @param xlab,ylab labels of x- and y-axis. +#' @param las style of axis labels. (0=parallel, 1=all horizontal, 2=all perpendicular to axis, 3=all vertical) +#' @param ... additional parameters to \code{\link[graphics]{plot}}, such as type, main, sub, pch, col. +#' @import graphics +#' @import stats +#' @import scales +#' @import investr +#' @export +#' @details The linear models (line2P, line3P, log2P) in this package are estimated by \code{\link[stats]{lm}} function, \cr while the nonlinear models (exp2P, exp3P, power2P, power3P) are estimated by \code{\link[stats]{nls}} function (i.e., least-squares method).\cr\cr The argument 'Pvalue.corrected' is workful for non-linear regression only.\cr\cr If "Pvalue.corrected = TRUE", the P-value is calculated by using "Residual Sum of Squares" and "Corrected Total Sum of Squares (i.e. sum((y-mean(y))^2))".\cr If "Pvalue.corrected = TRUE", the P-value is calculated by using "Residual Sum of Squares" and "Uncorrected Total Sum of Squares (i.e. sum(y^2))". +#' @note +#' Confidence intervals for nonlinear regression (i.e., objects of class +#' \code{nls}) are based on the linear approximation described in Bates & Watts (2007) and Greenwell & Schubert-Kabban (2014). +#' +#' @references +#' Bates, D. M., and Watts, D. G. (2007) +#' \emph{Nonlinear Regression Analysis and its Applications}. Wiley. +#' +#' Greenwell B. M., and Schubert-Kabban, C. M. (2014) +#' \emph{investr: An R Package for Inverse Estimation}. The R Journal, 6(1), 90-100. +#' @return NULL +#' @examples +#' library(basicTrendline) +#' x <- c(1, 3, 6, 9, 13, 17) +#' y <- c(5, 8, 11, 13, 13.2, 13.5) +#' +#' # [case 1] default +#' trendline(x, y, model="line2P", ePos.x = "topleft", summary=TRUE, eDigit=5) + +#' # [case 2] draw lines of confidenc interval only (set CI.fill = FALSE) +#' trendline(x, y, model="line3P", CI.fill = FALSE, CI.color = "black", CI.lty = 2, linecolor = "blue") +#' +#' # [case 3] draw trendliine only (set CI.color = NA) +#' trendline(x, y, model="log2P", ePos.x= "top", linecolor = "red", CI.color = NA) +#' +#' # [case 4] show regression equation only (set show.Rpvalue = FALSE) +#' trendline(x, y, model="exp2P", show.equation = TRUE, show.Rpvalue = FALSE) +#' +#' # [case 5] specify the name of parameters in equation +#' # see Arguments c('xname', 'yname', 'yhat', 'Rname', 'Pname'). +#' trendline(x, y, model="exp3P", xname="T", yname=paste(delta^15,"N"), +#' yhat=FALSE, Rname=1, Pname=0, ePos.x = "bottom") +#' +#' # [case 6] change the digits, font size, and color of equation. +#' trendline(x, y, model="power2P", eDigit = 3, eSize = 1.4, text.col = "blue") +#' +#' # [case 7] don't show equation (set ePos.x = NA) +#' trendline(x, y, model="power3P", ePos.x = NA) +#' +#' ### NOT RUN +#' # [case 8] set graphical parameters by par {graphics} +#' +#' par(mgp=c(1.5,0.4,0), mar=c(3,3,1,1), tck=-0.01, cex.axis=0.9) +#' +#' trendline(x, y) +#' +#' dev.off() +#' +#' ### END (NOT RUN) +#' +#' @author Weiping Mei, Guangchuang Yu +#' @seealso \code{\link{trendline}}, \code{\link{SSexp3P}}, \code{\link{SSpower3P}}, \code{\link[stats]{nls}}, \code{\link[stats]{selfStart}}, \code{\link[investr]{plotFit}} + +trendline <- function(x, y, model="line2P", Pvalue.corrected = TRUE, + linecolor = "blue", lty = 1, lwd = 1, + show.equation = TRUE, show.Rpvalue = TRUE, + Rname = 1, Pname = 0, xname = "x", yname = "y", yhat = FALSE, + summary = TRUE, + ePos.x = NULL, ePos.y = NULL, text.col="black", eDigit = 5, eSize = 1, + CI.fill = TRUE, CI.level = 0.95, CI.color = "grey", CI.alpha = 1, CI.lty = 1, CI.lwd = 1, + las = 1, xlab=NULL, ylab=NULL, ...) +{ + model=model + if(is.null(xlab)) xlab = deparse(substitute(x)) else xlab = xlab + if(is.null(ylab)) ylab = deparse(substitute(y)) else ylab = ylab + if(Rname==0) Rname = "r" else Rname = "R" + if(Pname==0) Pname = "p" else Pname = "P" + xname = substitute(xname) + if(yhat == TRUE) yname = substitute(hat(yname)) else yname = substitute(yname) + + OK <- complete.cases(x, y) + x <- x[OK] + y <- y[OK] + z<-data.frame(x,y) + + return <- trendline_summary(x=x, y=y, model=model, Pvalue.corrected=Pvalue.corrected, summary = FALSE, eDigit = eDigit) + a = return$parameter$a + b = return$parameter$b + if (is.null(return$parameter$c)==FALSE){ + c = return$parameter$c + }else{} + if (return$p.value >= 0.0001){ + pval <- return$p.value + pval <- paste("=" , unname(pval)) + }else{ + pval <- "< 0.0001" + } + r2 <- return$R.squared + adjr2<- return$adj.R.squared + + +# 1) model="line2P" +if (model== c("line2P")) + { Pvalue.corrected=TRUE + formula = 'y = a*x + b' + fitting <- lm(y~x) + + if (summary==TRUE){ + trendline_summary(x,y,"line2P", Pvalue.corrected=Pvalue.corrected, eDigit=eDigit) + }else{} + + aa = abs(a) + bb = abs(b) + aa = format(aa, digits = eDigit) + bb = format(bb, digits = eDigit) + + param <- vector('expression',2) + if (aa==1){aa=c("")} + + if (a>0) + { + if (b>=0) + {param[1] <- substitute(expression(italic(yname) == aa~italic(xname) + bb))[2] + }else{param[1] <- substitute(expression(italic(yname) == aa~italic(xname) - bb))[2] + } + }else{ + if (b>=0) + {param[1] <- substitute(expression(italic(yname) == -aa~italic(xname) + bb))[2] + }else{param[1] <- substitute(expression(italic(yname) == -aa~italic(xname) - bb))[2] + } + } + + param[2] <- substitute(expression(italic(Rname)^2 == r2*","~~italic(Pname)~~pval))[2] + + } + +# 2) model="line3P" + if (model== c("line3P")) + { Pvalue.corrected=TRUE + formula = 'y = a*x^2 + b*x + c' + fitting <- lm(y~I(x^2)+x) + + if (summary==TRUE){ + trendline_summary(x,y,"line3P", Pvalue.corrected=Pvalue.corrected, eDigit=eDigit) + }else{} + + + aa = abs(a) + bb = abs(b) + cc = abs(c) + aa = format(aa, digits = eDigit) + bb = format(bb, digits = eDigit) + cc = format(cc, digits = eDigit) + + param <- vector('expression',2) + + if (aa==1){aa=c("")} + if (bb==1){bb=c("")} + + if (a>0) + { + if (b>=0) + { + if(c>=0) + {param[1] <- substitute(expression(italic(yname) == aa~italic(xname)^2 + bb~italic(xname) +cc))[2] + }else{param[1] <- substitute(expression(italic(yname) == aa~italic(xname)^2 + bb~italic(xname) -cc))[2] + } + }else{ + if(c>=0) + {param[1] <- substitute(expression(italic(yname) == aa~italic(xname)^2 - bb~italic(xname) +cc))[2] + }else{param[1] <- substitute(expression(italic(yname) == aa~italic(xname)^2 - bb~italic(xname) -cc))[2] + } + } + + }else{ + if (b>=0) + { + if(c>=0) + {param[1] <- substitute(expression(italic(yname) == -aa~italic(xname)^2 + bb~italic(xname) +cc))[2] + }else{param[1] <- substitute(expression(italic(yname) == -aa~italic(xname)^2 + bb~italic(xname) -cc))[2] + } + }else{ + if(c>=0) + {param[1] <- substitute(expression(italic(yname) == -aa~italic(xname)^2 - bb~italic(xname) +cc))[2] + }else{param[1] <- substitute(expression(italic(yname) == -aa~italic(xname)^2 - bb~italic(xname) -cc))[2] + } + } + + } + + param[2] <- substitute(expression(italic(Rname)^2 == r2*","~~italic(Pname)~~pval*" "))[2] + + } + +# 3) model="log2P" +if (model== c("log2P")) + { + Pvalue.corrected=TRUE + formula = 'y = a*ln(x) + b' + if (summary==TRUE){ + trendline_summary(x,y,"log2P", Pvalue.corrected=Pvalue.corrected, eDigit=eDigit) + }else{} + + if (min(x)>0) + { + fitting <- lm(y~log(x)) + + aa = abs(a) + bb = abs(b) + + param <- vector('expression',2) + + aa = format(aa, digits = eDigit) + bb = format(bb, digits = eDigit) + + if (aa==1){aa=c("")} + + if (a>0) + { + if (b>=0) + { + param[1] <- substitute(expression(italic(yname) == aa~"ln(x)" + bb))[2] + }else{ + param[1] <- substitute(expression(italic(yname) == aa~"ln(x)" - bb))[2] + } + + }else{ + if (b>=0) + { + param[1] <- substitute(expression(italic(yname) == -aa~"ln(x)" + bb))[2] + }else{ + param[1] <- substitute(expression(italic(yname) == -aa~"ln(x)" - bb))[2] + } + } + + param[2] <- substitute(expression(italic(Rname)^2 == r2*","~~italic(Pname)~~pval))[2] + + }else{ + stop(" +'log2P' model need ALL x values greater than 0. Try other models.") + } +} + +# 4.2) model="exp2P" + if (model== "exp2P") + { + formula = 'y = a*exp(b*x)' + fitting <- nls(y~SSexp2P(x,a,b),data=z) + + if (summary==TRUE){ + trendline_summary(x, y, "exp2P", Pvalue.corrected=Pvalue.corrected, eDigit=eDigit) + }else{} + + aa= abs(a) + bb= abs(b) + param <- vector('expression',2) + + aa = format(aa, digits = eDigit) + bb = format(bb, digits = eDigit) + + if (aa==1){aa=c("")} + if (bb==1){bb=c("")} + + if (a>=0) + { + if (b>=0) + { + param[1] <- substitute(expression(italic(yname) == aa~"e"^{bb~italic(xname)}))[2] + }else{ + param[1] <- substitute(expression(italic(yname) == aa~"e"^{-bb~italic(xname)}))[2] + } + + }else{ + if (b>=0) + { + param[1] <- substitute(expression(italic(yname) == -aa~"e"^{bb~italic(xname)}))[2] + }else{ + param[1] <- substitute(expression(italic(yname) == -aa~"e"^{-bb~italic(xname)}))[2] + } + } + + param[2] <- substitute(expression(italic(Rname)^2 == r2*","~~italic(Pname)~~pval))[2] + + } + + +# 4.3) model="exp3P" + if (model== "exp3P") + { + formula = 'y = a*exp(b*x) + c' + fitting <- nls(y~SSexp3P(x,a,b,c),data=z) + + if (summary==TRUE){ + trendline_summary(x, y, "exp3P", Pvalue.corrected=Pvalue.corrected, eDigit=eDigit) + }else{} + + aa= abs(a) + bb= abs(b) + cc= abs(c) + + param <- vector('expression',2) + + aa = format(aa, digits = eDigit) + bb = format(bb, digits = eDigit) + cc = format(cc, digits = eDigit) + + if (aa==1){aa=c("")} + if (bb==1){bb=c("")} + + if (a>=0) + { + if (b>=0) + { + if (c>=0){ + param[1] <- substitute(expression(italic(yname) == aa~"e"^{bb~italic(xname)}~+cc))[2] + }else{ + param[1] <- substitute(expression(italic(yname) == aa~"e"^{bb~italic(xname)}~-cc))[2] + } + }else{ + if (c>=0){ + param[1] <- substitute(expression(italic(yname) == aa~"e"^{-bb~italic(xname)}~+cc))[2] + }else{ + param[1] <- substitute(expression(italic(yname) == aa~"e"^{-bb~italic(xname)}~-cc))[2] + } + } + + }else{ + if (b>=0) + { + if (c>=0){ + param[1] <- substitute(expression(italic(yname) == -aa~"e"^{bb~italic(xname)}~+cc))[2] + }else{ + param[1] <- substitute(expression(italic(yname) == -aa~"e"^{bb~italic(xname)}~-cc))[2] + } + }else{ + if (c>=0){ + param[1] <- substitute(expression(italic(yname) == -aa~"e"^{-bb~italic(xname)}~+cc))[2] + }else{ + param[1] <- substitute(expression(italic(yname) == -aa~"e"^{-bb~italic(xname)}~-cc))[2] + } + } +} + param[2] <- substitute(expression(italic(Rname)^2 == r2*","~~italic(Pname)~~pval))[2] + } + + +# 5.2) model="power2P" +if (model== "power2P") + {formula = 'y = a*x^b' + + if (summary==TRUE){ + trendline_summary(x, y, "power2P", Pvalue.corrected=Pvalue.corrected, eDigit=eDigit) + }else{} + + if (min(x)>0){ + fitting <- nls(y~SSpower2P(x,a,b),data=z) + + aa<-abs(a) + + param <- vector('expression',2) + + aa = format(aa, digits = eDigit) + + if (aa==1){aa=c("")} + + if (a>=0) + { + param[1] <- substitute(expression(italic(yname) == aa~italic(xname)^b))[2] + }else{ + param[1] <- substitute(expression(italic(yname) == -aa~italic(xname)^b))[2] + } + param[2] <- substitute(expression(italic(Rname)^2 == r2*","~~italic(Pname)~~pval))[2] + + }else{ + stop(" + 'power2P' model need ALL x values greater than 0. Try other models.") + } +} + + + # 5.3) model="power3P" +if (model== "power3P") + {formula = 'y = a*x^b + c' + + if (summary==TRUE){ + trendline_summary(x,y,"power3P", Pvalue.corrected=Pvalue.corrected, eDigit=eDigit) + }else{} + + if (min(x)>0){ + fitting <- nls(y~SSpower3P(x,a,b,c),data=z) + + aa<-abs(a) + cc<-abs(c) + + param <- vector('expression',2) + + aa = format(aa, digits = eDigit) + cc = format(cc, digits = eDigit) + + if (aa==1){aa=c("")} + + if (a>=0) + { + if (c>=0){ + param[1] <- substitute(expression(italic(yname) == aa~italic(xname)^b ~ + cc))[2] + }else{ + param[1] <- substitute(expression(italic(yname) == aa~italic(xname)^b ~ - cc))[2] + } + + }else{ + if (c>=0){ + param[1] <- substitute(expression(italic(yname) == -aa~italic(xname)^b ~ + cc))[2] + }else{ + param[1] <- substitute(expression(italic(yname) == -aa~italic(xname)^b ~ - cc))[2] + } + } + param[2] <- substitute(expression(italic(Rname)^2 == r2*","~~italic(Pname)~~pval))[2] + + }else{ + stop(" +'power3P' model need ALL x values greater than 0. Try other models.") + } + +# 100) beyond the built-in models. + +}else{ + Check<-c("line2P","line3P","log2P","exp2P","exp3P","power2P","power3P") + if (!model %in% Check) + stop(" +\"model\" should be one of c(\"lin2P\",\"line3P\",\"log2P\",\"exp2P\",\"exp3P\",\"power2P\",\"power3P\".") +} + +### plot and draw trendline + if (requireNamespace(c("investr", "scales"), quietly = TRUE)){ + investr::plotFit(fitting, interval = "confidence", level = CI.level, data=z, + shade = CI.fill, + col.fit = linecolor, lty.fit = lty, lwd.fit = lwd, + col.conf = scales::alpha(CI.color, alpha = CI.alpha),lty.conf = CI.lty, lwd.conf = CI.lwd, + las = las, xlab = xlab, ylab = ylab, ...) + } + +### show legend + if (show.equation == TRUE) param[1] = param[1] else param[1]=NULL + if (show.Rpvalue == TRUE) param[2] = param[2] else param[2]=NULL + if (is.null(ePos.x)) ePos.x = "topleft" else ePos.x = ePos.x + legend(ePos.x, ePos.y, text.col = text.col, legend = param, cex = eSize, bty = 'n') + +} diff --git a/R/trendline_summary.R b/R/trendline_summary.R new file mode 100644 index 0000000..855fb35 --- /dev/null +++ b/R/trendline_summary.R @@ -0,0 +1,588 @@ +#' Summarized Results of Each Regression Model +#' +#' Summarizing the results of each regression model which built in the 'trendline()' function. The function includes the following models in the latest version: +#' "line2P" (formula as: y=a*x+b), "line3P" (y=a*x^2+b*x+c), "log2P" (y=a*ln(x)+b), "exp2P" (y=a*exp(b*x)),"exp3P" (y=a*exp(b*x)+c), "power2P" (y=a*x^b), and "power3P" (y=a*x^b+c). +#' +#' @param x,y the x and y arguments provide the x and y coordinates for the plot. Any reasonable way of defining the coordinates is acceptable. +#' @param model select which model to fit. Default is "line2P". The "model" should be one of c("line2P", "line3P", "log2P", "exp2P", "exp3P", "power2P", "power3P"), their formulas are as follows:\cr "line2P": y=a*x+b \cr "line3P": y=a*x^2+b*x+c \cr "log2P": y=a*ln(x)+b \cr "exp2P": y=a*exp(b*x) \cr "exp3P": y=a*exp(b*x)+c \cr "power2P": y=a*x^b \cr "power3P": y=a*x^b+c +#' @param Pvalue.corrected if P-value corrected or not, the vlaue is one of c("TRUE", "FALSE"). +#' @param summary summarizing the model fits. Default is TRUE. +#' @param eDigit the numbers of digits for summarized results. Default is 5. +#' @import stats +#' @export +#' @details The linear models (line2P, line3P, log2P) in this package are estimated by \code{\link[stats]{lm}} function, \cr while the nonlinear models (exp2P, exp3P, power2P, power3P) are estimated by \code{\link[stats]{nls}} function (i.e., least-squares method).\cr\cr The argument 'Pvalue.corrected' is workful for non-linear regression only.\cr\cr If "Pvalue.corrected = TRUE", the P-vlaue is calculated by using "Residual Sum of Squares" and "Corrected Total Sum of Squares (i.e. sum((y-mean(y))^2))".\cr If "Pvalue.corrected = TRUE", the P-vlaue is calculated by using "Residual Sum of Squares" and "Uncorrected Total Sum of Squares (i.e. sum(y^2))". + +#' @return R^2, indicates the R-Squared value of each regression model. +#' @return p, indicates the p-value of each regression model. +#' @return N, indicates the sample size. +#' @return AIC or BIC, indicate the Akaike's Information Criterion or Bayesian Information Criterion for fitted model. Click \code{\link[stats]{AIC}} for details. The smaller the AIC or BIC, the better the model. +#' @return RSS, indicate the value of "Residual Sum of Squares". +#' @examples +#' library(basicTrendline) +#' x1<-1:5 +#' x2<- -2:2 +#' x3<- c(101,105,140,200,660) +#' x4<- -5:-1 +#' x5<- c(1,30,90,180,360) +#' +#' y1<-c(2,14,18,19,20) # increasing convex trend +#' y2<- c(-2,-14,-18,-19,-20) # decreasing concave trend +#' y3<-c(2,4,16,38,89) # increasing concave trend +#' y4<-c(-2,-4,-16,-38,-89) # decreasing convex trend +#' y5<- c(600002,600014,600018,600019,600020) # high y values with low range. +#' +#' trendline_summary(x1,y1,model="line2P",summary=TRUE,eDigit=10) +#' trendline_summary(x2,y2,model="line3P",summary=FALSE) +#' trendline_summary(x3,y3,model="log2P") +#' trendline_summary(x4,y4,model="exp3P") +#' trendline_summary(x5,y5,model="power3P") +#' +#' @author Weiping Mei, Guangchuang Yu +#' @seealso \code{\link{trendline}}, \code{\link{SSexp3P}}, \code{\link{SSpower3P}}, \code{\link[stats]{nls}}, \code{\link[stats]{selfStart}} + +trendline_summary <- function(x,y,model="line2P", Pvalue.corrected=TRUE, summary=TRUE, eDigit=5) +{ + model=model + + OK <- complete.cases(x, y) + x <- x[OK] + y <- y[OK] + z<-data.frame(x,y) + z<-na.omit(z) + nrow = nrow(z) + + # 1) model="line2P" + if (model== c("line2P")) + { + Pvalue.corrected=TRUE + + formula = 'y = a*x + b' + + fit<- lm(y~x) + sum.line2P <- summary(fit) + ss.res<-sum((residuals(fit))^2) # Residual Sum of Squares, DF= n-k + + if (summary==TRUE){ + print(sum.line2P,digits=eDigit) + }else{} + + coeff<-sum.line2P$coefficients + a<-coeff[2,1] # slope + b<-coeff[1,1] # intercept + + n<-length(x) + pval <- coeff[2,4] # p-value of parameter "a", indicates the p-value of whole model. + pval<-unname(pval) + r2 <- sum.line2P$r.squared + adjr2<- sum.line2P$adj.r.squared + + a = format(a, digits = eDigit) + b = format(b, digits = eDigit) + r2 = format(r2, digits = eDigit) + adjr2 = format(adjr2, digits = eDigit) + pval = format(pval, digits = eDigit) + + a=as.numeric(a) + b=as.numeric(b) + param.out<- c(list("a"=a,"b"=b)) # for return values + + } + + # 2) model="line3P" + if (model== c("line3P")) + { + Pvalue.corrected=TRUE + + formula = 'y = a*x^2 + b*x + c' + + fit<-lm(y~I(x^2)+x) + + sum.line3P <- summary(fit) + ss.res<-sum((residuals(fit))^2) # Residual Sum of Squares, DF= n-k + + if (summary==TRUE){ + print(sum.line3P,digits=eDigit) + }else{} + + coeff<-coef(sum.line3P) + a<-coeff[2,1] # slope of x.square + b<-coeff[3,1] # slope of x + c<-coeff[1,1] # intercept c + + n<-length(x) + r2<-sum.line3P$r.squared + adjr2 <- sum.line3P$adj.r.squared + + fstat<-sum.line3P$fstatistic + pval<-pf(fstat[1], fstat[2], fstat[3], lower.tail=FALSE) #p-value of whole model. + pval<-unname(pval) + + a = format(a, digits = eDigit) + b = format(b, digits = eDigit) + c = format(c, digits = eDigit) + r2 = format(r2, digits = eDigit) + adjr2 = format(adjr2, digits = eDigit) + pval = format(pval, digits = eDigit) + + a=as.numeric(a) + b=as.numeric(b) + c=as.numeric(c) + param.out<- c(list("a"=a,"b"=b,"c"=c)) + } + + # 3) model="log2P" + if (model== c("log2P")) + { + Pvalue.corrected=TRUE + + formula = 'y = a*ln(x) + b' + + yadj<-y-min(y) #adjust + + if (min(x)>0) + { + if (summary==TRUE){ + fit0<-lm(y~log(x)) + sum.log0<-summary(fit0) + ss.res<-sum((residuals(fit0))^2) # Residual Sum of Squares, DF= n-k + print(sum.log0, digits = eDigit) + }else{} + + fit<-lm(yadj~log(x)) # adjusted y used + sum.log<-summary(fit) + ss.res<-sum((residuals(fit))^2) # Residual Sum of Squares, DF= n-k + a<-sum.log$coefficients[2,1] # slope + b<-sum.log$coefficients[1,1] # intercept + b=b+min(y) #re-adjust + + fstat<-sum.log$fstatistic + pval<-pf(fstat[1], fstat[2], fstat[3], lower.tail=FALSE) #p-value of whole model. + pval<-unname(pval) + + n<-length(x) + r2<-sum.log$r.squared + adjr2 <- sum.log$adj.r.squared + + a = format(a, digits = eDigit) + b = format(b, digits = eDigit) + r2 = format(r2, digits = eDigit) + adjr2 = format(adjr2, digits = eDigit) + pval= format(pval, digits = eDigit) + + a=as.numeric(a) + b=as.numeric(b) + param.out<- c(list("a"=a,"b"=b)) + + }else{ + stop(" + 'log2P' model need ALL x values greater than 0. Try other models.") + } + } + + + # 4.2) model="exp2P" + if (model== c("exp2P")) + { + formula = 'y = a*exp(b*x)' + + n=length(x) + k = 2 # k means the count numbers of parameters(i.e., 'a', 'b' and 'c' in this case) + + fit<-nls(y~SSexp2P(x,a,b),data=z) + sum.exp2P <- summary(fit) # Get the exact value of each parameter. + + ### calculate the F-statistic and p-value for model + ss.res<-sum((residuals(fit))^2) # Residual Sum of Squares, DF= n-k + ss.total.uncor<-sum(y^2) # Uncorrected Total Sum of Squares, DF=n + ss.total.cor<-sum((y-mean(y))^2) # Corrected Total Sum of Squares, DF=n-1 + + if (Pvalue.corrected==TRUE){ + ss.reg <- ss.total.cor - ss.res # Regression Sum of Squares, DF= (n-1)-(n-k) = k-1 in this case + dfR= k-1 + }else{ + ss.reg <- ss.total.uncor - ss.res # Regression Sum of Squares, DF= n-(n-k) = k in this case + dfR= k + } + + dfE= n-k # degrees of freedom for Error (or Residuals) + + Fval=(ss.reg/dfR)/(ss.res/dfE) + pval=pf(Fval,dfR,dfE,lower.tail = F) + pval<-unname(pval) + + RSE<-sum.exp2P$sigma # Residual standard error, type ?summary.nls in R for more detials. + SSE<-(RSE^2)*(n-1) # Sum of Squares for Error, not equal to 'ss.res'. + + adjr2 <- 1-SSE/((var(y))*(n-1)) + r2<-1-(1-adjr2)*((n-k)/(n-1)) + + if (summary==TRUE){ + ### Start print step by step + coeff = coef(sum.exp2P) + # print + cat("\nNonlinear regression model\n") + cat("\nFormula: y = a*exp(b*x)","\n") + df <- sum.exp2P$df + rdf <- df[2L] + cat("\nParameters:\n") + printCoefmat(coeff, digits = eDigit) + cat("\nResidual standard error:", + format(sum.exp2P$sigma, digits = eDigit), "on", rdf, "degrees of freedom","\n") + + convInfo = fit$convInfo + iterations<-convInfo$finIter + tolerance<-convInfo$finTol + + cat("\nNumber of iterations to convergence:", + format(iterations, digits = eDigit)) + cat("\nAchieved convergence tolerance:",format(tolerance, digits = eDigit),"\n") + + cat("\nMultiple R-squared:", + format(r2, digits = eDigit), ", Adjusted R-squared: ", + format(adjr2, digits = eDigit)) + cat("\nF-statistic:", + format(Fval, digits = eDigit), "on", dfR, "and", dfE, "DF, ", "p-value:", format(pval, digits = eDigit), "\n") + ### finished print + }else{} + + coeffs<-sum.exp2P$coefficients + a<-coeffs[1,1] + b<-coeffs[2,1] + + a = format(a, digits = eDigit) + b = format(b, digits = eDigit) + r2 = format(r2, digits = eDigit) + adjr2 = format(adjr2, digits = eDigit) + pval= format(pval, digits = eDigit) + + a=as.numeric(a) + b=as.numeric(b) + param.out<- c(list("a"=a,"b"=b)) + + } + + # 4.3) model="exp3P" + if (model== c("exp3P")) + { + formula = 'y = a*exp(b*x) + c' + + yadj<-y-min(y)+1 + zzz<-data.frame(x,yadj) + + n=length(x) + k = 3 # k means the count numbers of parameters(i.e., 'a', 'b' and 'c' in this case) + + # use selfStart function 'SSexp3P' for y = a *exp(b*x)+ c + # fit model + fit<-nls(yadj~SSexp3P(x,a,b,c),data=zzz) # use 'yadj', in case of extreme high y-values with low range, such as y= c(600002,600014,600018,600019,600020). + sum.exp3P <- summary(fit) # Get the exact value of each parameter. + + ### calculate the F-statistic and p-value for model + ss.res<-sum((residuals(fit))^2) # Residual Sum of Squares, DF= n-k + ss.total.uncor<-sum(y^2) # Uncorrected Total Sum of Squares, DF=n + ss.total.cor<-sum((y-mean(y))^2) # Corrected Total Sum of Squares, DF=n-1 + + if (Pvalue.corrected==TRUE){ + ss.reg <- ss.total.cor - ss.res # Regression Sum of Squares, DF= (n-1)-(n-k) = k-1 in this case + dfR= k-1 + }else{ + ss.reg <- ss.total.uncor - ss.res # Regression Sum of Squares, DF= n-(n-k) = k in this case + dfR= k + } + + dfE= n-k # degrees of freedom for Error (or Residuals) + + Fval=(ss.reg/dfR)/(ss.res/dfE) + pval=pf(Fval,dfR,dfE,lower.tail = F) + pval<-unname(pval) + + RSE<-sum.exp3P$sigma # Residual standard error, type ?summary.nls in R for more detials. + SSE<-(RSE^2)*(n-1) # Sum of Squares for Error, not equal to 'ss.res'. + + adjr2 <- 1-SSE/((var(y))*(n-1)) + r2<-1-(1-adjr2)*((n-k)/(n-1)) + + if (summary==TRUE){ + ### Start print step by step + #re-adjust the output of coefficients. + coeffadj = coef(sum.exp3P) + ab.param<-coeffadj[1:2,] + # re-adjust the Estimate value of parameter c + c.param<-coeffadj[3,] + c.p1<-c.param[1] + c.p1 = c.p1 + min(y)-1 # re-adjust 'Estimate' value + c.se<-c.param[2] # Std.Error value + c.tval<-c.p1/c.se #re-adjust 't-value' + c.pval<-2 * pt(abs(c.tval), n-k, lower.tail = FALSE) #re-adjust 'p-value' + c<-c(c.p1,c.se,c.tval,c.pval) # re-adjust + coeff.re.adj<- rbind(ab.param,c) + + # print + cat("\nNonlinear regression model\n") + cat("\nFormula: y = a*exp(b*x) + c","\n") + df <- sum.exp3P$df + rdf <- df[2L] + cat("\nParameters:\n") + printCoefmat(coeff.re.adj, digits = eDigit) + cat("\nResidual standard error:", + format(sum.exp3P$sigma, digits = eDigit), "on", rdf, "degrees of freedom","\n") + + convInfo = fit$convInfo + iterations<-convInfo$finIter + tolerance<-convInfo$finTol + + cat("\nNumber of iterations to convergence:", + format(iterations, digits = eDigit)) + cat("\nAchieved convergence tolerance:",format(tolerance, digits = eDigit),"\n") + + cat("\nMultiple R-squared:", + format(r2, digits = eDigit), ", Adjusted R-squared: ", + format(adjr2, digits = eDigit)) + cat("\nF-statistic:", + format(Fval, digits = eDigit), "on", dfR, "and", dfE, "DF, ", "p-value:", format(pval, digits = eDigit), "\n") + ### finished print + }else{} + + coeff<-sum.exp3P$coefficients + a<-coeff[1,1] + b<-coeff[2,1] + c<-coeff[3,1]+min(y)-1 + + a = format(a, digits = eDigit) + b = format(b, digits = eDigit) + c = format(c, digits = eDigit) + r2 = format(r2, digits = eDigit) + adjr2 = format(adjr2, digits = eDigit) + pval= format(pval, digits = eDigit) + + a=as.numeric(a) + b=as.numeric(b) + c=as.numeric(c) + param.out<- c(list("a"=a,"b"=b,"c"=c)) + + } + + + # 5.2) model="power2P" + if (model== c("power2P")) + { + formula = 'y = a*x^b' + + n<-length(x) + k = 2 # k means the count numbers of parameters (i.e., a, b and c in this case) + + if (min(x)>0){ + # use selfStart function 'SSpower2P' for y = a *x^b + # trendline model + fit<-nls(y ~ SSpower2P(x,a,b),data=z) + sum.power2P <- summary(fit) + + ### calculate the F-statistic and p-value for model + ss.res<-sum((residuals(fit))^2) # Residual Sum of Squares, DF= n-k + ss.total.uncor<-sum(y^2) # Uncorrected Total Sum of Squares, DF=n + ss.total.cor<-sum((y-mean(y))^2) # Corrected Total Sum of Squares, DF=n-1 + + if (Pvalue.corrected==TRUE){ + ss.reg <- ss.total.cor - ss.res # Regression Sum of Squares, DF= (n-1)-(n-k) = k-1 in this case + dfR= k-1 + }else{ + ss.reg <- ss.total.uncor - ss.res # Regression Sum of Squares, DF= n-(n-k) = k in this case + dfR= k + } + + dfE = n-k # degrees of freedom for Error (or Residuals) + + Fval = (ss.reg/dfR)/(ss.res/dfE) + pval = pf(Fval,dfR,dfE,lower.tail = F) + pval <- unname(pval) + + RSE<-sum.power2P$sigma # Residual standard error, type ?summary.nls in R for more detials. + SSE<-(RSE^2)*(n-1) # Sum of Squares for Error, not equal to 'ss.res'. + + adjr2 <- 1-SSE/((var(y))*(n-1)) + r2 <- 1-(1-adjr2)*((n-k)/(n-1)) + + if (summary==TRUE){ + + ### Start print step by step + coeff = coef(sum.power2P) + ab.param<-coeff[1:2,] + + # print + cat("\nNonlinear regression model\n") + cat("\nFormula: y = a*x^b","\n") + df <- sum.power2P$df + rdf <- df[2L] + cat("\nParameters:\n") + printCoefmat(coeff, digits = eDigit) + cat("\nResidual standard error:", + format(sum.power2P$sigma, digits = eDigit), "on", rdf, "degrees of freedom","\n") + + convInfo = fit$convInfo + iterations<-convInfo$finIter + tolerance<-convInfo$finTol + + cat("\nNumber of iterations to convergence:", + format(iterations, digits = eDigit)) + cat("\nAchieved convergence tolerance:",format(tolerance, digits = eDigit),"\n") + + cat("\nMultiple R-squared:", + format(r2, digits = eDigit), ", Adjusted R-squared: ", + format(adjr2, digits = eDigit)) + cat("\nF-statistic:", + format(Fval, digits = eDigit), "on", dfR, "and", dfE, "DF, ", "p-value:", format(pval, digits = eDigit), "\n") + ### finished print + }else{} + + coeffs<-sum.power2P$coefficients + a<-coeffs[1,1] + b<-coeffs[2,1] + + a = format(a, digits = eDigit) + b = format(b, digits = eDigit) + + r2 = format(r2, digits = eDigit) + adjr2 = format(adjr2, digits = eDigit) + pval = format(pval, digits = eDigit) + + a=as.numeric(a) + b=as.numeric(b) + param.out<- c(list("a"=a,"b"=b)) + }else{ + stop(" + 'power2P' model need ALL x values greater than 0. Try other models.") + } + } + + + # 5.3) model="power3P" + if (model== c("power3P")) + { + formula = 'y = a*x^b + c' + + yadj<-y-min(y)+1 + zzz<-data.frame(x,yadj) + + n<-length(x) + k = 3 # k means the count numbers of parameters (i.e., a, b and c in this case) + + if (min(x)>0){ + # use selfStart function 'SSpower3P' for y = a *x^b+ c + # trendline model + fit<-nls(yadj~SSpower3P(x,a,b,c),data=zzz) # use 'yadj', in case of extreme high y-values with low range. + sum.power3P <- summary(fit) + + ### calculate the F-statistic and p-value for model + ss.res<-sum((residuals(fit))^2) # Residual Sum of Squares, DF= n-k + ss.total.uncor<-sum(y^2) # Uncorrected Total Sum of Squares, DF=n + ss.total.cor<-sum((y-mean(y))^2) # Corrected Total Sum of Squares, DF=n-1 + + if (Pvalue.corrected==TRUE){ + ss.reg <- ss.total.cor - ss.res # Regression Sum of Squares, DF= (n-1)-(n-k) = k-1 in this case + dfR= k-1 + }else{ + ss.reg <- ss.total.uncor - ss.res # Regression Sum of Squares, DF= n-(n-k) = k in this case + dfR= k + } + + dfE= n-k # degrees of freedom for Error (or Residuals) + + Fval=(ss.reg/dfR)/(ss.res/dfE) + pval=pf(Fval,dfR,dfE,lower.tail = F) + pval<-unname(pval) + + RSE<-sum.power3P$sigma # Residual standard error, type ?summary.nls in R for more detials. + SSE<-(RSE^2)*(n-1) # Sum of Squares for Error, not equal to 'ss.res'. + + adjr2 <- 1-SSE/((var(y))*(n-1)) + r2<-1-(1-adjr2)*((n-k)/(n-1)) + + if (summary==TRUE){ + + ### Start print step by step + #re-adjust the output of coefficients. + coeffadj = coef(sum.power3P) + ab.param<-coeffadj[1:2,] + # re-adjust the 'Estimate\ value of parameter 'c' + c.param<-coeffadj[3,] + c.p1<-c.param[1] + c.p1 = c.p1 + min(y)-1 # re-adjust + c.se<-c.param[2] # Std.Error value + c.tval<-c.p1/c.se #re-adjust 't-value' + c.pval<-2 * pt(abs(c.tval), n-k, lower.tail = FALSE) #re-adjust 'p-value' + c<-c(c.p1,c.se,c.tval,c.pval) # re-adjust + coeff.re.adj<- rbind(ab.param,c) + + # print + cat("\nNonlinear regression model\n") + cat("\nFormula: y = a*x^b + c","\n") + df <- sum.power3P$df + rdf <- df[2L] + cat("\nParameters:\n") + printCoefmat(coeff.re.adj, digits = eDigit) + cat("\nResidual standard error:", + format(sum.power3P$sigma, digits = eDigit), "on", rdf, "degrees of freedom","\n") + + convInfo = fit$convInfo + iterations<-convInfo$finIter + tolerance<-convInfo$finTol + + cat("\nNumber of iterations to convergence:", + format(iterations, digits = eDigit)) + cat("\nAchieved convergence tolerance:",format(tolerance, digits = eDigit),"\n") + + cat("\nMultiple R-squared:", + format(r2, digits = eDigit), ", Adjusted R-squared: ", + format(adjr2, digits = eDigit)) + cat("\nF-statistic:", + format(Fval, digits = eDigit), "on", dfR, "and", dfE, "DF, ", "p-value:", format(pval, digits = eDigit), "\n") + ### finished print + }else{} + + coeff<-sum.power3P$coefficients + a<-coeff[1,1] + b<-coeff[2,1] + c<-coeff[3,1] + c<-c+min(y)-1 #re-adjust + + a = format(a, digits = eDigit) + b = format(b, digits = eDigit) + c = format(c, digits = eDigit) + r2 = format(r2, digits = eDigit) + adjr2 = format(adjr2, digits = eDigit) + pval = format(pval, digits = eDigit) + + a=as.numeric(a) + b=as.numeric(b) + c=as.numeric(c) + param.out<- c(list("a"=a,"b"=b,"c"=c)) + }else{ + stop(" + 'power3P' model need ALL x values greater than 0. Try other models.") + } + + # 100) beyond the built-in models. + + }else{ + Check<-c("line2P","line3P","log2P","exp2P","exp3P","power2P","power3P") + if (!model %in% Check) + stop(" + \"model\" should be one of c(\"lin2P\",\"line3P\",\"log2P\",\"exp2P\",\"exp3P\",\"power2P\",\"power3P\".") + } + + nrow = as.numeric(nrow) + r2=as.numeric(r2) + adjr2=as.numeric(adjr2) + pval=as.numeric(pval) + AIC = as.numeric(format(AIC(fit), digits = eDigit)) + BIC = as.numeric(format(BIC(fit), digits = eDigit)) + ss.res=as.numeric(format(ss.res, digits = eDigit)) + + if (summary==TRUE){ + ##print N, AIC, BIC and RSS + cat("\nN:", nrow, ", AIC:", AIC, ", BIC: ", BIC, "\nResidual Sum of Squares: ", ss.res,"\n") + }else{} + + invisible(list(formula=formula, parameter=param.out, R.squared=r2, adj.R.squared=adjr2, p.value = pval, + N = nrow, AIC=AIC, BIC=BIC, RSS=ss.res)) + } diff --git a/README.md b/README.md new file mode 100644 index 0000000..e3d747b --- /dev/null +++ b/README.md @@ -0,0 +1,298 @@ +# basicTrendline: an R package for adding trendline of basic regression models to plot + + + +[![cran version](http://www.r-pkg.org/badges/version/basicTrendline)](http://cran.rstudio.com/web/packages/basicTrendline) +[![rstudio mirror downloads](http://cranlogs.r-pkg.org/badges/grand-total/basicTrendline)](https://github.com/metacran/cranlogs.app) +[![rstudio mirror downloads](http://cranlogs.r-pkg.org/badges/basicTrendline)](https://github.com/metacran/cranlogs.app) +[![HitCount](http://hits.dwyl.io/PhDMeiwp/basicTrendline.svg)](http://hits.dwyl.io/PhDMeiwp/basicTrendline) + + +## Authors + + + +Weiping MEI https://PhDMeiwp.github.io + + +Graduate School of Fisheries and Environmental Sciences, Nagasaki University + +## Citation + +Mei W, Yu G, Lai J, Rao Q, Umezawa Y (2018) basicTrendline: Add Trendline and Confidence Interval of Basic Regression Models to Plot. R package version 2.0.3. http://CRAN.R-project.org/package=basicTrendline + +## Installation + +Get the released version from CRAN: + + install.packages("basicTrendline") + +Or the development version from github: + + install.packages("devtools") + devtools::install_github("PhDMeiwp/basicTrendline@master", force = TRUE) + + +## Changes in version 2.0.3 + +- add several arguments to `trendline()` function, including show.equation, show.Rpvalue, Rname, Pname, xname, yname, yhat, CI.fill, CI.level, CI.alpha, CI.color, CI.lty, CI.lwd, ePos.x, ePos.y, las. +- enable to draw confidence interval for regression models (arguments CI.fill, CI.level, etc.) +- add 'show.equation' and show.Rpvale' arguments to enable to choose which parameter to show +- add 'Rname' and 'Pname' arguments to specify the character of R-square and P-vlaue (i.e. R^2 or r^2; P or p) +- add 'xname' and 'ynameto' arguments to specify the character of 'x' and 'y' in the equation +- add 'yhat' argument to enable to add a hat symbol on the top of 'y' in the equation +- add 'ePos.x' and 'ePos.y' arguments to specify the x and y co-ordinates of equation's position +- deleted the 'ePos' argument +- add "Residual Sum of Squares" to the output of 'trendline_summary()' function + +## Changes in version 1.2.0 + +- change the expression for `model` of `exp3P` using a supscript +- add `trendline_summary()` function +- add `SSexp2P()` function +- add `SSpower2P` function +- add `Pvalue.corrected` argument in `trendline()` and `trendline_summary()` , for P-vlaue calculation for non-linear regression. +- add `Details` in `trendline()` and `trendline_summary()` +- add `...` argument in `trendline()` as same as those in `plot()` + +--- + +# Examples + + library(basicTrendline) + x <- c(1, 3, 6, 9, 13, 17) + y <- c(5, 8, 11, 13, 13.2, 13.5) + +# [case 1] default + trendline(x, y, model="line2P", ePos.x = "topleft", summary=TRUE, eDigit=5) + + + +# [case 2] draw lines of confidenc interval only (set CI.fill = FALSE) + trendline(x, y, model="line3P", CI.fill = FALSE, CI.color = "black", CI.lty = 2, linecolor = "blue") + + + +# [case 3] draw trendliine only (set CI.color = NA) + trendline(x, y, model="log2P", ePos.x= "top", linecolor = "red", CI.color = NA) + + + +# [case 4] show regression equation only (set show.Rpvalue = FALSE) + trendline(x, y, model="exp2P", show.equation = TRUE, show.Rpvalue = FALSE) + + + +# [case 5] specify the name of parameters in equation +** see Arguments c('xname', 'yname', 'yhat', 'Rname', 'Pname') ** + trendline(x, y, model="exp3P", xname="T", yname=paste(delta^15,"N"), + yhat=FALSE, Rname=1, Pname=0, ePos.x = "bottom") + + + +# [case 6] change the digits, font size, and color of equation. + trendline(x, y, model="power2P", eDigit = 3, eSize = 1.4, text.col = "blue") + + + +# [case 7] don't show equation (set ePos.x = NA) + trendline(x, y, model="power3P", ePos.x = NA) + + + +# [case 8] set graphical parameters by par {graphics} + ### NOT RUN + par(mgp=c(1.5,0.4,0), mar=c(3,3,1,1), tck=-0.01, cex.axis=0.9) + + trendline(x, y) + + dev.off() + + ### END (NOT RUN) + + + +--- + +## Description + +Plot, draw regression line and confidence interval, and show regression equation, R-square and P-value, as simple as possible, + +by using different models built in the 'trendline()' function. The function includes the following models in the latest version: + +"line2P" (formula as: y=a\*x+b), + +"line3P" (y=a\*x2+b\*x+c), + +"log2P" (y=a\*ln(x)+b), + +"exp2P" (y=a\*eb\*x), + +"exp3P" (y=a\*eb\*x+c), + +"power2P" (y=a\*xb), + +"power3P" (y=a\*xb+c). + +Besides, the summarized results of each fitted model are also output by default. + +## Usage + + trendline(x, y, model = "line2P", Pvalue.corrected = TRUE, + linecolor = "blue", lty = 1, lwd = 1, + show.equation = TRUE, show.Rpvalue = TRUE, + Rname = 1, Pname = 0, xname = "x", yname = "y", + yhat = FALSE, + summary = TRUE, + ePos.x = NULL, ePos.y = NULL, text.col = "black", eDigit = 5, eSize = 1, + CI.fill = TRUE, CI.level = 0.95, CI.color = "grey", CI.alpha = 1, CI.lty = 1, CI.lwd = 1, + las = 1, xlab = NULL, ylab = NULL, ...) + +## Arguments + +
**x, y**
+the x and y arguments provide the x and y coordinates for the plot. Any reasonable way of defining the coordinates is acceptable. + +
**model**
+select which model to fit. Default is "line2P". The "model" should be one of c("line2P", "line3P", "log2P", "exp3P", "power3P"), their formulas are as follows: +
"line2P": y=a\*x+b +
"line3P": y=a\*x2+b\*x+c +
"log2P": y=a\*ln(x)+b +
"exp2P": y=a\*eb\*x +
"exp3P": y=a\*eb\*x+c +
"power2P": y=a\*xb +
"power3P": y=a\*xb+c + +
**Pvalue.corrected**
+if P-value corrected or not, the vlaue is one of c("TRUE", "FALSE"). + +
**linecolor**
+color of regression line. + +
**lty**
+line type. lty can be specified using either text c("blank","solid","dashed","dotted","dotdash","longdash","twodash") or number c(0, 1, 2, 3, 4, 5, 6). Note that lty = "solid" is identical to lty=1. + +
**lwd**
+line width. Default is 1. + +
**show.equation**
+whether to show the regression equation, the value is one of c("TRUE", "FALSE"). + +
**show.Rpvalue**
+whether to show the R-square and P-value, the value is one of c("TRUE", "FALSE"). + +
**Rname**
+to specify the character of R-square, the value is one of c(0, 1), corresponding to c(r^2, R^2). + +
**Pname**
+to specify the character of P-value, the value is one of c(0, 1), corresponding to c(p, P). + +
**xname**
+to specify the character of "x" in equation, see Examples [case 5]. + +
**yname**
+to specify the character of "y" in equation, see Examples [case 5]. + +
**yhat**
+whether to add a hat symbol (^) on the top of "y" in equation. Default is FALSE. + +
**summary**
+summarizing the model fits. Default is TRUE. + +
**ePos.x, ePos.y**
+equation position. Default as ePos.x = "topleft". If no need to show equation, set ePos.x = NA. It's same as those in legend. + +
**text.col**
+the color used for the legend text. + +
**eDigit**
+the numbers of digits for equation parameters. Default is 5. + +
**eSize**
+font size in percentage of equation. Default is 1. + +
**CI.fill**
+fill the confidance interval? (TRUE by default, see 'CI.level' to control) + +
**CI.level**
+level of confidence interval to use (0.95 by default) + +
**CI.color**
+line or fill color of confidence interval. + +
**CI.alpha**
+alpha value of fill color of confidence interval. + +
**CI.lty**
+line type of confidence interval. + +
**CI.lwd**
+line width of confidence interval. + +
**las**
+style of axis labels. (0=parallel, 1=all horizontal, 2=all perpendicular to axis, 3=all vertical) + +
**xlab, ylab**
+labels of x- and y-axis. + +**...**
+additional parameters to plot,such as type, main, sub, xlab, ylab, col. + +## Details + +The linear models (line2P, line3P, log2P) in this package are estimated by **lm** function, while the **nonlinear models (exp2P, exp3P, power2P, power3P)** are estimated by **nls** function (i.e., **least-squares method**). +
The argument 'Pvalue.corrected' is workful for non-linear regression only. +
If "Pvalue.corrected = TRUE", the P-vlaue is calculated by using "Residual Sum of Squares" and "Corrected Total Sum of Squares (i.e. sum((y-mean(y))^2))". +
If "Pvalue.corrected = TRUE", the P-vlaue is calculated by using "Residual Sum of Squares" and "Uncorrected Total Sum of Squares (i.e. sum(y^2))". + +## Note + +Confidence intervals for nonlinear regression (i.e., objects of class nls) are based on the linear approximation described in Bates & Watts (2007). + +## References + +Bates, D. M., and Watts, D. G. (2007) *Nonlinear Regression Analysis and its Applications*. Wiley. + +Greenwell B. M., and Schubert-Kabban, C. M. (2014) *investr: An R Package for Inverse Estimation*. The R Journal, 6(1), 90-100. + +## Value + +R2, indicates the R-Squared value of each regression model. + +p, indicates the p-value of each regression model. + +AIC or BIC, indicate the Akaike's Information Criterion or Bayesian Information Criterion for fitted model. Click AIC for details. The smaller the AIC or BIC, the better the model. + +RSS, indicates the "Residual Sum of Squares” of regression model. + +--- + +To see examples on how to use "basicTrendline" in R software, you can run the following R code if you have the "basicTrendline" package installed: + + library(basicTrendline) + ?trendline() + + +## Acknowledgements + +We would like to express my special thanks to **Uwe Ligges, Swetlana Herbrandt, and CRAN team** for their very valuable comments to the 'basicTrendline' package. +Our thanks also go to those who contributed R codes by: + +- adding conficende interval for both lm and nls objects: https://github.com/bgreenwell/investr +- adding-regression-line-equation-and-r2-on-graph-1: http://blog.sciencenet.cn/blog-267448-1021594.html +- adding-regression-line-equation-and-r2-on-graph-2: https://stackoverflow.com/questions/7549694/adding-regression-line-equation-and-r2-on-graph +- What is non-linear regression?: https://datascienceplus.com/first-steps-with-non-linear-regression-in-r/ +- adding regression line for nonlinear regression: http://blog.sciencenet.cn/blog-651374-1014133.html +- R codes for 'print.summary.nls of exp3P and power3P' cite from https://github.com/SurajGupta/r-source/blob/master/src/library/stats/R/nls.R + +## Contact + +- If you have any question or comment to this package, tell me at [here](http://meiweiping.cn/en/basicTrendline-an-R-package-for-adding-trendline-of-basic-regression-models-to-plot/). + +- Bugs and feature requests can be filed to https://github.com/PhDMeiwp/basicTrendline/issues. BTW, [Pull requests](https://github.com/PhDMeiwp/basicTrendline/pulls) are also welcome. + +## Appendix + +The **PDF file** of this R package is available at https://cran.r-project.org/web/packages/basicTrendline/index.html + +> 点击进入 [basicTrendline函数包中文介绍入口](http://meiweiping.cn/%E7%94%A8%E4%BA%8E%E5%B8%B8%E8%A7%84%E7%BA%BF%E6%80%A7%E9%9D%9E%E7%BA%BF%E6%80%A7%E6%8B%9F%E5%90%88%E7%9A%84R%E5%87%BD%E6%95%B0%E5%8C%85%EF%BC%88basicTrendline%EF%BC%89%E4%BB%8B%E7%BB%8D/) \ No newline at end of file diff --git a/basicTrendline.Rproj b/basicTrendline.Rproj new file mode 100644 index 0000000..cba1b6b --- /dev/null +++ b/basicTrendline.Rproj @@ -0,0 +1,21 @@ +Version: 1.0 + +RestoreWorkspace: No +SaveWorkspace: No +AlwaysSaveHistory: Default + +EnableCodeIndexing: Yes +UseSpacesForTab: Yes +NumSpacesForTab: 2 +Encoding: UTF-8 + +RnwWeave: Sweave +LaTeX: pdfLaTeX + +AutoAppendNewline: Yes +StripTrailingWhitespace: Yes + +BuildType: Package +PackageUseDevtools: Yes +PackageInstallArgs: --no-multiarch --with-keep.source +PackageRoxygenize: rd,collate,namespace diff --git a/docs/images/basicTrendline.hex.png b/docs/images/basicTrendline.hex.png new file mode 100644 index 0000000..5f91a38 Binary files /dev/null and b/docs/images/basicTrendline.hex.png differ diff --git a/docs/images/case1.png b/docs/images/case1.png new file mode 100644 index 0000000..c5e2961 Binary files /dev/null and b/docs/images/case1.png differ diff --git a/docs/images/case2.png b/docs/images/case2.png new file mode 100644 index 0000000..20e09a6 Binary files /dev/null and b/docs/images/case2.png differ diff --git a/docs/images/case3.png b/docs/images/case3.png new file mode 100644 index 0000000..21e784e Binary files /dev/null and b/docs/images/case3.png differ diff --git a/docs/images/case4.png b/docs/images/case4.png new file mode 100644 index 0000000..4ef93fb Binary files /dev/null and b/docs/images/case4.png differ diff --git a/docs/images/case5.png b/docs/images/case5.png new file mode 100644 index 0000000..40eb46e Binary files /dev/null and b/docs/images/case5.png differ diff --git a/docs/images/case6.png b/docs/images/case6.png new file mode 100644 index 0000000..de842f2 Binary files /dev/null and b/docs/images/case6.png differ diff --git a/docs/images/case7.png b/docs/images/case7.png new file mode 100644 index 0000000..a15f7f2 Binary files /dev/null and b/docs/images/case7.png differ diff --git a/docs/images/case8.png b/docs/images/case8.png new file mode 100644 index 0000000..f64e4dd Binary files /dev/null and b/docs/images/case8.png differ diff --git a/man/SSexp2P.Rd b/man/SSexp2P.Rd new file mode 100644 index 0000000..af45c4f --- /dev/null +++ b/man/SSexp2P.Rd @@ -0,0 +1,34 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/SSexp2P.R +\name{SSexp2P} +\alias{SSexp2P} +\title{Self-Starting Nls 'exp2P' Regression Model} +\usage{ +SSexp2P(predictor, a, b) +} +\arguments{ +\item{predictor}{a numeric vector of values at which to evaluate the model.} + +\item{a, b}{The numeric parameters responsing to the exp2P model.} +} +\description{ +This selfStart model evaluates the power regression function (formula as: y=a*exp(b*x)). It has an initial attribute that will evaluate initial estimates of the parameters 'a' and 'b' for a given set of data. +} +\examples{ +library(basicTrendline) +x<-1:5 +y<-c(2,4,8,20,25) +xy<-data.frame(x,y) +getInitial(y ~ SSexp2P(x,a,b), data = xy) +## Initial values are in fact the converged values + +fitexp2P <- nls(y~SSexp2P(x,a,b), data=xy) +summary(fitexp2P) + +} +\seealso{ +\code{\link{trendline}}, \code{\link{SSexp3P}}, \code{\link{SSpower3P}}, \code{\link[stats]{nls}}, \code{\link[stats]{selfStart}} +} +\author{ +Weiping Mei \email{meiweipingg@163.com} +} diff --git a/man/SSexp3P.Rd b/man/SSexp3P.Rd new file mode 100644 index 0000000..72e8de8 --- /dev/null +++ b/man/SSexp3P.Rd @@ -0,0 +1,34 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/SSexp3P.R +\name{SSexp3P} +\alias{SSexp3P} +\title{Self-Starting Nls 'exp3P' Regression Model} +\usage{ +SSexp3P(predictor, a, b, c) +} +\arguments{ +\item{predictor}{a numeric vector of values at which to evaluate the model.} + +\item{a, b, c}{Three numeric parameters responsing to the exp3P model.} +} +\description{ +This selfStart model evaluates the exponential regression function (formula as: y=a*exp(b*x)+c). It has an initial attribute that will evaluate initial estimates of the parameters a, b, and c for a given set of data. +} +\examples{ +library(basicTrendline) +x<-1:5 +y<-c(2,4,8,16,28) +xy<-data.frame(x,y) +getInitial(y ~ SSexp3P(x,a,b,c), data = xy) +## Initial values are in fact the converged values + +fitexp3P <- nls(y~SSexp3P(x,a,b,c), data=xy) +summary(fitexp3P) + +} +\seealso{ +\code{\link{trendline}}, \code{\link{SSexp3P}}, \code{\link{SSpower3P}}, \code{\link[stats]{nls}}, \code{\link[stats]{selfStart}} +} +\author{ +Weiping Mei \email{meiweipingg@163.com} +} diff --git a/man/SSpower2P.Rd b/man/SSpower2P.Rd new file mode 100644 index 0000000..146b439 --- /dev/null +++ b/man/SSpower2P.Rd @@ -0,0 +1,34 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/SSpower2P.R +\name{SSpower2P} +\alias{SSpower2P} +\title{Self-Starting Nls 'power2P' Regression Model} +\usage{ +SSpower2P(predictor, a, b) +} +\arguments{ +\item{predictor}{a numeric vector of values at which to evaluate the model.} + +\item{a, b}{The numeric parameters responsing to the exp2P model.} +} +\description{ +This selfStart model evaluates the power regression function (formula as: y=a*x^b). It has an initial attribute that will evaluate initial estimates of the parameters 'a' and 'b' for a given set of data. +} +\examples{ +library(basicTrendline) +x<-1:5 +y<-c(2,4,8,20,25) +xy<-data.frame(x,y) +getInitial(y ~ SSpower2P(x,a,b), data = xy) +## Initial values are in fact the converged values + +fitpower2P <- nls(y~SSpower2P(x,a,b), data=xy) +summary(fitpower2P) + +} +\seealso{ +\code{\link{trendline}}, \code{\link{SSexp3P}}, \code{\link{SSpower3P}}, \code{\link[stats]{nls}}, \code{\link[stats]{selfStart}} +} +\author{ +Weiping Mei \email{meiweipingg@163.com} +} diff --git a/man/SSpower3P.Rd b/man/SSpower3P.Rd new file mode 100644 index 0000000..d300c9f --- /dev/null +++ b/man/SSpower3P.Rd @@ -0,0 +1,34 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/SSpower3P.R +\name{SSpower3P} +\alias{SSpower3P} +\title{Self-Starting Nls 'power3P' Regression Model} +\usage{ +SSpower3P(predictor, a, b, c) +} +\arguments{ +\item{predictor}{a numeric vector of values at which to evaluate the model.} + +\item{a, b, c}{Three numeric parameters responsing to the exp3P model.} +} +\description{ +This selfStart model evaluates the power regression function (formula as: y=a*x^b+c). It has an initial attribute that will evaluate initial estimates of the parameters a, b, and c for a given set of data. +} +\examples{ +library(basicTrendline) +x<-1:5 +y<-c(2,4,8,20,25) +xy<-data.frame(x,y) +getInitial(y ~ SSpower3P(x,a,b,c), data = xy) +## Initial values are in fact the converged values + +fitpower3P <- nls(y~SSpower3P(x,a,b,c), data=xy) +summary(fitpower3P) + +} +\seealso{ +\code{\link{trendline}}, \code{\link{SSexp3P}}, \code{\link{SSpower3P}}, \code{\link[stats]{nls}}, \code{\link[stats]{selfStart}} +} +\author{ +Weiping Mei \email{meiweipingg@163.com} +} diff --git a/man/trendline.Rd b/man/trendline.Rd new file mode 100644 index 0000000..b07fa20 --- /dev/null +++ b/man/trendline.Rd @@ -0,0 +1,134 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/trendline.R +\name{trendline} +\alias{trendline} +\title{Add Trendline and Show Equation to Plot} +\usage{ +trendline(x, y, model = "line2P", Pvalue.corrected = TRUE, + linecolor = "blue", lty = 1, lwd = 1, show.equation = TRUE, + show.Rpvalue = TRUE, Rname = 1, Pname = 0, xname = "x", yname = "y", + yhat = FALSE, summary = TRUE, ePos.x = NULL, ePos.y = NULL, + text.col = "black", eDigit = 5, eSize = 1, CI.fill = TRUE, + CI.level = 0.95, CI.color = "grey", CI.alpha = 1, CI.lty = 1, + CI.lwd = 1, las = 1, xlab = NULL, ylab = NULL, ...) +} +\arguments{ +\item{x, y}{the x and y arguments provide the x and y coordinates for the plot. Any reasonable way of defining the coordinates is acceptable.} + +\item{model}{select which model to fit. Default is "line2P". The "model" should be one of c("line2P", "line3P", "log2P", "exp2P", "exp3P", "power2P", "power3P"), their formulas are as follows:\cr "line2P": y=a*x+b \cr "line3P": y=a*x^2+b*x+c \cr "log2P": y=a*ln(x)+b \cr "exp2P": y=a*exp(b*x) \cr "exp3P": y=a*exp(b*x)+c \cr "power2P": y=a*x^b \cr "power3P": y=a*x^b+c} + +\item{Pvalue.corrected}{if P-value corrected or not, the value is one of c("TRUE", "FALSE").} + +\item{linecolor}{color of regression line.} + +\item{lty}{line type. lty can be specified using either text c("blank","solid","dashed","dotted","dotdash","longdash","twodash") or number c(0, 1, 2, 3, 4, 5, 6). Note that lty = "solid" is identical to lty=1.} + +\item{lwd}{line width. Default is 1.} + +\item{show.equation}{whether to show the regression equation, the value is one of c("TRUE", "FALSE").} + +\item{show.Rpvalue}{whether to show the R-square and P-value, the value is one of c("TRUE", "FALSE").} + +\item{Rname}{to specify the character of R-square, the value is one of c(0, 1), corresponding to c(r^2, R^2).} + +\item{Pname}{to specify the character of P-value, the value is one of c(0, 1), corresponding to c(p, P).} + +\item{xname}{to specify the character of "x" in equation, see Examples [case 5].} + +\item{yname}{to specify the character of "y" in equation, see Examples [case 5].} + +\item{yhat}{whether to add a hat symbol (^) on the top of "y" in equation. Default is FALSE.} + +\item{summary}{summarizing the model fits. Default is TRUE.} + +\item{ePos.x, ePos.y}{equation position. Default as ePos.x = "topleft". If no need to show equation, set ePos.x = NA. It's same as those in \code{\link[graphics]{legend}}.} + +\item{text.col}{the color used for the equation text.} + +\item{eDigit}{the numbers of digits for equation parameters. Default is 5.} + +\item{eSize}{font size in percentage of equation. Default is 1.} + +\item{CI.fill}{fill the confidance interval? (TRUE by default, see 'CI.level' to control)} + +\item{CI.level}{level of confidence interval to use (0.95 by default)} + +\item{CI.color}{line or fill color of confidence interval.} + +\item{CI.alpha}{alpha value of fill color of confidence interval.} + +\item{CI.lty}{line type of confidence interval.} + +\item{CI.lwd}{line width of confidence interval.} + +\item{las}{style of axis labels. (0=parallel, 1=all horizontal, 2=all perpendicular to axis, 3=all vertical)} + +\item{xlab, ylab}{labels of x- and y-axis.} + +\item{...}{additional parameters to \code{\link[graphics]{plot}}, such as type, main, sub, pch, col.} +} +\description{ +Plot, draw regression line and confidence interval, and show regression equation, R-square and P-value, as simple as possible, +by using different models built in the 'trendline()' function. The function includes the following models in the latest version: +"line2P" (formula as: y=a*x+b), "line3P" (y=a*x^2+b*x+c), "log2P" (y=a*ln(x)+b), "exp2P" (y=a*exp(b*x)),"exp3P" (y=a*exp(b*x)+c), "power2P" (y=a*x^b), and "power3P" (y=a*x^b+c). +Besides, the summarized result of each fitted model is also output by default. +} +\details{ +The linear models (line2P, line3P, log2P) in this package are estimated by \code{\link[stats]{lm}} function, \cr while the nonlinear models (exp2P, exp3P, power2P, power3P) are estimated by \code{\link[stats]{nls}} function (i.e., least-squares method).\cr\cr The argument 'Pvalue.corrected' is workful for non-linear regression only.\cr\cr If "Pvalue.corrected = TRUE", the P-value is calculated by using "Residual Sum of Squares" and "Corrected Total Sum of Squares (i.e. sum((y-mean(y))^2))".\cr If "Pvalue.corrected = TRUE", the P-value is calculated by using "Residual Sum of Squares" and "Uncorrected Total Sum of Squares (i.e. sum(y^2))". +} +\note{ +Confidence intervals for nonlinear regression (i.e., objects of class +\code{nls}) are based on the linear approximation described in Bates & Watts (2007) and Greenwell & Schubert-Kabban (2014). +} +\examples{ +library(basicTrendline) +x <- c(1, 3, 6, 9, 13, 17) +y <- c(5, 8, 11, 13, 13.2, 13.5) + +# [case 1] default +trendline(x, y, model="line2P", ePos.x = "topleft", summary=TRUE, eDigit=5) +# [case 2] draw lines of confidenc interval only (set CI.fill = FALSE) +trendline(x, y, model="line3P", CI.fill = FALSE, CI.color = "black", CI.lty = 2, linecolor = "blue") + +# [case 3] draw trendliine only (set CI.color = NA) +trendline(x, y, model="log2P", ePos.x= "top", linecolor = "red", CI.color = NA) + +# [case 4] show regression equation only (set show.Rpvalue = FALSE) +trendline(x, y, model="exp2P", show.equation = TRUE, show.Rpvalue = FALSE) + +# [case 5] specify the name of parameters in equation +# see Arguments c('xname', 'yname', 'yhat', 'Rname', 'Pname'). +trendline(x, y, model="exp3P", xname="T", yname=paste(delta^15,"N"), + yhat=FALSE, Rname=1, Pname=0, ePos.x = "bottom") + +# [case 6] change the digits, font size, and color of equation. +trendline(x, y, model="power2P", eDigit = 3, eSize = 1.4, text.col = "blue") + +# [case 7] don't show equation (set ePos.x = NA) +trendline(x, y, model="power3P", ePos.x = NA) + +### NOT RUN +# [case 8] set graphical parameters by par {graphics} + +par(mgp=c(1.5,0.4,0), mar=c(3,3,1,1), tck=-0.01, cex.axis=0.9) + +trendline(x, y) + +dev.off() + +### END (NOT RUN) + +} +\references{ +Bates, D. M., and Watts, D. G. (2007) +\emph{Nonlinear Regression Analysis and its Applications}. Wiley. + +Greenwell B. M., and Schubert-Kabban, C. M. (2014) +\emph{investr: An R Package for Inverse Estimation}. The R Journal, 6(1), 90-100. +} +\seealso{ +\code{\link{trendline}}, \code{\link{SSexp3P}}, \code{\link{SSpower3P}}, \code{\link[stats]{nls}}, \code{\link[stats]{selfStart}}, \code{\link[investr]{plotFit}} +} +\author{ +Weiping Mei, Guangchuang Yu +} diff --git a/man/trendline_summary.Rd b/man/trendline_summary.Rd new file mode 100644 index 0000000..d7fb69a --- /dev/null +++ b/man/trendline_summary.Rd @@ -0,0 +1,65 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/trendline_summary.R +\name{trendline_summary} +\alias{trendline_summary} +\title{Summarized Results of Each Regression Model} +\usage{ +trendline_summary(x, y, model = "line2P", Pvalue.corrected = TRUE, + summary = TRUE, eDigit = 5) +} +\arguments{ +\item{x, y}{the x and y arguments provide the x and y coordinates for the plot. Any reasonable way of defining the coordinates is acceptable.} + +\item{model}{select which model to fit. Default is "line2P". The "model" should be one of c("line2P", "line3P", "log2P", "exp2P", "exp3P", "power2P", "power3P"), their formulas are as follows:\cr "line2P": y=a*x+b \cr "line3P": y=a*x^2+b*x+c \cr "log2P": y=a*ln(x)+b \cr "exp2P": y=a*exp(b*x) \cr "exp3P": y=a*exp(b*x)+c \cr "power2P": y=a*x^b \cr "power3P": y=a*x^b+c} + +\item{Pvalue.corrected}{if P-value corrected or not, the vlaue is one of c("TRUE", "FALSE").} + +\item{summary}{summarizing the model fits. Default is TRUE.} + +\item{eDigit}{the numbers of digits for summarized results. Default is 5.} +} +\value{ +R^2, indicates the R-Squared value of each regression model. + +p, indicates the p-value of each regression model. + +N, indicates the sample size. + +AIC or BIC, indicate the Akaike's Information Criterion or Bayesian Information Criterion for fitted model. Click \code{\link[stats]{AIC}} for details. The smaller the AIC or BIC, the better the model. + +RSS, indicate the value of "Residual Sum of Squares". +} +\description{ +Summarizing the results of each regression model which built in the 'trendline()' function. The function includes the following models in the latest version: +"line2P" (formula as: y=a*x+b), "line3P" (y=a*x^2+b*x+c), "log2P" (y=a*ln(x)+b), "exp2P" (y=a*exp(b*x)),"exp3P" (y=a*exp(b*x)+c), "power2P" (y=a*x^b), and "power3P" (y=a*x^b+c). +} +\details{ +The linear models (line2P, line3P, log2P) in this package are estimated by \code{\link[stats]{lm}} function, \cr while the nonlinear models (exp2P, exp3P, power2P, power3P) are estimated by \code{\link[stats]{nls}} function (i.e., least-squares method).\cr\cr The argument 'Pvalue.corrected' is workful for non-linear regression only.\cr\cr If "Pvalue.corrected = TRUE", the P-vlaue is calculated by using "Residual Sum of Squares" and "Corrected Total Sum of Squares (i.e. sum((y-mean(y))^2))".\cr If "Pvalue.corrected = TRUE", the P-vlaue is calculated by using "Residual Sum of Squares" and "Uncorrected Total Sum of Squares (i.e. sum(y^2))". +} +\examples{ +library(basicTrendline) +x1<-1:5 +x2<- -2:2 +x3<- c(101,105,140,200,660) +x4<- -5:-1 +x5<- c(1,30,90,180,360) + +y1<-c(2,14,18,19,20) # increasing convex trend +y2<- c(-2,-14,-18,-19,-20) # decreasing concave trend +y3<-c(2,4,16,38,89) # increasing concave trend +y4<-c(-2,-4,-16,-38,-89) # decreasing convex trend +y5<- c(600002,600014,600018,600019,600020) # high y values with low range. + +trendline_summary(x1,y1,model="line2P",summary=TRUE,eDigit=10) +trendline_summary(x2,y2,model="line3P",summary=FALSE) +trendline_summary(x3,y3,model="log2P") +trendline_summary(x4,y4,model="exp3P") +trendline_summary(x5,y5,model="power3P") + +} +\seealso{ +\code{\link{trendline}}, \code{\link{SSexp3P}}, \code{\link{SSpower3P}}, \code{\link[stats]{nls}}, \code{\link[stats]{selfStart}} +} +\author{ +Weiping Mei, Guangchuang Yu +}