diff --git a/DESCRIPTION b/DESCRIPTION index 4070f88..ab8db35 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -46,9 +46,11 @@ Collate: 'SCAPars-class.R' 'SCAMCMC-class.R' 'a4aFit-class.R' + 'a4aFits-class.R' 'a4aFitSA-class.R' 'a4aFitSAs-class.R' 'a4aFitMCMC-class.R' + 'a4aFitMCMCs-class.R' 'a4aFitresiduals-class.R' 'coef-methods.R' 'vcov-methods.R' diff --git a/NAMESPACE b/NAMESPACE index be260d8..cb556c7 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -50,10 +50,12 @@ export(sep.sa) #export(plotIters) exportClasses(a4aFit) +exportClasses(a4aFits) exportClasses(a4aFitResiduals) exportClasses(a4aFitSA) exportClasses(a4aFitSAs) exportClasses(a4aFitMCMC) +exportClasses(a4aFitMCMCs) exportClasses(a4aFitCatchDiagn) exportClasses(a4aGr) exportClasses(a4aM) @@ -80,9 +82,11 @@ exportMethods(coefficients) exportMethods("coefficients<-") exportMethods(as.mcmc) exportMethods(a4aFit) +exportMethods(a4aFits) exportMethods(a4aFitSA) exportMethods(a4aFitSAs) exportMethods(a4aFitMCMC) +exportMethods(a4aFitMCMCs) exportMethods(a4aGr) exportMethods(a4aM) exportMethods(breakpts) diff --git a/NEWS b/NEWS index 7a69cb4..6b0ab5f 100644 --- a/NEWS +++ b/NEWS @@ -4,6 +4,8 @@ - plot for residuals by age - multisca method to run several sca - plot for multiple fits +- new classes for plural fit objects + BUG FIXES: - sca use of covariates fixed diff --git a/R/a4aFitMCMCs-class.R b/R/a4aFitMCMCs-class.R new file mode 100644 index 0000000..a6cbac8 --- /dev/null +++ b/R/a4aFitMCMCs-class.R @@ -0,0 +1,111 @@ +#==================================================================== +# plural class for a4aFitMCMC (used for model averaging) +#==================================================================== + +#' @rdname a4aFitMCMC-class +#' @aliases a4aFitMCMCs-class + +setClass("a4aFitMCMCs", + contains="a4aFitSAs" +) + +setValidity("a4aFitMCMCs", + function(object) { + if(!all(sapply(object, is, 'a4aFitMCMC'))) { + "Components must be a4aFitMCMC" + } else { + TRUE + } +}) + + +#' @rdname a4aFitMCMC-class +#' @aliases a4aFitMCMCs a4aFitMCMCs-methods +setGeneric("a4aFitMCMCs", function(object, ...) standardGeneric("a4aFitMCMCs")) + +#' @rdname a4aFitMCMC-class +setMethod("a4aFitMCMCs", signature(object="list"), + function(object, ...) { + args <- list(...) + + # names in args, ... + if("names" %in% names(args)) { + names <- args[['names']] + } else { + # ... or in object, + if(!is.null(names(object))) { + names <- names(object) + # ... or in elements, ... + } else { + names <- unlist(lapply(object, name)) + # ... or 1:n + idx <- names == "NA" | names == "" + if(any(idx)) + names[idx] <- as.character(length(names))[idx] + } + } + + # desc & lock + args <- c(list(Class="a4aFitMCMCs", .Data=object, names=names), + args[!names(args)%in%'names']) + + return( + do.call('new', args) + ) + +}) + +#' @rdname a4aFitMCMC-class +setMethod("a4aFitMCMCs", signature(object="a4aFitMCMC"), function(object, ...) { + lst <- c(object, list(...)) + a4aFitMCMCs(lst) +}) + +#' @rdname a4aFitMCMC-class +setMethod("a4aFitMCMCs", signature(object="missing"), + function(...) { + # empty + if(missing(...)){ + new("a4aFitMCMCs") + # or not + } else { + args <- list(...) + object <- args[!names(args)%in%c('names', 'desc', 'lock')] + args <- args[!names(args)%in%names(object)] + do.call('a4aFitMCMCs', c(list(object=object), args)) + } + } +) + +#' @title Plot of metrics of multiple fits +#' @name plot metrics of multiple fits +#' @docType methods +#' @rdname plot-mmcfits +#' @aliases plot,a4aFitMCMCs, missing-method +#' @description Method to plot fitting statistics of multiple fits, useful to compare fits. +#' @param x an \code{a4aFitMCMCs} object with multiple fits +#' @param y ignored +#' @param ... additional argument list that might never be used +#' @return a \code{plot} with fitting statistics +#' @examples +#' data(ple4) +#' data(ple4.index) +#' qmods <- list(list(~s(age, k=6))) +#' fmods = list() +#' for(i in 1:6) { +#' fmods[[paste0(i)]] <- as.formula(paste0("~te(age, year, k = c(6,", i+14,"), bs = 'tp') + s(age, k = 6)")) +#' } +#' myFits <- FLa4a:::multisca(FLStocks(ple4), list(FLIndices(ple4.index)), fmodel = fmods, qmodel=qmods, fit="MCMC", mcmc=SCAMCMC()) +#' plot(myFits) + +setMethod("plot", c("a4aFitMCMCs", "missing"), function(x, y=missing, ...){ + args <- list() + par(mar=c(5, 4, 4, 4) + 0.1) + accrate = lapply(x,function(x) fitSumm(x)['accrate',]) + df <- data.frame(unlist(accrate)) + df$fit <- as.numeric(gsub("fit", "",names(accrate))) + names(df) <- c("accrate","fit") + plot(df$fit, df$accrate, type = "b", col = "blue", ylab = "accrate", xlab = "fit", main="Analysis of fit metrics", ylim=c(0.1,0.5)) + abline(h=0.3, col = "blue",lty = 2) +}) + diff --git a/R/a4aFitSAs-class.R b/R/a4aFitSAs-class.R index 80f8fef..798f9b7 100644 --- a/R/a4aFitSAs-class.R +++ b/R/a4aFitSAs-class.R @@ -6,7 +6,7 @@ #' @aliases a4aFitSAs-class setClass("a4aFitSAs", - contains="FLComps" + contains="a4aFits" ) setValidity("a4aFitSAs", @@ -18,7 +18,6 @@ setValidity("a4aFitSAs", } }) - #' @rdname a4aFitSA-class #' @aliases a4aFitSAs a4aFitSAs-methods setGeneric("a4aFitSAs", function(object, ...) standardGeneric("a4aFitSAs")) @@ -78,43 +77,3 @@ setMethod("a4aFitSAs", signature(object="missing"), ) -#' @title Plot of metrics of multiple fits -#' @name plot metrics of multiple fits -#' @docType methods -#' @rdname plot-mfits -#' @aliases plot,a4aFitSAs, missing-method -#' @description Method to plot fitting statistics of multiple fits, useful to compare fits. -#' @param x an \code{a4aFitSAs} object with multiple fits -#' @param y ignored -#' @param ... additional argument list that might never be used -#' @return a \code{plot} with fitting statistics -#' @examples -#' data(ple4) -#' data(ple4.index) -#' qmods <- list(list(~s(age, k=6))) -#' fmods = list() -#' for(i in 1:6) { -#' fmods[[paste0(i)]] <- as.formula(paste0("~te(age, year, k = c(6,", i+14,"), bs = 'tp') + s(age, k = 6)")) -#' } -#' myFits <- FLa4a:::multisca(FLStocks(ple4), list(FLIndices(ple4.index)), fmodel = fmods, qmodel=qmods) -#' plot(myFits) - -setMethod("plot", c("a4aFitSAs", "missing"), function(x, y=missing, ...){ - args <- list() - par(mar=c(5, 4, 4, 4) + 0.1) - gcv = lapply(x,function(x) fitSumm(x)['gcv',]) - bic = lapply(x, function(x) BIC(x)) - df <- data.frame(unlist(gcv), unlist(bic)) - df$fit <- as.numeric(gsub("fit", "",names(gcv))) - names(df) <- c("GCV","BIC","fit") - df <- df[complete.cases(df),] - plot(df$fit, df$GCV, type = "b", col = "blue", ylab = "GCV", xlab = "fit", main="Analysis of fit metrics") - par(new = TRUE) - plot(df$fit, df$BIC, type = "b", col = "red", axes = FALSE, xlab = "", ylab = "") - axis(4) - mtext("BIC", side=4, line=3) - abline(v=df[min(df$GCV)==df$GCV,]$fit, col = "blue",lty = 2) - abline(v=df[min(df$BIC)==df$BIC,]$fit, col = "red",lty = 2) - legend("topleft", legend = c("GCV", "BIC"), col = c("blue", "red"), lty = 1, bg="white") -}) - diff --git a/R/a4aFits-class.R b/R/a4aFits-class.R new file mode 100644 index 0000000..a78e41b --- /dev/null +++ b/R/a4aFits-class.R @@ -0,0 +1,120 @@ +#==================================================================== +# plural class for a4aFit (used for model averaging) +#==================================================================== + +#' @rdname a4aFit-class +#' @aliases a4aFits-class + +setClass("a4aFits", + contains="FLComps" +) + +setValidity("a4aFits", + function(object) { + if(!all(sapply(object, is, 'a4aFit'))) { + "Components must be a4aFit" + } else { + TRUE + } +}) + + +#' @rdname a4aFit-class +#' @aliases a4aFits a4aFits-methods +setGeneric("a4aFits", function(object, ...) standardGeneric("a4aFits")) + +#' @rdname a4aFit-class +setMethod("a4aFits", signature(object="list"), + function(object, ...) { + args <- list(...) + + # names in args, ... + if("names" %in% names(args)) { + names <- args[['names']] + } else { + # ... or in object, + if(!is.null(names(object))) { + names <- names(object) + # ... or in elements, ... + } else { + names <- unlist(lapply(object, name)) + # ... or 1:n + idx <- names == "NA" | names == "" + if(any(idx)) + names[idx] <- as.character(length(names))[idx] + } + } + + # desc & lock + args <- c(list(Class="a4aFits", .Data=object, names=names), + args[!names(args)%in%'names']) + + return( + do.call('new', args) + ) + +}) + +#' @rdname a4aFit-class +setMethod("a4aFits", signature(object="a4aFit"), function(object, ...) { + lst <- c(object, list(...)) + a4aFits(lst) +}) + +#' @rdname a4aFit-class +setMethod("a4aFits", signature(object="missing"), + function(...) { + # empty + if(missing(...)){ + new("a4aFits") + # or not + } else { + args <- list(...) + object <- args[!names(args)%in%c('names', 'desc', 'lock')] + args <- args[!names(args)%in%names(object)] + do.call('a4aFits', c(list(object=object), args)) + } + } +) + + +#' @title Plot of metrics of multiple fits +#' @name plot metrics of multiple fits +#' @docType methods +#' @rdname plot-mfits +#' @aliases plot,a4aFits, missing-method +#' @description Method to plot fitting statistics of multiple fits, useful to compare fits. +#' @param x an \code{a4aFits} object with multiple fits +#' @param y ignored +#' @param ... additional argument list that might never be used +#' @return a \code{plot} with fitting statistics +#' @examples +#' data(ple4) +#' data(ple4.index) +#' qmods <- list(list(~s(age, k=6))) +#' fmods = list() +#' for(i in 1:6) { +#' fmods[[paste0(i)]] <- as.formula(paste0("~te(age, year, k = c(6,", i+14,"), bs = 'tp') + s(age, k = 6)")) +#' } +#' myFits <- FLa4a:::multisca(FLStocks(ple4), list(FLIndices(ple4.index)), fmodel = fmods, qmodel=qmods, fit="MP") +#' plot(myFits) + +setMethod("plot", c("a4aFits", "missing"), function(x, y=missing, ...){ + args <- list() + par(mar=c(5, 4, 4, 4) + 0.1) + gcv = lapply(x,function(x) fitSumm(x)['gcv',]) + bic = lapply(x, function(x) BIC(x)) + df <- data.frame(unlist(gcv), unlist(bic)) + df$fit <- as.numeric(gsub("fit", "",names(gcv))) + names(df) <- c("GCV","BIC","fit") + df <- df[complete.cases(df),] + plot(df$fit, df$GCV, type = "b", col = "blue", ylab = "GCV", xlab = "fit", main="Analysis of fit metrics") + par(new = TRUE) + plot(df$fit, df$BIC, type = "b", col = "red", axes = FALSE, xlab = "", ylab = "") + axis(4) + mtext("BIC", side=4, line=3) + abline(v=df[min(df$GCV)==df$GCV,]$fit, col = "blue",lty = 2) + abline(v=df[min(df$BIC)==df$BIC,]$fit, col = "red",lty = 2) + legend("topleft", legend = c("GCV", "BIC"), col = c("blue", "red"), lty = 1, bg="white") +}) + diff --git a/R/fittingFunctions.R b/R/fittingFunctions.R index 74c1fde..bc7ddf3 100644 --- a/R/fittingFunctions.R +++ b/R/fittingFunctions.R @@ -1499,7 +1499,7 @@ fitADMB <- function(fit, wkdir, df.data, stock, indices, full.df, #' @param srmodel a list of \code{srmodel} objects, each with a formula object depicting the model for log recruitment #' @param n1model a list of \code{n1model} objects, each with a formula object depicting the model for the first year of catch data #' @param vmodel a list of \code{vmodel} objects, each with a list of formula objects depicting the models for log survey and log fishing mortality variance -#' @param stock an \code{FLStocks} object, each component with a \ciode{FLStock} object containing catch and stock information +#' @param stock an \code{FLStocks} object, each component with a \code{FLStock} object containing catch and stock information #' @param combination.all bolean parameter (default is FALSE) to define if a full factorial across all stocks, indices, and submodel is run or just a sequence of runs. #' @param ... all other arguments to be passed to \code{sca} #' @return an \code{a4aFits} or \code{a4aFitSAs} or \code{a4aFitMCMCs} depending on the argument \code{fit} @@ -1577,8 +1577,12 @@ multisca <- function(stocks, indicess, fmodel = missing, qmodel = missing, srmod do.call("sca", args) }) - names(fits) <- paste0("fit", c(1:length(fits))) - return(a4aFitSAs(fits)) + names(fits) <- paste0("fit", c(1:length(fits))) + # the sequqnce of the following commands matter + if(is(fits[[1]], "a4aFitMCMC")) return(a4aFitMCMCs(fits)) + if(is(fits[[1]], "a4aFitSA")) return(a4aFitSAs(fits)) + if(is(fits[[1]], "a4aFit")) return(a4aFits(fits)) + } diff --git a/tests/addition.Rout b/tests/addition.Rout deleted file mode 100644 index 81e3d70..0000000 --- a/tests/addition.Rout +++ /dev/null @@ -1,328 +0,0 @@ - -R version 3.5.2 (2018-12-20) -- "Eggshell Igloo" -Copyright (C) 2018 The R Foundation for Statistical Computing -Platform: x86_64-pc-linux-gnu (64-bit) - -R is free software and comes with ABSOLUTELY NO WARRANTY. -You are welcome to redistribute it under certain conditions. -Type 'license()' or 'licence()' for distribution details. - - Natural language support but running in an English locale - -R is a collaborative project with many contributors. -Type 'contributors()' for more information and -'citation()' on how to cite R or R packages in publications. - -Type 'demo()' for some demos, 'help()' for on-line help, or -'help.start()' for an HTML browser interface to help. -Type 'q()' to quit R. - -> #==================================================================== -> # tests for predict -> #==================================================================== -> library(FLa4a) -Loading required package: FLCore -Loading required package: lattice -Loading required package: iterators -FLCore (Version 2.6.12.9002, packaged: 2019-03-21 15:32:42 UTC) -Loading required package: triangle -Loading required package: copula -Loading required package: coda -This is FLa4a 1.6.8. For overview type 'help(package="FLa4a")' - -> data(ple4) -> data(ple4.indices) -> ple4.indices <- ple4.indices[c(3,4,5)] -> data(ple4.index) -> nits <- 2 -> -> #==================================================================== -> # abundance indices -> #==================================================================== -> #-------------------------------------------------------------------- -> # both assessment types -> #-------------------------------------------------------------------- -> # 1 x 1 -> fit0 <- sca(ple4, FLIndices(ple4.index), qmodel=list(~s(age, k=4))) -> stk0 <- ple4 + fit0 -> fit <- sca(ple4, FLIndices(ple4.index), qmodel=list(~s(age, k=4)), fit="assessment") -> stk <- ple4 + fit -> all.equal(stk, stk0) -[1] TRUE -> -> #-------------------------------------------------------------------- -> # N -> #-------------------------------------------------------------------- -> # N x 1 -> set.seed(1234) -> stk0 <- propagate(ple4, nits) + fit -> dims(stk0)$iter==nits -[1] TRUE -> -> # are the right slots being copyed -> identical(c(stock.n(stk0)), c(stock.n(fit)[,,,,,rep(1,nits)])) -[1] TRUE -> identical(c(catch.n(stk0)), c(catch.n(fit)[,,,,,rep(1,nits)])) -[1] TRUE -> identical(c(harvest(stk0)), c(harvest(fit)[,,,,,rep(1,nits)])) -[1] TRUE -> -> set.seed(1234) -> stk <- propagate(ple4, nits) * fit -> dims(stk)$iter==nits -[1] TRUE -> -> # ref sim for comparison -> fits <- simulate(fit, nits, 1234) -> -> # are the right slots being copyed -> identical(stock.n(stk), stock.n(fits)) -[1] TRUE -> identical(catch.n(stk), catch.n(fits)) -[1] TRUE -> identical(harvest(stk), harvest(fits)) -[1] TRUE -> -> # must be different because there's no simulation involved in stk0 -> !identical(stk, stk0) -[1] TRUE -> -> # 1 x N -> stk0 <- ple4 + fits -> dims(stk0)$iter==nits -[1] TRUE -> -> # are the right slots being copyed -> identical(stock.n(stk0), stock.n(fits)) -[1] TRUE -> identical(catch.n(stk0), catch.n(fits)) -[1] TRUE -> identical(harvest(stk0), harvest(fits)) -[1] TRUE -> -> set.seed(1234) -> stk <- ple4 * propagate(fit, nits) -> dims(stk)$iter==nits -[1] TRUE -> -> # are the right slots being copyed -> identical(stock.n(stk), stock.n(fits)) -[1] TRUE -> identical(catch.n(stk), catch.n(fits)) -[1] TRUE -> identical(harvest(stk), harvest(fits)) -[1] TRUE -> -> # must be equal because simulations are controled in both cases -> identical(stk, stk0) -[1] TRUE -> -> # N x N -> stk0 <- propagate(ple4, 2) + fits -> dims(stk0)$iter==nits -[1] TRUE -> -> # are the right slots being copyed -> identical(stock.n(stk0), stock.n(fits)) -[1] TRUE -> identical(catch.n(stk0), catch.n(fits)) -[1] TRUE -> identical(harvest(stk0), harvest(fits)) -[1] TRUE -> -> set.seed(1234) -> stk <- propagate(ple4, 2) * propagate(fit, 2) -> dims(stk)$iter==nits -[1] TRUE -> -> # are the right slots being copyed -> identical(stock.n(stk), stock.n(fits)) -[1] TRUE -> identical(catch.n(stk), catch.n(fits)) -[1] TRUE -> identical(harvest(stk), harvest(fits)) -[1] TRUE -> -> # must be equal because simulations are controled in both cases -> identical(stk, stk0) -[1] TRUE -> -> #==================================================================== -> # biomass indices -> #==================================================================== -> -> # creating idx 1 -> bioidx <- FLIndexBiomass(FLQuant(NA, dimnames=list(age="all", year=range(ple4)["minyear"]:range(ple4)["maxyear"])), name="bioidx") -> index(bioidx) <- stock(ple4)*0.001 -> index(bioidx) <- index(bioidx)*exp(rnorm(index(bioidx), sd=0.1)) -> range(bioidx)[c("startf","endf")] <- c(0,0) -> -> #-------------------------------------------------------------------- -> # both assessment types -> #-------------------------------------------------------------------- -> fit0 <- sca(ple4, FLIndices(bioidx), qmodel=list(~1)) -> stk0 <- ple4 + fit0 -> fit <- sca(ple4, FLIndices(bioidx), qmodel=list(~1), fit="assessment") -> stk <- ple4 + fit -> all.equal(stk, stk0) -[1] TRUE -> -> #-------------------------------------------------------------------- -> # N -> #-------------------------------------------------------------------- -> # N x 1 -> set.seed(1234) -> stk0 <- propagate(ple4, nits) + fit -> dims(stk0)$iter==nits -[1] TRUE -> -> # are the right slots being copyed -> identical(c(stock.n(stk0)), c(stock.n(fit)[,,,,,rep(1,nits)])) -[1] TRUE -> identical(c(catch.n(stk0)), c(catch.n(fit)[,,,,,rep(1,nits)])) -[1] TRUE -> identical(c(harvest(stk0)), c(harvest(fit)[,,,,,rep(1,nits)])) -[1] TRUE -> -> # ref sim for comparison -> fits <- simulate(fit, nits, 1234) -> -> # 1 x N -> stk0 <- ple4 + fits -> dims(stk0)$iter==nits -[1] TRUE -> -> # are the right slots being copyed -> identical(stock.n(stk0), stock.n(fits)) -[1] TRUE -> identical(catch.n(stk0), catch.n(fits)) -[1] TRUE -> identical(harvest(stk0), harvest(fits)) -[1] TRUE -> -> # N x N -> stk0 <- propagate(ple4, 2) + fits -> dims(stk0)$iter==nits -[1] TRUE -> -> # are the right slots being copyed -> identical(stock.n(stk0), stock.n(fits)) -[1] TRUE -> identical(catch.n(stk0), catch.n(fits)) -[1] TRUE -> identical(harvest(stk0), harvest(fits)) -[1] TRUE -> -> ## must fail -> ##stk0 <- propagate(ple4, 2) + fit0 -> -> ## N x 1 -> #set.seed(1234) -> #stk0 <- propagate(ple4, 2) + fit -> #set.seed(1234) -> #stk <- propagate(ple4, 2) * fit -> #all.equal(stk, stk0) -> -> ## 1 x N -> #stk0 <- ple4 + simulate(fit, 2, 1234) -> #stk <- ple4 * simulate(fit, 2, 1234) -> #all.equal(stk, stk0) -> -> #stk0 <- ple4 + simulate(fit, 2, 1234, ple4) -> #stk <- ple4 * simulate(fit, 2, 1234, ple4) -> #all.equal(stk, stk0) -> -> ## N x N -> #stk0 <- propagate(ple4, 2) + simulate(fit, 2, 1234) -> #stk <- propagate(ple4, 2) * simulate(fit, 2, 1234) -> #all.equal(stk, stk0) -> -> #==================================================================== -> # both indices -> #==================================================================== -> -> #-------------------------------------------------------------------- -> # both assessment types -> #-------------------------------------------------------------------- -> fit0 <- sca(ple4, FLIndices(bioidx, ple4.index), qmodel=list(~1, ~s(age, k=4))) -> stk0 <- ple4 + fit0 -> fit <- sca(ple4, FLIndices(bioidx, ple4.index), qmodel=list(~1, ~s(age, k=4)), fit="assessment") -> stk <- ple4 + fit -> all.equal(stk, stk0) -[1] TRUE -> -> #-------------------------------------------------------------------- -> # N -> #-------------------------------------------------------------------- -> # N x 1 -> set.seed(1234) -> stk0 <- propagate(ple4, nits) + fit -> dims(stk0)$iter==nits -[1] TRUE -> -> # are the right slots being copyed -> identical(c(stock.n(stk0)), c(stock.n(fit)[,,,,,rep(1,nits)])) -[1] TRUE -> identical(c(catch.n(stk0)), c(catch.n(fit)[,,,,,rep(1,nits)])) -[1] TRUE -> identical(c(harvest(stk0)), c(harvest(fit)[,,,,,rep(1,nits)])) -[1] TRUE -> -> # ref sim for comparison -> fits <- simulate(fit, nits, 1234) -> -> # 1 x N -> stk0 <- ple4 + fits -> dims(stk0)$iter==nits -[1] TRUE -> -> # are the right slots being copyed -> identical(stock.n(stk0), stock.n(fits)) -[1] TRUE -> identical(catch.n(stk0), catch.n(fits)) -[1] TRUE -> identical(harvest(stk0), harvest(fits)) -[1] TRUE -> -> # N x N -> stk0 <- propagate(ple4, 2) + fits -> dims(stk0)$iter==nits -[1] TRUE -> -> # are the right slots being copyed -> identical(stock.n(stk0), stock.n(fits)) -[1] TRUE -> identical(catch.n(stk0), catch.n(fits)) -[1] TRUE -> identical(harvest(stk0), harvest(fits)) -[1] TRUE -> -> ## must fail -> ##stk0 <- propagate(ple4, 2) + fit0 -> -> ## N x 1 -> #set.seed(1234) -> #stk0 <- propagate(ple4, 2) + fit -> #set.seed(1234) -> #stk <- propagate(ple4, 2) * fit -> #all.equal(stk, stk0) -> -> ## 1 x N -> #stk0 <- ple4 + simulate(fit, 2, 1234) -> #stk <- ple4 * simulate(fit, 2, 1234) -> #all.equal(stk, stk0) -> -> #stk0 <- ple4 + simulate(fit, 2, 1234, ple4) -> #stk <- ple4 * simulate(fit, 2, 1234, ple4) -> #all.equal(stk, stk0) -> -> ## N x N -> #stk0 <- propagate(ple4, 2) + simulate(fit, 2, 1234) -> #stk <- propagate(ple4, 2) * simulate(fit, 2, 1234) -> #all.equal(stk, stk0) -> -> -> -> proc.time() - user system elapsed - 39.721 2.658 42.435 diff --git a/tests/multisca.R b/tests/multisca.R new file mode 100644 index 0000000..f2d4509 --- /dev/null +++ b/tests/multisca.R @@ -0,0 +1,46 @@ +#==================================================================== +# tests for multiple sca +#==================================================================== +library(FLa4a) +data(ple4) +data(ple4.index) + +#==================================================================== +# run fits sca MP +#==================================================================== + +qmods <- list(list(~s(age, k=6))) +fmods = list() +for(i in 1:6) { + fmods[[paste0(i)]] <- as.formula(paste0("~te(age, year, k = c(6,", i+14,"), bs = 'tp') + s(age, k = 6)")) +} + +#-------------------------------------------------------------------- +# MP +#-------------------------------------------------------------------- + +fits <- FLa4a:::multisca(FLStocks(ple4), list(FLIndices(ple4.index)), fmodel = fmods, qmodel=qmods, fit="MP") + +is(fits, "a4aFits") +is(fits[[1]], "a4aFit") + + +#-------------------------------------------------------------------- +# SA +#-------------------------------------------------------------------- + +fits <- FLa4a:::multisca(FLStocks(ple4), list(FLIndices(ple4.index)), fmodel = fmods, qmodel=qmods, fit="assessment") + +is(fits, "a4aFitSAs") +is(fits[[1]], "a4aFitSA") + +#-------------------------------------------------------------------- +# MCMC +#-------------------------------------------------------------------- + +fits <- FLa4a:::multisca(FLStocks(ple4), list(FLIndices(ple4.index)), fmodel = fmods, qmodel=qmods, fit="MCMC", mcmc=SCAMCMC()) + +is(fits, "a4aFitMCMCs") +is(fits[[1]], "a4aFitMCMC") + + diff --git a/tests/sca.R b/tests/sca.R index 4b31586..8cdbd6b 100644 --- a/tests/sca.R +++ b/tests/sca.R @@ -519,4 +519,3 @@ fit2 <- sca(ple4, ple4.indices, center=1) identical(fitSumm(fit0)["nlogl",], fitSumm(fit1)["nlogl",]) identical(fitSumm(fit0)["nlogl",], fitSumm(fit2)["nlogl",]) identical(fitSumm(fit1)["nlogl",], fitSumm(fit2)["nlogl",]) - diff --git a/tests/testthat.R b/tests/testthat.R deleted file mode 100644 index 9cfca68..0000000 --- a/tests/testthat.R +++ /dev/null @@ -1,12 +0,0 @@ -# This file is part of the standard setup for testthat. -# It is recommended that you do not modify it. -# -# Where should you do additional test configuration? -# Learn more about the roles of various files in: -# * https://r-pkgs.org/testing-design.html#sec-tests-files-overview -# * https://testthat.r-lib.org/articles/special-files.html - -library(testthat) -library(FLa4a) - -test_check("FLa4a")