Hello,
I am attempting to fit a non-linear model (Von Bertalanffy growth model)
to fish length-at-age data with the purpose of determining if any of the
three parameters differ between male and female fish. I believe that I
can successfully accomplish this goal assuming an additive error
structure (illustrated in section 1 below). My trouble begins when I
attempt this analysis using a model with a multiplicative error
structure. I believe that I can fit the multiplicative error model to
data without attempting to test between males and females (section 2
below), but I am not sure how to declare the model formula if I want to
test between males and females (section 3 below). In addition, I
assumed what I did in section 2 was correct, separated the data into a
dataframe of males and a dataframe of females, and then attempted to fit
separate models to each group (Section 4). This was "successful" for
males but not for females. The problem with the female group appears to
be that my "logging" creates NaN cells. When I put a trace on the nls
function it appears to do one iteration and then fails. My assumption
is that the nls function must be attempting different parameter
estimates that return a NaN but I cannot test this as I'm not sure how
to see the new parameters that it is trying. [I tried different
starting values but did not find any that corrected this problem.]
So (1) how do I declare the model formula for using multiplicative
errors (i.e., is what I did in Section 2 correct), (2) how do I declare
the model formula for comparing two groups and using multiplicative
errors, and (3) any suggestions for how to find the "issue" leading to
the error for just the female model in section 4.
Sys.info information is in Section 4 and I am using R 2.6.1. I have not
included the data as it is a rather large file.
Thank you very much for reading this long message and for any help that
you can offer.
### Section 1 ###
> library(xlsReadWrite)
> fwd <- read.xls("FWDrum_SDF.xls",sheet=1)
> fwd$Gm <- fwd$Gf <- rep(0,length(fwd$age))
> fwd$Gf[fwd$sex=="female"] <- 1
> fwd$Gm[fwd$sex=="male"] <- 1
> str(fwd)
'data.frame': 719 obs. of 5 variables:
$ age: num 1.27 2.25 2.25 2.25 2.25 ...
$ tl : num 226 208 226 226 234 ...
$ sex: Factor w/ 2 levels "female","male": 1 1 1 1 1 1 1 1
1 1 ...
$ Gf : num 1 1 1 1 1 1 1 1 1 1 ...
$ Gm : num 0 0 0 0 0 0 0 0 0 0 ...
# starting values
> sLinf <- 405
> sK <- 0.11
> sto <- -5.2
> svb <- list(Linf=sLinf,K=sK,to=sto)
# Fit the additive error structure model to both groups combined
> vbla <- nls(tl~Linf*(1-exp(-K*(age-to))),start=svb,data=fwd)
> summary(vbl1)
Formula: tl ~ Linf * (1 - exp(-K * (age - to)))
Parameters:
Estimate Std. Error t value Pr(>|t|)
Linf 520.751048 19.560552 26.623 < 2e-16 ***
K 0.061097 0.008506 7.183 1.72e-12 ***
to -6.950346 1.112839 -6.246 7.24e-10 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05
'.' 0.1 ' ' 1
Residual standard error: 83.75 on 716 degrees of freedom
Number of iterations to convergence: 9
Achieved convergence tolerance: 7.525e-06
# Fit the additive error structure model to the separate groups
# compare this to vbla to see if any parameters differ (fit other
models to see which parameter)
> vbla.gen<- nls(tl~Gf*Linff*(1-exp(-Kf*(age-tof))) +
Gm*Linfm*(1-exp(-Km*(age-tom))),
+
start=list(Linff=sLinf,Kf=sK,tof=sto,Linfm=sLinf,Km=sK,tom=sto),data=fwd
)
> summary(vbla.gen)
Formula: tl ~ Gf * Linff * (1 - exp(-Kf * (age - tof))) + Gm * Linfm *
(1 - exp(-Km * (age - tom)))
Parameters:
Estimate Std. Error t value Pr(>|t|)
Linff 692.719734 47.191732 14.679 < 2e-16 ***
Kf 0.038532 0.006469 5.957 4.04e-09 ***
tof -7.589941 1.235644 -6.142 1.35e-09 ***
Linfm 372.224132 13.268491 28.053 < 2e-16 ***
Km 0.095814 0.024983 3.835 0.000137 ***
tom -8.058334 2.445654 -3.295 0.001033 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05
'.' 0.1 ' ' 1
Residual standard error: 68.15 on 713 degrees of freedom
Number of iterations to convergence: 8
Achieved convergence tolerance: 6.589e-06
### Section 2 ###
> vblm <- nls(log(tl)~log(Linf*(1-exp(-K*(age-to)))),start=svb,data=fwd)
> summary(vblm)
Formula: log(tl) ~ log(Linf * (1 - exp(-K * (age - to))))
Parameters:
Estimate Std. Error t value Pr(>|t|)
Linf 515.414166 21.239627 24.267 < 2e-16 ***
K 0.055431 0.007486 7.404 3.71e-13 ***
to -8.289902 1.002402 -8.270 6.51e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05
'.' 0.1 ' ' 1
Residual standard error: 0.1907 on 716 degrees of freedom
Number of iterations to convergence: 8
Achieved convergence tolerance: 2.236e-06
### Section 3 ###
# FIRST TRY -- Fit the multiplicative error structure model to the
separate groups
> vblm.gen1<- nls(log(tl)~log(Gf*Linff*(1-exp(-Kf*(age-tof))) +
Gm*Linfm*(1-exp(-Km*(age-tom)))),
start=list(Linff=sLinf,Kf=sK,tof=sto,Linfm=sLinf,Km=sK,tom=sto),data=fwd
)
Error in numericDeriv(form[[3]], names(ind), env) :
Missing value or an infinity produced when evaluating the model
In addition: Warning message:
In log(Gf * Linff * (1 - exp(-Kf * (age - tof))) + Gm * Linfm * :
NaNs produced
# SECOND TRY -- Fit the multiplicative error structure model to the
separate groups
> vblm.gen2<- nls(log(tl)~log(Gf*Linff*(1-exp(-Kf*(age-tof)))) +
log(Gm*Linfm*(1-exp(-Km*(age-tom)))),
start=list(Linff=sLinf,Kf=sK,tof=sto,Linfm=sLinf,Km=sK,tom=sto),data=fwd
)
Error in numericDeriv(form[[3]], names(ind), env) :
Missing value or an infinity produced when evaluating the model
> THIRD TRY -- Fit the multiplicative error structure model to the
separate groups
> vblm.gen3<- nls(log(tl)~Gf*log(Linff*(1-exp(-Kf*(age-tof)))) +
Gm*log(Linfm*(1-exp(-Km*(age-tom)))),
start=list(Linff=sLinf,Kf=sK,tof=sto,Linfm=sLinf,Km=sK,tom=sto),data=fwd
)
Error in numericDeriv(form[[3]], names(ind), env) :
Missing value or an infinity produced when evaluating the model
In addition: Warning message:
In log(Linff * (1 - exp(-Kf * (age - tof)))) : NaNs produced
### Section 4 ###
# Just males
> fwd.m <- subset(fwd,sex=="male")
> vblm.m <-
nls(log(tl)~log(Linf*(1-exp(-K*(age-to)))),start=svb,data=fwd.m)
> summary(vblm.m)
Formula: log(tl) ~ log(Linf * (1 - exp(-K * (age - to))))
Parameters:
Estimate Std. Error t value Pr(>|t|)
Linf 363.07557 4.88130 74.381 < 2e-16 ***
K 0.11765 0.01184 9.936 < 2e-16 ***
to -6.29861 0.75724 -8.318 4.15e-15 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05
'.' 0.1 ' ' 1
Residual standard error: 0.09035 on 276 degrees of freedom
Number of iterations to convergence: 5
Achieved convergence tolerance: 7.262e-06
# Just females
> fwd.f <- subset(fwd,sex=="female")
> vblm.f <-
nls(log(tl)~log(Linf*(1-exp(-K*(age-to)))),start=svb,data=fwd.f)
26.89721 : 405.00 0.11 -5.20
Error in numericDeriv(form[[3]], names(ind), env) :
Missing value or an infinity produced when evaluating the model
In addition: Warning message:
In log(Linf * (1 - exp(-K * (age - to)))) : NaNs produced
### Section 5 ###
> Sys.info()
sysname release
version
"Windows" "NT 5.1"
"(build
2600) Service Pack 2"
nodename machine
login
"CSE229-001" "x86"
[[alternative HTML version deleted]]
