importCol <- function(res.file, Dev=FALSE, CPUE=FALSE, Survey=FALSE, CAc=FALSE,
CAs=FALSE, CLc=FALSE, CLs=FALSE, LA=FALSE, quiet=TRUE)
{
## Implementation notes
## read* functions: Vector, Matrix
## get* functions: N, B, Sel, Dev, CPUE, Survey, CAc, CAs, CLc, CLs, LA
## Global objects (used inside get* functions) are *.file, *.vector, and quiet
## 1 Define functions
readVector <- function(keyword, same.line=TRUE, file=res.file,
vector=res.vector)
{
## Extract white-space delimited numeric vector followed by keyword
line <- match(keyword, substring(vector,1,nchar(keyword)))
v <- if(same.line)
as.numeric(scan(file, what="", skip=line-1, nlines=1, quiet=TRUE)[-1])
else
as.numeric(scan(file, what="", skip=line, nlines=1, quiet=TRUE))
if(!quiet) cat("vector...")
v
}
readMatrix <- function(keyword, nrow, header=FALSE,
stripe=c("no","left","right","upper","lower"),
file=res.file, vector=res.vector)
{
## Extract white-space delimited numeric matrix followed by keyword
## Striped matrices, like N@A in Coleraine: Female|Male|Female|Male|...
stripe <- match.arg(stripe)
line <- match(keyword,substring(vector,1,nchar(keyword))) +
as.numeric(header)
m <- scan(file, skip=line, nlines=nrow, quiet=TRUE)
m <- matrix(m, byrow=TRUE, nrow=nrow)
m <- switch(stripe,
left=m[,seq(1,ncol(m)/2)],
right=m[,seq(ncol(m)/2+1,ncol(m))],
upper=m[seq(1,nrow(m)-1,by=2),],
lower=m[seq(2,nrow(m),by=2),], m)
if(!quiet) cat("matrix...")
m
}
getN <- function(sexes, years, ages)
{
## Look for "\"Numbers_at_age_by_Year,sex_and_age\" in colera31,
## or "Numbers_at_age_by_Year,sex_and_age" in colera32
if(!quiet) cat("N ")
nsexes <- length(sexes)
nyears <- length(years)
nages <- length(ages)
if(nsexes == 1)
{
Nu <- readMatrix("Numbers_at_age_by_Year,sex_and_age", nrow=nyears*nsexes)
N <- data.frame(Sex=rep(sexes,nyears*nages), Year=rep(years,each=nages),
Age=rep(ages,nyears), N=as.vector(t(Nu)),
stringsAsFactors=FALSE)
}
if(nsexes == 2)
{
Nf <- readMatrix("Numbers_at_age_by_Year,sex_and_age",
nrow=nyears*nsexes, stripe="upper")
Nm <- readMatrix("Numbers_at_age_by_Year,sex_and_age",
nrow=nyears*nsexes, stripe="lower")
N <- data.frame(Sex=rep(sexes,each=nyears*nages),
Year=rep(rep(years,each=nages),2), Age=rep(ages,2*nyears),
N=as.vector(t(rbind(Nf,Nm))), stringsAsFactors=FALSE)
}
if(!quiet) cat("OK\n")
N
}
getB <- function(years, gears)
{
ngears <- length(gears)
if(!quiet) cat("B ")
vb <- readMatrix("Vulnerable_Biomass_by_Method_and_Year", nrow=ngears)
sb <- readVector("Spawning_Biomass_by_Year", same.line=FALSE)
y <- c(readVector("Total_Catch_by_Method_and_Year", same.line=FALSE), NA)
B <- data.frame(years=years, vb=t(vb), sb=sb, y=y)
names(B) <- if(ngears==1) c("Year", "VB", "SB", "Y")
else c("Year", paste("VB",gears,sep="."), "SB", "Y")
if(!quiet) cat("OK\n")
B
}
getSel <- function(gears, surveys, years, sexes, ages)
{
if(!quiet) cat("Sel ")
ngears <- length(gears)
nsurveys <- length(surveys)
nyears <- length(years)
nsexes <- length(sexes)
nages <- length(ages)
com <- readMatrix(
"Commercial_age-specific_selectivity_by_method,Year,sex_and_age",
nrow=ngears*nyears*nsexes)
## First year: selectivity
com <- com[seq(1, to=ngears*nyears*nsexes, by=nyears),]
srv <- readMatrix(
"Survey_age-specific_selectivity_by_survey,Year,sex_and_age",
nrow=nsurveys*nsexes)
fecundity <- readVector("Fecundity_by_year_and_age", same.line=FALSE)
weight <- readVector("Weight_by_year,sex_and_age", same.line=FALSE)
## Sexes have same maturity
mat <- rep(ifelse(weight>0,fecundity/weight,0), nsexes)
if(is.numeric(gears))
gears <- paste("Gear", gears)
if(is.numeric(surveys))
surveys <- paste("Survey", surveys)
Sel <- data.frame(Series=
c(rep(gears,each=nsexes*nages),
rep(surveys,each=nsexes*nages),
rep("Maturity",nsexes*nages)),
Sex=rep(rep(sexes,each=nages),ngears+nsurveys+1),
Age=rep(ages,(ngears+nsurveys+1)*nsexes),
P=c(t(com),t(srv),mat),
stringsAsFactors=FALSE)
if(!quiet) cat("OK\n")
Sel
}
getDev <- function(ages, years)
{
if(!quiet) cat("Dev ")
Dev <- list()
Dev$sigmaR <- c(readVector("p_log_InitialDev",same.line=TRUE)[6],
readVector("p_log_RecDev",same.line=TRUE)[6])
names(Dev$sigmaR) <- c("Initial", "Annual")
Dev$Initial <- readVector("log_InitialDev", same.line=TRUE)
names(Dev$Initial) <- ages[-c(1,length(ages))] # exclude first and last age
Dev$Annual <- readVector("log_RecDev", same.line=TRUE)
names(Dev$Annual) <- years[-length(years)] # exclude last year
if(!quiet) cat("OK\n")
Dev
}
getCPUE <- function(gears, years)
{
if(!quiet) cat("CPUE ")
nseries <- readVector("NCPUEindex")
ngears <- length(gears)
nyears <- length(years)
obs <- readMatrix("indexmethodyearvaluecv",
nrow=readVector("Number_of_CPUE_data",same.line=FALSE))
obs <- data.frame(Series=obs[,1], Gear=obs[,2], Year=obs[,3], Obs=obs[,4],
CV=obs[,5], stringsAsFactors=FALSE)
fit <- readMatrix("CPUE_Index_Trajectories", nrow=nseries)
fit <- data.frame(Series=rep(1:nseries,each=nyears),
Year=rep(years,nseries), Fit=as.vector(t(fit)),
stringsAsFactors=FALSE)
## Merge without looking at gear
CPUE <- merge(obs[,names(obs)!="Gear"], fit, all=TRUE)
sgkey <- unique(obs[,c("Series","Gear")])
CPUE <- merge(sgkey, CPUE) # add gear column
CPUE <- data.frame(Series=paste("Series ",CPUE$Series,"-",CPUE$Gear,sep=""),
Year=as.integer(CPUE$Year), Obs=CPUE$Obs, CV=CPUE$CV,
Fit=CPUE$Fit, stringsAsFactors=FALSE)
if(!quiet) cat("OK\n")
CPUE
}
getSurvey <- function(years)
{
if(!quiet) cat("Survey ")
nyears <- length(years)
nseries <- readVector("Nsurveyindex")
obs <- readMatrix("indexyearvaluecv",
nrow=readVector("Number_of_survey_data",same.line=FALSE))
obs <- data.frame(Series=obs[,1], Year=obs[,2], Obs=obs[,3], CV=obs[,4])
fit <- readMatrix("Survey_Index_Trajectories", nrow=nseries)
fit <- data.frame(Series=rep(1:nseries,each=nyears),
Year=rep(years,nseries), Fit=as.vector(t(fit)),
stringsAsFactors=FALSE)
Survey <- merge(obs, fit, all=TRUE)
Survey$Series <- as.integer(Survey$Series)
Survey$Year <- as.integer(Survey$Year)
if(!quiet) cat("OK\n")
Survey
}
getCAc <- function(sexes, ages)
{
if(!quiet) cat("CAc ")
nsexes <- length(sexes)
nages <- length(ages)
nobs <- readVector("Number_of_Commercial_C@A", same.line=FALSE)
obs <- readMatrix("methodyearsamplesizesex1a1sex1a2sex1a3", nrow=nobs)
fit <- readMatrix("methodyearsamplesizesex1a1sex1a2sex1a3", nrow=nobs,
header=2*(nobs+1))
CAc <- data.frame(Series=rep(obs[,1],each=nsexes*nages),
Year=rep(obs[,2],each=nsexes*nages),
SS=rep(obs[,3],each=nsexes*nages),
Sex=rep(rep(sexes,each=nages),nobs),
Age=rep(ages,nsexes*nobs),
Obs=as.vector(t(obs[,-(1:3)])),
Fit=as.vector(t(fit)),
stringsAsFactors=FALSE)
CAc$Series <- as.integer(CAc$Series)
CAc$Year <- as.integer(CAc$Year)
CAc$Age <- as.integer(CAc$Age)
if(!quiet) cat("OK\n")
CAc
}
getCAs <- function(sexes, ages)
{
if(!quiet) cat("CAs ")
nsexes <- length(sexes)
nages <- length(ages)
nobs <- readVector("Number_of_survey_C@A",same.line=FALSE)
obs <- readMatrix("surveyyearsamplesizesex1a1sex1a2sex1a3", nrow=nobs)
fit <- readMatrix("surveyyearsamplesizesex1a1sex1a2sex1a3", nrow=nobs,
header=2*(nobs+1))
CAs <- data.frame(Series=rep(obs[,1],each=nsexes*nages),
Year=rep(obs[,2],each=nsexes*nages),
SS=rep(obs[,3],each=nsexes*nages),
Sex=rep(rep(sexes,each=nages),nobs),
Age=rep(ages,nsexes*nobs),
Obs=as.vector(t(obs[,-(1:3)])),
Fit=as.vector(t(fit)),
stringsAsFactors=FALSE)
CAs$Series <- as.integer(CAs$Series)
CAs$Year <- as.integer(CAs$Year)
CAs$Age <- as.integer(CAs$Age)
if(!quiet) cat("OK\n")
CAs
}
getCLc <- function(sexes, lengths)
{
if(!quiet) cat("CLc ")
nsexes <- length(sexes)
nlengths <- length(lengths)
nobs <- readVector("Number_of_Commercial_C@L", same.line=FALSE)
obs <- readMatrix("methodyearsamplesizesex1l1sex1l2sex1l3", nrow=nobs)
fit <- readMatrix("methodyearsamplesizesex1l1sex1l2sex1l3", nrow=nobs,
header=nobs+1)
CLc <- data.frame(Series=rep(obs[,1],each=nsexes*nlengths),
Year=rep(obs[,2],each=nsexes*nlengths),
SS=rep(obs[,3],each=nsexes*nlengths),
Sex=rep(rep(sexes,each=nlengths),nobs),
Length=rep(lengths,nsexes*nobs),
Obs=as.vector(t(obs[,-(1:3)])),
Fit=as.vector(t(fit)),
stringsAsFactors=FALSE)
CLc$Series <- as.integer(CLc$Series)
CLc$Year <- as.integer(CLc$Year)
CLc$Length <- as.integer(CLc$Length)
if(!quiet) cat("OK\n")
CLc
}
getCLs <- function(sexes, lengths)
{
if(!quiet) cat("CLs ")
nsexes <- length(sexes)
nlengths <- length(lengths)
nobs <- readVector("Number_of_surveyC@L",same.line=FALSE)
obs <- readMatrix("surveyyearsamplesizesex1l1sex1l2sex1l3", nrow=nobs)
fit <- readMatrix("surveyyearsamplesizesex1l1sex1l2sex1l3", nrow=nobs,
header=2*(nobs+1))
CLs <- data.frame(Series=rep(obs[,1],each=nsexes*nlengths),
Year=rep(obs[,2],each=nsexes*nlengths),
SS=rep(obs[,3],each=nsexes*nlengths),
Sex=rep(rep(sexes,each=nlengths),nobs),
Length=rep(lengths,nsexes*nobs),
Obs=as.vector(t(obs[,-(1:3)])),
Fit=as.vector(t(fit)),
stringsAsFactors=FALSE)
CLs$Series <- as.integer(CLs$Series)
CLs$Year <- as.integer(CLs$Year)
CLs$Length <- as.integer(CLs$Length)
if(!quiet) cat("OK\n")
CLs
}
getLA <- function(sexes, ages)
{
## Sex | Age | Obs | Fit | CV
if(!quiet) cat("LA ")
nsexes <- length(sexes)
nages <- length(ages)
nobs <- readVector("#femalesmales", same.line=FALSE, file=latage.file,
vector=latage.vector) # two elements
obs <- readMatrix("VonBertalanfy--Lenght-at-agefit--Likelihood",
nrow=sum(nobs), header=8)
obs <- data.frame(Sex=rep(sexes,nobs), Age=obs[,1], Obs=obs[,2],
stringsAsFactors=FALSE)
Linf <- suppressWarnings(readVector("VonBeratalanfy:Linf")[-(1:3)])
K <- suppressWarnings(readVector("VonBeratalanfy:k")[-(1:3)])
t0 <- suppressWarnings(readVector("VonBeratalanfy:to")[-(1:3)])
CV1 <- suppressWarnings(readVector("cvoftheFitbysex")[-(1:5)])
CVratio <- suppressWarnings(
readVector("ratioofcv(L_an)/cv(L_a1)oftheFitbysex")[-(1:7)])
sigmaLA <- readVector(
"#LinearrelationshipofsigmaL@A(1=age;2=length---ignoreifW@Aissupplied)",
same.line=FALSE, file=txt.file, vector=txt.vector)[1]
max.age <- c(max(obs$Age[obs$Sex==sexes[1]]),
max(obs$Age[obs$Sex==sexes[2]])) # max ages in LA data
fit <- data.frame(Sex=rep(sexes,max.age), Age=c(1:max.age[1],1:max.age[2]),
stringsAsFactors=FALSE)
fit$Fit[fit$Sex==sexes[1]] <- Linf[1]*
(1-exp(-K[1]*(fit$Age[fit$Sex==sexes[1]]-t0[1])))
fit$Fit[fit$Sex==sexes[2]] <- Linf[2]*
(1-exp(-K[2]*(fit$Age[fit$Sex==sexes[2]]-t0[2])))
if(sigmaLA == 1) # CV ~ age
{
A <- rep(max(ages), 2)
a <- cbind(fit$Age[fit$Sex==sexes[1]], fit$Age[fit$Sex==sexes[2]])
fit$CV[fit$Sex==sexes[1]] <- CV1[1] +
CV1[1]*(CVratio[1]-1)/(A[1]-1)*(a[,1]-1)
fit$CV[fit$Sex==sexes[2]] <- CV1[2] +
CV1[2]*(CVratio[2]-1)/(A[2]-1)*(a[,2]-1)
}
if(sigmaLA == 2) # CV ~ length
{
L1 <- Linf*(1-exp(-K*(1-t0)))
Ln <- Linf*(1-exp(-K*(max(ages)-t0)))
fit$CV[fit$Sex==sexes[1]] <- CV1[1] +
CV1[1]*(CVratio[1]-1)/(Ln[1]-L1[1])*(fit$Fit[fit$Sex==sexes[1]]-L1[1])
fit$CV[fit$Sex==sexes[2]] <- CV1[2] +
CV1[2]*(CVratio[2]-1)/(Ln[2]-L1[2])*(fit$Fit[fit$Sex==sexes[2]]-L1[2])
}
LA <- merge(obs, fit, by=c("Sex","Age"), all=TRUE)
LA$Age <- as.integer(LA$Age)
LA$Fit <- LA$Fit
LA$CV <- LA$CV
if(!quiet) cat("OK\n")
LA
}
## 2 Parse args
if(!file.exists(res.file))
stop("file ", res.file, " not found; use / or \\\\ separators")
## 3 Read dimensions
res.vector <- readLines(res.file) # string vector, one element being one line
## Remove white space and quotes
res.vector <- gsub("\"","",
gsub("\t","",gsub(" ","",res.vector)))
if(!quiet) cat("\nParsing text file ", res.file, ":\n\nPreamble ", sep="")
sexes <- if(readVector("Nsexes")==1) "Unisex" else c("Female","Male")
## Vectors 'gears' and 'surveys' are possibly of length zero
gears <- seq(1, length.out=readVector("Nmethods"))
surveys <- seq(1, length.out=readVector("Nsurveyindex"))
years <- seq(from=readVector("StartYear"), to=readVector("EndYear")+1)
ages <- seq(from=1, to=readVector("Nages"))
lengths <- seq(from=readVector("First_length"),
by=readVector("Length_class_increment"),
length.out=readVector("Number_of_length_classes"))
if(!quiet) cat("OK\n")
## 4 Read N, B, R, Sel
model <- list()
model$N <- getN(sexes, years, ages)
model$B <- getB(years, gears) # Recruits:
rec <- model$N[model$N$Age==1,] # by year and sex
rec <- tapply(rec$N, rec$Year, sum) # combine sexes
model$B$R <- c(rec[-1], NA) # align with spawning stock
model$Sel <- getSel(gears, surveys, years, sexes, ages)
## 5 Read Dev, CPUE, Survey, CAc, CAs, CLc, CLs, LA
if(Dev) model$Dev <- getDev(ages, years)
if(CPUE) model$CPUE <- getCPUE(gears, years)
if(Survey) model$Survey <- getSurvey(years)
if(CAc) model$CAc <- getCAc(sexes, ages)
if(CAs) model$CAs <- getCAs(sexes, ages)
if(CLc) model$CLc <- getCLc(sexes, lengths)
if(CLs) model$CLs <- getCLs(sexes, lengths)
if(LA)
{
latage.file <- paste(dirname(res.file),"l_at_age.dat",sep="/")
if(!file.exists(latage.file))
stop("file ", latage.file, " not found; use / or \\\\ separators")
latage.vector <- readLines(latage.file) # LAnobs for each sex
latage.vector <- gsub("\"","", gsub("\t","",gsub(" ","",latage.vector)))
txt.file <- gsub("\\.res", "\\.txt", res.file)
if(!file.exists(txt.file))
stop("file ", txt.file, " not found; use / or \\\\ separators")
txt.vector <- readLines(txt.file) # sigmaLA~x
txt.vector <- gsub("\"","", gsub("\t","",gsub(" ","",txt.vector)))
model$LA <- getLA(sexes, ages)
}
if(!quiet) cat("\n")
## 6 Create attributes
attr(model,"call") <- match.call()
class(model) <- "scape"
model
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.