MigratingCellCounting.m

clear % clear workspace to avoid problems with global variables
%open czi image (stitiched or unstitched)

cziImages = bfopen();

%determine how many dimensions the czi has, i.e. how many series and their
%sizes

sizeCziImages = size(cziImages);
listOfImageDimensions = [];

for i = 1:sizeCziImages(1)
    layer = cziImages{i}; % 144 x 2 matrix
    imgDimensions = size(layer{1, 1});
    listOfImageDimensions = [listOfImageDimensions, imgDimensions];
end

%find the image with the dimensions closest to 1050 - this is because the
%original bead finding code was written for images of size 1050x1050
%pixels, so everything is normalized to this size
listOfImageDimensions2 = listOfImageDimensions(1:2:end);
n=1050;
[val,idx]=min(abs(listOfImageDimensions2-n));
imagePixelSize=listOfImageDimensions2(idx);

% extract channels from series of the largest image - this will give us the
% best resolution for seeing the migratory cells and sprouts. The first
% image series is always the largest one
biggestImage = cziImages(1);
biggestImage = biggestImage{1,1};
biggestImageSize = [listOfImageDimensions(1), listOfImageDimensions(2)];
pyramidBase = listOfImageDimensions2(1)/listOfImageDimensions2(idx);

channel_1 = biggestImage(1:3:end, 1); % series 2 channel 1 is BF
channel_2 = biggestImage(2:3:end, 1); % series 2 channel 2 is green fluorescent
channel_3 = biggestImage(3:3:end, 1); % series 2 channel 3 is red fruorescent

% Extract green fluorescent channel from largest image series
newImage = zeros(biggestImageSize);
for i = 1:size(channel_2) 
   concatImage = cat(3, channel_2{i}, newImage);
   newImage = concatImage;
end

% find max projection of green channel stack using concat matrix
[C,I] = max(newImage, [], 3); % c is the max projection, I is the index

% Extract red fluorescent from largest image series
newImageC2 = zeros(biggestImageSize);
for i = 1:size(channel_3)
   concatImage2 = cat(3, channel_3{i}, newImageC2);
   newImageC2 = concatImage2;
end

% find the max projection of the red stack 
[C2,I2] = max(newImageC2, [], 3); % c2 is the max projection, I2 is the index

% make and RGB color map
redChannel = C2;
blueChannel = zeros(biggestImageSize);
greenChannel = C;  % we don't use the blue channel so this is all zeros

rgbImage = cat(3, redChannel, greenChannel, blueChannel);


% Get higher quality BF img using series 2 (2100x2100) 
bfImage = round(size(channel_1)/2);
bfImage = channel_1{bfImage};

% Overlay brighfield image and RGB fluorescent images
bf3D = cat(3, bfImage, bfImage, bfImage); % make bf 3x3 by adding zeros to make same dimensions as rgbImg fluorescent channels
AllChannelsImg = imadd(bf3D, rgbImage); % composite bf and RGB together 

% Loop through every tile so that it is zoomed in enough for the user to
% select the migratory cells... when the 'return' key is pressed, it will
% switch to the next tile

sizebf = size(bfImage);
axisLimitNumber = sizebf(1);

%create a list with 'M' and 'F' to denote which sex the user will be
%looking for when picking out migratory cells
MFList = ['A', 'B'];
sizeMFList = size(MFList);

%create empty lists to append locations of selected migratory cells, these
%will be used to obtain the cell count later
FemaleCellsX = [];
MaleCellsX = [];
FemaleCellsY = [];
MaleCellsY = [];

%loop through everything twice: once to select male cells and once to
%select female cells
for L = 1:sizeMFList(2)
    clear global %clear globals becasue they don't autonmatically reset each loop
    
    global ENTER_IS_PRESSED
    ENTER_IS_PRESSED = 0;

    global migratoryCellCoords

    figure = imshow(AllChannelsImg); %figure used is overlayed image with bf, red, and green

    set(gcf, 'KeyPressFcn', @EnterKeyPressFcnCallback)

    migratoryCellIndiciesList = []; %blank list to store indices of migratory selected cells

    %set axis limits for panning across image
    X1 = 0;
    X2 = axisLimitNumber/9; %divide the size of the image into 9 tiles - all images are stitched into 9 tiles
    Y1 = axisLimitNumber - axisLimitNumber/9;
    Y2 = axisLimitNumber;
    for b  = 1:81 %81 is total number of tiles in the image
        title("Tile # " + b + " " + MFList(L) + " ");
        xlim([X1 X2]);
        ylim([Y1 Y2]);

        figure.ButtonDownFcn = @(~,~) migratoryCellLocationButtonPressCallback(figure.Parent);

        while ~ENTER_IS_PRESSED
           drawnow

        end

        ENTER_IS_PRESSED = 0;


        s = size(migratoryCellCoords);
        s = s(2);

        migratoryCellIndicesListSize = size(migratoryCellIndiciesList); % this gets the size of the migratoryCellIndicesList
        migratoryCellIndicesListSize = migratoryCellIndicesListSize(2); % This gets second dimension of migratoryCellIndicesListSize
        for numberOfIndicies = 1:((s/2)-migratoryCellIndicesListSize)
            migratoryCellIndiciesList = [migratoryCellIndiciesList, b];
        end

        if rem(b,9) == 0 %every time b is divisible by 9, so once one row of 9 tiles is complete reset x limits and movce to a different row
            X1 = 0;
            X2 = axisLimitNumber/9;
            Y1 = Y1-axisLimitNumber/9;
            Y2 = Y2-axisLimitNumber/9;
        else % if not at end of row, stay in same row and move x limits to new tile
            X1 = X1+axisLimitNumber/9;
            X2 = X2+axisLimitNumber/9;
        end
    end
    migratoryCellX = migratoryCellCoords(1:2:end);
    migratoryCellY = migratoryCellCoords(2:2:end);

    % Display the RGB/BF image with the migratory cells shown as red * -
    % this doesnt really do anything, except keep the image from stalling
    % when switching from M to F
    imshow(AllChannelsImg);
    hold on
    scatter(migratoryCellX, migratoryCellY, 'r*');
    
    %store info in table form
    MigratoryCellTable = table(migratoryCellIndiciesList.', migratoryCellX.', migratoryCellY.', 'VariableNames', {'Bead#', 'migratoryCell_x_coord', 'migratoryCell_y_coord'});

    % This code deletes right clicked migratoryCells from the list of all migratoryCells and
    % then removes them from the migratoryCellTable
    global deletedmigratoryCells

    deletedsX = deletedmigratoryCells(1:2:end);
    deletedsY = deletedmigratoryCells(2:2:end);

    if ~ isempty(deletedsX)
        for i = 1:size(deletedsX,2)
            deletedIndexS = find(MigratoryCellTable.migratoryCell_x_coord==deletedsX(i));
            MigratoryCellTable(deletedIndexS, :) = [];
        end
    end
    
    %if the loop is looking for female cells, append all x and y values to
    %the female lists, otherwise apend to the male lists
    if strcmpi(MFList(L), 'B')
        FemaleCellsX = [FemaleCellsX, MigratoryCellTable.migratoryCell_x_coord];
        FemaleCellsY = [FemaleCellsY, MigratoryCellTable.migratoryCell_y_coord];
    else
        MaleCellsX = [MaleCellsX, MigratoryCellTable.migratoryCell_x_coord];
        MaleCellsY = [MaleCellsY, MigratoryCellTable.migratoryCell_y_coord];
    end
    
end

close %close the figure

%calculate number of migratory cells from their coordinates
numberOfFMigratoryCells = size(FemaleCellsX);
FMigratoryCellCount = numberOfFMigratoryCells(1);
numberOfMMigratoryCells = size(MaleCellsX);
MMigratoryCellCount = numberOfMMigratoryCells(1);

%store cell counts in table form
OutputTable = table(FMigratoryCellCount, MMigratoryCellCount, 'VariableNames', {'Population_A_count', 'Population_B_count'});

% export data in excel
filename = 'MigratoryCellCount.xlsx';
writetable(OutputTable,filename,'Sheet',1);

clear % clear workspace to avoid problems with global variables
%Functions needed
function EnterKeyPressFcnCallback(~ , eventdata) % record the sex of the bead when either an 'M' or 'F' key is pressed... m for male and f for female
     global ENTER_IS_PRESSED
     
%      global beadSex
     if strcmpi(eventdata.Key, 'Return')
         ENTER_IS_PRESSED = 1;
%          beadSex = [beadSex, 'M'];
%      end
%      if strcmpi(eventdata.Key, 'F')
%          SEX_IS_CHOSEN = 1;
%          beadSex = [beadSex, 'F'];
     end
end

function deletedmigratoryCellsCallback(~,eventData) %when a migratoryCell is right clicked, this deletes it from the image and adds it to a list of migratoryCells that should be deleted completely
    global deletedmigratoryCells
    if strcmp(eventData.SelectionType, 'right') %if the roi is right clicked, delete the clicked selection
        set(eventData.Source,'visible','off');
        
        deletedmigratoryCell = eventData.Source.Position;
        deletedmigratoryCells = [deletedmigratoryCells, deletedmigratoryCell];

    end
end

function migratoryCellLocationButtonPressCallback(idk) % records the location of the migratoryCell when a spot is clicked, and displays a point ROI on that spot
    global migratoryCellCoords

   coordinates = idk.CurrentPoint; 
   coordinates = [coordinates(1,1) coordinates(1,2)];
   
   h = drawpoint('Position',coordinates);
   
   addlistener(h,'ROIClicked',@deletedmigratoryCellsCallback);
   
   coordinates = coordinates(1,1:2);
   migratoryCellCoords = [migratoryCellCoords, coordinates];

end