deletemodelgenes.m

预测得到最优代谢通量分布而进行最适单基因敲除的例程

function [model,hasEffect,constrRxnNames] = deleteModelGenes(model,geneList,downRegFraction)
%deleteModelGenes Delete one or more genes and constrain the reactions
%affected to zero
%
% [model,hasEffect] = deleteModelGenes(model,geneList,downRegFraction)
%
% model             Constraint-based model with the appropriate constrains for
%                   a particular condition
% geneList          List of genes to be deleted (opt, default all genes in
%                   model)
% downRegFraction   Fraction of the original bounds that the reactions
%                   corresponding to downregulated genes will be assigned
%                   (opt, default = 0 corresponding to a full deletion)
%
% model             Model with the reactions deleted
% hasEffect         True if the gene deletion has an effect on the model
% rxnList           reactions that are associated to the genes in geneList
%
% Markus Herrgard 8/28/06

if (nargin < 2)
    geneList = model.genes;
end

if (nargin < 3)
    downRegFraction = 0;
end

if (~iscell(geneList))
    geneName = geneList;
    clear geneList;
    geneList{1} = geneName;
end

if (~isfield(model,'genes'))
    error('Gene-reaction associations not included with the model');
end

hasEffect = false;
constrRxnNames = {};

% Find gene indices in model
[isInModel,geneInd] = ismember(geneList,model.genes);
if (any(isInModel))
    % Find rxns associated with this gene
    rxnInd = find(any(model.rxnGeneMat(:,geneInd),2));
    if (~isempty(rxnInd))
        x = true(size(model.genes));
        x(geneInd) = false;
        constrainRxn = false(length(rxnInd),1);
        % Figure out if any of the reaction states is changed
        for j = 1:length(rxnInd)
            if (~eval(model.rules{rxnInd(j)}))
                constrainRxn(j) = true;
            end
        end
        % Constrain flux through the reactions associated with these genes
        if (any(constrainRxn))
            constrRxnNames = model.rxns(rxnInd(constrainRxn));
            if (nargin > 2)
                model = changeRxnBounds(model,constrRxnNames,downRegFraction*model.lb(findRxnIDs(model,constrRxnNames)),'l');
                model = changeRxnBounds(model,constrRxnNames,downRegFraction*model.ub(findRxnIDs(model,constrRxnNames)),'u');
            else
                % Full deletion
                model = changeRxnBounds(model,constrRxnNames,0,'b');
            end
            hasEffect = true;
        end
    end
end