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

[David Winsemius]

On May 4, 2010, at 5:25 PM, ivo welch 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.

There would of course be try()

But ... Again, No error on my device:
 > xy <- scan()
1:  1  5.017
3:  2  7.993
5:  3 11.014
7:  4 13.998
9:  5 17.003
11:  6 19.977
13:  7 23.011
15:  8 25.991
17:  9 29.003
19:  10 32.014
21:  11 31.995
23:  12 32.004
25:  13 32.012
27:  14 31.994
29:  15 31.998
31:  16 32.000
33:  17 32.009
35:  18 31.995
37:  19 32.000
39:  20 31.982
41:
Read 40 items
 > xym <-matrix(xy, ncol=2, byrow=TRUE)
 > colnames(xym)<-c("x","y")
 > r1= nls( y~ a+b*pmin(c,x), start=list(a=2, b=3, c=10), trace=TRUE )
0.001505770 :   2  3 10
0.001361737 :   2.003063  2.999924 10.000135

Plotting predict(r1) confirmed a very close fit.

 > sessionInfo()
R version 2.10.1 RC (2009-12-09 r50695)
x86_64-apple-darwin9.8.0

locale:
[1] en_US.UTF-8/en_US.UTF-8/C/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base

other attached packages:
[1] misc3d_0.7-0      rgl_0.91          spatstat_1.18-3   deldir_0.0-12
[5] mgcv_1.6-1        ks_1.6.12         mvtnorm_0.9-9      
KernSmooth_2.23-3
[9] lattice_0.18-3

loaded via a namespace (and not attached):
[1] grid_2.10.1        Matrix_0.999375-38 nlme_3.1-96         
tools_2.10.1


>
>> 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
>>
>>

David Winsemius, MD
West Hartford, CT