C @(#)knn.for 17.1.1.1 (ES0-DMD) 01/25/02 17:16:38 C=========================================================================== C Copyright (C) 1995 European Southern Observatory (ESO) C C This program is free software; you can redistribute it and/or C modify it under the terms of the GNU General Public License as C published by the Free Software Foundation; either version 2 of C the License, or (at your option) any later version. C C This program is distributed in the hope that it will be useful, C but WITHOUT ANY WARRANTY; without even the implied warranty of C MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the C GNU General Public License for more details. C C You should have received a copy of the GNU General Public C License along with this program; if not, write to the Free C Software Foundation, Inc., 675 Massachusetss Ave, Cambridge, C MA 02139, USA. C C Corresponding concerning ESO-MIDAS should be addressed as follows: C Internet e-mail: midas@eso.org C Postal address: European Southern Observatory C Data Management Division C Karl-Schwarzschild-Strasse 2 C D 85748 Garching bei Muenchen C GERMANY C=========================================================================== C C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++C C C C Carry out a K-NN DISCRIMINANT ANALYSIS, C C with two groups defined in a training set, and with C C assignment of members of a test set. C C C C Parameters: C C C C TRAIN(N,M) training set, where first N1 rows relate to the C C first group, and the next N2 rows to the second C C group. Must have N1 + N2 = N. C C TEST(N3,M) test set; C C K number of nearest neighbours to consider; C C KLIST(K), DK(K), KPOP(K) are used for storing the K NNs, C C their distances to the object under consider- C C ation, and the group to which the NNs belong. C C C C NOTE: this routine is written in VAX-11 Fortran and not in C C Fortran-77. C C It also contains MIDAS STTPUT calls. C C C C F. Murtagh, ESA/ESO/STECF, Garching. February 1986. C C C C-----------------------------------------------------------------C SUBROUTINE KNN(TRAIN,N1,N2,N,TEST,N3,M,K,KLIST,DK,KPOP,IERR) REAL TRAIN(N,M),TEST(N3,M),DK(K),D INTEGER KLIST(K),KPOP(K) INTEGER I,IX,N1,N2,N3,N,M,K,IERR,ISTAT INTEGER IND,I2,IK,IK2,NUM1,NUM2 CHARACTER*80 LINE C WRITE(LINE,20) CALL STTPUT(LINE,ISTAT) 20 FORMAT(' Object --> group with probability') WRITE(LINE,30) CALL STTPUT(LINE,ISTAT) 30 FORMAT(' ') C DO I = 1, N3 DO IX = 1, K KLIST(IX) = 0 DK(IX) = 1.E+15 KPOP(IX) = 0 ENDDO IND = 0 DO I2 = 1, N1 CALL DISTK(I2,TRAIN,I,TEST,N,N3,M,D) IF (IND.GT.0) THEN DO IK = 1, K IF (D.LT.DK(IK)) THEN IF (IK.LT.K) THEN DO IK2 = K, IK+1, -1 DK(IK2) = DK(IK2-1) KLIST(IK2) = KLIST(IK2-1) KPOP(IK2) = KPOP(IK2-1) ENDDO ENDIF DK(IK) = D KLIST(IK) = I2 KPOP(IK) = 1 GOTO 100 ENDIF ENDDO ELSE IND = IND + 1 DK(IND) = D KLIST(IND) = I2 KPOP(IND) = 1 ENDIF 100 CONTINUE ENDDO C C C DO I2 = N1+1, N CALL DISTK(I2,TRAIN,I,TEST,N,N3,M,D) IF (IND.GT.0) THEN DO IK = 1, K IF (D.LT.DK(IK)) THEN IF (IK.LT.K) THEN DO IK2 = K, IK+1, -1 DK(IK2) = DK(IK2-1) KLIST(IK2) = KLIST(IK2-1) KPOP(IK2) = KPOP(IK2-1) ENDDO ENDIF DK(IK) = D KLIST(IK) = I2 KPOP(IK) = 2 GOTO 200 ENDIF ENDDO ELSE IND = IND + 1 DK(IND) = D KLIST(IND) = I2 KPOP(IND) = 2 ENDIF 200 CONTINUE ENDDO C NUM1 = 0 NUM2 = 0 DO IX = 1, K IF (KPOP(IX).EQ.1) NUM1 = NUM1 + 1 IF (KPOP(IX).EQ.2) NUM2 = NUM2 + 1 ENDDO C (Error check:) IF ((NUM1+NUM2).NE.K) THEN LINE(1:) = ' Total assignments to gp. 1 and to gp. 2 ' X //'dont equal K; is K too large ? ' CALL STTPUT(LINE,ISTAT) IERR = 1 RETURN ENDIF IF (NUM1.GT.NUM2) THEN WRITE (LINE,500) I,FLOAT(NUM1)*100./FLOAT(K) CALL STTPUT(LINE,ISTAT) ENDIF IF (NUM2.GT.NUM1) THEN WRITE (LINE,525) I,FLOAT(NUM2)*100./FLOAT(K) CALL STTPUT(LINE,ISTAT) ENDIF IF (NUM1.EQ.NUM2) THEN WRITE (LINE,550) I,FLOAT(NUM1)*100./FLOAT(K) CALL STTPUT(LINE,ISTAT) ENDIF 500 FORMAT(I6,6X,' 1 ',F8.2,'%') 525 FORMAT(I6,6X,' 2 ',F8.2,'%') 550 FORMAT(I6,6X,' 1 or 2 ',F8.2,'% (equiprobable).') ENDDO C RETURN END C SUBROUTINE DISTK(I,ARR1,J,ARR2,NR1,NR2,NC,D) INTEGER NC,NR1,NR2,K,J,I REAL ARR1(NR1,NC),ARR2(NR2,NC),D C D = 0.0 DO K = 1, NC D = D + (ARR1(I,K)-ARR2(J,K))**2 ENDDO C RETURN END