/*
** TESTING MEAN BY BOOTSTRAP 
** PREPARED BY S.C. AHN      
** DESIGNED TO EXAMINE 
** FINITE-SAMPLE SIZE PROPERTIES OF MEAN TESTS
*/


@ CLEAN UP MEMORY @

    new ;

@ OPEN OUTPUT FILE @

    output file = bmean.out reset ;

@ DATA GENERATION @

    sd    = 4    ; @ seed number for random number generation @
    tt    = 10   ; @ number of observation @
    nboot = 1000 ; @ number of bootstrap samples @

@ THE VALUE OF POPULATION MEAN TO TEST @

    mu = 5   ; @ H_0: mu = 5 @

@ X IID FROM CHI(KK) @

    kk = 5                        ; @ degrees of freedom @
    xx = sumc(rndns(kk,tt,sd)^2)  ; @ tt*1 data vector ~ chi^2(5) @

@ COMPUTING SAMPLE MEAN AND SAMPLE VARIANCE @

    xb  = meanc(xx)   ; @ sample mean     @
    xs2 = stdc(xx)^2  ; @ sample variance @

@ ASYMPTOTIC CHI_SQUARE WALD TEST FOR H_0: MU = TRUE MU @

    awt  = tt*(xb-mu)^2/xs2  ; @ chi^2 Wald stat.  @
	pawt = cdfchic(awt,1)    ; @ p-value           @


/* BOOTSTRAP WITH PROB = 1/T FOR ALL OBSERVATIONS */

bootxb = zeros(nboot,1) ;
bootwt = zeros(nboot,1) ;

i = 1 ; do while i <= nboot ;

    cboot = ceil( tt*rndus(tt,1,sd) ) ; @ choosing bootstrap data points @
    bsam  = xx[cboot,.] ;

    bxb   = meanc(bsam)  ; @ boot-sample mean     @
    bxs2  = stdc(bsam)^2 ; @ boot-sample variance @

    bwt       = tt*(bxb-xb)^2/bxs2 ; @ boot-chi^2 Wald stat. @
    bootxb[i] = bxb                ; @ saving boot-means     @
    bootwt[i] = bwt                ; @ saving boot-stat      @

i = i + 1 ; endo ;

@ BIAS @

    bias = meanc(bootxb) - xb ;

@ BIAS-CORRECTED ESTIMATE @

    bcxb = xb - bias ;

@ BOOTSTRAP CRITICAL VALUES @

    bootwt = sortc(bootwt,1) ;
    tc     = nboot           ;
    bootc  = bootwt[ceil(tc*0.99)]|bootwt[ceil(tc*0.95)]|bootwt[ceil(tc*0.9)];

/* REPORTING RESULTS */

format /rd 12,4 ;

"POINT ESTIMATES " ;
"" ;
"   SAMPLE MEAN  BIAS-COR. MEAN" ;
xb bcxb;
"";
"ASYMPTOTIC WALD TEST FOR H_0: MU = " mu ;
"";
"    STATISTICS      P-VAL";
awt pawt ;
"" ;
"BOOTSTRAP WALD TEST FOR H_0: MU =  " mu ;
"";
"    STATISTICS    C-VAL-1%     C-VAL-5%     C-VAL-10% ";
awt bootc' ;

output off ;