[R] R formula language---a min and max function?

Gabor Grothendieck ggrothendieck at gmail.com
Wed May 5 01:53:20 CEST 2010


On Tue, May 4, 2010 at 7:39 PM, Gabor Grothendieck
<ggrothendieck at gmail.com> wrote:
> For fixed c, y is linear in pmin(c, x) so the first two statements
> find c and the next solves the remaining linear problem
>

Here is a slightly shorter f:

f <- function(c) sum(resid(lm(y ~ pmin(c, x)))^2)

> f <- function(c) sum((y - fitted(lm(y ~ pmin(c, x))))^2)
> fit.c <- optimize(f, c(0, max(DF$x))); fit.c
>
> lm(y ~ pmin(fit.c$minimum, x))
>
>
> as a function of c
>
> On Tue, May 4, 2010 at 5:25 PM, ivo welch <ivo.welch at gmail.com> wrote:
>> thank you, david and gabor.  very much appreciated.  I should have
>> thought of setting the seed.  this was only an example, of course.
>> alas, such intermittent errors could still be of concern to me,
>> because I need to simulate this nls() to find out its properties under
>> the NULL, so I can't easily tolerate errors.  fortunately, I had the
>> window still open, so getting my y's out was easy, and the rounded
>> figures produce the same nls error.
>>
>>> cbind(x,round(y,3))
>>       x      y
>>  [1,]  1  5.017
>>  [2,]  2  7.993
>>  [3,]  3 11.014
>>  [4,]  4 13.998
>>  [5,]  5 17.003
>>  [6,]  6 19.977
>>  [7,]  7 23.011
>>  [8,]  8 25.991
>>  [9,]  9 29.003
>> [10,] 10 32.014
>> [11,] 11 31.995
>> [12,] 12 32.004
>> [13,] 13 32.012
>> [14,] 14 31.994
>> [15,] 15 31.998
>> [16,] 16 32.000
>> [17,] 17 32.009
>> [18,] 18 31.995
>> [19,] 19 32.000
>> [20,] 20 31.982
>>
>>> r1= nls( y~ a+b*pmin(c,x), start=list(a=2, b=3, c=10), trace=TRUE )
>> 0.002138 :   2  3 10
>> 0.002117 :   2.004  3.000  9.999
>> 0.002113 :   2.006  2.999 10.001
>> 0.002082 :   2.005  2.999 10.000
>> 0.002077 :   2.005  2.999 10.000
>> 0.002077 :   2.005  2.999 10.000
>> Error in nls(y ~ a + b * pmin(c, x), start = list(a = 2, b = 3, c = 10),  :
>>  step factor 0.000488281 reduced below 'minFactor' of 0.000976562
>>
>> I really don't care about this example, of course---only about
>> learning how to avoid nls() from dying on me.  so, any advice would be
>> appreciated.
>>
>> regards,
>>
>> /iaw
>>
>>
>>
>> ----
>> Ivo Welch (ivo.welch at brown.edu, ivo.welch at gmail.com)
>>
>>
>>
>> On Tue, May 4, 2010 at 3:59 PM, David Winsemius <dwinsemius at comcast.net> wrote:
>>>
>>> On May 4, 2010, at 3:52 PM, ivo welch wrote:
>>>
>>>> thank you, david.  indeed.  works great (almost).  an example for
>>>> anyone else googling this in the future:
>>>>
>>>>> x=1:20
>>>>> y= 2+3*ifelse(x>10, 10, x)+rnorm(20,0,0.01)
>>>>> r1= nls( y~ a+b*pmin(c,x), start=list(a=2, b=3, c=10), trace=TRUE )
>>>>
>>>> 0.002142 :   2  3 10
>>>> 0.002115 :   2.004  3.000 10.000
>>>> 0.002114 :   2.006  2.999 10.001
>>>> 0.002084 :   2.005  2.999 10.000
>>>> ...
>>>> 0.002079 :   2.005  2.999 10.000
>>>> Error in nls(y ~ a + b * pmin(c, x), start = list(a = 2, b = 3, c = 10),
>>>>  :
>>>>  step factor 0.000488281 reduced below 'minFactor' of 0.000976562
>>>>
>>>> strange error, but unrelated to my question.  will figure this one out
>>>> next.
>>>
>>> I get no error. May be difficult to sort out unless you can reproduce after
>>> setting a random seed.
>>>
>>>> x=1:20
>>>> y= 2+3*ifelse(x>10, 10, x)+rnorm(20,0,0.01)
>>>> r1= nls( y~ a+b*pmin(c,x), start=list(a=2, b=3, c=10), trace=TRUE )
>>> 0.001560045 :   2  3 10
>>> 0.001161253 :   2.003824  2.998973 10.000388
>>> 0.001161253 :   2.003824  2.998973 10.000388
>>>
>>> --
>>> David.
>>>
>>>>
>>>> regards,
>>>>
>>>> /iaw
>>>>
>>>>
>>>> On Tue, May 4, 2010 at 3:40 PM, David Winsemius <dwinsemius at comcast.net>
>>>> wrote:
>>>>>
>>>>> On May 4, 2010, at 3:33 PM, ivo welch wrote:
>>>>>
>>>>>> Dear R experts---I would like to estimate a non-linear least squares
>>>>>> expression that looks something like
>>>>>>
>>>>>>  y ~ a+b*min(c,x)
>>>>>>
>>>>>> where a, b, and c are the three parameters.  how do I define a min
>>>>>> function in the formula language of R?  advice appreciated.
>>>>>
>>>>> ?pmin
>>>>>
>>>>>>
>>>>>> sincerely,
>>>>>>
>>>>>> /iaw
>>>>>>
>>>>>> ______________________________________________
>>>>>> 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.
>>>>>
>>>>> David Winsemius, MD
>>>>> West Hartford, CT
>>>>>
>>>>>
>>>
>>> David Winsemius, MD
>>> West Hartford, CT
>>>
>>>
>>
>> ______________________________________________
>> 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.
>>
>



More information about the R-help mailing list