Hello Gabor, It is very kind of you to reply and give suggestion so rapid. I will try to learn and use it.
Thanks very much for your help! Best regards, Jianling On 22 September 2015 at 06:45, Gabor Grothendieck <ggrothendi...@gmail.com> wrote: > Or if you really can't bear to write out 20 terms have R do it for you: > > # number of terms is the number of unique values in ref column > nterms <- length(unique(dproot$ref)) > > dproot2 <- do.call(data.frame, transform(dproot, ref = outer(dproot$ref, > seq(nterms), "==") + 0)) > > # construct the formula as a string > terms <- paste( sprintf("Rm%d*ref.%d", 1:nterms, 1:nterms), collapse = "+") > fo <- sprintf("den ~ (%s)/(1+(depth/d50)^c)", terms) > > library(nlmrt) > fm <- nlxb(fo, data = dproot2, masked = "Rm6", > start = c(Rm1=1.01, Rm2=1.01, Rm3=1.01, Rm4=6.65, Rm5=1.01, Rm6=1, > d50=20, c=-1)) > > > On Tue, Sep 22, 2015 at 7:04 AM, Gabor Grothendieck > <ggrothendi...@gmail.com> wrote: >> >> Just write out the 20 terms. >> >> On Mon, Sep 21, 2015 at 10:26 PM, Jianling Fan <fanjianl...@gmail.com> >> wrote: >>> >>> Hello, Gabor, >>> >>> Thanks again for your suggestion. And now I am trying to improve the >>> code by adding a function to replace the express "Rm1 * ref.1 + Rm2 * >>> ref.2 + Rm3 * ref.3 + Rm4 * ref.4 + Rm5 * ref.5 + Rm6 * ref.6" because >>> I have some other dataset need to fitted to the same model but with >>> more groups (>20). >>> >>> I tried to add the function as: >>> >>> denfun<-function(i){ >>> for(i in 1:6){ >>> Rm<-sum(Rm[i]*ref.i) >>> return(Rm)} >>> } >>> >>> but I got another error when I incorporate this function into my >>> regression: >>> >>> >fitdp1<-nlxb(den ~ denfun(6)/(1+(depth/d50)^c), >>> data = dproot2, >>> start = c(Rm1=1.01, Rm2=1.01, Rm3=1.01, Rm4=6.65, >>> Rm5=1.01, Rm6=1, d50=20, c=-1), >>> masked = "Rm6") >>> >>> Error in deriv.default(parse(text = resexp), names(start)) : >>> Function 'denfun' is not in the derivatives table >>> >>> I think there must be something wrong with my function. I tried some >>> times but am not sure how to improve it because I am quite new to R. >>> >>> Could anyone please give me some suggestion. >>> >>> Thanks a lot! >>> >>> >>> Jianling >>> >>> >>> On 22 September 2015 at 00:43, Gabor Grothendieck >>> <ggrothendi...@gmail.com> wrote: >>> > Express the formula in terms of simple operations like this: >>> > >>> > # add 0/1 columns ref.1, ref.2, ..., ref.6 >>> > dproot2 <- do.call(data.frame, transform(dproot, ref = >>> > outer(dproot$ref, >>> > seq(6), "==") + 0)) >>> > >>> > # now express the formula in terms of the new columns >>> > library(nlmrt) >>> > fitdp1<-nlxb(den ~ (Rm1 * ref.1 + Rm2 * ref.2 + Rm3 * ref.3 + Rm4 * >>> > ref.4 + >>> > Rm5 * ref.5 + Rm6 * ref.6)/(1+(depth/d50)^c), >>> > data = dproot2, >>> > start = c(Rm1=1.01, Rm2=1.01, Rm3=1.01, Rm4=6.65, Rm5=1.01, >>> > Rm6=1, >>> > d50=20, c=-1), >>> > masked = "Rm6") >>> > >>> > where we used this input: >>> > >>> > Lines <- " depth den ref >>> > 1 20 0.5730000 1 >>> > 2 40 0.7800000 1 >>> > 3 60 0.9470000 1 >>> > 4 80 0.9900000 1 >>> > 5 100 1.0000000 1 >>> > 6 10 0.6000000 2 >>> > 7 20 0.8200000 2 >>> > 8 30 0.9300000 2 >>> > 9 40 1.0000000 2 >>> > 10 20 0.4800000 3 >>> > 11 40 0.7340000 3 >>> > 12 60 0.9610000 3 >>> > 13 80 0.9980000 3 >>> > 14 100 1.0000000 3 >>> > 15 20 3.2083491 4 >>> > 16 40 4.9683383 4 >>> > 17 60 6.2381133 4 >>> > 18 80 6.5322348 4 >>> > 19 100 6.5780660 4 >>> > 20 120 6.6032064 4 >>> > 21 20 0.6140000 5 >>> > 22 40 0.8270000 5 >>> > 23 60 0.9500000 5 >>> > 24 80 0.9950000 5 >>> > 25 100 1.0000000 5 >>> > 26 20 0.4345774 6 >>> > 27 40 0.6654726 6 >>> > 28 60 0.8480684 6 >>> > 29 80 0.9268951 6 >>> > 30 100 0.9723207 6 >>> > 31 120 0.9939966 6 >>> > 32 140 0.9992400 6" >>> > >>> > dproot <- read.table(text = Lines, header = TRUE) >>> > >>> > >>> > >>> > On Mon, Sep 21, 2015 at 12:22 PM, Jianling Fan <fanjianl...@gmail.com> >>> > wrote: >>> >> >>> >> Thanks Prof. Nash, >>> >> >>> >> Sorry for late reply. I am learning and trying to use your nlmrt >>> >> package since I got your email. It works good to mask a parameter in >>> >> regression but seems does work for my equation. I think the problem is >>> >> that the parameter I want to mask is a group-specific parameter and I >>> >> have a "[]" syntax in my equation. However, I don't have your 2014 >>> >> book on hand and couldn't find it in our library. So I am wondering if >>> >> nlxb works for group data? >>> >> Thanks a lot! >>> >> >>> >> following is my code and I got a error form it. >>> >> >>> >> > fitdp1<-nlxb(den~Rm[ref]/(1+(depth/d50)^c),data=dproot, >>> >> + start =c(Rm1=1.01, Rm2=1.01, Rm3=1.01, Rm4=6.65, >>> >> Rm5=1.01, Rm6=1, d50=20, c=-1), >>> >> + masked=c("Rm6")) >>> >> >>> >> Error in deriv.default(parse(text = resexp), names(start)) : >>> >> Function '`[`' is not in the derivatives table >>> >> >>> >> >>> >> Best regards, >>> >> >>> >> Jianling >>> >> >>> >> >>> >> On 20 September 2015 at 12:56, ProfJCNash <profjcn...@gmail.com> >>> >> wrote: >>> >> > I posted a suggestion to use nlmrt package (function nlxb to be >>> >> > precise), >>> >> > which has masked (fixed) parameters. Examples in my 2014 book on >>> >> > Nonlinear >>> >> > parameter optimization with R tools. However, I'm travelling just >>> >> > now, >>> >> > or >>> >> > would consider giving this a try. >>> >> > >>> >> > JN >>> >> > >>> >> > >>> >> > On 15-09-20 01:19 PM, Jianling Fan wrote: >>> >> >> >>> >> >> no, I am doing a regression with 6 group data with 2 shared >>> >> >> parameters >>> >> >> and 1 different parameter for each group data. the parameter I want >>> >> >> to >>> >> >> coerce is for one group. I don't know how to do it. Any suggestion? >>> >> >> >>> >> >> Thanks! >>> >> >> >>> >> >> On 19 September 2015 at 13:33, Jeff Newmiller >>> >> >> <jdnew...@dcn.davis.ca.us> >>> >> >> wrote: >>> >> >>> >>> >> >>> Why not rewrite the function so that value is not a parameter? >>> >> >>> >>> >> >>> >>> >> >>> >>> >> >>> --------------------------------------------------------------------------- >>> >> >>> Jeff Newmiller The ..... ..... >>> >> >>> Go >>> >> >>> Live... >>> >> >>> DCN:<jdnew...@dcn.davis.ca.us> Basics: ##.#. ##.#. >>> >> >>> Live >>> >> >>> Go... >>> >> >>> Live: OO#.. Dead: OO#.. >>> >> >>> Playing >>> >> >>> Research Engineer (Solar/Batteries O.O#. #.O#. >>> >> >>> with >>> >> >>> /Software/Embedded Controllers) .OO#. .OO#. >>> >> >>> rocks...1k >>> >> >>> >>> >> >>> >>> >> >>> >>> >> >>> --------------------------------------------------------------------------- >>> >> >>> Sent from my phone. Please excuse my brevity. >>> >> >>> >>> >> >>> On September 18, 2015 9:54:54 PM PDT, Jianling Fan >>> >> >>> <fanjianl...@gmail.com> wrote: >>> >> >>>> >>> >> >>>> Hello, everyone, >>> >> >>>> >>> >> >>>> I am using a nls regression with 6 groups data. I am trying to >>> >> >>>> coerce >>> >> >>>> a parameter to 1 by using a upper and lower statement. but I >>> >> >>>> always >>> >> >>>> get an error like below: >>> >> >>>> >>> >> >>>> Error in ifelse(internalPars < upper, 1, -1) : >>> >> >>>> (list) object cannot be coerced to type 'double' >>> >> >>>> >>> >> >>>> does anyone know how to fix it? >>> >> >>>> >>> >> >>>> thanks in advance! >>> >> >>>> >>> >> >>>> My code is below: >>> >> >>>> >>> >> >>>> >>> >> >>>> >>> >> >>>>> dproot >>> >> >>>> >>> >> >>>> depth den ref >>> >> >>>> 1 20 0.5730000 1 >>> >> >>>> 2 40 0.7800000 1 >>> >> >>>> 3 60 0.9470000 1 >>> >> >>>> 4 80 0.9900000 1 >>> >> >>>> 5 100 1.0000000 1 >>> >> >>>> 6 10 0.6000000 2 >>> >> >>>> 7 20 0.8200000 2 >>> >> >>>> 8 30 0.9300000 2 >>> >> >>>> 9 40 1.0000000 2 >>> >> >>>> 10 20 0.4800000 3 >>> >> >>>> 11 40 0.7340000 3 >>> >> >>>> 12 60 0.9610000 3 >>> >> >>>> 13 80 0.9980000 3 >>> >> >>>> 14 100 1.0000000 3 >>> >> >>>> 15 20 3.2083491 4 >>> >> >>>> 16 40 4.9683383 4 >>> >> >>>> 17 60 6.2381133 4 >>> >> >>>> 18 80 6.5322348 4 >>> >> >>>> 19 100 6.5780660 4 >>> >> >>>> 20 120 6.6032064 4 >>> >> >>>> 21 20 0.6140000 5 >>> >> >>>> 22 40 0.8270000 5 >>> >> >>>> 23 60 0.9500000 5 >>> >> >>>> 24 80 0.9950000 5 >>> >> >>>> 25 100 1.0000000 5 >>> >> >>>> 26 20 0.4345774 6 >>> >> >>>> 27 40 0.6654726 6 >>> >> >>>> 28 60 0.8480684 6 >>> >> >>>> 29 80 0.9268951 6 >>> >> >>>> 30 100 0.9723207 6 >>> >> >>>> 31 120 0.9939966 6 >>> >> >>>> 32 140 0.9992400 6 >>> >> >>>> >>> >> >>>>> fitdp<-nls(den~Rm[ref]/(1+(depth/d50)^c),data=dproot, >>> >> >>>> >>> >> >>>> + start = list(Rm=c(1.01, 1.01, 1.01, 6.65,1.01,1), d50=20, >>> >> >>>> c=-1)) >>> >> >>>>> >>> >> >>>>> summary(fitdp) >>> >> >>>> >>> >> >>>> >>> >> >>>> Formula: den ~ Rm[ref]/(1 + (depth/d50)^c) >>> >> >>>> >>> >> >>>> Parameters: >>> >> >>>> Estimate Std. Error t value Pr(>|t|) >>> >> >>>> Rm1 1.12560 0.07156 15.73 3.84e-14 *** >>> >> >>>> Rm2 1.57643 0.11722 13.45 1.14e-12 *** >>> >> >>>> Rm3 1.10697 0.07130 15.53 5.11e-14 *** >>> >> >>>> Rm4 7.23925 0.20788 34.83 < 2e-16 *** >>> >> >>>> Rm5 1.14516 0.07184 15.94 2.87e-14 *** >>> >> >>>> Rm6 1.03658 0.05664 18.30 1.33e-15 *** >>> >> >>>> d50 22.69426 1.03855 21.85 < 2e-16 *** >>> >> >>>> c -1.59796 0.15589 -10.25 3.02e-10 *** >>> >> >>>> --- >>> >> >>>> Signif. codes: 0 ?**?0.001 ?*?0.01 ??0.05 ??0.1 ??1 >>> >> >>>> >>> >> >>>> Residual standard error: 0.1094 on 24 degrees of freedom >>> >> >>>> >>> >> >>>> Number of iterations to convergence: 8 >>> >> >>>> Achieved convergence tolerance: 9.374e-06 >>> >> >>>> >>> >> >>>>> fitdp1<-nls(den~Rm[ref]/(1+(depth/d50)^c),data=dproot, >>> >> >>>> >>> >> >>>> algorithm="port", >>> >> >>>> + start = list(Rm=c(1.01, 1.01, 1.01, 6.65, 1.01, 1), d50=20, >>> >> >>>> c=-1), >>> >> >>>> + lower = list(Rm=c(1.01, 1.01, 1.01, 6.65, 1.01, 1), d50=20, >>> >> >>>> c=-1), >>> >> >>>> + upper = list(Rm=c(2.1, 2.2, 2.12, 12.5, 2.3, 1), d50=50, c=1)) >>> >> >>>> >>> >> >>>> Error in ifelse(internalPars < upper, 1, -1) : >>> >> >>>> (list) object cannot be coerced to type 'double' >>> >> >>>> >>> >> >>>> ______________________________________________ >>> >> >>>> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see >>> >> >>>> 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. >>> >> >>> >>> >> >>> >>> >> >> >>> >> >> >>> >> >> >>> >> > >>> >> > ______________________________________________ >>> >> > R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see >>> >> > 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. >>> >> >>> >> >>> >> >>> >> -- >>> >> Jianling Fan >>> >> 樊建凌 >>> >> >>> >> ______________________________________________ >>> >> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see >>> >> 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. >>> > >>> > >>> > >>> > >>> > -- >>> > Statistics & Software Consulting >>> > GKX Group, GKX Associates Inc. >>> > tel: 1-877-GKX-GROUP >>> > email: ggrothendieck at gmail.com >>> >>> >>> >>> -- >>> Jianling Fan >>> 樊建凌 >> >> >> >> >> -- >> Statistics & Software Consulting >> GKX Group, GKX Associates Inc. >> tel: 1-877-GKX-GROUP >> email: ggrothendieck at gmail.com > > > > > -- > Statistics & Software Consulting > GKX Group, GKX Associates Inc. > tel: 1-877-GKX-GROUP > email: ggrothendieck at gmail.com -- Jianling Fan 樊建凌 ______________________________________________ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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.