spmboot {MQMF}

spmboot conducts a bootstrap analysis on a spm model

Description

spmboot conducts a bootstrap analysis on a spm model. It does this by saving the original fishery data, estimating the cpue residuals, and multiplying the optimum predicted CPUE by a bootstrap sample of the log-normal residuals (Haddon, 2011, p311; and the On Uncertainty chapter in the URMQMF book, p195). This bootstrap sample of CPUE replaces the original fish[,"cpue"] and the model is re-fitted. This is repeated iter times and the outputs reported ready for the derivation of percentile confidence intervals. The optimum solution is used as the first bootstrap replicate (it is standard practice to include the original fit in the bootstrap analysis). If 1000 replicates are run this procedure can take a couple of minutes on a reasonably fast computer. A comparison of the mean with the median should provide some notion of any bias in the mean estimate.

Usage

spmboot(optpar, fishery, iter = 1000, schaefer = TRUE)

Arguments

optpar

The optimum model parameters from fitting a surplus production model

fishery

fishery data containing the original observed cpue values

iter

the number of bootstrap replicates to be run, default=1000

schaefer

default=TRUE, should a Schaefer or a Fox model be run

Value

a list of two matrices. One containing the bootstrap parameters and the other containing some of the dynamics, including the ModelB, the bootstrap CPUE sample, the Depletion, and annual harvest rate.

Examples

 data(dataspm); fish <- as.matrix(dataspm)
 pars <- log(c(r=0.24,K=5150,Binit=2800,0.15))
 ans <- fitSPM(pars,dataspm,schaefer=TRUE,maxiter=1000)
 boots <- spmboot(ans$estimate,fishery=fish,iter=5,schaefer=TRUE) 
 dynam <- boots$dynam
 bootpar <- boots$bootpar
 rows <- colnames(bootpar)
 columns <- c(c(0.025,0.05,0.5,0.95,0.975),"Mean")
 bootCI <- matrix(NA,nrow=length(rows),ncol=length(columns),
                dimnames=list(rows,columns))
 for (i in 1:length(rows)) {
   tmp <- sort(bootpar[,i])
   qtil <- quantile(tmp,probs=c(0.025,0.05,0.5,0.95,0.975),na.rm=TRUE)
   bootCI[i,] <- c(qtil,mean(tmp,na.rm=TRUE))
 }
 round(bootCI,3) # we used only 5 bootstraps for a speedy example, normally 
 pairs(bootpar[,c("r","K","Binit","MSY")]) # use 1000, 2000, or more