Add files via upload

NOBIAS_preparedata_angle.m

@@ -0,0 +1,63 @@
+function data=NOBIAS_preparedata_angle(AllTracks)
+% input All tracks is a N by 1 cell array, where each element of the cell
+% array is a tracectories with at least 4 colmns, where the 2nd columns
+% should correspond to the frame of the tracks (gaps are allowed and 
+% should be pre-filtered in tracking step to avoid too huge gaps), the 3rd
+% and 4th columns should be the rows and columns of the 2D tracks
+% coordinates, pay attention to x-y and row-col difference.
+
+TrID=[];
+All_steps={};
+All_corrstep={};
+All_corrstep_angle = {};
+min_length=3;
+for i=1:length(AllTracks)
+    temptr=AllTracks{i};
+    fixedTrack = nan(max(temptr(:,2)),size(temptr,2));
+    fixedTrack(temptr(:,2),:) = temptr;
+    fixedTrack(1:find(all(isnan(fixedTrack),2)==0,1,'first')-1,:)=[];
+    tempstep=fixedTrack(2:end,[4,3])-fixedTrack(1:end-1,[4,3]);
+    tempcorrstep=[tempstep(1:end-1,:).*tempstep(2:end,:); nan, nan];
+    temp_angle = get_steps_angle(tempstep);
+    gapsID=(~isnan(tempstep(:,1)))&(~isnan(tempstep(:,2)));
+    tempstep=tempstep(~isnan(tempstep(:,1)),:);
+    tempstep=tempstep(~isnan(tempstep(:,2)),:);
+    tempcorrstep=tempcorrstep(gapsID,:);
+    temp_angle = temp_angle(gapsID);
+    if size(tempstep,1)>min_length
+        All_steps{end+1}=tempstep;
+        All_corrstep{end+1}=tempcorrstep;
+        All_corrstep_angle{end+1} = temp_angle;
+    end
+end
+for i=1:length(All_steps)
+    TrID=[TrID, ones(1,size(All_steps{i},1))*i];
+end
+data.obs=cat(1,All_steps{:})';
+data.TrID=TrID;
+data.obs_corr=cat(1,All_corrstep{:})'; %not needed if motion blur correction not wanted
+data.obs_angle = cat(1,All_corrstep_angle{:})';
+% scale data to certain varience
+data = NOBIAS_scale_data(data);
+end
+
+function angle = get_steps_angle(vectors)
+    angle = zeros(length(vectors),1);
+    for i=1:length(vectors)-1
+        vector1=vectors(i,:);
+        vector2=vectors(i+1,:);
+        dot_product = dot(vector1, vector2);
+        magnitude_product = norm(vector1) * norm(vector2);
+    if abs(dot_product / magnitude_product) > 1
+        angle(i) = NaN;  % Invalid input, return NaN
+    else
+        % Calculate the angle in radians
+        angle_rad = acos(dot_product / magnitude_product);
+
+        % Convert the angle to degrees
+        angle(i) = rad2deg(angle_rad);
+
+    end
+    end
+    angle(end)=nan;
+end

NOBIAS_stackbar.m

@@ -0,0 +1,74 @@
+function NOBIAS_stackbar(out, data)
+
+
+
+c=colormap('lines');
+pixel_size = 49;
+raw_obs = data.obs * data.scale_factor * pixel_size ;
+raw_displacements = sqrt(raw_obs(1,:).^2 + raw_obs(2,:).^2);
+used_state = unique(out.reord_stateSeq);
+figure;
+hold on;
+
+numBins = 20;
+
+% Get the bin edges
+edges = linspace(min(raw_displacements), max(raw_displacements), numBins+1);
+
+% Initialize the matrix to store counts for each class
+counts = zeros(numBins, 3);
+
+% Count the number of elements in each bin for each class
+for i = 1:length(used_state)
+    state_disp= raw_displacements(out.reord_stateSeq == i);
+    counts(:, i) = histcounts(state_disp, edges);
+end
+
+h = bar(edges(1:end-1) + diff(edges)/2, counts, 'stacked');
+
+for i = 1:length(used_state)
+    h(i).FaceColor = c(i, :);
+end
+legend('State 1', 'State 2', 'State 3');
+title('Cumulative step displacement histogram by state');
+xlabel('Displacement, nm');
+ylabel('Occurance');
+hold off
+
+
+% Calculate the size of each batch
+[uniqueTrIDs, ~, trIndices] =  unique(data.TrID);
+trackSizes = accumarray(trIndices, 1);
+
+% Find unique batch sizes
+uniquetrackSizes = unique(trackSizes);
+
+% Initialize the matrix to store counts for each class in each batch size
+numtrackSizes = length(uniquetrackSizes);
+numStates = length(unique(out.reord_stateSeq));
+counts = zeros(numtrackSizes, numStates);
+
+% Count the number of elements in each batch size for each class
+for i = 1:numtrackSizes
+    trackSize = uniquetrackSizes(i);
+    trackIDsWithThisSize = uniqueTrIDs(trackSizes == trackSize);
+    for j = 1:numStates
+        counts(i, j) = sum(ismember(data.TrID, trackIDsWithThisSize) & (out.reord_stateSeq == j));
+%         weightedCounts(i, j) = sum(ismember(data.TrID, trackIDsWithThisSize) & (out.reord_stateSeq == j)) * trackSize;
+    end
+end
+
+% Plot the cumulative bar graph (stacked histogram)
+figure;
+h = bar(uniquetrackSizes, counts, 'stacked');
+title('Cumulative Bar Graph by Track Length');
+xlabel('Track Length');
+ylabel('Occurance');
+
+
+for i = 1:numStates
+    h(i).FaceColor = c(i, :);
+end
+legend('State 1', 'State 2', 'State 3');
+
+end

NOBIAS_stepscatter.m

@@ -0,0 +1,27 @@
+function NOBIAS_stepscatter(out, data)
+
+c=colormap('lines');
+pixel_size = 49;
+raw_obs = data.obs * data.scale_factor * pixel_size ;
+raw_displacements = sqrt(raw_obs(1,:).^2 + raw_obs(2,:).^2);
+used_state = unique(out.reord_stateSeq);
+figure;
+hold on;
+jitterAmount = 0.5; % Adjust as needed
+jitter = (rand(size(out.reord_stateSeq)) - 0.5) * jitterAmount;
+
+figure;
+hold on
+for i=1:length(used_state)
+    index = out.reord_stateSeq == i;
+    scatter(out.reord_stateSeq(index) + jitter(index), raw_displacements(index), 5,c(i,:),'filled',...
+        'MarkerFaceAlpha', sum(index)/length(raw_displacements));
+end
+title('Step Displacement by diffusive states');
+xticks(used_state);
+xticklabels(arrayfun(@num2str , used_state, 'UniformOutput', false));
+ylabel('Displacement, nm')
+set(gca,'yscale','log')
+ylim([10 1000])
+
+end