options nocenter;
%macro fallback(dataset);
/*
Computes adjusted p-values for the fallback procedure.

Written by Alex Dmitrienko. 

Inputs:
 
DATASET = Data set with raw p-values and weights.
*/
proc iml;    
    * Compute the Bonferroni p-value for an intersection hypothesis;
    start bonf(p,u);
        bonfp=min(p[loc(u)]/u[loc(u)]);
        return(bonfp);
    finish bonf;
    use &dataset;
    read all var {weight raw_p} into data;        
    data=t(data);
    w=data[1,];
    p=data[2,];    
    nhyps=ncol(p); 
    nints=2**nhyps-1;
    h=j(nints,nhyps,0);
    v=j(nints,nhyps,0);
    hyp=j(nints,nhyps,0); 
    adjp=j(nhyps,1,0);
    window=j(nhyps,nhyps,0); 
    do i=1 to nhyps;
        do j=0 to nints-1;
            k=floor(j/2**(nhyps-i));
            if k/2=floor(k/2) then h[j+1,i]=1;
        end;
    end;        
    do i=1 to nhyps;
        do j=1 to nhyps;            
            if i>=j then window[i,j]=1;
        end;
    end;            
    do i=1 to nints;
        do j=1 to nhyps;
            if h[i,j]=1 then v[i,j]=sum(w#window[j,])-sum(v[i,]); else v[i,j]=0;
        end;
        hyp[i,]=h[i,]*bonf(p,v[i,]);
    end;
    do i=1 to nhyps; 
        adjp[i]=max(hyp[,i]); 
    end;        
    create adjp from adjp[colname={adjp}];
    append from adjp;        
    quit;
data fallback;
    merge &dataset adjp;
proc print data=fallback noobs;
    title "Fallback procedure: Adjusted p-values";
    format raw_p adjp 6.4;
    run;
%mend fallback;
* Example (problem with four hypotheses);
data scenario;
    input weight raw_p;
    datalines;
    0.25 0.0391
    0.25 0.0242
    0.25 0.0225
    0.25 0.0103
    ;  
%fallback(scenario);