I did not understand what "proportion of harvestable biomass" exactly mean when specifying argument 'effort'. I have temporarily chosen it as the ratio between catch/biomass observed in the last year, but the forward predictions that I get appear to be not reliable....
I have the same problem to specify 'ratio' from argument 'fleets'.

Any suggestion is greatly appreciated.

Gadget models have three types of output files:

• Model output (e.g. likelihood prediction, biomass by year, step, age etc)
• Likelihood output (values of parameters and likelihood components)
• parameter files (both in wide and long format)

A model output file is created with a specific printfile and is of the form preamble, tabular data and in the case of survey indices a postamble

; Gadget version 2.2.00-BETA running on ComputerName Wed Apr 27 22:03:35 2016
; Likelihood output file for the likelihood component si.F.20-50
; year-step-area-label-number
1997    4          1      len10    6687762.6
....
2015    4          1      len10     30062.25
; Regression information for area 1
; len10 intercept 6.805886 slope -5.3238671e-10 sse 1e+10

Here you want the tabular data to be read as a data.frame, names of columns read from the preamble and the information from the postamble and joined with the table. The postamble can include one or more regression results, i.e. the length of the postamble is 2*unique(label) lines.

The likelihood output files has three components, first a listing where the variable occurs, second the likelihood components names, type and weights, and third table of the likelihood function values (total and each component) along with the parameter values. The file is created using:
gadget -s -i params.in -o lik.out

; Gadget version 2.2.00-BETA running on ComputerName Sat Jun 11 17:00:34 2016
; Listing of the switches used in the current Gadget run
Linf    codimmgrowergrowthcalcHmodelvariable   
....      
;
;
; Listing of the likelihood components used in the current Gadget run
;
; Component     Type    Weight
matp.igfs       6       48.095
alkeys.aut      3       33599.6
alkeys.comm     3       16937.3
alkeys.igfs     3       35086.4
ldist.aut       3       35162.7
ldist.comm      3       23701.9
ldist.igfs      3       30961.6
si.gp1  1       7.1061
si.gp1a 1       1.08626
si.gp2  1       29.1428
si.gp2a 1       2.40657
si.gp3  1       85.0689
si.gp3a 1       100
understocking   2       100
bounds  8       10
;
; Listing of the output from the likelihood components for the current Gadget run
;
0              160         90          6          1       1000         35         25         15          7          7          5          5         70         70          1          1         50         50         40         30         30         30         30         30         30         30         30         30         30         30         30         30         30         30         30         30         30         30         30         30         30         30         30         30         30         30         25         10         10        -10         19         -4          9         -4          9               328.9    10.02    17.34    27.55    5.653    38.11     22.4    24.14     15.1    7.866    1e+10    2.352    2.783 1.22e+18        0          1.21991794478447e+20

The third format is gadget parameter files which has format very similar to the gadget output files.
These outputs files are not currently handled correctly by read.gadget.file, as they are mainly tabular. Currently these files are read in using read.printfiles, read.gadget.lik.out and read.gadget.parameters but ideally read.gadget.file should digest all of these.

Thanks Bjarki, the new gadget.forward deals now nicely with the argument custom.print. Only drawback is that the format of gadget.forward object is changed with each species into separate elements of the list.
I tend to prefer the previous format, so I've a small code to arrange the new output into the old one. I assume that the custom printfiles do not contain the words "lw", "recruitment", "num.trials", "stochastic", "sim.begin"

arrange.gadget.forward <- function(x){
    out <- vector("list",7)
    names(out) <- c("lw","catch","custom.print","recruitment","num.trials","stochastic","sim.begin")
# $catch
    tmp <- x[grepl("catch",names(x)) & grepl("lw",names(x))]
    for(i in 1:length(tmp)){
        stk <- names(tmp)[i]
        stk <- gsub("catch.","",stk)
        stk <- gsub(".lw","",stk)
        tmp[[i]]$stock <- stk
    }
    out[["catch"]] <- bind_rows(tmp)
# $lw
    tmp <- x[grepl("lw",names(x)) & !(grepl("catch",names(x)))]
    for(i in 1:length(tmp)){
        stk <- names(tmp)[i]
        stk <- gsub(".lw","",stk)
        tmp[[i]]$stock <- stk
    }
    out[["lw"]] <- bind_rows(tmp)
# $custom.print
    out[["custom.print"]] <- x[!(grepl("lw",names(x))         |
                                 grepl("sim.begin",names(x))  |
                                 grepl("stochastic",names(x)) |
                                 grepl("num.trials",names(x)) |
                                 grepl("recruitment",names(x)))]
# $recruitment, $num.trials, $stochastic, $sim.begin
    out[["recruitment"]] <- x[["recruitment"]]
    out[["num.trials"]] <- x[["num.trials"]]
    out[["stochastic"]] <- x[["stochastic"]]
    out[["sim.begin"]] <- x[["sim.begin"]]
    out
}

I use the function below as follows:
stockDistr.plot(data=read.printfiles("WGTS/out.fit/")$sprat.std, ages=1:8, plusGroup=T, lines=F)

Plot the estimated age structure of the stock

stockDistr.plot <- function(data, ages, plusGroup=T, step=1, area=1, lines=F){

data - data.frame with the columns 'number', 'age', 'year'

ages - vector with ages to plot

step - specify the step to plot (only one at a time), default = 1

area - vector with the areas to sum across, default = 1

plusGroup - if the last year is a plusGroup

lines - add lines to follow individual cohorts

require(RColorBrewer)

data <- data[data$step==step & data$area %in% area,]
adist <- tapply(data$number, list(data$year, data$age), sum)

if(plusGroup==TRUE){
    tmp <- adist[,which(colnames(adist) %in% ages[length(ages)]):dim(adist)[2]]
    tmp <- apply(tmp,1,sum)
    adist <- adist[,which(colnames(adist) %in% ages)]
    adist <- cbind(adist[,-dim(adist)[2]], tmp)
    colnames(adist) <- ages
    colnames(adist)[dim(adist)[2]] <- paste(colnames(adist)[dim(adist)[2]], "+",sep="")}
if(plusGroup==FALSE){
    adist <- adist[,which(colnames(adist) %in% ages)]} else {}

col.pal <- brewer.pal(11,"Spectral")
col.pal <- rep(col.pal,4)

age <- colnames(adist)
yr <- rownames(adist)
xx <- seq(0.75, 0.75+(1+0.2)*(length(yr)-1), 1+0.2) # x unit needed in barplot to account for inter-bar space

par(mfrow=c(length(age),1), mai=c(0.1,0,0,0),omi=c(0.5,0.5,0.1,0.1))
for(a in 1:length(age)){
    barplot(adist[,a],xaxt="n",col=col.pal[10:31-a], space=0.2)
    if(lines==TRUE){
        for(i in 1:length(xx)){
            lines(c(xx[i]-0.5,xx[i]+0.5),c(max(adist[,a]),0), col="grey")
        }
    } else {}
    legend("topright",paste("Age",age[a]), bty="n",cex=1.5)
}
axis(1,xx,yr)
mtext("Year", 1, line=2, outer=T)
mtext("Number of fish", 2, line=2.5, outer=T)

}

I wonder what does "mortality" represent in fit.forward$catch, where fit.forward is an object fitted using gadget.forward().

I stumbled upon an error when importing gadget fleet files using read.gadget.file. A typical fleet file looks like this

[fleetcomponent]
totalfleet              comm
livesonareas            1
multiplicative          1.0
suitability
codimm                  function exponential    #acomm (* 0.01 #bcomm)  0 1
codmat                  function exponential    #acomm (* 0.01 #bcomm)  0 1
amount                  Data/cod.fleet.data

Whereas the output from gadget.file is:

[fleetcomponent]
totalfleet  comm
livesonareas    1
multiplicative  1
codimm  function exponential    #acomm (* 0.01 #bcomm)  0 1
codmat  function exponential    #acomm (* 0.01 #bcomm)  0 1
amount  Data/cod.fleet.data

Note the missing suitability line.

Hello !

When I run:
fit <- gadget.fit()
for whiting I get the error:
Error in data.frame(name = x, pos = pos, ord = regexpr(x, text[pos])[1], :
arguments imply differing number of rows: 1, 0

So, I can't figure out what is the issue.

Thanks very much !

it would be very useful if 'gadget.forward(..., compact=F)' could return also projected catches by 'fleet' and 'length'

I noticed that the values of fishing mortality in the historical (not forecasted) part extracted by gadget.forward() are different from those estimated by gadget.fit(). This does not happen when I provide a print.file disaggregated by age, but in this case there are problems with the forecasted future time series that shows a very little random variation.

Gadget formulae can mention 2 file types:

• Time Variables: A data.frame of year, step, value. Read in the data on read, when eval'ing give a year/step to choose the value
• Stock Variables: See https://github.com/Hafro/gadget/blob/master/src/stockvariable.cc - Contains a list of stock names. Read it in, worry about how it's actually eval'ed later.

Length distribution from survey are only predicted in one area but not the others (despite of present observation data and catch data in that areas for the survey fleet): column predicted in those are all NAs, number.y and total.pred are all 0s.

Is there a method in rgadget for calculating confidence intervals around gadget estimates (e.g. CI around the estimated catch fitted to a time-series, and for estimates of recruitment and F)?

When running multi-area model with gadget.iterative I got an error:
*** Error ingadget': free(): invalid next size (fast): 0x0000000000d5d5f0 ***`
Run continued though after the error, but at the end all files in WGTS folder starting with lik. or params. were empty.

Is there a reason why a recsigma object is created only when the minage of a stock is >0 (see code below)? I always seem to have problems with these lines of code anyways, but wondering why it is not well supported to have a minage of 0 when doesrenew == 1?

if (length(x$sigma) == 1) {
            x$sigma <- rep(x$sigma, length(x$minage:x$maxage))
        }
        else if (length(x$sigma) >= x$maxage & x$minage > 0) {
            recsigma <- x$sigma[1]
            x$sigma <- x$sigma[x$minage:x$maxage]
        }
        else if (length(x$sigma) > length(x$minage:x$maxage)) {
            x$sigma <- x$sigma[1:length(x$minage:x$maxage)]
        }
        else {
            stop("Error -- initial sigma useless")
        }

I'm trying 'gadget.forward' on a 3 species, 4 fleets model with some default values as follows:
fwd <- gadget.forward(params.file="WGTS/params.final", years=10, fleets=data.frame(fleet=c('comm.act','comm.pas','comm.her','comm.spr'), ratio=c(1,1,1,1)), num.trials=10, effort=0.5, rec.window=2000:2012, stochastic=F)

regardless the values I get the following error:
Error in read.table(sprintf("%s/out/%s.lw", pre, x), comment.char = ";") : no lines available in input In addition: Warning messages: 1: In initialize(value, ...) : NAs introduced by coercion 2: running command 'gadget -s -i ./PRE/params.forward -main ./PRE/main.pre 2>/dev/null' had status 1

Thankful for any suggestion

simple scatterplot of obs vs. pred survey indices.

This how I use it:
surveyIndex.plot(data=fit, yearLab=TRUE, cex=0.8)

surveyIndex.plot <- function(data, h=TRUE, yearLab=FALSE, ...){
## data: gadget.fit object returned from 'fit.gadget()'
## h: horizontal orientation of the multi-panel plot
## yearLab: replace dots with year labels

surv <- sort(unique(data$sidat$name))
par(mfrow=panel.size(length(surv),h=h), oma=c(1,1,0.2,0.2)) # see additional function below

# not sure why $acoustic and $number.x are distinct (this merges the two columns)
data$sidat$number.x <- ifelse(is.na(data$sidat$number.x), data$sidat$acoustic, data$sidat$number.x)

for(i in 1:length(surv)){
    tmp <- data$sidat[data$sidat$name==surv[i],]
    plot(tmp$number.x, tmp$predict, type="n", xlab="Observed", ylab="Predicted", main=surv[i])
    if(yearLab==FALSE){
        points(tmp$number.x, tmp$predict, ...)}
    if(yearLab==TRUE){
        text(tmp$number.x, tmp$predict, substring(tmp$year,3,4), ...)}
}}

automatically define multi-panel size (nx,ny)

panel.size <- function(n, h=TRUE){

n = number of panels to plot

h = T (nx<=ny)

if(h==TRUE){
    out <- c(ceiling(n/sqrt(2_n)), floor(sqrt(2_n)))
} else {
    out <- c(floor(sqrt(2_n)), ceiling(n/sqrt(2_n)))}
if(out[1] \* out[2] < n){
    out <- c(ceiling(n/sqrt(2_n)), ceiling(sqrt(2_n)))
} else {NULL}
out  

}

I noticed that whether I select fixedslopeloglinear or loglinear fit to survey indices will influence my model a lot (age distribution, length distribution, etc.). How then to choose a fitting function? If the slope is estimated (in loglinear function), the values of it are very different. Some are close to one, while other are 0.5 or even 6.5. What those values really tell us? How should we select a value of the slope when choosing fixedslopeloglinear function?

I tried to read a stock file with rgadget, but unfortunately got an error
Error: could not find function "read.gadget.stockfiles"

TimeVariables are not handled by Rgadget, ideally eval.gadget.formula should do something intelligent when asked to deal with timevariables.

After fitting a model with gadget.fit() I got a file out$likelihoodsummary. Weights of likelihood component there are always 1. What does it mean? Should they always be 1 at that file?

Hi, the function gadget.ypr is only for a normal parametric distribution in renewal. It doesn't work for the other distribution types.

I've noticed that when 'linearfleet' is used, fit$fleet.info$amount reports the effort and not the catch as expected. I suspect that also other output ie fit$catchdist.fleets may have issues.

When running gadget.forward (with last Rgadget version installed) on a model assuming seasonal recruitments I get the following error:

fit.forward <- gadget.forward(years = 11, num.trials = 101,
fleets = data.frame(fleet = c("commBMT","TUNcommBMT","survey"),
ratio = c(1,1,1)), effort = 0.62, stochastic = TRUE,
rec.window = 2002:2014)

Errore: wrong result size (44), expected 176 or 1

This does not happen for a model assuming annual recruitments.

Gadget.fit currently assumes that all lengths are in centimeters, making illustrations of fit for a model of shrimp less than ideal:)

I have a 3-area model, and recruitment is in each of the areas (in recruitment file each year-area combination has it's own parameter value). However something went wrong after using gadget.fit():

head(fit$res.by.year[fit$res.by.year$stock=="codim",])

year stock area catch num.catch F total.number total.biomass
4 1974 codim area1 29463112 30475039 0.5114301 366729520 91762669
5 1974 codim area1 29463112 30475039 0.5114301 366729520 91762669
6 1974 codim area1 29463112 30475039 0.5114301 366729520 91762669
7 1974 codim area2 21673890 22514621 0.3737798 360020534 94085918
8 1974 codim area2 21673890 22514621 0.3737798 360020534 94085918
9 1974 codim area2 21673890 22514621 0.3737798 360020534 94085918
harv.biomass ssb .id recruitment
4 0 0 codim 138727350
5 0 0 codim 278095920
6 0 0 codim 269194550
7 0 0 codim 138727350
8 0 0 codim 278095920
9 0 0 codim 269194550

Since, number of rows for immature cod got tripled and recruitment value is not unique for year-area combination, but rather for the year and row position.

Will need to define a standard plot of the model fit to stomach data.

The only thing that ought be a "Depends" for any ~modern R package is the R versions you support, everything else ought be in "Imports".

The reason for this is that it adds the packages into the session's search path, potentially breaking other code (although I'm guessing it's making your life easier by loading all the packages you're using at the R prompt).

Hi, I think that there is a typo in line 149 of gadget.fit function, the line says:
f.age.range <- c(max(out[[sprintf("%s.prey", x)]]$age), max(out[[sprintf("%s.prey", x)]]$age))
and perhaps it should be:
f.age.range <- c(min(out[[sprintf("%s.prey", x)]]$age), max(out[[sprintf("%s.prey", x)]]$age)),

Hi,

I updated Rgadget to version 0.5. Now I get the error with gadget.iterative:
Error: Each variable must be named. Problem variables: 5

traceback()
35: stop(..., call. = FALSE, domain = NA)
34: stopc(problem, ".\n", "Problem variables: ", format_n(vars))
33: invalid_df("Each variable must be named", x, which(bad_name))
32: check_tibble(x)
31: list_to_tibble(x, validate, raw_rownames(x))
30: as_data_frame.data.frame(data)
29: as_data_frame(data)
28: tbl_df(.data)
27: filter_(tbl_df(.data), .dots = dots)
26: as.data.frame(filter_(tbl_df(.data), .dots = dots))
25: filter_.data.frame(.data, .dots = lazyeval::lazy_dots(...))
24: filter_(.data, .dots = lazyeval::lazy_dots(...))
23: dplyr::filter(., year %in% year_range)
22: function_list[k]
21: withVisible(function_list[k])
20: freduce(value, _function_list)
19: _fseq(_lhs)
18: eval(expr, envir, enclos)
17: eval(quote(_fseq(_lhs)), env, env)
16: withVisible(eval(quote(_fseq(_lhs)), env, env))
15: dat %>% dplyr::filter(year %in% year_range)
14: .fun(piece, ...)
13: (function (i)
{
piece <- pieces[[i]]
if (.inform) {
res <- try(.fun(piece, ...))
if (inherits(res, "try-error")) {
piece <- paste(utils::capture.output(print(piece)),
collapse = "\n")
stop("with piece ", i, ": \n", piece, call. = FALSE)
}
}
else {
res <- .fun(piece, ...)
}
progress$step()
res
})(1L)
12: .Call(loop_apply_, as.integer(n), f, env)
11: loop_apply(n, do.ply)
10: llply(.data = pieces, .fun = .fun, ..., .progress = .progress,
.inform = .inform, .parallel = .parallel, .paropts = .paropts)
9: dlply(x, "name", read.func)
8: .fun(piece, ...)
7: (function (i)
{
piece <- pieces[[i]]
if (.inform) {
res <- try(.fun(piece, ...))
if (inherits(res, "try-error")) {
piece <- paste(utils::capture.output(print(piece)),
collapse = "\n")
stop("with piece ", i, ": \n", piece, call. = FALSE)
}
}
else {
res <- .fun(piece, ...)
}
progress$step()
res
})(2L)
6: .Call(loop_apply_, as.integer(n), f, env)
5: loop_apply(n, do.ply)
4: llply(.data = pieces, .fun = .fun, ..., .progress = .progress,
.inform = .inform, .parallel = .parallel, .paropts = .paropts)
3: dlply(subset(likelihood$weights, !(type %in% c("penalty", "understocking",
"migrationpenalty"))), "type", function(x) dlply(x, "name",
read.func))
2: read.gadget.data(likelihood, year_range = time[[1]]$firstyear:time[[1]]$lastyear)
1: gadget.iterative(rew.sI = TRUE, grouping = list(sind.1.1 = grep("bits1",
lik$surveyindices$name, value = TRUE), sind.4 = grep("bits4",
lik$surveyindices$name, value = TRUE), surveyQ1 = c("ldist.bits1",
"alk.bits1"), surveyQ4 = c("ldist.bits4", "alk.bits4"), catch.a = c("ldist.comm.active"),
catch.p = c("ldist.comm.passive"), catch = c("adist.comm")),
params.file = "params.base", wgts = "WGTS")

I tried few different models and all of them got an error.

There are two groups of features which I haven't been able to find in 'gadget.forward' but which would be relevant.

1. specify the seasonal pattern in the effort of the fleet. Ideally, the most flexible setting would be to give the user the freedom to specify the fleet.ratio for all combinations of fleet, year, timestep

2. arguments like rec.window and mat.par should accept stock specific input

When reading a fleet file of the form:

[fleetcomponent]
totalfleet      igfs
livesonareas    1
multiplicative  1
suitability
lingimm function        exponentiall50  #ling.igfs.alpha        #ling.igfs.l50
lingmat function        exponentiall50  #ling.igfs.alpha        #ling.igfs.l50
amount  three_fleets/Data/fleet.igfs.data

writing this out results in:

[fleetcomponent]
totalfleet      igfs
livesonareas    1
multiplicative  1
suitability
suitability     
lingimm function        exponentiall50  #ling.igfs.alpha        #ling.igfs.l50
lingmat function        exponentiall50  #ling.igfs.alpha        #ling.igfs.l50
amount  three_fleets/Data/fleet.igfs.data

I.e. the line containing suitabilityis duplicated.

In most of the models the fit$res.by.year$F is zero for the first 2-3 years, despite there are catches.
Also the "mortality" reported in WGTS/out.fit/stockName.prey (which refers to fishing mortality as far as I understand) is different from zero.
Any suggestion?

Hi,

I run gadget.iterative for my multi-species model. Everything works smoothly. Then I run gadget.fit(compile.fleet.info=FALSE) (since I'm using linearfleet for one species; for other two it's totalfleet). I got an error here: Error in names(len.agg)[1:3] <- c("length", "lower", "upper") :
attempt to set an attribute on NULL. However if I close the file with script, reopen it and only run gadget.fit(compile.fleet.info=FALSE) ( with weights being ready by previous step) I don't receive this error anymore. Why can this happen? Thank you a lot in advance!

P.S. I have separate fleet files for three species.

With best wishes,
Nataliia

I've a model implemented with a linearfleet and effort data, and a 'catchinkilos' likelihood component.

The gadget.iterative always returns weight = 1 for the 'catchinkilos' likelihood component.
In practice I do not get this component weighted, but checking the model I cannot really figure out why (no warning or error is reported and all files seem in place).

Could it be related to gadget.iterative?

read.printfiles reads the first 7 lines in from output file and tries to extrapolate the header information from the last line with ';' . This can be problematic for small survey index output files as they include regression information in the last few lines.

There are two separate issues here:

1. fit$fleet.info does not account for the argument multiplicative in the fleetfile

2. when I use a linearfleet the catch amount in fit$fleet.info does not correspond to the catch amount in fit$res.by.year

Recursively read everything in a model (i.e. so reading a mainfile will read everything).

Don't worry about file pointers, if the whole thing doesn't fit in memory then you can't run the Gadget model anyway.

The function gadget.fit currently only allow for length and age based indices, other survey indices currently crash the function.

This function bind different gadget runs ("most of the output") for easy comparison plot.
Gfit is a list of 2 or more gadget.fit-class objects:

gadget.bind <- function(Gfit, ...){
    out <- Gfit[[1]]
    out[] <- NA
    for(i in 1:length(out)){
        if(sum(class(Gfit[[1]][[names(out)[i]]]) == "list") == 1 &
           length(Gfit[[1]][[names(out)[i]]]) > 1){ NULL } else {
               tmp <- sapply(Gfit,function(x){x[names(out)[i]]})
               names(tmp) <- names(Gfit) # replace with readable run names
               out[[i]] <- bind_rows(tmp, .id="run")
        }}
    return(out)
}

This is an example of how to use it:

fitL <- gadget.bind(list(r16=fit16, r17=fit17, r18=fit18))
ggplot(fitL$res.by.year,aes(year,total.biomass)) +
    geom_line(aes(color=run)) +
    theme_bw() +
    ylab("TSB") +
    xlab('Year')

gadget.forward gives an error when the renewal is of the type normal condition distribution, to fix it temporarily it's necessary to remove two lines around lines 204 and add a new one as follows:

[email protected] <- rbind.fill(subset([email protected],
year < sim.begin), data.frame(year = rec.years,
step = tail([email protected]$step, 1), area = tail([email protected]$area,
1), age = tail([email protected]$age, 1), number = sprintf("(* 0.0001 #rec%s )",
rec.years), mean = tail([email protected]$mean,
1), stddev = tail([email protected]$stddev,
1), * #relcond = tail([email protected]$relcond, 1),
#alpha = tail(stocks$[email protected]$alpha, 1),
#beta = tail([email protected]$beta, 1),*
stringsAsFactors = FALSE))

and also in the stockfile inside folder PRE, to replace by hand normalparamfile by normalcondfile.

I am attempting to set up a gadget model with spawning instead of renewal for the recruitment; however, gadget.fit does not want to work with a spawning model. I think the issue is on line(s) 180-181 where bio.by.year gets merged with stock.recruitment.

# code above this creates the bio.by.year data.frame (lines 162-179)
   }) %>% dplyr::bind_rows() %>% dplyr::left_join(stock.recruitment %>% 
   dplyr::mutate(area = paste0("area", area), year = as.numeric(year)))

I'd like to rectify this issue by avoiding the left join with stock.recruitment (since there is no stock.recruitment data.frame in this instance). Is there a particular way that you'd like to see it done? or just any way that works.

Thanks,
Paul

Minor nit: gadget.forward has an empty @return line, which Roxygen moans about so you can't compile the .Rd files.

I believe there's a small bug in the suitability function.

The code is:
else if (tolower(type) == "gamma") {
S <- array(rep((l/((params[1] - 1) * params[2] * params[3]))^(params[1] -
1) * exp(params[1] - 1 - 1/(params[2] * params[3])),
each = length(L)), c(length(L), length(l)))
}

and I think it should be changed with:
else if (tolower(type) == "gamma") {
S <- array(rep((l/((params[1] - 1) * params[2] * params[3]))^(params[1] -
1) * exp(params[1] - 1 - l/(params[2] * params[3])),
each = length(L)), c(length(L), length(l)))
}

I'm using Rgadget_0.1

Does anyone have a suggestion on how to do forecasts based on predefined F levels (ie Fmsy). For instance in the ICES area this is relevant if we would like to test the current single species management into our multispecies models.

As I see it a first step would be to extend the function gadget.forward so that the arguments effort and/or fleet(ratio) could be specified for each year and possibly also timestep of the projection. This links to issue #26

I'm currently quite stuck with this which is an urgent matter. Grateful for any help.

After I run gadget.iterative():
`> tmp <- gadget.iterative(rew.sI=TRUE, grouping=list(SI=c('cpue','cpue2', 'biomass','biomass2')))

Joining by: comp
[1] "Running final"
[1] "Comp final completed"`
,

I run gadget.fit():
`> fit <- gadget.fit()

Joining, by = c(".id", "length")
Joining, by = c("year", "area", "fleet")
Error in fix.by(by.y, y) : 'by' must specify uniquely valid columns`

So, why do I get that error ?

When I run fit <- gadget.fit() and then inspect fit$fleet.info I notice that harv.rate can assume values
greater than 1. Is it possible?

fit$fleet.info
Source: local data frame [70 x 7]
Groups: year, step, area, fleet [70]

year  step  area                 fleet       harv.bio         amount            harv.rate

1 2002 1 1 commBMT 1625700.0 2379448.4 1.463645e+00
2 2002 1 1 commDP 10874676.2 107551.6 9.890097e-03
3 2002 1 1 survey 1436912.2 93.3 6.493090e-05
4 2002 1 1 TUNcommBMT 1020173.6 1538000.0 1.507587e+00
5 2002 1 1 TUNcommDP 10874676.2 0.0 0.000000e+00
6 2003 1 1 commBMT 1337630.3 1925937.1 1.439813e+00
7 2003 1 1 commDP 8104348.2 87052.9 1.074151e-02
8 2003 1 1 survey 1356842.0 106.6 7.856478e-05
9 2003 1 1 TUNcommBMT 850947.9 1202000.0 1.412542e+00
10 2003 1 1 TUNcommDP 8104348.2 0.0 0.000000e+00

... with 60 more rows

I have updated Rgadget and I'm having problem in running 'gadget.forward' on previous models and codes:

Joining, by = "stock"
Error in if (!is_readable(full_path)) return(NULL) : 
  argument is of length zero
In addition: Warning message:
In left_join_impl(x, y, by$x, by$y, suffix$x, suffix$y) :
  joining factor and character vector, coercing into character vector```

I have problem to run 'gadget.iterative' more than once in the same session and I start to get the following error:

Error in file(con, "r") : all connections are in use

Then a new R session is required.

I use this LH for weight-at-age data.

> fit <- gadget.fit()
Joining by: c(".id", "length")
Joining by: c("year", "area", "fleet")
Warning message:
In FUN(X[[i]], ...) :
  Error in read.printfile -- Header could not be read from file WGTS/out.fit/wgtAtAge.comm

I think the reason is that the header made by Gadget in WGTS/out.fit/wgtAtAge.comm is:
; year-step-area-age-number-mean[-stddev]
and '-' within the last column gets confused in read.printfile

Somehow simulated gadget fleet is not catching any fish at a certain period (despite of existent observational data; and despite if suitability parameters are or aren't optimized)