trackingManual4
Previous (Tracking Neutrophils with Matlab 3) / Index / Next (Tracking Neutrophils with Matlab 5)

Tracking Neutrophils with Matlab: 4

Contents

Exploring Stats a bit further

There are some metrics that depend on the wound region, in particular the size of the wound region can alter the time at which a neutrophil arrives, or the quantity that arrive. There are some other metrics that do not depend on this manually delineated region

For the following metrics, when a neutrophil is moving towards the wound area is considered to be in translation, and when it has reached the wound, it is considered to be in exploration


Metrics

  • 1 absolute number of Neutrophils in Exploration
  • 2 relative number of Neutrophils in Exploration
  • 3 forwardRatio Total
  • 4 forwardRatio Exp
  • 5 forwardRatio Trans
  • 6 inExplorationRatio
  • 7 inWoundRatio Tot
  • 8 IdleWoundRatio
  • 9 backwardRatioTot
  • 10 LeaveWoundRatio
  • 11 LeaveWoundRatio2
  • 12 LeaveWoundRatio3
  • 13 staticRatio
  • 14 returnRatio
  • 15 returnRatio2

To calculate these metrics, for all movements of cells between frames, the orientation (angle) was calculated and a histogram revealed two main orientations, one that corresponded *towards* the wound, and another opposite. For every hop, the distance and angle were calculated:

  • absDistPerHop Absolute distance per hop
  • anglePerHop angle per hop
  • avDisHop average of all hops

Using the peak of towards, the X,Y coordinates were rotated so that the wound was on the right and X dimension indicated movement to/from wound and Y indicated lateral movements. From this rotated coordinate system we get several measurements:

  • maxDispl maximum displacement towards the wound (extreme X)
  • oriDistPerHop oriented distance per hop, projection to X axis, which is always lower or equal to the absolute distance and may change direction
  • latDistPerHop lateral distance per hop, projection to Y axis
  • furthestPoint 0 for all points BEFORE cell reaches maxDispl, 1 afterwards
  • furthestPoint2 1 for all points that are within n*average Hop from maxDispl, 0 otherwise
  • effDistPerHop Effective distance per hop = oriented distance / absolute distance, will vary between [-1,+1].

We also need a record if the cell is inside or outside the wound (inWound), 1 inside, 0 outside and changes as the cell moves inside or outside and another variable that records just if the cell reached the wound (networkReachWound) 0 before the cell reaches the wound and 1 afterwards.

With those metrics, we then calculate for every track, the measurements that summarise the behaviour:

  • 3 forwardRatioTot= sum(effDistPerHop>0.6)./numHops;
  • 4 forwardRatioE = sum((inWound.*effDistPerHop)>0.6)./ numHopsExploration;
  • 5 forwardRatioT = sum(((posFinalNetwork-networkReachWound).*effDistPerHop)>0.6)./numHopsTranslation;
  • 6 inExplorationRatio = numHopsInWound./ numHops;
  • 7 inWoundRatio= numHopsInWound./ numHopsExploration;
  • 8 IdleWoundRatio = sum(inWound.*(absDistPerHop<0.7))./numHopsInWound;
  • 9 backwardRatioTot = sum(oriDistPerHop<-0.6)./ numHops;
  • 10 LeaveWoundRatio= sum(networkReachWound.*(oriDistPerHop<-1) ) ./ numHopsExploration;
  • 11 LeaveWoundRatio2 = sum(networkReachWound.*(effDistPerHop<-0) ) ./ numHopsExploration;
  • 12 LeaveWoundRatio3 = sum(networkReachWound.*(oriDistPerHop<-1) ) ./ sum(networkReachWound.*(absDistPerHop>1) );
  • 13 staticRatio = (sum(furthestPoint2)-1) ./ numHops;
  • 14 returnRatio = sum(furthestPoint.*(oriDistPerHop<0) ) ./ sum(furthestPoint);
  • 15 returnRatio2= sum(furthestPoint.*(oriDistPerHop<-(3.5*avDisHop)) ) ./ sum(furthestPoint);
  • 17 average of absolute velocity per track
  • 18 average of absolute velocity per track in Exploration
  • 19 average of absolute velocity per track in Translation
  • 20 average of oriented velocity per Track
  • 21 average of oriented velocity per Track in Exploration
  • 22 average of oriented velocity per Track in Translation
  • 23 average of lateral velocity per Track
  • 24 average of lateral velocity per Track in Exploration
  • 25 average of lateral velocity per Track in Translation
  • 26 how many tracks enter the wound and leave
  • 27 how many tracks enter the wound and leave / total number of tracks

To access these metrics (from where you can copy to excel if needed!) you need to open the handles by double clicking on that variable on the workspace:









Then you can either read the values stored there (rows and columns of the data, for instance), or double click on the fields of interest. In this case we are interested on handles.distMaps and handles.distanceNetwork:









The fields there will contain either:

  • (a) the average value for the field,
  • (b) an average value for the track, in these cases you will see that the value is something like <1 x 222 double>, that is, there were 222 tracks registered and there is one double value for each of them, or
  • (c) there are values for each hop of a neutrophil along the movement in the fish. In these cases, you have something like <324 x 222 double>. This means that the longest track was 324 hops, and there were 222 tracks and for each hop and every track there is a value. Since no all tracks will be equally long, the rest of the column for each track will be filled with zeros.

If there is a metric for which a certain track does not qualify (for instance, it did not reach the wound and the metric is related to the neutrophil reaching the wound) it will have a value of NaN (not a number). These are to be ignored. Look at the following examples:













Tracks characteristics as "maps"

A different way to visualise the data is through "maps" in which the data calculated per region and (if the suitable data is available) overlaid on the data. These maps require more fields that are generated if the "wound" can be determined. The wound allows the calculation of the distance projected towards the line of direct distance, that is the oriented distance, or the perpendicular movement, the lateral distance.

load NeutropTest_mat_Re/T00001.mat

imagesc(dataR(:,:,8))
woundRegion                                 = roipoly();

roipoly selects a region of interest (roi) with a polygon. With the mouse one can select the region by clicking points, and a double-click finishes the selection.


Once this wound region has been selected it is passed as a variable to the following functions

handles       = effectiveDistance(handles,woundRegion);

disp (handles)
             numFrames: 9
                  rows: 250
                  cols: 250
                  levs: 24
               minBlob: 10
            thresLevel: [2.1054e+03 3.1856e+03]
           nodeNetwork: [132x31 double]
          neighNetwork: [9x15 double]
        reducedNetwork: [9x15 double]
         linkedNetwork: [9x15 double]
          finalNetwork: [9x15 double]
       distanceNetwork: [1x1 struct]
            finalLabel: [1 2 3 4 5 6 7 8 9 10 11 12 13 14 15]
    explorationNetwork: [9x15 double]
    translationNetwork: [9x15 double]
        inWoundNetwork: [9x15 logical]



To understand better the differences between the absolute, lateral and oriented velocities, it is useful to think of the images as a Cartesian plane, where one of the axis is determined by the orientation of the wound that has been previously selected, and the other axis, perpendicular to the wound should be parallel to the body of the fish as shown in the figure below.

Once this orientation has been determined, the movement of the neutrophils is described in terms of the movement towards the wound, that is the oriented movement and the movement parallel to the wound, the lateral movement. Together these two components will form a total or absolute movement. From these movements the corresponding velocities are derived. In the example below, the first hop in a series of displacements of a Neutrophil has a larger oriented than lateral component. In the second hop, the movement is almost completely oriented towards the wound. The third hop moves away from the wound with a negative oriented velocity and a large lateral velocity.




The calculation of effectiveDistance adds handles.translationNetwork (those movements before the neutrophil reaches the wound), handles.explorationNetwork (those movements after the neutrophil reaches the wound) and handles.inWoundNetwork (a register if the neutrophil is inside or outside the wound). The following fields are created in handles.distanceNetwork:


% numHopsTranslation,       = number of hops in Translation / 
% numHopsExploration        = number of hops in Exploration 
% perHopT, perHopE,         = distance of each hop in Translation / Exploration        
% perHopOrientedT,          = oriented distance in Translation
% perHopEffectiveT,         = effective distance in Translation
% perHopOrientedE,          = oriented distance in Exploration
% perHopEffectiveE,         = effective distance in Exploration
% forwardRatioE,            = Forward ratio in Exploration
% forwardRatioT,            = Forward ratio in Translation
% forwardRatioTot,          = Forward ratio in Total



To calculate the maps, use:

handles      = effectiveTracks(handles,woundRegion);

disp (handles)
Warning: Divide by zero.
Warning: Divide by zero.
Warning: Divide by zero.
Warning: Divide by zero.
Warning: Divide by zero.
Warning: Divide by zero.
Warning: Divide by zero.
Warning: Divide by zero.
Warning: Divide by zero.
Warning: Divide by zero.
Warning: Divide by zero.
Warning: Divide by zero.
Warning: Divide by zero.
Warning: Divide by zero.
             numFrames: 9
                  rows: 250
                  cols: 250
                  levs: 24
               minBlob: 10
            thresLevel: [2.1054e+03 3.1856e+03]
           nodeNetwork: [132x31 double]
          neighNetwork: [9x15 double]
        reducedNetwork: [9x15 double]
         linkedNetwork: [9x15 double]
          finalNetwork: [9x15 double]
       distanceNetwork: [1x1 struct]
            finalLabel: [1 2 3 4 5 6 7 8 9 10 11 12 13 14 15]
    explorationNetwork: [9x15 double]
    translationNetwork: [9x15 double]
        inWoundNetwork: [9x15 logical]
               inWound: [9x15 double]
              distMaps: [1x1 struct]

Do not worry about the warnings, these indicate that in some cases there are not enough data to calculate a division. The maps are inside handles.distMaps

disp (handles.distMaps)
     absDistPerHop: [9x15 double]
       anglePerHop: [9x15 double]
     oriDistPerHop: [9x15 double]
     furthestPoint: [9x15 double]
    furthestPoint2: [9x15 double]
     latDistPerHop: [9x15 double]
     effDistPerHop: [9x15 double]
       nodeNetwork: [132x37 double]
        absDistMap: [250x250 double]
        oriDistMap: [250x250 double]
          angleMap: [250x250 double]
        latDistMap: [250x250 double]
        effDistMap: [250x250 double]

The maps can be visualised with the command "mesh":

mesh(handles.distMaps.oriDistMap)

For this example, since there are not many time frames, the result does not look very revealing, but here are other examples:

Intensity image and maps combined.

Finally, if you want to display the maps combined with the intensity images of the data, you can follow this code to place the map slightly below the image (good when the map has only positive values)

load NeutropTest_mat_Re/T00001.mat

[rowsP,colsP,levsP]   = size(dataR);
meshdataR       = zeros(rowsP,colsP);
[RR,CC]         = meshgrid(1:colsP,1:rowsP);

dataToSurf=dataR(:,:,6);
dataToSurf=dataToSurf-min(dataToSurf(:));

stepMeshMap=1;
stepMesh=6;

dataToMesh = handles.distMaps.absDistMap;

hold off
surf(RR(1:stepMeshMap:end,1:stepMeshMap:end),CC(1:stepMeshMap:end,1:stepMeshMap:end),min(dataToMesh(:))+0.01+meshdataR(1:stepMeshMap:end,1:stepMeshMap:end), dataToSurf(1:stepMeshMap:end,1:stepMeshMap:end,1),'edgecolor','interp')
hold on

mesh(RR(1:stepMesh:end,1:stepMesh:end),CC(1:stepMesh:end,1:stepMesh:end),dataToMesh(1:stepMesh:end,1:stepMesh:end),'marker','none','LineStyle','-','facecolor','none','edgecolor','b')
axis tight
axis ij
rotate3d on
view(45,40)
title('Absolute velocity')

Or to place the map above the image (good when the map has negative values) use this:

dataToMesh = handles.distMaps.oriDistMap;

hold off
surf(RR(1:stepMeshMap:end,1:stepMeshMap:end),CC(1:stepMeshMap:end,1:stepMeshMap:end),min(dataToMesh(:))+0.01+meshdataR(1:stepMeshMap:end,1:stepMeshMap:end), dataToSurf(1:stepMeshMap:end,1:stepMeshMap:end,1),'edgecolor','interp')
hold on
mesh(RR(1:stepMesh:end,1:stepMesh:end),CC(1:stepMesh:end,1:stepMesh:end),dataToMesh(1:stepMesh:end,1:stepMesh:end),'marker','none','LineStyle','-','facecolor','none','edgecolor','b')
axis tight
axis ij
rotate3d on
view(45,40)
title('Effective velocity')

Previous (Tracking Neutrophils with Matlab 3) / Index / Next (Tracking Neutrophils with Matlab 5)