options crt;  ? Multivariate ARCH via ML
?  Example by Clint Cummins 9/2/92
?  p=neq=3 equations
?  n=2 lags  (pure ARCH, not GARCH)
?
?  basic multivariate normal likelihood function;
?  define the components below
?
frml march logl = -.5*log(dets) -.5*eie/dets + con*neq;
set con = -.5*log(6.2831853);   ? constant term -- log(2*pi)
?
set neq = 3;
?
?  s.. = sigma matrix
?  a.. = constant part of sigma matrix
?  b.. = coefficients of first lag
?  c.. = coefficients of second lag
?  e12 = equation 1, lag 2 -- define this below
?  e21 = equation 2, lag 1, etc.
?
frml s11 a11 + b11*e11*e11 + c11*e12*e12;
frml s21 a21 + b21*e21*e11 + c21*e22*e12;
frml s31 a31 + b31*e31*e11 + c31*e32*e12;
frml s22 a22 + b22*e21*e21 + c22*e22*e22;
frml s32 a32 + b32*e31*e21 + c32*e32*e22;
frml s33 a33 + b33*e31*e31 + c33*e32*e32;
list se s11 s21 s31 s22 s32 s33;
list ap a11 a21 a31 a22 a32 a33;
list bp b11 b21 b31 b22 b32 b33;
list cp c11 c21 c31 c22 c32 c33;
list allp ap bp cp;
?
?  i.. = sigma-inverse (times dets) = cofactors
?
frml i11 s22*s33 - s32*s32;
frml i21 -(s21*s33 - s31*s32);
frml i31 s21*s32 - s22*s31;
frml i22 s11*s33 - s31*s31;
frml i32 -(s11*s32 - s21*s31);
frml i33 s11*s22 - s21*s21;
list ie i11 i21 i31 i22 i32 i33;
?
?  dets = determinant of sigma
?
frml dets s11*i11 + s21*i21 + s31*i31;
?
?  eie = e'(sigma-inverse)*e (times dets)
?
frml eie e1*(e1*i11 + 2*e2*i21 + 2*e3*i31)
       + e2*(e2*i22 + 2*e3*i32)
       + e3*e3*i33;
?
? simple example of structural equations and their lags
?  e12 = equation 1, lag 2 -- define this below
?
frml e1 y1-a1-b1*x;
frml e11 y1(-1)-a1-b1*x(-1);
frml e12 y1(-2)-a1-b1*x(-2);
frml e2 y2-a2;
frml e21 y2(-1)-a2;
frml e22 y2(-2)-a2;
frml e3 y3-a3;
frml e31 y3(-1)-a3;
frml e32 y3(-2)-a3;
list elag0 e1-e3;
list e1lags e11-e12;
list e2lags e21-e22;
list e3lags e31-e32;
param a1-a3 b1;
set T = 100;
smpl 1,T;   ? draw some data for example
random(seed=491479) x;
random y1; y1 = 1+x+y1;
random(mean=2) y2;
random(mean=3) y3;
?
smpl 3,T;
lsq e1-e3;   ? estimate with LSQ to get Sigma starting values
?
param allp;  ? sets all covariance parameters to zero
eqsub march dets eie ie se elag0 e1lags-e3lags;
? initialize covariance parameters
set a11 = @covu(1,1);
set a21 = @covu(2,1);
set a31 = @covu(3,1);
set a22 = @covu(2,2);
set a32 = @covu(3,2);
set a33 = @covu(3,3);
options limnum=9999,limwarn=0;
const bp cp;  ? hold arch params at zero to start with
ml march;
param bp;  ? allow arch params on first lag to be non-zero
if (@ifconv); then; ml march;
param cp;  ? allow arch params on second lag to be non-zero
           ? (does not converge on this data, even in 100 iterations)
           ? (but note that the data-generating parameters are zero)
if (@ifconv); then; ml march;