#' Plot tag fits
#'
#' Tag plots when there are 2 tag files
#'
#' @export
#'
Tag <- function(stock, source.dir, target.dir = source.dir)
{
# Read files
allFiles <- list.files(source.dir)
i <- grep("TagResids.out", allFiles)
Files <- allFiles[i]
nFiles <- length(Files)
# Combine data
tag <- NULL
for (I in 1:nFiles)
{
tmp <- read.table(paste(source.dir, "/", Files[I], sep = ""), header = TRUE, as.is = TRUE)
tmp <- data.frame(file = I, tmp)
tag <- rbind(tag, tmp)
}
tag2 <- 1
nmoult <- (tag$perrec - tag$perrel)
nmoult[tag$sex==2] <- 1+floor(nmoult[tag$sex==2]/2)-(tag$perrec[tag$sex==2]+1)%%2
nmoult[tag$sex==2 & (tag$perrel+1) %% 2==(tag$perrec+1)%%2] <- (tag$perrec[tag$sex==2 & (tag$perrel+1)%%2==(tag$perrec+1)%%2]-tag$perrel[tag$sex==2 & (tag$perrel+1)%%2==(tag$perrec+1)%%2])/2
#write.table(cbind(tag,nmoult),"moulttest.csv",sep=",",row.names=F)
#tag$sex[tag$sex == 1] <- "Male"
#tag$sex[tag$sex == 2] <- "Female"
# ==============================================================================
# StdResid vs Pred
# ==============================================================================
PlotType(paste(target.dir, "/", stock, "TagResid", sep = ""),
width = 2 * PlotOptions$plotsize[1],
height = 10 + nFiles * PlotOptions$plotsize[2])
#ggplot(data = tag, aes(x = size2Pred, y = StdRes)) +
# geom_point() + facet_grid(.~sex) + labs()
par(mfrow = c(nFiles,2), oma = c(2,2,1,1), mar=c(4,2,1,1), mgp=c(4,1,0), cex = 3)
for ( ff in 1:nFiles )
{
for( sex in 1:2 )
{
if ( sex==1 )
{
xlim <- TagOptions$xlim.male
} else {
xlim <- TagOptions$xlim.female
}
plot(StdRes ~ size2Pred, data = tag[tag$sex==sex & tag$file==ff,],
pch = 1, ylab = "", xlab = "", xaxt = "s", cex = 0.7, xlim = xlim, las = 1)
abline(h = 0, lty = 2)
mtext(c("Male","Female")[sex], side = 3, line = 0.1, las = 1, adj = 0, cex = 2)
}
}
mtext("Predicted", side = 1, line = -1, outer = TRUE, cex = 2)
mtext("Standardised residual", side = 2, line = 1, outer = TRUE, cex = 2)
if ( PlotOptions$Captions) mtext(paste(source.dir," ",stock,": Tags Plot"),
side = 1, line = 1, outer = TRUE, cex = 2)
dev.off()
# ==============================================================================
# Plot obs vs pred
# ==============================================================================
PlotType(paste(target.dir, "/", stock, "TagPredObs", sep = ""),
width = 2 * PlotOptions$plotsize[1],
height = 10 + nFiles * PlotOptions$plotsize[2])
par(mfrow = c(nFiles,2), oma = c(2,2,1,1), mar=c(4,3,1,1), mgp=c(4,1,0), cex = 2)
for ( ff in 1:nFiles )
{
for ( sex in 1:2 )
{
if ( sex == 1 )
{
xlim <- TagOptions$xlim.male
} else {
xlim <- TagOptions$xlim.female
}
plot(size2Pred~size2Obs,data=tag[tag$sex==sex & tag$file==ff,],
pch = 1, ylab = "", xlab = "", xaxt = "s", xlim = xlim, cex = 2, las = 1)
abline(a=0,b=1)
mtext(c("Male","Female")[sex], side = 3, line = 0.1, las = 1, adj = 0, cex = 3)
}
}
mtext("Observed", side = 1, line = -1, outer = TRUE, cex = 3)
mtext("Predicted", side = 2, line = -0.4, outer = TRUE, cex = 3)
if ( PlotOptions$Captions ) mtext(paste(source.dir," ",stock,": Increment vs size at release plot"),
side = 1, line = 1, outer = TRUE, cex = 2)
dev.off()
# ==============================================================================
# QQ normal plots by sex
# ==============================================================================
PlotType(paste(target.dir, "/", stock, "TagQQ", sep = ""),
width = 2 * PlotOptions$plotsize[1],
height = 10 + nFiles * PlotOptions$plotsize[2])
par(mfrow = c(nFiles,2), oma = c(2,2,1,1), mar=c(4,2,1,1), mgp=c(4,1,0), cex = 2)
for ( ff in 1:nFiles )
{
for ( sex in 1:2 )
{
qqnorm(tag$StdRes[tag$sex==sex & tag$file==ff],xlab="",ylab="",main="",las=1, cex = 2)
abline(a = 0, b = 1)
abline(h=quantile(tag$StdRes,p=c(5,10,25,50,75,90,95)/100),lty=2)#Plot quntiles of the stdevs
mtext(c("Male","Female")[sex],side=3,line=0.1,las=1,adj=0, cex = 2.5)
}
}
mtext("Theoretical quantiles", side = 1, line = -1, outer = TRUE, cex = 2)
mtext("Standardised residual", side = 2, line = 1, outer = TRUE, cex = 2)
if ( PlotOptions$Captions) mtext(paste(source.dir," ",stock,": Tag resid Plot"),
side = 1, line = 1, outer = TRUE, cex = 2)
dev.off()
# ==============================================================================
# Size
# ==============================================================================
tag$sc <- floor(tag$size2Pred / 4) * 4 + 2
if ( TagOptions$StdResSize )
{
PlotType(paste(target.dir, "/", stock, "TagResSize", sep = ""),
width = 2 * PlotOptions$plotsize[1],
height = 10 + nFiles * PlotOptions$plotsize[2])
par(mfrow = c(nFiles,2), oma = c(2,2,1,1), mar=c(4,3,1,1), mgp=c(4,1,0), cex = 2.5)
for (ff in 1:nFiles)
{
for ( sex in 1:2 )
{
nn <- nrow(tag[tag$sex == sex & tag$file == ff,])
boxplot(StdRes ~ sc, data = tag[tag$sex==sex & tag$file==ff,],
ylab = "", xlab = "", line = 2, medlwd = 1, las = 1,
varwidth = TRUE)
abline(h = 0, lty = 2)
legend("topright", legend = paste("N: ", nn, sep = ""), bty = "n", cex = 1)
mtext(c("Male","Female")[sex], side = 3, line = 0.1,
las = 1, adj = 0, cex = 2.5)
}
}
mtext("Size class", side = 1, line = -1, outer = TRUE, cex = 3)
mtext("Standardised residual", side = 2, line = -0.4, outer = TRUE, cex = 3)
if ( PlotOptions$Captions ) mtext(paste(source.dir," ",stock,": Standardised residuals by 4mm size classes for each sex"),
side = 1, line = 1, outer = TRUE, cex = 2)
dev.off()
}
# ==============================================================================
# Release Count
# ==============================================================================
if ( TagOptions$StdResRel )
{
PlotType(paste(target.dir, "/", stock, "TagResRel", sep = ""),
width = 2 * PlotOptions$plotsize[1],
height = 10 + nFiles * PlotOptions$plotsize[2])
par(mfrow = c(nFiles,2), oma = c(2,2,1,1), mar=c(4,3,1,1), mgp=c(4,1,0), cex = 2.5)
for (ff in 1:nFiles)
{
for(sex in 1:2)
{
nn <- nrow(tag[tag$sex == sex & tag$file == ff,])
boxplot(StdRes~release,data=tag[tag$sex==sex & tag$file==ff,],
ylab="",xlab="",line=2, medlwd = 1, varwidth = TRUE, las = 1)
abline(h=0,lty=2)
legend("topright", legend = paste("N: ", nn, sep = ""), bty = "n", cex = 1.5)
mtext(c("Male","Female")[sex],side=3,line=0.1,las=1,adj=0, cex = 2.5)
}
mtext("Release", side = 1, line = -1, outer = TRUE, cex = 3)
mtext("Standardised residual", side = 2, line = -0.4, outer = TRUE, cex = 3)
}
if(PlotOptions$Captions) mtext(paste(source.dir," ",stock,": Standardised residuals by release number for each sex"),
side = 1, line = 1, outer = TRUE, cex = 2)
dev.off()
}
# ==============================================================================
# Area
# ==============================================================================
if( TagOptions$StdResArea )
{
PlotType(paste(target.dir, "/", stock, "TagResArea", sep = ""),
width = 2 * PlotOptions$plotsize[1],
height = 10 + nFiles * PlotOptions$plotsize[2])
par(mfrow = c(nFiles,2), oma = c(2,2,1,1), mar=c(4,3,1,1), mgp=c(4,1,0), cex = 2)
for (ff in 1:nFiles)
{
for( sex in 1:2 )
{
nn <- nrow(tag[tag$sex == sex & tag$file == ff,])
boxplot(StdRes ~ area, data = tag[tag$sex == sex & tag$file == ff,],
ylab = "", xlab = "", line = 2,
medlwd = 1, las = 1, varwidth = TRUE)
abline(h = 0, lty = 2)
legend("top", legend = paste("N: ", nn, sep = ""), bty = "n", cex = 1)
mtext(c("Male","Female")[sex], side = 3,
line = 0.1, las = 1, adj = 0, cex = 3)
}
}
mtext("Area", side = 1, line = -1, outer = TRUE, cex = 3)
mtext("Standardised residual", side = 2, line = -0.4, outer = TRUE, cex = 3)
if ( PlotOptions$Captions ) mtext(paste(source.dir," ",stock,": Standardised residuals by area for each sex"),
side = 1, line = 1, outer = TRUE, cex = 2)
dev.off()
}
#Condition
#if(TagOptions$StdResCond==TRUE){
# par(las=1,mfrow=c(2,2),oma=c(1,1,1,1),mar=c(3,3,1,1),mgp=c(2,1,0))
# for(sex in 1:2){
# boxplot(StdRes~cond,data=tag[tag$sex==sex & tag$file==1,],ylab="Standardised residual",xlab="Condition",line=2,medlwd=1)
# abline(h=0,lty=2)
# mtext(c("Male (File 1)","Female (File 1)")[sex],side=3,line=0.1,las=1,adj=0)
# }
#for(sex in 1:2){
# boxplot(StdRes~cond,data=tag[tag$sex==sex & tag$file==2,],ylab="Standardised residual",xlab="Condition",line=2,medlwd=1)
# abline(h=0,lty=2)
# mtext(c("Male (File 2)","Female (File 2)")[sex],side=3,line=0.1,las=1,adj=0)
# }
# if(Captions) mtext(paste(source.dir," ",stock,": Standardised residuals by condition code for each sex"),side=1,line=0,outer=TRUE,cex=0.7)
# savePlot(filename=paste(target.dir,"/",stock,"TagResCond.wmf",sep=""),type=plottype,device=dev.cur())
#}
#Type
#if(TagOptions$StdResType==TRUE){
# par(las=1,mfrow=c(2,2),oma=c(1,1,1,1),mar=c(3.5,3,1,1),mgp=c(2,1,0))
# for(sex in 1:2){
# type<-factor(c("NULL","WRL","SPYR","TP","HALL","Disc")[tag$type[tag$sex==sex]])
# boxplot(tag$StdRes[tag$sex==sex & tag$file==1]~type,ylab="Standardised residual",xlab="Type",line=2,medlwd=1)
# abline(h=0,lty=2)
# mtext(c("Male (File 1)","Female (File 1)")[sex],side=3,line=0.1,las=1,adj=0)
# }
# for(sex in 1:2){
# type<-factor(c("NULL","WRL","SPYR","TP","HALL","Disc")[tag$type[tag$sex==sex]])
# boxplot(tag$StdRes[tag$sex==sex & tag$file==2]~type,ylab="Standardised residual",xlab="Type",line=2,medlwd=1)
# abline(h=0,lty=2)
# mtext(c("Male (File 2)","Female (File 2)")[sex],side=3,line=0.1,las=1,adj=0)
# }
# if(Captions) mtext(paste(source.dir," ",stock,": Standardised residuals by tag type for each sex"),side=1,line=0,outer=TRUE,cex=0.7)
# savePlot(filename=paste(target.dir,"/",stock,"TagResType.wmf",sep=""),type=plottype,device=dev.cur())
#}
# ==============================================================================
# Histograms by sex
# ==============================================================================
if(TagOptions$Reshist)
{
PlotType(paste(target.dir, "/", stock, "TagBHistSex", sep = ""),
width = 2 * PlotOptions$plotsize[1],
height = 10 + nFiles * PlotOptions$plotsize[2])
par(mfrow = c(nFiles,2), oma = c(2,2,1,1), mar=c(4,2,1,1), mgp=c(4,1,0), cex = 3)
for ( ff in 1:nFiles )
{
for( sex in 1:2 )
{
sub<-tag[tag$sex==sex & tag$file==ff,]
brmin<-floor(min(sub$StdRes))
brmax<-ceiling(max(sub$StdRes))
by<-0.1
br<-seq(brmin-by/2,brmax+by/2,by=by)
aa<-hist(sub$StdRes,br=br,right=F,plot=F)
plot(aa$mids,aa$counts/sum(aa$counts),xlim=c(brmin,brmax),
ylim=range(aa$counts/sum(aa$counts)),xlab="",ylab="", las = 1)
abline(v=0,lty=2)
x<-seq(brmin,(brmax+by),by=by)
#a<-pnorm(x, mean=0, sd=1, lower.tail = TRUE, log.p = FALSE)
a<-pnorm(x, mean=mean(sub$StdRes), sd=sqrt(var(sub$StdRes)), lower.tail = TRUE, log.p = FALSE)
b<-rep(0,length(a))
b[2:length(a)]<-a[2:length(a)]-a[1:(length(a)-1)]
lines(br,b,col="blue",lwd=2)
mtext(c("Male","Female")[sex], side = 3, line = 0.1, las = 1, adj = 0, cex = 2.5)
}
}
mtext("Standardised residual", side = 1, line = -1, outer = TRUE, cex = 3)
mtext("Frequency", side = 2, line = 1, outer = TRUE, cex = 3)
if ( PlotOptions$Captions ) mtext(paste(source.dir," ",stock,": Histograms of residuals for tag by sex"),
side = 1, line = 1, outer = TRUE, cex = 3)
dev.off()
}
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.