cplexFastSL/cplexFastSL.m at master · omkarmohite310/cplexFastSL · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
function cplexFastSL(model,cutoff,order,eliList,atpm)
%% fastSL(model,cutoff,order,eliList,atpm)
% Requires the openCOBRA toolbox
% http://opencobra.sourceforge.net/openCOBRA/Welcome.html
%
% INPUT
% model (the following fields are required - others can be supplied)
%   S            Stoichiometric matrix
%   b            Right hand side = dx/dt
%   c            Objective coefficients
%   lb           Lower bounds
%   ub           Upper bounds
%   rxns         Reaction Names
%OPTIONAL
% cutoff         cutoff percentage value for lethality.Default is 0.01.
% order          Order of SLs required.Default order is 2. Max value 3.
% eliList        List of reactions to be ignored for lethality
% analysis:Exchange Reactions, ATPM etc.
% atpm           ATPM Reaction Id in model.rxns if other than 'ATPM'
%
%
% Aditya Pratapa      6/26/14.
%%
profile on;

if exist('cutoff', 'var')
    if isempty(cutoff)
        cutoff = 0.01;
    end
else
    cutoff = 0.01;
end

if exist('order', 'var')
    if isempty(order)
        order = 2;
    else
        if (order>3)
        err = MException('ResultChk:OutOfRange', ...
        'Resulting value is outside expected range. Maximum value is 3.');
         throw(err)
        end
    end
else
    order = 2;
end


%Please change this according to your model
if exist('atpm', 'var')
    if isempty(atpm)
        atpm = 'ATPM'; %Reaction Id of ATP maintenance reaction- by default it takes 'ATPM'
    end
else
    atpm = 'ATPM';
end


if exist('eliList', 'var')
    if isempty(eliList)
        eliList = model.rxns(ismember(model.rxns,atpm)); %To eliminate ATPM.
    end
else
    eliList = model.rxns(ismember(model.rxns,atpm));
end

fname=strcat(model.description,'Cplex_Rxn_lethals.mat');
% fname=strcat('Ecoli_lethals.mat');


%%

switch order
    case 1
        [Jsl]=cplexSingleSL(model,cutoff,eliList,atpm);

        fprintf('\n Saving Single Lethal Reactions List...\n');
        save(fname,'Jsl');
        fprintf('Done. \n');
    case 2
        [Jsl,Jdl]=cplexDoubleSL(model,cutoff,eliList,atpm);

        fprintf('\n Saving Single and Double Lethal Reactions List...\n');
        save(fname,'Jsl');
        save(fname,'Jdl','-append');
        fprintf('Done. \n');
    case 3
        [Jsl,Jdl,Jtl]=cplexTripleSL(model,cutoff,eliList,atpm);

        fprintf('\n Saving Single, Double and Triple Lethal Reactions List...\n');
        save(fname,'Jsl');
        save(fname,'Jdl','-append');
        save(fname,'Jtl','-append');
        fprintf('Done. \n');
end


profview