PROGRAM NMODEL ****************************************************************************** * * * ** ** ****** ** ** ******* *** ** * * *** ** ** ** ** ** ** ** ** ** * * ** ** ** ****** ** ** ** *** *** ** ** * * ** *** ** ** ** ** ** ********** ** * * ** *** ****** *** ** *** ** ** ******** * * * * * ******* SIMPLE NEURAL NETWORK (BACK-PROPAGATION) FORECASTING PROGRAM ******* * * * Geraint Johnes * * * * * * O --- O * * X > O * * O --- O * * * * * * * * This program is designed to produce forecasts by way of a simple * * back-propagation algorithm. A single hidden layer feedforward neural * * network provides the framework. It forecasts values for 4 differenced * * variables using information about (4-period) lagged values of those * * same (differenced) variables, and about the (1 period) lagged levels * * of these four variables. It is therefore a nonlinear error correction * * process with a vector of variables. * * * * The program enters data on the following: 12th difference in log * * industrial production; 12th difference in log prices; log unemployment * * rate; log interest rate. You retain some observations at the end of * * your series in order to evaluate forecasting performance of the * * neural network. You set the forecast length (FCLONG), and * * the program output gives data on rate of change of industrial * * production, rate of change of prices, unemployment rate and interest * * rate (all in per cent) which are forecast FCLONG periods ahead by * * the neural network. The number of such forecasts is equal to the * * number of retained observations minus FCLONG. * * * * You may use this program (or simple variants) freely so long as you * * both (a) acknowledge its source and (b) absolve its author of any * * responsibilities which might arise. Although I believe this program * * to be bug-free, I can offer no guarantees. * * * ****************************************************************************** * DOUBLE PRECISION DATA1(1000,50),INP(1000,20),DESOUT(1000,20), 1TOTERR(20),WEIGHT(20),W(20,20),SIGNAL(20),WEIGHT1(20),V(20,20), 2OUTRES(20),ERR1(20),INCERR(20),FINERR(20),COM(20),CON(20), 3TT1(20),TTT(20),ERRSQ(1000,20),WDELTA1(20,20),WDELTA2(20,20), 4TT11(20),outres1(20),DATA2(1000,50),data3(1000,50),zz1(4,1000), 5zz2(4,4),yak(1000,4),pc(1000),vvv(4,4),ZZ4(1000,20), 6data4(1000,20),yak1(1000,4),yaki(1000,4) DOUBLE PRECISION BETA,ERRMAX1,errmax2,errmax3,errmax4, 1ALPHA,SEED,pcc,pccc,errmax INTEGER HIDDEN,INPUT,OUTPUT,NPASS,NOBS,NPATTS,tt,I,J,K,II, 1nnobs,nobstot,nmax,ifail,ia,ivr,ivi,iaa,ig,ik,jg,jk,ii3,j33 parameter (nmax=4,ia=nmax,ivr=nmax,ivi=nmax) double precision a(ia,nmax),p(nmax),ri(nmax),rr(nmax), 1vi(ivi,nmax),vr(ivr,nmax) integer intger(nmax) double precision g05caf EXTERNAL G05CAF EXTERNAL G05CBF external f02agf external f01aaf CALL G05CBF(0) OPEN(5,FILE='neumodel.dat') OPEN(6,FILE='nmodel.out') * * set number of neurodes in hidden layer * HIDDEN=10 * * set number of neurodes in input layer * INPUT=20 * * set number of neurodes in output layer * OUTPUT=4 * * set value of the learning constant * BETA=.1 * * set maximum acceptable error * errmax=0.0001 errmax1=0.00005 errmax2=0.00005 errmax3=0.0005 errmax4=0.0005 * * set maximum number of passes * NPASS=10000 * * set length of forecast required * FCLONG=6 * * set momentum * ALPHA=.3 * * set number of observations to be used in training * NOBS=408 * * set number of observations in total (including those withheld for * forecast evaluation) * nobstot=436 * * choose criterion for algorithm exit (1 => errmax is maximum value * of RMSE, 0 => errmax is maximum value for the worst fit during last pass; * 2 => errmax is maximum value of RMSE variable by variable) * EXIT=2 * * In this case the inputs comprise 4 lags of each of the 4 endogenous * variables in differenced form plus lagged levels of the 4 endogenous * variables (to allow, in an unrestricted manner, for cointegration). * The next few lines read in the data and sort them into patterns. Note * that all variables must be transformed to ensure that they lie * within the unit interval. The data read in at this stage should include * all patterns to be used in training, but should not include * out-of-sample patterns which are retained by the forecaster to assess * forecasting performance * pccc=10**10 DO 1 I=1,OUTPUT DO 1 J=1,HIDDEN WDELTA1(I,J)=0 1 WDELTA2(I,J)=0 * * read in the basic data. Column 1 is the 12th difference of the logged * industrial production index. Column 2 is the 12th difference of the * logged consumer price index. Column 3 is the logged rate of unemployment * and column 4 is the logged interest rate * read(5,9980)((data1(i,j),j=1,4),i=1,nobstot) * write(6,9980)((data1(i,j),j=1,4),i=1,nobstot) * * transform data so that all observations lie in the unit interval * NPATTS=NOBS-5 DO 77 I=2,NOBStot DO 77 J=1,OUTPUT data3(i,j)=data1(i,j) 77 DATA1(I,J+4)=DATA1(I,J)-DATA1(I-1,J) DO 78 I=3,NOBStot DO 78 J=1,OUTPUT 78 DATA1(I,J+8)=DATA1(I-1,J+4) DO 79 I=4,NOBStot DO 79 J=1,OUTPUT 79 DATA1(I,J+12)=DATA1(I-1,J+8) DO 76 I=5,NOBStot DO 76 J=1,OUTPUT 76 DATA1(I,J+16)=DATA1(I-1,J+12) DO 75 I=6,NOBStot DO 75 J=1,OUTPUT 75 DATA1(I,J+20)=DATA1(I-1,J+16) DO 74 I=2,NOBStot DO 74 J=1,OUTPUT 74 DATA1(I,J+24)=DATA1(I-1,J) do 16 i=1,nobstot do 16 j=1,28 16 data1(i,j)=exp(data1(i,j))/(1+exp(data1(i,j))) DO 72 I=1,NOBSTOT DO 72 J=1,28 72 DATA2(I,J)=DATA1(I,J) DO 71 I=NOBS+1,NOBSTOT DO 71 J=1,28 71 DATA1(I,J)=0 * DO 6 I=1,NPATTS DO 6 J=1,INPUT 6 INP(I,J)=DATA1(I+5,J+8) do 608 i=1,npatts DO 608 K=1,OUTPUT 608 DESOUT(I,K)=DATA1(I+5,K+4) 7777 format(10f8.3) * **** train network **** train network **** train network **** * DO 14 J=1,20 DO 14 K=1,20 * * initialise weights at random * SEED=G05CAF(SEED) V(J,K)=2*(SEED-.5) SEED=G05CAF(SEED) 14 W(J,K)=2*(SEED-.5) * DO 10 I=1,NPASS DO 11 II3=1,NPATTS DO 11 J33=1,OUTPUT 11 ERRSQ(II3,J33)=0 DO 15 J=1,OUTPUT TT1(J)=0 15 TOTERR(J)=0 TT=0 DO 59 J=1,OUTPUT 59 TT11(J)=0 DO 12 II=5,NPATTS * * compute signal from each hidden layer neurode * DO 20 J=1,HIDDEN WEIGHT(J)=0 DO 30 K=1,INPUT 30 WEIGHT(J)=WEIGHT(J)+W(J,K)*INP(II,K) 20 SIGNAL(J)=1/(1+EXP(-WEIGHT(J))) * * compute signal from each output layer neurode * DO 40 J=1,OUTPUT WEIGHT1(J)=0 DO 50 K=1,HIDDEN 50 WEIGHT1(J)=WEIGHT1(J)+V(J,K)*SIGNAL(K) 40 OUTRES(J)=1/(1+EXP(-WEIGHT1(J))) * * compute error * DO 60 J=1,OUTPUT ERR1(J)=DESOUT(II,J)-OUTRES(J) IF(ABS(ERR1(J)).GT.TT11(J))TT11(J)=ABS(ERR1(J)) 60 ERRSQ(II,J)=ERR1(J)**2+ERRSQ(II,J) * * if(i.eq.277)then * if(i.eq.1)write(6,6441) * write(6,1212)(errsq(ii,j),j=1,output) * ENDIF * * * backpropagate errors from output layer to the hidden layer * DO 70 J=1,HIDDEN INCERR(J)=0 DO 80 K=1,OUTPUT 80 INCERR(J)=INCERR(J)+ERR1(K)*V(K,J) 70 FINERR(J)=INCERR(J)*SIGNAL(J)*(1-SIGNAL(J)) 1218 FORMAT(2F20.12) * * adjust weights used by output layer * DO 90 J=1,OUTPUT DO 90 K=1,HIDDEN WDELTA1(J,K)=BETA*ERR1(J)*SIGNAL(K)+ALPHA*WDELTA1(J,K) 90 V(J,K)=V(J,K)+WDELTA1(J,K) * * adjust weights used by hidden layer * DO 100 J=1,HIDDEN DO 100 K=1,INPUT WDELTA2(J,K)=BETA*FINERR(J)*INP(II,K)+ALPHA*WDELTA2(J,K) 100 W(J,K)=W(J,K)+WDELTA2(J,K) 12 CONTINUE 1212 format(4f15.10) 1213 format(f10.4) * * evaluate whether error is sufficiently small to exit training algorithm * 7236 format(7h aloha ) do 1133 ig=1,npatts do 1133 jg=1,output 1133 ZZ4(ig,jg)=ERRSQ(ig,jg) 7776 format(7h hello ) 8888 format(1x,a,i5) 7775 format(4f10.4) 1138 format(f10.4) 207 continue IF(EXIT.ne.1)GOTO108 DO 106 J=1,OUTPUT TTT(J)=ERRSQ(NPATTS,J) 106 TTT(J)=(TTT(J)/NPATTS)**0.5 TT=0 write(6,1212)(ttt(j),j=1,4) DO 107 J=1,OUTPUT 107 IF(TTT(J).LE.ERRMAX)TT=TT+1 IF(TT.EQ.OUTPUT)GOTO110 108 continue if(exit.ne.0)goto9444 TT=0 DO 111 J=1,OUTPUT 111 IF(TT11(J).GT.ERRMAX)TT=TT+1 IF(TT.EQ.0)GOTO110 9444 if(errsq(npatts,1).le.errmax1.and. 1errsq(npatts,2).le.errmax2.and. 2errsq(npatts,3).le.errmax3.and. 3errsq(npatts,4).le.errmax4)goto110 1010 do 1139 ig=1,npatts do 1139 jg=1,output 1139 errsq(ig,jg)=0 10 CONTINUE 110 CONTINUE * **** end training **** end training **** end training **** * * display weights resulting from trained network * WRITE(6,9990)I WRITE(6,9960) WRITE(6,9970)((W(J,K),K=1,INPUT),J=1,HIDDEN) WRITE(6,9960) WRITE(6,9970)((V(J,K),K=1,HIDDEN),J=1,OUTPUT) * * use trained network to forecast * do 210 nnobs=nobs,nobstot-fclong * do 210 nnobs=nobs+fclong,nobstot * do 432 ig=1,nobstot do 432 jg=1,output 432 yaki(ig,jg)=0 * DO 140 I=NnobS,NnobS+FCLONG do 141 ii=1,i do 141 jj=1,28 141 data1(ii,jj)=data2(ii,jj) * * put forecasts into data1 where outturns not available * if(i.le.nnobs)then goto181 else do 143 jjj=1,4 data1(i,jjj+4)=outres(jjj) data1(i,jjj+8)=data1(i-1,jjj+4) data1(i,jjj+12)=data1(i-1,jjj+8) data1(i,jjj+16)=data1(i-1,jjj+12) data1(i,jjj+20)=data1(i-1,jjj+16) data1(i,jjj+24)= 1exp(log(data1(i-1,jjj+24)/(1-data1(i-1,jjj+24)))+ 2log(data1(i-1,jjj+4)/(1-data1(i-1,jjj+4))))/(1+exp( 3log(data1(i-1,jjj+24)/(1-data1(i-1,jjj+24)))+ 4log(data1(i-1,jjj+4)/(1-data1(i-1,jjj+4))))) 143 continue endif 181 continue do 146 jjj=1,input inp(i,jjj)=data1(i,jjj+8) 146 continue DO 155 K1=1,HIDDEN 155 CON(K1)=0 DO 160 K1=1,HIDDEN DO 160 K2=1,INPUT 5555 format(f12.6) 160 CON(K1)=CON(K1)+W(K1,K2)*INP(I,K2) DO 170 K1=1,HIDDEN 170 SIGNAL(K1)=1/(1+EXP(-CON(K1))) DO 175 K1=1,OUTPUT 175 COM(K1)=0 DO 180 K1=1,OUTPUT DO 180 K2=1,HIDDEN 180 COM(K1)=COM(K1)+V(K1,K2)*SIGNAL(K2) DO 190 K1=1,OUTPUT OUTRES(K1)=1/(1+EXP(-COM(K1))) yaki(i,k1)=yaki(i,k1)+log(outres(k1)/(1-outres(k1))) if(i.eq.nnobs+fclong)then yak(i,k1)=yaki(i,k1) * write(6,1212)yak(i,k1) endif 190 continue 140 CONTINUE ********************************************************************** * * display forecasts * do 611 ig=1,nobstot do 611 jg=1,output data4(ig,jg+4)=log(data2(ig,jg+4)/(1-data2(ig,jg+4))) 611 continue 210 continue do 246 nnobs=nobs,nobstot-fclong do 248 ig=nnobs+fclong,nobstot do 247 jg=1,2 247 yak1(ig,jg)=100*data3(ig-fclong,jg) do 248 jg=3,4 248 yak1(ig,jg)=exp(data3(ig-fclong,jg)) ***************************************************************** do 243 ig=nnobs+fclong,nnobs+fclong do 243 k1=1,2 yak1(ig,k1)=yak1(ig,k1)+100*yak(ig,k1) 243 continue do 244 ig=nnobs+fclong,nnobs+fclong do 244 k1=3,4 1313 format(2f20.10) 244 yak1(ig,k1)=(exp(log(yak1(ig,k1))+yak(ig,k1))) ***************************************************************** 246 continue write(6,9950) write(6,9951) * write(6,6332) * do 211 nnobs=nobs,nobstot-fclong * 211 if(nnobs.gt.nobs) * 1write(6,9870)(data4(nnobs,j+4), * 2yak(nnobs+fclong,j),j=1,output) write(6,6333) * do 213 nnobs=nobs,nobstot-fclong do 213 j=3,4 213 data3(nnobs+fclong,j)=exp(data3(nnobs+fclong,j)) do 214 nnobs=nobs,nobstot-fclong do 214 j=1,2 214 data3(nnobs+fclong,j)=100*(data3(nnobs+fclong,j)) do 212 nnobs=nobs,nobstot-fclong 212 if(nnobs.gt.nobs) 1write(6,9870)(data3(nnobs+fclong,j),yak1(nnobs+fclong,j), 2j=1,output) write(6,6334) write(6,6335) write(6,6336) write(6,6337) write(6,6338) 9990 FORMAT(17H NUMBER OF PASSES ,I10/) 9980 FORMAT(4F15.6) 9970 FORMAT(2F10.4) 9870 format(8f9.5) 9960 FORMAT(/8H WEIGHTS ) 9950 FORMAT(/10H FORECASTS ) 9951 format(49h In pairs of columns, actual first then forecast ) 6332 format(/45h FORECASTS OF DIFFERENCED (LOGGED) VARIABLES ) 6339 format(/25h WITHIN-SAMPLE FORECASTS ) 6333 format(/41h FORECASTS OF LEVEL (UNLOGGED) VARIABLES ) 6334 format(/50h These 4 pairs of columns represent, respectively ) 6335 format(/45h 12th difference of logged indl prodn, x 100 ) 6336 format(38h 12th difference of logged CPI, x 100 ) 6337 format(19h unemployment rate ) 6338 format(15h interest rate ) 6441 format(/16h SQUARED ERRORS ) STOP END