0001 function [val] = optGeneFitness(rxn_vector_matrix, model, targetRxn, rxnListInput, isGeneList)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017 global MaxKnockOuts
0018
0019
0020 popsize = size(rxn_vector_matrix,1);
0021 val = zeros(1,popsize);
0022
0023 for i = 1:popsize
0024 rxn_vector = rxn_vector_matrix(i,:);
0025 rxnList = rxnListInput(logical(rxn_vector));
0026
0027
0028
0029 val_temp = memoize(rxn_vector);
0030 if ~ isempty(val_temp)
0031 val(i) = val_temp;
0032 continue;
0033 end
0034
0035
0036 nummutations = sum(rxn_vector);
0037 if nummutations > MaxKnockOuts
0038 continue;
0039 end
0040
0041
0042 if isGeneList
0043 modelKO = deleteModelGenes(model, rxnList);
0044 else
0045 [isValidRxn,removeInd] = ismember(rxnList,model.rxns);
0046 removeInd = removeInd(isValidRxn);
0047 modelKO = model;
0048 modelKO.ub(removeInd) = 0;
0049 modelKO.lb(removeInd) = 0;
0050 end
0051
0052
0053
0054
0055 if exist('LPBasis', 'var')
0056 modelKO.LPBasis = LPBasis;
0057 end
0058
0059 [slnKO, LPOUT] = solveCobraLPCPLEX(modelKO, 0,1);
0060 LPBasis = LPOUT.LPBasis;
0061 growthrate = slnKO.obj;
0062
0063
0064
0065 if growthrate < .10
0066 continue;
0067 end
0068
0069
0070
0071
0072
0073
0074
0075 gamma = 0.01;
0076
0077
0078 indBOF = find(modelKO.c);
0079 indTar = findRxnIDs(modelKO, targetRxn);
0080
0081 modelKOsetGR = modelKO;
0082 modelKOsetGR.lb(indBOF) = growthrate - gamma;
0083 modelKOsetGR.c = zeros(size(modelKO.c));
0084 modelKOsetGR.c(indTar) = -1;
0085
0086
0087
0088
0089
0090
0091 if exist('LPBasis2', 'var')
0092 modelKOsetGR.LPBasis = LPBasis2;
0093 end
0094
0095 [slnKOsetGR, LPOUT2] = solveCobraLPCPLEX(modelKOsetGR, 0,1);
0096 LPBasis2 = LPOUT2.LPBasis;
0097 minProdAtSetGR = -slnKOsetGR.obj;
0098
0099
0100
0101
0102 val(i) = -minProdAtSetGR;
0103
0104
0105
0106
0107
0108
0109 if val(i) > -1e-3
0110 val(i) = 0;
0111 end
0112
0113 memoize(rxn_vector, val(i));
0114 end
0115
0116 return;
0117
0118
0119
0120
0121
0122
0123
0124
0125
0126 function [value] = memoize(gene_vector, value)
0127 global HTABLE
0128 hashkey = num2str(gene_vector);
0129 hashkey = strrep(hashkey,' ','');
0130
0131 if nargin == 1
0132 value = HTABLE.get(hashkey);
0133 return;
0134 else
0135 if HTABLE.size() > 10000
0136 HTABLE = java.util.Hashtable;
0137 end
0138 HTABLE.put(hashkey, value);
0139 value = [];
0140 return;
0141 end
0142 return