[R] function in nls argument -- robust estimation
Fernando Moyano
nanomoyano at yahoo.com
Tue May 20 12:39:25 CEST 2008
Hi Kate and others,
thanks for the info.
Btw, you sent the different
methods to analyze the data: nls, nls.lm and nlrob. Comparing the
results visually nlrob performed better then nls, but nls.lm (using the
0.9 quantile of residuals) was still better than nlrob. My data may
have a rather large amount of contamination, so that an M-estimator
with a higher breakdown point should be used (least trimmed squares?).
I haven't found this in R and wouldn't know how to implement it. But I
can live with my results. Then remains the question of obtaining the
parameter st. errors. Jackknife was suggested. Is there an R function I
could use for that?
cheers,
Fernando
Katharine Mullen wrote:
>
> Dear Martin,
>
> Thanks for the ideas regarding the relation of what Fernando is doing with
> robust regression. Indeed, it's an important point that he can't consider
> the standard error estimates on his parameters correct.
>
> I know from discussion off-list that he's happy with the results he has
> now; nevertheless the robust regression route may be an interesting
> alternative. I'm posting a scipt to R-SIG-robust now that compares the 3
> ways (nls, nlrob and nls.lm w/residuals above a certain quantile set to
> zero).
>
> best,
> Kate
>
> On Sat, 10 May 2008, Martin Maechler wrote:
>
>> Hi Kate and Fernando,
>>
>> I'm late into this thread,
>> but from reading it I get the impression that Fernando really
>> wants to do *robust* (as opposed to least-squares) non-linear
>> model fitting. His proposal to set residuals to zero when they
>> are outside a given bound is a very special case of an
>> M-estimator, namely (if I'm not mistaken) the so-called "Huber
>> skipped-mean", an M-estimator with psi-function
>> psi <- function(x, k) ifelse(abs(x) <= k, x, 0)
>> It is known that this can be far from optimal, and either using
>> Huber-psi or "a redescender" such as Tukey's biweight can be
>> considerably better.
>> Also note that the standard inference (std.errors, P-values, ...)
>> that you'd get from summary(nlsfit) or anova(nls1, nl2) is
>> *invalid* here, since you are effectively using *random* weighting.
>>
>> The nlrob() function in package 'robustbase'
>> implements M-estimation of nonlinear models directly.
>> Unfortunately, how to do correct inference in this situation
>> is a hard problem, probably even an open research question in
>> parts. I would expect that "the" bootstrap should work if you only
>> have a few outliers.
>>
>> I don't have time at the moment to look at the example data and
>> the model, and show you how to use it for nlrob();
>> if you find a way to you it for nls() , then the same should
>> work for nlrob().
>>
>> I'm CCing this to the specialists for "Robust Stats with R"
>> mailing list, R-SIG-robust.
>>
>> Best regards,
>> Martin Maechler
>> ETH Zurich
>>
>> >>>>> "KateM" == Katharine Mullen <kate at few.vu.nl>
>> >>>>> on Fri, 9 May 2008 15:50:08 +0200 (CEST) writes:
>>
>> KateM> You can take minpack.lm_1.1-0 (source code and MS Windows
>> build,
>> KateM> respectively) from here:
>>
>> KateM> http://www.nat.vu.nl/~kate/minpack.lm_1.1-0.tar.gz
>> KateM> http://www.nat.vu.nl/~kate/minpack.lm_1.1-0.zip
>>
>> KateM> The bug that occurs when nprint = 0 is fixed. Also fixed is
>> another
>> KateM> problem suggested your example: when the argument par is a
>> list, calling
>> KateM> summary on the output of nls.lm was not working.
>>
>> KateM> I'll submit the new version to CRAN soon.
>>
>> KateM> This disscusion has been fruitful - thanks for it.
>>
>> KateM> On Fri, 9 May 2008, Katharine Mullen wrote:
>>
>> >> You indeed found a bug. I can reproduce it (which I should have
>> tried to
>> >> do on other examples in the first place!). Thanks for finding it.
>> >>
>> >> It will be fixed in version 1.1-0 which I will submit to CRAN
>> soon.
>> >>
>> >> On Fri, 9 May 2008, elnano wrote:
>> >>
>> >> >
>> >> > Find the data (data_nls.lm_moyano.txt) here:
>> >> > ftp://ftp.bgc-jena.mpg.de/pub/outgoing/fmoyano
>> >> >
>> >> >
>> >> >
>> >> > Katharine Mullen wrote:
>> >> > >
>> >> > > Thanks for the details - it sounds like a bug. You can either
>> send me the
>> >> > > data in an email off-list or make it available on-line
>> somewhere, so that
>> >> > > I and other people can download it.
>> >> > >
>> >> > >
>> >> > > ______________________________________________
>> >> > > 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.
>> >> > >
>> >> > >
>> >> >
>> >> > --
>> >> > View this message in context:
>> http://www.nabble.com/function-in-nls-argument-tp17108100p17146812.html
>> >> > Sent from the R help mailing list archive at Nabble.com.
>> >> >
>> >> > ______________________________________________
>> >> > 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.
>> >> >
>> >>
>> >> ______________________________________________
>> >> 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.
>> >>
>>
>> KateM> ______________________________________________
>> KateM> R-help at r-project.org mailing list
>> KateM> https://stat.ethz.ch/mailman/listinfo/r-help
>> KateM> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> KateM> and provide commented, minimal, self-contained, reproducible
>> code.
>>
>
> ______________________________________________
> 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.
>
>
--
View this message in context: http://www.nabble.com/function-in-nls-argument-tp17108100p17337520.html
Sent from the R help mailing list archive at Nabble.com.
More information about the R-help
mailing list