Home > SuperSegger > viz > makeLineage.m

makeLineage

PURPOSE ^

makeLineage : creates a lineage tree for the cells with ID_.

SYNOPSIS ^

function data = makeLineage( clist, ID_, min_width )

DESCRIPTION ^

 makeLineage : creates a lineage tree for the cells with ID_.

 INPUT :
   clist : list of cells and characteristics
   ID_ : vector of cell ID to include. If empty all the cells are used
   min_width : min number of cells in a tree

 Copyright (C) 2016 Wiggins Lab
 Written by Paul Wiggins.
 University of Washington, 2016
 This file is part of SuperSegger.

 SuperSegger is free software: you can redistribute it and/or modify
 it under the terms of the GNU General Public License as published by
 the Free Software Foundation, either version 3 of the License, or
 (at your option) any later version.

 SuperSegger is distributed in the hope that it will be useful,
 but WITHOUT ANY WARRANTY; without even the implied warranty of
 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 GNU General Public License for more details.

 You should have received a copy of the GNU General Public License
 along with SuperSegger.  If not, see <http://www.gnu.org/licenses/>.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

SOURCE CODE ^

0001 function data = makeLineage( clist, ID_, min_width )
0002 % makeLineage : creates a lineage tree for the cells with ID_.
0003 %
0004 % INPUT :
0005 %   clist : list of cells and characteristics
0006 %   ID_ : vector of cell ID to include. If empty all the cells are used
0007 %   min_width : min number of cells in a tree
0008 %
0009 % Copyright (C) 2016 Wiggins Lab
0010 % Written by Paul Wiggins.
0011 % University of Washington, 2016
0012 % This file is part of SuperSegger.
0013 %
0014 % SuperSegger is free software: you can redistribute it and/or modify
0015 % it under the terms of the GNU General Public License as published by
0016 % the Free Software Foundation, either version 3 of the License, or
0017 % (at your option) any later version.
0018 %
0019 % SuperSegger is distributed in the hope that it will be useful,
0020 % but WITHOUT ANY WARRANTY; without even the implied warranty of
0021 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
0022 % GNU General Public License for more details.
0023 %
0024 % You should have received a copy of the GNU General Public License
0025 % along with SuperSegger.  If not, see <http://www.gnu.org/licenses/>.
0026 
0027 global daughter1_index;
0028 global daughter2_index;
0029 global cell_birth_index;
0030 global cell_div_index;
0031 global stat_0_index;
0032 global cell_error_index;
0033 global cell_age_index;
0034 global pole_age_index;
0035 
0036 figure(1)
0037 clf;
0038 figure(2);
0039 clf;
0040 figure(3);
0041 clf;
0042 
0043 
0044 clist_ = gate( clist );
0045 gated  = clist_.data(:,1); 
0046 
0047 
0048 
0049 if isempty(clist) 
0050     disp ('No clist found');
0051     return;
0052 end
0053 
0054 if ~exist( 'ID_', 'var' );
0055     ID_ = [];
0056 else
0057     ID_ = unique(ID_);
0058 end
0059 
0060 if ~exist( 'min_width' ) || isempty( min_width )
0061     min_width = 0;
0062 end
0063 
0064 
0065 clf;
0066 num = size( clist.data, 1 );
0067 
0068 
0069 starter = [];
0070 hh = [];
0071 hh2 = [];
0072 
0073 data.n0 = [];
0074 data.n  = {};
0075 data.t0 = [];
0076 data.t  = {};
0077 
0078 data.n1_max = [];
0079 data.n2_max = [];
0080 
0081 
0082 legend_text = {};
0083 count = 0;
0084 
0085 mother_index = grabClistIndex(clist, 'mother id');
0086 daughter1_index = grabClistIndex(clist, 'daughter1 id');
0087 daughter2_index = grabClistIndex(clist, 'daughter2 id');
0088 cell_div_index = grabClistIndex(clist, 'Cell Death Time');
0089 stat_0_index = grabClistIndex(clist,'stat0');
0090 cell_birth_index = grabClistIndex(clist,'Cell Birth Time');
0091 cell_error_index = grabClistIndex(clist,'Error Frame');
0092 cell_age_index = grabClistIndex(clist,'Cell Age' );
0093 pole_age_index = grabClistIndex(clist,'Old Pole Age' );
0094 
0095 if isempty(cell_error_index)
0096     disp ('no error frame found in clist')
0097 end
0098 
0099 if isempty( ID_ )
0100     for ii = 1: num
0101         
0102         
0103         if clist.data(ii,mother_index) == 0
0104             % Make new lineage
0105             ID = clist.data(ii,1);
0106             
0107             [width,list,time,death,stat0,starter,error,age,gen,pole] ...
0108                 = intGetWidth(ID, clist, starter,1);
0109             
0110             if (width >= min_width) || time(1) == 1
0111                 [hh, hh2, data] = intDoDraw( clist, list, time, death, ...
0112                     stat0, starter, hh, hh2, error, gated, data );
0113                 count = count + 1;
0114                 
0115                 data.width = width;
0116                 data.ID    = list;
0117                 data.birth = time;
0118                 data.death = death;
0119                 data.stat0 = stat0;
0120                 data.error = error;
0121                 data.age   = age;
0122                 data.gen   = gen;
0123                 data.pole  = pole;
0124                 
0125                 legend_text{count} = ['Cell ',num2str( ID )];
0126             end
0127             
0128         end
0129         
0130     end
0131 else
0132     for ii = 1:numel(ID_)
0133         
0134         ID = ID_(ii);
0135         [width,list,time,death,stat0,starter,error,age,gen,pole] = ...
0136             intGetWidth(ID, clist, starter, 1);
0137         
0138         [hh,hh2,data] = intDoDraw( clist, list, time, death, stat0, starter, ...
0139             hh, hh2, error, gated, data );
0140         
0141         data.width = width;
0142         data.ID    = list;
0143         data.birth = time;
0144         data.death = death;
0145         data.stat0 = stat0;
0146         data.error = error;
0147         data.age   = age;
0148         data.gen   = gen;
0149         data.pole  = pole;
0150         
0151         legend_text{ii} = ['Cell ',num2str( ID )];
0152     end
0153 end
0154 
0155 
0156 figure(1);
0157 set( gca, 'YDir', 'Reverse'  );
0158 
0159 ylabel( 'Time (Frames)' );
0160 xlabel( 'Cells' );
0161 title( 'Cell Lineage' );
0162 
0163 width1 = numel(starter);
0164 xlim( [-.1,1.1]*width1 );
0165 
0166 height1 = max(clist.data(:,5));
0167 ylim( [-.1,1.1]*height1 );
0168 
0169 
0170 figure(2);
0171 ylabel( 'Number of cells' );
0172 xlabel( 'time (frames)' );
0173 legend( hh, legend_text, 'Location' , 'NorthWest' );
0174 set( gca, 'YScale', 'log'  );
0175 title( 'Cummulative  Number of Cells' );
0176 
0177 try
0178 ylim_ = data.n1_max;
0179 ylim( [0.5,2*ylim_]);
0180 xlim( [-.1,1.1]*height1 );
0181 end
0182 
0183 figure(3);
0184 ylabel( 'Log_2 Number of cells' );
0185 xlabel( 'Time (frames)' );
0186 legend( hh, legend_text, 'Location' , 'NorthWest' );
0187 %set( gca, 'YScale', 'log'  );
0188 title( 'Number of Cells' );
0189 
0190 legend( hh2, legend_text, 'Location' , 'NorthWest' );
0191 
0192 try
0193 ylim_ = data.n2_max;
0194 ylim( ceil(log([0.5,2*ylim_(end)])/log(2)) );
0195 xlim( [-.1,1.1]*height1 );
0196 end
0197 
0198 set( gca, 'YGrid', 'on')
0199 
0200 end
0201 
0202 function [width,list,time,death,stat0,starter,error,age,gen,pole] = ...
0203     intGetWidth(ID, clist, starter, gen__)
0204 global cell_birth_index;
0205 global cell_div_index;
0206 global cell_age_index;
0207 global stat_0_index;
0208 global cell_error_index;
0209 global pole_age_index;
0210 
0211 end_time = max( clist.data(:,cell_div_index) );
0212 ind = find( clist.data(:,1)==ID );
0213 
0214 if isempty(ind) || (ID == 0) || isnan(ID)
0215     width = 0;
0216     list  = [];
0217     time  = [];
0218     death = [];
0219     stat0 = [];
0220     error = [];
0221     age   = [];
0222     gen   = [];
0223     pole  = [];
0224 else
0225     [ID1,ID2] = intGetDaughter( ID, clist );
0226     
0227     error_ = clist.data(ind,cell_error_index);
0228     death_ = clist.data(ind,cell_div_index);
0229     age_   = clist.data(ind,cell_age_index);
0230     gen_   = gen__;
0231     
0232     if isempty(error_)
0233         error_ = nan;
0234     end
0235 
0236     if isnan(ID1) || isnan(ID2) || (ID1==0) || (ID2==0) || ~isnan(error_)
0237         starter = [starter,ID];
0238         
0239         if death_ ~= end_time && isnan(error_);
0240             error_ = death_;
0241         end
0242             
0243     end
0244     
0245     [w1,l1,t1,d1,s1,starter,e1,a1,g1,p1] = intGetWidth( ID1, clist, starter,gen__ + 1 );
0246     [w2,l2,t2,d2,s2,starter,e2,a2,g2,p2] = intGetWidth( ID2, clist, starter,gen__ + 1  );
0247     
0248     stat0_ = clist.data(ind,stat_0_index);
0249     time_  = clist.data(ind,cell_birth_index);
0250     pole_  = clist.data(ind,pole_age_index);
0251     %death_ = death_ + time_;
0252     
0253     width =  1 + w1 + w2;
0254     list  = [ID,l1,l2];
0255     time  = [time_,t1,t2];
0256     death = [death_,d1,d2];
0257     stat0 = [stat0_,s1,s2];
0258     error = [error_,e1,e2];
0259     age   = [age_,a1,a2];
0260     gen   = [gen_,g1,g2];
0261     pole  = [pole_,p1,p2];
0262 end
0263 end
0264 
0265 
0266 function [hh,hh2,data] = intDoDraw( clist, list, time, death, stat0, starter,...
0267     hh, hh2, error, gated, data  )
0268 
0269 %% Show the lineages
0270 figure(1);
0271 
0272 num = numel( list );
0273 for ii = 1:num
0274     
0275     if stat0(ii) == 2
0276         cc = 'b';
0277     elseif stat0(ii) == 1
0278         cc = 'g';
0279     else
0280         cc = 'g';
0281     end
0282     
0283     ID = list(ii);
0284     
0285     pos = intGetPos( ID, clist, starter );
0286     
0287     if ismember( ID, gated )
0288         style = '-';
0289     else
0290         style = ':';    
0291     end
0292     
0293     plot( pos+[0,0], [time(ii)-1,death(ii)], style, 'Color', cc, 'LineWidth', 1);
0294     hold on;
0295     
0296     if ~isnan(error(ii))
0297             plot( pos, error(ii), '.', 'Color', 'r', 'MarkerSize', 20 );
0298     end
0299     
0300     
0301     [ID1,ID2] = intGetDaughter( ID, clist );
0302     
0303     if isnan(ID1) || (ID1==0)
0304         pos1 = pos;
0305     else
0306         pos1 = intGetPos( ID1, clist, starter );
0307     end
0308     
0309     if isnan(ID1) || (ID1==0)
0310         pos2 = pos;
0311     else
0312         pos2 = intGetPos( ID2, clist, starter );
0313     end
0314     
0315     plot( [min(pos1),max(pos2)], [0,0] + death(ii), style, 'Color', 'b', ...
0316         'LineWidth', .5 );
0317     
0318     
0319     h = text( pos(1), time(ii)-1, [' ',num2str( ID )], ...
0320         'VerticalAlignment', 'middle', 'HorizontalAlignment', 'left' );
0321     
0322     set( h, 'rotation', 90 );
0323     
0324 end
0325 
0326 %% Show the cell number;
0327 figure(2);
0328 
0329 tt = sort( time );
0330 
0331 nn = 1:numel(tt);
0332 
0333 h_ = stairs( tt, nn );
0334 
0335 data.n1_max = max( [data.n1_max, max(nn) ] );
0336 
0337 
0338 hh = [hh, h_];
0339 
0340 hold on;
0341 
0342 cc = get(h_,'Color' );
0343 
0344 first_death_time = min(error);
0345 
0346 ind = sum(first_death_time>tt);
0347 if ~ind
0348     ind = 1;
0349 end
0350 
0351 plot( first_death_time, nn(ind), '.', ...
0352     'MarkerSize', 20, 'Color', cc );
0353 
0354 data.t = { data.t{:}, tt };
0355 data.n = { data.n{:}, nn };
0356 
0357 data.t0 = [ data.t0, first_death_time];
0358 data.n0 = [ data.n0, nn(ind) ];
0359 
0360 
0361 
0362 figure(3);
0363 tt   = sort( time );
0364 dtt  = sort( death+1 );
0365 ntt  = ones(size(tt));
0366 ndtt = -ones(size(dtt));
0367 
0368 ntt = [ntt,ndtt];
0369 tt  = [tt,dtt]; 
0370 [ttt,ord] = sort( tt );
0371 ntt = ntt(ord);
0372 nn = cumsum( ntt );
0373 
0374 
0375 [ttt,ord] = unique( ttt );
0376 nn = nn(ord);
0377 
0378 h_ = stairs( ttt, log(nn)/log(2), '-', 'Color', cc );
0379 hold on;
0380 
0381 data.n2_max = max( [data.n2_max, max(nn) ] );
0382 
0383 ind = sum(first_death_time>ttt);
0384 if ~ind
0385     ind = 1;
0386 end
0387 
0388 plot( first_death_time, log(nn(ind))/log(2), '.', ...
0389     'MarkerSize', 20, 'Color', cc );
0390 
0391 hh2 = [hh2, h_];
0392 
0393 
0394 end
0395 
0396 
0397 function  pos = intGetPos( ID, clist, starter )
0398 % intGetPos : gets the position of where the cell should be drawn.
0399 
0400 [ID1,ID2] = intGetDaughter( ID, clist );
0401 
0402 flag1 = isnan(ID1) + (ID1==0);
0403 flag2 = isnan(ID2) + (ID2==0);
0404 
0405 pos = find( starter== ID );
0406 
0407 if (~flag1) && (flag2)
0408     pos = (pos + intGetPos( ID1, clist, starter ))/2;
0409 elseif (flag1) && (~flag2)
0410     pos = (pos + intGetPos( ID2, clist, starter ))/2;
0411 elseif (~flag1) && (~flag2)
0412     pos = (intGetPos( ID1, clist, starter ) + ...
0413         intGetPos( ID2, clist, starter ) )/2;
0414 end
0415 
0416 end
0417 
0418 
0419 function [ID1,ID2] = intGetDaughter( ID, clist )
0420 global daughter1_index
0421 global daughter2_index
0422 ind = find( clist.data(:,1)==ID );
0423 
0424 if isempty(ind)
0425     ID1 = nan;
0426     ID2 = nan;
0427 else
0428     ID1 = clist.data(ind,daughter1_index);
0429     ID2 = clist.data(ind,daughter2_index);
0430 end
0431 
0432 end
0433 
0434 
0435 
0436 
0437

Generated on Thu 19-Jan-2017 13:55:21 by m2html © 2005