Howdy,

Last week I got some great help on why I was getting an error code when trying to run this model, thanks everyone!  I was able to get the code up and running beautifully for several data sets.  Now I am getting different errors with this new data set.  I can't figure out why, I have more data points with this species, and it is ordered exactly the same as the other species I have been working with.  

The two errors I get are for these specific lines of code: 

This error I actually got with all the data sets but its just assumption checking so I did this manually. 

> fitGen <- nls(vbGen,data=BLG,start=svGen)
Error in nlsModel(formula, mf, start, wts) : 
  singular gradient matrix at initial parameter estimates

and 

This code is only not working with this data set. 

> fitK <- nls(vbK,data=BLG,start=svK)
Error in numericDeriv(form[[3L]], names(ind), env) : 
  Missing value or an infinity produced when evaluating the model

Thanks,
April


Here is all my code and my data set: 

Data set: 

structure(list(Length = c(60.96, 60.96, 63.5, 63.5, 66.04, 66.04, 
71.12, 71.12, 96.52, 96.52, 96.52, 101.6, 101.6, 106.68, 111.76, 
116.84, 116.84, 119.38, 124.46, 124.46, 134.62, 137.16, 144.78, 
147.32, 149.86, 154.94, 157.48, 167.64, 172.72, 175.26, 175.26, 
180.34, 180.34, 182.88, 190.5, 193.04, 193.04, 193.04, 195.58, 
198.12, 203.2, 210.82), Age = c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 
4L, 3L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 6L, 3L, 5L, 5L, 6L, 5L, 7L, 
6L, 6L), Cohort = c(2007L, 2009L, 2004L, 2005L, 2008L, 2011L, 
2006L, 2010L, 2006L, 2007L, 2008L, 2004L, 2009L, 2003L, 2010L, 
2005L, 2007L, 2008L, 2005L, 2006L, 2003L, 2009L, 2006L, 2008L, 
2007L, 2004L, 2002L, 2005L, 2004L, 2003L, 2007L, 2005L, 2006L, 
2006L, 2002L, 2002L, 2004L, 2005L, 2003L, 2005L, 2004L, 2002L
), Year = c(2008L, 2010L, 2005L, 2006L, 2009L, 2012L, 2007L, 
2011L, 2008L, 2009L, 2010L, 2006L, 2011L, 2005L, 2012L, 2007L, 
2010L, 2011L, 2008L, 2009L, 2006L, 2012L, 2010L, 2012L, 2011L, 
2007L, 2006L, 2009L, 2008L, 2007L, 2012L, 2010L, 2011L, 2012L, 
2005L, 2007L, 2009L, 2011L, 2008L, 2012L, 2010L, 2008L)), .Names = c("Length", 
"Age", "Cohort", "Year"), class = "data.frame", row.names = c(NA, 
-42L))


Code: 

BLG=read.csv("BLGa.csv", header=T)
attach(BLG)
names(BLG)
dput(BLG,"")

$ Von Bertalanffy Growth Model$

source("http://www.rforge.net/NCStats/InstallNCStats.R")
source("http://www.rforge.net/FSA/InstallFSA.R")

library(FSA)
library(NCStats)
library(nlstools)

# provide starting values
svCom <- list(Linf=199.8132, K=0.5631745, t0=-1.960147)
svCom <- vbStarts(Length~Age,data=BLG)

# repeat the starting values for each Cohort
svGen <- lapply(svCom, rep, length(unique(Cohort)))

# define a new variable, uniquely numbered for each Cohort
Cohortf <- as.factor(Cohort)

# fit a single model with different parameter estimates for each Cohort
vbGen <- Length~Linf[Cohortf]*(1-exp(-K[Cohortf]*(Age-t0[Cohortf])))
fitGen <- nls(vbGen,data=BLG,start=svGen)

#fit a single model with the same parameter estimates using all the Cohorts data#
vbCom <- Length~Linf*(1-exp(-K*(Age-t0)))
fitCom <- nls(vbCom,data=BLG,start=svCom)
overview(fitCom)

# Checking assumptions#
   #want residual plot to look like a shot gun blast, no funnels, no trends, no arches, (ie Unbiased and homoscedastic)#
residPlot(fitCom)
   #Want histogram to be normally distributed#
hist(residuals(fitCom),main="")

# Create a model that holds Linf and t0 steady but varies by Cohort#
vbK=Length~Linf*(1-exp(-K[Cohortf]*(Age-t0)))
svK=mapply(rep,svCom,c(1,8,1))

#Fit model to the data#
fitK <- nls(vbK,data=BLG,start=svK)

#Compare the confidence intervals for each Cohort's estimates of K, remember Linf, and t0 were found for entire dataset#
overview(fitK)
overview(fitCom)

#fitted plots of Length at Age with K varying by Cohort#
plot(Length~Age,data=BLG,subset=Cohortf=="2005",pch=7,xlab="Age (yrs)",ylab="Mean Length at Capture (Length)", xlim=c(1,5), ylim=c(0,250))
points(Length~Age,data=BLG,subset=Cohortf=="2006",pch=0,col="black")
points(Length~Age,data=BLG,subset=Cohortf=="2007",pch=1,col="black")
points(Length~Age,data=BLG,subset=Cohortf=="2008",pch=2,col="black")
points(Length~Age,data=BLG,subset=Cohortf=="2009",pch=3,col="black")
points(Length~Age,data=BLG,subset=Cohortf=="20010",pch=4,col="black")
points(Length~Age,data=BLG,subset=Cohortf=="20011",pch=5,col="black")
points(Length~Age,data=BLG,subset=Cohortf=="20012",pch=6,col="black")

#Add growth curves to data set#
vbTypical <- vbFuns("typical")
curve(vbTypical(x,Linf=coef(fitK)[-c(3,4,5,6,7,8,9)]),from=1,to=5,lwd=2,lty=2,add=TRUE,xlim=c(1,5), ylim=c(0,250))
curve(vbTypical(x,Linf=coef(fitK)[-c(2,4,5,6,7,8,9)]),from=1,to=5,lwd=2,lty=3,add=TRUE)
curve(vbTypical(x,Linf=coef(fitK)[-c(2,3,5,6,7,8,9)]),from=1,to=5,lwd=2,lty=4,add=TRUE)
curve(vbTypical(x,Linf=coef(fitK)[-c(2,3,4,6,7,8,9)]),from=1,to=5,lwd=2,lty=5,add=TRUE)
curve(vbTypical(x,Linf=coef(fitK)[-c(2,3,4,5,7,8,9)]),from=1,to=5,lwd=2,lty=6,add=TRUE)
curve(vbTypical(x,Linf=coef(fitK)[-c(2,3,4,5,6,8,9)]),from=1,to=5,lwd=2,lty=99,add=TRUE)
curve(vbTypical(x,Linf=coef(fitK)[-c(2,3,4,5,6,7,9)]),from=1,to=5,lwd=2,lty=4194,add=TRUE)
curve(vbTypical(x,Linf=coef(fitK)[-c(2,3,4,5,6,7,9)]),from=1,to=5,lwd=2,lty=5995,add=TRUE)

#and a legend for the points#
legend("right",legend=c("2005","2006", "2007","2008","2009","2010","2011","2012"),pch=c(7,0,1,2,3,4,5,6),lwd=1,lty=1,cex=0.75)

#legend for the curves#
legend("bottomright",legend=c("2005","2006", "2007","2008","2009","2010","2011","2012", "Overall"),lty=c(2,3,4,5,6,99,4194,5995, 1),lwd=2,cex=0.75)

#all points on 1 graph, with 1 curve#
vbCom <- Length~Linf*(1-exp(-K*(Age-t0)))
fitCom <- nls(vbCom,data=BLG,start=svCom)
plot(Length~Age,data=BLG,pch=21,xlab="Age (yrs)",ylab="Mean Length at Capture (mm)", xlim=c(1,5), ylim=c(0,250))
curve(vbTypical(x,Linf=coef(fitCom)),from=1,to=5,lwd=4,add=TRUE)
	[[alternative HTML version deleted]]

