OPTIONS CRT;  ? PDL DONE MANUALLY WITH LSQ -- complex example
              ? with FAR constraint and AR(2) error term
FREQ Q;
SMPL 53:1 67:4;
TREND T;
LOAD;
SMPL 55:3,67:4;
?FRML EQ1 Y = A0 + BT*T + B0*X + B1*X(-1) + B2*X(-2) + B3*X(-3) + B4*X(-4) +
?                               B5*X(-5) + B6*X(-6) + B7*X(-7) + B8*X(-8);
? AR(2) ERROR:  Y = X*B + U;  U = E + P1*U(-1) + P2*U(-2);
?       SO      Y = X*B + P1*[Y(-1)-X(-1)*B] + P2*[Y(-2)-X(-2)*B] + E;
FRML EQ2 Y = P1*Y(-1) + P2*Y(-2) + A0*(1-P1-P2) + BT*[T-P1*(T(-1)-P2*T(-2))] +
 B0*X + (B1-P1*B0)*X(-1) + (B2-P1*B1-P2*B0)*X(-2) + (B3-P1*B2-P2*B1)*X(-3) +
  (B4-P1*B3-P2*B2)*X(-4) + (B5-P1*B4-P2*B3)*X(-5) + (B6-P1*B5-P2*B4)*X(-6) +
  (B7-P1*B6-P2*B5)*X(-7) + (B8-P1*B7-P2*B6)*X(-8) +
    (-P1*B8-P2*B7)*X(-9) +         (-P2*B8)*X(-10);
AMAT 9,4,NZVAR,2,A;
PRINT A;
SET I = 1;
DOT 0-8;
  FRML EB. B. = C1*A.1 + C2*A.2 + C3*A.3;
  ? FRML EB0 B0 = C1*A01 + C2*A02 + C3*A03;
  ? ... ETC. ...
  ? FRML EB8 B8 = C1*A81 + C2*A82 + C3*A83;
  SET J = 1;
  DOT 1-3;
    SET IJ = I + (J-1)*9;
    SET A.. = A(IJ);  ? SETS A02 = A(1,2), ETC.
    SET J = J+1;
  ENDDOT;
  SET I = I+1;
ENDDOT;
EQSUB EQ2 EB0-EB8;
PARAM A0 BT C1-C3;   ? NOTE:  C1-C3 ARE THE "SCRAMBLED" COEFFICIENTS
? DO COMPARISONS FOR AR(0), AR(1), AND AR(2)
? AR(0)
OLSQ Y C T X(4,9,FAR);
CONST P1 P2;  ? HOLDS AR COEFFICIENTS AT ZERO (I.E. OLS ESTIMATION)
LSQ EQ2;
ANALYZ EB0-EB8;   ? GET "UNSCRAMBLED" COEFFICIENTS AND STANDARD ERRORS
? AR(1)
SMPL 55:2,67:4;  ? SPECIAL SAMPLE FOR AR1(CORC), WHICH DROPS FIRST OBS.
AR1(METHOD=CORC) Y C T X(4,9,FAR);
SMPL 55:3,67:4;
PARAM P1;  ? ESTIMATE FIRST AR COEFFICIENT
LSQ EQ2;
ANALYZ EB0-EB8;   ? GET "UNSCRAMBLED" COEFFICIENTS AND STANDARD ERRORS
? AR(2)
PARAM P2;  ? ESTIMATE BOTH AR COEFFICIENTS
LSQ EQ2;
ANALYZ EB0-EB8;   ? GET "UNSCRAMBLED" COEFFICIENTS AND STANDARD ERRORS
?
PROC AMAT P,N,NZVAR,JZERO,A;
?
? THIS MAKES THE LAGRANGIAN INTERPOLATION POLYNOMIALS USING COOPER'S FORMULA.
? CODED DIRECTLY FROM TSP SUBROUTINE AMAT.
?  P = NUMBER OF LAGS
?  N = DEGREE OF POLYNOMIAL
?  NZVAR = N - (NUMBER OF ZERO RESTRICTIONS)
?  JZERO = 1  BOTH NEAR AND FAR RESTRICTIONS
?          2  FAR ONLY
?          3  NEAR ONLY
?          4  NO RESTRICTIONS
?
LOCAL TAUJ,TAU0,TAUN,SPACNG,J,J1,DENOM,I,TAU,PHINUM,TJ;
SET NZVAR = N - 2*(JZERO=1) - 1*(JZERO=2 .OR. JZERO.EQ.3);
MFORM(NROW=P,NCOL=NZVAR) A=0;
MFORM(NROW=N,NCOL=1) TAUJ=0;
IF (N.GT.1); THEN; DO;
  SET TAU0 = 0;
  IF (JZERO.EQ.1 .OR. JZERO.EQ.3); THEN; SET TAU0 = -1;
  SET TAUN = P-1;
  IF (JZERO.LE.2); THEN; SET TAUN = P;
  SET SPACNG = (TAUN-TAU0)/(N-1);
  SET TAUJ(1) = TAU0;
  DO J=2,N;
    SET J1 = J-1;
    SET TAUJ(J) = TAUJ(J1) + SPACNG;
  ENDDO;
ENDDO;
DO J=1,NZVAR;
  SET JJ=J;
  IF (JZERO.EQ.1 .OR. JZERO.EQ.3); THEN; SET JJ = J+1;
  SET DENOM=1;
  DO I=1,N;
    IF (I.NE.JJ); THEN; SET DENOM = DENOM*(TAUJ(JJ)-TAUJ(I));
  ENDDO;
  DO TAU=1,P;
    SET PHINUM=1;
    DO I=1,N;
      IF (I.NE.JJ); THEN; SET PHINUM = PHINUM*(TAU-1-TAUJ(I));
    ENDDO;
    SET TJ = TAU + (J-1)*P;  ? FOR A(TAU,J) SUBSCRIPT
    SET A(TJ) = PHINUM/DENOM;
  ENDDO;
ENDDO;
ENDPROC;
END;
?
NOPRINT;
LOAD Y X;  ? ORIGINAL ALMON DATA FROM MADDALA, P. 370
2072 1660 2077 1926 2078 2181 2043 1897 2062 1695 2067 1705 1964 1731 1981 2151
1914 2556 1991 3152 2129 3763 2309 3903 2614 3912 2896 3571 3058 3199 3309 3262
3446 3476 3466 2993 3435 2262 3183 2011 2697 1511 2338 1631 2140 1990 2012 1993
2071 2520 2192 2804 2240 2919 2421 3024 2639 2725 2733 2321 2721 2131 2640 2552
2513 2234 2448 2282 2429 2533 2516 2517 2534 2772 2494 2380 2596 2568 2572 2944
2601 2629 2648 3133 2840 3449 2937 3764 3136 3983 3299 4381 3514 4786 3815 4094
4093 4870 4262 5344 4531 5433 4825 5911 5160 6109 5319 6542 5574 5785 5749 5707
5715 5412 5637 5465 5383 5550 5467 5465;