FLash: Forecasting and Projection of Fish Stocks

# fwd.R
# FLash/R/fwd.R
# Copyright 2003-2007 FLR Team. Distributed under the GPL 2 or later
# Maintainer: Finlay Scott, Cefas
# Last Change: 06 Mar 2009 19:17
# $Id: fwdControl.R 366 2009-10-26 09:37:13Z lauriekell $

trgtNms    <-function() return(c("year","min","val","max","quantity","season","area","unit","spp","fleet","metier","rel.year","rel.season","rel.area","rel.unit"))
effNms     <-function() return(c("year","min","val","max","fleet","metier","rel.year","rel.fleet","rel.metier","rel.bound"))
quantityNms<-function() return(c("ssb","biomass","catch","landings","discards","f","z","f.landings","f.discards","effort","costs","revenue","profit","mnsz"))

validFwdControl <- function(object){
	return(TRUE)

  if (dim(object@target)[1]!=dim(object@trgtArray)[1]){
     warning("rows in target & trgtArray don't match")
     return(FALSE)}

  if (any(object@target[,"quantity"] %in% quantityNms)){
     warning("quantity not recognised")
     return(FALSE)}

  if (length(slot(object, 'effort'))>0){
     if (dim(object@effort)[1]!=dim(object@effArray)[1]){
        warning("rows in effort & effArray don't match")
        return(FALSE)}

     if (dim(object@target)[1]!=dim(object@effort)[1]){
        warning("rows in target & effort don't match")
        return(FALSE)}

     if (dim(object@trgtArray)[3]!=dim(object@effArray)[3]){
        warning("iter in trgtArray & effArray don't match")
        return(FALSE)}
     }

	# Everything is fine
	return(TRUE)
  }

setClass("fwdControl",
	representation(
		target   ="data.frame",
		effort   ="data.frame",
    trgtArray="array",
		effArray ="array",
		block    ="numeric"), ## specifies if mulitple rows are done together
	prototype=prototype(
		target   =data.frame(NULL),
		effort   =data.frame(NULL),
    trgtArray=array(),
		effArray =array(),
		block    =numeric()),
	validity=validFwdControl
  )

if (!isGeneric("fwdControl")) {
	setGeneric("fwdControl", function(object, ...){
		value  <-  standardGeneric("fwdControl")
		value
	})}

setMethod("fwdControl", signature(object="data.frame"),
fwdControl.<-function(object,effort=NULL,trgtArray=NULL,effArray=NULL,...){

    ##### Internal Functions ###################################################
    setArray<-function(x,nrws,nits=NULL,type="trgtArray"){
       if (is(x,"list") & any(names(x) %in% c("min","val","max"))){
         if (!all(lapply(x,class) %in% c("array","matrix","numeric")))
            stop(paste(type,": elements of list neither 'array', 'matrix' or 'numeric'"))

         if (is.null(nits))
            if      (is(x[[1]],"numeric"))                     nits<-length(x[[1]])
            else if (is(x[[1]],"array") | is(x[[1]],"matrix")) nits<-dim(x[[1]])[length(dim(x[[1]]))]
            else stop("")

         res<-array(NA,dim=c(nrws,3,nits),dimnames=list(1:nrws,c("min","val","max"),iter=1:nits))
         if ("val" %in% names(x)){
            if (is.vector(x$val)) x$val<-array(x$val,dim=c(1,length(x$val)))
            if (nits == dim(x$val)[2])
               res[,"val",]<-x$val
            }
         if ("min" %in% names(x)){
            if (is.vector(x$min)) x$min<-array(x$min,dim=c(1,length(x$min)))
            if (nits == dim(x$min)[2])
               res[,"min",]<-x$min
            }
         if ("max" %in% names(x)) {
            if (is.vector(x$max)) x$max<-array(x$max,dim=c(1,length(x$max)))
            if (nits == dim(x$max)[2])
               res[,"max",]<-x$max}
            }
       else if (is(x,"array") & (length(dim(x))==3)){
          if (is.null(nits))
             nits<-dim(x)[3]

          res<-array(NA,dim=c(nrws,3,nits),dimnames=list(1:nrws,c("min","val","max"),iter=1:nits))

          res[dimnames(x)[[1]],dimnames(x)[[2]],]<-x
          }
       else stop("Has to be either a 3D array or list with 'min', 'max' or 'val' vectors")

       return(res)
       }
	# Creates data.frame with desired column names (nms) and no. of rows (no. yrs)
	# Used for creating target and effort dataframes
    df<-function(yrs,nms){
      df<-NULL
      for (i in nms)
         df<-cbind(df,rep(NA,length(yrs)))
      dimnames(df)<-list(1:length(yrs),nms)
      return(data.frame(df))
      }

    checkMinMax<-function(object)
        {
        # check that if max or min specified then no target & vice versa
        if (any((!is.na(object[,"min"]) | !is.na(object[,"max"])) & !is.na(object[,"val"]))) {
           cat("Can't specify val and both a min or max values")
           return(FALSE)}
        else if (any((!is.na(object[,"min"]) & !is.na(object[,"max"])) & object[,"max"]<=object[,"min"])){
           cat("max less than than min value")
           return(FALSE)}
        else
           return(TRUE)
        }
    ##### End Internal Functions ###############################################

    if (!is(object,"data.frame"))
       stop("target not data.frame")

    if (!("year" %in% names(object)))
       stop("year not specified in object")
    yrs<-object[,"year"]

    res<-new("fwdControl")

    ##Targets ##################################################################
    ## Create complete target data frame
    res@target<-df(yrs,trgtNms())
    res@target[,dimnames(object)[[2]]]<-object[,dimnames(object)[[2]]]
    if (!checkTarget(res@target))
       stop("target not valid")

    if (!is.null(trgtArray)){
       res@trgtArray<-setArray(trgtArray,length(yrs),type="trgtArray")
       if (length(dim(res@trgtArray[,1,]))==2){
          res@target[,"min"]<-apply(res@trgtArray[,"min",],1,median)
          res@target[,"max"]<-apply(res@trgtArray[,"max",],1,median)
          res@target[,"val"]<-apply(res@trgtArray[,"val",],1,median)}
      else{
          res@target[,"min"]<-median(res@trgtArray[,"min",])
          res@target[,"max"]<-median(res@trgtArray[,"max",])
          res@target[,"val"]<-median(res@trgtArray[,"val",])}}
    else{
       res@trgtArray<-array(as.numeric(NA),dim=c(length(res@target[,1]),3,1),dimnames=list(1:length(res@target[,1]),c("min","val","max"),iter=1))}

    res@target[,"quantity"]<-factor(res@target[,"quantity"],levels=c("ssb","biomass","catch","landings","discards","f","z","f.landings","f.discards","effort","costs","revenue","profit","mnsz"))

    for (i in 1:length(res@target[,1])){
       if (any(is.na(res@trgtArray[i,"min",]))) res@trgtArray[i,"min",]<-res@target[i,"min"]
       if (any(is.na(res@trgtArray[i,"val",]))) res@trgtArray[i,"val",]<-res@target[i,"val"]
       if (any(is.na(res@trgtArray[i,"max",]))) res@trgtArray[i,"max",]<-res@target[i,"max"]}

    if (!checkMinMax(res@target)) {
       cat(" in target\n")
       stop()}

    ##Effort ###################################################################
    if (!is.null(effort)){
      res@effort<-df(yrs,effNms())
      res@effort[ ,dimnames(effort)[[2]]]<-effort[,dimnames(effort)[[2]]]
      if (!is.null(effArray))
         res@effArray<-setArray(effArray,length(yrs),type="effArray")

      if (!is.null(effArray)){
         res@effArray<-setArray(effArray,length(yrs),type="effArray")
         if (length(dim(res@effArray[,1,]))==2){
            res@effort[,"min"]<-apply(res@effArray[,"min",],1,median)
            res@effort[,"max"]<-apply(res@effArray[,"max",],1,median)
            res@effort[,"val"]<-apply(res@effArray[,"val",],1,median)}
         else{
            res@effort[,"min"]<-median(res@effArray[,"min",])
            res@effort[,"max"]<-median(res@effArray[,"max",])
            res@effort[,"val"]<-median(res@effArray[,"val",])}}
      else
         res@effArray<-array(as.numeric(NA),dim=c(length(res@effort[,1]),3,1),dimnames=list(1:length(res@effort[,1]),c("min","val","max"),iter=1))

    for (i in 1:length(res@effort[,1])){
       if (any(is.na(res@effArray[i,"min",]))) res@effArray[i,"min",]<-res@effort[i,"min"]
       if (any(is.na(res@effArray[i,"val",]))) res@effArray[i,"val",]<-res@effort[i,"val"]
       if (any(is.na(res@effArray[i,"max",]))) res@effArray[i,"max",]<-res@effort[i,"max"]}

    if (!checkMinMax(res@effort)){
       cat(" in effort\n")
       stop()}}

   return(res)
   })

showArray<-function(object){
    if(dim(object)[3] > 1){
		  v1 <- apply(object, 1:2, median, na.rm=TRUE)
  		v2 <- apply(object, 1:2, mad,    na.rm=TRUE)
      v3 <- paste(format(v1,digits=5),"(", format(v2, digits=3), ")", sep="")}
    else
      v3 <- paste(format(apply(object, 1:2, median, na.rm=TRUE),digits=5))

    print(array(v3, dim=dim(object)[1:2], dimnames=dimnames(object)[1:2]), quote=FALSE)

		if(dim(object)[3] != 1)
			cat("iter: ", dim(object)[3],"\n\n")}

setMethod('show', signature(object='fwdControl'),
  function(object){

  showDFTarget<-function(object){

      nm      <-names(object@target)
      optional<-c("season","area","unit","rel.year","rel.season","rel.area","rel.unit")
      flag    <-apply(as.matrix(!is.na(object@target[,optional])),2,any)
      
      print(object@target[,c("year","quantity","min","val","max",names(flag[flag]))])

      cat("\n")}

  showDFEffort<-function(object){

      nm      <-names(object@effort)
      optional<-c("fleet","metier","rel.year","rel.fleet","rel.metier","rel.bound")
      flag    <-apply(as.matrix(!is.na(object@target[,optional])),2,any)

      print(object@effort[,c("year","min","val","max",names(flag[flag]))])

      cat("\n")}

  cat("\nTarget\n")
  showDFTarget(object)
  if (any(!is.na(object@trgtArray)))
     showArray(object@trgtArray)

  if (length(slot(object, 'effort'))>0){
     cat("\n\nEffort\n")
     showDFEffort(object)
     if (any(!is.na(object@effArray)))
        showArray(object@effArray)}
  })

chkFwdControl<-function(ctrl,sr,x,y=NULL){
   if (is(x,"FLStock")){

      return(ctrl)
      }
   else if (is(x,"FLBiol")){
      return(ctrl)
      }
   }

checkTarget<-function(target)
    {
    # check that if max or min specified then no target & vice versa
    if (any((!is.na(target[,"min"]) | !is.na(target[,"max"])) & !is.na(target[,"val"]))) {
       warning("Can't specify a val and a min or max values")
       return(FALSE)}

    if (any((!is.na(target[,"min"]) & !is.na(target[,"max"])) & target[,"max"]<=target[,"min"])){
       warning("max less than than min value")
       return(FALSE)}

	# Should also check quantity

    return(TRUE)
    }

matrixTarget <- function(target)
    {
    #reorder columns for C code (???)
    target <- target[,trgtNms()]
    for(i in names(target))
        target[,i] <- as.double(target[,i])

    return(matrix(unlist(target),dim(target)))
    }

checkTarget<-function(target)
    {
    # check that if max or min specified then no target & vice versa
    if (any((!is.na(target[,"min"]) | !is.na(target[,"max"])) & !is.na(target[,"val"]))) {
       warning("Can't specify a val and a min or max values")
       return(FALSE)}

    if (any((!is.na(target[,"min"]) & !is.na(target[,"max"])) & target[,"max"]<=target[,"min"])){
       warning("max less than than min value")
       return(FALSE)}

	  # Should also check quantity

    return(TRUE)
    }

matrixEffort <- function(effort)
    matrix(apply(effort,2,as.double),dim(effort))


chkTrgtArrayIters <- function(object,trgtArray,sr)
    {
    if (is(object,'FLlst')) object <- object[[1]]
    # get iterations from trgtArray, stock, SR parameters and SR residuals
    its<-sort(unique(c(length(dimnames(trgtArray)$iter), dims(object)$iter, length(dimnames(sr$params[[1]])$iter), length(dimnames(sr$residuals[[1]])$iter))))
    if (length(its)>2 | (length(its)>1 & its[1]!=1)) stop("iter not 1 or n")
    if (length(its)==2 & length(dimnames(trgtArray)$iter == 1)){
        dmns<-dimnames(trgtArray)
        dmns$iter<-1:its[2]
        trgtArray<-array(trgtArray,dim=unlist(lapply(dmns,length)),dimnames=dmns)}

    return(trgtArray)
    }

# check target quantity is factor and that it is currently implemented
chkTargetQuantity <- function(target,object)
    {
    ordDmn<-function(dmn,val){
      tmp       <-1:length(dmn)
      names(tmp)<-dmn

      return(tmp[ac(val)])
      }

    if (!is(target[,"quantity"],"factor"))
        target[,"quantity"]<-factor(target[,"quantity"],quantityNms())
    if (!all(as.character(target[,"quantity"]) %in% quantityNms()))
        stop("invalid quantity in control target")
    if (any(as.character(target[,"quantity"]) %in% c("effort","costs","revenue","profit")))
        stop("fwd not yet implemented for 'effort','costs','revenue' or 'profit'")

    if (!is.numeric(target[,"season"])) target[,"season"]<-ordDmn(dimnames(m(object))$season,target[,"season"])
    if (!is.numeric(target[,"unit"]  )) target[,"unit"]  <-ordDmn(dimnames(m(object))$unit,  target[,"unit"])
    if (!is.numeric(target[,"area"]  )) target[,"area"]  <-ordDmn(dimnames(m(object))$area,  target[,"area"])

	  return(target)
    }

# FLQuants -> fwdControl {{{

setAs("FLQuants", "fwdControl",
  function(from) {
		

    # CONVERT
    target <- as.data.frame(from)[,c('year', 'iter', 'data', 'qname')]
    names(target)[3:4] <- c('val', 'quantity')

    # dims and check
    nits <- length(unique(target$iter))

		# ITERS
		if(nits == 1) {

      target <- cbind(target[,-2], fishery=as.numeric(NA), catch=as.numeric(NA),
        biol=1)

      return(fwdControl(target))
		} else {

      dft <- target[target$iter == target$iter[1],][,c('year', 'val', 'quantity')]

			arrt <- array(NA, dim=c(dim(dft)[1], 3, nits),
				dimnames=list(seq(dim(dft)[1]), c('min', 'val', 'max'), iter=unique(target$iter)))
			arrt[,'val',] <- target$val
      
			return(fwdControl(dft, trgtArray=arrt))
		}
	stop('Conversion unsucessful')
} ) # }}}