Hi,
I am not sure this is a bug but I can repeat it, The functions and data
are below.
I know this is nasty data, and it is very questionable whether a 4pl
model
is appropriate, but it is data fed to an automated tool and I would
have hoped for an error. Does this repeat for anyone else?
My details:
> version
_
platform i686-pc-linux-gnu
arch i686
os linux-gnu
system i686, linux-gnu
status
major 2
minor 6.0
year 2007
month 10
day 03
svn rev 43063
language R
version.string R version 2.6.0 (2007-10-03)
start=c(-1.5, 9.5, 0.09, 10.25)
names(start)<-c("A","B","xmid","scal")
gnls(response~SSllogis(conc,A,B,xmid,scal),tdat,start=start,weights=varPower(),verbose=TRUE)
**Iteration 1
GLS step: Objective: NULLvarStruct parameters:
power
0.3373199
NLS step: RSS = 0
model parameters:-0.799941 8.99983 -0.522623 212.314
iterations: 2
Convergence:
params varStruct
1.172208 1.000000
### After about 10 min hit cntrl C
Warning messages:
1: In log(xmid) : NaNs produced
2: In log(xmid) : NaNs produced
3: In log(xmid) : NaNs produced
#####################################################################
#### Data
tdat<-data.frame(conc=c(0.00203,0.0061,0.0183,0.0549,0.165,0.494,1.48,4.44,13.3,40),
response=c(12,-4,19,11,-5,-3,1,6,0,-8))
#### Self start function
SSllogis <- selfStart(~ A + (B-A)/(1 + exp(scal*(log(x)-log(xmid)))),
function(mCall, data, LHS)
{
#browser()
xy <- sortedXyData(mCall[["x"]], LHS, data)
if (nrow(xy) < 5) {
stop("too few distinct input values to fit a four-parameter
logistic")
}
rng <- range(xy$y)
drng <- diff(rng)
B <- rng[2] + 0.001
A <- rng[1] - 0.001
xy$prop <- log((B-xy$y)/(xy$y-A+0.001))
#(xy$y - rng[1] + 0.05 * drng)/(1.1 * drng)
ir <- as.vector(coef(lm(prop ~ log(x), data = xy)))
scal <- ir[2]
xmid <- exp(-ir[1]/ir[2])
#pars <- as.vector(coef(nls(y ~ cbind(1, 1/(1 + exp((xmid -
# x)/exp(lscal)))), data = xy, start = list(xmid = ir[1],
# lscal = log(abs(ir[2]))), algorithm = "plinear")))
value <- c(A, B, xmid, scal)
names(value) <- mCall[c("A", "B", "xmid", "scal")]
value
}, c("A", "B", "xmid", "scal"))
______________________________________________
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
