mapCompartments

 mapCompartments
   Maps compartments in the geneScoreStructure. This is used if you do not
   want a models that uses all of the compartment from the predictor. This
   function will then let you define rules on how the compartments should
   be merged.

   Any number of rules could be defined as consecutive strings or in a cell array.
   'comp1'             comp1 should be kept in the structure

0001 function geneScoreStructure=mapCompartments(geneScoreStructure,varargin)
0002 % mapCompartments
0003 %   Maps compartments in the geneScoreStructure. This is used if you do not
0004 %   want a models that uses all of the compartment from the predictor. This
0005 %   function will then let you define rules on how the compartments should
0006 %   be merged.
0007 %
0008 %   Any number of rules could be defined as consecutive strings or in a cell array.
0009 %   'comp1'             comp1 should be kept in the structure
0010 
0011 %   'comp1=comp2'       The scores in comp2 are merged to comp1 and comp2 is
0012 %                       removed from the structure. This automatically
0013 %                       keeps comp1 in the structure
0014 
0015 %   'comp1=comp2 comp3' The scores in comp2 and comp3 are merged to comp1
0016 %                       and comp2 & comp2 are removed from the structure.
0017 %                       This automatically keeps comp1 in the structure
0018 
0019 %   'comp1 comp2=comp3' The scores in comp3 are split between comp1 and
0020 %                       comp2. This automatically keeps comp1 and comp2 in
0021 %                       the structure
0022 
0023 %   'comp1=other'       The scores in any compartment not included are
0024 %                       merged to comp1. This is applied after all other
0025 %                       rules.
0026 %
0027 %   When one compartment is merged to another the resulting scores will be
0028 %   the best for each gene in either of the compartments. In the case where
0029 %   one compartment is split among several, the scores for the compartment
0030 %   to be merged is weighted with the number of compartments to split to.
0031 %
0032 %   Example: The predictor you use give prediction for Extracellular,
0033 %   Cytosol, Nucleus, Peroxisome, Mitochondria, ER, and Lysosome. You want to
0034 %   have a model with Extracellular, Cytosol, Mitochondria, and Peroxisome
0035 %   where Lysosome is merged with Peroxisome and all other compartments
0036 %   are merged to the Cytosol.
0037 %
0038 %   GSS=mapCompartments(GSS,'Extracellular','Mitochondria','Peroxisome=Lyso
0039 %   some','Cytosol=other');
0040 %
0041 %   geneScoreStructure  a structure to be used in predictLocalization
0042 %
0043 %   Usage: geneScoreStructure=mapCompartments(geneScoreStructure,varargin)
0044 
0045 varargin=upper(varargin);
0046 
0047 %First find the compartment that will end up in the final structure. They
0048 %are the ones that stand alone or are to the left of some '='
0049 toKeep={};
0050 toMerge={};
0051 I=regexp(varargin,'=','split');
0052 for i=1:numel(varargin)
0053     if numel(I{i})==1
0054         toKeep=[toKeep;I{i}];
0055     else
0056         J=regexp(I{i}(1),' ','split');
0057         K=regexp(I{i}(2),' ','split');
0058         toKeep=[toKeep;J{1}(:)];
0059         toMerge=[toMerge;K{1}(:)];
0060     end
0061 end
0062 
0063 %Check that there are no compartment that should both be merged and kept
0064 if ~isempty(intersect(toKeep,toMerge))
0065     EM='There are inconsistencies where one or more compartment(s) should be both kept and merged to another';
0066     dispEM(EM);
0067 end
0068 
0069 %Check that there are no compartments in the rules that are not in the
0070 %geneScoreStructure.
0071 uComps=upper(geneScoreStructure.compartments);
0072 J=[uComps,{'OTHER'}];
0073 
0074 if ~isempty(setdiff([toKeep;toMerge],J))
0075     EM='There are compartment in the rules that are not in geneScoreStructure.compartments';
0076     dispEM(EM);
0077 end
0078 
0079 %Loop through it again and do the mapping
0080 otherIndex=[]; %This stores the rule which maps 'other'.
0081 
0082 for i=1:numel(I)
0083     if numel(I{i})>1
0084         %Get the compartment indexes that should be mapped
0085         J=regexp(I{i}(2),' ','split');
0086         if strcmpi(J{1},'other')
0087             otherIndex=i;
0088             continue;
0089         end
0090         [k, K]=ismember(J{1},uComps);
0091         
0092         %And to where they should be mapped
0093         J=regexp(I{i}(1),' ','split');
0094         [l, L]=ismember(J{1},uComps);
0095         
0096         %It's not allowed to have rules like A B=C D
0097         if numel(K)>1 && numel(L)>1
0098             EM='It is not allowed to have rules like "A B=C D" (map more than one compartment to more than one compartment)';
0099             dispEM(EM);
0100         end
0101         
0102         if ~all(k) || ~all(l)
0103             EM='Error in mapping. This most likely means that some compartment(s) are mapped to different compartments in different rules. Use A B=C if you want to map C to several compartments';
0104             dispEM(EM);
0105         end
0106         
0107         %Get the sum of the scores for the compartments that should be
0108         %merged to something else
0109         S=max(geneScoreStructure.scores(:,K),[],2);
0110         for j=1:numel(L)
0111             %If the scores are mapped to several different compartments
0112             %then split the scores between them
0113             geneScoreStructure.scores(:,L(j))=max(geneScoreStructure.scores(:,L(j)),S./numel(L));
0114         end
0115         
0116         %Remove the comparement that were merged
0117         geneScoreStructure.compartments(K)=[];
0118         geneScoreStructure.scores(:,K)=[];
0119         uComps(K)=[];
0120     end
0121 end
0122 
0123 %Then check if there are remaining compartments that should be removed or
0124 %mapped as 'other'
0125 J=find(~ismember(uComps,toKeep));
0126 if any(J)
0127     if any(otherIndex)
0128         K=regexp(I{otherIndex}(1),' ','split');
0129         [l, L]=ismember(K{1},uComps);
0130         if l==1 && numel(l)==1
0131             S=max(geneScoreStructure.scores(:,J),[],2);
0132             geneScoreStructure.scores(:,L)=max(geneScoreStructure.scores(:,L),S);
0133         else
0134             EM='Could not map "other" to more than one compartment';
0135             dispEM(EM);
0136         end
0137     else
0138         EM='There are compartments that are not defined if they should be kept or removed. Use "A=other" or define more rules if you do not want them to be deleted';
0139         dispEM(EM,false);
0140     end
0141     
0142     %Remove the comparement that were merged
0143     geneScoreStructure.compartments(J)=[];
0144     geneScoreStructure.scores(:,J)=[];
0145 end
0146 
0147 %Renormalize
0148 I=max(geneScoreStructure.scores,[],2);
0149 geneScoreStructure.scores=bsxfun(@times, geneScoreStructure.scores, 1./I);
0150 
0151 %If there are genes that have score 0 in all compartments, remove them and
0152 %print a warning.
0153 I=find(isnan(geneScoreStructure.scores(:,1))); %Only looks a the first colum as it will be the same for the other ones
0154 if any(I)
0155     EM='The following genes had score 0.0 in all compartments. They have been removed from the structure. Consider using more rules or "A=other" in order to prevent this:';
0156     dispEM(EM,false,geneScoreStructure.genes(I));
0157     geneScoreStructure.scores(I,:)=[];
0158     geneScoreStructure.genes(I)=[];
0159 end
0160 end