2015-09-02 5 views
0

У меня есть кадр данных, как этотОшибка в NLS единственном числе градиента R

time  WT WTIC ATIC RHIC esaIC   k uIC 
1 0.00 25.191 25.191 21.4 67.925 25.49 0.06152572 3.53 
2 0.05 25.186 25.191 21.4 67.925 25.49 0.06152572 3.53 
3 0.10 25.179 25.191 21.4 67.925 25.49 0.06152572 3.53 
4 0.15 25.168 25.191 21.4 67.925 25.49 0.06152572 3.53 
5 0.20 25.158 25.191 21.4 67.925 25.49 0.06152572 3.53 
6 0.25 25.147 25.191 21.4 67.925 25.49 0.06152572 3.53 

который я хотел бы, чтобы соответствовать, используя эту нелинейную функцию

f <- function(x,a1,a2,a3,a4,a5,a6,par1,par2,par3) { 
     Tinf <- a2 - (par2*(1-a3/100)*a4)/(1+par2*a5*a4) 
     kC <-par1*sqrt(a6) 
     V <- kC + par3 
     tau <- 1/(V*(1+par2*a5*a4)) 

     func <- a1 -(a1-Tinf)*(1-exp(-x/tau)) 
     return(func) 
} 

Однако при использовании nls

nls(WT~f(time,WTIC,ATIC,RHIC,esaIC,k,uIC,par1,par2,par3), data=df, start=c(par1=1, par2=1,par3=1)) 

Я получаю эту ошибку

Error in nlsModel(formula, mf, start, wts) : 
    singular gradient matrix at initial parameter estimates 

Я попытался изменить начальные значения параметров, но по-прежнему получаю ту же ошибку. Любая помощь?

+0

Сколько данных вы пытаетесь установить? Возможно, вы переуплотываете свою модель. В этом случае вы должны получить больше данных или указать меньше параметров. Если проблема связана только с начальными значениями, вы можете использовать поиск грубой силы, чтобы найти их. Вы видели эту тему? http://r.789695.n4.nabble.com/Error-quot-singular-gradient-matrix-at-initial-parameter-estimates-quot-in-nls-td1745042.html –

+0

У меня возникли проблемы с воспроизведением вашей проблемы , Является ли переменная «a» в знаменателе «tau» опечаткой или она определена в другом месте вашего кода? Аналогично для переменной «r» в kC. Кроме того, kC, по-видимому, определяется в пределах «f», но тогда он не используется нигде в вычислении возвращаемого значения «func». – frank2165

+0

@ RomanLuštrik в каждом наборе данных есть 300 баллов. – user3036416

ответ

2

Эта ошибка «сингулярного градиента» является недостатком функции «nls». "Оптим" ведет себя лучше:

df <- read.table(text = 
' time  WT WTIC ATIC RHIC esaIC   k uIC 
    "0.00" "25.191" "25.191" "21.4" "67.925" "25.49" "0.06152572" "3.53" 
    "0.05" "25.186" "25.191" "21.4" "67.925" "25.49" "0.06152572" "3.53" 
    "0.10" "25.179" "25.191" "21.4" "67.925" "25.49" "0.06152572" "3.53" 
    "0.15" "25.168" "25.191" "21.4" "67.925" "25.49" "0.06152572" "3.53" 
    "0.20" "25.158" "25.191" "21.4" "67.925" "25.49" "0.06152572" "3.53" 
    "0.25" "25.147" "25.191" "21.4" "67.925" "25.49" "0.06152572" "3.53"', 
    header = TRUE) 


f <- function(x,a1,a2,a3,a4,a5,a6,par1,par2,par3) { 
    Tinf <- a2 - (par2*(1-a3/100)*a4)/(1+par2*a5*a4) 
    kC <-par1*sqrt(a6) 
    V <- kC + par3 
    tau <- 1/(V*(1+par2*a5*a4)) 

    func <- a1 -(a1-Tinf)*(1-exp(-x/tau)) 
    return(func) 
} 

#---------------------------------------------------------------------- 
# Essentially the same as f: 

g <- function(v){f(v[1],v[2],v[3],v[4],v[5],v[6],v[7],v[8],v[9],v[10])} 

#---------------------------------------------------------------------- 
# The function we want to minimize: 

squaredError <- function(par) 
{ 
    sum((df$"WT"-apply(cbind(df[,-2],par[1],par[2],par[3]),1,g))^2) 
} 

#---------------------------------------------------------------------- 
# Optimization of the parameters: 

opt <- optim(par = c(1,1,1), 
       fn  = squaredError, 
       method = "BFGS") 

#---------------------------------------------------------------- 
# Result: 

opt 

squaredError(opt$par + c(1, 0, 0)*1e-3) 
squaredError(opt$par + c(-1, 0, 0)*1e-3) 
squaredError(opt$par + c(0, 1, 0)*1e-3) 
squaredError(opt$par + c(0,-1, 0)*1e-3) 
squaredError(opt$par + c(0, 0, 1)*1e-3) 
squaredError(opt$par + c(0, 0,-1)*1e-3) 

.

> opt 
$par 
[1] -0.04261273 -0.23600921 0.44504195 

$value 
[1] 4.572781e-05 

$counts 
function gradient 
    137  100 

$convergence 
[1] 1 

$message 
NULL 


> squaredError(opt$par + c(1, 0, 0)*1e-3) 
[1] 4.581051e-05 

> squaredError(opt$par + c(-1, 0, 0)*1e-3) 
[1] 4.583096e-05 

> squaredError(opt$par + c(0, 1, 0)*1e-3) 
[1] 4.900303e-05 

> squaredError(opt$par + c(0,-1, 0)*1e-3) 
[1] 4.939846e-05 

> squaredError(opt$par + c(0, 0, 1)*1e-3) 
[1] 4.57487e-05 

> squaredError(opt$par + c(0, 0,-1)*1e-3) 
[1] 4.575957e-05 
> 
Смежные вопросы