monotonicity now "on" by default for adjustment models

extra documentation and fixed tests for the above

bug fix in check.bins -- missing $distance not caught in tests

Author: David Lawrence Miller

R/create.bins.R

@@ -16,7 +16,10 @@ create.bins <- function(data,cutpoints){
   cp <- cutpoints
 
   # remove distances outside bins
-  in.cp.ind <- data>=cp[1] & data<=cp[length(cp)]
+  in.cp.ind <- data$distance>=cp[1] & data$distance<=cp[length(cp)]
+  if(!all(in.cp.ind)){
+    warning("Some distances were outside bins and have been removed.")
+  }
   data <- data[in.cp.ind,]
 
   # pull out the distances (removing the NAs for now)
@@ -33,6 +36,12 @@ create.bins <- function(data,cutpoints){
     distbegin[ind] <- cp[i]
     distend[ind]   <- cp[i+1]
   }
+  # last cutpoint, include those observations AT the truncation point
+  ind <- which(d>=cp[i] & d<=cp[i+1])
+
+  distbegin[ind] <- cp[i]
+  distend[ind]   <- cp[i+1]
+
 
   # handle NA distances, that we need to preserve
   distbegin.na <- rep(NA,length(data$distance))

R/ds.R

@@ -35,10 +35,7 @@
 #'        \code{distend} then these will be used as bins if \code{cutpoints}
 #'        is not specified. If both are specified, \code{cutpoints} has
 #'        precedence.
-#' @param monotonicity should the detection function be constrained for
-#'        monotonicity weakly ("weak"), strictly ("strict") or not at all
-#'        ("none" or \code{FALSE}). See Montonicity, below. (Default
-#'        \code{FALSE}).
+#' @param monotonicity should the detection function be constrained for monotonicity weakly (\code{"weak"}), strictly (\code{"strict"}) or not at all (\code{"none"} or \code{FALSE}). See Montonicity, below. (Default \code{"strict"}).
 #' @param dht.group should density abundance estimates consider all groups to be
 #'        size 1 (abundance of groups) \code{dht.group=TRUE} or should the
 #'        abundance of individuals (group size is taken into account),
@@ -147,7 +144,7 @@
 #' library(Distance)
 #' data(book.tee.data)
 #' tee.data<-book.tee.data$book.tee.dataframe[book.tee.data$book.tee.dataframe$observer==1,]
-#' ds.model<-ds(tee.data,4,monotonicity="strict")
+#' ds.model<-ds(tee.data,4)
 #' summary(ds.model)
 #' plot(ds.model)
 #'
@@ -158,17 +155,20 @@
 #' samples<-book.tee.data$book.tee.samples
 #' obs<-book.tee.data$book.tee.obs
 #'
-#' ds.dht.model<-ds(tee.data,4,region.table=region,monotonicity="strict",
+#' ds.dht.model<-ds(tee.data,4,region.table=region,
 #'              sample.table=samples,obs.table=obs)
 #' summary(ds.dht.model)
 #'
 #' # specify order 2 cosine adjustments
 #' ds.model.cos2<-ds(tee.data,4,adjustment="cos",order=2)
 #' summary(ds.model.cos2)
 #'
-#' # specify order 2 and 3 cosine adjustments - LOTS of non-monotonicity!
-#' ds.model.cos24<-ds(tee.data,4,adjustment="cos",order=c(2,3))
-#' summary(ds.model.cos24)
+#' # specify order 2 and 3 cosine adjustments, turning monotonicity
+#' # constraints off
+#' ds.model.cos24<-ds(tee.data,4,adjustment="cos",order=c(2,3),
+#'                    monotonicity=FALSE)
+#' # check for non-monotonicity -- actually no problems
+#' check.mono(ds.model.cos24$ddf,plot=TRUE,n.pts=100)
 #'
 #' # truncate the largest 10% of the data and fit only a hazard-rate
 #' # detection function
@@ -185,7 +185,7 @@ ds<-function(data, truncation=ifelse(is.null(cutpoints),
              formula=~1, key=c("hn","hr","unif"),
              adjustment=c("cos","herm","poly"),
              order=NULL, scale=c("width","scale"),
-             cutpoints=NULL, monotonicity=FALSE, dht.group=FALSE,
+             cutpoints=NULL, monotonicity="strict", dht.group=FALSE,
              region.table=NULL, sample.table=NULL, obs.table=NULL,
              convert.units=1, method="nlminb", quiet=FALSE, debug.level=0,
              initial.values=NULL){
@@ -488,6 +488,9 @@ ds<-function(data, truncation=ifelse(is.null(cutpoints),
                             " with ", adj.name,"(",
                             paste(order[1:i],collapse=","),
                             ") adjustments", sep="")
+    }else{
+      # if we have only the key function, turn off monotonicity
+      meta.data$mono <- meta.data$mono.strict <- FALSE
     }
 
     model.formula<-paste(model.formula,")",sep="")

man/ds.Rd

@@ -6,7 +6,7 @@ ds(data, truncation = ifelse(is.null(cutpoints), ifelse(is.null(data$distend),
   max(data$distance), max(data$distend)), max(cutpoints)),
   transect = c("line", "point"), formula = ~1, key = c("hn", "hr",
   "unif"), adjustment = c("cos", "herm", "poly"), order = NULL,
-  scale = c("width", "scale"), cutpoints = NULL, monotonicity = FALSE,
+  scale = c("width", "scale"), cutpoints = NULL, monotonicity = "strict",
   dht.group = FALSE, region.table = NULL, sample.table = NULL,
   obs.table = NULL, convert.units = 1, method = "nlminb", quiet = FALSE,
   debug.level = 0, initial.values = NULL)
@@ -75,9 +75,10 @@ ds(data, truncation = ifelse(is.null(cutpoints), ifelse(is.null(data$distend),
   \code{cutpoints} has precedence.}
 
   \item{monotonicity}{should the detection function be
-  constrained for monotonicity weakly ("weak"), strictly
-  ("strict") or not at all ("none" or \code{FALSE}). See
-  Montonicity, below. (Default \code{FALSE}).}
+  constrained for monotonicity weakly (\code{"weak"}),
+  strictly (\code{"strict"}) or not at all (\code{"none"}
+  or \code{FALSE}). See Montonicity, below. (Default
+  \code{"strict"}).}
 
   \item{dht.group}{should density abundance estimates
   consider all groups to be size 1 (abundance of groups)
@@ -261,7 +262,7 @@ using \code{ds()}.
 library(Distance)
 data(book.tee.data)
 tee.data<-book.tee.data$book.tee.dataframe[book.tee.data$book.tee.dataframe$observer==1,]
-ds.model<-ds(tee.data,4,monotonicity="strict")
+ds.model<-ds(tee.data,4)
 summary(ds.model)
 plot(ds.model)
 
@@ -272,17 +273,20 @@ region<-book.tee.data$book.tee.region
 samples<-book.tee.data$book.tee.samples
 obs<-book.tee.data$book.tee.obs
 
-ds.dht.model<-ds(tee.data,4,region.table=region,monotonicity="strict",
+ds.dht.model<-ds(tee.data,4,region.table=region,
              sample.table=samples,obs.table=obs)
 summary(ds.dht.model)
 
 # specify order 2 cosine adjustments
 ds.model.cos2<-ds(tee.data,4,adjustment="cos",order=2)
 summary(ds.model.cos2)
 
-# specify order 2 and 3 cosine adjustments - LOTS of non-monotonicity!
-ds.model.cos24<-ds(tee.data,4,adjustment="cos",order=c(2,3))
-summary(ds.model.cos24)
+# specify order 2 and 3 cosine adjustments, turning monotonicity
+# constraints off
+ds.model.cos24<-ds(tee.data,4,adjustment="cos",order=c(2,3),
+                   monotonicity=FALSE)
+# check for non-monotonicity -- actually no problems
+check.mono(ds.model.cos24$ddf,plot=TRUE,n.pts=100)
 
 # truncate the largest 10\% of the data and fit only a hazard-rate
 # detection function

tests/testthat/test_ds.R

@@ -45,7 +45,7 @@ test_that("Simple models work",{
 
   # specify order 2 cosine adjustments
   ds.model.cos2<-ds(egdata,4,adjustment="cos",order=2, region.table=region,
-                    sample.table=samples,obs.table=obs)
+                    sample.table=samples,obs.table=obs,monotonicity=FALSE)
   # pars and lnl
   #result <- ddf(dsmodel=~mcds(key="hn", formula=~1, adj.series="cos",
   #                            adj.order=2), data=egdata, method="ds",
@@ -58,7 +58,8 @@ test_that("Simple models work",{
 
   # specify order 2 and 4 cosine adjustments
   ds.model.cos24<-ds(egdata,4,adjustment="cos",order=c(2,4),
-                     region.table=region, sample.table=samples, obs.table=obs)
+                     region.table=region, sample.table=samples, obs.table=obs,
+                     monotonicity=FALSE)
   tp <- c(0.92121582, -0.03712634, -0.03495348)
   names(tp) <- c("X.Intercept.","V2","V3")
   expect_equal(ds.model.cos24$ddf$par, tp)