GFLS: Generalized Flexible Least Squares
for Approximately Linear Models

Last Updated: 15 April 2017

Site maintained by:
Leigh Tesfatsion
Professor of Econ, Math, and ECpE
Heady Hall 375
Iowa State University
Ames, Iowa 50011-1070
FAX: 515-294-0221
http://www2.econ.iastate.edu/tesfatsi/
tesfatsi AT iastate.edu

FLS Homepage
GFLS Help Page

About GFLS:

The Fortran program GFLS, developed by Robert E. Kalaba and Leigh Tesfatsion, implements a Generalized Flexible Least Squares (GFLS) method for the estimation of nonlinear systems described by approximately linear dynamic and measurement relations.

The basic GFLS objective is to determine the Cost-Efficient Frontier (CEF) for such systems in analogy to the determination of Pareto-efficient frontiers for multicriteria decision problems. The CEF is the collection of all state trajectory estimates that are efficient for the system at hand in the sense that they are minimally incompatible with the specified dynamic and measurement relations.

The GFLS method was proposed by Robert Kalaba and Leigh Tesfatsion, "Flexible Least Squares for Approximately Linear Systems" (pdf,1.2MB), IEEE Transactions on Systems, Man, and Cybernetics, Vol. 20, No. 5 (1990), pp. 978-989. The published article is available from IEEE Xplore.

The GFLS program has been incorporated into the statistical package GAUSS - TSM (Time Series Methods).

A GFLS Help Page is available. See, also, the FLS Homepage for additional information about GFLS.

The GFLS Program:

Software Release Disclaimer:

The software below is provided as-is, without warranty of any kind. It is released here as freeware by Robert E. Kalaba and Leigh Tesfatsion (the copyright holders) under the terms of the Artistic License Agreement.
 

// EXEC FORT1CLG,REGION=512K                                            00000010
//*                                                                     00000020
/*JOBPARM COPIES=4                                                      00000030
C                                                                       00000040
C     PROGRAM GFLS:                                                 00000050
C        EQUATION NUMBERS IN COMMENT STATEMENTS REFER TO                00000060
C        ''FLEXIBLE LEAST SQUARES FOR APPROXIMATELY LINEAR SYSTEMS''    00000070
C        BY R. KALABA AND L. TESFATSION, IEEE TRANSACTIONS ON SYSTEMS,  00000080
C        MAN, AND CYBERNETICS 20 (1990), PAGES 978-989.                 00000090
C     LAST UPDATED: 14 JUNE 1994                                        00000100
C                                                                       00000110
      IMPLICIT REAL*8(A-H,O-Z)                                          00000120
      INTEGER T,TCAP,TCAP1                                              00000130
      REAL*8 M                                                          00000140
C                                                                       00000150
C     THIS PROGRAM IS CURRENTLY DIMENSIONED FOR A MAXIMUM OF TCAP=110   00000160
C     OBSERVATION VECTORS Y OF MAXIMUM DIMENSION MOBS = 15 WITH STATE   00000170
C     VECTORS X OF MAXIMUM DIMENSION N = 15.                            00000180
C                                                                       00000190
      DIMENSION QO(15,15),PO(15,15),RO(15,15),QZERO(15,15)              00000200
      DIMENSION PZERO(15,15),RZERO(15,15)                               00000210
      DIMENSION F(15,15),A(15,15),H(15,15),B(15,15),D(15,15)            00000220
      DIMENSION M(15,15),Y(15,15),TRUEX(15,110),YY(15,110)              00000230
      DIMENSION HT(15,15),U(15,15),C(15,15),W(15,15),V(15,15)           00000240
      DIMENSION E(15,15),Z(15,15),G(15,15),QNEW(15,15),PNEW(15,15)      00000250
      DIMENSION S(15,15),RNEW(15,15),XTCAP(15,15),X(15,110)             00000260
      DIMENSION AA(15,15),BB(15,15),CC(15,15),DD(15,15),EE(15,15)       00000270
      DIMENSION FF(15,15),HH(15,15),OO(15,15),PP(15,15),QQ(15,15)       00000280
      DIMENSION RR(15,15),TT(15,15)                                     00000290
C                                                                       00000300
C     ADDITIONAL ARRAYS IF SMOOTHED ESTIMATES ARE TO BE CALCULATED      00000310
C     FOR INTERMEDIATE X VALUES (I.E., IF IFLAGS IS SET AT 1)           00000320
C                                                                       00000330
      DIMENSION GG(15,15,110),SS(15,110)                                00000340
C                                                                       00000350
C     THE FOLLOWING SUBROUTINE INITIALIZES THE PENALTY WEIGHT AMU,      00000360
C     THE NUMBER OF OBSERVATIONS TCAP, THE STATE VECTOR DIMENSION       00000370
C     N, THE OBSERVATION VECTOR DIMENSION MOBS, AND THE INITIAL COST    00000380
C     FUNCTION CHARACTERISTICS QZERO, PZERO, AND RZERO.  IT ALSO SETS   00000390
C     THE VALUE FOR A FLAG "IFLAGR" TO DETERMINE IF RNEW IS TO BE       00000400
C     CALCULATED (1) OR NOT (0) AND A FLAG "IFLAGS" TO DETERMINE IF     00000410
C     SMOOTHED ESTIMATES FOR INTERMEDIATE X VALUES ARE TO BE CALCULATED 00000420
C     (1) OR NOT (0).                                                   00000430
C                                                                       00000440
      CALL INPUT(AMU,TCAP,N,MOBS,QZERO,PZERO,RZERO,IFLAGR,IFLAGS)       00000450
      CALL SHIFT(N,N,QZERO,QO)                                          00000460
      CALL SHIFT(N,1,PZERO,PO)                                          00000470
      CALL SHIFT(1,1,RZERO,RO)                                          00000480
C                                                                       00000490
C     ENTERING THE MAIN DO LOOP FOR GENERATING Q,P,AND R FOR            00000500
C     SUCCESSIVE TIMES T = 1,TCAP USING EQS.(24), (26), AND (27).       00000510
C                                                                       00000520
      DO 50 T=1,TCAP                                                    00000530
      CALL MODEL(T,F,A,H,B,D,M,Y,TRUEX)                                 00000540
      DO 5 I=1,MOBS                                                     00000550
      YY(I,T) = Y(I,1)                                                  00000560
   5  CONTINUE                                                          00000570
C                                                                       00000580
C     GETTING U=HT*M*H + QO IN EQ.(28)                                  00000590
C                                                                       00000600
      CALL MUL(MOBS,MOBS,N,M,H,AA)                                      00000610
      CALL TRANS(MOBS,N,H,HT)                                           00000620
      CALL MUL(N,MOBS,N,HT,AA,BB)                                       00000630
      CALL ADD(N,N,BB,QO,U)                                             00000640
C                                                                       00000650
C     GETTING C=FT*D                                                    00000660
C                                                                       00000670
      CALL TRANS(N,N,F,AA)                                              00000680
      CALL MUL(N,N,N,AA,D,C)                                            00000690
C                                                                       00000700
C     GETTING W=AMU*C*F+U                                               00000710
C                                                                       00000720
      CALL MUL(N,N,N,C,F,AA)                                            00000730
      CALL MULCON(N,N,AMU,AA,BB)                                        00000740
      CALL ADD(N,N,BB,U,W)                                              00000750
C                                                                       00000760
C     GETTING V=WINV IN EQ.(17)                                         00000770
C                                                                       00000780
      CALL INV(N,W,V)                                                   00000790
C                                                                       00000800
C     GETTING E = (Y-B)                                                 00000810
C                                                                       00000820
      CALL SUB(MOBS,1,Y,B,E)                                            00000830
C                                                                       00000840
C     GETTING Z = HT*M*E + PO IN EQ.(29)                                00000850
C                                                                       00000860
      CALL MUL(MOBS,MOBS,1,M,E,AA)                                      00000870
      CALL MUL(N,MOBS,1,HT,AA,BB)                                       00000880
      CALL ADD(N,1,BB,PO,Z)                                             00000890
C                                                                       00000900
C     GETTING G = AMU*V*C IN EQ.(20)                                    00000910
C                                                                       00000920
      CALL MUL(N,N,N,V,C,AA)                                            00000930
      CALL MULCON(N,N,AMU,AA,G)                                         00000940
      IF(IFLAGS.EQ.0) GO TO 110                                         00000950
C                                                                       00000960
C     STORE G FOR CALCULATION OF SMOOTHED ESTIMATES                     00000970
C                                                                       00000980
      DO 10 I=1,N                                                       00000990
      DO 20 J=1,N                                                       00001000
      GG(I,J,T)=G(I,J)                                                  00001010
  20  CONTINUE                                                          00001020
  10  CONTINUE                                                          00001030
 110  CONTINUE                                                          00001040
C                                                                       00001050
C     GETTING QNEW = AMU*D*(I-F*G) IN EQ.(24)                           00001060
C                                                                       00001070
      CALL MUL(N,N,N,F,G,AA)                                            00001080
      CALL IDEN(N,BB)                                                   00001090
      CALL SUB(N,N,BB,AA,CC)                                            00001100
      CALL MUL(N,N,N,D,CC,DD)                                           00001110
      CALL MULCON(N,N,AMU,DD,QNEW)                                      00001120
C                                                                       00001130
C     GETTING PNEW = GT*Z+QNEWT*A IN EQ.(26)                            00001140
C                                                                       00001150
      CALL TRANS(N,N,G,AA)                                              00001160
      CALL MUL(N,N,1,AA,Z,BB)                                           00001170
      CALL TRANS(N,N,QNEW,CC)                                           00001180
      CALL MUL(N,N,1,CC,A,DD)                                           00001190
      CALL ADD(N,1,BB,DD,PNEW)                                          00001200
C                                                                       00001210
C     GETTING S = V*(Z - AMU*C*A) IN EQ.(19)                            00001220
C                                                                       00001230
      CALL MUL(N,N,1,C,A,BB)                                            00001240
      CALL MULCON(N,1,AMU,BB,CC)                                        00001250
      CALL SUB(N,1,Z,CC,DD)                                             00001260
      CALL MUL(N,N,1,V,DD,S)                                            00001270
      IF(IFLAGS.EQ.0) GO TO 210                                         00001280
C                                                                       00001290
C     STORE S FOR CALCULATION OF SMOOTHED ESTIMATES                     00001300
C                                                                       00001310
      DO 30 I=1,N                                                       00001320
      SS(I,T)=S(I,1)                                                    00001330
  30  CONTINUE                                                          00001340
 210  CONTINUE                                                          00001350
      IF(IFLAGR.EQ.0) GO TO 310                                         00001360
C                                                                       00001370
C     GETTING RNEW = RO + ET*M*E + AMU*AT*D*A - ST*W*S IN EQ.(27)       00001380
C                                                                       00001390
      CALL MUL(MOBS,MOBS,1,M,E,AA)                                      00001400
      CALL TRANS(MOBS,1,E,BB)                                           00001410
      CALL MUL(1,MOBS,1,BB,AA,CC)                                       00001420
      CALL ADD(1,1,RO,CC,DD)                                            00001430
      CALL MUL(N,N,1,D,A,EE)                                            00001440
      CALL TRANS(N,1,A,FF)                                              00001450
      CALL MUL(1,N,1,FF,EE,HH)                                          00001460
      CALL MULCON(1,1,AMU,HH,OO)                                        00001470
      CALL ADD(1,1,DD,OO,PP)                                            00001480
      CALL MUL(N,N,1,W,S,QQ)                                            00001490
      CALL TRANS(N,1,S,RR)                                              00001500
      CALL MUL(1,N,1,RR,QQ,TT)                                          00001510
      CALL SUB(1,1,PP,TT,RNEW)                                          00001520
 310  CONTINUE                                                          00001530
      IF(T.EQ.TCAP) GO TO 50                                            00001540
C                                                                       00001550
C     UPDATING QO,PO, AND RO                                            00001560
C                                                                       00001570
      CALL SHIFT(N,N,QNEW,QO)                                           00001580
      CALL SHIFT(N,1,PNEW,PO)                                           00001590
      IF(IFLAGR.EQ.0) GO TO 50                                          00001600
      CALL SHIFT(1,1,RNEW,RO)                                           00001610
  50  CONTINUE                                                          00001620
C                                                                       00001630
C     GETTING THE FLS FILTER ESTIMATE FOR XTCAP = UINV*Z IN EQ.(30)     00001640
C                                                                       00001650
      CALL INV(N,U,AA)                                                  00001660
      CALL MUL(N,N,1,AA,Z,XTCAP)                                        00001670
      DO 65 I=1,N                                                       00001680
      X(I,TCAP)=XTCAP(I,1)                                              00001690
  65  CONTINUE                                                          00001700
      IF (IFLAGS.EQ.1) GOTO 410                                         00001710
C                                                                       00001720
C     PRINTING OUT THE FLS FILTER ESTIMATE FOR XTCAP                    00001730
C                                                                       00001740
      CALL OUTPUT(TCAP,N,X,TRUEX)                                       00001750
      IF(IFLAGS.EQ.0) GOTO 510                                          00001760
 410  CONTINUE                                                          00001770
C                                                                       00001780
C     GETTING SMOOTHED ESTIMATES FOR X1,..., XTCAP-1 IN EQS.(33A)       00001790
C                                                                       00001800
      TCAP1=TCAP-1                                                      00001810
      DO 70 T=1,TCAP1                                                   00001820
      L=TCAP-T                                                          00001830
      DO 80 I=1,N                                                       00001840
      X(I,L)=SS(I,L)                                                    00001850
      DO 90 J=1,N                                                       00001860
      X(I,L)=X(I,L)+GG(I,J,L)*X(J,L+1)                                  00001870
  90  CONTINUE                                                          00001880
  80  CONTINUE                                                          00001890
  70  CONTINUE                                                          00001900
C                                                                       00001910
C     PRINTING OUT THE FLS ESTIMATES FOR X1,...,XTCAP                   00001920
C                                                                       00001930
      DO 150 T=1,TCAP                                                   00001940
      CALL OUTPUT(T,N,X,TRUEX)                                          00001950
 150  CONTINUE                                                          00001960
C     VALIDATION TEST: HOW WELL DO THE FLS ESTIMATES SATISFY THE        00001970
C     FIRST-ORDER CONDITIONS FOR THE COST MINIMIZATION PROBLEM (5)      00001980
      CALL FOCTST(X,YY)                                                 00001990
 510  CONTINUE                                                          00002000
      STOP                                                              00002010
      END                                                               00002020
C                                                                       00002030
C     MATRIX SUBROUTINES FOR ADDITION, MULTIPLICATION, TRANSPOSITION,   00002040
C     SUBTRACTION, INVERSION, MULTIPLICATION BY A SCALAR, SHIFT, AND    00002050
C     FORMATION OF AN IDENTITY MATRIX                                   00002060
C                                                                       00002070
C     OBTAINING THE SUM C=A+B OF TWO NROW X MCOL MATRICES A AND B       00002080
C                                                                       00002090
      SUBROUTINE ADD(NROW,MCOL,A,B,C)                                   00002100
      IMPLICIT REAL*8(A-H,O-Z)                                          00002110
      DIMENSION A(15,15),B(15,15),C(15,15)                              00002120
      DO 10 I=1,NROW                                                    00002130
      DO 20 J=1,MCOL                                                    00002140
      C(I,J)=A(I,J)+B(I,J)                                              00002150
  20  CONTINUE                                                          00002160
  10  CONTINUE                                                          00002170
      RETURN                                                            00002180
      END                                                               00002190
C                                                                       00002200
C     OBTAINING THE PRODUCT C=A*B OF AN NROW X L MATRIX A AND AN        00002210
C     L X MCOL MATRIX B                                                 00002220
C                                                                       00002230
      SUBROUTINE MUL(NROW,L,MCOL,A,B,C)                                 00002240
      IMPLICIT REAL*8(A-H,O-Z)                                          00002250
      DIMENSION A(15,15),B(15,15),C(15,15)                              00002260
      DO 10 I=1,NROW                                                    00002270
      DO 20 J=1,MCOL                                                    00002280
      SUM=0.0D+00                                                       00002290
      DO 30 K=1,L                                                       00002300
      SUM=SUM+A(I,K)*B(K,J)                                             00002310
  30  CONTINUE                                                          00002320
      C(I,J)=SUM                                                        00002330
  20  CONTINUE                                                          00002340
  10  CONTINUE                                                          00002350
      RETURN                                                            00002360
      END                                                               00002370
C                                                                       00002380
C     OBTAINING THE TRANSPOSE B OF AN NROW X MCOL MATRIX A              00002390
C                                                                       00002400
      SUBROUTINE TRANS(NROW,MCOL,A,B)                                   00002410
      IMPLICIT REAL*8(A-H,O-Z)                                          00002420
      DIMENSION A(15,15),B(15,15)                                       00002430
      DO 10 I=1,NROW                                                    00002440
      DO 20 J=1,MCOL                                                    00002450
      B(J,I)=A(I,J)                                                     00002460
  20  CONTINUE                                                          00002470
  10  CONTINUE                                                          00002480
      RETURN                                                            00002490
      END                                                               00002500
C                                                                       00002510
C     OBTAINING THE DIFFERENCE C=A-B BETWEEN NROW X MCOL MATRICES       00002520
C     A AND B                                                           00002530
C                                                                       00002540
      SUBROUTINE SUB(NROW,MCOL,A,B,C)                                   00002550
      IMPLICIT REAL*8(A-H,O-Z)                                          00002560
      DIMENSION A(15,15),B(15,15),C(15,15)                              00002570
      DO 10 I=1,NROW                                                    00002580
      DO 20 J=1,MCOL                                                    00002590
      C(I,J)=A(I,J)-B(I,J)                                              00002600
  20  CONTINUE                                                          00002610
  10  CONTINUE                                                          00002620
      RETURN                                                            00002630
      END                                                               00002640
C                                                                       00002650
C     OBTAINING THE INVERSE C OF A K X K MATRIX A                       00002660
C                                                                       00002670
      SUBROUTINE INV(K,A,C)                                             00002680
      IMPLICIT REAL*8(A-H,O-Z)                                          00002690
      DIMENSION A(15,15),B(15,30),C(15,15)                              00002700
      DO 5 J=1,K                                                        00002710
      DO 6 I=1,K                                                        00002720
      B(I,J)=A(I,J)                                                     00002730
   6  CONTINUE                                                          00002740
   5  CONTINUE                                                          00002750
      K2=K*2                                                            00002760
      DO 7 J=1,K                                                        00002770
      DO 8 I=1,K                                                        00002780
      B(I,K+J)=0.0D+00                                                  00002790
      IF(I.EQ.J) B(I,K+J)=1.0D+00                                       00002800
   8  CONTINUE                                                          00002810
   7  CONTINUE                                                          00002820
C                                                                       00002830
C     THE PIVOT OPERATION STARTS HERE                                   00002840
C                                                                       00002850
      DO 9 L=1,K                                                        00002860
      PIVOT = B(L,L)                                                    00002870
      DO 13 J=L,K2                                                      00002880
      B(L,J)=B(L,J)/PIVOT                                               00002890
  13  CONTINUE                                                          00002900
C                                                                       00002910
C     TO IMPROVE THE ROWS                                               00002920
C                                                                       00002930
      DO 14 I=1,K                                                       00002940
      IF(I.EQ.L) GO TO 14                                               00002950
      AIL=B(I,L)                                                        00002960
      DO 15 J=L,K2                                                      00002970
      B(I,J)=B(I,J)-AIL*B(L,J)                                          00002980
  15  CONTINUE                                                          00002990
  14  CONTINUE                                                          00003000
   9  CONTINUE                                                          00003010
      DO 45 I=1,K                                                       00003020
      DO 46 J=1,K                                                       00003030
      C(I,J)=B(I,K+J)                                                   00003040
  46  CONTINUE                                                          00003050
  45  CONTINUE                                                          00003060
      RETURN                                                            00003070
      END                                                               00003080
C                                                                       00003090
C     OBTAINING THE PRODUCT C*A OF A SCALAR C AND AN NROW X MCOL        00003100
C     MATRIX A                                                          00003110
C                                                                       00003120
      SUBROUTINE MULCON(NROW,MCOL,C,A,CA)                               00003130
      IMPLICIT REAL*8(A-H,O-Z)                                          00003140
      DIMENSION A(15,15),CA(15,15)                                      00003150
      DO 10 I=1,NROW                                                    00003160
      DO 20 J=1,MCOL                                                    00003170
      CA(I,J)=C*A(I,J)                                                  00003180
  20  CONTINUE                                                          00003190
  10  CONTINUE                                                          00003200
      RETURN                                                            00003210
      END                                                               00003220
C                                                                       00003230
C     PUTTING AN NROW X MCOL MATRIX A INTO AN NROW X MCOL MATRIX B      00003240
C                                                                       00003250
      SUBROUTINE SHIFT(NROW,MCOL,A,B)                                   00003260
      IMPLICIT REAL*8(A-H,O-Z)                                          00003270
      DIMENSION A(15,15),B(15,15)                                       00003280
      DO 10 I=1,NROW                                                    00003290
      DO 20 J=1,MCOL                                                    00003300
      B(I,J)=A(I,J)                                                     00003310
  20  CONTINUE                                                          00003320
  10  CONTINUE                                                          00003330
      RETURN                                                            00003340
      END                                                               00003350
C                                                                       00003360
C     FORMING THE N X N IDENTITY MATRIX E                               00003370
C                                                                       00003380
      SUBROUTINE IDEN(N,E)                                              00003390
      IMPLICIT REAL*8(A-H,O-Z)                                          00003400
      DIMENSION E(15,15)                                                00003410
      ZERO=0.0D+00                                                      00003420
      ONE=1.0D+00                                                       00003430
      DO 10 I=1,N                                                       00003440
      DO 20 J=1,N                                                       00003450
      E(I,J)=ZERO                                                       00003460
  20  CONTINUE                                                          00003470
  10  CONTINUE                                                          00003480
      DO 30 L=1,N                                                       00003490
      E(L,L)=ONE                                                        00003500
  30  CONTINUE                                                          00003510
      RETURN                                                            00003520
      END                                                               00003530
C                                                                       00003540
      SUBROUTINE INPUT(AMU,TCAP,N,MOBS,QZERO,PZERO,RZERO,IFLAGR,IFLAGS) 00003550
      IMPLICIT REAL*8(A-H,O-Z)                                          00003560
      INTEGER TCAP                                                      00003570
      DIMENSION QZERO(15,15),PZERO(15,15),RZERO(15,15)                  00003580
      AMU = 1.0D+00                                                     00003590
      TCAP = 30                                                         00003600
      N = 2                                                             00003610
      MOBS = 1                                                          00003620
      DO 10 J = 1,N                                                     00003630
      DO 20 I = 1,N                                                     00003640
      QZERO(I,J) = 0.0D+00                                              00003650
      PZERO(I,J) = 0.0D+00                                              00003660
      RZERO(I,J) = 0.0D+00                                              00003670
  20  CONTINUE                                                          00003680
  10  CONTINUE                                                          00003690
      IFLAGR=1                                                          00003700
      IFLAGS=1                                                          00003710
      RETURN                                                            00003720
      END                                                               00003730
C                                                                       00003740
C                                                                       00003750
      SUBROUTINE MODEL(T,F,A,H,B,D,M,Y,TRUEX)                           00003760
      IMPLICIT REAL*8(A-H,O-Z)                                          00003770
      REAL*8 M                                                          00003780
      INTEGER T,TCAP                                                    00003790
      DIMENSION F(15,15),A(15,15),H(15,15),B(15,15),D(15,15),M(15,15)   00003800
      DIMENSION Y(15,15),TRUEX(15,110),ZERO(15,15)                      00003810
      DIMENSION QZERO(15,15),PZERO(15,15),RZERO(15,15)                  00003820
C                                                                       00003830
C     REGIME SHIFT EXPERIMENT WITH POSSIBLY NOISY OBSERVATIONS          00003840
C                                                                       00003850
      CALL INPUT(AMU,TCAP,N,MOBS,QZERO,PZERO,RZERO,IFLAGR,IFLAGS)       00003860
      SIGMA = 0.00D+00                                                  00003870
      DO 10 I=1,15                                                      00003880
      DO 20 J=1,15                                                      00003890
      ZERO(I,J) = 0.0D+00                                               00003900
  20  CONTINUE                                                          00003910
  10  CONTINUE                                                          00003920
      CALL IDEN(N,F)                                                    00003930
      CALL SHIFT(N,1,ZERO,A)                                            00003940
      H(1,1)=1.0D+00                                                    00003950
      H(1,2)=1.0D+00                                                    00003960
      AT=DFLOAT(T)                                                      00003970
      IF(T.EQ.1) GO TO 200                                              00003980
      H(1,1)=DSIN(10.0D+00+(AT))+0.01D+00                               00003990
      H(1,2)=DCOS(10.0D+00+(AT))                                        00004000
 200  CONTINUE                                                          00004010
      CALL SHIFT(MOBS,1,ZERO,B)                                         00004020
      CALL IDEN(N,D)                                                    00004030
      CALL IDEN(MOBS,M)                                                 00004040
      IF (T.GT.15) GOTO 150                                             00004050
      TRUEX(1,T) = 2.0D+00                                              00004060
      TRUEX(2,T) = 3.0D+00                                              00004070
      GOTO 175                                                          00004080
 150  TRUEX(1,T) = 4.0D+00                                              00004090
      TRUEX(2,T) = 5.0D+00                                              00004100
 175  CONTINUE                                                          00004110
C     THE DISCREPANCY TERM UU IN THE NEXT LINE IS SET TO ZERO;
C     IT CAN BE REPLACED BY A CALL TO A PSEUDO-RANDOM NUMBER GENERATOR
C     TO TEST GFLS IN THE PRESENCE OF NOISE IN THE OBSERVATIONS
      UU = 0.0D+00                                                      00004120
      Y(1,1)=H(1,1)*TRUEX(1,T) + H(1,2)*TRUEX(2,T) + UU                 00004130
      RETURN                                                            00004140
      END                                                               00004150
C                                                                       00004160
      SUBROUTINE OUTPUT(T,N,X,TRUEX)                                    00004170
      IMPLICIT REAL*8(A-H,O-Z)                                          00004180
      INTEGER T                                                         00004190
      DIMENSION X(15,110),TRUEX(15,110)                                 00004200
      L = T                                                             00004210
      WRITE(6,100) L,(X(I,L),I=1,N)                                     00004220
 100  FORMAT(1H0,'TIME EQUALS',I3/1X,'FLS ESTIMATES',7X,2D25.10)        00004230
      WRITE(6,200) (TRUEX(I,L),I=1,N)                                   00004240
 200  FORMAT(1X,'TRUE X VALUES',7X,2D25.10)                             00004250
      RETURN                                                            00004260
      END                                                               00004270
C                                                                       00004280
C     VALIDATION TEST: HOW WELL DO THE FLS ESTIMATES SATISFY THE        00004290
C     FIRST-ORDER CONDITIONS FOR THE COST MINIMIZATION PROBLEM (5)      00004300
C                                                                       00004310
      SUBROUTINE FOCTST(X,YY)                                           00004320
      IMPLICIT REAL*8(A-H,O-Z)                                          00004330
      INTEGER T,TP1,TCAP,TCAP1                                          00004340
      REAL*8 M,MH                                                       00004350
      DIMENSION QZERO(15,15),PZERO(15,15),RZERO(15,15)                  00004360
      DIMENSION XT(15,15),X(15,110),XTT(15,15),E(15,15)                 00004370
      DIMENSION PZEROT(15,15),EE(15,15),CO(15,15),YT(15,15),YY(15,110)  00004380
      DIMENSION F(15,15),A(15,15),H(15,15),B(15,15),D(15,15),M(15,15)   00004390
      DIMENSION Y(15,15),TRUEX(15,110)                                  00004400
      DIMENSION MH(15,15),EM(15,15),EMT(15,15),W(15,15),XTP1(15,15)     00004410
      DIMENSION ED(15,15),EDT(15,15),U(15,15),V(15,15),FOCD(15,15)      00004420
      C = -1.0D+00                                                      00004430
C     FORM THE STATE VECTOR FOR TIME T = 1                              00004440
      CALL INPUT(AMU,TCAP,N,MOBS,QZERO,PZERO,RZERO,IFLAGR,IFLAGS)       00004450
      DO 100 I=1,N                                                      00004460
      XT(I,1) = X(I,1)                                                  00004470
 100  CONTINUE                                                          00004480
C     FORM THE INITIAL INCREMENTAL COST CO = -(X1'QO - PO')             00004490
      CALL TRANS(N,1,XT,XTT)                                            00004500
      CALL MUL(1,N,N,XTT,QZERO,E)                                       00004510
      CALL TRANS(N,1,PZERO,PZEROT)                                      00004520
      CALL SUB(1,N,E,PZEROT,EE)                                         00004530
      CALL MULCON(1,N,C,EE,CO)                                          00004540
C     DO LOOP FOR THE SEQUENTIAL CHECK OF THE FOC FOR T=1,TCAP          00004550
      DO 200 T=1,TCAP                                                   00004560
C     FORM THE TIME-T STATE VECTOR XT                                   00004570
      DO 300 I=1,N                                                      00004580
      XT(I,1) = X(I,T)                                                  00004590
 300  CONTINUE                                                          00004600
C     FORM THE TIME-T OBSERVATION VECTOR YT                             00004610
      DO 400 J=1,MOBS                                                   00004620
      YT(J,1) = YY(J,T)                                                 00004630
 400  CONTINUE                                                          00004640
      CALL MODEL(T,F,A,H,B,D,M,Y,TRUEX)                                 00004650
C     FORM W = (YT - H(T)XT - B(T))'M(T)H(T)                            00004660
      CALL MUL(MOBS,MOBS,N,M,H,MH)                                      00004670
      CALL RME(N,MOBS,YT,XT,H,B,EM)                                     00004680
      CALL TRANS(MOBS,1,EM,EMT)                                         00004690
      CALL MUL(1,MOBS,N,EMT,MH,W)                                       00004700
      IF(T.EQ.TCAP) GOTO 600                                            00004710
C     FORM THE TIME-T+1 STATE VECTOR XTP1                               00004720
      TP1 = T + 1                                                       00004730
      DO 500 I=1,N                                                      00004740
      XTP1(I,1) = X(I,TP1)                                              00004750
 500  CONTINUE                                                          00004760
C     FORM U = AMU*(XTP1 - F(T)XT - A(T))'*D(T)                         00004770
      CALL RDE(N,XTP1,XT,F,A,ED)                                        00004780
      CALL TRANS(N,1,ED,EDT)                                            00004790
      CALL MUL(1,N,N,EDT,D,E)                                           00004800
      CALL MULCON(1,N,AMU,E,U)                                          00004810
C     FORM V = U*F                                                      00004820
      CALL MUL(1,N,N,U,F,V)                                             00004830
      GOTO 800                                                          00004840
 600  CONTINUE                                                          00004850
      DO 700 I=1,N                                                      00004860
      V(1,I) = 0.0D+00                                                  00004870
 700  CONTINUE                                                          00004880
 800  CONTINUE                                                          00004890
C     DETERMINE THE FOC DISCREPANCIES FOR TIME T                        00004900
C     GIVEN BY FOCD = CO + V + W                                        00004910
      CALL ADD(1,N,CO,V,E)                                              00004920
      CALL ADD(1,N,E,W,FOCD)                                            00004930
C     PRINT OUT THE FOC DISCREPANCIES FOCD FOR TIME T                   00004940
      WRITE (6,36) T                                                    00004950
  36  FORMAT(1H0,'FOC DISCREPANCIES FOR TIME',I3)                       00004960
      WRITE (6,37) (FOCD(1,I),I=1,N)                                    00004970
  37  FORMAT(1X,13D10.2)                                                00004980
C     UPDATE THE INITIAL INCREMENTAL COST CO                            00004990
      CALL MULCON(1,N,C,U,CO)                                           00005000
 200  CONTINUE                                                          00005010
      RETURN                                                            00005020
      END                                                               00005030
C                                                                       00005040
C     SUBROUTINE FOR EVALUATING THE MEASUREMENT SPECIFICATION ERROR     00005050
C     EM = (YT - H(T)XT - B(T)) FOR TIME T                              00005060
C                                                                       00005070
      SUBROUTINE RME(N,MOBS,YT,XT,H,B,EM)                               00005080
      IMPLICIT REAL*8(A-H,O-Z)                                          00005090
      DIMENSION YT(15,15),XT(15,15),H(15,15),B(15,15),EM(15,15)         00005100
      DIMENSION HX(15,15),HXPB(15,15)                                   00005110
      CALL MUL(MOBS,N,1,H,XT,HX)                                        00005120
      CALL ADD(MOBS,1,HX,B,HXPB)                                        00005130
      CALL SUB(MOBS,1,YT,HXPB,EM)                                       00005140
      RETURN                                                            00005150
      END                                                               00005160
C                                                                       00005170
C     SUBROUTINE FOR EVALUATING THE DYNAMIC SPECIFICATION ERROR         00005180
C     ED = (XTP1 - F(T)XT - A(T)) FOR TIME T                            00005190
C                                                                       00005200
      SUBROUTINE RDE(N,XTP1,XT,F,A,ED)                                  00005210
      IMPLICIT REAL*8(A-H,O-Z)                                          00005220
      DIMENSION XTP1(15,15),XT(15,15),F(15,15),A(15,15),ED(15,15)       00005230
      DIMENSION FXT(15,15),FXTPA(15,15)                                 00005240
      CALL MUL(N,N,1,F,XT,FXT)                                          00005250
      CALL ADD(N,1,FXT,A,FXTPA)                                         00005260
      CALL SUB(N,1,XTP1,FXTPA,ED)                                       00005270
      RETURN                                                            00005280
      END                                                               00005290

Copyright © Leigh Tesfatsion. All Rights Reserved.