options crt ;
name hhex 'Minimum Distance Estimation of a Dynamic Factor Model' ;
?
?  TSP program to attempt an estimate of 2-factor q,ds,dr,dk
?  model using LSQ.  This model is described in Hall and Hayashi,
?  NBER Working Paper No. 2973.
?  Version: n=3, tau=2, lagged q. 
?
?  TSP 4.1 (1990?)
?  Author: Bronwyn Hall
?
smpl 1 2052 ;
read (file='hhdif85a.dat')
    ardsic 
?   Current value and 4 lags of Stock market ROR, Log Sales, 
?   Log R&D, Log Capital.     
    q0 q1 q2 q3 q4
    s0 s1 s2 s3 s4
    r0 r1 r2 r3 r4
    k0 k1 k2 k3 k4
?    
?   Products of the same variables. 
    q4s0 q3s0 q2s0 q1s0 q0s0 q0s1 q0s2 q0s3 q0s4
    q4r0 q3r0 q2r0 q1r0 q0r0 q0r1 q0r2 q0r3 q0r4
    q4k0 q3k0 q2k0 q1k0 q0k0 q0k1 q0k2 q0k3 q0k4
    s4r0 s3r0 s2r0 s1r0 s0r0 s0r1 s0r2 s0r3 s0r4
    s4k0 s3k0 s2k0 s1k0 s0k0 s0k1 s0k2 s0k3 s0k4
    r4k0 r3k0 r2k0 r1k0 r0k0 r0k1 r0k2 r0k3 r0k4
;
?
set nob = @nob ; set factor = 1/nob ;
const zero 0 ;
frml eqq0 q0 = sx*sx+sy*sy+sz*sz ; ? Current q has 3 factors. 
                                   ? z is idiosyncratic.
dot 1-4 ;
  frml eqq. q. = zero ;            ? Model for the covariance of 
  frml eqq.s0 q.s0 = zero ;        ? q with current variables. 
  frml eqq.r0 q.r0 = zero ;        ? (strong market efficiency).
  frml eqq.k0 q.k0 = zero ;
enddot ;
dot 0 1 2 3 4;
  frml eqq0s. q0s. = a.*sx + h.*sy ;  ? x and y factors allowed to enter
  frml eqq0r. q0r. = l.*sx+g.*sy ;    ? equations for s, r, k. 
  frml eqq0k. q0k. = b.*sx+d.*sy ;
enddot ;
?
?   General specification of the equations for the covariance of the
?   observed variables.  
frml eqs0 s0 = (a0*a0+a1*a1+a2*a2+a3*a3+a4*a4)
             + (h0*h0+h1*h1+h2*h2+h3*h3+h4*h4) ;
frml eqs1 s1 = (a0*a1+a1*a2+a2*a3+a3*a4)
             + (h0*h1+h1*h2+h2*h3+h3*h4) ;
frml eqs2 s2 = (a0*a2+a1*a3+a2*a4)
             + (h0*h2+h1*h3+h2*h4) ;
frml eqs3 s3 = (a0*a3+a1*a4) + (h0*h3+h1*h4) ;
frml eqs4 s4 = (a0*a4) + (h0*h4) ;
 
frml eqr0 r0 = (l0*l0+l1*l1+l2*l2+l3*l3+l4*l4)
             + (g0*g0+g1*g1+g2*g2+g3*g3+g4*g4) ;
frml eqr1 r1 = (l0*l1+l1*l2+l2*l3+l3*l4)
             + (g0*g1+g1*g2+g2*g3+g3*g4) ;
frml eqr2 r2 = (l0*l2+l1*l3+l2*l4)
             + (g0*g2+g1*g3+g2*g4) ;
frml eqr3 r3 = (l0*l3+l1*l4) + (g0*g3+g1*g4) ;
frml eqr4 r4 = (l0*l4) + (g0*g4) ;

frml eqk0 k0 = (b0*b0+b1*b1+b2*b2+b3*b3+b4*b4)
             + (d0*d0+d1*d1+d2*d2+d3*d3+d4*d4) +sige*sige ;
frml eqk1 k1 = (b0*b1+b1*b2+b2*b3+b3*b4)
             + (d0*d1+d1*d2+d2*d3+d3*d4) ;
frml eqk2 k2 = (b0*b2+b1*b3+b2*b4)
             + (d0*d2+d1*d3+d2*d4) ;
frml eqk3 k3 = (b0*b3+b1*b4) + (d0*d3+d1*d4) ;
frml eqk4 k4 = (b0*b4) + (d0*d4) ;
 
  
frml eqs4r0 s4r0 = (a4*l0)
                 + (h4*g0) ;
frml eqs3r0 s3r0 = (a3*l0+a4*l1)
                 + (h3*g0+h4*g1) ;
frml eqs2r0 s2r0 = (a2*l0+a3*l1+a4*l2)
                 + (h2*g0+h3*g1+h4*g2) ;
frml eqs1r0 s1r0 = (a1*l0+a2*l1+a3*l2+a4*l3)
                 + (h1*g0+h2*g1+h3*g2+h4*g3) ;
frml eqs0r0 s0r0 = (a0*l0+a1*l1+a2*l2+a3*l3+a4*l4)
                 + (h0*g0+h1*g1+h2*g2+h3*g3+h4*g4) ;
frml eqs0r1 s0r1 = (a0*l1+a1*l2+a2*l3+a3*l4)
                 + (h0*g1+h1*g2+h2*g3+h3*g4) ;
frml eqs0r2 s0r2 = (a0*l2+a1*l3+a2*l4)
                 + (h0*g2+h1*g3+h2*g4) ;
frml eqs0r3 s0r3 = (a0*l3+a1*l4)
                 + (h0*g3+h1*g4) ;
frml eqs0r4 s0r4 = (a0*l4)
                 + (h0*g4) ;
  
frml eqs4k0 s4k0 = (a4*b0)
                 + (h4*d0) ;
frml eqs3k0 s3k0 = (a3*b0+a4*b1)
                 + (h3*d0+h4*d1) ;
frml eqs2k0 s2k0 = (a2*b0+a3*b1+a4*b2)
                 + (h2*d0+h3*d1+h4*d2) ;
frml eqs1k0 s1k0 = (a1*b0+a2*b1+a3*b2+a4*b3)
                 + (h1*d0+h2*d1+h3*d2+h4*d3) ;
frml eqs0k0 s0k0 = (a0*b0+a1*b1+a2*b2+a3*b3+a4*b4)
                 + (h0*d0+h1*d1+h2*d2+h3*d3+h4*d4) ;
frml eqs0k1 s0k1 = (a0*b1+a1*b2+a2*b3+a3*b4)
                 + (h0*d1+h1*d2+h2*d3+h3*d4) ;
frml eqs0k2 s0k2 = (a0*b2+a1*b3+a2*b4)
                 + (h0*d2+h1*d3+h2*d4) ;
frml eqs0k3 s0k3 = (a0*b3+a1*b4)
                 + (h0*d3+h1*d4) ;
frml eqs0k4 s0k4 = (a0*b4)
                 + (h0*d4) ;
 
frml eqr4k0 r4k0 = (l4*b0)
                 + (g4*d0) ;
frml eqr3k0 r3k0 = (l3*b0+l4*b1)
                 + (g3*d0+g4*d1) ;
frml eqr2k0 r2k0 = (l2*b0+l3*b1+l4*b2)
                 + (g2*d0+g3*d1+g4*d2) ;
frml eqr1k0 r1k0 = (l1*b0+l2*b1+l3*b2+l4*b3)
                 + (g1*d0+g2*d1+g3*d2+g4*d3) ;
frml eqr0k0 r0k0 = (l0*b0+l1*b1+l2*b2+l3*b3+l4*b4)
                 + (g0*d0+g1*d1+g2*d2+g3*d3+g4*d4) ;
frml eqr0k1 r0k1 = (l0*b1+l1*b2+l2*b3+l3*b4)
                 + (g0*d1+g1*d2+g2*d3+g3*d4) ;
frml eqr0k2 r0k2 = (l0*b2+l1*b3+l2*b4)
                 + (g0*d2+g1*d3+g2*d4) ;
frml eqr0k3 r0k3 = (l0*b3+l1*b4)
                 + (g0*d3+g1*d4) ;
frml eqr0k4 r0k4 = (l0*b4)
                 + (g0*d4) ;
 
list eqlist eqq0-eqq4
            eqq4s0 eqq3s0 eqq2s0 eqq1s0
            eqq0s0 eqq0s1 eqq0s2 eqq0s3 eqq0s4
            eqq4r0 eqq3r0 eqq2r0 eqq1r0
            eqq0r0 eqq0r1 eqq0r2 eqq0r3 eqq0r4
            eqq4k0 eqq3k0 eqq2k0 eqq1k0
            eqq0k0 eqq0k1 eqq0k2 eqq0k3 eqq0k4
            eqs0-eqs4
            eqs0r4 eqs0r3 eqs0r2 eqs0r1 eqs0r0
            eqs1r0 eqs2r0 eqs3r0 eqs4r0
            eqs0k4 eqs0k3 eqs0k2 eqs0k1 eqs0k0
            eqs1k0 eqs2k0 eqs3k0 eqs4k0
            eqr0-eqr4 
            eqk0-eqk4
            eqr0k4 eqr0k3 eqr0k2 eqr0k1 eqr0k0
            eqr1k0 eqr2k0 eqr3k0 eqr4k0
            ;
list covlist q0-q4
            q4s0 q3s0 q2s0 q1s0
            q0s0 q0s1 q0s2 q0s3 q0s4
            q4r0 q3r0 q2r0 q1r0
            q0r0 q0r1 q0r2 q0r3 q0r4
            q4k0 q3k0 q2k0 q1k0
            q0k0 q0k1 q0k2 q0k3 q0k4
            s0-s4
            s0r4 s0r3 s0r2 s0r1 s0r0
            s1r0 s2r0 s3r0 s4r0
            s0k4 s0k3 s0k2 s0k1 s0k0
            s1k0 s2k0 s3k0 s4k0
            r0-r4 
            k0-k4
            r0k4 r0k3 r0k2 r0k1 r0k0
            r1k0 r2k0 r3k0 r4k0
            ;
msd covlist ;
cova (noprint) covlist ;   ? Compute the matrix of fourth moments. 
mform fourth = @cova ;     ? Save it. 
?
?   Correct for degrees of freedom in the eventual SUR estimation
?
mmult fourth factor fourth ;
smpl 1 1 ;
unmake @mean covlist ;     ? Save the values of the second moments in
                           ? the first observation of the data; from
                           ? now on, we do not use the rest of the data,
                           ? since these are sufficient statistics. 
?
?    Substitute in the constraints implied by the structural model.
?
?         x factor in ds equation.
frml eqa0 a0 = lam0 ;
frml eqa1 a1 = lam1-lam0 ;
frml eqa2 a2 = lam2-lam1+beta1+sigma*eta2 ;
frml eqa3 a3 = -lam2+beta2-beta1+sigma*eta3 ;
frml eqa4 a4 = beta3-beta2+sigma*eta4 ;
frml eqa5 a5 = beta4-beta3+sigma*eta5 ;
frml eqa6 a6 = -beta4 ;
?         x factor in dk equation.
frml eqb0 b0 = beta1+sigma*eta2 ;
frml eqb1 b1 = beta2-beta1+sigma*eta3 ;
frml eqb2 b2 = beta3-beta2+sigma*eta4 ;
frml eqb3 b3 = beta4-beta3+sigma*eta5 ;
frml eqb4 b4 = -beta4 ;
?         x factor in dr equation.
frml eql0 l0 = eta2 ;
frml eql1 l1 = eta3+(sigma-1)*eta2 ;
frml eql2 l2 = eta4+(sigma-1)*eta3 ;
frml eql3 l3 = eta5+(sigma-1)*eta4 ;
frml eql4 l4 = (sigma-1)*eta5 ;
?         y factor in ds equation.
frml eqh0 h0 = zero ;
frml eqh1 h1 = zero ;
frml eqh2 h2 = delt0+sigma*gam0 ;
frml eqh3 h3 = delt1+sigma*gam1 ;
frml eqh4 h4 = delt2+sigma*gam2 ;
?        y factor in dk equation.
frml eqd0 d0 = delt0+sigma*gam0 ;
frml eqd1 d1 = delt1+sigma*gam1 ;
frml eqd2 d2 = delt2+sigma*gam2 ;
frml eqd3 d3 = delt3+sigma*gam3 ;
frml eqd4 d4 = delt4+sigma*gam4 ;
?        y factor in dr equation.
frml eqg0 g0 = gam0 ;
frml eqg1 g1 = delt0+(sigma-1)*gam0+gam1 ;
frml eqg2 g2 = delt1+(sigma-1)*gam1+gam2 ;
frml eqg3 g3 = delt2+(sigma-1)*gam2+gam3 ;
frml eqg4 g4 = delt3+(sigma-1)*gam3+gam4 ;
dot covlist ;
  eqsub eq. eqa0-eqa4 eqb0-eqb4 eql0-eql4
            eqh0-eqh4 eqd0-eqd4 eqg0-eqg4 ;
enddot ;
?
?      Estimate the model in stages. 
?
param sige .064 lam0 .002  lam1 .036 lam2 .0014
 beta1 .02348  beta2 .03302
 eta2 -.00386 eta3 .0144 ;
const gam0-gam4 delt0-delt4 ;
?
?    Constrain n to be 3 (the order of x, or demand, process).
const beta3 0 beta4 0 eta4 0 eta5 0 ;
const sigma 1 ;
const sy 0 ; param sx .064 sz .1944;
supres covu covt w regout ;
search sigma .3 .3 .2 ;     ?  Estimate by first searching over sigma, 
                            ?  the technology accumulation parameter.
param delt0-delt4 gam0-gam4 ;
search sy .005 .05 .005 ;   ?  And then search over the technology factor. 
?
? ========================================================================
?
proc search param minval maxval incr ;
?
?    Procedure to minimize the function by searching over one
?  difficult parameter and choosing the best value as a starting
?  value.
?
  local mintr minpar ; const param mintr minpar ;
  set mintr = 999999. ;set minpar = minval ;
  do param = minval to maxval by incr ;
    sur (silent,tol=.01,maxit=80,wname=fourth) eqlist ;
    if mintr>@tr ; then ; do ;
      set mintr = @tr ; set minpar = param ;
    enddo ;
    print @tr param ;
  enddo ;
  set param = minpar ; param param ;
  sur (tol=.001,maxit=80,wname=fourth) eqlist ;
endproc search ;
?
? ======================================================================
?
stop ;end ;