Difference between revisions of "Calculating MSD and Diffusion Coefficients"
From Course Wiki
Juliesutton (Talk | contribs) (Created page with "Category:20.309 Category:Matlab {{Template:20.309}} Example MSD calculation of sum and difference particles using MATLAB. T...") |
MAXINE JONAS (Talk | contribs) (→Diffusion measurements from particle tracking data) |
||
(7 intermediate revisions by one user not shown) | |||
Line 5: | Line 5: | ||
[[Image:stabilityFigure.jpg|thumb|right|Example MSD calculation of sum and difference particles using MATLAB.]] | [[Image:stabilityFigure.jpg|thumb|right|Example MSD calculation of sum and difference particles using MATLAB.]] | ||
− | This page | + | This page contains example MATLAB code for calculating the mean squared displacement (MSD) and diffusion coefficient from a movie of fluorescent microspheres. |
+ | <onlyinclude> | ||
==Diffusion measurements from particle tracking data== | ==Diffusion measurements from particle tracking data== | ||
− | + | The function <tt>CalculateDiffusionCoefficientFromTrajectories</tt> can be used to plot the MSD as a function of time for fluorescent particles and calculate their average diffusion coefficient. The input argument <tt>Centroids</tt> can be generated by the function <tt>track</tt>, which creates a 4 column matrix listing the x- and y- coordinates, the frame number and the particle number. Additional input arguments to <tt>CalculateDiffusionCoefficientFromTrajectories</tt> include the pixel size, the movie frame rate, and an optional figure handle for plotting the MSD. | |
− | + | <pre>function out = CalculateDiffusionCoefficientFromTrajectories(Centroids, PixelSize, FrameRate, FigureHandle) | |
− | <pre> | + | |
− | function out = CalculateDiffusionCoefficientFromTrajectories(Centroids, PixelSize, FrameRate, FigureHandle) | + | |
Centroids( :, 1:2 ) = Centroids( :, 1:2 ) * PixelSize; | Centroids( :, 1:2 ) = Centroids( :, 1:2 ) * PixelSize; | ||
Line 21: | Line 20: | ||
longestTrack = 0; | longestTrack = 0; | ||
+ | particlesToSkip = []; | ||
+ | |||
if( nargin < 4 || isempty( FigureHandle )) | if( nargin < 4 || isempty( FigureHandle )) | ||
FigureHandle = figure; | FigureHandle = figure; | ||
end | end | ||
− | figure(FigureHandle) | + | figure( FigureHandle ) |
+ | |||
for ii = 1:numberOfParticles | for ii = 1:numberOfParticles | ||
% get the centroids for a single particle | % get the centroids for a single particle | ||
Line 38: | Line 40: | ||
particles{ii}.number = numberOfMsds; | particles{ii}.number = numberOfMsds; | ||
− | for tau = 1:numberOfMsds | + | if numberOfMsds==0 |
− | + | particlesToSkip = [particlesToSkip ii]; | |
− | + | else | |
− | + | ||
− | + | for tau = 1:numberOfMsds | |
− | + | dx = particle(1:1:(end-tau), 1) - particle((tau+1):1:end, 1); | |
− | + | dy = particle(1:1:(end-tau), 2) - particle((tau+1):1:end, 2); | |
− | + | drSquared = dx.^2+dy.^2; | |
+ | msd = mean( drSquared ); | ||
+ | |||
+ | particles{ii}.msd(tau) = msd; | ||
+ | particles{ii}.msdUncertainty(tau) = std( drSquared ) / sqrt(length( drSquared )); | ||
+ | end | ||
+ | |||
+ | particles{ii}.diffusionCoefficient = particles{ii}.msd(1) / (2 * 2 * FrameRate^-1); | ||
+ | |||
+ | longestTrack = max( longestTrack, numberOfSamples ); | ||
+ | loglog( (1:numberOfMsds)/FrameRate, particles{ii}.msd(1:numberOfMsds), 'color', rand(1,3) ); | ||
+ | hold on; | ||
end | end | ||
− | |||
− | |||
− | |||
− | |||
− | |||
− | |||
end | end | ||
Line 61: | Line 68: | ||
for ii = 1:numberOfParticles | for ii = 1:numberOfParticles | ||
− | numberOfMsds = length(particles{ii}.msd); | + | if sum(particlesToSkip==ii)==0 |
− | + | numberOfMsds = length(particles{ii}.msd); | |
− | + | d(ii) = particles{ii}.diffusionCoefficient; | |
− | + | averageMsd(1:numberOfMsds) = averageMsd(1:numberOfMsds) + particles{ii}.msd; | |
+ | averageMsdCount(1:numberOfMsds) = averageMsdCount(1:numberOfMsds) + 1; | ||
+ | end | ||
end | end | ||
− | figure(FigureHandle); | + | figure( FigureHandle ); |
averageMsd = averageMsd ./ averageMsdCount; | averageMsd = averageMsd ./ averageMsdCount; | ||
loglog((1:longestTrack)/FrameRate, averageMsd, 'k', 'linewidth', 3); | loglog((1:longestTrack)/FrameRate, averageMsd, 'k', 'linewidth', 3); | ||
− | ylabel('MSD ( | + | ylabel('MSD (m^2)'); |
xlabel('Time (s)'); | xlabel('Time (s)'); | ||
title('MSD versus Time Interval'); | title('MSD versus Time Interval'); | ||
Line 78: | Line 87: | ||
out.diffusionCoefficient = mean(d); | out.diffusionCoefficient = mean(d); | ||
out.diffusionCoefficientUncertainty = std(d)/sqrt(length(d)); | out.diffusionCoefficientUncertainty = std(d)/sqrt(length(d)); | ||
− | out.particles = particles; | + | out.particles = particles;</pre> |
− | </pre> | + | </onlyinclude> |
Latest revision as of 12:21, 10 March 2020
This page contains example MATLAB code for calculating the mean squared displacement (MSD) and diffusion coefficient from a movie of fluorescent microspheres.
Diffusion measurements from particle tracking data
The function CalculateDiffusionCoefficientFromTrajectories can be used to plot the MSD as a function of time for fluorescent particles and calculate their average diffusion coefficient. The input argument Centroids can be generated by the function track, which creates a 4 column matrix listing the x- and y- coordinates, the frame number and the particle number. Additional input arguments to CalculateDiffusionCoefficientFromTrajectories include the pixel size, the movie frame rate, and an optional figure handle for plotting the MSD.
function out = CalculateDiffusionCoefficientFromTrajectories(Centroids, PixelSize, FrameRate, FigureHandle) Centroids( :, 1:2 ) = Centroids( :, 1:2 ) * PixelSize; numberOfParticles = max( Centroids(:,4) ); particles = cell( 1, numberOfParticles ); longestTrack = 0; particlesToSkip = []; if( nargin < 4 || isempty( FigureHandle )) FigureHandle = figure; end figure( FigureHandle ) for ii = 1:numberOfParticles % get the centroids for a single particle particle = Centroids(Centroids(:,4) == ii, :); numberOfSamples = size(particle,1); % number of samples in trajectory particles{ii} = struct(); numberOfMsds = round(numberOfSamples-1);% / 8); particles{ii}.msd = zeros(1, numberOfMsds); particles{ii}.msdUncertainty = zeros(1, numberOfMsds); particles{ii}.number = numberOfMsds; if numberOfMsds==0 particlesToSkip = [particlesToSkip ii]; else for tau = 1:numberOfMsds dx = particle(1:1:(end-tau), 1) - particle((tau+1):1:end, 1); dy = particle(1:1:(end-tau), 2) - particle((tau+1):1:end, 2); drSquared = dx.^2+dy.^2; msd = mean( drSquared ); particles{ii}.msd(tau) = msd; particles{ii}.msdUncertainty(tau) = std( drSquared ) / sqrt(length( drSquared )); end particles{ii}.diffusionCoefficient = particles{ii}.msd(1) / (2 * 2 * FrameRate^-1); longestTrack = max( longestTrack, numberOfSamples ); loglog( (1:numberOfMsds)/FrameRate, particles{ii}.msd(1:numberOfMsds), 'color', rand(1,3) ); hold on; end end % compute the average value of D and the ensemble average of msd d = zeros(1,numberOfParticles); averageMsd = zeros( 1, longestTrack ); averageMsdCount = zeros( 1, longestTrack); for ii = 1:numberOfParticles if sum(particlesToSkip==ii)==0 numberOfMsds = length(particles{ii}.msd); d(ii) = particles{ii}.diffusionCoefficient; averageMsd(1:numberOfMsds) = averageMsd(1:numberOfMsds) + particles{ii}.msd; averageMsdCount(1:numberOfMsds) = averageMsdCount(1:numberOfMsds) + 1; end end figure( FigureHandle ); averageMsd = averageMsd ./ averageMsdCount; loglog((1:longestTrack)/FrameRate, averageMsd, 'k', 'linewidth', 3); ylabel('MSD (m^2)'); xlabel('Time (s)'); title('MSD versus Time Interval'); %hold off; out = struct(); out.diffusionCoefficient = mean(d); out.diffusionCoefficientUncertainty = std(d)/sqrt(length(d)); out.particles = particles;