function [Lhat,Vhat]=get_initial_L(Y,X,omega,Smat,resid,N1,Qparam)
[ystar,sig,cm]=transformH1(Y(:,1:N1),omega,resid(:,N1+1:end));
Z=X';
XX=vec(ystar');
iomega=invpd(sig);
Zomega=kron(Z,iomega);
ZZomega=kron(Z*Z',iomega);
z=XX-(kron(Z',eye(N1)))*Qparam;
Lhat=invpd(Smat'*(ZZomega)*Smat)*(Smat'*(Zomega)*z);
Vhat=invpd(Smat'*(ZZomega)*Smat);