From 7134418edb2caccaca0ae6901abf616e31ca1401 Mon Sep 17 00:00:00 2001
From: Mark Winter <mwinter@drexel.edu>
Date: Wed, 20 Jan 2021 15:20:56 +0100
Subject: [PATCH] Fixes for vesicle stats pull and only use ~6 sec/frame in
tracking statistics
---
+Stats/collate_tracks.m | 34 +++++---
+Stats/export_intensities.m | 19 +++++
+Stats/load_ectostats.m | 145 ++++++++++++++++++++++++++--------
+Stats/plot_ecto_stats.m | 97 +++++++++++++++--------
+Stats/plot_vel_histograms.m | 4 +-
+Stats/region_avg_intensity.m | 128 ++++++++++++++++++++++++++++++
6 files changed, 344 insertions(+), 83 deletions(-)
create mode 100644 +Stats/export_intensities.m
create mode 100644 +Stats/region_avg_intensity.m
diff --git a/+Stats/collate_tracks.m b/+Stats/collate_tracks.m
index eba30b8..3407f52 100644
--- a/+Stats/collate_tracks.m
+++ b/+Stats/collate_tracks.m
@@ -25,6 +25,12 @@ function track_data = collate_tracks(ROOT, filename)
%% Estimate parallelogram mask for first-frame image
imagePath = CONSTANTS.imageData.imageDir;
+ chkimd = MicroscopeData.ReadMetadataFile(fullfile(imagePath,filesep));
+ if ( isempty(chkimd) )
+ warning(['Unable to load image data for: ' datasetName]);
+ return;
+ end
+
im = MicroscopeData.Reader(imagePath, 'timeRange',[1,1]);
parallelMask = Segment.llsm_mask(im);
@@ -52,9 +58,10 @@ function track_data = collate_tracks(ROOT, filename)
return;
end
- %% Skip if movie is >= 7fps
+ %% Skip if movie is >= 7fps or < 5
medTsmp = median(diff(tsmp));
- if ( medTsmp >= 7.0 )
+ if ( medTsmp >= 7.0 || medTsmp < 5 )
+% if ( medTsmp >= 7.0 )
return;
end
@@ -86,15 +93,17 @@ function track_data = collate_tracks(ROOT, filename)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Read per-frame vesicle info
- bEcto = false(0,1);
- bSkel = false(0,1);
- bEndo = false(0,1);
+ cid_max = max([trackInfo{:}]);
+
+ bEcto = false(cid_max,1);
+ bSkel = false(cid_max,1);
+ bEndo = false(cid_max,1);
- tids = zeros(0,1);
- coms = zeros(0,3);
- ctimes = zeros(0,1);
- csizes = zeros(0,1);
- cframes = zeros(0,1);
+ tids = zeros(cid_max,1);
+ coms = zeros(cid_max,3);
+ ctimes = zeros(cid_max,1);
+ csizes = zeros(cid_max,1);
+ cframes = zeros(cid_max,1);
imSize = [size(im,1), size(im,2), size(im,3)];
for t=1:keeplen
@@ -131,11 +140,10 @@ function track_data = collate_tracks(ROOT, filename)
tlens = cellfun(@(x)(length(x)), trackInfo);
bKeep = (tlens >= minTrackLen);
- lastId = length(csizes);
trackInfo = trackInfo(bKeep);
for i=1:length(trackInfo)
- bValid = (trackInfo{i} <= lastId);
- trackInfo{i} = trackInfo{i}(bValid);
+ bValidHulls = ((cframes(trackInfo{i}) > 0 ) & cframes(trackInfo{i}) <= keeplen);
+ trackInfo{i} = trackInfo{i}(bValidHulls);
end
bValidInfo = cellfun(@(x)(length(x) >= 5), trackInfo);
diff --git a/+Stats/export_intensities.m b/+Stats/export_intensities.m
new file mode 100644
index 0000000..590a8e1
--- /dev/null
+++ b/+Stats/export_intensities.m
@@ -0,0 +1,19 @@
+function export_intensities()
+ s = load('ds_intensities.mat');
+
+ fid = fopen('ds_intensities.csv','wt');
+ fprintf(fid, 'DatasetName,cal_ff_skel,cal_ff_ecto,cal_skel,cal_ecto,fm_ff_skel,fm_ff_ecto,fm_skel,fm_ecto\n');
+ for i=1:length(s.slist)
+ fprintf(fid, '%s,%f,%f,%f,%f,%f,%f,%f,%f\n',...
+ s.slist(i).datasetName,...
+ s.slist(i).means(2).ff_skel_mean,...
+ s.slist(i).means(2).ff_ecto_mean,...
+ s.slist(i).means(2).skel_mean,...
+ s.slist(i).means(2).ecto_mean,...
+ s.slist(i).means(1).ff_skel_mean,...
+ s.slist(i).means(1).ff_ecto_mean,...
+ s.slist(i).means(1).skel_mean,...
+ s.slist(i).means(1).ecto_mean);
+ end
+ fclose(fid);
+end
diff --git a/+Stats/load_ectostats.m b/+Stats/load_ectostats.m
index 9b22f6d..09fd7f2 100644
--- a/+Stats/load_ectostats.m
+++ b/+Stats/load_ectostats.m
@@ -123,14 +123,17 @@ function [size_stats, track_stats] = load_ectostats(ROOT, filename)
%% Read per-frame vesicle info
- bEcto = false(0,1);
- bSkel = false(0,1);
- bEndo = false(0,1);
+ cid_max = max([trackInfo{:}]);
- coms = zeros(0,3);
- ctimes = zeros(0,1);
- csizes = zeros(0,1);
- cframes = zeros(0,1);
+ bEcto = false(cid_max,1);
+ bSkel = false(cid_max,1);
+ bEndo = false(cid_max,1);
+
+ tids = zeros(cid_max,1);
+ coms = zeros(cid_max,3);
+ ctimes = zeros(cid_max,1);
+ csizes = zeros(cid_max,1);
+ cframes = zeros(cid_max,1);
imSize = [size(im,1), size(im,2), size(im,3)];
for t=1:keeplen
@@ -175,9 +178,10 @@ function [size_stats, track_stats] = load_ectostats(ROOT, filename)
size_stats.endo.com = coms(bEndo,:);
%% Tracking stats
- %% Skip if movie is >= 7fps
+ %% Skip if movie is >= 7fps or < 5
medTsmp = median(diff(tsmp));
- if ( medTsmp >= 7.0 )
+ if ( medTsmp >= 7.0 || medTsmp < 5 )
+% if ( medTsmp >= 7.0 )
return;
end
@@ -187,11 +191,10 @@ function [size_stats, track_stats] = load_ectostats(ROOT, filename)
tlens = cellfun(@(x)(length(x)), trackInfo);
bKeep = (tlens >= minTrackLen);
- lastId = length(csizes);
trackInfo = trackInfo(bKeep);
for i=1:length(trackInfo)
- bValid = (trackInfo{i} <= lastId);
- trackInfo{i} = trackInfo{i}(bValid);
+ bValidHulls = ((cframes(trackInfo{i}) > 0 ) & cframes(trackInfo{i}) <= keeplen);
+ trackInfo{i} = trackInfo{i}(bValidHulls);
end
bValidInfo = cellfun(@(x)(length(x) >= 5), trackInfo);
@@ -245,12 +248,109 @@ function [size_stats, track_stats] = load_ectostats(ROOT, filename)
%% Compute per-track info
cids = trackInfo{i};
+ tframe = cframes(cids);
ttime = ctimes(cids);
tsize = csizes(cids);
tpos = coms(cids,:);
tdt = diff(ttime);
+ tdfr = diff(tframe);
+
+ tdx = diff(tpos,1,1);
+ tsd = sum(tdx.^2, 2);
+ tssd = cumsum(tsd);
+
+ n_msd = floor(length(tdfr)/2);
+ mx_msd = sum(tdfr);
+% ndf = (1:mx_msd).';
+% cdf = zeros(mx_msd,1);
+% csd = zeros(mx_msd,1);
+% for k=1:n_msd
+% idxk = (1:ceil(length(tdt)/k)).';
+% idxsum = reshape(repmat(idxk,1,k).',[],1);
+% idxsum = idxsum(1:length(tdt));
+%
+% chk_df = accumarray(idxsum,tdfr);
+% chk_sd = accumarray(idxsum,tsd);
+%
+% cdf = cdf + accumarray(chk_df,ones(length(chk_sd),1), [mx_msd,1]);
+% csd = csd + accumarray(chk_df,chk_sd, [mx_msd,1]);
+% end
+%
+% bvdf = (cdf > 0);
+% ndf = ndf(bvdf);
+% cdf = cdf(bvdf);
+% csd = csd(bvdf);
+%
+% if ( length(csd) >= 6 )
+% ndf = ndf(1:ceil(end/2));
+% cdf = cdf(1:ceil(end/2));
+% csd = csd(1:ceil(end/2));
+% end
+%
+% % Compute Diffusion coefficient (with bias)
+% tmsd = (csd ./ cdf);
+% tmdt = median(tdt ./ tdfr);
+%
+% x1_rd = ndf * tmdt;
+% x2_rd = ones(length(ndf),1);
+%
+% [b_rd,~,~,~,rstat_rd] = regress(tmsd, [x1_rd x2_rd]);
+% D_rd = b_rd(1)/(2*msd_gamma);
+
+ % Use single delta-t diffusoin computation if (D_rd < 0) this
+ % generally indicates a tracking error on a short track
+% if ( D_rd < 0 )
+ x1 = cumsum(tdt);
+ x2 = ones(size(tdt,1),1);
+ [b,~,~,~,rstat] = regress(tssd, [x1 x2]);
+
+ track_msd_D(i) = b(1)/(2*msd_gamma);
+ track_msd_sig(i) = b(2)/(2*msd_gamma);
+ track_msd_R(i) = sqrt(rstat(1));
+
+% msde = mean(tssd ./ x1);
+% rsqe = 1 - mean((tssd - msde*x1).^2) / var(tssd);
+%
+% track_msd_D(i) = msde/(2*msd_gamma);
+% track_msd_R(i) = sqrt(rsqe);
+% else
+% track_msd_D(i) = D_rd;
+% track_msd_sig(i) = b_rd(2)/(2*msd_gamma);
+% track_msd_R(i) = sqrt(rstat_rd(1));
+% end
+
+% if ( length(tssd) > 30 && rstat(1) > 0.9 )
+% xl = [0,max(x1)+1];
+% yl = [0,max(tssd)+(max(tssd)/length(tssd))];
+% xx = [xl(1) xl(2)+1];
+% yy = b(1)*xx + b(2);
+%
+% m2 = mean(tssd ./ x1);
+% r2 = 1 - mean((tssd - m2*x1).^2) / var(tssd);
+%
+% yy = m2*xx;
+%
+% cmap = lines(7);
+% figure;hold on
+% plot(x1,tssd,'.', 'Color',cmap(7,:), 'MarkerSize',12);
+% plot(xx,yy, '-', 'Color',cmap(1,:), 'LineWidth',2);
+% xlim(xl);
+% ylim(yl);
+% end
+
+
+% tmsd = cumsum(sum(tdx.^2, 2));
+
+ tmot = sqrt(sum(tdx.^2, 2));
+ tvel = tmot ./ tdt;
+
+ nrmtdx = tdx ./ sqrt(tdt);
+ brntvel = sqrt(sum(nrmtdx.^2, 2));
+% brntvel = tmot ./ sqrt(dt);
+
+ % Handle spicule stats if has defined spicule
if ( ~isempty(sp_struct) )
T = sp_struct.root .* pxSize;
tds = spicule_distances(cs,Rs,rs,als, T,tpos);
@@ -264,29 +364,10 @@ function [size_stats, track_stats] = load_ectostats(ROOT, filename)
track_sp_ivel(i) = mean(tds_delta ./ tdt);
end
- tdx = diff(tpos,1,1);
- tmsd = cumsum(sum(tdx.^2, 2));
-
- tmot = sqrt(sum(tdx.^2, 2));
- tvel = tmot ./ tdt;
-
- nrmtdx = tdx ./ sqrt(tdt);
- brntvel = sqrt(sum(nrmtdx.^2, 2));
-% brntvel = tmot ./ sqrt(dt);
-
- % Compute Diffusion coefficient (with bias)
- x1 = cumsum(tdt);
- x2 = ones(size(tdt,1),1);
- [b,~,~,~,rstat] = regress(tmsd, [x1 x2]);
-
- track_msd_D(i) = b(1)/(2*msd_gamma);
- track_msd_sig(i) = b(2)/(2*msd_gamma);
- track_msd_R(i) = sqrt(rstat(1));
-
% Compute non-ideal diffusion coeffient
x1 = log(cumsum(tdt));
x2 = ones(size(tdt,1),1);
- [b,~,~,~,rstat] = regress(log(tmsd), [x1 x2]);
+ [b,~,~,~,rstat] = regress(log(tssd), [x1 x2]);
track_nimsd_alpha(i) = b(1);
track_nimsd_D(i) = exp(b(2))/(2*msd_gamma);
diff --git a/+Stats/plot_ecto_stats.m b/+Stats/plot_ecto_stats.m
index de36c14..f67325b 100644
--- a/+Stats/plot_ecto_stats.m
+++ b/+Stats/plot_ecto_stats.m
@@ -242,37 +242,37 @@ bBB94_skel = cellfun(@(x)(strcmpi(x,'BB94')), exp_skel_grp);
exp_skel_grp = exp_skel_grp(~bBB94_skel);
skel_wdiridx = skel_wdiridx(~bBB94_skel);
-figure;
-boxplot(skel_wdiridx, exp_skel_grp, 'notch','on');
-xtickangle(30);
-zoom('reset');
-% ylim([0.3,0.6]);
-title('Directionality Index - Mesoderm');
-
-%
-exp_endo_grp_cell = arrayfun(@(x,y)(repmat(x,size(y.mean_sizes,1),1)), exp_types_track, endo_track_info, 'UniformOutput',false);
-exp_endo_grp = vertcat(exp_endo_grp_cell{:});
-
-endo_wdir_cell = vertcat(endo_track_info.avg_wdir_index);
-endo_wdiridx = cellfun(@(x)(x(didx_idx)), endo_wdir_cell);
-
-bBB94_ecto = cellfun(@(x)(strcmpi(x,'BB94')), exp_endo_grp);
-exp_endo_grp = exp_endo_grp(~bBB94_ecto);
-endo_wdiridx = endo_wdiridx(~bBB94_ecto);
-
-figure;
-boxplot(endo_wdiridx, exp_endo_grp, 'notch','on');
-xtickangle(30);
-zoom('reset');
+% figure;
+% boxplot(skel_wdiridx, exp_skel_grp, 'notch','on');
+% xtickangle(30);
+% zoom('reset');
% ylim([0.3,0.6]);
-title('Directionality Index - Endoderm');
+% title('Directionality Index - Mesoderm');
+%
+% %
+% exp_endo_grp_cell = arrayfun(@(x,y)(repmat(x,size(y.mean_sizes,1),1)), exp_types_track, endo_track_info, 'UniformOutput',false);
+% exp_endo_grp = vertcat(exp_endo_grp_cell{:});
+%
+% endo_wdir_cell = vertcat(endo_track_info.avg_wdir_index);
+% endo_wdiridx = cellfun(@(x)(x(didx_idx)), endo_wdir_cell);
+%
+% bBB94_ecto = cellfun(@(x)(strcmpi(x,'BB94')), exp_endo_grp);
+% exp_endo_grp = exp_endo_grp(~bBB94_ecto);
+% endo_wdiridx = endo_wdiridx(~bBB94_ecto);
+%
+% figure;
+% boxplot(endo_wdiridx, exp_endo_grp, 'notch','on');
+% xtickangle(30);
+% zoom('reset');
+% % ylim([0.3,0.6]);
+% title('Directionality Index - Endoderm');
%%
exp_ecto_grp_cell = arrayfun(@(x,y)(repmat(x,size(y.mean_sizes,1),1)), exp_types_track, ecto_track_info, 'UniformOutput',false);
exp_ecto_grp = vertcat(exp_ecto_grp_cell{:});
ecto_msdD = vertcat(ecto_track_info.msd_D);
-ecto_msdR = vertcat(ecto_track_info.msd_R);
+ecto_msdR = real(vertcat(ecto_track_info.msd_R));
ecto_msdsig = vertcat(ecto_track_info.msd_sig);
bBB94_ecto = cellfun(@(x)(strcmpi(x,'BB94')), exp_ecto_grp);
@@ -394,6 +394,7 @@ didx_idx = 8;
all_mean_sizes = [vertcat(ecto_track_info.mean_sizes); vertcat(skel_track_info.mean_sizes)];
all_mean_ivel = [vertcat(ecto_track_info.mean_inst_vel); vertcat(skel_track_info.mean_inst_vel)];
all_diff_coeff = [vertcat(ecto_track_info.msd_D); vertcat(skel_track_info.msd_D)];
+all_diff_R = [vertcat(ecto_track_info.msd_R); vertcat(skel_track_info.msd_R)];
all_mean_didx = [vertcat(ecto_track_info.avg_wdir_index); vertcat(skel_track_info.avg_wdir_index)];
all_mean_didx = cellfun(@(x)(x(didx_idx)), all_mean_didx);
@@ -410,7 +411,7 @@ all_direct_len = [vertcat(ecto_track_info.direct_len); vertcat(skel_track_info.d
all_path_ratio = all_direct_len ./ all_path_len;
%%
-[~,ctr] = hist(all_mean_ivel(bSkelDMSO), ctr);
+[~,ctr] = hist(all_mean_ivel(bSkelDMSO), 100);
figure;
hist(all_mean_ivel(bSkelDMSO), ctr);
xlim([0,0.3]);
@@ -466,11 +467,17 @@ ylim([0,0.2]);
%
plot(all_mean_sizes(bSkelAxt),all_mean_ivel(bSkelAxt), '.', 'Color',cmap(1,:));
-[b,~,r,~,rst] = regress(all_mean_ivel(bSkelAxt|bSkelDMSO), [ones(size(all_mean_ivel(bSkelAxt|bSkelDMSO))), all_mean_sizes(bSkelAxt|bSkelDMSO)]);
+[b,~,r,~,rstDMSO] = regress(all_mean_ivel(bSkelDMSO), [ones(size(all_mean_ivel(bSkelDMSO))), all_mean_sizes(bSkelDMSO)]);
+xx = 0:0.1:1.3;
+yy = b(1) + xx*b(2);
+plot(xx,yy, '-', 'Color',cmap(2,:), 'LineWidth',2);
+
+[b,~,r,~,rstAxt] = regress(all_mean_ivel(bSkelAxt), [ones(size(all_mean_ivel(bSkelAxt))), all_mean_sizes(bSkelAxt)]);
xx = 0:0.1:1.3;
yy = b(1) + xx*b(2);
-plot(xx,yy, '-k', 'LineWidth',2);
-legend('Mesoderm: Control','Mesoderm: VEGFR Inhib.');
+plot(xx,yy, '-', 'Color',cmap(1,:),'LineWidth',2);
+
+legend('Skeletogenic Cells: Control','Skeletogenic Cells: VEGFR Inhib.');
saveas(gcf,'C:\Users\mwinter\Desktop\figs\sup2_meso_speed_v_size.png');
% figure;plot(all_mean_sizes(bEctoDMSO|bSkelDMSO), r, '.r');
@@ -483,10 +490,16 @@ ylim([0,0.2]);
%
plot(all_mean_sizes(bEctoAxt),all_mean_ivel(bEctoAxt), '.', 'Color',cmap(1,:));
-[b,~,r,~,rst] = regress(all_mean_ivel(bEctoAxt|bEctoDMSO), [ones(size(all_mean_ivel(bEctoAxt|bEctoDMSO))), all_mean_sizes(bEctoAxt|bEctoDMSO)]);
+[b,~,r,~,rstDMSO] = regress(all_mean_ivel(bEctoDMSO), [ones(size(all_mean_ivel(bEctoDMSO))), all_mean_sizes(bEctoDMSO)]);
+xx = 0:0.1:1.3;
+yy = b(1) + xx*b(2);
+plot(xx,yy, '-', 'Color',cmap(2,:), 'LineWidth',2);
+
+[b,~,r,~,rstAxt] = regress(all_mean_ivel(bEctoAxt), [ones(size(all_mean_ivel(bEctoAxt))), all_mean_sizes(bEctoAxt)]);
xx = 0:0.1:1.3;
yy = b(1) + xx*b(2);
-plot(xx,yy, '-k', 'LineWidth',2);
+plot(xx,yy, '-', 'Color',cmap(1,:),'LineWidth',2);
+
legend('Ectoderm: Control','Ectoderm: VEGFR Inhib.');
saveas(gcf,'C:\Users\mwinter\Desktop\figs\sup2_ecto_speed_v_size.png');
@@ -500,11 +513,18 @@ ylim([0,0.05]);
%
plot(all_mean_sizes(bSkelAxt),all_diff_coeff(bSkelAxt), '.', 'Color',cmap(1,:));
-[b,~,r,~,rst] = regress(all_diff_coeff(bSkelAxt|bSkelDMSO), [ones(size(all_diff_coeff(bSkelAxt|bSkelDMSO))), all_mean_sizes(bSkelAxt|bSkelDMSO)]);
+% [b,~,r,~,rst] = regress(all_diff_coeff(bSkelDMSO|bSkelAxt), [ones(size(all_diff_coeff(bSkelDMSO|bSkelAxt))), all_mean_sizes(bSkelDMSO|bSkelAxt)]);
+[b,~,r,~,rstDMSO] = regress(all_diff_coeff(bSkelDMSO), [ones(size(all_diff_coeff(bSkelDMSO))), all_mean_sizes(bSkelDMSO)]);
+xx = 0:0.1:1.3;
+yy = b(1) + xx*b(2);
+plot(xx,yy, '-', 'Color',cmap(2,:), 'LineWidth',2);
+
+[b,~,r,~,rstAxt] = regress(all_diff_coeff(bSkelAxt), [ones(size(all_diff_coeff(bSkelAxt))), all_mean_sizes(bSkelAxt)]);
xx = 0:0.1:1.3;
yy = b(1) + xx*b(2);
-plot(xx,yy, '-k', 'LineWidth',2);
-legend('Mesoderm: Control','Mesoderm: VEGFR Inhib.');
+plot(xx,yy, '-', 'Color',cmap(1,:),'LineWidth',2);
+
+legend('Skeletogenic Cells: Control','Skeletogenic Cells: VEGFR Inhib.');
saveas(gcf,'C:\Users\mwinter\Desktop\figs\sup3_meso_diff_v_size.png');
% figure;plot(all_mean_sizes(bEctoDMSO|bSkelDMSO), r, '.r');
@@ -517,10 +537,17 @@ ylim([0,0.05]);
%
plot(all_mean_sizes(bEctoAxt),all_diff_coeff(bEctoAxt), '.', 'Color',cmap(1,:));
-[b,~,r,~,rst] = regress(all_diff_coeff(bEctoAxt|bEctoDMSO), [ones(size(all_diff_coeff(bEctoAxt|bEctoDMSO))), all_mean_sizes(bEctoAxt|bEctoDMSO)]);
+% [b,~,r,~,rst] = regress(all_diff_coeff(bEctoAxt|bEctoDMSO), [ones(size(all_diff_coeff(bEctoAxt|bEctoDMSO))), all_mean_sizes(bEctoAxt|bEctoDMSO)]);
+[b,~,r,~,rstDMSO] = regress(all_diff_coeff(bEctoDMSO), [ones(size(all_diff_coeff(bEctoDMSO))), all_mean_sizes(bEctoDMSO)]);
xx = 0:0.1:1.3;
yy = b(1) + xx*b(2);
-plot(xx,yy, '-k', 'LineWidth',2);
+plot(xx,yy, '-', 'Color',cmap(2,:), 'LineWidth',2);
+
+[b,~,r,~,rstAxt] = regress(all_diff_coeff(bEctoAxt), [ones(size(all_diff_coeff(bEctoAxt))), all_mean_sizes(bEctoAxt)]);
+xx = 0:0.1:1.3;
+yy = b(1) + xx*b(2);
+plot(xx,yy, '-', 'Color',cmap(1,:),'LineWidth',2);
+
legend('Ectoderm: Control','Ectoderm: VEGFR Inhib.');
saveas(gcf,'C:\Users\mwinter\Desktop\figs\sup3_ecto_diff_v_size.png');
diff --git a/+Stats/plot_vel_histograms.m b/+Stats/plot_vel_histograms.m
index 5d0dcc0..55e1264 100644
--- a/+Stats/plot_vel_histograms.m
+++ b/+Stats/plot_vel_histograms.m
@@ -69,7 +69,6 @@ xlabel('Vesicle Speed [\mum/s]')
subplot(2,2,3)
hist(log10(axt_skel_ivel), nbins);
-hist(log10(dmso_ecto_ivel), nbins);
oldL = get(gca, 'XTickLabel');
newL = cellfun(@(x)(['10^{' x '}']), oldL,'UniformOutput',false);
set(gca, 'XTickLabel',newL);
@@ -78,12 +77,11 @@ xlabel('Vesicle Speed [\mum/s]')
subplot(2,2,4)
hist(log10(axt_ecto_ivel), nbins);
-hist(log10(dmso_ecto_ivel), nbins);
oldL = get(gca, 'XTickLabel');
newL = cellfun(@(x)(['10^{' x '}']), oldL,'UniformOutput',false);
set(gca, 'XTickLabel',newL);
title('Instantaneous Speed - VEGFR Inhib. - Ectoderm')
xlabel('Vesicle Speed [\mum/s]')
-saveas(gcf, 'C:\Users\mwinter\Desktop\figs\sup_figure2.png');
+saveas(gcf, 'C:\Users\mwinter\Desktop\figs\sup_inst_speed_hist_fix.png');
diff --git a/+Stats/region_avg_intensity.m b/+Stats/region_avg_intensity.m
new file mode 100644
index 0000000..922d488
--- /dev/null
+++ b/+Stats/region_avg_intensity.m
@@ -0,0 +1,128 @@
+function region_avg_intensity()
+ dslist = {'23-06-2016_TimeLapse1_DMSO_20hrs_Calcein_FM464'
+ '22-06-2016_TimeLapse1_DMSO_24hrs_Calcein_FM464'
+ '22-06-2016_TimeLapse3_DMSO_20hrs_Calcein_FM464'
+ '22-06-2016_TimeLapse2_DMSO_24hrs_Calcein_FM464'
+ '01-07-2016_TimeLapse4_DMSO_21hrs_Calcein_FM464'
+ '01-07-2016_TimeLapse3_DMSO_21hrs_Calcein_FM464'
+ '01-07-2016_TimeLapse5_DMSO_23hrs_Calcein_FM464'
+ '28-06-2016_TimeLapse1_Pos1_DMSO_18hrs_Calcein_FM464'
+ '28-06-2016_TimeLapse1_Pos2_DMSO_18hrs_Calcein_FM464'
+ '24-06-2016_TimeLapse1_Axtinib_150_19hrs_Calcein_FM464'
+ '24-06-2016_TimeLapse2_Axtinib_150_20hrs_Calcein_FM464'
+ '24-06-2016_TimeLapse3_Axtinib_150_21hrs_Calcein_FM464'
+ '23-06-2016_Timelapse_3_Axitinib_36hr_Calcein_FM464'
+ '30-06-2016_TimeLapse1_Pos0_Axtinib_150_19hrs_Calcein_FM464'
+ '30-06-2016_TimeLapse1_Pos1_Axtinib_150_19hrs_Calcein_FM464'
+ '30-06-2016_TimeLapse1_Pos2_Axtinib_150_19hrs_Calcein_FM464'
+ '01-07-2016_TimeLapse1_Pos0_Axtinib_150_18hrs_Calcein_FM464'
+ '01-07-2016_TimeLapse1_Pos1_Axtinib_150_18hrs_Calcein_FM464'
+ '01-07-2016_TimeLapse2_Pos2_Axtinib_150_19hrs_Calcein_FM464'
+ '01-07-2016_TimeLapse2_Pos3_Axtinib_150_19hrs_Calcein_FM464'};
+
+ ROOT = 'E:\Smadar\LeverJS';
+
+ slist = [];
+ for i=1:length(dslist)
+ leverFile = fullfile(ROOT,[dslist{i} '.LEVER']);
+ [datasetName,int_means] = get_region_intensities(leverFile);
+ if ( isempty(int_means) )
+ continue;
+ end
+
+ s = struct('datasetName',{datasetName}, 'means',{int_means});
+ slist = [slist; s];
+ end
+
+ save('ds_intensities.mat', 'slist');
+end
+
+function [datasetName,int_means] = get_region_intensities(leverFile)
+ int_means = [];
+
+ %% Load constants
+ conn = database(leverFile, '','', 'org.sqlite.JDBC', 'jdbc:sqlite:');
+ CONSTANTS = Read.getConstants(conn);
+
+ %% Estimate parallelogram mask for first-frame image
+ datasetName = CONSTANTS.imageData.DatasetName;
+
+ imagePath = CONSTANTS.imageData.imageDir;
+ chkimd = MicroscopeData.ReadMetadataFile(fullfile(imagePath,filesep));
+ if ( isempty(chkimd) )
+ warning(['Unable to load image data for: ' datasetName]);
+ return;
+ end
+
+ im = MicroscopeData.Reader(imagePath, 'timeRange',[1,1]);
+
+ parallelMask = Segment.llsm_mask(im);
+ % HACK: Drop messy top/bottom slices in mask
+ parallelMask(:,:,1) = false;
+ parallelMask(:,:,end) = false;
+
+ %% Get ecto/endo regions
+ dataPath = 'Data';
+ maskPath = fullfile(dataPath, 'ectoMasks');
+
+ ectoPath = fullfile(maskPath, [datasetName '_ectoregion.mat']);
+ endoPath = fullfile(maskPath, [datasetName '_endoregion.mat']);
+
+ bHasEcto = exist(ectoPath, 'file');
+ bHasEndo = exist(endoPath, 'file');
+
+ endoInfo = [];
+ bwEndo = false(size(parallelMask));
+
+ %%%%%%%%%%%%%%%%%%%%
+ %% Exit if there's not ectoderm region info
+ if ( ~bHasEcto )
+ warning(['Unable to load ectoderm region info for: ' datasetName]);
+ return;
+ end
+
+ %% Load ectoderm and endoderm region info
+ ectoInfo = load(ectoPath);
+ bwEcto = (parallelMask & ectoInfo.regionInfo.bwMask(:,:,:,1));
+
+ if ( bHasEndo )
+ endoInfo = load(endoPath);
+ bwEndo = (parallelMask & endoInfo.regionInfo.bwMask(:,:,:,1));
+ end
+
+ bwSkel = (parallelMask & ~bwEcto & ~bwEndo);
+
+ %% Get average intensities for first and first 50 frames
+ int_means = struct('ff_ecto_mean',{},'ff_skel_mean',{},'ecto_mean',{},'skel_mean',{});
+
+ mxFrame = min(50, CONSTANTS.imageData.NumberOfFrames);
+
+ ecto_m = zeros(mxFrame,size(im,4));
+ skel_m = zeros(mxFrame,size(im,4));
+ parfor t=1:min(mxFrame)
+ im = MicroscopeData.Reader(imagePath, 'timeRange',[t,t]);
+ % For simplicity we use mean of means for intensity
+
+ ecto_mt = zeros(1,size(im,4));
+ skel_mt = zeros(1,size(im,4));
+ for c = 1:size(im,4)
+ imC = im(:,:,:,c);
+
+ ecto_mean = mean(imC(bwEcto));
+ skel_mean = mean(imC(bwSkel));
+
+ ecto_mt(1,c) = ecto_mean;
+ skel_mt(1,c) = skel_mean;
+ end
+
+ ecto_m(t,:) = ecto_mt;
+ skel_m(t,:) = skel_mt;
+ end
+
+ for c=1:size(im,4)
+ int_means(c).ff_ecto_mean = ecto_m(1,c);
+ int_means(c).ff_skel_mean = skel_m(1,c);
+ int_means(c).ecto_mean = mean(ecto_m(:,c));
+ int_means(c).skel_mean = mean(skel_m(:,c));
+ end
+end
--
GitLab