From 9a0e2e73b806e57c1fe54abc4f6f4499d99d61cd Mon Sep 17 00:00:00 2001 From: Andre Chalom Date: Fri, 22 May 2015 02:08:41 -0300 Subject: [PATCH 1/3] updatesad now works for models with fixed params --- R/updatesad.R | 14 ++++++++++---- 1 file changed, 10 insertions(+), 4 deletions(-) diff --git a/R/updatesad.R b/R/updatesad.R index eb607cb..1eaad8a 100644 --- a/R/updatesad.R +++ b/R/updatesad.R @@ -1,7 +1,7 @@ ## updates a fitsad/fitrad object by running the optimizer again starting ## on better fit returned by profile updatesad <- function(object, ...) { - dots <- as.list(...) + dots <- list(...) prof <- profile(object) if(class(prof) == "profile.mle2") stop("Cannot update, profile did not find a better fit!") newcall <- as.list(object@call) @@ -9,9 +9,15 @@ updatesad <- function(object, ...) { newcall[["control"]] <- NULL # removes the "control" slot # merges the original call with arguments supplied by ... for (v in names(dots)) newcall[[v]] <- dots[[v]] - names(newcall$start) -> name # bugfix? profile returns coefficient "prob" name as "prob.prob" - newcall$start <- list(prof@fullcoef) # TODO: might be coef instead? how to deal with fixed? - names(newcall$start) <- name + # bugfix? profile returns coefficient "prob" name as "prob.prob", sometimes changes parameters + # from start to fixed... + names(newcall$start) -> name + for (v in name) { + if (v %in% names(prof@coef)) + newcall$start[[v]] <- as.numeric(prof@coef[[v]]) + else # then it must be in fixed + newcall$start[[v]] <- as.numeric(prof@call.orig$fixed[v]) + } newobj <- do.call("mle2", newcall) if(class(object) == "fitsad") return (new("fitsad", newobj, sad=object@sad, distr=object@distr, trunc=object@trunc)) From b04f6148229a84e86d43b520984de87d6bd562f8 Mon Sep 17 00:00:00 2001 From: Andre Chalom Date: Sat, 23 May 2015 18:28:15 -0300 Subject: [PATCH 2/3] improved values for fitmand --- R/fitmand.R | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/R/fitmand.R b/R/fitmand.R index 931e730..0b4b011 100644 --- a/R/fitmand.R +++ b/R/fitmand.R @@ -8,7 +8,9 @@ fitmand <- function(x, trunc, start.value, ...){ if (min(y)<=trunc) stop("truncation point should be lower than the lowest data value") } if(missing(start.value)){ - shat <- 2 + x75 <- rad.tab[(floor(dim(rad.tab)[1]/4)):dim(rad.tab)[1],] + x75 <- x75[x75$abund > 1, ] + shat <- - coef(lm(log(abund)~log(rank), data=x75))[2] vhat <- 30 } else{ From abcd6020480d626bc4aa2ea7373297f73a7a0e1f Mon Sep 17 00:00:00 2001 From: Andre Chalom Date: Thu, 28 May 2015 19:06:46 -0300 Subject: [PATCH 3/3] Quickfix for #52 plotprofmle was issuing an error when the entire profile was lower than "ratio", i.e., no blue line was to be drawn. Fixed adding an "if" to check for that condition --- R/plotprofmle.R | 1 + 1 file changed, 1 insertion(+) diff --git a/R/plotprofmle.R b/R/plotprofmle.R index 6807b2b..40530ce 100644 --- a/R/plotprofmle.R +++ b/R/plotprofmle.R @@ -41,6 +41,7 @@ plotprofmle <- function(profobj, nseg=20, ratio=log(8), which=1:length(profobj@p change <- (interpol$y - ratio)[2:l] * (interpol$y - ratio)[1:(l-1)] endpoints <- which(change < 0) corr <- (interpol$x[2]-interpol$x[1])/2 + if(length(endpoints) > 0) for (j in 1:(length(endpoints)/2)) { lower <-interpol$x[endpoints[(2*j)-1]]+corr upper <- interpol$x[endpoints[2*j]]+corr