% Author: Elsa W. Birch, ewbirch@stanford.edu, Covert Lab - Stanford Bioengineering
% Date: 4/24/2011
% Version: associated with production of the figures for Birch et al. 2012 
%
% - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 
%
% Summary: Parses a regulatory rule to check taht the dependancies exist and is properly formatted for error free evaluation
% 
%
% Inputs: rawRegulationString - the definition of the rule in the flat file
%		  elementName - regulation rule associated with 
%		  organismModelStruct - organism structure being parsed in 
%
% Outputs: 	regulationRulesEvalString - string of rule converted to a form including properly indexed dependancies
%			regulationDependancies - the referenced gnees/tfs/concentrationsstiuli etc
%			regulationErrors - erros to display, all collected at once for more convivnet model fixing
%
%
% - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 
%
% Notes: 
%
%
% - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

function [regulationRulesEvalString, regulationDependancies, regulationErrors ] = ...
    FBA_EB_regulationRuleParse( rawRegulationString, elementName, organismModelStruct )

% initialize regulation evaluation statement
regulationRulesEvalString = '';
regulationDependancies = {};
regulationErrors = {};

% create dummy logical arrays to use in testing rules for errors later
geneProductActivity = true(organismModelStruct.nGene, 1 );
geneActivity = true(organismModelStruct.nGene, 1 );
tfActivity = true(organismModelStruct.nTFs, 1 );
currentMedia = true(organismModelStruct.nMediaExchange, 1 );
currentStimuli = organismModelStruct.stimuliValuesDefault; %  
currentFlux = ones(organismModelStruct.nRxns, 1);

% initialize flag for continuation of rule parsing
continueRuleParse = 1;


% Basic reasons rule parsing would be unnecesary: ovbious errors or no rule is being
% specified.
% % . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
 
% % test for regulation empty, only whitespace characters, or unexpected characters
if isempty(rawRegulationString) || isempty( regexp(rawRegulationString, '\S','once' ) )

    % set flag to skip other rule parsing steps
    continueRuleParse = 0;
    % add to regulation evaluation statment that this reaction is always allowed
    regulationRulesEvalString = 'true';
    regulationDependancies = {''};

else

    % find any unexpected strings in rule
    invalidRuleMatchExp = '[^\w\s\(\)\>\<\=\[\]\.]'; % anything other than alphanumeric, underscores, whitespace and parenthases
    [~, ~, ~, unexpectedStrings, ~, ~, ~] = ...
                            regexp(rawRegulationString, invalidRuleMatchExp);
    if ~isempty( cat(2, unexpectedStrings{:}) )
        % set flag to skip other rule parsing steps
        continueRuleParse = 0;
        regulationErrors = [ regulationErrors; ...
                                ['unexpected characters in : ''', rawRegulationString,...
                                ''' associated with ', elementName] ];

    end % if unexpected charecters in regulation rule string

end % regulation string found to be empty or whitespace

% check for on/off true/false
[~, ~, ~, trueMatchStrings, ~, ~, trueMatchLeftovers] = ...
                            regexp(rawRegulationString, '((?i)\<true\>|((?i)\<on\>))');
[~, ~, ~, falseMatchStrings, ~, ~, falseMatchLeftovers] = ...
                            regexp(rawRegulationString, '((?i)\<false\>|((?i)\<off\>))');
if ~isempty(trueMatchStrings)
    
    % constant on/true found
    regulationRulesEvalString = 'true'; 
    regulationDependancies = {};
    
    % no need to continue parsing whether or not error found later
    continueRuleParse = 0;
    
    % check for leftovers
    
    [~, ~, ~, unexpectedLeftovers, ~, ~, ~] = regexp( trueMatchLeftovers, '\S' );
    
    
    if ~isempty( cat(2,unexpectedLeftovers{:}) )
        regulationErrors = [ regulationErrors; ...
                                {['unexpected characters found with true/on : ''', rawRegulationString,...
                                ''' associated with ', elementName]} ];
    end % if leftovers non-whitespace
    
elseif ~isempty(falseMatchStrings)
    
    % constant off/false found
    regulationRulesEvalString = 'false'; 
    regulationDependancies = {};
    
    % no need to continue parsing whether or not error found later
    continueRuleParse = 0;
    
    % check for leftovers
    [~, ~, ~, unexpectedLeftovers, ~, ~, ~] = regexp( falseMatchLeftovers, '\S' );
    
    if ~isempty( cat(2,unexpectedLeftovers{:}) )
        regulationErrors = [ regulationErrors; ...
                                {['unexpected characters found with false/off : ''', rawRegulationString,...
                                ''' associated with ', elementName]} ];
    end % if leftovers non-whitespace
    
end % if ture/false rule

% if no failure so far, find genes, etc
if continueRuleParse
    
    % replace logical operator strings (use regexprep so can ignore case) with spaces to
    % prevent identification as TFs durring mathching
    keywordsRemovedRegulationString = regexprep(rawRegulationString,'( or )|( and )|(not)',' ','ignorecase');
    
    % find using regular expressions all gene names of valid form: b folowed by four digits as whole word
    [~, ~, ~, geneNameMatchStrings, ~, ~, geneNameMatchLeftovers] = ...
                            regexp(keywordsRemovedRegulationString, '\<(b\d{4})\>');
    % number of valid gene strings found
    nRegulationGenes = length(geneNameMatchStrings);
    
    % find using regular expressions reactions names of valid form: rxnname[x] where x is
    % any valid single lowercase letter compartment designation
    [~, ~, ~, rxnNameMatchStrings, ~, ~, rxnNameMatchLeftovers] = ...
                            regexp( cat(2, geneNameMatchLeftovers{:}) , '\<([a-zA-Z]\w*(\[)[a-z]{1}(\]))\>');
    % number of valid gene strings found
    nRegulationRxn = length(rxnNameMatchStrings);
    
    % find using regular expressions media component names of valid form: exchangeMetName[media] 
    [~, ~, ~, mediaNameMatchStrings, ~, ~, mediaNameMatchLeftovers] = ...
                            regexp(cat(2, rxnNameMatchLeftovers{:}), '\<([0-9a-zA-Z]\w*(\[)media(\]))\>');
    % number of valid gene strings found
    nRegulationMedia = length(mediaNameMatchStrings);
    
    % find using regular expressions stimuli names of valid form: stimName[stim] 
    [~, ~, ~, stimNameMatchStrings, ~, ~, stimNameMatchLeftovers] = ...
                            regexp(cat(2, mediaNameMatchLeftovers{:}), '\<[a-zA-Z]\w*\[stim\]\>');
    % number of valid gene strings found
    nRegulationStim = length(stimNameMatchStrings);
    
    % find using regular expressions TF names of valid form: exchangeMetName[media] 
    [~, ~, ~, tfNameMatchStrings, ~, ~, tfNameMatchLeftovers] = ...
                            regexp(cat(2, stimNameMatchLeftovers{:}), '\<([a-zA-Z]\w*)\>');
    % number of valid gene strings found
    nRegulationTF = length(tfNameMatchStrings);
    
    % check leftovers for unexpected characters
    validLeftovers = '(AND)|(NOT)|(OR)|(and)|(or)|(not)|\s|\(|\)|[\<]|\=|(\d+\.\d+)|(\d)|[\>]';
    [~, ~, ~, ~, ~, ~, unexpectedLeftovers] = ...
                            regexp(cat(2, tfNameMatchLeftovers{:}), validLeftovers);
    % throw error if unexpected strings leftover
    if ~isempty(cat(2,unexpectedLeftovers{:}))
        % leftovers after all valid matches and operators removed, add error
        unexpectedLeftovers;
        continueRuleParse = 0;
        regulationErrors = [ regulationErrors; ...
                                ['formatting of : ''', rawRegulationString,''' does not match required formatting'] ];
        
    end % if unexpected leftovers after string matches collected

end % no initial failures found, continue with rule parsing by finding gene matches 



% check and update strings for genes in rule, and replace logical operator strings
if continueRuleParse

    % initialize modified rule string
    modifiedRegulationString = rawRegulationString;

    % for each gene match verify that it is found in model, and get index
    if nRegulationGenes > 0
    for iGeneMatch = 1 : nRegulationGenes

        geneMatch = geneNameMatchStrings{iGeneMatch};

        if isfield(organismModelStruct.genes, geneMatch)

            % gene is found in model, replace occurance in reaction regulation rule
            % with reference to corresponding index of logical gene activity array
            idxGeneMatch = organismModelStruct.genes.(geneMatch).idxGene;
            if organismModelStruct.TFgeneActivityInstant
                % flag for using gene activities, not gene product activities in TF rules
                geneMatchEvalString = [' geneActivity(',num2str(idxGeneMatch),') '];
            else
                geneMatchEvalString = [' geneProductActivity(',num2str(idxGeneMatch),') '];
            end
            modifiedRegulationString = strrep(modifiedRegulationString, geneMatch, geneMatchEvalString);

        else
            % gene match from rule not found in model
            continueRuleParse = 0;
            regulationErrors = [ regulationErrors; ...
                                ['unknown gene, ', geneMatch, '( found in regulation associated with: ',  elementName,')' ] ];
        end % if gene found in model

    end % for iGeneMatch
    end % if any regulation genes
    
    % for each media component match verify that it is found in model, and get index
    if nRegulationMedia > 0
        
        if nRegulationMedia > 1
            mediaMatchLengths = zeros(nRegulationMedia,1);
            for iSort = 1 : nRegulationMedia
                mediaMatchLengths(iSort) = length(mediaNameMatchStrings{iSort});
            end % for iSort
            
            [~,sortIdx] = sort(mediaMatchLengths,'descend');
            mediaNameMatchStrings = mediaNameMatchStrings(sortIdx);
            
        end % if nRegulationMedia > 1
        
        for iMediaMatch = 1 : nRegulationMedia

            mediaMatchRaw = mediaNameMatchStrings{iMediaMatch};
            mediaMatch = ['e',upper(mediaMatchRaw(1 : end - 7))];

            if isfield(organismModelStruct.mediaExchange, mediaMatch)

                % media is found in model, replace occurance in reaction regulation rule
                % with reference to corresponding index of logical media concentration activity array
                idxMediaMatch = organismModelStruct.mediaExchange.(mediaMatch).idxMedia;
                mediaMatchEvalString = [' currentMedia(',num2str(idxMediaMatch),') '];
                modifiedRegulationString = strrep(modifiedRegulationString, mediaMatchRaw, mediaMatchEvalString);

            else
                % media match from rule not found in model
                continueRuleParse = 0;
                regulationErrors = [ regulationErrors; ...
                                    ['media, ', mediaMatch, ' found in regulation associated with: ',  elementName,')' ] ];
            end % if media found in model

        end % for iMediaMatch
    end % if nRegualtionMedia > 0
    
    % for each stimuli match verify that it is found in model, and get index
    if nRegulationStim > 0
    for iStimMatch = 1 : nRegulationStim

        stimMatchRaw = stimNameMatchStrings{iStimMatch};
        stimMatch = upper(stimMatchRaw(1 : end - 6));
        
        if isfield(organismModelStruct.stimuli, stimMatch)
            
            % stimuli is found in model, replace occurance in reaction regulation rule
            % with reference to corresponding index of logical stimuli activity array
            idxStimMatch = organismModelStruct.stimuli.(stimMatch).idxStim;
            stimMatchEvalString = [' currentStimuli{',num2str(idxStimMatch),'} '];
            modifiedRegulationString = strrep(modifiedRegulationString, stimMatchRaw, stimMatchEvalString);

        else
            % stim match from rule not found in model
            continueRuleParse = 0;
            regulationErrors = [ regulationErrors; ...
                                {['stimuli, ', stimMatch, ' found in regulation associated with: ',  elementName,')' ]} ];
        end % if stim found in model

    end % for iStimMatch
    end % nRegulationStim > 0
    
    % for each rxn match verify that it is found in model, and get index
    if nRegulationRxn > 0
    for iRxnMatch = 1 : nRegulationRxn

        rxnMatchRaw = rxnNameMatchStrings{iRxnMatch};
        rxnMatch = [rxnMatchRaw(end - 1),upper(rxnMatchRaw(1 : end - 3))];

        if isfield(organismModelStruct.rxns, rxnMatch)

            % rxn is found in model, replace occurance in reaction regulation rule
            % with reference to corresponding index of flux activity array
            idxRxnMatch = organismModelStruct.rxns.(rxnMatch).idxRxn;
            rxnMatchEvalString = [' currentFlux(',num2str(idxRxnMatch),') '];
            modifiedRegulationString = strrep(modifiedRegulationString, rxnMatchRaw, rxnMatchEvalString);

        else
            % rxn match from rule not found in model
            continueRuleParse = 0;
            regulationErrors = [ regulationErrors; ...
                                ['rxn, ', rxnMatch, ' found in regulation associated with: ',  elementName,')' ] ];
        end % if rxn found in model

    end % for iRxnMatch
    end % if nRegulationRxn > 0
    
    % for each TF match verify that it is found in model, and get index
    if nRegulationTF > 0
    for iTFmatch = 1 : nRegulationTF
        
        tfMatch = tfNameMatchStrings{iTFmatch};

        if isfield(organismModelStruct.TFs, tfMatch)

            % TF is found in model, replace occurance in reaction regulation rule
            % with reference to corresponding index of logical TF activity array
            idxTFmatch = organismModelStruct.TFs.(tfMatch).idxTF;
            tfMatchEvalString = [' tfActivity(',num2str(idxTFmatch),') '];
            modifiedRegulationString = strrep(modifiedRegulationString, tfMatch, tfMatchEvalString);

        else
            % TFs match from rule not found in model
            continueRuleParse = 0;
            regulationErrors = [ regulationErrors; ...
                                ['unknown TF, ', tfMatch, ' found in regulation associated with: ',  elementName,')' ] ];
        end % if TF found in model

    end % for iTFmatch
    end % if nRegulationTf > 0
    
    % replace logical operator strings (use regexprep so can ignore case)
    modifiedRegulationString = regexprep(modifiedRegulationString,' or ','|','ignorecase');
    modifiedRegulationString = regexprep(modifiedRegulationString,' and ','&','ignorecase');
    modifiedRegulationString = regexprep(modifiedRegulationString,'not',' ~ ','ignorecase');     

end % if continue parsing (check and update strings for genes in rule, and replace logical operator strings)



% if continueRuleParse, test validity of expression (for error producing evaluation issues)
if continueRuleParse

    % try evaluating rule, if error is produced
    try
        regulationStringEvalTest = islogical( eval(modifiedRegulationString) );

        if regulationStringEvalTest

            % produced a logical value, store for output
            regulationRulesEvalString = modifiedRegulationString;
            regulationDependancies = unique(geneNameMatchStrings);

        else

            % did not produce a logical value, fail
             continueRuleParse = 0;
             regulationErrors = [ regulationErrors ;
                           ['non-logical values produced upon eval of ''', modifiedRegulationString,...
                            ]]; %''' which was orrigionally ''',rawRegulationString, ...
                            % ''' associated with: ', elementName] ] ;

        end % logical value producted by rule evaluation test

    catch errorOnRuleEval % figure out how to use this value?

        % error was encountered in evaluation, flag and store info
%         continueRuleParse = 0;
        regulationErrors = [ regulationErrors ;
                           ['error was producted in evaluation of ''', modifiedRegulationString,...
                            ''' which was orrigionally ''',rawRegulationString, ...
                            ''' associated with: ', elementName] ] ;

    end % try/catch evaluation of logical

end % if continueRuleParse, test validity of expression (for error producing evaluation issues)
    
