PROGRAM BT2_8
      PARAMETER (NMAX=100,MMAX=6)
	DIMENSION X(NMAX,MMAX)
	DIMENSION TB(MMAX),DX(MMAX),SX(MMAX),A(MMAX),E(MMAX)
	REAL TRIMEAN(MMAX),TBHC(MMAX),IQR(MMAX),MAD(MMAX)
	REAL SHC(MMAX),YULKED(MMAX) 
	REAL MEDIAN(MMAX)
	DIMENSION R(MMAX,MMAX),RC(MMAX,MMAX)
      DIMENSION XX(NMAX*MMAX)

	OPEN(1,FILE='C:\MASTER\G_TR\BAITAP~1\NEW\BANG8.TXT')
	READ(1,*)
	READ(1,*)
	READ(1,*) N,M
	DO I=1,N
         READ(1,*) TMP, (X(I,J),J=1,M)
	ENDDO
	K=0
	DO J=1,M
	DO I=1,N
	   K=K+1
	   XX(K)=X(I,J)
	ENDDO
	ENDDO
	ALFA=0.10
	CALL STATISTICS(XX,N,M,ALFA,TB,DX,SX,A,E,MEDIAN,
     &            TRIMEAN,TBHC,IQR,MAD,SHC,YULKED)
	
	WRITE(*,*)'    CAC DAC TRUNG THONG KE CO BAN '
	WRITE(*,*)' TT ','   TB ','   DX ','   SX ','    A ',
     &	'    E ','MEDIAN',' TRIM ',' TBHC ','  IQR ',
     &    '  MAD ','  SHC ','YULKED'
	DO J=1,M
	    WRITE(*,'(I4,12F6.2)')J,TB(J),DX(J),SX(J),A(J),E(J),
     &   MEDIAN(J),TRIMEAN(J),TBHC(J),IQR(J),MAD(J),SHC(J),YULKED(J)
	ENDDO

	END

	SUBROUTINE STATISTICS(XX,N,M,ALFA,TB,DX,SX,A,E,MEDIAN,
     &             TRIMEAN,TBHC,IQR,MAD,SHC,YULKED)
! CHUONG TRINH NAY TINH CAC DAC TRUNG THONG	KE CUA HE CAC 
! DAI LUONG NGAU NHIEN MA X LA MA TRAN SO LIEU BAN DAU
! INPUT: + X(N,M) MA TRAN SO LIEU N HANG (N QUAN TRAC), M COT
!          (M BIEN)
!        + N DUNG LUONG MAU
!        + M BIEN
!        + ALFA (%) SO PHAN TRAM CAC THANH PHAN BI CAT BO
! OUTPUT: + TB,DX,SX,A,E LA CAC MANG DO DAI M TUONG UNG CHUA
!           TRUNG BINH, PHUONG SAI, DO LECH CHUAN, DO BAT DOI
!           XUNG VA DO NHON CUA CAC BIEN
!         + TRIMEAN,TBHC,IQR,MAD,SHC,YULKED LA CAC MANG DO DAI
!           M TUONG UNG CHUA TRIMEAN, TRUNG BINH HIEU CHINH,
!           BIEN DO PHAN TU, MAD (MEDIAN ABSOLUTE DEVIATION),
!           PHUONG SAI HIEU CHINH VA CHI SO YULE-KENDALL CUA
!           CAC BIEN
!         + MEDIAN MANG DO DAI M CHUA TRUNG VI
!
	  DIMENSION X(N,M),XX(N*M)
	  INTEGER N,M,I,J,K
	  DIMENSION TB(M),DX(M),SX(M),A(M),E(M)
	  REAL TRIMEAN(M),TBHC(M),IQR(M),MAD(M),SHC(M),YULKED(M)
	  REAL MEDIAN(M)
	  DIMENSION XTMP(N)
	K=0
	DO J=1,M
	  DO I=1,N
	     K=K+1
		 X(I,J)=XX(K)
	  ENDDO
	ENDDO
	  DO 10 J=1,M
	     DO I=1,N
		    XTMP(I)=X(I,J)
		 ENDDO
		 TB(J)=AMMGOC(XTMP,N,1)
		 DX(J)=AMMTAM(XTMP,N,2)
		 SX(J)=SQRT(DX(J))
		  A(J)=AMMTAM(XTMP,N,3)/(DX(J)*SX(J))
		  E(J)=AMMTAM(XTMP,N,4)/(DX(J)*DX(J))-3.0
	     CALL ANOVA(XTMP,N,ALFA,TB(J),MEDIAN(J),TRIMEAN(J),TBHC(J),
     &             DX(J),SX(J),IQR(J),MAD(J),SHC(J),A(J),YULKED(J))
  10  CONTINUE
	  RETURN
	  END

!+++++++++++++++++++++++++++++++++++++++++++++++++
      SUBROUTINE SORT(X,N,X1,CHISO,N1)
!
!  CHUONG TRINH NAY SAP XEP CHUOI THANH CHUOI TRINH TU
!  VA CHUOI XEP HANG
!
! INPUT: + MANG X DO DAI N CHUA CHUOI SO LIEU BAN DAU
!        + N DUNG LUONG MAU
! OUTPUT: + MANG X DO DAI N CHUA CHUOI TRINH TU
!         + MANG X1 DO DAI N1 CHUA CHUOI XEP HANG
!         + MANG CHISO DO DAI N1 CHUA CHI SO CUA CHUOI XEP HANG
!         + N1 SO THANH PHAN TRONG CHUOI XEP HANG
!
	  DIMENSION X(1),X1(1),CHISO(1)
	  INTEGER N,N1,I,J,K
	  REAL TMP
! SAP XEP THANH CHUOI TRINH TU TANG DAN
	  DO I=1,N-1
	  DO J=I+1,N
	     IF (X(J).LT.X(I)) THEN
	        TMP=X(J)
			X(J)=X(I)
			X(I)=TMP
		 ENDIF
	  ENDDO
	  ENDDO
! XAC DINH CAC THANH PHAN CUA CHUOI XEP HANG
	  N1=1
	  I=1
	  X1(N1)=X(I)
	  CHISO(N1)=REAL(I)
  10  I=I+1
      IF (I.GT.N) GOTO 100
	  IF (X(I).NE.X(I-1)) THEN
	     N1=N1+1
		 X1(N1)=X(I)
	  ENDIF
      GOTO 10
 100  CONTINUE
!  TINH CHI SO CUA CAC THANH PHAN TRONG CHUOI XEP HANG
	  DO J=1,N1
         K=0
		 CHISO(J)=0.0
		 DO I=1,N
            IF (X(I).EQ.X1(J)) THEN
			   K=K+1
			   CHISO(J)=CHISO(J)+REAL(I)
			ENDIF
		 ENDDO
		 CHISO(J)=CHISO(J)/REAL(K)
	  ENDDO
      RETURN
	  END

	  FUNCTION Q05(X,N)
!
! HAM NAY TINH TRUNG VI CUA CHUOI X
! INPUT: + X LA MANG SO LIEU DA SAP XEP THANH CHUOI TRINH TU
!        + N LA DO DAI CUA X (DUNG LUONG MAU)
! OUTPUT: + TRUNG VI CUA CHUOI
!
      DIMENSION X(1)
	  INTEGER N,I,NC
	  NC=N/2
	  IF (MOD(N,2).EQ.0) THEN
		 Q05=(X(NC)+X(NC+1))/2.0
	  ELSE
	     Q05=X(NC+1)
	  ENDIF
	  RETURN
	  END

	  SUBROUTINE QUANTILES(X,N,Q50,Q25,Q75)
!
! CHUONG TRINH NAY TINH CAC PHAN VI CHINH CUA CHUOI X DO DAI N
!
! INPUT: + X MANG CHUA CHUOI SO LIEU DA SAP XEP TRINH TU
!        + N DUNG LUONG MAU
! OUTPUT: + TRUNG VI Q50
!         + CAC PHAN VI DUOI VA TREN Q25, Q75
! FUNCTION/SUBROUTINE DUOC GOI TOI: Q05(X,N)
!
      DIMENSION X(1), X1(N)
	  INTEGER N,I,NC1,NC2,NT,J
	  REAL Q50,Q25,Q75,Q125,Q375,Q625,Q875
	  Q50=Q05(X,N)
	  NT=N
	  NC1=NT/2
	  IF (MOD(NT,2).EQ.0) THEN
	     NC2=NC1+1
	  ELSE
	     NC1=NC1+1
		 NC2=NC1
	  ENDIF
	  DO I=1,NC1
	     X1(I)=X(I)
	  ENDDO
	  Q25=Q05(X1,NC1)
	  J=0
	  DO I=NC2,N
	     J=J+1
	     X1(J)=X(I)
	  ENDDO
	  Q75=Q05(X1,NC1)
	  RETURN
	  END

	  FUNCTION AMMGOC(X,N,K)
!
! HAM NAY TINH MOMEN GOC BAC K CUA CHUOI SO LIEU
! INPUT: + X MANG SO LIEU BAN DAU DO DAI N
!        + N DUNG LUONG MAU
!        + K BAC CUA MOMEN
! OUTPUT: MOMEN GOC BAC K
!
      DIMENSION X(N)
	  INTEGER N,K,I
	  REAL TMP
	  TMP=0.0
	  DO I=1,N
	     TMP=TMP+X(I)**K
	  ENDDO
	  AMMGOC=TMP/REAL(N)
	  RETURN
	  END

	  FUNCTION FAC(N)
! TINH N! (N GIAI THUA)
	  REAL TMP
	  IF (N.LT.0) THEN
	     WRITE(*,*)' INVALID NUMERIC INPUT IN FAC FUNCTION'
		 STOP
	  ELSE IF (N.EQ.0.OR.N.EQ.1) THEN
	     FAC=1.0
		 RETURN
	  ELSE
	     TMP=1.0
		 DO I=2,N
	        TMP=TMP*REAL(I)
		 ENDDO
		 FAC=TMP
		 RETURN
	  ENDIF
	  END

      FUNCTION COMBIN(N,K)
! TINH TO HOP CHAP K CUA N
! FUNCTION/SUBROUTINE DUOC GOI TOI: FAC(N)
	  IF (N.LT.0.OR.K.LT.0.OR.N.LT.K) THEN
	     WRITE(*,*)' INVALID NUMERIC INPUT IN COMBIN FUNCTION'
		 STOP
      ELSE 
         COMBIN=FAC(N)/FAC(K)/FAC(N-K)
	     RETURN
	  ENDIF
	  END

	  FUNCTION AMMTAM(X,N,K)
! HAM NAY TINH MOMEN TRUNG TAM BAC K TU CHUOI SO LIEU X
! INPUT: + X MANG SO LIEU BAN DAU DO DAI N
!        + N DUNG LUONG MAU
!        + K BAC CUA MOMEN
! OUTPUT: MOMEN TRUNG TAM BAC K
! FUNCTION/SUBROUTINE DUOC GOI TOI: AMMGOC(X,N,K), COMBIN(N,K)
	  DIMENSION X(N)
	  INTEGER N,K,I
	  REAL TMP,TB
	  TB=AMMGOC(X,N,1)
	  TMP=0
	  DO I=0,K
	     TMP=TMP+(-1)**I*COMBIN(K,I)*TB**I*AMMGOC(X,N,K-I)
	  ENDDO
	  AMMTAM=TMP
	  RETURN
	  END

	  SUBROUTINE ANOVA(XX,N,ALFA,TB,MEDIAN,TRIMEAN,TBHC,
     &                   DX,SX,IQR,MAD,SHC,A,YULKED)
! CHUONG TRINH NAY TINH CAC DAC TRUNG THONG KE DON GIAN UNG DUNG
! TRONG PHAN TICH KHAO SAT SO LIEU
!
! INPUT: + X MANG SO LIEU BAN DAU DO DAI N
!        + N DUNG LUONG MAU
!        + ALFA () PHAN TRAM SO LIEU BI CAT BO
! OUTPUT: + TB TRUNG BINH SO HOC
!         + MEDIAN TRUNG VI
!         + TRIMEAN
!         + TBHC TRUNG BINH HIEU CHINH
!         + DX PHUONG SAI
!         + SX DO LECH CHUAN
!         + IQR CHI SO BIEN DO PHAN TU
!         + MAD CHI SO MAD (MEDIAN CUA DO LECH TUYET DOI SO VOI
!                            MEDIAN BAN DAU CUA CHUOI)
!         + SHC PHUONG SAI HIEU CHINH
!         + A DO BAT DOI XUNG
!         + YULKED CHI SO YULE-KENDALL
! FUNCTION/SUBROUTINE DUOC GOI TOI: SORT, QUANTILES
!
	  DIMENSION XX(N),X(N),X1(N),CHISO(N)
	  REAL ALFA,TB,MEDIAN,TRIMEAN,TBHC
	  REAL DX,SX,IQR,MAD,SHC,A,YULKED
	  REAL Q50,Q25,Q75,TMP
	  INTEGER N,I,N1,K
!
	  X=XX
	  TB=AMMGOC(X,N,1)
	  DX=AMMTAM(X,N,2)
	  SX=SQRT(DX)
	  A =AMMTAM(X,N,3)/(DX*SX)

	  CALL SORT(X,N,X1,CHISO,N1)
	  CALL QUANTILES(X,N,Q50,Q25,Q75)
!
	  MEDIAN=Q50
	  TRIMEAN=(Q25+2.0*Q50+Q75)/4.0
	  IQR=Q75-Q25
	  YULKED=(Q25-2.0*Q50+Q75)/IQR
!
	  K=INT(REAL(N)*ALFA)
	  TMP=0.0
	  DO I=K+1,N-K
	     TMP=TMP+X(I)
	  ENDDO
	  TBHC=TMP/REAL(N-2*K)
!
	  TMP=0.0
	  DO I=K+1,N-K
	     TMP=TMP+(X(I)-TBHC)**2
	  ENDDO
	  SHC=TMP/REAL(N-2*K)
!
	  DO I=1,N
	     X(I)=ABS(X(I)-Q50)
	  ENDDO
	  CALL QUANTILES(X,N,Q50,Q25,Q75)
	  MAD=Q50
	  RETURN
	  END