Hi there,
I don't have much experience with fitting at all and I'd like to get
some advice how to use the "weights"-argument with nls correctly.
I have created some data with a sigmoidal curve shape. Each y-Value was
generated by the mean of three values. A standard deviation was
calculated too.
Now, I'd like to weight the data points respective to its standard
deviation. If this is very high, the datapoint should be considered less
important than a datapoint with a low standard deviation.
I've tested it with the following code but I'd like to know whether I
should change anything applying the weights. I've only guessed what
could give a nice result...
Further, I'd like to know, why I get an error message, if I change the
bottom-parameter to min(yData)??? Can anybody explain, why it ends up with:
"singular gradient"
8< ------------------------------------------------------- >8
sigmoid <- function(x, bottom, top, slope, logec50) { bottom + (top -
bottom) / (1 + 10^((logec50-x)*slope)) }
xData <- seq(-10,10,1.2)
yData1 <- sigmoid(xData,10,1,-0.2,3) + rnorm(1:length(xData), sd = 0.5)
yData2 <- sigmoid(xData,10,1,-0.2,3) + rnorm(1:length(xData), sd = 0.2)
yData3 <- sigmoid(xData,10,1,-0.2,3) + rnorm(1:length(xData), sd = 1.1)
yMat <- rbind(yData1, yData2, yData3)
yData <- apply(yMat, 2, mean)
yDataSd <- apply(yMat, 2, sd)
plot(yData ~ xData, ylim = c(min(yData - yDataSd), max(yData + yDataSd)))
arrows(xData, yData - yDataSd, xData, yData + yDataSd, angle = 90,
length = 0.1, code = 3)
# without weights
model <- nls(yData ~ sigmoid(xData,bottom, top, slope, logec50), start =
list(bottom = yData[1], top = max(yData), slope = 0.5, logec50 =
median(xData)), trace = TRUE)
lines(xData, predict(model), col ="blue")
# with weights
model <- nls(yData ~ sigmoid(xData,bottom, top, slope, logec50), start =
list(bottom = yData[1], top = max(yData), slope = 0.5, logec50 =
median(xData)), weights = 1 / (yDataSd^2), trace = TRUE)
lines(xData, predict(model), col ="green")
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.