% Function: This is the program that calculates the alpha cutoff to  
% determine consistent hits. To do this we use a formula to calculate the 
% p-value required for a desired consistency, which is dependent on the
% number of experiments. This value is then adjusted based on the empirical
% false positive rate determined by comparing experimental replicates.

% Note this function is currently only compatible for experiments with 8 or
% less experimental replicates. 

% Inputs:
% MHpvalue - this is the multiple hypothesis corrected p-value cutoff being 
%    used, usually 0.05 divided by the number of IP/Probe combinations 
%    (i.e. bonferroni post hoc correction)
% RepPValueCell - this is a cell array, where each entry contains the cell 
%    array with the replicate comparison p-values
% NRepValues - the number of replicates across the entire RepPValueCell,
%    entered as a numerical entry 

% Outputs:
% EmpiricalAlphaCell - a cell array that has the p-value cutoff for each of
%    the consistency cutoffs
% PerAgreeCell - a cell array that matches the EmpiricalAlphaCell and lists
%    the number of experiments needed to have a given consistency level
% RepList4Fit - a list of all the p-values between the replicates, which is
%    used by the program to perform the empirical false positive rate
%    correction
% BadProbeCell - a list of probes wells that were eliminated in the false 
%    positive rate correction process because of obvious technical errors,
%    defined when the median error in a well is less than MHpvalue (if
%    filtering is selected, default is hard coded as yes)

function [EmpiricalAlphaCell,PerAgreeCell,RepList4Fit,BadProbeCell] = ...
    EmpAdjPValueList(MHpvalue,RepPValueCell,NRepValues)

BinomCoefCell = {[1,0,0;0,0,1];...
    [2,0,0;0,0,2;1,1,0;1,0,1;0,1,1];...
    [3,0,0;0,0,3;2,1,0;2,0,1;0,1,2;1,0,2;1,2,0;1,1,1;0,2,1];...
    [4,0,0;0,0,4;3,1,0;3,0,1;0,1,3;1,0,3;2,2,0;2,1,1;2,0,2;0,2,2;1,1,2;...
    1,3,0;1,2,1;0,3,1];...
    [5,0,0;0,0,5;4,1,0;4,0,1;0,1,4;1,0,4;3,2,0;3,1,1;3,0,2;0,2,3;1,1,3;...
    2,0,3;2,3,0;2,2,1;2,1,2;1,2,2;0,3,2;1,4,0;1,3,1;0,4,1];...
    [6,0,0;0,0,6;5,1,0;5,0,1;0,1,5;1,0,5;4,2,0;4,1,1;4,0,2;0,2,4;1,1,4;...
    2,0,4;3,3,0;3,2,1;3,1,2;3,0,3;0,3,3;1,2,3;2,1,3;2,4,0;2,3,1;2,2,2;...
    0,4,2;1,3,2;1,5,0;1,4,1;0,5,1];...
    [7,0,0;0,0,7;6,1,0;6,0,1;0,1,6;1,0,6;5,2,0;5,1,1;5,0,2;0,2,5;1,1,5;...
    2,0,5;4,3,0;4,2,1;4,1,2;4,0,3;0,3,4;1,2,4;2,1,4;3,0,4;3,4,0;3,3,1;...
    3,2,2;3,1,3;0,4,3;1,3,3;2,2,3;2,5,0;2,4,1;2,3,2;0,5,2;1,4,2;1,6,0;...
    1,5,1;0,6,1];...
    [8,0,0;0,0,8;7,1,0;7,0,1;0,1,7;1,0,7;6,2,0;6,1,1;6,0,2;0,2,6;1,1,6;...
    2,0,6;5,3,0;5,2,1;5,1,2;5,0,3;0,3,5;1,2,5;2,1,5;3,0,5;4,4,0;4,3,1;...
    4,2,2;4,1,3;4,0,4;1,3,4;2,2,4;3,1,4;3,5,0;3,4,1;3,3,2;3,2,3;0,5,3;...
    1,4,3;2,6,0;2,5,1;2,4,2;0,6,2;1,5,2;1,7,0;1,6,1;0,7,1]};

PerAgreeCell = {{1},{2,[2,1]},{3,[3,2],[3,2,1]},...
    {4,[4,3],[4,3,2],[4,3,2,1]},{5,[5,4],[5,4,3],[5,4,3,2],[5,4,3,2,1]},...
    {6,[6,5],[6,5,4],[6,5,4,3],[6,5,4,3,2],[6,5,4,3,2,1]},{7,[7,6],...
    [7,6,5],[7,6,5,4],[7,6,5,4,3],[7,6,5,4,3,2],[7,6,5,4,3,2,1]},...
    {8,[8,7],[8,7,6],[8,7,6,5],[8,7,6,5,4],[8,7,6,5,4,3],...
    [8,7,6,5,4,3,2],[8,7,6,5,4,3,2,1]}};

POutputCell = cell(length(RepPValueCell),1);

for alpha = 1:length(RepPValueCell)
    
    POutputCell{alpha} = BinomPValue(BinomCoefCell{alpha}, MHpvalue);
    %Function - Folder/Statistics/MetaConsistency
    
end


EmpiricalAlphaCell = cell(length(POutputCell),1);

for i = 1:length(POutputCell)
    
    for j = 1:length(POutputCell{i})
        
        EmpPOutputCellj = zeros(length(POutputCell{i}{j}),1);
        
        for k = 1:length(POutputCell{i}{j})
        
            [EmpPOutputCellj(k), RepList4Fit,BadProbeCell] = ....
                EmpiricalAlpha(RepPValueCell,POutputCell{i}{j}(k),...
                NRepValues,'Y'); 
            %Function - Folder/Statistics/MetaConsistency
        end
        
        EmpiricalAlphaCell{i}{j} = EmpPOutputCellj;
        
    end
end

