GetKegg

 GetKegg
 reads flux values obtained from COBRA optimization, and uses the 
 input BiGG reactions to retreive the equivalent KEGG IDs obtained from 
 Bigg2Kegg function. 
 Inputs:
 Bigg = A cell array of BiGG reaction IDs (can be obtained from readCbModel
        function of COBRA toolbox) for which flux values to be displayed
        on KEGG maps.
 Inflx = A mXn matrix containing input flux data, which m corresponds to
         BiGG IDs and n corresponds to time-series values.
 ModelName = Model description according to BiGG database (e.g. RECON1)
 Output:
 MapChoice = A cell array of KEGG pathway identifiers chosen by the user
             via a simple GUI. These set of KEGG pathway IDs will be used
             as the input of KeggDraw function.
 RxnCds = KEGG reaction identifiers equivalent to input Bigg IDs.
 Outflx = Similar to Inflx, with the difference that m corresponds to
          RxnCds (Note: if KEGG rxn IDs for some Bigg elements cannot be 
          found, these IDs should be trimmed).
 
 O. Jamialahmadi
 TMU, Chem. Eng. Dept., Biotech. Group 
 Nov. 2015

0001 function [MapChoice, Outflx, RxnCds] = GetKegg(Bigg, Inflx, ModelName)
0002 % GetKegg
0003 % reads flux values obtained from COBRA optimization, and uses the
0004 % input BiGG reactions to retreive the equivalent KEGG IDs obtained from
0005 % Bigg2Kegg function.
0006 % Inputs:
0007 % Bigg = A cell array of BiGG reaction IDs (can be obtained from readCbModel
0008 %        function of COBRA toolbox) for which flux values to be displayed
0009 %        on KEGG maps.
0010 % Inflx = A mXn matrix containing input flux data, which m corresponds to
0011 %         BiGG IDs and n corresponds to time-series values.
0012 % ModelName = Model description according to BiGG database (e.g. RECON1)
0013 % Output:
0014 % MapChoice = A cell array of KEGG pathway identifiers chosen by the user
0015 %             via a simple GUI. These set of KEGG pathway IDs will be used
0016 %             as the input of KeggDraw function.
0017 % RxnCds = KEGG reaction identifiers equivalent to input Bigg IDs.
0018 % Outflx = Similar to Inflx, with the difference that m corresponds to
0019 %          RxnCds (Note: if KEGG rxn IDs for some Bigg elements cannot be
0020 %          found, these IDs should be trimmed).
0021 %
0022 % O. Jamialahmadi
0023 % TMU, Chem. Eng. Dept., Biotech. Group
0024 % Nov. 2015
0025 
0026 % Retrieve BiGG2KEGG rxn identifiers from Bigg2Kegg function---------------
0027 Pth1 = which('GetKegg.m');
0028 tind = find(Pth1=='\',1,'last');
0029 Pth2 = Pth1(1:tind-1);
0030 Nme = [ModelName,'KEGG.mat'];
0031 TargetPth = fullfile(Pth2,'BIGG2KEGG',Nme);
0032 load(TargetPth);
0033 Brxns1 = B2Kegg.B; Krxns1 = B2Kegg.K;
0034 % -------------------------------------------------------------------------
0035 % Find common rxns between Brxns and Bigg
0036 % First, seperate multi-step rxns------------------------------------------
0037 MultiB1 = Brxns1(find(strcmp(Brxns1,'MULTIR'))+1:end);
0038 MultiK1 = Krxns1(find(strcmp(Brxns1,'MULTIR'))+1:end);
0039 MultiK = MultiK1(ismember(MultiB1, Bigg));
0040 MultiB = MultiB1(ismember(MultiB1, Bigg));
0041 [~,Fx22] = ismember(MultiB1, Bigg);
0042 Fx22(Fx22==0) = [];
0043 MultiFlx = Inflx(Fx22,:);
0044 setappdata(0,'MultiB',MultiB); % For later use in KeggDraw
0045 setappdata(0,'MultiK',MultiK);
0046 setappdata(0,'MultiFlx',MultiFlx);
0047 % -------------------------------------------------------------------------
0048 Brxns = Brxns1(1:find(strcmp(Brxns1,'MULTIR'))-1);
0049 Krxns = Krxns1(1:find(strcmp(Brxns1,'MULTIR'))-1);
0050 RxnCdsI = Krxns(ismember(Brxns, Bigg));
0051 BiGGI = Brxns(ismember(Brxns, Bigg));
0052 [~,Fx2] = ismember(Brxns, Bigg);
0053 Fx2(Fx2==0) = [];
0054 FlxI = Inflx(Fx2,:);
0055 % If any, ask user to fill missed KEGG rxn IDs
0056 if any(~ismember(Bigg,Brxns))
0057     disp('Some KEGG rxns cannot be found for input BiGG IDs')
0058     disp('Insert them manually')
0059     disp('NOTE: Exchange reactions are not displayed!')
0060     NotData = Bigg(~ismember(Bigg,Brxns));
0061     % Exclude Exchange reactions : EX_w\*
0062     NoDataIdf = regexp(NotData,'\<EX_\w*');
0063     NoDataIdf1 = ~cellfun('isempty', NoDataIdf);
0064     NotData(NoDataIdf1) = [];
0065     %------------------------------------
0066     setappdata(0,'NotData',NotData)
0067     waitfor(GetKeggTable)
0068     ModfdTab = getappdata(0,'Final');
0069     FindUndfnd = ModfdTab(:,2);
0070     BiGGUndfnd =  ModfdTab(:,1);
0071     DfinedLoci = ~strcmpi(FindUndfnd,'Unknown');
0072     RxnCdsII = FindUndfnd(DfinedLoci);
0073     BiGGII = BiGGUndfnd(DfinedLoci);
0074     Flxs = Inflx(~ismember(Bigg,Brxns),:);
0075     FlxII = Flxs(DfinedLoci,:);
0076     RxnCds = [RxnCdsI';RxnCdsII];
0077     BiGG4KeggDraw = [BiGGI';BiGGII];
0078     setappdata(0,'BiGG4KeggDraw',BiGG4KeggDraw);
0079     Outflx = [FlxI;FlxII];
0080 else
0081     RxnCds = RxnCdsI;
0082     Outflx = FlxI;
0083 end
0084 %--------------------------------------------------------------------------
0085 % Find all KEGG pathways for RxnCds
0086 Fileid1 = fopen('rxn2map.txt','r');
0087 rxn2map = textscan(Fileid1,'%s %s');
0088 RawRxns = rxn2map{1};
0089 RawMaps = rxn2map{2};
0090 [Loci1,~] = regexp(RawMaps, 'path:rn');
0091 Loci2 = ~cellfun('isempty', Loci1);
0092 RawMaps(Loci2) = [];
0093 RawRxns(Loci2) = [];
0094 Maps = strrep(RawMaps,'path:map','');
0095 Rxns = strrep(RawRxns,'rn:','');
0096 RxnTemp = Rxns(ismember(Rxns,RxnCds));
0097 MapTemp = Maps(ismember(Rxns,RxnCds));
0098 count=1;
0099 RxnList = {0}; MapList = {0}; SortTemp = (0);
0100 RmvMaps = {'01100','01110','01120','01130','01200','01210','01212',...
0101     '01220','01230','00121'};
0102 while numel(MapTemp)
0103     Loci3 = ismember(MapTemp,MapTemp(1));
0104     if ~ismember(MapTemp(1),RmvMaps)
0105         RxnList{count} = RxnTemp(Loci3);
0106         SortTemp(count) = numel(RxnList{count});
0107         MapList(count) = MapTemp(1);
0108         count=count+1;
0109     end
0110     RxnTemp(Loci3)=[];
0111     MapTemp(Loci3)=[];
0112 end
0113 % Sort MapList based on rxns in RxnList
0114 [~,ST] = sort(SortTemp,'descend');
0115 RxnList = RxnList(ST);
0116 MapList = MapList(ST);
0117 % Add pathway names to map codes of MapList
0118 Fileid2 = fopen('KEGGmaps.txt','r');
0119 TempMaps = textscan(Fileid2,'%s %[^\n]');
0120 TempMaps1 = cell(size(TempMaps{1},1),1);
0121 for count = 1:size(TempMaps{1},1)
0122     TempMaps1{count}=[TempMaps{1}{count},' ',TempMaps{2}{count}];
0123 end
0124 TempMaps1 = strrep(TempMaps1,'path:map','');
0125 MapNames2 = cell(numel(MapList),1);
0126 for count = 1:numel(MapList)
0127     Loci4 = regexp(TempMaps1,MapList{count});
0128     MapNames2{count} = TempMaps1{~cellfun('isempty', Loci4)};
0129 end
0130 
0131 setappdata(0,'MapNames2',MapNames2)
0132 setappdata(0,'RxnList',RxnList)
0133 waitfor(GetKeggMaps)
0134 MapChoiceLoci = getappdata(0,'AllChoices');
0135 MapChoice = MapList(MapChoiceLoci);
0136 % Remove previous data in appdata------------------------------------------
0137 AT=getappdata(0);
0138 AT1 = fieldnames(AT);
0139 for co = 1:length(AT1)
0140     if strcmp(AT1{co},'BiGGModel') || strcmp(AT1{co},'Mdl') || ...
0141             strcmp(AT1{co},'CRxns')
0142         rmappdata(0,AT1{co});
0143     end
0144 end