我有一个包含自变量和一组因变量的数据集.我想使用自举非线性最小二乘过程为每组自变量拟合一个函数.在某些情况下,自变量是“良好质量”,即合理地拟合函数.在其他情况下,他们很吵.
在所有情况下,我都可以使用nls()来估算参数.但是,当数据有噪声时,引导程序会在nls(…)中抛出错误Error：奇异梯度.我能理解为什么适合噪声数据的nls会失败,例如通过在迭代次数过多后无法收敛,但我不明白为什么它是一个奇异的梯度误差,为什么我只得到它重新采样的质量差的数据集.
码：
require(ggplot2)
require(plyr)
require(boot)
# Data are in long form: columns are 'enzyme','x',and 'y'
enz <- read.table("http://dl.dropBox.com/s/ts3ruh91kpr47sj/SE.txt",header=TRUE)
# Nonlinear formula to fit to data
mmFormula <- formula(y ~ (x*Vmax) / (x + Km))
nls完全能够拟合数据(即使在某些情况下,就像a,我怀疑模型是否适合数据.
# Use nls to fit mmFormula to the data - this works well enough fitDf <- ddply(enz,.(enzyme),function(x) coefficients(nls(mmFormula,x,start=list(Km=100,Vmax=0.5)))) # Create points to plot for the simulated fits xGrid <- 0:200 simFits <- dlply(fitDf,function(x) data.frame(x=xGrid,y=(xGrid * x$Vmax)/(xGrid + x$Km))) simFits <- ldply(simFits,identity) ggplot() + geom_point(data=enz,aes(x=x,y=y)) + geom_line(data=simFits,y=y)) + facet_wrap(~enzyme,scales="free_y") + aes(ymin=0)
Bootstrapping适用于高质量数据：
# Function to pass to bootstrap; returns coefficients of nls fit to formula
nlsCoef <- function(df,i) {
KmGuess <- median(df$x)
VmaxGuess <- max(df$y)
dfSamp <- df[i,]
nlsCoef <- coefficients(nls(mmFormula,dfSamp,Vmax=0.5)))
}
eBoot <- boot(subset(enz,enzyme=="e"),nlsCoef,R=1000) #No error
但不是因为质量差的数据
dBoot <- boot(subset(enz,enzyme=="d"),R=10) > Error in nls(mmFormula,start = list(Km = KmGuess,Vmax = VmaxGuess)) : singular gradient
导致此错误的原因是什么？考虑到我想使用plyr同时执行大量的bootstrap模拟,我应该怎么做呢？
解决方法
这允许您检查发生的情况：
#modified function
#returns NAs if fit is not sucessfull
#does global assignment to store bootstrap permutations
nlsCoef <- function(df,]
fit <- NULL
try(fit <- nls(mmFormula,Vmax=0.5)))
if(!is.null(fit)){
res <- coef(fit)
} else{
res <- c(Km=NA,Vmax=NA)
}
istore[k,] <<- i
k <<- k+1
res
}
n <- 100
istore <- matrix(nrow=n+1,ncol=9)
k <- 1
dat <- subset(enz,enzyme=="d")
dBoot <- boot(dat,R=n)
#permutations that create samples that cannot be fitted
nais <- istore[-1,][is.na(dBoot$t[,1]),]
#look at first bootstrap sample
#that could not be fitted
samp <- dat[nais[1,],]
plot(y~x,data=samp)
fit <- nls(mmFormula,samp,Vmax=0.5))
#error
您还可以使用自启动模型：
try(fit <- nls(y ~ SSmicmen(x,Vmax,Km),data = dfSamp))
too few distinct input values to fit a Michaelis-Menten model
这意味着,一些自助样本含有少于三个不同的浓度.
但也有一些其他错误：
step factor 0.000488281 reduced below 'minFactor' of 0.000976562
您可以通过减少minFactor来避免.
以下是令人讨厌的.你可以尝试不同的拟合算法或起始值：
singular gradient matrix at initial parameter estimates singular gradient