Difference between revisions of "Framework for estimation of membrane thickness"
(13 intermediate revisions by the same user not shown) | |||
Line 3: | Line 3: | ||
== Reading the particles == | == Reading the particles == | ||
− | The easiest way to bring particles aligned by a project is through the use of the <tt>dpkdata.Reader</tt> object. | + | The easiest way to bring particles aligned by a project is through the use of the <tt>dpkdata.Reader</tt> object. First store the [[table]] and [[data folder]] of your project in variables <tt>myTable</tt> and <tt>myData</tt>. This can be easily done by browsing into the folder structure, or just by using <tt>ddb</tt> to extract items from the project. Option <tt>-d<tt> will just display the file or folder location. |
+ | |||
+ | <nowiki>ddb myProject:rt -d # prints in screen the last available refined table | ||
+ | ddb myProject:data -d % prints on screen the data folder linked to this project;</nowiki> | ||
+ | |||
+ | In any case, | ||
<nowiki >r = dpkdata.Reader(); | <nowiki >r = dpkdata.Reader(); | ||
− | r.setSource('table', myTable,'data',myData); | + | r.setSource('table', myTable,'data',myData);</nowiki> |
− | </nowiki> | ||
you can access the particles throught the <tt>getParticle</tt> method: | you can access the particles throught the <tt>getParticle</tt> method: | ||
+ | |||
<tt>p = r.getParticle(12,'align',1);</tt> | <tt>p = r.getParticle(12,'align',1);</tt> | ||
gets the particle with tag 12 as a variable <tt>p</tt> in the workspace. You can also read all the particles as a [[matlab cell array]] through | gets the particle with tag 12 as a variable <tt>p</tt> in the workspace. You can also read all the particles as a [[matlab cell array]] through | ||
Line 31: | Line 36: | ||
</nowiki> | </nowiki> | ||
− | + | You can create a compact version of this cell array: | |
<tt>pac = mbg.cellCompact(pa);</tt> | <tt>pac = mbg.cellCompact(pa);</tt> | ||
Line 50: | Line 55: | ||
</nowiki> | </nowiki> | ||
− | == Creating | + | === Checking the read particles=== |
+ | |||
+ | We can only measure the membrane distance on aligned particles. At this point we make certain that everything is ok by quickly adding all the particles in the cell array and checking visually the result. | ||
+ | <nowiki> | ||
+ | myAverage = dynamo_sumarray(pac); | ||
+ | dview(myAverage);</nowiki> | ||
+ | |||
+ | [[File:MembraneThicknessAverage.png|thumb|center|500px| Y-view of the averaged membrane]] | ||
+ | |||
+ | == Creating distance filters== | ||
We first create a set of filters that model the distance between the two layers of the membrane. | We first create a set of filters that model the distance between the two layers of the membrane. | ||
− | <tt> t = dpkfw.mdm.templateSeries('size',16,'thickness',1,'distance',[3:7]); </tt> | + | <tt> t = dpkfw.mdm.templateSeries('size',[16,16,16],'thickness',1,'distance',[3:7]); </tt> |
+ | |||
+ | Here, <tt>t</tt> is a [[matlab cell array | cell array]], containing 5 volumes of 16x16x16 pixels each. Each volume contains a "model" of a membrane. The different models are coded with the flag <tt>distance</tt>. In this case, the model number one contains a membrane which a distance between layers of 3 pixels. In the last model, the distance is 7. | ||
+ | |||
+ | You can visualize each membrane distance filter with: | ||
+ | |||
+ | <tt>figure;dslices(t,'projy','*','dim',[1,5]);</tt> | ||
− | + | [[File:MembraneDistanceFilters.png|thumb|center|400px| Y projections of the five filters, with distances from 3 to 7 pixels]] | |
== Applying the filters onto a particle == | == Applying the filters onto a particle == | ||
If you have your particle in the Matlab/''Dynamo'' workspace with name <tt>p</tt> you can type. | If you have your particle in the Matlab/''Dynamo'' workspace with name <tt>p</tt> you can type. | ||
− | < | + | <tt>o = dpkfw.mdm.bestFilteringDistance(p,'templates',t);</tt> |
The output structure <tt>o</tt> contains several fields with different output aspects. | The output structure <tt>o</tt> contains several fields with different output aspects. | ||
Line 77: | Line 97: | ||
</nowiki> | </nowiki> | ||
− | In particular, <tt>ccmax</tt> is telling you the response of the particle to each filter. In this case, the largest response is to the third filter in <tt>t</tt>, corresponding to an space of 5 pixels. | + | In particular, <tt>ccmax</tt> is telling you the response of the particle to each filter. In this case, the largest response is to the third filter in <tt>t</tt>, corresponding to an space of 5 pixels. Note that in this particle the maximum of the cc is significantly bigger for one filter than for the others. Ideally, most of your particles should behave this way. |
− | |||
== Loop on particles and distance histogram== | == Loop on particles and distance histogram== | ||
Now you can put all these elements together and compute the largest response for each particle. | Now you can put all these elements together and compute the largest response for each particle. | ||
+ | |||
+ | |||
+ | <nowiki>for i=1:length(pac) | ||
+ | o = dpkfw.mdm.bestFilteringDistance(pa{i},'templates',t); | ||
+ | [maxcc,myIndex] = max(o.ccmax); | ||
+ | indexList(i) = myIndex(1); % in case there are two maxima | ||
+ | disp(i); | ||
+ | end | ||
+ | </nowiki> | ||
+ | Now, in the position <tt>i</tt> of numeric array <tt>indexList</tt> you have the index of the filter in <tt>t</tt> that responded best to the particle <tt>pa{i}</tt>. This can be depicted with the <tt>hist</tt> command in Matlab. | ||
+ | |||
+ | <nowiki> figure;hist(indexList,[1:5]); | ||
+ | >> xlabel('Best filter'); | ||
+ | >> ylabel('# Particles'); | ||
+ | </nowiki> | ||
+ | [[File:MembraneDistanceHistogram.png|thumb|center|300px| Histogram of best responding membrane filters]] | ||
+ | You'll need to convert each index to the actual number of pixels of the intermembrane space, and of course to physical units. |
Latest revision as of 11:09, 9 May 2017
This framework is more useful in plastic sections. The premise is that an alignment project has been created to fix the position an orientation of a set of points of interest inside a membrane. In a later stage, we want to get an estimation on the distance between the membrane layers on each particle. This task can be performed through a multireference analysis, but it can also be approached with an analysis a posteriori as sketeched here.
Contents
Reading the particles
The easiest way to bring particles aligned by a project is through the use of the dpkdata.Reader object. First store the table and data folder of your project in variables myTable and myData. This can be easily done by browsing into the folder structure, or just by using ddb to extract items from the project. Option -d will just display the file or folder location.
ddb myProject:rt -d # prints in screen the last available refined table ddb myProject:data -d % prints on screen the data folder linked to this project;
In any case,
r = dpkdata.Reader(); r.setSource('table', myTable,'data',myData);
you can access the particles throught the getParticle method:
p = r.getParticle(12,'align',1);
gets the particle with tag 12 as a variable p in the workspace. You can also read all the particles as a matlab cell array through
pa = r.getParticle('*','align',1);
Note that pa may have empty spaces: only the positions with an actual particle tag in the data folder will be occupied. For instance, if you don't have a particle labelled wiht 4 in your data folder, you will see it as.
>> pa pa = 1×194 cell array Columns 1 through 5 [96×96×96 double] [96×96×96 double] [96×96×96 double] [] [96×96×96 double] Columns 6 through 9 [96×96×96 double] [96×96×96 double] [96×96×96 double] [96×96×96 double]
You can create a compact version of this cell array:
pac = mbg.cellCompact(pa);
>> pac us = 1×193 cell array Columns 1 through 4 [96×96×96 double] [96×96×96 double] [96×96×96 double] [96×96×96 double] Columns 5 through 8 [96×96×96 double] [96×96×96 double] [96×96×96 double] [96×96×96 double]
Checking the read particles
We can only measure the membrane distance on aligned particles. At this point we make certain that everything is ok by quickly adding all the particles in the cell array and checking visually the result.
myAverage = dynamo_sumarray(pac); dview(myAverage);
Creating distance filters
We first create a set of filters that model the distance between the two layers of the membrane.
t = dpkfw.mdm.templateSeries('size',[16,16,16],'thickness',1,'distance',[3:7]);
Here, t is a cell array, containing 5 volumes of 16x16x16 pixels each. Each volume contains a "model" of a membrane. The different models are coded with the flag distance. In this case, the model number one contains a membrane which a distance between layers of 3 pixels. In the last model, the distance is 7.
You can visualize each membrane distance filter with:
figure;dslices(t,'projy','*','dim',[1,5]);
Applying the filters onto a particle
If you have your particle in the Matlab/Dynamo workspace with name p you can type.
o = dpkfw.mdm.bestFilteringDistance(p,'templates',t);
The output structure o contains several fields with different output aspects.
>> o o = output with properties: ok: 0 report: {} cc: {[1×1 struct] [1×1 struct] [1×1 struct] [1×1 struct] [1×1 struct]} ccmax: [0.3385 0.3321 0.7576 0.3305 0.3393]
In particular, ccmax is telling you the response of the particle to each filter. In this case, the largest response is to the third filter in t, corresponding to an space of 5 pixels. Note that in this particle the maximum of the cc is significantly bigger for one filter than for the others. Ideally, most of your particles should behave this way.
Loop on particles and distance histogram
Now you can put all these elements together and compute the largest response for each particle.
for i=1:length(pac) o = dpkfw.mdm.bestFilteringDistance(pa{i},'templates',t); [maxcc,myIndex] = max(o.ccmax); indexList(i) = myIndex(1); % in case there are two maxima disp(i); end
Now, in the position i of numeric array indexList you have the index of the filter in t that responded best to the particle pa{i}. This can be depicted with the hist command in Matlab.
figure;hist(indexList,[1:5]); >> xlabel('Best filter'); >> ylabel('# Particles');
You'll need to convert each index to the actual number of pixels of the intermembrane space, and of course to physical units.