Search All of the Math Forum:
Views expressed in these public forums are not endorsed by
Drexel University or The Math Forum.
|
|
Torsten
Posts:
1,128
Registered:
11/8/10
|
|
Re: Maximum Likelihood Estimation and Confidence Intervals
Posted:
Jan 18, 2013 5:23 AM
|
|
"Topi Kaaresoja" wrote in message <kdb1a9$75m$1@newscl01ah.mathworks.com>...
> Hi,
>
> I have conducted an experiment where users needed to judge where two events were simultaneous or not. The time (milliseconds) between the two events were varied like this:
>
> t=[0, 10, 20, 30, 50, 70, 100, 150, 300];
>
> The question was repeated n=96 times for each value of t. The proportion of "YES" answers is in a response vector
>
> y=[0.875 0.89583 0.89583 0.86458 0.69792 0.63542 0.34375 0.1875 0];
> x=n*y;
>
> I wanted to apply a Gaussian model to that data using maximum likelihood estimation (MLE) similarly than in an MLE Tutorial (http://people.physics.anu.edu.au/~tas110/Teaching/Lectures/L3/Material/Myung03.pdf).
>
> So the PDF of one x is binomial distribution, and my model is
> a*exp(-.5*(((t-mu)/sigma).^2);
>
> So, after some calculations, the negative log-likelihood function becomes:
>
> function loglik = gauss_mle_bin(w,t,y,n)
> % gauss_mle The log-likelihood function of the gaussian model
> % based on binary distribution on separate responses (Myung MLE tutorial)
> % t = vector or m latencies
> % y = vector of m proportion of ?simultaneous? responses
> % Estimates: w(1) = mu, w(2) = sigma, w(3) = a
> x=y*n;
> mu = w(1);
> sigma = w(2);
> a = w(3);
> h = -.5*(((t-mu)/sigma).^2);
> loglik = (-1) * (sum(log(a*exp(h)).*x) + sum(log(1-a*exp(h)).*(n-x)));
>
> I get nice estimates for mu, sigma and a with fminsearch or fminunc.
> 2.4895 77.6230 0.8956
>
> My problem is: How to get confidence intervals for these parameters? Function mle does them automatically, but I have no idea how to pass my data to mle, because my data is not equally spaced.
>
> I know that I can get Hessian matrix from fminunc and by computing an inverse (variance-covariance matrix) and taking square root I will get approximations of standard errors for the parameters. From them I can calculate the 95% confidence intervals as follows:
>
> hessian = hessian returned by fminunc
> varcov = inv(hessian);
> stderrs = sqrt(diag(varcov));
> c= stderrs*sqrt(chi2inv(0.95,1));
>
> I would like to double check if these are correct with mle function (or by other means) if possible.
>
> Best,
>
> Topi Kaaresoja
> Nokia Research Center
> Finland
help nlparci
help nlpredci
Best wishes
Torsten.
|
|
|
|
