version 1.4

cran · Mar 27, 2019 · c363720 · c363720
1 parent a8a9c84
commit c363720
Show file tree

Hide file tree

Showing 51 changed files with 5,010 additions and 313 deletions.
diff --git a/DESCRIPTION b/DESCRIPTION
@@ -1,6 +1,6 @@
 Package: QPot
-Version: 1.1
-Date: 2016-04-03
+Version: 1.4
+Date: 2019-03-24
 Title: Quasi-Potential Analysis for Stochastic Differential Equations
 Authors@R: c(
 	person("Christopher", "Moore", role = c("aut"), email = "life.dispersing@gmail.com"),
@@ -12,16 +12,18 @@ Authors@R: c(
 	)
 Depends: R (>= 3.0.2)
 Imports: MASS
+Suggests: R.rsp, R.devices,phaseR,plot3D,viridis,markdown
+VignetteBuilder: R.rsp
 Maintainer: Christopher Stieha <stieha@hotmail.com>
 Description: Tools to 1) simulate and visualize stochastic differential
     equations and 2) determine stability of equilibria using the ordered-upwind
     method to compute the quasi-potential.
 License: GPL-2
-URL: http://www.r-project.org, https://github.com/bmarkslash7/QPot
+URL: https://www.r-project.org, https://github.com/bmarkslash7/QPot
 BugReports: https://github.com/bmarkslash7/QPot/issues
 RoxygenNote: 5.0.0
 NeedsCompilation: yes
-Packaged: 2016-04-03 20:35:00 UTC; stieha
+Packaged: 2019-03-27 02:50:58 UTC; stieha
 Author: Christopher Moore [aut],
   Christopher Stieha [aut, cre],
   Ben Nolting [aut],
@@ -30,4 +32,4 @@ Author: Christopher Moore [aut],
   James Gregson [cph] (author of expression_parser library:
     https://github.com/jamesgregson/expression_parser)
 Repository: CRAN
-Date/Publication: 2016-04-04 01:07:18
+Date/Publication: 2019-03-27 10:00:20 UTC
diff --git a/MD5 b/MD5
@@ -1,40 +1,67 @@
-2a260c22eb95c9fc698e1997afd9441f *DESCRIPTION
-f2aab4d8f56b75ecbc1374060aea63e4 *NAMESPACE
-f76207f3f595bc7121e1504a94a81dd6 *NEWS
-e9d2ef44b2e0cb94ef9bc12395563b41 *R/Model2String.R
-2f57f4f8e9898f139f970eaf1f40d6bb *R/QPContour.R
-c48318359eeda0f6a430d74669621bd7 *R/QPGlobal.R
+2bcc74299996b18f050e6d7cde90bd47 *DESCRIPTION
+431aba3d96c814a6c2c221228e18a176 *NAMESPACE
+50308d8b808c8adb3fb136946379b012 *NEWS
+731e322dc1f6a1a614a0d65196c310d6 *R/Model2String.R
+bd3019a30fc82903fa732a339c901fce *R/QPContour.R
+6aa53b010d5580b5c6d23e5dc911b561 *R/QPGlobal.R
 d2d85e310c7fe29bd9a1effea8e6dfd4 *R/QPInterp.R
-0690c15bda55eba6aadbd9c595572af8 *R/QPotential.R
-c6e629e3b3a842ea2a46c4a0abe575ee *R/TSDensity.R
-1befa650d2ec75ddbc9cbe02a26c42ea *R/TSPlot.R
-054df490227be3cb2a8f333944ee5956 *R/TSTraj.R
+fd44efe68b8d6e327ff25538768abfde *R/QPotential.R
+946a3aee165f87b58dfb4e48215211b1 *R/TSDensity.R
+8734c3d0367c27b91ab6f4cf61f9c31e *R/TSPlot.R
+a95384b299ff340b2202a033993c552b *R/TSTraj.R
 bb373ebb9588b477bfd6c31c2f56c251 *R/VecDecomAll.R
 46075a3bffe1208a55673dac4558d4c8 *R/VecDecomGrad.R
-b2224cabd16a940d6ad2bc1ced7e9fb8 *R/VecDecomPlot.R
+c64c824a4209be9617c917a22785359f *R/VecDecomPlot.R
 466fe85e70ea7873a65dd39a77f4d688 *R/VecDecomRem.R
 02dbf45099e60b3d404092a22594620b *R/VecDecomVec.R
-2d3995a70b02de03e92a581a6e76922f *README.md
+3590644b98f1b1c245c8aacfb132abaf *README.md
+ad009431debbfc678b5a15384264b7bb *build/vignette.rds
 f038e36263be3e8501717d54d112627f *demo/00Index
-aaec4bec9812614d74ac4369e2009cef *demo/Ex1_test.R
-9e4d61e14bcc679a92eb290fa5a99780 *demo/Ex2_test.R
-2454bc11e4acb976da29ea8216cdb090 *demo/Ex3_test.R
-31ee8d0078f0cd237547f337363e3b43 *man/Model2String.Rd
-5b6189a092fbe2b0b09ddac3f28368cd *man/QPContour.Rd
+f85d17a0c2d4c8454ae311773d67b6fb *demo/Ex1_test.R
+0ff8c061dc757a3bd8ab21514fee261a *demo/Ex2_test.R
+49227ef52a6bd9fc1fe001d2f4a23fb0 *demo/Ex3_test.R
+89af1889d3be855e30f2ae4d733ef600 *inst/CITATION
+9f266c138ea1bd3ad6aec3b7101ed9a5 *inst/doc/QPot_Example1.R
+83d2a424eb93b9ade35c5cb3e7fed9b6 *inst/doc/QPot_Example1.html
+ced148d27dd28de8d221e0fb131dd315 *inst/doc/QPot_Example1.md.rsp
+8951794eb1ca74174648eff4b4833067 *inst/doc/QPot_Example2.R
+6b3c5a0906bcfe524c3aee428152946c *inst/doc/QPot_Example2.html
+3bddacbc5c685a79f3960785b3d6f52d *inst/doc/QPot_Example2.md.rsp
+b81af4f3e4404a6aaa6e35a656419ac1 *inst/doc/QPot_Rarticle.pdf
+6cd529b8faf6ac2753dcc22389199e23 *inst/doc/QPot_Rarticle.pdf.asis
+43a39e0d522e8aa1d6d070ba96d4940a *man/Model2String.Rd
+70a3950a11b949b0821e23e3abddb02c *man/QPContour.Rd
 5f9c14f28311635a96640a1691f7bba4 *man/QPGlobal.Rd
 329a6d792faa3a12e0e9d7c155f2fd59 *man/QPInterp.Rd
-62318ca6beffa6de894928552d3bed9b *man/QPotential.Rd
-0a56faa17f3f9a16215f09b3856aa32e *man/TSDensity.Rd
-f439ec2dc7981da54bde8520aa1ebcc5 *man/TSPlot.Rd
-8f845957abe8e93ac6b07512dfc91e00 *man/TSTraj.Rd
+90b43f02da14c4b9053227618e57cc16 *man/QPotential.Rd
+a1a5fd501beea7f26d8e0d1c8974cb12 *man/TSDensity.Rd
+9730702ae8605298d9d1cc231d66cd6e *man/TSPlot.Rd
+c6cef26dcf63b727ffa10fc8e124a0f1 *man/TSTraj.Rd
 540725662fc3b2a7b67509eec857763a *man/VecDecomAll.Rd
 c769c52bbf3c781da28538e2e556f0a4 *man/VecDecomGrad.Rd
-9734e0acfc790db09301eea4b5dd9443 *man/VecDecomPlot.Rd
+44859ef8aaaf128c2c76744eccc4d364 *man/VecDecomPlot.Rd
 47b0908dcce791cff9a0216c4a343eef *man/VecDecomRem.Rd
 f4e7ac43678d1a7b24f0836e5403475b *man/VecDecomVec.Rd
 4884e5a52bce9acf6488c38b8976ebd7 *man/figures/Example3.png
 f349d3fa85852bfefe11ef0abbb8bd11 *src/Makevars
 6ced1da2b745886b816f9ff206f2edf3 *src/Makevars.win
 7f1c9d34bf1a028f338345cf522c3a34 *src/expression_parser.c
 7022188b5b6811cc0d3512c34ab8ddda *src/expression_parser.h
-b176e7c9d494b36fd4ebc6609fe84d02 *src/upwindorderedv4.c
+e7bbf9590ab2bb0ad570fed58a0948a9 *src/upwindorderedv4.c
+aded4b6e87090a7bd79e15db5adef62a *vignettes/0_Introduction_to_quasi-potential_analysis.Rmd
+14155298f48eae7f1ddd4f31c8b357f2 *vignettes/0_Introduction_to_quasi-potential_analysis.html
+ef00ffac171cb17bdf0396905f63909a *vignettes/1_Analyzing_the_deterministic_skeleton.Rmd
+848119eb6286bdf38265a775be1af37d *vignettes/1_Analyzing_the_deterministic_skeleton.html
+a0b9621b3141209c70e0003f3c60751b *vignettes/2_Stochastic_simulation.Rmd
+d338df125849026ffd88f145680dfb98 *vignettes/2_Stochastic_simulation.html
+da744328333828918951918f55a90d7f *vignettes/3_Local_quasi-potential_calculation.Rmd
+44463f3a6c35c3a277cec782036865c9 *vignettes/3_Local_quasi-potential_calculation.html
+ab60b5620ae289d417a51b57d9abb80f *vignettes/4_Global_quasi-potential_calculation.Rmd
+f4f0bb8b4c3c74ab6df3ea732dca9191 *vignettes/4_Global_quasi-potential_calculation.html
+f19fde7b69eaaf230172e4f4a7a64da3 *vignettes/5_Global_quasi-potential_visualization.Rmd
+1027d29ede5cfa5f0eaa3b6a884e7293 *vignettes/5_Global_quasi-potential_visualization.html
+a42c3ba019cd159d79cd8b33d7501332 *vignettes/6_Vector_field_decomposition.Rmd
+7f6977453f718c921f605e2593144a60 *vignettes/6_Vector_field_decomposition.html
+ced148d27dd28de8d221e0fb131dd315 *vignettes/QPot_Example1.md.rsp
+3bddacbc5c685a79f3960785b3d6f52d *vignettes/QPot_Example2.md.rsp
+6cd529b8faf6ac2753dcc22389199e23 *vignettes/QPot_Rarticle.pdf.asis
diff --git a/NAMESPACE b/NAMESPACE
@@ -5,6 +5,7 @@ importFrom("graphics", ".filled.contour", "abline", "arrows", "axis",
 			"contour", "lines", "par", "plot", "polygon")
 importFrom("stats", "density", "rnorm")
 importFrom("utils", "read.table")
+importFrom("graphics", "mtext")
 
 export(Model2String)
 export(QPContour)

diff --git a/NEWS b/NEWS
@@ -6,22 +6,3 @@ R CMD build QPot
 #copy tar.gz to /home to install
 #can't install across partitions
 R CMD INSTALL QPot_
-
-#################################
-examples taken from functions in R/
-##################################
-
-QPContour
-
-# #' @examples
-# #' # First, use a surface (example from QPGlobal)
-# #' global.qp <- QPGlobal(list(local.1,local.2),c(0,4),c(0,4),c(-1,5),c(-1,5))
-# #'
-# #' # Second, input that surface into QPContour
-# #' QPContour(surface=global.qp, density=c(100,100), xlim=c(-0.5,20), y.lim=c(-0.5,20), n.filled.contour=20, n.contour.lines=20, col=c("red", "white", "blue"), contour.lines = TRUE)
-
-
-QPGlobal
-# #' @export
-# #' @examples 
-# #' QPGlobal(list(local.1,local.2),c(0,4),c(0,4),c(-1,5),c(-1,5))
diff --git a/R/Model2String.R b/R/Model2String.R
@@ -1,6 +1,6 @@
 #' Inserts parameter values into equations
 #'
-#' Converts differential equations from string-format (or function-format) with parameters (e.g. "a*x+b) to string-format with parameter values (e.g. 2*x+3).  Specifically, Model2String reads in the equations, searches for the differential equations within the function (if required), and replaces the parameters with numerical values given by the user.  Returns an array of strings containing the differential equations.
+#' Converts differential equations from string-format (or function-format) with parameters (e.g. "a*x+b) to string-format with parameter values (e.g. 2*x+3).  Specifically, Model2String reads in the equations, searches for the differential equations within the function (if required), and replaces the parameters with numerical values given by the user.  Returns an array of strings containing the differential equations.  This code is specifically given so that any problems can be found in R as opposed to receiving terse errors from the C code if this was built into the QPotential() function.
 #' 
 #' @param model contains the differential equations as given to \code{\link{TSTraj}}.  Can either be a string or a function used by the package \code{deSolve} (see third example).
 #' @param parms a named vector of paramters and their respective values for the deterministic equations.  If inputing a function and parms is empty, Model2String will return the equation as a string.
@@ -42,10 +42,11 @@
 #'  x.lhs.term = 'dx', y.lhs.term = 'dy') 
 
 
-Model2String <- function(model, parms = 'NULL', deSolve.form = FALSE, x.lhs.term = 'dx', y.lhs.term = 'dy', supress.print = FALSE, width.cutoff = 500) {
+Model2String <- function(model = NULL, parms = NULL, deSolve.form = FALSE, x.lhs.term = 'dx', y.lhs.term = 'dy', supress.print = FALSE, width.cutoff = 500) {
 	if (!supress.print) {
-		print("Note: This function is supplied as duct tape.  Long equations, equations spanning multiple lines, equations with strange notation, etc, may not work.  Always check the output.")
+		message("Note: This function is supplied to help convert equations to strings.  Long equations, equations spanning multiple lines, equations with strange notation, etc, may not work.  Always check the output.")
 	}
+	if (is.null(model)) {stop("No equation supplied to function Model2String")}
 
 	if (deSolve.form == TRUE) {
 
@@ -59,7 +60,8 @@ Model2String <- function(model, parms = 'NULL', deSolve.form = FALSE, x.lhs.term
 		foundy = 0	#flag for making sure dy is only found once
 
 	#remove the lhs and return the rhs
-		for (i in 1:length(temp)) {
+#		for (i in 1:length(temp)) { 
+		for (i in seq_along(temp)) {
 		#when searching, first look for the lhs defining whether the derivative is for x or y
 		#once found, look inside the string and use either '<-' or '=' to separate
 		# the lhs from the rhs
@@ -97,7 +99,8 @@ Model2String <- function(model, parms = 'NULL', deSolve.form = FALSE, x.lhs.term
 	}
 
 #if parameters are not declared, then we do not have to replace anything
-	if (!(parms[1] == 'NULL')){
+	if (!is.null(parms)) {
+#	if (!(parms[1] == 'NULL')){
 #replace the parameter names in the equations with their values
 		allnames <-names(parms)
 		for (i in 1:length(parms)) {

diff --git a/R/QPContour.R b/R/QPContour.R
@@ -5,13 +5,16 @@
 #' @param dens vector respectively for the number of \code{x} and \code{y} points to be plotted.
 #' @param x.bound a two-element vector with the minimum and maximum x values used for computing the quasi-potential.
 #' @param y.bound a two-element vector with the minimum and maximum y values used for computing the quasi-potential.
-#' @param xlim numeric vectors of length 2, giving the x coordinate range.
-#' @param ylim numeric vectors of length 2, giving the y coordinates range.
+#' @param xlim numeric vectors of length 2, giving the x coordinate range. Default \code{= NULL} automatically sizes plot window.
+#' @param ylim numeric vectors of length 2, giving the y coordinate range. Default \code{= NULL} automatically sizes plot window.
 #' @param n.filled.contour numeric value for the nubmber of breaks in the filled contour.
 #' @param n.contour.lines numeric value for the number of breaks in the contour lines.
 #' @param c.parm contour line adjustment (see details).
 #' @param col.contour colors to interpolate; must be a valid argument to \code{\link{col2rgb}}.
 #' @param contour.lines if \code{TRUE}, then contour lines plotted over filled contour; vice versa if \code{FALSE}.
+#' @param xlab a title for the x axis.  Default is 'X'
+#' @param ylab a title for the y axis.  Default is 'Y'
+#' @param contour.lwd line width of contour lines.
 #' @param ... passes arguments to \code{\link{plot}}.
 #' @details Because, in general, capturing the topological features of a surface can be subtle, we implemented a feature in \code{\link{QPContour}} to keep the filled contour region while changing the contour lines.  Specifically, \code{\link{filled.contour}} takes the range of the surface values (\eqn{\phi}), divides by the number of the specified contours (i.e., \code{n.filled.contour}), and creates a contour at each break, which happenes to be equal across the range.  But because visualizing some topology may (i) require looking between contour breaks and (ii) adding contour lines would overload the plot with lines, we use an equation to modify the distribution of contour lines.  Namely, adjusting the \code{c} argument in the \code{\link{QPContour}} function adjusts the \eqn{c} paramter in the following equation: \deqn{max_\phi \times \left(\frac{x}{n-1}\right)^c.}  This allows the user to keep the same number of contour lines (i.e., specified with \code{n.contour.lines}), but focus them toward the troughs or peaks of the surfaces. At \eqn{c=1}, the contour lines correspond to the filled.contour breaks.  If \eqn{c > 1}, then the contour lines become more concentrated towards the trough.  Similarly, if \eqn{c < 1}, then the contour lines are more focused at the peaks of the surface.  As an example, we change \eqn{c} : \cr \figure{Example3.png}.
 #'
@@ -51,7 +54,7 @@
 #' 	QPContour(ex1.global, dens = c(100,100), x.bound = xbounds, 
 #'		y.bound = ybounds, c.parm = 5)
 
-QPContour <- function(surface, dens, x.bound, y.bound, xlim = 'NULL', ylim = 'NULL', n.filled.contour=25, n.contour.lines=25, c.parm=1, col.contour, contour.lines = TRUE, ...){
+QPContour <- function(surface, dens, x.bound, y.bound, xlim = NULL, ylim = NULL, n.filled.contour = 25, n.contour.lines = 25, c.parm = 1, col.contour, contour.lines = TRUE, xlab = "X", ylab = "Y", contour.lwd = 1, ...){
 	x.range <- max(x.bound)-min(x.bound)
 	y.range <- max(y.bound)-min(y.bound)
 
@@ -62,14 +65,17 @@ QPContour <- function(surface, dens, x.bound, y.bound, xlim = 'NULL', ylim = 'NU
 	row.max <- max(which(surface != 0 , arr.ind = T)[,1])
 	col.min <- min(which(surface != 0 , arr.ind = T)[,2])
 	col.max <- max(which(surface != 0 , arr.ind = T)[,2])
-
-	x.min <- ((row.min-1)/row.range)*x.range + min(x.bound)
-	x.max <- ((row.max-1)/row.range)*x.range + min(x.bound)
-	y.min <- ((col.min-1)/col.range)*y.range + min(y.bound)
-	y.max <- ((col.max-1)/col.range)*y.range + min(y.bound)
 
-	if(missing(xlim)) {xlim <- c(x.min,x.max)}
-	if(missing(ylim)) {ylim <- c(y.min,y.max)}
+	if(missing(xlim)) {
+		x.min <- ((row.min-1)/row.range)*x.range + min(x.bound)
+		x.max <- ((row.max-1)/row.range)*x.range + min(x.bound)
+		xlim <- c(x.min,x.max)
+		}
+	if(missing(ylim)) {
+		y.min <- ((col.min-1)/col.range)*y.range + min(y.bound)
+		y.max <- ((col.max-1)/col.range)*y.range + min(y.bound)		
+		ylim <- c(y.min,y.max)
+		}
 
  	sub.x <- seq(row.min, row.max, length.out=dens[1])
 	sub.y <- seq(col.min, col.max, length.out=dens[2])
@@ -79,13 +85,13 @@ QPContour <- function(surface, dens, x.bound, y.bound, xlim = 'NULL', ylim = 'NU
 
 	eq.sub <- surface[sub.x, sub.y]
 
-	plot(0 , type = "n" , xlim = xlim , ylim = ylim , las = 1, ...)
+	plot(0 , type = "n" , xlim = xlim , ylim = ylim, xlab = xlab, ylab = ylab, ...)
 	min.eq.sub <- min(eq.sub , na.rm = T)
 	max.eq.sub <- max(eq.sub , na.rm = T)
 	contour.breaks <- seq(min.eq.sub , max.eq.sub , length = n.filled.contour)
 	eq.max <- max(surface, na.rm = T)
 	line.contour.breaks <- (eq.max)*(((0:n.contour.lines)/(n.contour.lines-1)))^c.parm
 	myRamp <- if(missing(col.contour)){colorRampPalette(c("#FDE725FF","#E3E418FF","#C7E020FF","#ABDC32FF","#8FD744FF","#75D054FF","#5DC963FF","#47C06FFF","#35B779FF","#28AE80FF","#20A486FF","#1F9A8AFF","#21908CFF","#24868EFF","#287C8EFF","#2C728EFF" ,"#31688EFF","#355D8DFF","#3B528BFF","#404688FF","#443A83FF","#472D7BFF","#481F71FF","#471163FF","#440154FF"))(n.filled.contour)}else{colorRampPalette(col.contour)(n.filled.contour)}
 	.filled.contour(sub.x.val , sub.y.val , eq.sub , levels = contour.breaks , col = myRamp)
-	if(contour.lines==TRUE){contour(sub.x.val , sub.y.val , eq.sub , levels = line.contour.breaks, drawlabels = F ,  add = TRUE , col = "black" , ...)}
-	}
+	if(contour.lines==TRUE){contour(sub.x.val , sub.y.val , eq.sub , levels = line.contour.breaks, drawlabels = F ,  add = TRUE , col = "black", lwd = contour.lwd)}
+	}
diff --git a/R/QPGlobal.R b/R/QPGlobal.R
@@ -67,8 +67,8 @@ QPGlobal <- function(local.surfaces , unstable.eq.x , unstable.eq.y , x.bound ,
 		}
 	}
 
-	max.phi <- max(unstable.phi,na.rm=T)
-	max.phi.arr <- which(unstable.phi == max.phi , arr.ind=T)
+	max.phi <- max(unstable.phi, na.rm = TRUE)
+	max.phi.arr <- which(unstable.phi == max.phi , arr.ind = T)
 
 	if(nrow(max.phi.arr) != n.unstable.pts){ #if not all max(phi) are the same, then they need to be aligned
 		global.max <- unstable.phi[max.phi.arr[1,1],max.phi.arr[1,2]]
@@ -87,8 +87,8 @@ QPGlobal <- function(local.surfaces , unstable.eq.x , unstable.eq.y , x.bound ,
 	}
 	rm(eq.arr.unadj)
 	gc()
-	global.qp <- apply(eq.arr , c(1:n.surfaces) , min)
+	global.qp <- apply(eq.arr , c(1, 2) , min)
 	rm(eq.arr)
 	gc()
 	global.qp
-	}
+	}