addSSsummarize: Add a model to the list of models to compare
In r4ss: R code for Stock Synthesis

Description Usage Arguments Details Value Note Author(s) See Also Examples

Adds specified quantities from any model to the list of models returned from SSsummarize for further comparison.

1	addSSsummarize(origModels,newModels)

`origModels`	A list of models created by `SSsummarize`.
`newModels`	A list of models to add to the originals models list. Each new model is an element of the list, and is a list itself with possible components described in the details below.

The models to add do not need to contain all of the entries. Missing entries are given NA values, so are not plotted. Each model is a list with the following elements:

npars: the number of parameters in the model
maxgrad: the maximum gradient component (if used)
nsexes: the number of sexes
likelihoods: likelihoods from the model. A data.frame with the 2nd column as names, which matches on names from origModels. Names which do not match will be added SS uses the following names

TOTAL
Equil_catch
Survey
Length_comp
Age_comp
Recruitment
recast_Recruitment
Parm_priors
Parm_softbounds
Parm_devs
Crash_Pen
Size_at_age

likelambdas: NOT IMPLEMENTED
pars: NOT IMPLEMENTED YET FOR DIFFICULTY IN MATCHING PARAMETERS
parsSD: NOT IMPLEMENTED YET FOR DIFFICULTY IN MATCHING PARAMETERS
parsphases: NOT IMPLEMENTED YET FOR DIFFICULTY IN MATCHING PARAMETERS
SpawnBio: Spawning biomass matrix

1st column is year
2nd column is spawning biomass in same units as original models (SS reports female spawning biomass)
3rd column is the standard deviation of estimated spawning biomass
4th column is a lower bound of the confidence interval to be plotted (say from an MCMC)
5th column is an upper bound of the confidence interval to be plotted (say from an MCMC)

Bratio: Depletion matrix

1st column is year
2nd column is depletion
3rd column is the standard deviation of depletion (optional)
4th column is a lower bound of the confidence interval to be plotted (say from an MCMC)
5th column is an upper bound of the confidence interval to be plotted (say from an MCMC)

SPRratio: SPR ratio matrix

1st column is year
2nd column is depletion
3rd column is the standard deviation (optional)
4th column is a lower bound of the confidence interval to be plotted (say from an MCMC)
5th column is an upper bound of the confidence interval to be plotted (say from an MCMC)

recruits: Recruitment matrix

1st column is year
2nd column is recruitment as in original models (SS reports age-0 recruits)
3rd column is the standard deviation (optional)
4th column is a lower bound of the confidence interval to be plotted (say from an MCMC)
5th column is an upper bound of the confidence interval to be plotted (say from an MCMC)

recdevs: Recruitment deviate matrix

1st column is year
2nd column is deviate (matched with original models)

growth: NOT IMPLEMENTED
indices: Matrix of fits to indices

1st column is year
2nd column is observed index (data)
3rd column is expected index (prediction)
4th column is catchability coefficient (q)
5th column is standard error of index (total used in fitting)
6th column is a likelihood for this point, or enter any value to make sure it plots, or enter NA not to plot the estimate

InitAgeYrs: NOT IMPLEMENTED

Returns list as is returned from SSsummarize, but contains additions for the new models.

This function was made to compare TINSS results and SS results, and assumed that you would always start with a list of SS models output from SSsummarize. It has not been tested to see how it works when starting with an empty list.

Allan Hicks

SSsummarize SSplotComparisons

  ## Not run: 
  ######################################
  #DO NOT RUN
tinss1 <- list(npars=A$fit$npar,maxgrad=A$fit$maxgrad,nsexes=1,
               #note, there is an estimated parameter called sd_sbt,
               #      but it is a single value
               SpawnBio=data.frame(c(1964,1965,A$yrs),
                                   c(A$sbo,A$sbo,A$sbt)*1e6,0,
                                   qnorm(0.025,c(A$so,A$so,A$sbt)*1e6,0),
                                   qnorm(0.975,c(A$so,A$so,A$sbt)*1e6,0)),
               Bratio=data.frame(A$yrs,A$sbt/A$sbo,0,
                                 qnorm(0.025,A$sbt/A$sbo,0),
                                 qnorm(0.975,A$sbt/A$sbo,0)),
               SPRratio=data.frame(A$yr,A$spr,0,qnorm(0.025,A$spr,0),
                                   qnorm(0.975,A$spr,0)),
               recruits=data.frame(A$yr,A$nt[,1]*1e6,0,qnorm(0.025,A$nt[,1]*1e5,0),
                                   qnorm(0.975,A$nt[,1]*1e6,0)),
               #I'm not sure exactly what wt are,
               #   but it is important to line them up correctly
               recdevs=data.frame(A$recYrs,A$wt),  
               indices = data.frame(A$iyr,1e6*A$yt,1e6*A$qbt,
                                    rep(A$q,length(A$iyr)),rep(0.4,length(A$iyr)),
                                    rep(0,length(A$iyr)),rep(1,length(A$iyr)))
               )
  tinss <- list(tinss1,tinss1)   #can add more models here


  #add TINSS model to SS models already summarized                
  SSnTINSS <- addSSsummarize(models,tinss)
  mcmcInd <- seq(burnin+1,nrow(A$mc.sbt),thin)
  SSnTINSS$mcmc[[2]] <- data.frame(A$mc.sb0[mcmcInd],
                                   A$mc.sbt[mcmcInd,],
                                   A$mc.depl[mcmcInd,],
                                   A$mc.spr[mcmcInd,],
                                   A$mc.rt[mcmcInd,],
                                   log(A$mcmc[mcmcInd,"Ro"]*1e6),
                                   A$mcmc[mcmcInd,"msy"]*1e6)  
  names(SSnTINSS$mcmc[[2]]) <-
    c("SPB_Virgin",paste("SPB",A$yrs,sep="_"),
      paste("Bratio",A$yrs,sep="_"),
      paste("SPRratio",A$yr,sep="_"),
      paste("Recr",A$yr,sep="_"),"SR_R0","TotYield_MSY")
  modelnames <- c("SS", "TINSS","TINSS.MLE")
  
  SSplotComparisons(SSnTINSS, legendlabels=modelnames,
                    subplot=2,endyr=2011,mcmcVec=c(T,T,F))
  title(main="MCMC")
  SSplotComparisons(SSnTINSS, legendlabels=modelnames,
                    subplot=4,endyr=2011,mcmcVec=c(T,T,F))
  title(main="MCMC")
  ###############################################
  
## End(Not run)