# Heavy duty MLE

• Mar 17th 2010, 04:35 PM
Anonymous1
MLE challenge
I will be impressed if someone can get this. I will also be impressed if someone can explain why this problem is of great interest.

Two correlated random variables $\displaystyle x_i$ and $\displaystyle y_i$ are Bivariate normal

$\displaystyle \left[ \begin{array}{ccc} x_i \\ y_i\\ \end{array} \right]$ ~ $\displaystyle N(\left[ \begin{array}{ccc} 0 \\ 0\\ \end{array} \right], VC = \left[ \begin{array}{ccc} 1 & a\\ a & 1\\ \end{array} \right])$ Where $\displaystyle VC$ is the Variance Covariance matrix.

Given $\displaystyle n$ $\displaystyle i.i.d.$ samples the likelihood function is:

$\displaystyle lik= (2\pi)^{-n}|VC|^{-n/2} exp(-1/2\sum^n [x_i, y_i]VC ^{-1} \left[ \begin{array}{ccc} x_i \\ y_i\\ \end{array} \right]).$ (Assume $\displaystyle 1> a^2$)

Some useful information:

$\displaystyle |VC|= (1-a^2),$ $\displaystyle VC^{-1}= \frac{1}{1-a^2}\left[ \begin{array}{ccc} 1 & -a\\ -a & 1\\ \end{array} \right]$

Find the MLE and asymptotic variance.
• Mar 17th 2010, 06:44 PM
Anonymous1
I'm working out the answer now. If any one is curious I will post my solution.
• Mar 17th 2010, 10:17 PM
matheagle
I assume you want the MLE of a.
• Mar 18th 2010, 10:14 AM
Anonymous1
Yes I do. It's tedious to calculate, but can you see why this problem is interesting?
• Mar 18th 2010, 10:16 AM
Moo
Well, there are formulas to estimate a covariance matrix... so it should be possible to narrow it to a huh ?
• Mar 18th 2010, 10:35 AM
Anonymous1
Did I hear someone thinking maximization of a cubic?
• Mar 22nd 2010, 07:44 PM
mr fantastic
Quote:

Originally Posted by Anonymous1
Did I hear someone thinking maximization of a cubic?

Look, if you know the answer then say so and/or post it. If this was meant to be a challenge question to members then it should have been posted in the subforum whose link is given below and the rules of that subforum (plainly stated in the stickies) followed:

http://www.mathhelpforum.com/math-he...enge-problems/
• Apr 1st 2010, 09:20 PM
Anonymous1
Been working on this for a while now so I thought I'd share (Happy).

Here is a data correlation project that involves finding the solution to the above $\displaystyle MLE$ estimator, which entailed solving a cubic.

My instructor gave us this problem to point out that $\displaystyle MLE$ can run into trouble for "multiple hump" functions, $\displaystyle e.g.,$ cubics. This is because there may be various local maximums.

To make sure I found the proper solution to: $\displaystyle \frac{dl(\hat a)}{d} = 0$ $\displaystyle i.e.,$ the $\displaystyle global$ maximum, I used the naive estimator, $\displaystyle \hat a_0,$ as my initial guess $\displaystyle x_0 ,$ along with 10-step Newton's.

I hope someone finds this interesting. Cause I do!

Code:

u1= load('UNITED.txt'); u2= load('STATES.txt'); u3= load('HONG.txt'); u4= load('kONG.txt'); m1= sum(abs(u1.^2));    u1= u1./sqrt(m1); m2= sum(abs(u2.^2));    u2= u2./sqrt(m2); m3= sum(abs(u3.^2));    u3= u3./sqrt(m3); m4= sum(abs(u4.^2));    u4= u4./sqrt(m4); k= [20, 50, 100, 200, 300, 400]; a= sum(u1.*u2); aHatUS0= zeros(100,length(k));  aHatUS= zeros(100,length(k)); aHatHK0= zeros(100,length(k));  aHatHK= zeros(100,length(k)); for t= 1:100     for i= 1:length(k)         R= randn(length(u1),k(i));         v1= R'*u1;  v2= R'*u2; v3= R'*u3;  v4= R'*u4;         aHatUS0(t,i)= (1/k(i)).*sum(v1.*v2);         aHatHK0(t,i)= (1/k(i)).*sum(v3.*v4);                 xn1= aHatUS0;  xn2= aHatHK0;  n= length(v1);         for j=1:10             a1= xn1;    a2= xn2;                         f1= n*a1.^3 - a1.^2 .* sum (v1.*v2) + a1.*(-n+sum(v1.^2+v2.^2)) - sum (v1.*v2);             fp1= 3*n*a1.^2 - 2.*a1.* sum (v1.*v2) -n + sum(v1.^2+v2.^2);             f2= n*a2.^3 - a2.^2 .* sum (v3.*v4) + a2.*(-n+sum(v3.^2+v4.^2)) - sum (v3.*v4);             fp2= 3*n*a2.^2 - 2.*a2.* sum (v3.*v4) -n + sum(v3.^2+v4.^2);                         xn1= a1 - f1./fp1;  xn2= a2 - f2./fp2;         end         aHatUS(t,i)=xn1(1);         aHatHK(t,i)=xn2(1);     end end aHatUS0Var= (m1*m2+a^2)*k.^(-1);  aHatUSVar= (1-a^2)^2./((1+a^2)*k); aHatHK0Var= (m3*m4+a^2)*k.^(-1);  aHatHKVar= (1-a^2)^2./((1+a^2)*k); mSEaHatUS0= aHatUS0Var + (mean(aHatUS0-a)).^2; mSEaHatHK0= aHatHK0Var + (mean(aHatHK0-a)).^2; mSEaHatUS= aHatUSVar + (mean(aHatUS-a)).^2; mSEaHatHK= aHatHKVar + (mean(aHatHK-a)).^2; figure; hold on;    grid on; plot(k,mSEaHatUS0,'g'); plot(k,aHatUS0Var,'b'); plot(k,mSEaHatHK0,'r'); plot(k,aHatHK0Var,'y'); title('Figure 1: Empirical MSEs and Theoretical Variances') xlabel('k'); ylabel('MSE & Variance'); h=legend('Empirical MSE of aUS0_{hat}','Theoretical Variance of aUS0_{hat}','Empirical MSE of aHK0_{hat}','Theoretical Variance of aHK0_{hat}',4); set(h,'Location','NorthEast') hold off; figure; hold on; grid on; plot(k,mSEaHatUS,'g'); plot(k,aHatUSVar,'b'); plot(k,mSEaHatHK,'r'); plot(k,aHatHKVar,'y'); title('Figure 1: Empirical MSEs and Theoretical Variances') xlabel('k'); ylabel('MSE & Variance'); h=legend('Empirical MSE of aUS_{hat}','Theoretical Variance of aUS_{hat}','Empirical MSE of aHK_{hat}','Theoretical Variance of aHK_{hat}',4); set(h,'Location','NorthEast') hold off;
• Apr 2nd 2010, 10:59 AM
Anonymous1
Note that the final theoretical variance was:

$\displaystyle Var(\hat a)= \frac{(1-a^2)^2}{(1+a^2)n}$