'********************************************************** '* Calculate the mean and first third moments of a set of * '* data Y(i) * '* ------------------------------------------------------ * '* Description: * '* * '* Input data: * '* ndata: Integer size of data set * '* Y: Real Vector of ndata ordinates * '* Output data: * '* S1: ndata-Mean of ndata values * '* S2: Sum (Y(i)-S1)^2 for i=1..ndata * '* S3: Sum (Y(i)-S1)^3 for i=1..ndata * '* S4: Sum (Y(i)-S1)^4 for i=1..ndata * '* ------------------------------------------------------ * '* SAMPLE RUN: * '* Number of data: 5 * '* 1: 12 * '* 2: 9 * '* 3: 7 * '* 4: 15 * '* 5: 6 * '* * '* S1= 9.8 * '* S2= 54.8 * '* S3= 73.92001 * '* S4= 1024.976 * '* * '* Error code: 0 * '* * '* ------------------------------------------------------ * '* Ref.: "Journal of Applied Statistics (1972) vol.21, * '* page 226". * '* * '* Basic Release By J-P Moreau, Paris. * '* (www.jpmoreau.fr) * '********************************************************** 'PROGRAM TMOMNTS; DEFINT I-N NMAX = 25 CLS PRINT INPUT " Number of data: ", ndata DIM Y(ndata) FOR i = 1 TO ndata PRINT " "; i; ": "; : INPUT "", Y(i) NEXT i FOR i = 1 TO ndata X = Y(i): K = 4: N = i - 1: GOSUB 1000 'call MOMNTS(X,K,N,S1,S2,S3,S4,IFAULT) NEXT i PRINT PRINT " S1= "; S1 PRINT " S2= "; S2 PRINT " S3= "; S3 PRINT " S4= "; S4 PRINT PRINT " Error code: "; IFAULT PRINT END 1000 'Subroutine MOMNTS(X, K, N, S1, S2, S3, S4, IFAULT) ' ------------------------------------------------------------------ ' ALGORITHM AS 52 APPL. STATIST. (1972) VOL.21, P.226 ' ADDS A NEW VALUE, X, WHEN CALCULATING A MEAN, AND SUMS ' OF POWERS OF DEVIATIONS. N IS THE CURRENT NUMBER OF ' OBSERVATIONS, AND MUST BE SET TO ZERO BEFORE FIRST ENTRY ' ------------------------------------------------------------------ ' AN, AN1, DX, DX2, ZERO, ONE, TWO, THREE, FOUR, SIX: REAL ZERO = 0!: ONE = 1!: TWO = 2!: THREE = 3!: FOUR = 4!: SIX = 6! IF K > 0 AND K < 5 AND N >= 0 THEN GOTO 10 IFAULT = 1 RETURN 10 IFAULT = 0 N = N + 1 IF N > 1 THEN GOTO 20 ' FIRST ENTRY, SO INITIALISE S1 = X S2 = ZERO S3 = ZERO S4 = ZERO RETURN ' SUBSEQUENT ENTRY, SO UPDATE 20 AN = N AN1 = AN - ONE DX = (X - S1) / AN DX2 = DX * DX IF K = 1 THEN GOTO 60 'calculate only mean S1 ELSEIF K = 2 THEN GOTO 50 'calculate S1 and S2 ELSEIF K = 3 THEN GOTO 40 'calculate S1, S2 and S3 ELSEIF K = 4 THEN GOTO 30 'calculate S1..S4 END IF 30 S4 = S4 - DX * (FOUR * S3 - DX * (SIX * S2 + AN1 * (ONE + AN1 * AN1 * AN1) * DX2)) 40 S3 = S3 - DX * (THREE * S2 - AN * AN1 * (AN - TWO) * DX2) 50 S2 = S2 + AN * AN1 * DX2 60 S1 = S1 + DX RETURN 'end of file momnts.bas