disp('Running MATLAB script file NeilLung1.m') ;
%
%    For analysis of Neil Haye's Lung Cancer Data
%    So far data are those Yufeng Liu has been analyzing:
%        - Genes are filtered down to 2530 "expressed genes"
%        - Cases restricted to only:
%                 Carcinoid
%                 Colon
%                 Normal
%                 SmallCell
%
%    Does some visualizations
%
%    Initial input files are:
%        lung4group.nolabel.txt
%        lung4group.labelonly.txt
%
%


ipart = 332 ;    %  0 - Read data, check inputs and save as:
               %         NeilLung1.mat
               %          (everything else reads this)
               %  1 - PCA Plot of Combined Data  Matlab default colors)
               %  2 - PCA plot of combined data (all black)
               %  3 - PCA plot of combined data (Class colors)
               %  4 - PCA plot of combined data (Class colors and symbols)
               %  5 - PCA plot of combined data (Class symbols only)

               %  11 - Proj on MD Carc-Col, MD Norm-SC, & ortho PC 
               %  12 - Proj on MD Carc-Norm, MD Col-SC, & ortho PC 
               %  13 - Proj on MD Carc-SC, MD Norm-Col, & ortho PC 
               %  14 - Proj on MD (Carc&SC)-(Norm&Col), MD (Carc&Norm)-(SC&Col), & ortho PC 
               %  15 - Proj on MD (Carc&Norm)-(SC&Col), MD (Carc&Col)-(Norm&SC), & ortho PC 
               %  16 - Proj on MD (Carc&Col)-(SC&Norm), MD (Carc&SC)-(Norm&Col), & ortho PC 
               %  17 - Proj on MD Indiv vs. all

               %  21 - Proj on DWD Carc-Col, DWD Norm-SC, & ortho PC 
               %  22 - Proj on DWD Carc-Norm, DWD Col-SC, & ortho PC 
               %  23 - Proj on DWD Carc-SC, DWD Norm-Col, & ortho PC 
               %  24 - Proj on DWD (Carc&SC)-(Norm&Col), DWD (Carc&Norm)-(SC&Col), & ortho PC 
               %  25 - Proj on DWD (Carc&Norm)-(SC&Col), DWD (Carc&Col)-(Norm&SC), & ortho PC 
               %  26 - Proj on DWD (Carc&Col)-(SC&Norm), DWD (Carc&SC)-(Norm&Col), & ortho PC 
               %  27 - Proj on DWD Indiv vs. all

               %  31 - Proj on SVM Carc-Col, SVM Norm-SC, & ortho PC 
               %  32 - Proj on SVM Carc-Norm, SVM Col-SC, & ortho PC 
               %  33 - Proj on SVM Carc-SC, SVM Norm-Col, & ortho PC 
               %  34 - Proj on SVM (Carc&SC)-(Norm&Col), SVM (Carc&Norm)-(SC&Col), & ortho PC 
               %  35 - Proj on SVM (Carc&Norm)-(SC&Col), SVM (Carc&Col)-(Norm&SC), & ortho PC 
               %  36 - Proj on SVM (Carc&Col)-(SC&Norm), SVM (Carc&SC)-(Norm&Col), & ortho PC 
               %  37 - Proj on DWD Indiv vs. all

               %  41 - Tree1, level 0
               %            Proj on DWD (Carc&Norm)-(SC&Col), DWD Norm-Carc, DWD SC&Col, & ortho PC 
               %            All Black
               %  42 - Tree1, level 1
               %            Proj on DWD (Carc&Norm)-(SC&Col), DWD Norm-Carc, DWD SC&Col, & ortho PC 
               %            Level 1 Col
               %  43 - Tree1, level 2
               %            Proj on DWD (Carc&Norm)-(SC&Col), DWD Norm-Carc, DWD SC&Col, & ortho PC 
               %            All Col
               %  44 - Tree 2, level 0
               %            Proj on DWD (Carc&SC)-(Norm&Col), DWD Carc-SC, DWD Norm-Col, & ortho PC 
               %            All black
               %  45 - Tree 2, level 1
               %            Proj on DWD (Carc&SC)-(Norm&Col), DWD Carc-SC, DWD Norm-Col, & ortho PC 
               %            Level 1 Col
               %  46 - Tree 2, level 2
               %            Proj on DWD (Carc&SC)-(Norm&Col), DWD Carc-SC, DWD Norm-Col, & ortho PC 
               %            All Col
               %  47 - Tree 3, level 0
               %            Proj on DWD SC-(Carc&Norm&Col), DWD Carc-(Col&Norm), DWD Norm-Col, & ortho PC 
               %            All black
               %  48 - Tree 3, level 1
               %            Proj on DWD SC-(Carc&Norm&Col), DWD Carc-(Col&Norm), DWD Norm-Col, & ortho PC 
               %            Level 1
               %  49 - Tree 3, level 2
               %            Proj on DWD SC-(Carc&Norm&Col), DWD Carc-(Col&Norm), DWD Norm-Col, & ortho PC 
               %            Level 2
               %  50 - Tree 3, level 3
               %            Proj on DWD SC-(Carc&Norm&Col), DWD Carc-(Col&Norm), DWD Norm-Col, & ortho PC 
               %            Level 3 Col

               %  Using "new labelling", where one Normal is classified as Colon
               %  51 - Proj on DWD Carc-Col, DWD Norm-SC, & ortho PC 
               %  52 - Proj on DWD Carc-Norm, DWD Col-SC, & ortho PC 
               %  53 - Proj on DWD Carc-SC, DWD Norm-Col, & ortho PC 
               %  54 - Proj on DWD (Carc&SC)-(Norm&Col), DWD (Carc&Norm)-(SC&Col), & ortho PC 
               %  55 - Proj on DWD (Carc&Norm)-(SC&Col), DWD (Carc&Col)-(Norm&SC), & ortho PC 
               %  56 - Proj on DWD (Carc&Col)-(SC&Norm), DWD (Carc&SC)-(Norm&Col), & ortho PC 
               %  57 - Proj on DWD Indiv vs. all

               %  61 - Tree 3, based on Random Relabelling, Seed 1, Relab colors 
               %  62 - Tree 3, based on Random Relabelling, Seed 1, Original colors 
               %  63 - Tree 3, based on Random Relabelling, Seed 2, Relab colors 
               %  64 - Tree 3, based on Random Relabelling, Seed 2, Original colors 
               %  65 - Tree 3, based on Random Relabelling, Seed 3, Relab colors 
               %  66 - Tree 3, based on Random Relabelling, Seed 3, Original colors 

               %  71 - Tree 3, based on 80% subsample, full colors
               %  72 - Tree 3, based on 80% subsample, missings grayed out
               %  73 - Tree 3, based on 60% subsample, full colors
               %  74 - Tree 3, based on 60% subsample, missings grayed out
               %  75 - Tree 3, based on 40% subsample, full colors
               %  76 - Tree 3, based on 40% subsample, missings grayed out
               %  77 - Tree 3, based on 20% subsample, full colors
               %  78 - Tree 3, based on 20% subsample, missings grayed out

               %  101 - Tree 3, DWD, all 2530 genes unmasked, DWD SC vs. Rest
               %  102 - Tree 3, DWD, 1200 genes unmasked, DWD SC vs. Rest
               %  103 - Tree 3, DWD, 600 genes unmasked, DWD SC vs. Rest
               %  104 - Tree 3, DWD, 300 genes unmasked, DWD SC vs. Rest
               %  105 - Tree 3, DWD, 150 genes unmasked, DWD SC vs. Rest
               %  106 - Tree 3, DWD, 90 genes unmasked, DWD SC vs. Rest
               %  107 - Tree 3, DWD, 60 genes unmasked, DWD SC vs. Rest

               %  111 - Tree 3, DWD, all 2530 genes unmasked, DWD Carc vs. NormCol
               %  112 - Tree 3, DWD, 1200 genes unmasked, DWD Carc vs. NormCol
               %  113 - Tree 3, DWD, 600 genes unmasked, DWD Carc vs. NormCol
               %  114 - Tree 3, DWD, 300 genes unmasked, DWD Carc vs. NormCol
               %  115 - Tree 3, DWD, 150 genes unmasked, DWD Carc vs. NormCol
               %  116 - Tree 3, DWD, 90 genes unmasked, DWD Carc vs. NormCol
               %  117 - Tree 3, DWD, 60 genes unmasked, DWD Carc vs. NormCol

               %  121 - Tree 3, DWD, all 2530 genes unmasked, DWD Norm vs.Col
               %  122 - Tree 3, DWD, 1200 genes unmasked, DWD Norm vs.Col
               %  123 - Tree 3, DWD, 600 genes unmasked, DWD Norm vs.Col
               %  124 - Tree 3, DWD, 300 genes unmasked, DWD Norm vs.Col
               %  125 - Tree 3, DWD, 150 genes unmasked, DWD Norm vs.Col
               %  126 - Tree 3, DWD, 90 genes unmasked, DWD Norm vs.Col
               %  127 - Tree 3, DWD, 60 genes unmasked, DWD Norm vs.Col

               %  131 - Tree 3, DWD, 1024 genes unmasked, DWD Carc vs. NormCol
               %  132 - Tree 3, DWD, 512 genes unmasked, DWD Carc vs. NormCol
               %  133 - Tree 3, DWD, 256 genes unmasked, DWD Carc vs. NormCol
               %  134 - Tree 3, DWD, 128 genes unmasked, DWD Carc vs. NormCol
               %  135 - Tree 3, DWD, 64 genes unmasked, DWD Carc vs. NormCol
               %  136 - Tree 3, DWD, 32 genes unmasked, DWD Carc vs. NormCol
               %  137 - Tree 3, DWD, 16 genes unmasked, DWD Carc vs. NormCol
               %  138 - Tree 3, DWD, 8 genes unmasked, DWD Carc vs. NormCol
               %  139 - Tree 3, DWD, 4 genes unmasked, DWD Carc vs. NormCol

               %  140 - Tree 3, Gene Focussing, DWD Carc vs. NormCol, Compute Indices and save
               %  141 - Tree 3, Gene Focussing, DWD Carc vs. NormCol, Plots of Cluster Indices
               %  142 - Tree 3, Gene Focussing, DWD Carc vs. NormCol, Plots of Mean Diffs
               %  143 - Tree 3, Gene Focussing, DWD Carc vs. NormCol, Plots of MDPs

               %  151 - Tree 3, DWD, 1024 genes unmasked, DWD Carc vs. NormCol, recomputed axes
               %  152 - Tree 3, DWD, 512 genes unmasked, DWD Carc vs. NormCol, recomputed axes
               %  153 - Tree 3, DWD, 256 genes unmasked, DWD Carc vs. NormCol, recomputed axes
               %  154 - Tree 3, DWD, 128 genes unmasked, DWD Carc vs. NormCol, recomputed axes
               %  155 - Tree 3, DWD, 64 genes unmasked, DWD Carc vs. NormCol, recomputed axes
               %  156 - Tree 3, DWD, 32 genes unmasked, DWD Carc vs. NormCol, recomputed axes
               %  157 - Tree 3, DWD, 16 genes unmasked, DWD Carc vs. NormCol, recomputed axes
               %  158 - Tree 3, DWD, 8 genes unmasked, DWD Carc vs. NormCol, recomputed axes
               %  159 - Tree 3, DWD, 4 genes unmasked, DWD Carc vs. NormCol, recomputed axes

               %  201 - Indiv - All, DWD, 1024 genes unmasked, DWD Carc vs. Rest
               %  202 - Indiv - All, DWD, 512 genes unmasked, DWD Carc vs. Rest
               %  203 - Indiv - All, DWD, 256 genes unmasked, DWD Carc vs. Rest
               %  204 - Indiv - All, DWD, 128 genes unmasked, DWD Carc vs. Rest
               %  205 - Indiv - All, DWD, 64 genes unmasked, DWD Carc vs. Rest
               %  206 - Indiv - All, DWD, 32 genes unmasked, DWD Carc vs. Rest
               %  207 - Indiv - All, DWD, 16 genes unmasked, DWD Carc vs. Rest
               %  208 - Indiv - All, DWD, 8 genes unmasked, DWD Carc vs. Rest
               %  209 - Indiv - All, DWD, 4 genes unmasked, DWD Carc vs. Rest

               %  211 - Indiv - All, DWD, 1024 genes unmasked, DWD SC vs. Rest
               %  212 - Indiv - All, DWD, 512 genes unmasked, DWD SC vs. Rest
               %  213 - Indiv - All, DWD, 256 genes unmasked, DWD SC vs. Rest
               %  214 - Indiv - All, DWD, 128 genes unmasked, DWD SC vs. Rest
               %  215 - Indiv - All, DWD, 64 genes unmasked, DWD SC vs. Rest
               %  216 - Indiv - All, DWD, 32 genes unmasked, DWD SC vs. Rest
               %  217 - Indiv - All, DWD, 16 genes unmasked, DWD SC vs. Rest
               %  218 - Indiv - All, DWD, 8 genes unmasked, DWD SC vs. Rest
               %  219 - Indiv - All, DWD, 4 genes unmasked, DWD SC vs. Rest

               %  221 - Indiv - All, DWD, 1024 genes unmasked, DWD Norm vs. Rest
               %  222 - Indiv - All, DWD, 512 genes unmasked, DWD Norm vs. Rest
               %  223 - Indiv - All, DWD, 256 genes unmasked, DWD Norm vs. Rest
               %  224 - Indiv - All, DWD, 128 genes unmasked, DWD Norm vs. Rest
               %  225 - Indiv - All, DWD, 64 genes unmasked, DWD Norm vs. Rest
               %  226 - Indiv - All, DWD, 32 genes unmasked, DWD Norm vs. Rest
               %  227 - Indiv - All, DWD, 16 genes unmasked, DWD Norm vs. Rest
               %  228 - Indiv - All, DWD, 8 genes unmasked, DWD Norm vs. Rest
               %  229 - Indiv - All, DWD, 4 genes unmasked, DWD Norm vs. Rest

               %  231 - Indiv - All, DWD, 1024 genes unmasked, DWD Col vs. Rest
               %  232 - Indiv - All, DWD, 512 genes unmasked, DWD Col vs. Rest
               %  233 - Indiv - All, DWD, 256 genes unmasked, DWD Col vs. Rest
               %  234 - Indiv - All, DWD, 128 genes unmasked, DWD Col vs. Rest
               %  235 - Indiv - All, DWD, 64 genes unmasked, DWD Col vs. Rest
               %  236 - Indiv - All, DWD, 32 genes unmasked, DWD Col vs. Rest
               %  237 - Indiv - All, DWD, 16 genes unmasked, DWD Col vs. Rest
               %  238 - Indiv - All, DWD, 8 genes unmasked, DWD Col vs. Rest
               %  239 - Indiv - All, DWD, 4 genes unmasked, DWD Col vs. Rest

               %  240 - Indiv - All, Gene Focussing, DWD Carc vs. Rest, Compute Indices and save
               %  241 - Indiv - All, Gene Focussing, DWD Carc vs. Rest, Plots of Cluster Indices
               %  242 - Indiv - All, Gene Focussing, DWD Carc vs. Rest, Plots of Mean Diffs
               %  243 - Indiv - All, Gene Focussing, DWD Carc vs. Rest, Plots of MDPs

               %  250 - Indiv - All, Gene Focussing, DWD SC vs. Rest, Compute Indices and save
               %  251 - Indiv - All, Gene Focussing, DWD SC vs. Rest, Plots of Cluster Indices
               %  252 - Indiv - All, Gene Focussing, DWD SC vs. Rest, Plots of Mean Diffs
               %  253 - Indiv - All, Gene Focussing, DWD SC vs. Rest, Plots of MDPs

               %  260 - Indiv - All, Gene Focussing, DWD Norm vs. Rest, Compute Indices and save
               %  261 - Indiv - All, Gene Focussing, DWD Norm vs. Rest, Plots of Cluster Indices
               %  262 - Indiv - All, Gene Focussing, DWD Norm vs. Rest, Plots of Mean Diffs
               %  263 - Indiv - All, Gene Focussing, DWD Norm vs. Rest, Plots of MDPs

               %  270 - Indiv - All, Gene Focussing, DWD Col vs. Rest, Compute Indices and save
               %  271 - Indiv - All, Gene Focussing, DWD Col vs. Rest, Plots of Cluster Indices
               %  272 - Indiv - All, Gene Focussing, DWD Col vs. Rest, Plots of Mean Diffs
               %  273 - Indiv - All, Gene Focussing, DWD Col vs. Rest, Plots of MDPs

               %  301 - Remove Top 3 PCs, & PCA plot
               %  302 - Remove Top 3 PCs, & DWD Indiv vs. All plot

               %  311 - Remove Tree 1 subspace, & PCA plot
               %  312 - Remove Tree 1 subspace, & DWD Indiv vs. All plot

               %  321 - Remove Tree 2 subspace, & PCA plot
               %  322 - Remove Tree 2 subspace, & DWD Indiv vs. All plot

               %  331 - Remove Tree 3 subspace, & PCA plot
               %  332 - Remove Tree 3 subspace, & DWD Indiv vs. All plot



%  Check Matlab version
%
if version('-release') < 14 ;

  disp('!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!') ;
  disp('!!!   Warning from NeilLung1.m:        !!!') ;
  disp('!!!   A release of Matlab before 14    !!!') ;
  disp('!!!   is being used.  This can give    !!!') ;
  disp('!!!   unexpected results with data     !!!') ;
  disp('!!!   writing and reading              !!!') ;
  disp('!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!') ;
  pauseSM ;  

end ;



if ipart == 0 ;    %  Do main data read, check inputs, and save as .mat

  %  Read in Data Values
  %
  infilestr = 'lung4group-nolabel.txt' ;
  mdata = dlmread(infilestr) ;
  disp('  Finished Reading Main Data') ;

  d = size(mdata,1)
  n = size(mdata,2)


  %  Read in Header Values  (Case Names)
  %
  infilestr = 'lung4group-labelonly.txt' ;
  TextInput = textread(infilestr,'%s','delimiter','\t') ;
  disp('  Finished Reading Case Names') ;

  CaseNames = TextInput(1:n) ;
  ClassLabels = TextInput(n+1:2*n) ;


  %  Check Inputs
  %
  disp(' ') ;
  disp('Check upper left data value = 0.104613460185975') ;
  mdata(1,1)

  disp(' ') ;
  disp('Check upper right data value = -1.94102152862633') ;
  mdata(1,end)

  disp(' ') ;
  disp('Check lower left data value = 0.0141469045589287') ;
  mdata(end,1)

  disp(' ') ;
  disp('Check lower right data value = -0.156888709178854') ;
  mdata(end,end)


  disp(' ') ;
  disp('Check 1st Case Name = CL2001031609AA.CEL') ;
  CaseNames{1} 
  
  disp(' ') ;
  disp('Check last Case Name = CL2001032720AA.CEL') ;
  CaseNames{end} 


  disp(' ') ;
  disp('Check 1st CLass Label = Carcinoid') ;
  ClassLabels{1} 
  
  disp(' ') ;
  disp('Check last Class Label = SmallCell') ;
  ClassLabels{end} 


  disp(' ') ;
  disp('Check this is 0 (correct number of Case Names)') ;
  abs(length(CaseNames) - size(mdata,2))

  disp(' ') ;
  disp('Check this is 0 (correct number of Class Labels)') ;
  abs(length(ClassLabels) - size(mdata,2))




  %  Save results in .mat file
  %
  save NeilLung1.mat mdata CaseNames ClassLabels ;
      %  save these variables

  disp('Finished data save') ;




else ;    %  then do other parts of the program


  load NeilLung1 ;
      %  loads: 
      %      mdata 
      %      CaseNames
      %      ClassLabels
  disp('  Finished Data Load') ;

  d = size(mdata,1) ;
  n = size(mdata,2) ;


  figure(1) ;
  clf ;

  %  Create generic icolor & markerstr for Classes
  %
  icolor = [] ;
  markerstr = [] ;
  CarcinoidFlag = zeros(1,length(ClassLabels)) ;
  ColonFlag = zeros(1,length(ClassLabels)) ;
  NormalFlag = zeros(1,length(ClassLabels)) ;
  SmallCellFlag = zeros(1,length(ClassLabels)) ;
  for i = 1:length(ClassLabels) ;

    if strcmp(ClassLabels{i},'Carcinoid') ;
      icolor = [icolor; [1 0 1]] ;
      markerstr = strvcat(markerstr,'^') ;
      CarcinoidFlag(i) = 1 ;
    elseif strcmp(ClassLabels{i},'Colon') ;
      icolor = [icolor; [1 0 0]] ;
      markerstr = strvcat(markerstr,'s') ;
      ColonFlag(i) = 1 ;
    elseif strcmp(ClassLabels{i},'Normal') ;
      icolor = [icolor; [0 1 0]] ;
      markerstr = strvcat(markerstr,'o') ;
      NormalFlag(i) = 1 ;
    elseif strcmp(ClassLabels{i},'SmallCell') ;
      icolor = [icolor; [0 0 1]] ;
      markerstr = strvcat(markerstr,'x') ;
      SmallCellFlag(i) = 1 ;
    else ;
      disp('!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!') ;
      disp('!!!   Error from NeilLung1.m:         !!!') ;
      disp('!!!   Unmatched Class                 !!!')
      disp('!!!   Quit here and fix Class Labels  !!!')
      disp('!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!') ;
      pauseSM ;
    end ;

  end ;

  CarcinoidFlag = logical(CarcinoidFlag) ;
  ColonFlag = logical(ColonFlag) ;
  NormalFlag = logical(NormalFlag) ;
  SmallCellFlag = logical(SmallCellFlag) ;


  %  Set up generic legend
  %
  legendcellstr = {{'Carcinoid' 'Colon' 'Normal' 'SmallCell'}} ;
  mlegendcolor = [[1 0 1]; [1 0 0]; [0 1 0]; [0 0 1]] ;



  if ipart == 1 ;    %  PCA plot of combined data (Matlab default colors)

    titlecellstr = {{'Standard PCA View' 'Neil Hayes Lung data' 'Carc, Col, Norm, SC only' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'PC 1' 'PC 2' 'PC 3' 'PC 4'}} ;

    savestr = ['NeilLung1ip00' num2str(ipart) 'PCASPdefcol'] ;

    paramstruct = struct('npcadiradd',4, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;
 

    scatplotSM(mdata,[],paramstruct) ;


  elseif ipart == 2 ;    %  PCA plot of combined data (all black)

    titlecellstr = {{'Standard PCA View' 'Neil Hayes Lung data' 'Carc, Col, Norm, SC only' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'PC 1' 'PC 2' 'PC 3' 'PC 4'}} ;

    savestr = ['NeilLung1ip00' num2str(ipart) 'PCASPblack'] ;

    paramstruct = struct('icolor','k', ...
                         'npcadiradd',4, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;
 

    scatplotSM(mdata,[],paramstruct) ;


  elseif ipart == 3 ;    %  PCA plot of combined data (Class colors)

    titlecellstr = {{'Standard PCA View' 'Neil Hayes Lung data' 'Carc, Col, Norm, SC only' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'PC 1' 'PC 2' 'PC 3' 'PC 4'}} ;

    savestr = ['NeilLung1ip00' num2str(ipart) 'PCASPclasscol'] ;

    paramstruct = struct('icolor',icolor, ...
                         'npcadiradd',4, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;
 

    scatplotSM(mdata,[],paramstruct) ;


  elseif ipart == 4 ;    %  PCA plot of combined data (Class colors and symbols)

    titlecellstr = {{'Standard PCA View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'PC 1' 'PC 2' 'PC 3' 'PC 4'}} ;

    savestr = ['NeilLung1ip00' num2str(ipart) 'PCASPclasscolsym'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'npcadiradd',4, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;
 

    scatplotSM(mdata,[],paramstruct) ;


  elseif ipart == 5 ;    %  PCA plot of combined data (Class symbols only)

    titlecellstr = {{'Standard PCA View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'PC 1' 'PC 2' 'PC 3' 'PC 4'}} ;

    savestr = ['NeilLung1ip00' num2str(ipart) 'PCASPclasssym'] ;

    paramstruct = struct('icolor','k', ...
                         'markerstr',markerstr, ...
                         'npcadiradd',4, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;
 

    scatplotSM(mdata,[],paramstruct) ;




  elseif ipart == 11 ;     %  Proj on M-Carc - M-Col, M-Norm - M-SC, & ortho PC 

    dir1 = mean(mdata(:,CarcinoidFlag),2) - mean(mdata(:,ColonFlag),2) ;
    dir1 = dir1 / norm(dir1) ; 
    dir2 = mean(mdata(:,NormalFlag),2) - mean(mdata(:,SmallCellFlag),2) ;
    dir2 = dir2 / norm(dir2) ; 
    mdir = [dir1 dir2] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'MD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'MD Carc - Col' 'MD Norm-SC' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'MDCarc-ColMDNorm-SC'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;
 

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 12 ;     %  Proj on M-Carc - M-Norm, M-Col - M-SC, & ortho PC 

    dir1 = mean(mdata(:,CarcinoidFlag),2) - mean(mdata(:,NormalFlag),2) ;
    dir1 = dir1 / norm(dir1) ; 
    dir2 = mean(mdata(:,ColonFlag),2) - mean(mdata(:,SmallCellFlag),2) ;
    dir2 = dir2 / norm(dir2) ; 
    mdir = [dir1 dir2] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'MD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'MD Carc - Norm' 'MD Col-SC' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'MDCarc-NormMDCol-SC'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;
 

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 13 ;     %  Proj on M-Carc - M-SC, M-Norm - M-Col, & ortho PC 

    dir1 = mean(mdata(:,CarcinoidFlag),2) - mean(mdata(:,SmallCellFlag),2) ;
    dir1 = dir1 / norm(dir1) ; 
    dir2 = mean(mdata(:,NormalFlag),2) - mean(mdata(:,ColonFlag),2) ;
    dir2 = dir2 / norm(dir2) ; 
    mdir = [dir1 dir2] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'MD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'MD Carc - SC' 'MD Norm-Col' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'MDCarc-SCMDNorm-Col'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;
 

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 14 ;     %  Proj on M-(Carc&SC) - M-(Norm&Col), M-(Carc&Norm) - M-(SC&Col), & ortho PC 

    dir1 = mean(mdata(:,CarcinoidFlag|SmallCellFlag),2) - mean(mdata(:,NormalFlag|ColonFlag),2) ;
    dir1 = dir1 / norm(dir1) ; 
    dir2 = mean(mdata(:,CarcinoidFlag|NormalFlag),2) - mean(mdata(:,SmallCellFlag|ColonFlag),2) ;
    dir2 = dir2 / norm(dir2) ; 
    mdir = [dir1 dir2] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'MD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'MD CarcSC-NormCol ' 'MD CarcNorm-SCCol' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'MDCarcSC-NormColMDCarcNorm-SCCol'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;
 

  elseif ipart == 15 ;     %  Proj on M-(Carc&Norm) - M-(SC&Col), M-(Carc&Col) - M-(Norm&SC), & ortho PC 

    dir1 = mean(mdata(:,CarcinoidFlag|NormalFlag),2) - mean(mdata(:,SmallCellFlag|ColonFlag),2) ;
    dir1 = dir1 / norm(dir1) ; 
    dir2 = mean(mdata(:,CarcinoidFlag|ColonFlag),2) - mean(mdata(:,NormalFlag|SmallCellFlag),2) ;
    dir2 = dir2 / norm(dir2) ; 
    mdir = [dir1 dir2] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'MD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'MD CarcNorm-SCCol ' 'MD CarcCol-NormSC' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'MDCarcNorm-NormSCMDCarcNorm-SCCol'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 16 ;     %  Proj on M-(Carc&Col) - M-(SC&Norm), M-(Carc&SC) - M-(Norm&Col), & ortho PC 

    dir1 = mean(mdata(:,CarcinoidFlag|ColonFlag),2) - mean(mdata(:,SmallCellFlag|NormalFlag),2) ;
    dir1 = dir1 / norm(dir1) ; 
    dir2 = mean(mdata(:,CarcinoidFlag|SmallCellFlag),2) - mean(mdata(:,NormalFlag|ColonFlag),2) ;
    dir2 = dir2 / norm(dir2) ; 
    mdir = [dir1 dir2] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'MD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'MD CarcCol-SCNorm ' 'MD CarcSC-NormCol' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'MDCarcCol-SCNormMDCarcSC-NormCol'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 17 ;     %  Proj on MD Indiv vs. all 

    dir1 = mean(mdata(:,CarcinoidFlag),2) - mean(mdata(:,SmallCellFlag|NormalFlag|ColonFlag),2) ;
    dir1 = dir1 / norm(dir1) ; 
    dir2 = mean(mdata(:,SmallCellFlag),2) - mean(mdata(:,CarcinoidFlag|NormalFlag|ColonFlag),2) ;
    dir2 = dir2 / norm(dir2) ; 
    dir3 = mean(mdata(:,NormalFlag),2) - mean(mdata(:,CarcinoidFlag|SmallCellFlag|ColonFlag),2) ;
    dir3 = dir3 / norm(dir3) ; 
    dir4 = mean(mdata(:,ColonFlag),2) - mean(mdata(:,NormalFlag|CarcinoidFlag|SmallCellFlag),2) ;
    dir4 = dir4 / norm(dir4) ; 
    mdir = [dir1 dir2 dir3 dir4] ;

    titlecellstr = {{'MD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'MD Carc-all ' 'MD SC-all' 'MD Norm-all' 'MD Col-all'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'MDIndiv-all'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;



  elseif ipart == 21 ;     %  Proj on DWD Carc-Col, DWD Norm-SC, & ortho PC 

    mdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,SmallCellFlag))] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'DWD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD Carc - Col' 'DWD Norm-SC' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDCarc-ColDWDNorm-SC'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;
 

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 22 ;     %  Proj on DWD Carc-Norm, DWD Col-SC, & ortho PC 

    mdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,NormalFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,SmallCellFlag))] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'DWD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD Carc - Norm' 'DWD Col-SC' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDCarc-NormDWDCol-SC'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;
 

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 23 ;     %  Proj on DWD Carc-SC, DWD Norm-Col, & ortho PC 

    mdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,ColonFlag))] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'DWD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD Carc - SC' 'DWD Norm-Col' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDCarc-SCDWDNorm-Col'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;
 

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 24 ;     %  Proj on DWD (Carc&SC)-(Norm&Col), DWD (Carc&Norm)-(SC&Col), & ortho PC 

    mdir = DWD1SM(mdata(:,CarcinoidFlag|SmallCellFlag),mdata(:,NormalFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag|NormalFlag),mdata(:,SmallCellFlag|ColonFlag))] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'DWD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcSC-NormCol ' 'DWD CarcNorm-SCCol' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDCarcSC-NormColDWDCarcNorm-SCCol'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;
 

  elseif ipart == 25 ;     %  Proj on DWD (Carc&Norm)-(SC&Col), DWD (Carc&Col)-(Norm&SC), & ortho PC 

    mdir = DWD1SM(mdata(:,CarcinoidFlag|NormalFlag),mdata(:,SmallCellFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag|ColonFlag),mdata(:,NormalFlag|SmallCellFlag))] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'DWD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcNorm-SCCol ' 'DWD CarcCol-NormSC' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDCarcNorm-NormSCDWDCarcNorm-SCCol'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 26 ;     %  Proj on DWD (Carc&Col)-(SC&Norm), DWD (Carc&SC)-(Norm&Col), & ortho PC 

    mdir = DWD1SM(mdata(:,CarcinoidFlag|ColonFlag),mdata(:,SmallCellFlag|NormalFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag|SmallCellFlag),mdata(:,NormalFlag|ColonFlag))] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'DWD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcCol-SCNorm ' 'DWD CarcSC-NormCol' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDCarcCol-SCNormDWDCarcSC-NormCol'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 27 ;     %  Proj on DWD Indiv vs. all 

    mdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,CarcinoidFlag|SmallCellFlag|ColonFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,NormalFlag|CarcinoidFlag|SmallCellFlag))] ;


    titlecellstr = {{'DWD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD Carc-all ' 'DWD SC-all' 'DWD Norm-all' 'DWD Col-all'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDIndiv-all'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;



  elseif ipart == 31 ;     %  Proj on SVM Carc-Col, SVM Norm-SC, & ortho PC 

    mdir = SVM1SM(mdata(:,CarcinoidFlag),mdata(:,ColonFlag)) ;
    mdir = [mdir SVM1SM(mdata(:,NormalFlag),mdata(:,SmallCellFlag))] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'SVM View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'SVM Carc - Col' 'SVM Norm-SC' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'SVMCarc-ColSVMNorm-SC'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;
 

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 32 ;     %  Proj on SVM Carc-Norm, SVM Col-SC, & ortho PC 

    mdir = SVM1SM(mdata(:,CarcinoidFlag),mdata(:,NormalFlag)) ;
    mdir = [mdir SVM1SM(mdata(:,ColonFlag),mdata(:,SmallCellFlag))] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'SVM View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'SVM Carc - Norm' 'SVM Col-SC' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'SVMCarc-NormSVMCol-SC'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;
 

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 33 ;     %  Proj on SVM Carc-SC, SVM Norm-Col, & ortho PC 

    mdir = SVM1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag)) ;
    mdir = [mdir SVM1SM(mdata(:,NormalFlag),mdata(:,ColonFlag))] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'SVM View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'SVM Carc - SC' 'SVM Norm-Col' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'SVMCarc-SCSVMNorm-Col'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;
 

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 34 ;     %  Proj on SVM (Carc&SC)-(Norm&Col), SVM (Carc&Norm)-(SC&Col), & ortho PC 

    mdir = SVM1SM(mdata(:,CarcinoidFlag|SmallCellFlag),mdata(:,NormalFlag|ColonFlag)) ;
    mdir = [mdir SVM1SM(mdata(:,CarcinoidFlag|NormalFlag),mdata(:,SmallCellFlag|ColonFlag))] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'SVM View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'SVM CarcSC-NormCol ' 'SVM CarcNorm-SCCol' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'SVMCarcSC-NormColSVMCarcNorm-SCCol'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;
 

  elseif ipart == 35 ;     %  Proj on SVM (Carc&Norm)-(SC&Col), SVM (Carc&Col)-(Norm&SC), & ortho PC 

    mdir = SVM1SM(mdata(:,CarcinoidFlag|NormalFlag),mdata(:,SmallCellFlag|ColonFlag)) ;
    mdir = [mdir SVM1SM(mdata(:,CarcinoidFlag|ColonFlag),mdata(:,NormalFlag|SmallCellFlag))] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'SVM View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'SVM CarcNorm-SCCol ' 'SVM CarcCol-NormSC' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'SVMCarcNorm-NormSCSVMCarcNorm-SCCol'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 36 ;     %  Proj on SVM (Carc&Col)-(SC&Norm), SVM (Carc&SC)-(Norm&Col), & ortho PC 

    mdir = SVM1SM(mdata(:,CarcinoidFlag|ColonFlag),mdata(:,SmallCellFlag|NormalFlag)) ;
    mdir = [mdir SVM1SM(mdata(:,CarcinoidFlag|SmallCellFlag),mdata(:,NormalFlag|ColonFlag))] ;

    proj2md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 2-d subspace generated by mdir
    mdataortho = mdata - proj2md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj2md = [] ;
        %  to save space
    paramstruct = struct('npc',2,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;


    titlecellstr = {{'SVM View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'SVM CarcCol-SCNorm ' 'SVM CarcSC-NormCol' 'ortho PC 1' 'ortho PC 2'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'SVMCarcCol-SCNormSVMCarcSC-NormCol'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 37 ;     %  Proj on SVM Indiv vs. all 

    mdir = SVM1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
    mdir = [mdir SVM1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag))] ;
    mdir = [mdir SVM1SM(mdata(:,NormalFlag),mdata(:,CarcinoidFlag|SmallCellFlag|ColonFlag))] ;

    titlecellstr = {{'Tree I' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'SVM Carc-all ' 'SVM SC-all' 'SVM Norm-all' 'SVM Col-all'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'SVMIndiv-all'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 41 ;     %  Tree1, level 0
                           %  Proj on DWD (Carc&Norm)-(SC&Col), DWD Norm-Carc, DWD SC&Col, & ortho PC 
                           %  All Black

    mdir = DWD1SM(mdata(:,CarcinoidFlag|NormalFlag),mdata(:,SmallCellFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,CarcinoidFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,SmallCellFlag),mdata(:,ColonFlag))] ;

    proj3md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 3-d subspace generated by mdir
    mdataortho = mdata - proj3md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj3md = [] ;
        %  to save space
    paramstruct = struct('npc',1,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;

    titlecellstr = {{'DWD View, Tree I - 0' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcNorm-SCCol ' 'DWD Norm-Carc' 'DWD SC-Col' 'ortho PC 1'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDtreeI-lev0'] ;

    paramstruct = struct('icolor','k', ...
                         'markerstr',markerstr, ...
                         'isubpopkde',0, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 42 ;     %  Tree1, level 1
                           %  Proj on DWD (Carc&Norm)-(SC&Col), DWD Norm-Carc, DWD SC&Col, & ortho PC 
                           %  Level 1 Col

    mdir = DWD1SM(mdata(:,CarcinoidFlag|NormalFlag),mdata(:,SmallCellFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,CarcinoidFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,SmallCellFlag),mdata(:,ColonFlag))] ;

    proj3md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 3-d subspace generated by mdir
    mdataortho = mdata - proj3md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj3md = [] ;
        %  to save space
    paramstruct = struct('npc',1,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;

    icolornew = [] ;
    for i = 1:length(ClassLabels) ;

      if strcmp(ClassLabels{i},'Carcinoid') ;
        icolornew = [icolornew; [0.5 0.5 0.5]] ;
      elseif strcmp(ClassLabels{i},'Colon') ;
        icolornew = [icolornew; [0.5 0 0.5]] ;
      elseif strcmp(ClassLabels{i},'Normal') ;
        icolornew = [icolornew; [0.5 0.5 0.5]] ;
      elseif strcmp(ClassLabels{i},'SmallCell') ;
        icolornew = [icolornew; [0.5 0 0.5]] ;
      else ;
        disp('!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!') ;
        disp('!!!   Error from NeilLung1.m:         !!!') ;
        disp('!!!   Unmatched Class                 !!!')
        disp('!!!   Quit here and fix Class Labels  !!!')
        disp('!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!') ;
        pauseSM ;
      end ;
    end ;

    mlegendcolornew = [[0.5 0.5 0.5]; [0.5 0 0.5]; [0.5 0.5 0.5]; [0.5 0 0.5]] ;

    titlecellstr = {{'DWD View, Tree I - 1' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcNorm-SCCol ' 'DWD Norm-Carc' 'DWD SC-Col' 'ortho PC 1'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDtreeI-lev1'] ;

    paramstruct = struct('icolor',icolornew, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolornew, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 43 ;     %  Tree1, level 2
                           %  Proj on DWD (Carc&Norm)-(SC&Col), DWD Norm-Carc, DWD SC&Col, & ortho PC 
                           %  All Col

    mdir = DWD1SM(mdata(:,CarcinoidFlag|NormalFlag),mdata(:,SmallCellFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,CarcinoidFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,SmallCellFlag),mdata(:,ColonFlag))] ;

    proj3md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 3-d subspace generated by mdir
    mdataortho = mdata - proj3md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj3md = [] ;
        %  to save space
    paramstruct = struct('npc',1,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;

    titlecellstr = {{'DWD View, Tree I - 2' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcNorm-SCCol ' 'DWD Norm-Carc' 'DWD SC-Col' 'ortho PC 1'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDtreeI-lev2'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 44 ;     %  Tree2, level 0
                           %  Proj on DWD (Carc&SC)-(Norm&Col), DWD Carc-SC, DWD Norm-Col, & ortho PC 
                           %  All Black

    mdir = DWD1SM(mdata(:,CarcinoidFlag|SmallCellFlag),mdata(:,NormalFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,ColonFlag))] ;

    proj3md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 3-d subspace generated by mdir
    mdataortho = mdata - proj3md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj3md = [] ;
        %  to save space
    paramstruct = struct('npc',1,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;

    titlecellstr = {{'DWD View, Tree II - 0' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcSC-NormCol ' 'DWD Carc-SC' 'DWD Norm-Col' 'ortho PC 1'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDtreeII-lev0'] ;

    paramstruct = struct('icolor','k', ...
                         'markerstr',markerstr, ...
                         'isubpopkde',0, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 45 ;     %  Tree2, level 1
                           %  Proj on DWD (Carc&SC)-(Norm&Col), DWD Carc-SC, DWD Norm-Col, & ortho PC 
                           %  Level 1 Col

    mdir = DWD1SM(mdata(:,CarcinoidFlag|SmallCellFlag),mdata(:,NormalFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,ColonFlag))] ;

    proj3md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 3-d subspace generated by mdir
    mdataortho = mdata - proj3md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj3md = [] ;
        %  to save space
    paramstruct = struct('npc',1,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;

    icolornew = [] ;
    for i = 1:length(ClassLabels) ;

      if strcmp(ClassLabels{i},'Carcinoid') ;
        icolornew = [icolornew; [0.5 0 1]] ;
      elseif strcmp(ClassLabels{i},'Colon') ;
        icolornew = [icolornew; [0.5 0.5 0]] ;
      elseif strcmp(ClassLabels{i},'Normal') ;
        icolornew = [icolornew; [0.5 0.5 0]] ;
      elseif strcmp(ClassLabels{i},'SmallCell') ;
        icolornew = [icolornew; [0.5 0 1]] ;
      else ;
        disp('!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!') ;
        disp('!!!   Error from NeilLung1.m:         !!!') ;
        disp('!!!   Unmatched Class                 !!!')
        disp('!!!   Quit here and fix Class Labels  !!!')
        disp('!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!') ;
        pauseSM ;
      end ;
    end ;

    mlegendcolornew = [[0.5 0 1]; [0.5 0.5 0]; [0.5 0.5 0]; [0.5 0 1]] ;

    titlecellstr = {{'DWD View, Tree II - 1' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcSC-NormCol ' 'DWD Carc-SC' 'DWD Norm-Col' 'ortho PC 1'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDtreeII-lev1'] ;

    paramstruct = struct('icolor',icolornew, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolornew, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 46 ;     %  Tree2, level 2
                           %  Proj on DWD (Carc&SC)-(Norm&Col), DWD Carc-SC, DWD Norm-Col, & ortho PC 
                           %  All Col

    mdir = DWD1SM(mdata(:,CarcinoidFlag|SmallCellFlag),mdata(:,NormalFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,ColonFlag))] ;

    proj3md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 3-d subspace generated by mdir
    mdataortho = mdata - proj3md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj3md = [] ;
        %  to save space
    paramstruct = struct('npc',1,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;

    titlecellstr = {{'DWD View, Tree II - 2' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcSC-NormCol ' 'DWD Carc-SC' 'DWD Norm-Col' 'ortho PC 1'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDtreeII-lev2'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 47 ;     %  Tree3, level 0
                           %  Proj on DWD SC-(Carc&Norm&Col), DWD Carc-(Col&Norm), DWD Norm-Col, & ortho PC 
                           %  All Black

    mdir = DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag),mdata(:,ColonFlag|NormalFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,NormalFlag))] ;

    proj3md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 3-d subspace generated by mdir
    mdataortho = mdata - proj3md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj3md = [] ;
        %  to save space
    paramstruct = struct('npc',1,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;

    titlecellstr = {{'DWD View, Tree III - 0' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcNormCol-SC' 'DWD Carc-ColNorm' 'DWD Norm-Col' 'ortho PC 1'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDtreeIII-lev0'] ;

    paramstruct = struct('icolor','k', ...
                         'markerstr',markerstr, ...
                         'isubpopkde',0, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 48 ;     %  Tree3, level 1
                           %  Proj on DWD SC-(Carc&Norm&Col), DWD Carc-(Col&Norm), DWD Norm-Col, & ortho PC 
                           %  Level 1 Col

    mdir = DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag),mdata(:,ColonFlag|NormalFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,NormalFlag))] ;

    proj3md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 3-d subspace generated by mdir
    mdataortho = mdata - proj3md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj3md = [] ;
        %  to save space
    paramstruct = struct('npc',1,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;

    icolornew = [] ;
    for i = 1:length(ClassLabels) ;

      if strcmp(ClassLabels{i},'Carcinoid') ;
        icolornew = [icolornew; [0.67 0.33 0.33]] ;
      elseif strcmp(ClassLabels{i},'Colon') ;
        icolornew = [icolornew; [0.67 0.33 0.33]] ;
      elseif strcmp(ClassLabels{i},'Normal') ;
        icolornew = [icolornew; [0.67 0.33 0.33]] ;
      elseif strcmp(ClassLabels{i},'SmallCell') ;
        icolornew = [icolornew; [0 0 1]] ;
      else ;
        disp('!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!') ;
        disp('!!!   Error from NeilLung1.m:         !!!') ;
        disp('!!!   Unmatched Class                 !!!')
        disp('!!!   Quit here and fix Class Labels  !!!')
        disp('!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!') ;
        pauseSM ;
      end ;
    end ;

%    if strcmp(ClassLabels{i},'Carcinoid') ;
%      icolor = [icolor; [1 0 1]] ;
%    elseif strcmp(ClassLabels{i},'Colon') ;
%      icolor = [icolor; [1 0 0]] ;
%    elseif strcmp(ClassLabels{i},'Normal') ;
%      icolor = [icolor; [0 1 0]] ;
%    elseif strcmp(ClassLabels{i},'SmallCell') ;
%      icolor = [icolor; [0 0 1]] ;
%    end ;

    mlegendcolornew = [[0.67 0.33 0.33]; [0.67 0.33 0.33]; [0.67 0.33 0.33]; [0 0 1]] ;

    titlecellstr = {{'DWD View, Tree III - 1' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcNormCol-SC' 'DWD Carc-ColNorm' 'DWD Norm-Col' 'ortho PC 1'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDtreeIII-lev1'] ;

    paramstruct = struct('icolor',icolornew, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolornew, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 49 ;     %  Tree3, level 2
                           %  Proj on DWD SC-(Carc&Norm&Col), DWD Carc-(Col&Norm), DWD Norm-Col, & ortho PC 
                           %  Level 2 Col

    mdir = DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag),mdata(:,ColonFlag|NormalFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,NormalFlag))] ;

    proj3md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 3-d subspace generated by mdir
    mdataortho = mdata - proj3md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj3md = [] ;
        %  to save space
    paramstruct = struct('npc',1,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;

    icolornew = [] ;
    for i = 1:length(ClassLabels) ;

      if strcmp(ClassLabels{i},'Carcinoid') ;
        icolornew = [icolornew; [1 0 1]] ;
      elseif strcmp(ClassLabels{i},'Colon') ;
        icolornew = [icolornew; [0.5 0.5 0]] ;
      elseif strcmp(ClassLabels{i},'Normal') ;
        icolornew = [icolornew; [0.5 0.5 0]] ;
      elseif strcmp(ClassLabels{i},'SmallCell') ;
        icolornew = [icolornew; [0 0 1]] ;
      else ;
        disp('!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!') ;
        disp('!!!   Error from NeilLung1.m:         !!!') ;
        disp('!!!   Unmatched Class                 !!!')
        disp('!!!   Quit here and fix Class Labels  !!!')
        disp('!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!') ;
        pauseSM ;
      end ;
    end ;

%    if strcmp(ClassLabels{i},'Carcinoid') ;
%      icolor = [icolor; [1 0 1]] ;
%    elseif strcmp(ClassLabels{i},'Colon') ;
%      icolor = [icolor; [1 0 0]] ;
%    elseif strcmp(ClassLabels{i},'Normal') ;
%      icolor = [icolor; [0 1 0]] ;
%    elseif strcmp(ClassLabels{i},'SmallCell') ;
%      icolor = [icolor; [0 0 1]] ;
%    end ;

    mlegendcolornew = [[1 0 1]; [0.5 0.5 0]; [0.5 0.5 0]; [0 0 1]] ;

    titlecellstr = {{'DWD View, Tree III - 2' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcNormCol-SC' 'DWD Carc-ColNorm' 'DWD Norm-Col' 'ortho PC 1'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDtreeIII-lev2'] ;

    paramstruct = struct('icolor',icolornew, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolornew, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif ipart == 50 ;     %  Tree3, level 3
                           %  Proj on DWD SC-(Carc&Norm&Col), DWD Carc-(Col&Norm), DWD Norm-Col, & ortho PC 
                           %  Level 3 Col

    mdir = DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag),mdata(:,ColonFlag|NormalFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,NormalFlag))] ;

    proj3md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 3-d subspace generated by mdir
    mdataortho = mdata - proj3md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj3md = [] ;
        %  to save space
    paramstruct = struct('npc',1,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;

    titlecellstr = {{'DWD View, Tree III - 3' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcNormCol-SC' 'DWD Carc-ColNorm' 'DWD Norm-Col' 'ortho PC 1'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDtreeIII-lev3'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;

    subplot(4,4,1) ;
    vax = axis ;
    text(vax(1) + 0.05 * (vax(2) - vax(1)), ...
         vax(3) + 0.93 * (vax(4) - vax(3)), ...
         'SigClust, P = 7 x 10^{-12}') ;
             %  P-value from Yufeng

    subplot(4,4,6) ;
    vax = axis ;
    text(vax(1) + 0.05 * (vax(2) - vax(1)), ...
         vax(3) + 0.93 * (vax(4) - vax(3)), ...
         'SigClust, P = 4 x 10^{-17}') ;
             %  P-value from Yufeng

    subplot(4,4,11) ;
    vax = axis ;
    text(vax(1) + 0.05 * (vax(2) - vax(1)), ...
         vax(3) + 0.93 * (vax(4) - vax(3)), ...
         'SigClust, P = 3 x 10^{-16}') ;
             %  P-value from Yufeng


    orient landscape ;
    print('-dpsc2',savestr) ;


  elseif  ipart == 51  | ...
          ipart == 52  | ...
          ipart == 53  | ...
          ipart == 54  | ...
          ipart == 55  | ...
          ipart == 56  | ...
          ipart == 57  ;     %  Use "new labelling", where one Normal is classified as Colon


    NewColonFlag = ColonFlag ;
    NewNormalFlag = NormalFlag ;
    for i = 1:length(ClassLabels) ;

      if strcmp(CaseNames{i},'CL2001032256AA.CEL') ;
        NewColonFlag(i) = 1 ;
        NewNormalFlag(i) = 0 ;
        disp(['Changed labels for Case ' 'CL2001032256AA.CEL']) ;
      end ;

    end ;

    NewColonFlag = logical(NewColonFlag) ;
    NewNormalFlag = logical(NewNormalFlag) ;


    if ipart == 51 ;     %  Proj on DWD Carc-NewCol, DWD NewNorm-SC, & ortho PC 

      mdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,NewColonFlag)) ;
      mdir = [mdir DWD1SM(mdata(:,NewNormalFlag),mdata(:,SmallCellFlag))] ;

      proj2md = mdir * inv(mdir' * mdir) * mdir' ;
          %  Projection matrix onto 2-d subspace generated by mdir
      mdataortho = mdata - proj2md * mdata ;
          %  Projection of data onto subspace ortho to mdir
      proj2md = [] ;
          %  to save space
      paramstruct = struct('npc',2,...
                           'viout',[0 1]) ;
      outstruct = pcaSM(mdataortho,paramstruct) ;
      mpcadir = getfield(outstruct,'meigvec') ;
      mdir = [mdir, mpcadir] ;


      titlecellstr = {{'DWD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

      labelcellstr = {{'DWD Carc - NewCol' 'DWD NewNorm-SC' 'ortho PC 1' 'ortho PC 2'}} ;

      savestr = ['NeilLung1ip0' num2str(ipart) 'DWDCarc-NewColDWDNewNorm-SC'] ;

      paramstruct = struct('icolor',icolor, ...
                           'markerstr',markerstr, ...
                           'isubpopkde',1, ...
                           'legendcellstr',legendcellstr, ...
                           'mlegendcolor',mlegendcolor, ...
                           'titlecellstr',titlecellstr, ...
                           'labelcellstr',labelcellstr, ...
                           'savestr',savestr, ...
                           'iscreenwrite',1) ;
 

      scatplotSM(mdata,mdir,paramstruct) ;


    elseif ipart == 52 ;     %  Proj on DWD Carc-NewNorm, DWD NewCol-SC, & ortho PC 

      mdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,NewNormalFlag)) ;
      mdir = [mdir DWD1SM(mdata(:,NewColonFlag),mdata(:,SmallCellFlag))] ;

      proj2md = mdir * inv(mdir' * mdir) * mdir' ;
          %  Projection matrix onto 2-d subspace generated by mdir
      mdataortho = mdata - proj2md * mdata ;
          %  Projection of data onto subspace ortho to mdir
      proj2md = [] ;
          %  to save space
      paramstruct = struct('npc',2,...
                           'viout',[0 1]) ;
      outstruct = pcaSM(mdataortho,paramstruct) ;
      mpcadir = getfield(outstruct,'meigvec') ;
      mdir = [mdir, mpcadir] ;


      titlecellstr = {{'DWD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

      labelcellstr = {{'DWD Carc - NewNorm' 'DWD NewCol-SC' 'ortho PC 1' 'ortho PC 2'}} ;

      savestr = ['NeilLung1ip0' num2str(ipart) 'DWDCarc-NewNormDWDNewCol-SC'] ;

      paramstruct = struct('icolor',icolor, ...
                           'markerstr',markerstr, ...
                           'isubpopkde',1, ...
                           'legendcellstr',legendcellstr, ...
                           'mlegendcolor',mlegendcolor, ...
                           'titlecellstr',titlecellstr, ...
                           'labelcellstr',labelcellstr, ...
                           'savestr',savestr, ...
                           'iscreenwrite',1) ;
 

      scatplotSM(mdata,mdir,paramstruct) ;


    elseif ipart == 53 ;     %  Proj on DWD Carc-SC, DWD NewNorm-NewCol, & ortho PC 

      mdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag)) ;
      mdir = [mdir DWD1SM(mdata(:,NewNormalFlag),mdata(:,NewColonFlag))] ;

      proj2md = mdir * inv(mdir' * mdir) * mdir' ;
          %  Projection matrix onto 2-d subspace generated by mdir
      mdataortho = mdata - proj2md * mdata ;
          %  Projection of data onto subspace ortho to mdir
      proj2md = [] ;
          %  to save space
      paramstruct = struct('npc',2,...
                           'viout',[0 1]) ;
      outstruct = pcaSM(mdataortho,paramstruct) ;
      mpcadir = getfield(outstruct,'meigvec') ;
      mdir = [mdir, mpcadir] ;


      titlecellstr = {{'DWD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

      labelcellstr = {{'DWD Carc - SC' 'DWD NewNorm-NewCol' 'ortho PC 1' 'ortho PC 2'}} ;

      savestr = ['NeilLung1ip0' num2str(ipart) 'DWDCarc-SCDWDNewNorm-NewCol'] ;

      paramstruct = struct('icolor',icolor, ...
                           'markerstr',markerstr, ...
                           'isubpopkde',1, ...
                           'legendcellstr',legendcellstr, ...
                           'mlegendcolor',mlegendcolor, ...
                           'titlecellstr',titlecellstr, ...
                           'labelcellstr',labelcellstr, ...
                           'savestr',savestr, ...
                           'iscreenwrite',1) ;
 

      scatplotSM(mdata,mdir,paramstruct) ;


    elseif ipart == 54 ;     %  Proj on DWD (Carc&SC)-(NewNorm&NewCol), DWD (Carc&NewNorm)-(SC&NewCol), & ortho PC 

      mdir = DWD1SM(mdata(:,CarcinoidFlag|SmallCellFlag),mdata(:,NewNormalFlag|NewColonFlag)) ;
      mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag|NewNormalFlag),mdata(:,SmallCellFlag|NewColonFlag))] ;

      proj2md = mdir * inv(mdir' * mdir) * mdir' ;
          %  Projection matrix onto 2-d subspace generated by mdir
      mdataortho = mdata - proj2md * mdata ;
          %  Projection of data onto subspace ortho to mdir
      proj2md = [] ;
          %  to save space
      paramstruct = struct('npc',2,...
                           'viout',[0 1]) ;
      outstruct = pcaSM(mdataortho,paramstruct) ;
      mpcadir = getfield(outstruct,'meigvec') ;
      mdir = [mdir, mpcadir] ;


      titlecellstr = {{'DWD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

      labelcellstr = {{'DWD CarcSC-NewNormNewCol ' 'DWD CarcNewNorm-SCNewCol' 'ortho PC 1' 'ortho PC 2'}} ;

      savestr = ['NeilLung1ip0' num2str(ipart) 'DWDCarcSC-NewNormNewColDWDCarcNewNorm-SCNewCol'] ;

      paramstruct = struct('icolor',icolor, ...
                           'markerstr',markerstr, ...
                           'isubpopkde',1, ...
                           'legendcellstr',legendcellstr, ...
                           'mlegendcolor',mlegendcolor, ...
                           'titlecellstr',titlecellstr, ...
                           'labelcellstr',labelcellstr, ...
                           'savestr',savestr, ...
                           'iscreenwrite',1) ;

      scatplotSM(mdata,mdir,paramstruct) ;
 

    elseif ipart == 55 ;     %  Proj on DWD (Carc&NewNorm)-(SC&NewCol), DWD (Carc&NewCol)-(NewNorm&SC), & ortho PC 

      mdir = DWD1SM(mdata(:,CarcinoidFlag|NewNormalFlag),mdata(:,SmallCellFlag|NewColonFlag)) ;
      mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag|NewColonFlag),mdata(:,NewNormalFlag|SmallCellFlag))] ;

      proj2md = mdir * inv(mdir' * mdir) * mdir' ;
          %  Projection matrix onto 2-d subspace generated by mdir
      mdataortho = mdata - proj2md * mdata ;
          %  Projection of data onto subspace ortho to mdir
      proj2md = [] ;
          %  to save space
      paramstruct = struct('npc',2,...
                           'viout',[0 1]) ;
      outstruct = pcaSM(mdataortho,paramstruct) ;
      mpcadir = getfield(outstruct,'meigvec') ;
      mdir = [mdir, mpcadir] ;


      titlecellstr = {{'DWD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

      labelcellstr = {{'DWD CarcNewNorm-SCNewCol ' 'DWD CarcNewCol-NewNormSC' 'ortho PC 1' 'ortho PC 2'}} ;

      savestr = ['NeilLung1ip0' num2str(ipart) 'DWDCarcNewNorm-NewNormSCDWDCarcNewNorm-SCNewCol'] ;

      paramstruct = struct('icolor',icolor, ...
                           'markerstr',markerstr, ...
                           'isubpopkde',1, ...
                           'legendcellstr',legendcellstr, ...
                           'mlegendcolor',mlegendcolor, ...
                           'titlecellstr',titlecellstr, ...
                           'labelcellstr',labelcellstr, ...
                           'savestr',savestr, ...
                           'iscreenwrite',1) ;

      scatplotSM(mdata,mdir,paramstruct) ;


    elseif ipart == 56 ;     %  Proj on DWD (Carc&NewCol)-(SC&NewNorm), DWD (Carc&SC)-(NewNorm&NewCol), & ortho PC 

      mdir = DWD1SM(mdata(:,CarcinoidFlag|NewColonFlag),mdata(:,SmallCellFlag|NewNormalFlag)) ;
      mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag|SmallCellFlag),mdata(:,NewNormalFlag|NewColonFlag))] ;

      proj2md = mdir * inv(mdir' * mdir) * mdir' ;
          %  Projection matrix onto 2-d subspace generated by mdir
      mdataortho = mdata - proj2md * mdata ;
          %  Projection of data onto subspace ortho to mdir
      proj2md = [] ;
          %  to save space
      paramstruct = struct('npc',2,...
                           'viout',[0 1]) ;
      outstruct = pcaSM(mdataortho,paramstruct) ;
      mpcadir = getfield(outstruct,'meigvec') ;
      mdir = [mdir, mpcadir] ;


      titlecellstr = {{'DWD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

      labelcellstr = {{'DWD CarcNewCol-SCNewNorm ' 'DWD CarcSC-NewNormNewCol' 'ortho PC 1' 'ortho PC 2'}} ;

      savestr = ['NeilLung1ip0' num2str(ipart) 'DWDCarcNewCol-SCNewNormDWDCarcSC-NewNormNewCol'] ;

      paramstruct = struct('icolor',icolor, ...
                           'markerstr',markerstr, ...
                           'isubpopkde',1, ...
                           'legendcellstr',legendcellstr, ...
                           'mlegendcolor',mlegendcolor, ...
                           'titlecellstr',titlecellstr, ...
                           'labelcellstr',labelcellstr, ...
                           'savestr',savestr, ...
                           'iscreenwrite',1) ;

      scatplotSM(mdata,mdir,paramstruct) ;


    elseif ipart == 57 ;     %  Proj on DWD Indiv vs. all 

      mdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag|NewNormalFlag|NewColonFlag)) ;
      mdir = [mdir DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NewNormalFlag|NewColonFlag))] ;
      mdir = [mdir DWD1SM(mdata(:,NewNormalFlag),mdata(:,CarcinoidFlag|SmallCellFlag|NewColonFlag))] ;
      mdir = [mdir DWD1SM(mdata(:,NewColonFlag),mdata(:,NewNormalFlag|CarcinoidFlag|SmallCellFlag))] ;


      titlecellstr = {{'DWD View' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

      labelcellstr = {{'DWD Carc-all ' 'DWD SC-all' 'DWD NewNorm-all' 'DWD NewCol-all'}} ;

      savestr = ['NeilLung1ip0' num2str(ipart) 'DWDIndiv-all'] ;

      paramstruct = struct('icolor',icolor, ...
                           'markerstr',markerstr, ...
                           'isubpopkde',1, ...
                           'legendcellstr',legendcellstr, ...
                           'mlegendcolor',mlegendcolor, ...
                           'titlecellstr',titlecellstr, ...
                           'labelcellstr',labelcellstr, ...
                           'savestr',savestr, ...
                           'iscreenwrite',1) ;

      scatplotSM(mdata,mdir,paramstruct) ;


    end ;    %  of inner ipart if-block



  elseif  ipart == 61  | ...
          ipart == 62  | ...
          ipart == 63  | ...
          ipart == 64  | ...
          ipart == 65  | ...
          ipart == 66  ;     %  Tree 3, based on Random Relabelling
                             %  Proj on DWD SC-(Carc&Norm&Col), DWD Carc-(Col&Norm), DWD Norm-Col, & ortho PC 

    if ipart == 61 | ...
       ipart == 62 ;
      seed = 109874534 ;
      seedstr = '1' ;
    elseif ipart == 63 | ...
           ipart == 64 ;
      seed = 394750959 ;
      seedstr = '2' ;
    elseif ipart == 65 | ...
           ipart == 66 ;
      seed = 750298742 ;
      seedstr = '3' ;
    end ;

    rand('state',seed) ;
    vrand = randperm(n) ;
        %  random permutation of 1,...,n

    RRCarcinoidFlag = CarcinoidFlag(vrand) ;
    RRColonFlag = ColonFlag(vrand) ;
    RRNormalFlag = NormalFlag(vrand) ;
    RRSmallCellFlag = SmallCellFlag(vrand) ;
        %  Randomly Relabelled Flag Vectors (all with same permutation)


    if ipart == 61 | ...
       ipart == 63 | ...
       ipart == 65 ;    %  then use relabeled colors
      RRicolor = icolor(vrand,:) ;
      titcolstr = 'Relabelled Colors' ;
      savstrcol = 'RelCol' ;
    else ;    %  use original colors
      RRicolor = icolor ;
      titcolstr = 'Original Colors' ;
      savstrcol = 'OriCol' ;
    end ;


    mdir = DWD1SM(mdata(:,RRSmallCellFlag),mdata(:,RRCarcinoidFlag|RRNormalFlag|RRColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,RRCarcinoidFlag),mdata(:,RRColonFlag|RRNormalFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,RRColonFlag),mdata(:,RRNormalFlag))] ;

    proj3md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 3-d subspace generated by mdir
    mdataortho = mdata - proj3md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj3md = [] ;
        %  to save space
    paramstruct = struct('npc',1,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;

    titlecellstr = {{'DWD View, Tree III' ['Random Relabel Dir''ns ' seedstr] titcolstr ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcNormCol-SC' 'DWD ColNorm-Carc' 'DWD Norm-Col' 'ortho PC 1'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDtreeIIIRandRel' seedstr savstrcol] ;

    paramstruct = struct('icolor',RRicolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif  ipart == 71  | ...
          ipart == 72  | ...
          ipart == 73  | ...
          ipart == 74  | ...
          ipart == 75  | ...
          ipart == 76  | ...
          ipart == 77  | ...
          ipart == 78  ;     %  Tree 3, based on subsamples

    rand('state',349857797) ;
    vrand6 = randperm(6) ;
        %  random permutation of 1,...,6
    vrandrest = randperm(n-6) ;
        %  random permutation of 1,...,n-6

    SSSmallCellNum = zeros(1,n) ;
    SSSmallCellNum(SmallCellFlag) = vrand6 ;
    SSrestNum = zeros(1,n) ;
    SSrestNum(~SmallCellFlag) = vrandrest ;

    if ipart == 71 | ...
       ipart == 72 ;
      SSSmallCellFlag = (SSSmallCellNum > 1.5) ;
      SSprop = 0.8 ;
      SSrestFlag = (SSrestNum >= (1 - SSprop) * (n - 6)) ;
    elseif ipart == 73 | ...
           ipart == 74 ;
      SSSmallCellFlag = (SSSmallCellNum > 2.5) ;
      SSprop = 0.6 ;
      SSrestFlag = (SSrestNum >= (1 - SSprop) * (n - 6)) ;
    elseif ipart == 75 | ...
           ipart == 76 ;
      SSSmallCellFlag = (SSSmallCellNum > 3.5) ;
      SSprop = 0.4 ;
      SSrestFlag = (SSrestNum >= (1 - SSprop) * (n - 6)) ;
    elseif ipart == 77 | ...
           ipart == 78 ;
      SSSmallCellFlag = (SSSmallCellNum > 4.5) ;
      SSprop = 0.2 ;
      SSrestFlag = (SSrestNum >= (1 - SSprop) * (n - 6)) ;
    end ;

    SSCarcinoidFlag = CarcinoidFlag & SSrestFlag ;
    SSColonFlag = ColonFlag & SSrestFlag ;
    SSNormalFlag = NormalFlag & SSrestFlag ;
    SSoutFlag = ~(SSSmallCellFlag|SSCarcinoidFlag|SSColonFlag|SSNormalFlag) ;
        %  Randomly SubSampled Flag Vectors (all with same permutation)


    if ipart == 71 | ...
       ipart == 73 | ...
       ipart == 75 | ...
       ipart == 77 ;    %  then use relabeled colors
      SSicolor = icolor ;
      SSicolor(SSoutFlag,:) = ones(sum(SSoutFlag),1) * [0.5 0.5 0.5] ;
      titcolstr = 'Unused Points in Gray' ;
      savstrcol = 'MissGray' ;
    else ;    %  use full colors
      SSicolor = icolor ;
      titcolstr = 'Full Colors' ;
      savstrcol = 'FullCol' ;
    end ;


    mdir = DWD1SM(mdata(:,SSSmallCellFlag),mdata(:,SSCarcinoidFlag|SSNormalFlag|SSColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,SSCarcinoidFlag),mdata(:,SSColonFlag|SSNormalFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,SSColonFlag),mdata(:,SSNormalFlag))] ;

    proj3md = mdir * inv(mdir' * mdir) * mdir' ;
        %  Projection matrix onto 3-d subspace generated by mdir
    mdataortho = mdata - proj3md * mdata ;
        %  Projection of data onto subspace ortho to mdir
    proj3md = [] ;
        %  to save space
    paramstruct = struct('npc',1,...
                         'viout',[0 1]) ;
    outstruct = pcaSM(mdataortho,paramstruct) ;
    mpcadir = getfield(outstruct,'meigvec') ;
    mdir = [mdir, mpcadir] ;

    titlecellstr = {{'DWD View, Tree III' [num2str(100 * SSprop) '% Subampled Dir''ns '] titcolstr ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcNormCol-SC' 'DWD ColNorm-Carc' 'DWD Norm-Col' 'ortho PC 1'}} ;

    savestr = ['NeilLung1ip0' num2str(ipart) 'DWDtreeIIISubSamp' num2str(100 * SSprop) savstrcol] ;

    paramstruct = struct('icolor',SSicolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdir,paramstruct) ;


  elseif  ipart == 101  | ...
          ipart == 102  | ...
          ipart == 103  | ...
          ipart == 104  | ...
          ipart == 105  | ...
          ipart == 106  | ...
          ipart == 107  ;     %  Tree 3, gene masking, DWD SC vs. Rest


    mdir = DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
    vdir = mdir ;
    mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag),mdata(:,ColonFlag|NormalFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,NormalFlag))] ;

    maxlim = axisSM(mdir(:,1)' * mdata) ;
    maxlim = [maxlim; axisSM(mdir(:,2)' * mdata)] ;
    maxlim = [maxlim; axisSM(mdir(:,3)' * mdata)] ;
        %  Set this so have common axes for all masks

    if ipart == 101 ;
      ngunmask = d ;
    elseif ipart == 102 ;
      ngunmask = 1200 ;
    elseif ipart == 103 ;
      ngunmask = 600 ;
    elseif ipart == 104 ;
      ngunmask = 300 ;
    elseif ipart == 105 ;
      ngunmask = 150 ;
    elseif ipart == 106 ;
      ngunmask = 90 ;
    elseif ipart == 107 ;
      ngunmask = 60 ;
    end ;
    ngustr = num2str(ngunmask) ;

    [savdir,vi] = sort(abs(vdir)) ;
    viUM = vi((d - ngunmask + 1):d) ;
    flagUM = zeros(d,1) ;
    flagUM(viUM) = ones(ngunmask,1) ; 

    mdirUM = mdir .* (flagUM * ones(1,3)) ;
    mdirUM(:,1) = mdirUM(:,1) / norm(mdirUM(:,1)) ;
    mdirUM(:,2) = mdirUM(:,2) / norm(mdirUM(:,2)) ;
    mdirUM(:,3) = mdirUM(:,3) / norm(mdirUM(:,3)) ;

    titlecellstr = {{'4 Class Lung Data, Tree III' ['DWD SC-rest, ' ngustr ' Genes Unmasked'] ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcNormCol-SC' 'DWD Carc-ColNorm' 'DWD Norm-Col'}} ;

    savestr = ['NeilLung1ip' num2str(ipart) 'DWDSC-restTreeIII-GUmask' ngustr] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'maxlim',maxlim, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdirUM,paramstruct) ;


  elseif  ipart == 111  | ...
          ipart == 112  | ...
          ipart == 113  | ...
          ipart == 114  | ...
          ipart == 115  | ...
          ipart == 116  | ...
          ipart == 117  ;     %  Tree 3, gene masking, DWD Carc vs. NormCol


    mdir = DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
    vdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,ColonFlag|NormalFlag)) ;
    mdir = [mdir vdir] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,NormalFlag))] ;

    maxlim = axisSM(mdir(:,1)' * mdata) ;
    maxlim = [maxlim; axisSM(mdir(:,2)' * mdata)] ;
    maxlim = [maxlim; axisSM(mdir(:,3)' * mdata)] ;
        %  Set this so have common axes for all masks

    if ipart == 111 ;
      ngunmask = d ;
    elseif ipart == 112 ;
      ngunmask = 1200 ;
    elseif ipart == 113 ;
      ngunmask = 600 ;
    elseif ipart == 114 ;
      ngunmask = 300 ;
    elseif ipart == 115 ;
      ngunmask = 150 ;
    elseif ipart == 116 ;
      ngunmask = 90 ;
    elseif ipart == 117 ;
      ngunmask = 60 ;
    end ;
    ngustr = num2str(ngunmask) ;

    [savdir,vi] = sort(abs(vdir)) ;
    viUM = vi((d - ngunmask + 1):d) ;
    flagUM = zeros(d,1) ;
    flagUM(viUM) = ones(ngunmask,1) ; 

    mdirUM = mdir .* (flagUM * ones(1,3)) ;
    mdirUM(:,1) = mdirUM(:,1) / norm(mdirUM(:,1)) ;
    mdirUM(:,2) = mdirUM(:,2) / norm(mdirUM(:,2)) ;
    mdirUM(:,3) = mdirUM(:,3) / norm(mdirUM(:,3)) ;


    titlecellstr = {{'4 Class Lung Data, Tree III' ['DWD Carc vs. NormCol, ' ngustr ' Genes Unmasked'] ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcNormCol-SC' 'DWD Carc-ColNorm' 'DWD Norm-Col'}} ;

    savestr = ['NeilLung1ip' num2str(ipart) 'DWDCarc-NormColTreeIII-GUmask' ngustr] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'maxlim',maxlim, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdirUM,paramstruct) ;


  elseif  ipart == 121  | ...
          ipart == 122  | ...
          ipart == 123  | ...
          ipart == 124  | ...
          ipart == 125  | ...
          ipart == 126  | ...
          ipart == 127  ;     %  Tree 3, gene masking, DWD Norm vs. Col


    mdir = DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag),mdata(:,ColonFlag|NormalFlag))] ;
    vdir = DWD1SM(mdata(:,ColonFlag),mdata(:,NormalFlag)) ;
    mdir = [mdir vdir] ;

    maxlim = axisSM(mdir(:,1)' * mdata) ;
    maxlim = [maxlim; axisSM(mdir(:,2)' * mdata)] ;
    maxlim = [maxlim; axisSM(mdir(:,3)' * mdata)] ;
        %  Set this so have common axes for all masks

    if ipart == 121 ;
      ngunmask = d ;
    elseif ipart == 122 ;
      ngunmask = 1200 ;
    elseif ipart == 123 ;
      ngunmask = 600 ;
    elseif ipart == 124 ;
      ngunmask = 300 ;
    elseif ipart == 125 ;
      ngunmask = 150 ;
    elseif ipart == 126 ;
      ngunmask = 90 ;
    elseif ipart == 127 ;
      ngunmask = 60 ;
    end ;
    ngustr = num2str(ngunmask) ;

    [savdir,vi] = sort(abs(vdir)) ;
    viUM = vi((d - ngunmask + 1):d) ;
    flagUM = zeros(d,1) ;
    flagUM(viUM) = ones(ngunmask,1) ; 

    mdirUM = mdir .* (flagUM * ones(1,3)) ;
    mdirUM(:,1) = mdirUM(:,1) / norm(mdirUM(:,1)) ;
    mdirUM(:,2) = mdirUM(:,2) / norm(mdirUM(:,2)) ;
    mdirUM(:,3) = mdirUM(:,3) / norm(mdirUM(:,3)) ;


    titlecellstr = {{'4 Class Lung Data, Tree III' ['DWD Norm vs. Col, ' ngustr ' Genes Unmasked'] ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcNormCol-SC' 'DWD Carc-ColNorm' 'DWD Norm-Col'}} ;

    savestr = ['NeilLung1ip' num2str(ipart) 'DWDNorm-ColTreeIII-GUmask' ngustr] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'maxlim',maxlim, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdirUM,paramstruct) ;


  elseif  ipart == 131  | ...
          ipart == 132  | ...
          ipart == 133  | ...
          ipart == 134  | ...
          ipart == 135  | ...
          ipart == 136  | ...
          ipart == 137  | ...
          ipart == 138  | ...
          ipart == 139  ;     %  Tree 3, gene masking, DWD Carc vs. NormCol


    mdir = DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
    vdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,ColonFlag|NormalFlag)) ;
    mdir = [mdir vdir] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,NormalFlag))] ;

    maxlim = axisSM(mdir(:,1)' * mdata) ;
    maxlim = [maxlim; axisSM(mdir(:,2)' * mdata)] ;
    maxlim = [maxlim; axisSM(mdir(:,3)' * mdata)] ;
        %  Set this so have common axes for all masks

    if ipart == 131 ;
      ngunmask = 1024 ;
    elseif ipart == 132 ;
      ngunmask = 512 ;
    elseif ipart == 133 ;
      ngunmask = 256 ;
    elseif ipart == 134 ;
      ngunmask = 128 ;
    elseif ipart == 135 ;
      ngunmask = 64 ;
    elseif ipart == 136 ;
      ngunmask = 32 ;
    elseif ipart == 137 ;
      ngunmask = 16 ;
    elseif ipart == 138 ;
      ngunmask = 8 ;
    elseif ipart == 139 ;
      ngunmask = 4 ;
    end ;
    ngustr = num2str(ngunmask) ;

    [savdir,vi] = sort(abs(vdir)) ;
    viUM = vi((d - ngunmask + 1):d) ;
    flagUM = zeros(d,1) ;
    flagUM(viUM) = ones(ngunmask,1) ; 

    mdirUM = mdir .* (flagUM * ones(1,3)) ;
    mdirUM(:,1) = mdirUM(:,1) / norm(mdirUM(:,1)) ;
    mdirUM(:,2) = mdirUM(:,2) / norm(mdirUM(:,2)) ;
    mdirUM(:,3) = mdirUM(:,3) / norm(mdirUM(:,3)) ;


    titlecellstr = {{'4 Class Lung Data, Tree III' ['DWD Carc vs. NormCol, ' ngustr ' Genes Unmasked'] ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcNormCol-SC' 'DWD Carc-ColNorm' 'DWD Norm-Col'}} ;

    savestr = ['NeilLung1ip' num2str(ipart) 'DWDCarc-NormColTreeIII-GUmask' ngustr] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'maxlim',maxlim, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdirUM,paramstruct) ;


  elseif  ipart == 140  ;    %  Tree 3, Gene Focussing, DWD Carc vs. NormCol, Save indices

    mdir = DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
    vdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,ColonFlag|NormalFlag)) ;
    mdir = [mdir vdir] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,NormalFlag))] ;

    vngunm = [(1:10) (12:2:20) (25:5:50) (50:10:200) (220:20:400) (450:50:1000) (1100:100:2300)]' ;
        %  vector of numbers of genes to unmask 

    vCIf1 = [] ;
    vCIf2 = [] ;
    vCIf3 = [] ;
    vMD1 = [] ;
    vMD2 = [] ;
    vMD3 = [] ;
    vMDP1 = [] ;
    vMDP2 = [] ;
    vMDP3 = [] ;
    vCIp1 = [] ;
    vCIp2 = [] ;
    vCIp3 = [] ;
    for ngunm = vngunm' ;

      disp(['    Working on ngunm = ' num2str(ngunm) ' of ' num2str(vngunm(end))]) ;

      [savdir,vi] = sort(abs(vdir)) ;
      viUM = vi((d - ngunm + 1):d) ;
      flagUM = zeros(d,1) ;
      flagUM(viUM) = ones(ngunm,1) ; 

      mdataUM = mdata(viUM,:) ;

      %  Compute Full Space Cluster Indices
      %
      vCIf1 = [vCIf1; ClustIndSM(mdataUM,SmallCellFlag,CarcinoidFlag|NormalFlag|ColonFlag)] ;
      vCIf2 = [vCIf2; ClustIndSM(mdataUM,CarcinoidFlag,NormalFlag|ColonFlag)] ;
      vCIf3 = [vCIf3; ClustIndSM(mdataUM,NormalFlag,ColonFlag)] ;

      %  Compute Mean Differences
      %
      vMD1 = [vMD1; norm(mean(mdataUM(:,SmallCellFlag),2) - ...
                         mean(mdataUM(:,CarcinoidFlag|NormalFlag|ColonFlag),2))] ;
      vMD2 = [vMD2; norm(mean(mdataUM(:,CarcinoidFlag),2) - ...
                         mean(mdataUM(:,NormalFlag|ColonFlag),2))] ;
      vMD3 = [vMD3; norm(mean(mdataUM(:,NormalFlag),2) - ...
                         mean(mdataUM(:,ColonFlag),2))] ;

      %  Compute Maximal Data Piling Gap
      %
      if ngunm > n ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,SmallCellFlag),mdataUM(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
        vMDP1 = [vMDP1; MDPgap] ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,CarcinoidFlag),mdataUM(:,NormalFlag|ColonFlag)) ;
        vMDP2 = [vMDP2; MDPgap] ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,NormalFlag),mdataUM(:,ColonFlag)) ;
        vMDP3 = [vMDP3; MDPgap] ;
      else ;
        vMDP1 = [vMDP1; NaN] ;
        vMDP2 = [vMDP2; NaN] ;
        vMDP3 = [vMDP3; NaN] ;
      end ;

      %  Compute Projected Data
      %
      mdirUM = mdir .* (flagUM * ones(1,3)) ;
      mdirUM(:,1) = mdirUM(:,1) / norm(mdirUM(:,1)) ;
      mdirUM(:,2) = mdirUM(:,2) / norm(mdirUM(:,2)) ;
      mdirUM(:,3) = mdirUM(:,3) / norm(mdirUM(:,3)) ;
      mdataproj = mdirUM' * mdata ;
          %  projections of data, in 3 tree directions

      %  Compute Projected Cluster Indices
      %
      vCIp1 = [vCIp1; ClustIndSM(mdataproj,SmallCellFlag,CarcinoidFlag|NormalFlag|ColonFlag)] ;
      vCIp2 = [vCIp2; ClustIndSM(mdataproj,CarcinoidFlag,NormalFlag|ColonFlag)] ;
      vCIp3 = [vCIp3; ClustIndSM(mdataproj,NormalFlag,ColonFlag)] ;

    end ;

    %  Save results in .mat file
    %
    save NeilLung1ip140.mat vngunm vCIf1 vCIf2 vCIf3 vMD1 vMD2 vMD3 vMDP1 vMDP2 vMDP3 vCIp1 vCIp2 vCIp3 ;
        %  save these variables

    disp('Finished data save') ;


  elseif  ipart == 141  | ...
          ipart == 142  | ...
          ipart == 143  ;    %  Tree 3, Gene Focussing, DWD Carc vs. NormCol, Plots of Cluster Indices

    load NeilLung1ip140 ;
        %  loads: 
        %      vngunm
        %      vCIf1 vCIf2 vCIf3
        %      vMD1 vMD2 vMD3
        %      vMDP1 vMDP2 vMDP3
        %      vCIp1 vCIp2 vCIp3
    disp('  Finished Data Load') ;


    if ipart == 141 ;   % Tree 3, Gene Focussing, DWD Carc vs. NormCol, Plots of Cluster Indices

      bottom = 0 ;
      top = 1 ;
      left = 0 ; 
      right = 2300 ;

      subplot(2,1,1) ;
      plot(vngunm,vCIf1,'b-') ;
      axis([left right bottom top]) ;
      hold on ;
        plot(vngunm,vCIf2,'m-') ;
        plot(vngunm,vCIf3,'g-') ;
        plot(vngunm,vCIp1,'b--') ;
        plot(vngunm,vCIp2,'m--') ;
        plot(vngunm,vCIp3,'g--') ;
        plot([1024; 1024],[bottom top],'k-') ;
        plot([512; 512],[bottom top],'k-') ;
        plot([256; 256],[bottom top],'k-') ;
        plot([128; 128],[bottom top],'k-') ;
        plot([64; 64],[bottom top],'k-') ;
        plot([32; 32],[bottom top],'k-') ;
        plot([16; 16],[bottom top],'k-') ;
        plot([8; 8],[bottom top],'k-') ;
        plot([4; 4],[bottom top],'k-') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.8 * (top - bottom), ...
             'SmallCell vs. Carcinoid&Normal&Colon', ...
             'Color','b') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.6 * (top - bottom), ...
             'Carcinoid vs. Normal&Colon', ...
             'Color','m') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.4 * (top - bottom), ...
             'Normal vs. Colon', ...
             'Color','g') ;
      hold off ;
      title('Gene Focussing, based on DWD Carc vs. NormCol Direction') ;
      xlabel('Number of Genes Included') ;
      ylabel('Cluster Index (-- for Projected)') ;

      subplot(2,1,2) ;
      plot(log10(vngunm),vCIf1,'b-') ;
      right = log10(2300) ;
      axis([left right bottom top]) ;
      hold on ;
        plot(log10(vngunm),vCIf2,'m-') ;
        plot(log10(vngunm),vCIf3,'g-') ;
        plot(log10(vngunm),vCIp1,'b--') ;
        plot(log10(vngunm),vCIp2,'m--') ;
        plot(log10(vngunm),vCIp3,'g--') ;
        plot(log10([1024; 1024]),[bottom top],'k-') ;
        plot(log10([512; 512]),[bottom top],'k-') ;
        plot(log10([256; 256]),[bottom top],'k-') ;
        plot(log10([128; 128]),[bottom top],'k-') ;
        plot(log10([64; 64]),[bottom top],'k-') ;
        plot(log10([32; 32]),[bottom top],'k-') ;
        plot(log10([16; 16]),[bottom top],'k-') ;
        plot(log10([8; 8]),[bottom top],'k-') ;
        plot(log10([4; 4]),[bottom top],'k-') ;
      hold off ;
      title('Gene Focussing, based on DWD Carc vs. NormCol Direction') ;
      xlabel('log10(Number of Genes Included)') ;
      ylabel('Cluster Index (-- for Projected)') ;

      savestr = ['NeilLung1ip' num2str(ipart) 'GeneFocusDWDCarc-NormColClustInd'] ;
      orient landscape ;
      print('-dpsc2',savestr) ;


    elseif ipart == 142 ;   % Tree 3, Gene Focussing, DWD Carc vs. NormCol, Plots of Mean Diffs

      vax = axisSM([vMD1 vMD2 vMD3]) ;
      bottom = vax(1) ;
      top = vax(2) ;
      left = 0 ;
      right = 2300 ;

      subplot(2,1,1) ;
      plot(vngunm,vMD1,'b-') ;
      axis([left right bottom top]) ;
      hold on ;
        plot(vngunm,vMD2,'m-') ;
        plot(vngunm,vMD3,'g-') ;
        plot([1024; 1024],[bottom top],'k-') ;
        plot([512; 512],[bottom top],'k-') ;
        plot([256; 256],[bottom top],'k-') ;
        plot([128; 128],[bottom top],'k-') ;
        plot([64; 64],[bottom top],'k-') ;
        plot([32; 32],[bottom top],'k-') ;
        plot([16; 16],[bottom top],'k-') ;
        plot([8; 8],[bottom top],'k-') ;
        plot([4; 4],[bottom top],'k-') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.8 * (top - bottom), ...
             'SmallCell vs. Carcinoid&Normal&Colon', ...
             'Color','b') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.6 * (top - bottom), ...
             'Carcinoid vs. Normal&Colon', ...
             'Color','m') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.4 * (top - bottom), ...
             'Normal vs. Colon', ...
             'Color','g') ;
      hold off ;
      title('Gene Focussing, based on DWD Carc vs. NormCol Direction') ;
      xlabel('Number of Genes Included') ;
      ylabel('Mean Difference') ;

      subplot(2,1,2) ;
      plot(log10(vngunm),vMD1,'b-') ;
      right = log10(2300) ;
      axis([left right bottom top]) ;
      hold on ;
        plot(log10(vngunm),vMD2,'m-') ;
        plot(log10(vngunm),vMD3,'g-') ;
        plot(log10([1024; 1024]),[bottom top],'k-') ;
        plot(log10([512; 512]),[bottom top],'k-') ;
        plot(log10([256; 256]),[bottom top],'k-') ;
        plot(log10([128; 128]),[bottom top],'k-') ;
        plot(log10([64; 64]),[bottom top],'k-') ;
        plot(log10([32; 32]),[bottom top],'k-') ;
        plot(log10([16; 16]),[bottom top],'k-') ;
        plot(log10([8; 8]),[bottom top],'k-') ;
        plot(log10([4; 4]),[bottom top],'k-') ;
      hold off ;
      title('Gene Focussing, based on DWD Carc vs. NormCol Direction') ;
      xlabel('log10(Number of Genes Included)') ;
      ylabel('Mean Difference') ;

      savestr = ['NeilLung1ip' num2str(ipart) 'GeneFocusDWDCarc-NormColMeanDiff'] ;
      orient landscape ;
      print('-dpsc2',savestr) ;


    elseif ipart == 143 ;   % Tree 3, Gene Focussing, DWD Carc vs. NormCol, Plots of MDPs

      vax = axisSM([vMDP1 vMDP2 vMDP3]) ;
      bottom = vax(1) ;
      top = vax(2) ;
      left = 0 ;
      right = 2300 ;

      subplot(2,1,1) ;
      plot(vngunm,vMDP1,'b-') ;
      axis([left right bottom top]) ;
      hold on ;
        plot(vngunm,vMDP2,'m-') ;
        plot(vngunm,vMDP3,'g-') ;
        plot([1024; 1024],[bottom top],'k-') ;
        plot([512; 512],[bottom top],'k-') ;
        plot([256; 256],[bottom top],'k-') ;
        plot([128; 128],[bottom top],'k-') ;
        plot([64; 64],[bottom top],'k-') ;
        plot([32; 32],[bottom top],'k-') ;
        plot([16; 16],[bottom top],'k-') ;
        plot([8; 8],[bottom top],'k-') ;
        plot([4; 4],[bottom top],'k-') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.8 * (top - bottom), ...
             'SmallCell vs. Carcinoid&Normal&Colon', ...
             'Color','b') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.6 * (top - bottom), ...
             'Carcinoid vs. Normal&Colon', ...
             'Color','m') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.4 * (top - bottom), ...
             'Normal vs. Colon', ...
             'Color','g') ;
      hold off ;
      title('Gene Focussing, based on DWD Carc vs. NormCol Direction') ;
      xlabel('Number of Genes Included') ;
      ylabel('Maximal Data Piling Gap') ;

      subplot(2,1,2) ;
      plot(log10(vngunm),vMDP1,'b-') ;
      right = log10(2300) ;
      axis([left right bottom top]) ;
      hold on ;
        plot(log10(vngunm),vMDP2,'m-') ;
        plot(log10(vngunm),vMDP3,'g-') ;
        plot(log10([1024; 1024]),[bottom top],'k-') ;
        plot(log10([512; 512]),[bottom top],'k-') ;
        plot(log10([256; 256]),[bottom top],'k-') ;
        plot(log10([128; 128]),[bottom top],'k-') ;
        plot(log10([64; 64]),[bottom top],'k-') ;
        plot(log10([32; 32]),[bottom top],'k-') ;
        plot(log10([16; 16]),[bottom top],'k-') ;
        plot(log10([8; 8]),[bottom top],'k-') ;
        plot(log10([4; 4]),[bottom top],'k-') ;
      hold off ;
      title('Gene Focussing, based on DWD Carc vs. NormCol Direction') ;
      xlabel('log10(Number of Genes Included)') ;
      ylabel('Maximal Data Piling Gap') ;

      savestr = ['NeilLung1ip' num2str(ipart) 'GeneFocusDWDCarc-NormColMDPgap'] ;
      orient landscape ;
      print('-dpsc2',savestr) ;


    end ;


  elseif  ipart == 151  | ...
          ipart == 152  | ...
          ipart == 153  | ...
          ipart == 154  | ...
          ipart == 155  | ...
          ipart == 156  | ...
          ipart == 157  | ...
          ipart == 158  | ...
          ipart == 159  ;     %  Tree 3, gene masking, DWD Carc vs. NormCol


    mdir = DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
    vdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,ColonFlag|NormalFlag)) ;
    mdir = [mdir vdir] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,NormalFlag))] ;

    if ipart == 151 ;
      ngunmask = 1024 ;
    elseif ipart == 152 ;
      ngunmask = 512 ;
    elseif ipart == 153 ;
      ngunmask = 256 ;
    elseif ipart == 154 ;
      ngunmask = 128 ;
    elseif ipart == 155 ;
      ngunmask = 64 ;
    elseif ipart == 156 ;
      ngunmask = 32 ;
    elseif ipart == 157 ;
      ngunmask = 16 ;
    elseif ipart == 158 ;
      ngunmask = 8 ;
    elseif ipart == 159 ;
      ngunmask = 4 ;
    end ;
    ngustr = num2str(ngunmask) ;

    [savdir,vi] = sort(abs(vdir)) ;
    viUM = vi((d - ngunmask + 1):d) ;
    flagUM = zeros(d,1) ;
    flagUM(viUM) = ones(ngunmask,1) ; 

    mdirUM = mdir .* (flagUM * ones(1,3)) ;
    mdirUM(:,1) = mdirUM(:,1) / norm(mdirUM(:,1)) ;
    mdirUM(:,2) = mdirUM(:,2) / norm(mdirUM(:,2)) ;
    mdirUM(:,3) = mdirUM(:,3) / norm(mdirUM(:,3)) ;


    titlecellstr = {{'4 Class Lung Data, Tree III' ['DWD Carc vs. NormCol, ' ngustr ' Genes Unmasked'] ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD CarcNormCol-SC' 'DWD Carc-ColNorm' 'DWD Norm-Col'}} ;

    savestr = ['NeilLung1ip' num2str(ipart) 'DWDCarc-NormColTreeIII-GUmask' ngustr 'ChAx'] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdirUM,paramstruct) ;


  elseif  ipart == 201  | ...
          ipart == 202  | ...
          ipart == 203  | ...
          ipart == 204  | ...
          ipart == 205  | ...
          ipart == 206  | ...
          ipart == 207  | ...
          ipart == 208  | ...
          ipart == 209  ;     %  Indiv - All, gene masking, DWD Carc vs. Rest


    vdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
    mdir = vdir ;
    mdir = [mdir DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,CarcinoidFlag|SmallCellFlag|ColonFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,CarcinoidFlag|SmallCellFlag|NormalFlag))] ;

    if ipart == 201 ;
      ngunmask = 1024 ;
    elseif ipart == 202 ;
      ngunmask = 512 ;
    elseif ipart == 203 ;
      ngunmask = 256 ;
    elseif ipart == 204 ;
      ngunmask = 128 ;
    elseif ipart == 205 ;
      ngunmask = 64 ;
    elseif ipart == 206 ;
      ngunmask = 32 ;
    elseif ipart == 207 ;
      ngunmask = 16 ;
    elseif ipart == 208 ;
      ngunmask = 8 ;
    elseif ipart == 209 ;
      ngunmask = 4 ;
    end ;
    ngustr = num2str(ngunmask) ;

    [savdir,vi] = sort(abs(vdir)) ;
    viUM = vi((d - ngunmask + 1):d) ;
    flagUM = zeros(d,1) ;
    flagUM(viUM) = ones(ngunmask,1) ; 

    mdirUM = mdir .* (flagUM * ones(1,4)) ;
    mdirUM(:,1) = mdirUM(:,1) / norm(mdirUM(:,1)) ;
    mdirUM(:,2) = mdirUM(:,2) / norm(mdirUM(:,2)) ;
    mdirUM(:,3) = mdirUM(:,3) / norm(mdirUM(:,3)) ;
    mdirUM(:,4) = mdirUM(:,4) / norm(mdirUM(:,4)) ;


    titlecellstr = {{'4 Class Lung Data, Indiv - All' ['DWD Carc vs. Rest, ' ngustr ' Genes Unmasked'] ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD Carc-Rest' 'DWD SC-Rest' 'DWD Norm-Rest' 'DWD Col-Rest'}} ;

    savestr = ['NeilLung1ip' num2str(ipart) 'DWDCarc-RestIndivAll-GUmask' ngustr] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdirUM,paramstruct) ;


  elseif  ipart == 211  | ...
          ipart == 212  | ...
          ipart == 213  | ...
          ipart == 214  | ...
          ipart == 215  | ...
          ipart == 216  | ...
          ipart == 217  | ...
          ipart == 218  | ...
          ipart == 219  ;     %  Indiv - All, gene masking, DWD SC vs. Rest


    mdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
    vdir = DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
    mdir = [mdir vdir] ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,CarcinoidFlag|SmallCellFlag|ColonFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,CarcinoidFlag|SmallCellFlag|NormalFlag))] ;

    if ipart == 211 ;
      ngunmask = 1024 ;
    elseif ipart == 212 ;
      ngunmask = 512 ;
    elseif ipart == 213 ;
      ngunmask = 256 ;
    elseif ipart == 214 ;
      ngunmask = 128 ;
    elseif ipart == 215 ;
      ngunmask = 64 ;
    elseif ipart == 216 ;
      ngunmask = 32 ;
    elseif ipart == 217 ;
      ngunmask = 16 ;
    elseif ipart == 218 ;
      ngunmask = 8 ;
    elseif ipart == 219 ;
      ngunmask = 4 ;
    end ;
    ngustr = num2str(ngunmask) ;

    [savdir,vi] = sort(abs(vdir)) ;
    viUM = vi((d - ngunmask + 1):d) ;
    flagUM = zeros(d,1) ;
    flagUM(viUM) = ones(ngunmask,1) ; 

    mdirUM = mdir .* (flagUM * ones(1,4)) ;
    mdirUM(:,1) = mdirUM(:,1) / norm(mdirUM(:,1)) ;
    mdirUM(:,2) = mdirUM(:,2) / norm(mdirUM(:,2)) ;
    mdirUM(:,3) = mdirUM(:,3) / norm(mdirUM(:,3)) ;
    mdirUM(:,4) = mdirUM(:,4) / norm(mdirUM(:,4)) ;


    titlecellstr = {{'4 Class Lung Data, Indiv - All' ['DWD SC vs. Rest, ' ngustr ' Genes Unmasked'] ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD Carc-Rest' 'DWD SC-Rest' 'DWD Norm-Rest' 'DWD Col-Rest'}} ;

    savestr = ['NeilLung1ip' num2str(ipart) 'DWDSC-RestIndivAll-GUmask' ngustr] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdirUM,paramstruct) ;


  elseif  ipart == 221  | ...
          ipart == 222  | ...
          ipart == 223  | ...
          ipart == 224  | ...
          ipart == 225  | ...
          ipart == 226  | ...
          ipart == 227  | ...
          ipart == 228  | ...
          ipart == 229  ;     %  Indiv - All, gene masking, DWD Norm vs. Rest

    mdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag))] ;
    vdir = DWD1SM(mdata(:,NormalFlag),mdata(:,CarcinoidFlag|SmallCellFlag|ColonFlag)) ;
    mdir = [mdir vdir] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,CarcinoidFlag|SmallCellFlag|NormalFlag))] ;

    if ipart == 221 ;
      ngunmask = 1024 ;
    elseif ipart == 222 ;
      ngunmask = 512 ;
    elseif ipart == 223 ;
      ngunmask = 256 ;
    elseif ipart == 224 ;
      ngunmask = 128 ;
    elseif ipart == 225 ;
      ngunmask = 64 ;
    elseif ipart == 226 ;
      ngunmask = 32 ;
    elseif ipart == 227 ;
      ngunmask = 16 ;
    elseif ipart == 228 ;
      ngunmask = 8 ;
    elseif ipart == 229 ;
      ngunmask = 4 ;
    end ;
    ngustr = num2str(ngunmask) ;

    [savdir,vi] = sort(abs(vdir)) ;
    viUM = vi((d - ngunmask + 1):d) ;
    flagUM = zeros(d,1) ;
    flagUM(viUM) = ones(ngunmask,1) ; 

    mdirUM = mdir .* (flagUM * ones(1,4)) ;
    mdirUM(:,1) = mdirUM(:,1) / norm(mdirUM(:,1)) ;
    mdirUM(:,2) = mdirUM(:,2) / norm(mdirUM(:,2)) ;
    mdirUM(:,3) = mdirUM(:,3) / norm(mdirUM(:,3)) ;
    mdirUM(:,4) = mdirUM(:,4) / norm(mdirUM(:,4)) ;


    titlecellstr = {{'4 Class Lung Data, Indiv - All' ['DWD Norm vs. Rest, ' ngustr ' Genes Unmasked'] ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD Carc-Rest' 'DWD SC-Rest' 'DWD Norm-Rest' 'DWD Col-Rest'}} ;

    savestr = ['NeilLung1ip' num2str(ipart) 'DWDNorm-RestIndivAll-GUmask' ngustr] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdirUM,paramstruct) ;


  elseif  ipart == 231  | ...
          ipart == 232  | ...
          ipart == 233  | ...
          ipart == 234  | ...
          ipart == 235  | ...
          ipart == 236  | ...
          ipart == 237  | ...
          ipart == 238  | ...
          ipart == 239  ;     %  Indiv - All, gene masking, DWD Col vs. Rest

    mdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,CarcinoidFlag|SmallCellFlag|ColonFlag))] ;
    vdir = DWD1SM(mdata(:,ColonFlag),mdata(:,CarcinoidFlag|SmallCellFlag|NormalFlag)) ;
    mdir = [mdir vdir] ;

    if ipart == 231 ;
      ngunmask = 1024 ;
    elseif ipart == 232 ;
      ngunmask = 512 ;
    elseif ipart == 233 ;
      ngunmask = 256 ;
    elseif ipart == 234 ;
      ngunmask = 128 ;
    elseif ipart == 235 ;
      ngunmask = 64 ;
    elseif ipart == 236 ;
      ngunmask = 32 ;
    elseif ipart == 237 ;
      ngunmask = 16 ;
    elseif ipart == 238 ;
      ngunmask = 8 ;
    elseif ipart == 239 ;
      ngunmask = 4 ;
    end ;
    ngustr = num2str(ngunmask) ;

    [savdir,vi] = sort(abs(vdir)) ;
    viUM = vi((d - ngunmask + 1):d) ;
    flagUM = zeros(d,1) ;
    flagUM(viUM) = ones(ngunmask,1) ; 

    mdirUM = mdir .* (flagUM * ones(1,4)) ;
    mdirUM(:,1) = mdirUM(:,1) / norm(mdirUM(:,1)) ;
    mdirUM(:,2) = mdirUM(:,2) / norm(mdirUM(:,2)) ;
    mdirUM(:,3) = mdirUM(:,3) / norm(mdirUM(:,3)) ;
    mdirUM(:,4) = mdirUM(:,4) / norm(mdirUM(:,4)) ;


    titlecellstr = {{'4 Class Lung Data, Indiv - All' ['DWD Col vs. Rest, ' ngustr ' Genes Unmasked'] ['Genes filtered to ' num2str(d) ' exp''d']}} ;

    labelcellstr = {{'DWD Carc-Rest' 'DWD SC-Rest' 'DWD Norm-Rest' 'DWD Col-Rest'}} ;

    savestr = ['NeilLung1ip' num2str(ipart) 'DWDCol-RestIndivAll-GUmask' ngustr] ;

    paramstruct = struct('icolor',icolor, ...
                         'markerstr',markerstr, ...
                         'isubpopkde',1, ...
                         'legendcellstr',legendcellstr, ...
                         'mlegendcolor',mlegendcolor, ...
                         'titlecellstr',titlecellstr, ...
                         'labelcellstr',labelcellstr, ...
                         'savestr',savestr, ...
                         'iscreenwrite',1) ;

    scatplotSM(mdata,mdirUM,paramstruct) ;


  elseif  ipart == 240  ;    %  Indiv - All, Gene Focussing, DWD Carc vs. Rest, Save indices

    vdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
    mdir = vdir ;
    mdir = [mdir DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,CarcinoidFlag|SmallCellFlag|ColonFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,CarcinoidFlag|SmallCellFlag|NormalFlag))] ;

    vngunm = [(1:10) (12:2:20) (25:5:50) (50:10:200) (220:20:400) (450:50:1000) (1100:100:2300)]' ;
        %  vector of numbers of genes to unmask 

    vCIf1 = [] ;
    vCIf2 = [] ;
    vCIf3 = [] ;
    vCIf4 = [] ;
    vMD1 = [] ;
    vMD2 = [] ;
    vMD3 = [] ;
    vMD4 = [] ;
    vMDP1 = [] ;
    vMDP2 = [] ;
    vMDP3 = [] ;
    vMDP4 = [] ;
    vCIp1 = [] ;
    vCIp2 = [] ;
    vCIp3 = [] ;
    vCIp4 = [] ;
    for ngunm = vngunm' ;

      disp(['    Working on ngunm = ' num2str(ngunm) ' of ' num2str(vngunm(end))]) ;

      [savdir,vi] = sort(abs(vdir)) ;
      viUM = vi((d - ngunm + 1):d) ;
      flagUM = zeros(d,1) ;
      flagUM(viUM) = ones(ngunm,1) ; 

      mdataUM = mdata(viUM,:) ;

      %  Compute Full Space Cluster Indices
      %
      vCIf1 = [vCIf1; ClustIndSM(mdataUM,CarcinoidFlag,SmallCellFlag|NormalFlag|ColonFlag)] ;
      vCIf2 = [vCIf2; ClustIndSM(mdataUM,SmallCellFlag,CarcinoidFlag|NormalFlag|ColonFlag)] ;
      vCIf3 = [vCIf3; ClustIndSM(mdataUM,NormalFlag,CarcinoidFlag|SmallCellFlag|ColonFlag)] ;
      vCIf4 = [vCIf4; ClustIndSM(mdataUM,ColonFlag,CarcinoidFlag|SmallCellFlag|NormalFlag)] ;

      %  Compute Mean Differences
      %
      vMD1 = [vMD1; norm(mean(mdataUM(:,CarcinoidFlag),2) - ...
                         mean(mdataUM(:,SmallCellFlag|NormalFlag|ColonFlag),2))] ;
      vMD2 = [vMD2; norm(mean(mdataUM(:,SmallCellFlag),2) - ...
                         mean(mdataUM(:,CarcinoidFlag|NormalFlag|ColonFlag),2))] ;
      vMD3 = [vMD3; norm(mean(mdataUM(:,NormalFlag),2) - ...
                         mean(mdataUM(:,CarcinoidFlag|SmallCellFlag|ColonFlag),2))] ;
      vMD4 = [vMD4; norm(mean(mdataUM(:,ColonFlag),2) - ...
                         mean(mdataUM(:,CarcinoidFlag|NormalFlag|SmallCellFlag),2))] ;

      %  Compute Maximal Data Piling Gap
      %
      if ngunm > n ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,CarcinoidFlag),mdataUM(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
        vMDP1 = [vMDP1; MDPgap] ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,SmallCellFlag),mdataUM(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
        vMDP2 = [vMDP2; MDPgap] ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,NormalFlag),mdataUM(:,CarcinoidFlag|SmallCellFlag|ColonFlag)) ;
        vMDP3 = [vMDP3; MDPgap] ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,ColonFlag),mdataUM(:,CarcinoidFlag|SmallCellFlag|NormalFlag)) ;
        vMDP4 = [vMDP4; MDPgap] ;
      else ;
        vMDP1 = [vMDP1; NaN] ;
        vMDP2 = [vMDP2; NaN] ;
        vMDP3 = [vMDP3; NaN] ;
        vMDP4 = [vMDP4; NaN] ;
      end ;

      %  Compute Projected Data
      %
      mdirUM = mdir .* (flagUM * ones(1,4)) ;
      mdirUM(:,1) = mdirUM(:,1) / norm(mdirUM(:,1)) ;
      mdirUM(:,2) = mdirUM(:,2) / norm(mdirUM(:,2)) ;
      mdirUM(:,3) = mdirUM(:,3) / norm(mdirUM(:,3)) ;
      mdirUM(:,4) = mdirUM(:,4) / norm(mdirUM(:,4)) ;
      mdataproj = mdirUM' * mdata ;
          %  projections of data, in Indiv vs. Rest directions

      %  Compute Projected Cluster Indices
      %
      vCIp1 = [vCIp1; ClustIndSM(mdataproj,CarcinoidFlag,SmallCellFlag|NormalFlag|ColonFlag)] ;
      vCIp2 = [vCIp2; ClustIndSM(mdataproj,SmallCellFlag,CarcinoidFlag|NormalFlag|ColonFlag)] ;
      vCIp3 = [vCIp3; ClustIndSM(mdataproj,NormalFlag,CarcinoidFlag|SmallCellFlag|ColonFlag)] ;
      vCIp4 = [vCIp4; ClustIndSM(mdataproj,ColonFlag,CarcinoidFlag|SmallCellFlag|NormalFlag)] ;

    end ;

    %  Save results in .mat file
    %
    save NeilLung1ip240.mat vngunm vCIf1 vCIf2 vCIf3 vCIf4 vMD1 vMD2 vMD3 vMD4 vMDP1 vMDP2 vMDP3 vMDP4 vCIp1 vCIp2 vCIp3 vCIp4 ;
        %  save these variables

    disp('Finished data save') ;


  elseif  ipart == 241  | ...
          ipart == 242  | ...
          ipart == 243  ;    %  Indiv - All, Gene Focussing, DWD Carc vs. Rest, Plots of Cluster Indices

    load NeilLung1ip240 ;
        %  loads: 
        %      vngunm
        %      vCIf1 vCIf2 vCIf3 vCIf4
        %      vMD1 vMD2 vMD3 vMD4
        %      vMDP1 vMDP2 vMDP3 vMDP4
        %      vCIp1 vCIp2 vCIp3 vCIp4
    disp('  Finished Data Load') ;


    if ipart == 241 ;   % Indiv - All, Gene Focussing, DWD Carc vs. Rest, Plots of Cluster Indices

      bottom = 0 ;
      top = 1 ;
      left = 0 ; 
      right = 2300 ;

      subplot(2,1,1) ;
      plot(vngunm,vCIf1,'m-') ;
      axis([left right bottom top]) ;
      hold on ;
        plot(vngunm,vCIf2,'b-') ;
        plot(vngunm,vCIf3,'g-') ;
        plot(vngunm,vCIf4,'r-') ;
        plot(vngunm,vCIp1,'m--') ;
        plot(vngunm,vCIp2,'b--') ;
        plot(vngunm,vCIp3,'g--') ;
        plot(vngunm,vCIp4,'r--') ;
        plot([1024; 1024],[bottom top],'k-') ;
        plot([512; 512],[bottom top],'k-') ;
        plot([256; 256],[bottom top],'k-') ;
        plot([128; 128],[bottom top],'k-') ;
        plot([64; 64],[bottom top],'k-') ;
        plot([32; 32],[bottom top],'k-') ;
        plot([16; 16],[bottom top],'k-') ;
        plot([8; 8],[bottom top],'k-') ;
        plot([4; 4],[bottom top],'k-') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.8 * (top - bottom), ...
             'Carcinoid vs. SmallCell&Normal&Colon', ...
             'Color','m') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.6 * (top - bottom), ...
             'SmallCell vs. Carcinoid&Normal&Colon', ...
             'Color','b') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.4 * (top - bottom), ...
             'Normal vs. Carcinoid&SmallCell&Colon', ...
             'Color','g') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.2 * (top - bottom), ...
             'Colon vs. Carcinoid&SmallCell&Normal', ...
             'Color','r') ;
      hold off ;
      title('Gene Focussing, based on DWD Carc vs. Rest Direction') ;
      xlabel('Number of Genes Included') ;
      ylabel('Cluster Index (-- for Projected)') ;

      subplot(2,1,2) ;
      plot(log10(vngunm),vCIf1,'m-') ;
      right = log10(2300) ;
      axis([left right bottom top]) ;
      hold on ;
        plot(log10(vngunm),vCIf2,'b-') ;
        plot(log10(vngunm),vCIf3,'g-') ;
        plot(log10(vngunm),vCIf4,'r-') ;
        plot(log10(vngunm),vCIp1,'m--') ;
        plot(log10(vngunm),vCIp2,'b--') ;
        plot(log10(vngunm),vCIp3,'g--') ;
        plot(log10(vngunm),vCIp4,'r--') ;
        plot(log10([1024; 1024]),[bottom top],'k-') ;
        plot(log10([512; 512]),[bottom top],'k-') ;
        plot(log10([256; 256]),[bottom top],'k-') ;
        plot(log10([128; 128]),[bottom top],'k-') ;
        plot(log10([64; 64]),[bottom top],'k-') ;
        plot(log10([32; 32]),[bottom top],'k-') ;
        plot(log10([16; 16]),[bottom top],'k-') ;
        plot(log10([8; 8]),[bottom top],'k-') ;
        plot(log10([4; 4]),[bottom top],'k-') ;
      hold off ;
      title('Gene Focussing, based on DWD Carc vs. Rest Direction') ;
      xlabel('log10(Number of Genes Included)') ;
      ylabel('Cluster Index (-- for Projected)') ;

      savestr = ['NeilLung1ip' num2str(ipart) 'GeneFocusDWDCarc-RestClustInd'] ;
      orient landscape ;
      print('-dpsc2',savestr) ;


    elseif ipart == 242 ;   % Indiv - All, Gene Focussing, DWD Carc vs. Rest, Plots of Mean Diffs

      vax = axisSM([vMD1 vMD2 vMD3 vMD4]) ;
      bottom = vax(1) ;
      top = vax(2) ;
      left = 0 ;
      right = 2300 ;

      subplot(2,1,1) ;
      plot(vngunm,vMD1,'m-') ;
      axis([left right bottom top]) ;
      hold on ;
        plot(vngunm,vMD2,'b-') ;
        plot(vngunm,vMD3,'g-') ;
        plot(vngunm,vMD4,'r-') ;
        plot([1024; 1024],[bottom top],'k-') ;
        plot([512; 512],[bottom top],'k-') ;
        plot([256; 256],[bottom top],'k-') ;
        plot([128; 128],[bottom top],'k-') ;
        plot([64; 64],[bottom top],'k-') ;
        plot([32; 32],[bottom top],'k-') ;
        plot([16; 16],[bottom top],'k-') ;
        plot([8; 8],[bottom top],'k-') ;
        plot([4; 4],[bottom top],'k-') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.8 * (top - bottom), ...
             'Carcinoid vs. SmallCell&Normal&Colon', ...
             'Color','m') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.6 * (top - bottom), ...
             'SmallCell vs. Carcinoid&Normal&Colon', ...
             'Color','b') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.4 * (top - bottom), ...
             'Normal vs. Carcinoid&SmallCell&Colon', ...
             'Color','g') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.2 * (top - bottom), ...
             'Colon vs. Carcinoid&SmallCell&Normal', ...
             'Color','r') ;
      hold off ;
      title('Gene Focussing, based on DWD Carc vs. Rest Direction') ;
      xlabel('Number of Genes Included') ;
      ylabel('Mean Difference') ;

      subplot(2,1,2) ;
      plot(log10(vngunm),vMD1,'m-') ;
      right = log10(2300) ;
      axis([left right bottom top]) ;
      hold on ;
        plot(log10(vngunm),vMD2,'b-') ;
        plot(log10(vngunm),vMD3,'g-') ;
        plot(log10(vngunm),vMD4,'r-') ;
        plot(log10([1024; 1024]),[bottom top],'k-') ;
        plot(log10([512; 512]),[bottom top],'k-') ;
        plot(log10([256; 256]),[bottom top],'k-') ;
        plot(log10([128; 128]),[bottom top],'k-') ;
        plot(log10([64; 64]),[bottom top],'k-') ;
        plot(log10([32; 32]),[bottom top],'k-') ;
        plot(log10([16; 16]),[bottom top],'k-') ;
        plot(log10([8; 8]),[bottom top],'k-') ;
        plot(log10([4; 4]),[bottom top],'k-') ;
      hold off ;
      title('Gene Focussing, based on DWD Carc vs. Rest Direction') ;
      xlabel('log10(Number of Genes Included)') ;
      ylabel('Mean Difference') ;

      savestr = ['NeilLung1ip' num2str(ipart) 'GeneFocusDWDCarc-RestMeanDiff'] ;
      orient landscape ;
      print('-dpsc2',savestr) ;


    elseif ipart == 243 ;   % Indiv - All, Gene Focussing, DWD Carc vs. Rest, Plots of MDPs

      vax = axisSM([vMDP1 vMDP2 vMDP3 vMDP4]) ;
      bottom = vax(1) ;
      top = vax(2) ;
      left = 0 ;
      right = 2300 ;

      subplot(2,1,1) ;
      plot(vngunm,vMDP1,'m-') ;
      axis([left right bottom top]) ;
      hold on ;
        plot(vngunm,vMDP2,'b-') ;
        plot(vngunm,vMDP3,'g-') ;
        plot(vngunm,vMDP4,'r-') ;
        plot([1024; 1024],[bottom top],'k-') ;
        plot([512; 512],[bottom top],'k-') ;
        plot([256; 256],[bottom top],'k-') ;
        plot([128; 128],[bottom top],'k-') ;
        plot([64; 64],[bottom top],'k-') ;
        plot([32; 32],[bottom top],'k-') ;
        plot([16; 16],[bottom top],'k-') ;
        plot([8; 8],[bottom top],'k-') ;
        plot([4; 4],[bottom top],'k-') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.8 * (top - bottom), ...
             'SmallCell vs. Carcinoid&Normal&Colon', ...
             'Color','b') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.6 * (top - bottom), ...
             'Carcinoid vs. Normal&Colon', ...
             'Color','m') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.4 * (top - bottom), ...
             'Normal vs. Colon', ...
             'Color','g') ;
      hold off ;
      title('Gene Focussing, based on DWD Carc vs. Rest Direction') ;
      xlabel('Number of Genes Included') ;
      ylabel('Maximal Data Piling Gap') ;

      subplot(2,1,2) ;
      plot(log10(vngunm),vMDP1,'m-') ;
      right = log10(2300) ;
      axis([left right bottom top]) ;
      hold on ;
        plot(log10(vngunm),vMDP2,'b-') ;
        plot(log10(vngunm),vMDP3,'g-') ;
        plot(log10(vngunm),vMDP4,'r-') ;
        plot(log10([1024; 1024]),[bottom top],'k-') ;
        plot(log10([512; 512]),[bottom top],'k-') ;
        plot(log10([256; 256]),[bottom top],'k-') ;
        plot(log10([128; 128]),[bottom top],'k-') ;
        plot(log10([64; 64]),[bottom top],'k-') ;
        plot(log10([32; 32]),[bottom top],'k-') ;
        plot(log10([16; 16]),[bottom top],'k-') ;
        plot(log10([8; 8]),[bottom top],'k-') ;
        plot(log10([4; 4]),[bottom top],'k-') ;
      hold off ;
      title('Gene Focussing, based on DWD Carc vs. Rest Direction') ;
      xlabel('log10(Number of Genes Included)') ;
      ylabel('Maximal Data Piling Gap') ;

      savestr = ['NeilLung1ip' num2str(ipart) 'GeneFocusDWDCarc-RestMDPgap'] ;
      orient landscape ;
      print('-dpsc2',savestr) ;


    end ;


  elseif  ipart == 250  ;    %  Indiv - All, Gene Focussing, DWD SC vs. Rest, Save indices

    mdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
    vdir = DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
    mdir = [mdir vdir] ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,CarcinoidFlag|SmallCellFlag|ColonFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,CarcinoidFlag|SmallCellFlag|NormalFlag))] ;

    vngunm = [(1:10) (12:2:20) (25:5:50) (50:10:200) (220:20:400) (450:50:1000) (1100:100:2300)]' ;
        %  vector of numbers of genes to unmask 

    vCIf1 = [] ;
    vCIf2 = [] ;
    vCIf3 = [] ;
    vCIf4 = [] ;
    vMD1 = [] ;
    vMD2 = [] ;
    vMD3 = [] ;
    vMD4 = [] ;
    vMDP1 = [] ;
    vMDP2 = [] ;
    vMDP3 = [] ;
    vMDP4 = [] ;
    vCIp1 = [] ;
    vCIp2 = [] ;
    vCIp3 = [] ;
    vCIp4 = [] ;
    for ngunm = vngunm' ;

      disp(['    Working on ngunm = ' num2str(ngunm) ' of ' num2str(vngunm(end))]) ;

      [savdir,vi] = sort(abs(vdir)) ;
      viUM = vi((d - ngunm + 1):d) ;
      flagUM = zeros(d,1) ;
      flagUM(viUM) = ones(ngunm,1) ; 

      mdataUM = mdata(viUM,:) ;

      %  Compute Full Space Cluster Indices
      %
      vCIf1 = [vCIf1; ClustIndSM(mdataUM,CarcinoidFlag,SmallCellFlag|NormalFlag|ColonFlag)] ;
      vCIf2 = [vCIf2; ClustIndSM(mdataUM,SmallCellFlag,CarcinoidFlag|NormalFlag|ColonFlag)] ;
      vCIf3 = [vCIf3; ClustIndSM(mdataUM,NormalFlag,CarcinoidFlag|SmallCellFlag|ColonFlag)] ;
      vCIf4 = [vCIf4; ClustIndSM(mdataUM,ColonFlag,CarcinoidFlag|SmallCellFlag|NormalFlag)] ;

      %  Compute Mean Differences
      %
      vMD1 = [vMD1; norm(mean(mdataUM(:,CarcinoidFlag),2) - ...
                         mean(mdataUM(:,SmallCellFlag|NormalFlag|ColonFlag),2))] ;
      vMD2 = [vMD2; norm(mean(mdataUM(:,SmallCellFlag),2) - ...
                         mean(mdataUM(:,CarcinoidFlag|NormalFlag|ColonFlag),2))] ;
      vMD3 = [vMD3; norm(mean(mdataUM(:,NormalFlag),2) - ...
                         mean(mdataUM(:,CarcinoidFlag|SmallCellFlag|ColonFlag),2))] ;
      vMD4 = [vMD4; norm(mean(mdataUM(:,ColonFlag),2) - ...
                         mean(mdataUM(:,CarcinoidFlag|NormalFlag|SmallCellFlag),2))] ;

      %  Compute Maximal Data Piling Gap
      %
      if ngunm > n ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,CarcinoidFlag),mdataUM(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
        vMDP1 = [vMDP1; MDPgap] ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,SmallCellFlag),mdataUM(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
        vMDP2 = [vMDP2; MDPgap] ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,NormalFlag),mdataUM(:,CarcinoidFlag|SmallCellFlag|ColonFlag)) ;
        vMDP3 = [vMDP3; MDPgap] ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,ColonFlag),mdataUM(:,CarcinoidFlag|SmallCellFlag|NormalFlag)) ;
        vMDP4 = [vMDP4; MDPgap] ;
      else ;
        vMDP1 = [vMDP1; NaN] ;
        vMDP2 = [vMDP2; NaN] ;
        vMDP3 = [vMDP3; NaN] ;
        vMDP4 = [vMDP4; NaN] ;
      end ;

      %  Compute Projected Data
      %
      mdirUM = mdir .* (flagUM * ones(1,4)) ;
      mdirUM(:,1) = mdirUM(:,1) / norm(mdirUM(:,1)) ;
      mdirUM(:,2) = mdirUM(:,2) / norm(mdirUM(:,2)) ;
      mdirUM(:,3) = mdirUM(:,3) / norm(mdirUM(:,3)) ;
      mdirUM(:,4) = mdirUM(:,4) / norm(mdirUM(:,4)) ;
      mdataproj = mdirUM' * mdata ;
          %  projections of data, in Indiv vs. Rest directions

      %  Compute Projected Cluster Indices
      %
      vCIp1 = [vCIp1; ClustIndSM(mdataproj,CarcinoidFlag,SmallCellFlag|NormalFlag|ColonFlag)] ;
      vCIp2 = [vCIp2; ClustIndSM(mdataproj,SmallCellFlag,CarcinoidFlag|NormalFlag|ColonFlag)] ;
      vCIp3 = [vCIp3; ClustIndSM(mdataproj,NormalFlag,CarcinoidFlag|SmallCellFlag|ColonFlag)] ;
      vCIp4 = [vCIp4; ClustIndSM(mdataproj,ColonFlag,CarcinoidFlag|SmallCellFlag|NormalFlag)] ;

    end ;

    %  Save results in .mat file
    %
    save NeilLung1ip250.mat vngunm vCIf1 vCIf2 vCIf3 vCIf4 vMD1 vMD2 vMD3 vMD4 vMDP1 vMDP2 vMDP3 vMDP4 vCIp1 vCIp2 vCIp3 vCIp4 ;
        %  save these variables

    disp('Finished data save') ;


  elseif  ipart == 251  | ...
          ipart == 252  | ...
          ipart == 253  ;    %  Indiv - All, Gene Focussing, DWD SC vs. Rest, Plots of Cluster Indices

    load NeilLung1ip250 ;
        %  loads: 
        %      vngunm
        %      vCIf1 vCIf2 vCIf3 vCIf4
        %      vMD1 vMD2 vMD3 vMD4
        %      vMDP1 vMDP2 vMDP3 vMDP4
        %      vCIp1 vCIp2 vCIp3 vCIp4
    disp('  Finished Data Load') ;


    if ipart == 251 ;   % Indiv - All, Gene Focussing, DWD SC vs. Rest, Plots of Cluster Indices

      bottom = 0 ;
      top = 1 ;
      left = 0 ; 
      right = 2300 ;

      subplot(2,1,1) ;
      plot(vngunm,vCIf1,'m-') ;
      axis([left right bottom top]) ;
      hold on ;
        plot(vngunm,vCIf2,'b-') ;
        plot(vngunm,vCIf3,'g-') ;
        plot(vngunm,vCIf4,'r-') ;
        plot(vngunm,vCIp1,'m--') ;
        plot(vngunm,vCIp2,'b--') ;
        plot(vngunm,vCIp3,'g--') ;
        plot(vngunm,vCIp4,'r--') ;
        plot([1024; 1024],[bottom top],'k-') ;
        plot([512; 512],[bottom top],'k-') ;
        plot([256; 256],[bottom top],'k-') ;
        plot([128; 128],[bottom top],'k-') ;
        plot([64; 64],[bottom top],'k-') ;
        plot([32; 32],[bottom top],'k-') ;
        plot([16; 16],[bottom top],'k-') ;
        plot([8; 8],[bottom top],'k-') ;
        plot([4; 4],[bottom top],'k-') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.8 * (top - bottom), ...
             'Carcinoid vs. SmallCell&Normal&Colon', ...
             'Color','m') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.6 * (top - bottom), ...
             'SmallCell vs. Carcinoid&Normal&Colon', ...
             'Color','b') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.4 * (top - bottom), ...
             'Normal vs. Carcinoid&SmallCell&Colon', ...
             'Color','g') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.2 * (top - bottom), ...
             'Colon vs. Carcinoid&SmallCell&Normal', ...
             'Color','r') ;
      hold off ;
      title('Gene Focussing, based on DWD SC vs. Rest Direction') ;
      xlabel('Number of Genes Included') ;
      ylabel('Cluster Index (-- for Projected)') ;

      subplot(2,1,2) ;
      plot(log10(vngunm),vCIf1,'m-') ;
      right = log10(2300) ;
      axis([left right bottom top]) ;
      hold on ;
        plot(log10(vngunm),vCIf2,'b-') ;
        plot(log10(vngunm),vCIf3,'g-') ;
        plot(log10(vngunm),vCIf4,'r-') ;
        plot(log10(vngunm),vCIp1,'m--') ;
        plot(log10(vngunm),vCIp2,'b--') ;
        plot(log10(vngunm),vCIp3,'g--') ;
        plot(log10(vngunm),vCIp4,'r--') ;
        plot(log10([1024; 1024]),[bottom top],'k-') ;
        plot(log10([512; 512]),[bottom top],'k-') ;
        plot(log10([256; 256]),[bottom top],'k-') ;
        plot(log10([128; 128]),[bottom top],'k-') ;
        plot(log10([64; 64]),[bottom top],'k-') ;
        plot(log10([32; 32]),[bottom top],'k-') ;
        plot(log10([16; 16]),[bottom top],'k-') ;
        plot(log10([8; 8]),[bottom top],'k-') ;
        plot(log10([4; 4]),[bottom top],'k-') ;
      hold off ;
      title('Gene Focussing, based on DWD SC vs. Rest Direction') ;
      xlabel('log10(Number of Genes Included)') ;
      ylabel('Cluster Index (-- for Projected)') ;

      savestr = ['NeilLung1ip' num2str(ipart) 'GeneFocusDWDSS-RestClustInd'] ;
      orient landscape ;
      print('-dpsc2',savestr) ;


    elseif ipart == 252 ;   % Indiv - All, Gene Focussing, DWD SC vs. Rest, Plots of Mean Diffs

      vax = axisSM([vMD1 vMD2 vMD3 vMD4]) ;
      bottom = vax(1) ;
      top = vax(2) ;
      left = 0 ;
      right = 2300 ;

      subplot(2,1,1) ;
      plot(vngunm,vMD1,'m-') ;
      axis([left right bottom top]) ;
      hold on ;
        plot(vngunm,vMD2,'b-') ;
        plot(vngunm,vMD3,'g-') ;
        plot(vngunm,vMD4,'r-') ;
        plot([1024; 1024],[bottom top],'k-') ;
        plot([512; 512],[bottom top],'k-') ;
        plot([256; 256],[bottom top],'k-') ;
        plot([128; 128],[bottom top],'k-') ;
        plot([64; 64],[bottom top],'k-') ;
        plot([32; 32],[bottom top],'k-') ;
        plot([16; 16],[bottom top],'k-') ;
        plot([8; 8],[bottom top],'k-') ;
        plot([4; 4],[bottom top],'k-') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.8 * (top - bottom), ...
             'Carcinoid vs. SmallCell&Normal&Colon', ...
             'Color','m') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.6 * (top - bottom), ...
             'SmallCell vs. Carcinoid&Normal&Colon', ...
             'Color','b') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.4 * (top - bottom), ...
             'Normal vs. Carcinoid&SmallCell&Colon', ...
             'Color','g') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.2 * (top - bottom), ...
             'Colon vs. Carcinoid&SmallCell&Normal', ...
             'Color','r') ;
      hold off ;
      title('Gene Focussing, based on DWD SC vs. Rest Direction') ;
      xlabel('Number of Genes Included') ;
      ylabel('Mean Difference') ;

      subplot(2,1,2) ;
      plot(log10(vngunm),vMD1,'m-') ;
      right = log10(2300) ;
      axis([left right bottom top]) ;
      hold on ;
        plot(log10(vngunm),vMD2,'b-') ;
        plot(log10(vngunm),vMD3,'g-') ;
        plot(log10(vngunm),vMD4,'r-') ;
        plot(log10([1024; 1024]),[bottom top],'k-') ;
        plot(log10([512; 512]),[bottom top],'k-') ;
        plot(log10([256; 256]),[bottom top],'k-') ;
        plot(log10([128; 128]),[bottom top],'k-') ;
        plot(log10([64; 64]),[bottom top],'k-') ;
        plot(log10([32; 32]),[bottom top],'k-') ;
        plot(log10([16; 16]),[bottom top],'k-') ;
        plot(log10([8; 8]),[bottom top],'k-') ;
        plot(log10([4; 4]),[bottom top],'k-') ;
      hold off ;
      title('Gene Focussing, based on DWD SC vs. Rest Direction') ;
      xlabel('log10(Number of Genes Included)') ;
      ylabel('Mean Difference') ;

      savestr = ['NeilLung1ip' num2str(ipart) 'GeneFocusDWDSS-RestMeanDiff'] ;
      orient landscape ;
      print('-dpsc2',savestr) ;


    elseif ipart == 253 ;   % Indiv - All, Gene Focussing, DWD SC vs. Rest, Plots of MDPs

      vax = axisSM([vMDP1 vMDP2 vMDP3 vMDP4]) ;
      bottom = vax(1) ;
      top = vax(2) ;
      left = 0 ;
      right = 2300 ;

      subplot(2,1,1) ;
      plot(vngunm,vMDP1,'m-') ;
      axis([left right bottom top]) ;
      hold on ;
        plot(vngunm,vMDP2,'b-') ;
        plot(vngunm,vMDP3,'g-') ;
        plot(vngunm,vMDP4,'r-') ;
        plot([1024; 1024],[bottom top],'k-') ;
        plot([512; 512],[bottom top],'k-') ;
        plot([256; 256],[bottom top],'k-') ;
        plot([128; 128],[bottom top],'k-') ;
        plot([64; 64],[bottom top],'k-') ;
        plot([32; 32],[bottom top],'k-') ;
        plot([16; 16],[bottom top],'k-') ;
        plot([8; 8],[bottom top],'k-') ;
        plot([4; 4],[bottom top],'k-') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.8 * (top - bottom), ...
             'Carcinoid vs. SmallCell&Normal&Colon', ...
             'Color','m') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.6 * (top - bottom), ...
             'SmallCell vs. Carcinoid&Normal&Colon', ...
             'Color','b') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.4 * (top - bottom), ...
             'Normal vs. Carcinoid&SmallCell&Colon', ...
             'Color','g') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.2 * (top - bottom), ...
             'Colon vs. Carcinoid&SmallCell&Normal', ...
             'Color','r') ;
      hold off ;
      title('Gene Focussing, based on DWD SC vs. Rest Direction') ;
      xlabel('Number of Genes Included') ;
      ylabel('Maximal Data Piling Gap') ;

      subplot(2,1,2) ;
      plot(log10(vngunm),vMDP1,'m-') ;
      right = log10(2300) ;
      axis([left right bottom top]) ;
      hold on ;
        plot(log10(vngunm),vMDP2,'b-') ;
        plot(log10(vngunm),vMDP3,'g-') ;
        plot(log10(vngunm),vMDP4,'r-') ;
        plot(log10([1024; 1024]),[bottom top],'k-') ;
        plot(log10([512; 512]),[bottom top],'k-') ;
        plot(log10([256; 256]),[bottom top],'k-') ;
        plot(log10([128; 128]),[bottom top],'k-') ;
        plot(log10([64; 64]),[bottom top],'k-') ;
        plot(log10([32; 32]),[bottom top],'k-') ;
        plot(log10([16; 16]),[bottom top],'k-') ;
        plot(log10([8; 8]),[bottom top],'k-') ;
        plot(log10([4; 4]),[bottom top],'k-') ;
      hold off ;
      title('Gene Focussing, based on DWD SC vs. Rest Direction') ;
      xlabel('log10(Number of Genes Included)') ;
      ylabel('Maximal Data Piling Gap') ;

      savestr = ['NeilLung1ip' num2str(ipart) 'GeneFocusDWDSS-RestMDPgap'] ;
      orient landscape ;
      print('-dpsc2',savestr) ;


    end ;


  elseif  ipart == 260  ;    %  Indiv - All, Gene Focussing, DWD Norm vs. Rest, Save indices

    mdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag))] ;
    vdir = DWD1SM(mdata(:,NormalFlag),mdata(:,CarcinoidFlag|SmallCellFlag|ColonFlag)) ;
    mdir = [mdir vdir] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,CarcinoidFlag|SmallCellFlag|NormalFlag))] ;

    vngunm = [(1:10) (12:2:20) (25:5:50) (50:10:200) (220:20:400) (450:50:1000) (1100:100:2300)]' ;
        %  vector of numbers of genes to unmask 

    vCIf1 = [] ;
    vCIf2 = [] ;
    vCIf3 = [] ;
    vCIf4 = [] ;
    vMD1 = [] ;
    vMD2 = [] ;
    vMD3 = [] ;
    vMD4 = [] ;
    vMDP1 = [] ;
    vMDP2 = [] ;
    vMDP3 = [] ;
    vMDP4 = [] ;
    vCIp1 = [] ;
    vCIp2 = [] ;
    vCIp3 = [] ;
    vCIp4 = [] ;
    for ngunm = vngunm' ;

      disp(['    Working on ngunm = ' num2str(ngunm) ' of ' num2str(vngunm(end))]) ;

      [savdir,vi] = sort(abs(vdir)) ;
      viUM = vi((d - ngunm + 1):d) ;
      flagUM = zeros(d,1) ;
      flagUM(viUM) = ones(ngunm,1) ; 

      mdataUM = mdata(viUM,:) ;

      %  Compute Full Space Cluster Indices
      %
      vCIf1 = [vCIf1; ClustIndSM(mdataUM,CarcinoidFlag,SmallCellFlag|NormalFlag|ColonFlag)] ;
      vCIf2 = [vCIf2; ClustIndSM(mdataUM,SmallCellFlag,CarcinoidFlag|NormalFlag|ColonFlag)] ;
      vCIf3 = [vCIf3; ClustIndSM(mdataUM,NormalFlag,CarcinoidFlag|SmallCellFlag|ColonFlag)] ;
      vCIf4 = [vCIf4; ClustIndSM(mdataUM,ColonFlag,CarcinoidFlag|SmallCellFlag|NormalFlag)] ;

      %  Compute Mean Differences
      %
      vMD1 = [vMD1; norm(mean(mdataUM(:,CarcinoidFlag),2) - ...
                         mean(mdataUM(:,SmallCellFlag|NormalFlag|ColonFlag),2))] ;
      vMD2 = [vMD2; norm(mean(mdataUM(:,SmallCellFlag),2) - ...
                         mean(mdataUM(:,CarcinoidFlag|NormalFlag|ColonFlag),2))] ;
      vMD3 = [vMD3; norm(mean(mdataUM(:,NormalFlag),2) - ...
                         mean(mdataUM(:,CarcinoidFlag|SmallCellFlag|ColonFlag),2))] ;
      vMD4 = [vMD4; norm(mean(mdataUM(:,ColonFlag),2) - ...
                         mean(mdataUM(:,CarcinoidFlag|NormalFlag|SmallCellFlag),2))] ;

      %  Compute Maximal Data Piling Gap
      %
      if ngunm > n ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,CarcinoidFlag),mdataUM(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
        vMDP1 = [vMDP1; MDPgap] ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,SmallCellFlag),mdataUM(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
        vMDP2 = [vMDP2; MDPgap] ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,NormalFlag),mdataUM(:,CarcinoidFlag|SmallCellFlag|ColonFlag)) ;
        vMDP3 = [vMDP3; MDPgap] ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,ColonFlag),mdataUM(:,CarcinoidFlag|SmallCellFlag|NormalFlag)) ;
        vMDP4 = [vMDP4; MDPgap] ;
      else ;
        vMDP1 = [vMDP1; NaN] ;
        vMDP2 = [vMDP2; NaN] ;
        vMDP3 = [vMDP3; NaN] ;
        vMDP4 = [vMDP4; NaN] ;
      end ;

      %  Compute Projected Data
      %
      mdirUM = mdir .* (flagUM * ones(1,4)) ;
      mdirUM(:,1) = mdirUM(:,1) / norm(mdirUM(:,1)) ;
      mdirUM(:,2) = mdirUM(:,2) / norm(mdirUM(:,2)) ;
      mdirUM(:,3) = mdirUM(:,3) / norm(mdirUM(:,3)) ;
      mdirUM(:,4) = mdirUM(:,4) / norm(mdirUM(:,4)) ;
      mdataproj = mdirUM' * mdata ;
          %  projections of data, in Indiv vs. Rest directions

      %  Compute Projected Cluster Indices
      %
      vCIp1 = [vCIp1; ClustIndSM(mdataproj,CarcinoidFlag,SmallCellFlag|NormalFlag|ColonFlag)] ;
      vCIp2 = [vCIp2; ClustIndSM(mdataproj,SmallCellFlag,CarcinoidFlag|NormalFlag|ColonFlag)] ;
      vCIp3 = [vCIp3; ClustIndSM(mdataproj,NormalFlag,CarcinoidFlag|SmallCellFlag|ColonFlag)] ;
      vCIp4 = [vCIp4; ClustIndSM(mdataproj,ColonFlag,CarcinoidFlag|SmallCellFlag|NormalFlag)] ;

    end ;

    %  Save results in .mat file
    %
    save NeilLung1ip260.mat vngunm vCIf1 vCIf2 vCIf3 vCIf4 vMD1 vMD2 vMD3 vMD4 vMDP1 vMDP2 vMDP3 vMDP4 vCIp1 vCIp2 vCIp3 vCIp4 ;
        %  save these variables

    disp('Finished data save') ;


  elseif  ipart == 261  | ...
          ipart == 262  | ...
          ipart == 263  ;    %  Indiv - All, Gene Focussing, DWD Norm vs. Rest, Plots of Cluster Indices

    load NeilLung1ip260 ;
        %  loads: 
        %      vngunm
        %      vCIf1 vCIf2 vCIf3 vCIf4
        %      vMD1 vMD2 vMD3 vMD4
        %      vMDP1 vMDP2 vMDP3 vMDP4
        %      vCIp1 vCIp2 vCIp3 vCIp4
    disp('  Finished Data Load') ;


    if ipart == 261 ;   % Indiv - All, Gene Focussing, DWD Norm vs. Rest, Plots of Cluster Indices

      bottom = 0 ;
      top = 1 ;
      left = 0 ; 
      right = 2300 ;

      subplot(2,1,1) ;
      plot(vngunm,vCIf1,'m-') ;
      axis([left right bottom top]) ;
      hold on ;
        plot(vngunm,vCIf2,'b-') ;
        plot(vngunm,vCIf3,'g-') ;
        plot(vngunm,vCIf4,'r-') ;
        plot(vngunm,vCIp1,'m--') ;
        plot(vngunm,vCIp2,'b--') ;
        plot(vngunm,vCIp3,'g--') ;
        plot(vngunm,vCIp4,'r--') ;
        plot([1024; 1024],[bottom top],'k-') ;
        plot([512; 512],[bottom top],'k-') ;
        plot([256; 256],[bottom top],'k-') ;
        plot([128; 128],[bottom top],'k-') ;
        plot([64; 64],[bottom top],'k-') ;
        plot([32; 32],[bottom top],'k-') ;
        plot([16; 16],[bottom top],'k-') ;
        plot([8; 8],[bottom top],'k-') ;
        plot([4; 4],[bottom top],'k-') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.8 * (top - bottom), ...
             'Carcinoid vs. SmallCell&Normal&Colon', ...
             'Color','m') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.6 * (top - bottom), ...
             'SmallCell vs. Carcinoid&Normal&Colon', ...
             'Color','b') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.4 * (top - bottom), ...
             'Normal vs. Carcinoid&SmallCell&Colon', ...
             'Color','g') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.2 * (top - bottom), ...
             'Colon vs. Carcinoid&SmallCell&Normal', ...
             'Color','r') ;
      hold off ;
      title('Gene Focussing, based on DWD Norm vs. Rest Direction') ;
      xlabel('Number of Genes Included') ;
      ylabel('Cluster Index (-- for Projected)') ;

      subplot(2,1,2) ;
      plot(log10(vngunm),vCIf1,'m-') ;
      right = log10(2300) ;
      axis([left right bottom top]) ;
      hold on ;
        plot(log10(vngunm),vCIf2,'b-') ;
        plot(log10(vngunm),vCIf3,'g-') ;
        plot(log10(vngunm),vCIf4,'r-') ;
        plot(log10(vngunm),vCIp1,'m--') ;
        plot(log10(vngunm),vCIp2,'b--') ;
        plot(log10(vngunm),vCIp3,'g--') ;
        plot(log10(vngunm),vCIp4,'r--') ;
        plot(log10([1024; 1024]),[bottom top],'k-') ;
        plot(log10([512; 512]),[bottom top],'k-') ;
        plot(log10([256; 256]),[bottom top],'k-') ;
        plot(log10([128; 128]),[bottom top],'k-') ;
        plot(log10([64; 64]),[bottom top],'k-') ;
        plot(log10([32; 32]),[bottom top],'k-') ;
        plot(log10([16; 16]),[bottom top],'k-') ;
        plot(log10([8; 8]),[bottom top],'k-') ;
        plot(log10([4; 4]),[bottom top],'k-') ;
      hold off ;
      title('Gene Focussing, based on DWD Norm vs. Rest Direction') ;
      xlabel('log10(Number of Genes Included)') ;
      ylabel('Cluster Index (-- for Projected)') ;

      savestr = ['NeilLung1ip' num2str(ipart) 'GeneFocusDWDNorm-RestClustInd'] ;
      orient landscape ;
      print('-dpsc2',savestr) ;


    elseif ipart == 262 ;   % Indiv - All, Gene Focussing, DWD Norm vs. Rest, Plots of Mean Diffs

      vax = axisSM([vMD1 vMD2 vMD3 vMD4]) ;
      bottom = vax(1) ;
      top = vax(2) ;
      left = 0 ;
      right = 2300 ;

      subplot(2,1,1) ;
      plot(vngunm,vMD1,'m-') ;
      axis([left right bottom top]) ;
      hold on ;
        plot(vngunm,vMD2,'b-') ;
        plot(vngunm,vMD3,'g-') ;
        plot(vngunm,vMD4,'r-') ;
        plot([1024; 1024],[bottom top],'k-') ;
        plot([512; 512],[bottom top],'k-') ;
        plot([256; 256],[bottom top],'k-') ;
        plot([128; 128],[bottom top],'k-') ;
        plot([64; 64],[bottom top],'k-') ;
        plot([32; 32],[bottom top],'k-') ;
        plot([16; 16],[bottom top],'k-') ;
        plot([8; 8],[bottom top],'k-') ;
        plot([4; 4],[bottom top],'k-') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.8 * (top - bottom), ...
             'Carcinoid vs. SmallCell&Normal&Colon', ...
             'Color','m') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.6 * (top - bottom), ...
             'SmallCell vs. Carcinoid&Normal&Colon', ...
             'Color','b') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.4 * (top - bottom), ...
             'Normal vs. Carcinoid&SmallCell&Colon', ...
             'Color','g') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.2 * (top - bottom), ...
             'Colon vs. Carcinoid&SmallCell&Normal', ...
             'Color','r') ;
      hold off ;
      title('Gene Focussing, based on DWD Norm vs. Rest Direction') ;
      xlabel('Number of Genes Included') ;
      ylabel('Mean Difference') ;

      subplot(2,1,2) ;
      plot(log10(vngunm),vMD1,'m-') ;
      right = log10(2300) ;
      axis([left right bottom top]) ;
      hold on ;
        plot(log10(vngunm),vMD2,'b-') ;
        plot(log10(vngunm),vMD3,'g-') ;
        plot(log10(vngunm),vMD4,'r-') ;
        plot(log10([1024; 1024]),[bottom top],'k-') ;
        plot(log10([512; 512]),[bottom top],'k-') ;
        plot(log10([256; 256]),[bottom top],'k-') ;
        plot(log10([128; 128]),[bottom top],'k-') ;
        plot(log10([64; 64]),[bottom top],'k-') ;
        plot(log10([32; 32]),[bottom top],'k-') ;
        plot(log10([16; 16]),[bottom top],'k-') ;
        plot(log10([8; 8]),[bottom top],'k-') ;
        plot(log10([4; 4]),[bottom top],'k-') ;
      hold off ;
      title('Gene Focussing, based on DWD Norm vs. Rest Direction') ;
      xlabel('log10(Number of Genes Included)') ;
      ylabel('Mean Difference') ;

      savestr = ['NeilLung1ip' num2str(ipart) 'GeneFocusDWDNorm-RestMeanDiff'] ;
      orient landscape ;
      print('-dpsc2',savestr) ;


    elseif ipart == 263 ;   % Indiv - All, Gene Focussing, DWD Norm vs. Rest, Plots of MDPs

      vax = axisSM([vMDP1 vMDP2 vMDP3 vMDP4]) ;
      bottom = vax(1) ;
      top = vax(2) ;
      left = 0 ;
      right = 2300 ;

      subplot(2,1,1) ;
      plot(vngunm,vMDP1,'m-') ;
      axis([left right bottom top]) ;
      hold on ;
        plot(vngunm,vMDP2,'b-') ;
        plot(vngunm,vMDP3,'g-') ;
        plot(vngunm,vMDP4,'r-') ;
        plot([1024; 1024],[bottom top],'k-') ;
        plot([512; 512],[bottom top],'k-') ;
        plot([256; 256],[bottom top],'k-') ;
        plot([128; 128],[bottom top],'k-') ;
        plot([64; 64],[bottom top],'k-') ;
        plot([32; 32],[bottom top],'k-') ;
        plot([16; 16],[bottom top],'k-') ;
        plot([8; 8],[bottom top],'k-') ;
        plot([4; 4],[bottom top],'k-') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.8 * (top - bottom), ...
             'Carcinoid vs. SmallCell&Normal&Colon', ...
             'Color','m') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.6 * (top - bottom), ...
             'SmallCell vs. Carcinoid&Normal&Colon', ...
             'Color','b') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.4 * (top - bottom), ...
             'Normal vs. Carcinoid&SmallCell&Colon', ...
             'Color','g') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.2 * (top - bottom), ...
             'Colon vs. Carcinoid&SmallCell&Normal', ...
             'Color','r') ;
      hold off ;
      title('Gene Focussing, based on DWD Norm vs. Rest Direction') ;
      xlabel('Number of Genes Included') ;
      ylabel('Maximal Data Piling Gap') ;

      subplot(2,1,2) ;
      plot(log10(vngunm),vMDP1,'m-') ;
      right = log10(2300) ;
      axis([left right bottom top]) ;
      hold on ;
        plot(log10(vngunm),vMDP2,'b-') ;
        plot(log10(vngunm),vMDP3,'g-') ;
        plot(log10(vngunm),vMDP4,'r-') ;
        plot(log10([1024; 1024]),[bottom top],'k-') ;
        plot(log10([512; 512]),[bottom top],'k-') ;
        plot(log10([256; 256]),[bottom top],'k-') ;
        plot(log10([128; 128]),[bottom top],'k-') ;
        plot(log10([64; 64]),[bottom top],'k-') ;
        plot(log10([32; 32]),[bottom top],'k-') ;
        plot(log10([16; 16]),[bottom top],'k-') ;
        plot(log10([8; 8]),[bottom top],'k-') ;
        plot(log10([4; 4]),[bottom top],'k-') ;
      hold off ;
      title('Gene Focussing, based on DWD Norm vs. Rest Direction') ;
      xlabel('log10(Number of Genes Included)') ;
      ylabel('Maximal Data Piling Gap') ;

      savestr = ['NeilLung1ip' num2str(ipart) 'GeneFocusDWDNorm-RestMDPgap'] ;
      orient landscape ;
      print('-dpsc2',savestr) ;


    end ;


  elseif  ipart == 270  ;    %  Indiv - All, Gene Focussing, DWD Col vs. Rest, Save indices

    mdir = DWD1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,CarcinoidFlag|SmallCellFlag|ColonFlag))] ;
    vdir = DWD1SM(mdata(:,ColonFlag),mdata(:,CarcinoidFlag|SmallCellFlag|NormalFlag)) ;
    mdir = [mdir vdir] ;

    vngunm = [(1:10) (12:2:20) (25:5:50) (50:10:200) (220:20:400) (450:50:1000) (1100:100:2300)]' ;
        %  vector of numbers of genes to unmask 

    vCIf1 = [] ;
    vCIf2 = [] ;
    vCIf3 = [] ;
    vCIf4 = [] ;
    vMD1 = [] ;
    vMD2 = [] ;
    vMD3 = [] ;
    vMD4 = [] ;
    vMDP1 = [] ;
    vMDP2 = [] ;
    vMDP3 = [] ;
    vMDP4 = [] ;
    vCIp1 = [] ;
    vCIp2 = [] ;
    vCIp3 = [] ;
    vCIp4 = [] ;
    for ngunm = vngunm' ;

      disp(['    Working on ngunm = ' num2str(ngunm) ' of ' num2str(vngunm(end))]) ;

      [savdir,vi] = sort(abs(vdir)) ;
      viUM = vi((d - ngunm + 1):d) ;
      flagUM = zeros(d,1) ;
      flagUM(viUM) = ones(ngunm,1) ; 

      mdataUM = mdata(viUM,:) ;

      %  Compute Full Space Cluster Indices
      %
      vCIf1 = [vCIf1; ClustIndSM(mdataUM,CarcinoidFlag,SmallCellFlag|NormalFlag|ColonFlag)] ;
      vCIf2 = [vCIf2; ClustIndSM(mdataUM,SmallCellFlag,CarcinoidFlag|NormalFlag|ColonFlag)] ;
      vCIf3 = [vCIf3; ClustIndSM(mdataUM,NormalFlag,CarcinoidFlag|SmallCellFlag|ColonFlag)] ;
      vCIf4 = [vCIf4; ClustIndSM(mdataUM,ColonFlag,CarcinoidFlag|SmallCellFlag|NormalFlag)] ;

      %  Compute Mean Differences
      %
      vMD1 = [vMD1; norm(mean(mdataUM(:,CarcinoidFlag),2) - ...
                         mean(mdataUM(:,SmallCellFlag|NormalFlag|ColonFlag),2))] ;
      vMD2 = [vMD2; norm(mean(mdataUM(:,SmallCellFlag),2) - ...
                         mean(mdataUM(:,CarcinoidFlag|NormalFlag|ColonFlag),2))] ;
      vMD3 = [vMD3; norm(mean(mdataUM(:,NormalFlag),2) - ...
                         mean(mdataUM(:,CarcinoidFlag|SmallCellFlag|ColonFlag),2))] ;
      vMD4 = [vMD4; norm(mean(mdataUM(:,ColonFlag),2) - ...
                         mean(mdataUM(:,CarcinoidFlag|NormalFlag|SmallCellFlag),2))] ;

      %  Compute Maximal Data Piling Gap
      %
      if ngunm > n ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,CarcinoidFlag),mdataUM(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
        vMDP1 = [vMDP1; MDPgap] ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,SmallCellFlag),mdataUM(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
        vMDP2 = [vMDP2; MDPgap] ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,NormalFlag),mdataUM(:,CarcinoidFlag|SmallCellFlag|ColonFlag)) ;
        vMDP3 = [vMDP3; MDPgap] ;
        [MDPdir,MDPgap] = MaxDatPilJA(mdataUM(:,ColonFlag),mdataUM(:,CarcinoidFlag|SmallCellFlag|NormalFlag)) ;
        vMDP4 = [vMDP4; MDPgap] ;
      else ;
        vMDP1 = [vMDP1; NaN] ;
        vMDP2 = [vMDP2; NaN] ;
        vMDP3 = [vMDP3; NaN] ;
        vMDP4 = [vMDP4; NaN] ;
      end ;

      %  Compute Projected Data
      %
      mdirUM = mdir .* (flagUM * ones(1,4)) ;
      mdirUM(:,1) = mdirUM(:,1) / norm(mdirUM(:,1)) ;
      mdirUM(:,2) = mdirUM(:,2) / norm(mdirUM(:,2)) ;
      mdirUM(:,3) = mdirUM(:,3) / norm(mdirUM(:,3)) ;
      mdirUM(:,4) = mdirUM(:,4) / norm(mdirUM(:,4)) ;
      mdataproj = mdirUM' * mdata ;
          %  projections of data, in Indiv vs. Rest directions

      %  Compute Projected Cluster Indices
      %
      vCIp1 = [vCIp1; ClustIndSM(mdataproj,CarcinoidFlag,SmallCellFlag|NormalFlag|ColonFlag)] ;
      vCIp2 = [vCIp2; ClustIndSM(mdataproj,SmallCellFlag,CarcinoidFlag|NormalFlag|ColonFlag)] ;
      vCIp3 = [vCIp3; ClustIndSM(mdataproj,NormalFlag,CarcinoidFlag|SmallCellFlag|ColonFlag)] ;
      vCIp4 = [vCIp4; ClustIndSM(mdataproj,ColonFlag,CarcinoidFlag|SmallCellFlag|NormalFlag)] ;

    end ;

    %  Save results in .mat file
    %
    save NeilLung1ip270.mat vngunm vCIf1 vCIf2 vCIf3 vCIf4 vMD1 vMD2 vMD3 vMD4 vMDP1 vMDP2 vMDP3 vMDP4 vCIp1 vCIp2 vCIp3 vCIp4 ;
        %  save these variables

    disp('Finished data save') ;


  elseif  ipart == 271  | ...
          ipart == 272  | ...
          ipart == 273  ;    %  Indiv - All, Gene Focussing, DWD Col vs. Rest, Plots of Cluster Indices

    load NeilLung1ip270 ;
        %  loads: 
        %      vngunm
        %      vCIf1 vCIf2 vCIf3 vCIf4
        %      vMD1 vMD2 vMD3 vMD4
        %      vMDP1 vMDP2 vMDP3 vMDP4
        %      vCIp1 vCIp2 vCIp3 vCIp4
    disp('  Finished Data Load') ;


    if ipart == 271 ;   % Indiv - All, Gene Focussing, DWD Col vs. Rest, Plots of Cluster Indices

      bottom = 0 ;
      top = 1 ;
      left = 0 ; 
      right = 2300 ;

      subplot(2,1,1) ;
      plot(vngunm,vCIf1,'m-') ;
      axis([left right bottom top]) ;
      hold on ;
        plot(vngunm,vCIf2,'b-') ;
        plot(vngunm,vCIf3,'g-') ;
        plot(vngunm,vCIf4,'r-') ;
        plot(vngunm,vCIp1,'m--') ;
        plot(vngunm,vCIp2,'b--') ;
        plot(vngunm,vCIp3,'g--') ;
        plot(vngunm,vCIp4,'r--') ;
        plot([1024; 1024],[bottom top],'k-') ;
        plot([512; 512],[bottom top],'k-') ;
        plot([256; 256],[bottom top],'k-') ;
        plot([128; 128],[bottom top],'k-') ;
        plot([64; 64],[bottom top],'k-') ;
        plot([32; 32],[bottom top],'k-') ;
        plot([16; 16],[bottom top],'k-') ;
        plot([8; 8],[bottom top],'k-') ;
        plot([4; 4],[bottom top],'k-') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.8 * (top - bottom), ...
             'Carcinoid vs. SmallCell&Normal&Colon', ...
             'Color','m') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.6 * (top - bottom), ...
             'SmallCell vs. Carcinoid&Normal&Colon', ...
             'Color','b') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.4 * (top - bottom), ...
             'Normal vs. Carcinoid&SmallCell&Colon', ...
             'Color','g') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.2 * (top - bottom), ...
             'Colon vs. Carcinoid&SmallCell&Normal', ...
             'Color','r') ;
      hold off ;
      title('Gene Focussing, based on DWD Col vs. Rest Direction') ;
      xlabel('Number of Genes Included') ;
      ylabel('Cluster Index (-- for Projected)') ;

      subplot(2,1,2) ;
      plot(log10(vngunm),vCIf1,'m-') ;
      right = log10(2300) ;
      axis([left right bottom top]) ;
      hold on ;
        plot(log10(vngunm),vCIf2,'b-') ;
        plot(log10(vngunm),vCIf3,'g-') ;
        plot(log10(vngunm),vCIf4,'r-') ;
        plot(log10(vngunm),vCIp1,'m--') ;
        plot(log10(vngunm),vCIp2,'b--') ;
        plot(log10(vngunm),vCIp3,'g--') ;
        plot(log10(vngunm),vCIp4,'r--') ;
        plot(log10([1024; 1024]),[bottom top],'k-') ;
        plot(log10([512; 512]),[bottom top],'k-') ;
        plot(log10([256; 256]),[bottom top],'k-') ;
        plot(log10([128; 128]),[bottom top],'k-') ;
        plot(log10([64; 64]),[bottom top],'k-') ;
        plot(log10([32; 32]),[bottom top],'k-') ;
        plot(log10([16; 16]),[bottom top],'k-') ;
        plot(log10([8; 8]),[bottom top],'k-') ;
        plot(log10([4; 4]),[bottom top],'k-') ;
      hold off ;
      title('Gene Focussing, based on DWD Col vs. Rest Direction') ;
      xlabel('log10(Number of Genes Included)') ;
      ylabel('Cluster Index (-- for Projected)') ;

      savestr = ['NeilLung1ip' num2str(ipart) 'GeneFocusDWDCol-RestClustInd'] ;
      orient landscape ;
      print('-dpsc2',savestr) ;


    elseif ipart == 272 ;   % Indiv - All, Gene Focussing, DWD Col vs. Rest, Plots of Mean Diffs

      vax = axisSM([vMD1 vMD2 vMD3 vMD4]) ;
      bottom = vax(1) ;
      top = vax(2) ;
      left = 0 ;
      right = 2300 ;

      subplot(2,1,1) ;
      plot(vngunm,vMD1,'m-') ;
      axis([left right bottom top]) ;
      hold on ;
        plot(vngunm,vMD2,'b-') ;
        plot(vngunm,vMD3,'g-') ;
        plot(vngunm,vMD4,'r-') ;
        plot([1024; 1024],[bottom top],'k-') ;
        plot([512; 512],[bottom top],'k-') ;
        plot([256; 256],[bottom top],'k-') ;
        plot([128; 128],[bottom top],'k-') ;
        plot([64; 64],[bottom top],'k-') ;
        plot([32; 32],[bottom top],'k-') ;
        plot([16; 16],[bottom top],'k-') ;
        plot([8; 8],[bottom top],'k-') ;
        plot([4; 4],[bottom top],'k-') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.8 * (top - bottom), ...
             'Carcinoid vs. SmallCell&Normal&Colon', ...
             'Color','m') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.6 * (top - bottom), ...
             'SmallCell vs. Carcinoid&Normal&Colon', ...
             'Color','b') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.4 * (top - bottom), ...
             'Normal vs. Carcinoid&SmallCell&Colon', ...
             'Color','g') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.2 * (top - bottom), ...
             'Colon vs. Carcinoid&SmallCell&Normal', ...
             'Color','r') ;
      hold off ;
      title('Gene Focussing, based on DWD Col vs. Rest Direction') ;
      xlabel('Number of Genes Included') ;
      ylabel('Mean Difference') ;

      subplot(2,1,2) ;
      plot(log10(vngunm),vMD1,'m-') ;
      right = log10(2300) ;
      axis([left right bottom top]) ;
      hold on ;
        plot(log10(vngunm),vMD2,'b-') ;
        plot(log10(vngunm),vMD3,'g-') ;
        plot(log10(vngunm),vMD4,'r-') ;
        plot(log10([1024; 1024]),[bottom top],'k-') ;
        plot(log10([512; 512]),[bottom top],'k-') ;
        plot(log10([256; 256]),[bottom top],'k-') ;
        plot(log10([128; 128]),[bottom top],'k-') ;
        plot(log10([64; 64]),[bottom top],'k-') ;
        plot(log10([32; 32]),[bottom top],'k-') ;
        plot(log10([16; 16]),[bottom top],'k-') ;
        plot(log10([8; 8]),[bottom top],'k-') ;
        plot(log10([4; 4]),[bottom top],'k-') ;
      hold off ;
      title('Gene Focussing, based on DWD Col vs. Rest Direction') ;
      xlabel('log10(Number of Genes Included)') ;
      ylabel('Mean Difference') ;

      savestr = ['NeilLung1ip' num2str(ipart) 'GeneFocusDWDCol-RestMeanDiff'] ;
      orient landscape ;
      print('-dpsc2',savestr) ;


    elseif ipart == 273 ;   % Indiv - All, Gene Focussing, DWD Col vs. Rest, Plots of MDPs

      vax = axisSM([vMDP1 vMDP2 vMDP3 vMDP4]) ;
      bottom = vax(1) ;
      top = vax(2) ;
      left = 0 ;
      right = 2300 ;

      subplot(2,1,1) ;
      plot(vngunm,vMDP1,'m-') ;
      axis([left right bottom top]) ;
      hold on ;
        plot(vngunm,vMDP2,'b-') ;
        plot(vngunm,vMDP3,'g-') ;
        plot(vngunm,vMDP4,'r-') ;
        plot([1024; 1024],[bottom top],'k-') ;
        plot([512; 512],[bottom top],'k-') ;
        plot([256; 256],[bottom top],'k-') ;
        plot([128; 128],[bottom top],'k-') ;
        plot([64; 64],[bottom top],'k-') ;
        plot([32; 32],[bottom top],'k-') ;
        plot([16; 16],[bottom top],'k-') ;
        plot([8; 8],[bottom top],'k-') ;
        plot([4; 4],[bottom top],'k-') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.8 * (top - bottom), ...
             'Carcinoid vs. SmallCell&Normal&Colon', ...
             'Color','m') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.6 * (top - bottom), ...
             'SmallCell vs. Carcinoid&Normal&Colon', ...
             'Color','b') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.4 * (top - bottom), ...
             'Normal vs. Carcinoid&SmallCell&Colon', ...
             'Color','g') ;
        text(left + 0.6 * (right - left), ...
             bottom + 0.2 * (top - bottom), ...
             'Colon vs. Carcinoid&SmallCell&Normal', ...
             'Color','r') ;
      hold off ;
      title('Gene Focussing, based on DWD Col vs. Rest Direction') ;
      xlabel('Number of Genes Included') ;
      ylabel('Maximal Data Piling Gap') ;

      subplot(2,1,2) ;
      plot(log10(vngunm),vMDP1,'m-') ;
      right = log10(2300) ;
      axis([left right bottom top]) ;
      hold on ;
        plot(log10(vngunm),vMDP2,'b-') ;
        plot(log10(vngunm),vMDP3,'g-') ;
        plot(log10(vngunm),vMDP4,'r-') ;
        plot(log10([1024; 1024]),[bottom top],'k-') ;
        plot(log10([512; 512]),[bottom top],'k-') ;
        plot(log10([256; 256]),[bottom top],'k-') ;
        plot(log10([128; 128]),[bottom top],'k-') ;
        plot(log10([64; 64]),[bottom top],'k-') ;
        plot(log10([32; 32]),[bottom top],'k-') ;
        plot(log10([16; 16]),[bottom top],'k-') ;
        plot(log10([8; 8]),[bottom top],'k-') ;
        plot(log10([4; 4]),[bottom top],'k-') ;
      hold off ;
      title('Gene Focussing, based on DWD Col vs. Rest Direction') ;
      xlabel('log10(Number of Genes Included)') ;
      ylabel('Maximal Data Piling Gap') ;

      savestr = ['NeilLung1ip' num2str(ipart) 'GeneFocusDWDCol-RestMDPgap'] ;
      orient landscape ;
      print('-dpsc2',savestr) ;


    end ;


  end ;    %  of inner ipart if-block


  if  ipart == 301  | ...
      ipart == 302  ;    %  Remove Top 3 PCs


    PCAparamstruct = struct('npc',3,...
                            'viout',[0 1]) ;
    outstruct = pcaSM(mdata,PCAparamstruct) ;
    mdir = getfield(outstruct,'meigvec') ;

    mproj = mdir * pinv(mdir' * mdir) * mdir' ;
    mdatasub = mdata - mproj * mdata ;
        %  data projected on subspace ortho to top 3 PCs


    if  ipart == 301 ;    %  Make PCA plot

      titlecellstr = {{'PCA View, Ortho to PCs 123' 'Neil Hayes Lung data' 'Carc, Col, Norm, SC only' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

      labelcellstr = {{'PC 4' 'PC 5' 'PC 6' 'PC 7'}} ;

      savestr = ['NeilLung1ip' num2str(ipart) 'OrthoPC123-PCASP'] ;

      paramstruct = struct('icolor',icolor, ...
                           'markerstr',markerstr, ...
                           'isubpopkde',1, ...
                           'npcadiradd',4, ...
                           'legendcellstr',legendcellstr, ...
                           'mlegendcolor',mlegendcolor, ...
                           'titlecellstr',titlecellstr, ...
                           'labelcellstr',labelcellstr, ...
                           'savestr',savestr, ...
                           'iscreenwrite',1) ;

      scatplotSM(mdatasub,[],paramstruct) ;


    elseif  ipart == 302 ;    %  Proj on DWD Indiv. vs all plot

      mdir = DWD1SM(mdatasub(:,CarcinoidFlag),mdatasub(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
      mdir = [mdir DWD1SM(mdatasub(:,SmallCellFlag),mdatasub(:,CarcinoidFlag|NormalFlag|ColonFlag))] ;
      mdir = [mdir DWD1SM(mdatasub(:,NormalFlag),mdatasub(:,CarcinoidFlag|SmallCellFlag|ColonFlag))] ;
      mdir = [mdir DWD1SM(mdatasub(:,ColonFlag),mdatasub(:,NormalFlag|CarcinoidFlag|SmallCellFlag))] ;

      titlecellstr = {{'DWD View, Ortho to PCs 123' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

      labelcellstr = {{'Ortho DWD Carc-all ' 'Ortho DWD SC-all' 'Ortho DWD Norm-all' 'Ortho DWD Col-all'}} ;

      savestr = ['NeilLung1ip' num2str(ipart) 'OrthoPC123-DWDIndiv-all'] ;

      paramstruct = struct('icolor',icolor, ...
                           'markerstr',markerstr, ...
                           'isubpopkde',1, ...
                           'legendcellstr',legendcellstr, ...
                           'mlegendcolor',mlegendcolor, ...
                           'titlecellstr',titlecellstr, ...
                           'labelcellstr',labelcellstr, ...
                           'savestr',savestr, ...
                           'iscreenwrite',1) ;

      scatplotSM(mdatasub,mdir,paramstruct) ;


    end ;


  elseif  ipart == 311  | ...
          ipart == 312 ;    %  Remove Tree 1 subspace

    mdir = DWD1SM(mdata(:,CarcinoidFlag|NormalFlag),mdata(:,SmallCellFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,CarcinoidFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,SmallCellFlag),mdata(:,ColonFlag))] ;

    mproj = mdir * pinv(mdir' * mdir) * mdir' ;
    mdatasub = mdata - mproj * mdata ;
        %  data projected on subspace ortho to top 3 PCs


    if  ipart == 311 ;    %  Make PCA plot

      titlecellstr = {{'PCA View, Ortho to Tree 1' 'Neil Hayes Lung data' 'Carc, Col, Norm, SC only' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

      labelcellstr = {{'Ortho PC 1' 'Ortho PC 2' 'Ortho PC 3' 'Ortho PC 4'}} ;

      savestr = ['NeilLung1ip' num2str(ipart) 'OrthoTree1PCASP'] ;

      paramstruct = struct('icolor',icolor, ...
                           'markerstr',markerstr, ...
                           'isubpopkde',1, ...
                           'npcadiradd',4, ...
                           'legendcellstr',legendcellstr, ...
                           'mlegendcolor',mlegendcolor, ...
                           'titlecellstr',titlecellstr, ...
                           'labelcellstr',labelcellstr, ...
                           'savestr',savestr, ...
                           'iscreenwrite',1) ;

      scatplotSM(mdatasub,[],paramstruct) ;


    elseif  ipart == 312 ;    %  Proj on DWD Indiv. vs all plot

      mdir = DWD1SM(mdatasub(:,CarcinoidFlag),mdatasub(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
      mdir = [mdir DWD1SM(mdatasub(:,SmallCellFlag),mdatasub(:,CarcinoidFlag|NormalFlag|ColonFlag))] ;
      mdir = [mdir DWD1SM(mdatasub(:,NormalFlag),mdatasub(:,CarcinoidFlag|SmallCellFlag|ColonFlag))] ;
      mdir = [mdir DWD1SM(mdatasub(:,ColonFlag),mdatasub(:,NormalFlag|CarcinoidFlag|SmallCellFlag))] ;

      titlecellstr = {{'DWD View, Ortho to Tree 1' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

      labelcellstr = {{'Ortho DWD Carc-all ' 'Ortho DWD SC-all' 'Ortho DWD Norm-all' 'Ortho DWD Col-all'}} ;

      savestr = ['NeilLung1ip' num2str(ipart) 'OrthoTree1DWDIndiv-all'] ;

      paramstruct = struct('icolor',icolor, ...
                           'markerstr',markerstr, ...
                           'isubpopkde',1, ...
                           'legendcellstr',legendcellstr, ...
                           'mlegendcolor',mlegendcolor, ...
                           'titlecellstr',titlecellstr, ...
                           'labelcellstr',labelcellstr, ...
                           'savestr',savestr, ...
                           'iscreenwrite',1) ;

      scatplotSM(mdatasub,mdir,paramstruct) ;


    end ;


  elseif  ipart == 321  | ...
          ipart == 322  ;    %  Remove Tree 2 subspace

    mdir = DWD1SM(mdata(:,CarcinoidFlag|SmallCellFlag),mdata(:,NormalFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag),mdata(:,SmallCellFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,NormalFlag),mdata(:,ColonFlag))] ;

    mproj = mdir * pinv(mdir' * mdir) * mdir' ;
    mdatasub = mdata - mproj * mdata ;
        %  data projected on subspace ortho to top 3 PCs


    if  ipart == 321 ;    %  Make PCA plot

      titlecellstr = {{'PCA View, Ortho to Tree 2' 'Neil Hayes Lung data' 'Carc, Col, Norm, SC only' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

      labelcellstr = {{'Ortho PC 1' 'Ortho PC 2' 'Ortho PC 3' 'Ortho PC 4'}} ;

      savestr = ['NeilLung1ip' num2str(ipart) 'OrthoTree2PCASP'] ;

      paramstruct = struct('icolor',icolor, ...
                           'markerstr',markerstr, ...
                           'isubpopkde',1, ...
                           'npcadiradd',4, ...
                           'legendcellstr',legendcellstr, ...
                           'mlegendcolor',mlegendcolor, ...
                           'titlecellstr',titlecellstr, ...
                           'labelcellstr',labelcellstr, ...
                           'savestr',savestr, ...
                           'iscreenwrite',1) ;

      scatplotSM(mdatasub,[],paramstruct) ;


    elseif  ipart == 322 ;    %  Proj on DWD Indiv. vs all plot

      mdir = DWD1SM(mdatasub(:,CarcinoidFlag),mdatasub(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
      mdir = [mdir DWD1SM(mdatasub(:,SmallCellFlag),mdatasub(:,CarcinoidFlag|NormalFlag|ColonFlag))] ;
      mdir = [mdir DWD1SM(mdatasub(:,NormalFlag),mdatasub(:,CarcinoidFlag|SmallCellFlag|ColonFlag))] ;
      mdir = [mdir DWD1SM(mdatasub(:,ColonFlag),mdatasub(:,NormalFlag|CarcinoidFlag|SmallCellFlag))] ;

      titlecellstr = {{'DWD View, Ortho to Tree 2' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

      labelcellstr = {{'Ortho DWD Carc-all ' 'Ortho DWD SC-all' 'Ortho DWD Norm-all' 'Ortho DWD Col-all'}} ;

      savestr = ['NeilLung1ip' num2str(ipart) 'OrthoTree2DWDIndiv-all'] ;

      paramstruct = struct('icolor',icolor, ...
                           'markerstr',markerstr, ...
                           'isubpopkde',1, ...
                           'legendcellstr',legendcellstr, ...
                           'mlegendcolor',mlegendcolor, ...
                           'titlecellstr',titlecellstr, ...
                           'labelcellstr',labelcellstr, ...
                           'savestr',savestr, ...
                           'iscreenwrite',1) ;

      scatplotSM(mdatasub,mdir,paramstruct) ;


    end ;


  elseif  ipart == 331 | ...
          ipart == 332 ;    %  Remove Tree 3 subspace

    mdir = DWD1SM(mdata(:,SmallCellFlag),mdata(:,CarcinoidFlag|NormalFlag|ColonFlag)) ;
    mdir = [mdir DWD1SM(mdata(:,CarcinoidFlag),mdata(:,ColonFlag|NormalFlag))] ;
    mdir = [mdir DWD1SM(mdata(:,ColonFlag),mdata(:,NormalFlag))] ;

    mproj = mdir * pinv(mdir' * mdir) * mdir' ;
    mdatasub = mdata - mproj * mdata ;
        %  data projected on subspace ortho to top 3 PCs


    if  ipart == 331 ;    %  Make PCA plot

      titlecellstr = {{'PCA View, Ortho to Tree 3' 'Neil Hayes Lung data' 'Carc, Col, Norm, SC only' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

      labelcellstr = {{'Ortho PC 1' 'Ortho PC 2' 'Ortho PC 3' 'Ortho PC 4'}} ;

      savestr = ['NeilLung1ip' num2str(ipart) 'OrthoTree3PCASP'] ;

      paramstruct = struct('icolor',icolor, ...
                           'markerstr',markerstr, ...
                           'isubpopkde',1, ...
                           'npcadiradd',4, ...
                           'legendcellstr',legendcellstr, ...
                           'mlegendcolor',mlegendcolor, ...
                           'titlecellstr',titlecellstr, ...
                           'labelcellstr',labelcellstr, ...
                           'savestr',savestr, ...
                           'iscreenwrite',1) ;

      scatplotSM(mdatasub,[],paramstruct) ;


    elseif  ipart == 332 ;    %  Proj on DWD Indiv. vs all plot

      mdir = DWD1SM(mdatasub(:,CarcinoidFlag),mdatasub(:,SmallCellFlag|NormalFlag|ColonFlag)) ;
      mdir = [mdir DWD1SM(mdatasub(:,SmallCellFlag),mdatasub(:,CarcinoidFlag|NormalFlag|ColonFlag))] ;
      mdir = [mdir DWD1SM(mdatasub(:,NormalFlag),mdatasub(:,CarcinoidFlag|SmallCellFlag|ColonFlag))] ;
      mdir = [mdir DWD1SM(mdatasub(:,ColonFlag),mdatasub(:,NormalFlag|CarcinoidFlag|SmallCellFlag))] ;

      titlecellstr = {{'DWD View, Ortho to Tree 3' 'Neil Hayes Lung data' 'Carc Tri, Col Sq, Norm o, SC x' ['Genes filtered to ' num2str(d) ' exp''d']}} ;

      labelcellstr = {{'Ortho DWD Carc-all ' 'Ortho DWD SC-all' 'Ortho DWD Norm-all' 'Ortho DWD Col-all'}} ;

      savestr = ['NeilLung1ip' num2str(ipart) 'OrthoTree3DWDIndiv-all'] ;

      paramstruct = struct('icolor',icolor, ...
                           'markerstr',markerstr, ...
                           'isubpopkde',1, ...
                           'legendcellstr',legendcellstr, ...
                           'mlegendcolor',mlegendcolor, ...
                           'titlecellstr',titlecellstr, ...
                           'labelcellstr',labelcellstr, ...
                           'savestr',savestr, ...
                           'iscreenwrite',1) ;

      scatplotSM(mdatasub,mdir,paramstruct) ;


    end ;


  end ;    %  of inner ipart if-block
 

end ;    %  of outer ipart if-block




