2014-12-19 7 views
4

Я хотел бы подгонять кривые на ряде наборов данных, сгруппированных по обработке. Это очень хорошо работает с nlslist, но теперь я хотел бы ввести верхние оценки моих параметров.R Как использовать границы и алгоритм «port» в nlslist?

Представление границ очень хорошо работает, когда я поместил каждую группу отдельно с помощью nls, но, по-видимому, нет, когда я хочу ускорить работу (у меня есть много других наборов данных) с nlslist.

Может ли кто-нибудь помочь мне здесь с какой-либо идеей, как решить эту проблему?

Пример моего набора данных:

DF1<-data.frame(treatment = rep(c("mineral","residues"),4), 
      N_level = c(0,0,100,100,200,200,300,300), 
      yield = c(8,8.5,10,10.5,11,9.8,9.5,9.7)) 

DF1

treatment N_level yield 
1 mineral  0 8.0 
2 residues  0 8.5 
3 mineral  100 10.0 
4 residues  100 10.5 
5 mineral  200 11.0 
6 residues  200 9.8 
7 mineral  300 9.5 
8 residues  300 9.7 

Попытка соответствовать этот набор данных только с НЛС работает хорошо:

fit_mineral <- nls(formula = yield ~ a + b*0.99^N_level +c*N_level, 
       data=subset(DF1, subset = treatment == "mineral"), 
       algorithm = "port", start = list(a = 12, b = -8, c= -0.01), 
       upper = list(a=1000, b=-0.000001, c=-0.000001)) 

fit_mineral

Nonlinear regression model 
model: yield ~ a + b * 0.99^N_level + c * N_level 
data: subset(DF1, subset = treatment == "mineral") 
    a  b  c 
13.7882 -5.8685 -0.0126 
residual sum-of-squares: 0.4679 

Но как только я стараюсь сочетать вещи в nlslist, он не работает:

fit_mineral_and_residues <- nlsList(model = yield ~ a + b*0.99^N_level +c*N_level 
         | treatment, data=DF1, 
         algorithm = "port", start = list(a = 12, b = -8, c= -0.01), 
         upper = list(a=1000, b=-0.000001, c=-0.000001)) 

сообщение об ошибке:

Error in nlsList(model = yield ~ a + b * 0.99^N_level + c * N_level | : 
unused arguments (algorithm = "port", upper = list(a = 1000, b = -1e-06, c = -1e-06)) 
+0

'nlme :: nlsList' не поддерживает разные алгоритмы. – Roland

ответ

2

Я просто столкнулся с той же проблемой - это проблема, которую я считаю, действительно должно быть исправлено на уровне исходного кода!

В качестве временного решения можно, возможно, попробовать построить nlsList объект самостоятельно, то вдоль линий

library(nlme) 
DF1=data.frame(treatment = rep(c("mineral","residues"),4), 
       N_level = c(0,0,100,100,200,200,300,300), 
       yield = c(8,8.5,10,10.5,11,9.8,9.5,9.7)) 
nlslist=lapply(unique(DF1$treatment),function(i) {datasubs=DF1[DF1$treatment==i,]; 
                  nls(yield ~ a + b*0.99^N_level +c*N_level, 
                   data=datasubs, 
                   start = list(a = 12, b = -8, c= -0.01), 
                   upper = list(a=1000, b=-0.000001, c=-0.000001), 
                   algorithm="port", 
                   control=list(maxit=100000,tol=1e-10,warnOnly=T,minFactor=1e-10)) 
                }) 
names(nlslist)=unique(DF1$treatment) 
attr(nlslist, "dims")=list(N = nrow(DF1), M = length(nlslist)) 
attr(nlslist, "call")=NA # this line is not correct - should be fixed 
attr(nlslist,"groups")=names(nlslist) 
attr(nlslist, "origOrder")=1:length(unique(DF1$treatment)) 
attr(nlslist, "pool")=TRUE 
attr(nlslist, "groupsForm")=formula(~treatment) 
class(nlslist)=c("nlsList", "lmList") 

это почти заставляет меня там, за исключением того, что я не знаю, как правильно заполнить в "call" слота (в nlsList построен с использованием match.call() - любой, кто знает, как сделать это случайно

Если вы хотите, чтобы проверить правильную структуру, это может быть сделано, например, глядя на

test=nlsList(uptake ~ SSasympOff(conc, Asym, lrc, c0), 
       data = CO2, start = c(Asym = 30, lrc = -4.5, c0 = 52)) 
class(test)=NULL 
test 
+0

Привет, Том, Извините, что вы оказались в той же ситуации, но приятно слышать, что я не одинок в этом. Было бы хорошо, если бы это можно было решить на уровне исходного кода. В то же время: вы нашли работу? –

+0

Привет, Том, я скорректировал ваш код, чтобы использовать его в виде цикла и хранить коэффициенты в кадре данных (так что без проблем с выполнением функции nlslist - это было бы более приятным решением, я думаю, но не могу понять это) ,На данный момент работает достаточно хорошо для меня, хотя и не является долгосрочным решением. Спасибо за вашу помощь! Я могу проверить ваш ответ как решение, но, возможно, вы хотите приспособиться к работоспособному? –

+0

Hi Renske - Glad Я мог бы немного помочь - ну, может быть, оставил мой пост так, как есть, и немного подожду, чтобы принять ответ, возможно, кто-то рано или поздно опубликует решения, чтобы показать, как правильно построить рабочий объект 'nlsList', или исправить исходный код 'nlsList', не знаю ... –

0

Регулировка предложение Тома для создания обходного пути, я сейчас, используя следующий код:

DF2 <- lapply(unique(DF1$treatment),function(i) {datasubs=DF1[DF1$treatment==i,]; 
               coef(nls(yield ~ a + b*0.99^N_level +c*N_level, 
                data=datasubs, 
                start = list(a = 12, b = -8, c= -0.01), 
                upper = list(a=1000, b=-0.000001, c=-0.000001), 
                algorithm="port", 
                control=list(maxit=100000,tol=1e-10,warnOnly=T,minFactor=1e-10))) 
}) 

DF2 <- data.frame(DF2) 
names(DF2) <- levels(DF1$treatment) 
DF2 <- t(DF2) 
DF2 

Который дает мне коэффициенты припадков:

   a   b   c 
mineral 13.78825 -5.868506 -0.01260393 
residues 12.76157 -4.201832 -0.01006236 

работает достаточно хорошо для меня на момент, также для больших наборов данных.

 Смежные вопросы

  • Нет связанных вопросов^_^