[R-sig-finance] Distribution Fitting

[Lorenzo Isella]

On Mon, 2006-04-24 at 15:08 +0200, Martin Maechler wrote:
> >>>>> "Lorenzo" == Lorenzo Isella <lorenzo.isella at gmail.com>
> >>>>>     on Mon, 24 Apr 2006 00:10:41 +0200 writes:
> 
>     Lorenzo> Dear All,
>     Lorenzo> I am experiencing some problems in fitting a trimodal distribution (which
>     Lorenzo> should be thought as a sum of three Gaussian distributions) using the nls
>     Lorenzo> function for nonlinear fittings.
>     Lorenzo> As a test, just consider the very simple code:
> 
>     Lorenzo> rm(list=ls());
>     Lorenzo> mydata<-rnorm(10000,0,4);
>     Lorenzo> mydens<-density(mydata,kernel="gaussian");
>     Lorenzo> y1<-mydens$y;
>     Lorenzo> x1<-mydens$x;
>     Lorenzo> myfit<-nls(y1~A*exp(-x1^2/sig));
> 
> (it's slightly inefficient and completely unnecessary to append
>  an empty statement to every line -- which you 
>  do incidentally if you add superfluous ';' at end of lines )
> 
> 
>     Lorenzo> which I use to get the empirical density (as I would from real experimental
>     Lorenzo> data) and test it against a Gaussian ansatz.
>     Lorenzo> Well, either R always crashes for a segmentation fault and I have to restart
>     Lorenzo> it manually or I get this output:
> 
>     Lorenzo> Error in match.call(definition, call, expand.dots) :
>     Lorenzo> '.Primitiv�i�d�������������
...' is not a function
> 
> I cannot confirm any segmentation fault (and if you really get
> them, indeed your R *installation* is buggy (or outdated)), but
> indeed, you've revealed a buglet in nls() - which is still
> present in R 2.3.0 which has been released a few hours ago.
> Thank you for reporting it as a reproducible example!
> 
> I'm going to try to fix the bug.
> 
>     Lorenzo> Am I missing the obvious or is there some bug in my R build?
> 
> As others have said (implicitly): you're missing the facts
> 
> -  that it's rather bad idea to fit "data" to a density estimate.
> -  if your density has a closed from, you should rather use
>    maximum likelihood, typically most easily  via MASS::fitdistr
> 
> Martin

Thanks everybody for their advice.
Granted that I had better use fitdistr for my particular problem, how do
I feed a particular density function of my own choice into it?
E.g., consider the case of Gaussian-distributed data and suppose there
is no "normal" argument for the fitdistr().

mysample<-rnorm(10000,2,4)

You want to model the data as drawn from the probability distribution
A*exp(-(x-mu)^2/sig), where x is an experimental data and have A, mu and
sigma as parameters to fit.
I read the help(fitdistr) page, but there is not an example with a
user-provided probability density and I did not get very far with my
attempts.
Best regards and thanks in advance for your help

Lorenzo