Hi to all OK as I did not get any response and I really need some insight I try again with different subject line I have troubles with correct evaluating/structure of nls input Here is an example # data x <-1:10 y <-1/(.5-x)+rnorm(10)/100 # formula list form <- structure(list(a = list(quote(y ~ 1/(a - x)), "list(a=mean(y))")), .Names = "a") # This gives me an error due to not suitable default starting value fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y)) # This works and gives me a result fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y), start=list(a=mean(y))) *** How to organise list "form" and call to nls to enable to use other then default starting values***. I thought about something like fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y), start=get(form [[1]] [[2]])) ^^^^^^^^^^^^^^^^^^^ but this gives me an error so it is not correct syntax. (BTW I tried eval, assign, sustitute, evalq and maybe some other options but did not get it right. I know I can put starting values interactively but what if I want them computed by some easy way which is specified by second part of a list, like in above example. If it matters WXP, R2.9.0 devel. Regards Petr petr.pikal at precheza.cz
Petr PIKAL wrote:> Hi to all > > OK as I did not get any response and I really need some insight I try > again with different subject line > > I have troubles with correct evaluating/structure of nls input > > Here is an example > > # data > x <-1:10 > y <-1/(.5-x)+rnorm(10)/100 > > # formula list > form <- structure(list(a = list(quote(y ~ 1/(a - x)), "list(a=mean(y))")), > .Names = "a") > > # This gives me an error due to not suitable default starting value > > fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y)) > > # This works and gives me a result > > fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y), start=list(a=mean(y))) > > *** How to organise list "form" and call to nls to enable to use other > then default starting values***. > > I thought about something like > > fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y), start=get(form [[1]] > [[2]])) > ^^^^^^^^^^^^^^^^^^^ > but this gives me an error so it is not correct syntax. (BTW I tried eval, > assign, sustitute, evalq and maybe some other options but did not get it > right. > > I know I can put starting values interactively but what if I want them > computed by some easy way which is specified by second part of a list, > like in above example.If you really want to orgnize it that way, why not simpler as in: form <- list(y ~ 1/(a - x), a = mean(y)) fit <- nls(form[[1]], data.frame(x=x, y=y), start = form[2]) Uwe Ligges> If it matters > WXP, R2.9.0 devel. > > Regards > Petr > > petr.pikal at precheza.cz > > ______________________________________________ > R-help at 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.
It is not very clear what you are trying to do here, and> form <- structure(list(a = list(quote(y ~ 1/(a - x)), "list(a=mean(y))")), > .Names = "a")is using a historic anomaly (see the help page). I am gussing you want to give nls an object containing a formula and an expression for the starting value. It seems you are re-inventing self-starting nls models: see ?selfStart and MASS$ ca p. 216. One way to use them in your example is mod <- selfStart(~ 1/(a - x), function(mCall, data, LHS) { structure(mean(eval(LHS, data)), names="a") }, "a") nls(y ~ mod(x, a)) But if you want to follow ypur route, youer starting values would be better to be a list that you evaluate in an appropriate context (which y is this supposed to be?). nls() knows where it will find variables, but it is not so easy for you to replicate its logic without access to its evaluation frames. On Mon, 2 Mar 2009, Petr PIKAL wrote:> Hi to all > > OK as I did not get any response and I really need some insight I try > again with different subject line > > I have troubles with correct evaluating/structure of nls input > > Here is an example > > # data > x <-1:10 > y <-1/(.5-x)+rnorm(10)/100 > > # formula list > form <- structure(list(a = list(quote(y ~ 1/(a - x)), "list(a=mean(y))")), > .Names = "a") > > # This gives me an error due to not suitable default starting value > > fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y)) > > # This works and gives me a result > > fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y), start=list(a=mean(y))) > > *** How to organise list "form" and call to nls to enable to use other > then default starting values***. > > I thought about something like > > fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y), start=get(form [[1]] > [[2]])) > ^^^^^^^^^^^^^^^^^^^ > but this gives me an error so it is not correct syntax. (BTW I tried eval, > assign, sustitute, evalq and maybe some other options but did not get it > right. > > I know I can put starting values interactively but what if I want them > computed by some easy way which is specified by second part of a list, > like in above example. > > If it matters > WXP, R2.9.0 devel. > > Regards > Petr > > petr.pikal at precheza.cz-- Brian D. Ripley, ripley at stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UK Fax: +44 1865 272595
Thank you It was simplified version of my problem. I want to elaborate a function which can take predefined list of formulas, some data and evaluate which formulas can fit the data. I was inspired by some article in Chemical engineering in which some guy used excel solver for such task. I was curious if I can do it in R too. I am not sure if nls is appropriate tool for such task but I had to start somewhere. Here is a function which takes list of formulas and data and gives a result for each formula. modely <- function(formula, data, ...){ ll <- length(formula) #no of items in formula list result2 <- vector("list", ll) #prepare results result1 <- rep(NA, ll) for(i in 1:ll) { fit<-try(nls(formula[[i]], data)) if( class(fit)=="try-error") result1[i] <- NA else result1[i] <- sum(resid(fit)^2) if( class(fit)=="try-error") result2[[i]] <- NA else result2[[i]] <- coef(fit) } ooo<-order(result1) #order results according to residual sum #combine results into one list together with functions used result <- mapply(c, "sq.resid" = result1, result2) names(result) <- as.character(formula) # output result[ooo] } # data x <-1:10 y <-1/(.5-x)+rnorm(10)/100 # list of formulas fol <- structure(list(a = y ~ 1/(a - x), b = y ~ a * x^2 + b * log(x), c = y ~ x^a), .Names = c("a", "b", "c")) modely(fol, data.frame(x=x, y=y) does not use "correct" model because when using default start values it results in> nls(fol[[1]], data.frame(x=x, y=y))Error in numericDeriv(form[[3]], names(ind), env) : Missing value or an infinity produced when evaluating the model however nls(fol[[1]], data.frame(x=x, y=y), start=list(a=mean(y))) gives correct result. Therefore I started think about how to add a "better" starting value for some fits as a second part of my formula list to define structure like> list(a= formula1, start.formula1, b=formula2, start.formula2, ....) I wonder If you can push me to better direction. Thanks again Best regards Petr Uwe Ligges <ligges at statistik.tu-dortmund.de> napsal dne 02.03.2009 09:41:45:> Petr PIKAL wrote: > > Hi to all > > > > OK as I did not get any response and I really need some insight I try > > again with different subject line > > > > I have troubles with correct evaluating/structure of nls input > > > > Here is an example > > > > # data > > x <-1:10 > > y <-1/(.5-x)+rnorm(10)/100 > > > > # formula list > > form <- structure(list(a = list(quote(y ~ 1/(a - x)),"list(a=mean(y))")),> > .Names = "a") > > > > # This gives me an error due to not suitable default starting value > > > > fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y)) > > > > # This works and gives me a result > > > > fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y),start=list(a=mean(y)))> > > > *** How to organise list "form" and call to nls to enable to use other> > then default starting values***. > > > > I thought about something like > > > > fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y), start=get(form[[1]]> > [[2]])) > > ^^^^^^^^^^^^^^^^^^^ > > but this gives me an error so it is not correct syntax. (BTW I triedeval,> > assign, sustitute, evalq and maybe some other options but did not getit> > right. > > > > I know I can put starting values interactively but what if I want them> > computed by some easy way which is specified by second part of a list,> > like in above example. > > If you really want to orgnize it that way, why not simpler as in: > > form <- list(y ~ 1/(a - x), a = mean(y)) > fit <- nls(form[[1]], data.frame(x=x, y=y), start = form[2]) > > > Uwe Ligges > > > > If it matters > > WXP, R2.9.0 devel. > > > > Regards > > Petr > > > > petr.pikal at precheza.cz > > > > ______________________________________________ > > R-help at r-project.org mailing list > > https://stat.ethz.ch/mailman/listinfo/r-help > > PLEASE do read the posting guidehttp://www.R-project.org/posting-guide.html> > and provide commented, minimal, self-contained, reproducible code.
Petr PIKAL wrote:> Thank you > > It was simplified version of my problem. I want to elaborate a function > which can take predefined list of formulas, some data and evaluate which > formulas can fit the data. I was inspired by some article in Chemical > engineering in which some guy used excel solver for such task. I was > curious if I can do it in R too. I am not sure if nls is appropriate tool > for such task but I had to start somewhere. > > Here is a function which takes list of formulas and data and gives a > result for each formula. > > modely <- function(formula, data, ...){ > ll <- length(formula) #no of items in formula list > result2 <- vector("list", ll) #prepare results > result1 <- rep(NA, ll) > for(i in 1:ll) { > fit<-try(nls(formula[[i]], data)) > if( class(fit)=="try-error") result1[i] <- NA else result1[i] <- > sum(resid(fit)^2) > if( class(fit)=="try-error") result2[[i]] <- NA else result2[[i]] <- > coef(fit) > } > > ooo<-order(result1) #order results according to residual sum > > #combine results into one list together with functions used > > result <- mapply(c, "sq.resid" = result1, result2) > names(result) <- as.character(formula) > # output > result[ooo] > } > > # data > x <-1:10 > y <-1/(.5-x)+rnorm(10)/100 > > # list of formulas > fol <- structure(list(a = y ~ 1/(a - x), b = y ~ a * x^2 + b * log(x), > c = y ~ x^a), .Names = c("a", "b", "c")) > > modely(fol, data.frame(x=x, y=y) > > does not use "correct" model because when using default start values it > results in > >> nls(fol[[1]], data.frame(x=x, y=y)) > Error in numericDeriv(form[[3]], names(ind), env) : > Missing value or an infinity produced when evaluating the model > > however > > nls(fol[[1]], data.frame(x=x, y=y), start=list(a=mean(y))) > > gives correct result. Therefore I started think about how to add a > "better" starting value for some fits as a second part of my formula list > to define structure like> > > list(a= formula1, start.formula1, b=formula2, start.formula2, ....) > > I wonder If you can push me to better direction.You can make up a list of lists (each containing one formula and its starting values) or specify formulas in one list and starting values in a corresponding second list. You need just the corresponding subsetting in your call to nls such as in the simple case I suggested already. Best, Uwe> Thanks again > Best regards > Petr > > > > > Uwe Ligges <ligges at statistik.tu-dortmund.de> napsal dne 02.03.2009 > 09:41:45: > >> Petr PIKAL wrote: >>> Hi to all >>> >>> OK as I did not get any response and I really need some insight I try >>> again with different subject line >>> >>> I have troubles with correct evaluating/structure of nls input >>> >>> Here is an example >>> >>> # data >>> x <-1:10 >>> y <-1/(.5-x)+rnorm(10)/100 >>> >>> # formula list >>> form <- structure(list(a = list(quote(y ~ 1/(a - x)), > "list(a=mean(y))")), >>> .Names = "a") >>> >>> # This gives me an error due to not suitable default starting value >>> >>> fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y)) >>> >>> # This works and gives me a result >>> >>> fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y), > start=list(a=mean(y))) >>> *** How to organise list "form" and call to nls to enable to use other > >>> then default starting values***. >>> >>> I thought about something like >>> >>> fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y), start=get(form > [[1]] >>> [[2]])) >>> ^^^^^^^^^^^^^^^^^^^ >>> but this gives me an error so it is not correct syntax. (BTW I tried > eval, >>> assign, sustitute, evalq and maybe some other options but did not get > it >>> right. >>> >>> I know I can put starting values interactively but what if I want them > >>> computed by some easy way which is specified by second part of a list, > >>> like in above example. >> If you really want to orgnize it that way, why not simpler as in: >> >> form <- list(y ~ 1/(a - x), a = mean(y)) >> fit <- nls(form[[1]], data.frame(x=x, y=y), start = form[2]) >> >> >> Uwe Ligges >> >> >>> If it matters >>> WXP, R2.9.0 devel. >>> >>> Regards >>> Petr >>> >>> petr.pikal at precheza.cz >>> >>> ______________________________________________ >>> R-help at 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. >
Hi Prof Brian Ripley <ripley at stats.ox.ac.uk> napsal dne 02.03.2009 10:24:52:> It is not very clear what you are trying to do here, and > > > form <- structure(list(a = list(quote(y ~ 1/(a - x)),"list(a=mean(y))")),> > .Names = "a") > > is using a historic anomaly (see the help page). > > I am gussing you want to give nls an object containing a formula and > an expression for the starting value. It seems you are re-inventingYou are correct as usually.> self-starting nls models: see ?selfStart and MASS$ ca p. 216. > One way to use them in your example is > > mod <- selfStart(~ 1/(a - x), function(mCall, data, LHS) { > structure(mean(eval(LHS, data)), names="a") > }, "a") > > nls(y ~ mod(x, a)) > > But if you want to follow ypur route, youer starting values would be > better to be a list that you evaluate in an appropriate context > (which y is this supposed to be?). nls() knows where it will find > variables, but it is not so easy for you to replicate its logic > without access to its evaluation frames.It was simplified version of my problem. I want to elaborate a function which can take predefined list of formulas, some data and evaluate which formulas can fit the data. I was inspired by some article in Chemical engineering in which some guy used excel solver for such task. I was curious if I can do it in R too. I am not sure if nls is appropriate tool for such task but I had to start somewhere. Here is a function which takes list of formulas and data and gives a result for each formula. modely <- function(formula, data, ...){ ll <- length(formula) #no of items in formula list result2 <- vector("list", ll) #prepare results result1 <- rep(NA, ll) for(i in 1:ll) { fit<-try(nls(formula[[i]], data)) if( class(fit)=="try-error") result1[i] <- NA else result1[i] <- sum(resid(fit)^2) if( class(fit)=="try-error") result2[[i]] <- NA else result2[[i]] <- coef(fit) } ooo<-order(result1) #order results according to residual sum #combine results into one list together with functions used result <- mapply(c, "sq.resid" = result1, result2) names(result) <- as.character(formula) # output result[ooo] } # data x <-1:10 y <-1/(.5-x)+rnorm(10)/100 # list of formulas fol <- structure(list(a = y ~ 1/(a - x), b = y ~ a * x^2 + b * log(x), c = y ~ x^a), .Names = c("a", "b", "c")) modely(fol, data.frame(x=x, y=y) does not use "correct" model because when using default start values it results in> nls(fol[[1]], data.frame(x=x, y=y))Error in numericDeriv(form[[3]], names(ind), env) : Missing value or an infinity produced when evaluating the model I tried to establish such structure to get more appropriate starting values list(a= list(formula1, start.formula1), b=list(formula2, start.formula2), ....) But did not manage yet to get correct syntax for let say mean of response values. I try to look more closely what I can achieve with selfStart Thank you again Best regards Petr> > On Mon, 2 Mar 2009, Petr PIKAL wrote: > > > Hi to all > > > > OK as I did not get any response and I really need some insight I try > > again with different subject line > > > > I have troubles with correct evaluating/structure of nls input > > > > Here is an example > > > > # data > > x <-1:10 > > y <-1/(.5-x)+rnorm(10)/100 > > > > # formula list > > form <- structure(list(a = list(quote(y ~ 1/(a - x)),"list(a=mean(y))")),> > .Names = "a") > > > > # This gives me an error due to not suitable default starting value > > > > fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y)) > > > > # This works and gives me a result > > > > fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y),start=list(a=mean(y)))> > > > *** How to organise list "form" and call to nls to enable to use other > > then default starting values***. > > > > I thought about something like > > > > fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y), start=get(form[[1]]> > [[2]])) > > ^^^^^^^^^^^^^^^^^^^ > > but this gives me an error so it is not correct syntax. (BTW I triedeval,> > assign, sustitute, evalq and maybe some other options but did not getit> > right. > > > > I know I can put starting values interactively but what if I want them > > computed by some easy way which is specified by second part of a list, > > like in above example. > > > > If it matters > > WXP, R2.9.0 devel. > > > > Regards > > Petr > > > > petr.pikal at precheza.cz > > > -- > Brian D. Ripley, ripley at stats.ox.ac.uk > Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ > University of Oxford, Tel: +44 1865 272861 (self) > 1 South Parks Road, +44 1865 272866 (PA) > Oxford OX1 3TG, UK Fax: +44 1865 272595
Petr PIKAL
2009-Mar-03 10:55 UTC
[R] select nls starting values from list was: nls does not accept start values
Hallo Here is my "final" solution, but i still wonder if it could not be improved somehow. fit <- modely(fol, data.frame(x=x, y=y1) fit or fit <- modely(fol, data.frame(x=x, y=y2) fit shows that it arrives to some results. However selection starting values as means of response variable seems to me suboptimal. Is it possible to have some list of functions which will for each item in formula list compute starting values for nls? Or is it this mean guess as good as any other? Petr petr.pikal at precheza.cz # formulas fol <- list(y~1/(a-x), y~a*x^2+b*log(x)) set.seed(111) x<-1:10 y1 <- 1/(0.5-x) +rnorm(length(x))/10 y2 <- 2*x^2 + 3*log(x) + rnorm(length(x)) # here is a function modely <- function(formula, data, ...){ ll <- length(formula) # no of items in formula list # initiation of vectors result2 <- vector("list", ll) result1 <- rep(NA, ll) for(i in 1:ll) { # some start values ***this shall be improved somehow*** n <- length(all.vars(fol[[i]])) - 2 start <- as.list(rep(mean(data[,2]),n)) names(start) <- letters[1:n] # fit fit<-try(nls(formula[[i]], data, start)) if( class(fit)=="try-error") result1[i] <- NA else result1[i] <- sum(resid(fit)^2) if( class(fit)=="try-error") result2[[i]] <- NA else result2[[i]] <- coef(fit) } # ordering and combining results ooo<-order(result1) result <- mapply(c, "sq.resid" = result1, result2) names(result) <- as.character(formula) result[ooo] } Uwe Ligges <ligges at statistik.tu-dortmund.de> napsal dne 02.03.2009 10:54:16:> > > Petr PIKAL wrote: > > Thank you > > > > It was simplified version of my problem. I want to elaborate afunction> > which can take predefined list of formulas, some data and evaluatewhich> > formulas can fit the data. I was inspired by some article in Chemical > > engineering in which some guy used excel solver for such task. I was > > curious if I can do it in R too. I am not sure if nls is appropriatetool> > for such task but I had to start somewhere. > > > > Here is a function which takes list of formulas and data and gives a > > result for each formula. > > > > modely <- function(formula, data, ...){ > > ll <- length(formula) #no of items in formula list > > result2 <- vector("list", ll) #prepare results > > result1 <- rep(NA, ll) > > for(i in 1:ll) { > > fit<-try(nls(formula[[i]], data)) > > if( class(fit)=="try-error") result1[i] <- NA else result1[i] <- > > sum(resid(fit)^2) > > if( class(fit)=="try-error") result2[[i]] <- NA else result2[[i]] <- > > coef(fit) > > } > > > > ooo<-order(result1) #order results according to residual sum > > > > #combine results into one list together with functions used > > > > result <- mapply(c, "sq.resid" = result1, result2) > > names(result) <- as.character(formula) > > # output > > result[ooo] > > } > > > > # data > > x <-1:10 > > y <-1/(.5-x)+rnorm(10)/100 > > > > # list of formulas > > fol <- structure(list(a = y ~ 1/(a - x), b = y ~ a * x^2 + b * log(x),> > c = y ~ x^a), .Names = c("a", "b", "c")) > > > > modely(fol, data.frame(x=x, y=y) > > > > does not use "correct" model because when using default start valuesit> > results in > > > >> nls(fol[[1]], data.frame(x=x, y=y)) > > Error in numericDeriv(form[[3]], names(ind), env) : > > Missing value or an infinity produced when evaluating the model > > > > however > > > > nls(fol[[1]], data.frame(x=x, y=y), start=list(a=mean(y))) > > > > gives correct result. Therefore I started think about how to add a > > "better" starting value for some fits as a second part of my formulalist> > to define structure like> > > > > list(a= formula1, start.formula1, b=formula2, start.formula2, ....) > > > > I wonder If you can push me to better direction. > > > You can make up a list of lists (each containing one formula and its > starting values) or specify formulas in one list and starting values in > a corresponding second list. > You need just the corresponding subsetting in your call to nls such as > in the simple case I suggested already. > > Best, > Uwe > > > > > Thanks again > > Best regards > > Petr > > > > > > > > > > Uwe Ligges <ligges at statistik.tu-dortmund.de> napsal dne 02.03.2009 > > 09:41:45: > > > >> Petr PIKAL wrote: > >>> Hi to all > >>> > >>> OK as I did not get any response and I really need some insight Itry> >>> again with different subject line > >>> > >>> I have troubles with correct evaluating/structure of nls input > >>> > >>> Here is an example > >>> > >>> # data > >>> x <-1:10 > >>> y <-1/(.5-x)+rnorm(10)/100 > >>> > >>> # formula list > >>> form <- structure(list(a = list(quote(y ~ 1/(a - x)), > > "list(a=mean(y))")), > >>> .Names = "a") > >>> > >>> # This gives me an error due to not suitable default starting value > >>> > >>> fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y)) > >>> > >>> # This works and gives me a result > >>> > >>> fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y), > > start=list(a=mean(y))) > >>> *** How to organise list "form" and call to nls to enable to useother> > > >>> then default starting values***. > >>> > >>> I thought about something like > >>> > >>> fit <- nls(form [[1]] [[1]], data.frame(x=x, y=y), start=get(form > > [[1]] > >>> [[2]])) > >>> ^^^^^^^^^^^^^^^^^^^ > >>> but this gives me an error so it is not correct syntax. (BTW I tried> > eval, > >>> assign, sustitute, evalq and maybe some other options but did notget> > it > >>> right. > >>> > >>> I know I can put starting values interactively but what if I wantthem> > > >>> computed by some easy way which is specified by second part of alist,> > > >>> like in above example. > >> If you really want to orgnize it that way, why not simpler as in: > >> > >> form <- list(y ~ 1/(a - x), a = mean(y)) > >> fit <- nls(form[[1]], data.frame(x=x, y=y), start = form[2]) > >> > >> > >> Uwe Ligges > >> > >> > >>> If it matters > >>> WXP, R2.9.0 devel. > >>> > >>> Regards > >>> Petr > >>> > >>> petr.pikal at precheza.cz > >>> > >>> ______________________________________________ > >>> R-help at 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. > >