Image reconstruction from a set of localizations + uncertainty

[lechristophe]

Hi,

Note : I crosspost this on both ImageJ mailing list and Fiji-users
list, as I think both could contribute to useful advices.

I'd like to implement a script (or maybe a plugin ?) for
reconstructing images from "pointillist" super-resolution microscopy
experiments.

Basically, I already use a plugin that for detection and fitting of
the individual fluorophore localizations, and I end up with a table
for each localization, the XY coordinates (with a 0.01 nm numerical
precision) + a localization uncertainty around these coordinates (
typically ~10 nm).

From this localizations table, what I would like to do is to
reconstruct an image where each localization would be plotted as a
gaussian centered on its precise coordinates, with a width
corresponding to the localization uncertainty. Of course this would be
sampled to an arbitrary pixel size (let's say 5 nm), so each gaussian
would end up as a bunch of square pixels centered on the subpixel
center of the gaussian.

Anyone knows how I could do that in ImageJ? Notably, is there a simple
way to generate a sub-pixel centered gaussian directly on the
arbitrary-sized pixel array ?

Thanks for your help,

Christophe


--
Christophe Leterrier
Researcher
Axonal Domains Architecture Team
CRN2M CNRS UMR 7286
Aix Marseille University, France

[Stephan Preibisch]

Hi Christophe,

I have been implementing a similar method some time ago, called "Volume Injection". I have to look for the source code.

This was, however, just a first step in an inverse interpolation:
"Volume Reconstruction by Inverse Interpolation: Application to Interleaved MR Motion Correction"
http://www.ncbi.nlm.nih.gov/pubmed/18979819

The problem was that just placing the gaussians makes it somehow blurry. But I never finished to implementing the whole inverse interpolation thing.

Nice greetings,
Stephan

On Mar 8, 2012, at 5:05 , Christophe Leterrier wrote:

> Hi,
>
> Note : I crosspost this on both ImageJ mailing list and Fiji-users
> list, as I think both could contribute to useful advices.
>
> I'd like to implement a script (or maybe a plugin ?) for
> reconstructing images from "pointillist" super-resolution microscopy
> experiments.
>
> Basically, I already use a plugin that for detection and fitting of
> the individual fluorophore localizations, and I end up with a table
> for each localization, the XY coordinates (with a 0.01 nm numerical
> precision) + a localization uncertainty around these coordinates (
> typically ~10 nm).
>
> From this localizations table, what I would like to do is to
> reconstruct an image where each localization would be plotted as a
> gaussian centered on its precise coordinates, with a width
> corresponding to the localization uncertainty. Of course this would be
> sampled to an arbitrary pixel size (let's say 5 nm), so each gaussian
> would end up as a bunch of square pixels centered on the subpixel
> center of the gaussian.
>
> Anyone knows how I could do that in ImageJ? Notably, is there a simple
> way to generate a sub-pixel centered gaussian directly on the
> arbitrary-sized pixel array ?
>
> Thanks for your help,
>
> Christophe
>
>
> --
> Christophe Leterrier
> Researcher
> Axonal Domains Architecture Team
> CRN2M CNRS UMR 7286
> Aix Marseille University, France

[Robert Dougherty]

Christophe,

It may be slow, but you could use the explicit formula for the Gaussian.  Suppose your fluorophores are at (X[k], Y[k]), k = 0,...K-1, where X and Y are, say, 32-bit numbers scaled to the range of 0..w-1 and 0..h-1, where w and h are the width and height of the image in pixels. Suppose further that the uncertainty in the location is u[k] and the intensity is I[k].

ImageProcessor ipOut = new FloatProcessor(w,h);
float[] pixels = ipOut.getPixels();
for(j = 0, j < h; j++){
        for(i = 0; i < w; i++){
                for(k = 0; k < K; k++){
                        float r2 = (i - X[k])*(i - X[k]) + (j - Y[k])*(j - Y[k]);
                        float R2 = 2*u[k]*u[k];
                        float ratio = r2/R2;
                        if(ratio < 10){
                                pixels[i + w*j] += (float)(I[k]*Math.exp(-ratio));
                        }
                }
        }
}

Then you can convert ipOut to some other type, if you so desire, and display it.

If you also have fluorophore intensities, I[k], then you would multiply the exponential by I[k] before or after converting it to float.

Bob


On Mar 8, 2012, at 2:05 AM, Christophe Leterrier wrote:

> Hi,
>
> Note : I crosspost this on both ImageJ mailing list and Fiji-users
> list, as I think both could contribute to useful advices.
>
> I'd like to implement a script (or maybe a plugin ?) for
> reconstructing images from "pointillist" super-resolution microscopy
> experiments.
>
> Basically, I already use a plugin that for detection and fitting of
> the individual fluorophore localizations, and I end up with a table
> for each localization, the XY coordinates (with a 0.01 nm numerical
> precision) + a localization uncertainty around these coordinates (
> typically ~10 nm).
>
> From this localizations table, what I would like to do is to
> reconstruct an image where each localization would be plotted as a
> gaussian centered on its precise coordinates, with a width
> corresponding to the localization uncertainty. Of course this would be
> sampled to an arbitrary pixel size (let's say 5 nm), so each gaussian
> would end up as a bunch of square pixels centered on the subpixel
> center of the gaussian.
>
> Anyone knows how I could do that in ImageJ? Notably, is there a simple
> way to generate a sub-pixel centered gaussian directly on the
> arbitrary-sized pixel array ?
>
> Thanks for your help,
>
> Christophe
>
>
> --
> Christophe Leterrier
> Researcher
> Axonal Domains Architecture Team
> CRN2M CNRS UMR 7286
> Aix Marseille University, France

Robert Dougherty, Ph.D.
President, OptiNav, Inc.
1414 127th Place NE #106
Bellevue, WA 98005
[phone]
FAX [phone]
www.optinav.com
[hidden email]

[lechristophe]

Hi Bob,

Thanks a lot for this suggestion. Actually it helped me understand the
FloatProcessor class (I'm really learning the IJ API step by step !).
So it allows to manipulate images with no pixel sampling, and then
sample them by converting them to pixellated images? Maybe it could be
also possible to use the GaussianBlur.blurFloat() method to generate
individual gaussian sub-images and adding them to create the "big"
Floatprocessor where the reconstruction is done ?

Christophe

On Thu, Mar 8, 2012 at 17:41, Robert Dougherty <[hidden email]> wrote:

> Christophe,
>
> It may be slow, but you could use the explicit formula for the Gaussian.  Suppose your fluorophores are at (X[k], Y[k]), k = 0,...K-1, where X and Y are, say, 32-bit numbers scaled to the range of 0..w-1 and 0..h-1, where w and h are the width and height of the image in pixels. Suppose further that the uncertainty in the location is u[k] and the intensity is I[k].
>
> ImageProcessor ipOut = new FloatProcessor(w,h);
> float[] pixels = ipOut.getPixels();
> for(j = 0, j < h; j++){
>        for(i = 0; i < w; i++){
>                for(k = 0; k < K; k++){
>                        float r2 = (i - X[k])*(i - X[k]) + (j - Y[k])*(j - Y[k]);
>                        float R2 = 2*u[k]*u[k];
>                        float ratio = r2/R2;
>                        if(ratio < 10){
>                                pixels[i + w*j] += (float)(I[k]*Math.exp(-ratio));
>                        }
>                }
>        }
> }
>
> Then you can convert ipOut to some other type, if you so desire, and display it.
>
> If you also have fluorophore intensities, I[k], then you would multiply the exponential by I[k] before or after converting it to float.
>
> Bob
>
>
> On Mar 8, 2012, at 2:05 AM, Christophe Leterrier wrote:
>
>> Hi,
>>
>> Note : I crosspost this on both ImageJ mailing list and Fiji-users
>> list, as I think both could contribute to useful advices.
>>
>> I'd like to implement a script (or maybe a plugin ?) for
>> reconstructing images from "pointillist" super-resolution microscopy
>> experiments.
>>
>> Basically, I already use a plugin that for detection and fitting of
>> the individual fluorophore localizations, and I end up with a table
>> for each localization, the XY coordinates (with a 0.01 nm numerical
>> precision) + a localization uncertainty around these coordinates (
>> typically ~10 nm).
>>
>> From this localizations table, what I would like to do is to
>> reconstruct an image where each localization would be plotted as a
>> gaussian centered on its precise coordinates, with a width
>> corresponding to the localization uncertainty. Of course this would be
>> sampled to an arbitrary pixel size (let's say 5 nm), so each gaussian
>> would end up as a bunch of square pixels centered on the subpixel
>> center of the gaussian.
>>
>> Anyone knows how I could do that in ImageJ? Notably, is there a simple
>> way to generate a sub-pixel centered gaussian directly on the
>> arbitrary-sized pixel array ?
>>
>> Thanks for your help,
>>
>> Christophe
>>
>>
>> --
>> Christophe Leterrier
>> Researcher
>> Axonal Domains Architecture Team
>> CRN2M CNRS UMR 7286
>> Aix Marseille University, France
>
> Robert Dougherty, Ph.D.
> President, OptiNav, Inc.
> 1414 127th Place NE #106
> Bellevue, WA 98005
> [phone]
> FAX [phone]
> www.optinav.com
> [hidden email]

[Robert Dougherty]

Christophe,

Here is a more or less fully formed example code. Creating the FloatProcessor allocates its pixel array.  There are other ways to do it.  For example, you can just allocate and compute the pixel array  and later use it in a different constructor for a FloatProcessor.  I'm not sure how you could use a GaussianBlur, since the image you would blur with it would be forced to have values only at the pixel locations, not in between.  You could make a scaled-up blur to improve the resolution and then scale it down before adding to the "big" image, but this would be somewhat complicated.  You could also make the big image even bigger, put in the spots, blur that, and then scale down, but then you would not have the flexibility of applying different uncertainties to the different spots.

Bob



import ij.*;
import ij.plugin.PlugIn;
import ij.process.*;
import ij.gui.*;
import java.awt.*;

/** Make Gaussian Spots
 */
public class Gaussian_Points implements PlugIn{
        //Get parameters and save them in the preferences for the next run.
        public void run(String arg) {
                int w = (int)Prefs.get("gaussian_points.w", 300);
                int h = (int)Prefs.get("gaussian_points.h", 500);
                float X1 = (float)Prefs.get("gaussian_points.X1", 20.5);
                float Y1 = (float)Prefs.get("gaussian_points.Y1", 30.2);
                float u1 = (float)Prefs.get("gaussian_points.u1", 3.5);
                float I1 = (float)Prefs.get("gaussian_points.I1", 52);
                float X2 = (float)Prefs.get("gaussian_points.X2", 236.5);
                float Y2 = (float)Prefs.get("gaussian_points.Y2", 347.2);
                float u2 = (float)Prefs.get("gaussian_points.u2", 7.2);
                float I2 = (float)Prefs.get("gaussian_points.I2", 50);
                GenericDialog gd = new GenericDialog("Gaussian_Points...", IJ.getInstance());
                gd.addNumericField("w image width, pixels", w, 0);
                gd.addNumericField("h image height, pixels", h, 0);
                gd.addMessage("First Spot");
                gd.addNumericField("X1",X1,3);
                gd.addNumericField("Y1",Y1,3);
                gd.addNumericField("u1 uncertainty",u1,2);
                gd.addNumericField("I1 intensity",I1,2);
                gd.addMessage("Second Spot");
                gd.addNumericField("X2",X2,3);
                gd.addNumericField("Y2",Y2,3);
                gd.addNumericField("u2 uncertainty",u2,2);
                gd.addNumericField("I2 intensity",I2,2);
                gd.showDialog();
                if (gd.wasCanceled()) {
                        return;
                }
                w = (int)gd.getNextNumber();
                h = (int)gd.getNextNumber();
                X1 = (float)gd.getNextNumber();
                Y1 = (float)gd.getNextNumber();
                u1 = (float)gd.getNextNumber();
                I1 = (float)gd.getNextNumber();
                X2 = (float)gd.getNextNumber();
                Y2 = (float)gd.getNextNumber();
                u2 = (float)gd.getNextNumber();
                I2 = (float)gd.getNextNumber();
                Prefs.set("gaussian_points.w", w);
                Prefs.set("gaussian_points.h", h);
                Prefs.set("gaussian_points.X1", X1);
                Prefs.set("gaussian_points.Y1", Y1);
                Prefs.set("gaussian_points.u1", u1);
                Prefs.set("gaussian_points.I1", I1);
                Prefs.set("gaussian_points.h", h);
                Prefs.set("gaussian_points.X2", X2);
                Prefs.set("gaussian_points.Y2", Y2);
                Prefs.set("gaussian_points.u2", u2);
                Prefs.set("gaussian_points.I2", I2);
                //Make arrays
                float[] X = new float[]{X1,X2};
                float[] Y = new float[]{Y1,Y2};
                float[] u = new float[]{u1,u2};
                float[] I = new float[]{I1,I2};
                int K = 2;
                //Create FloatProcessor.  This allocates the pixel array.
                ImageProcessor ipOut = new FloatProcessor(w,h);
                //Get the pixel array to abel to change it.
                float[] pixels = (float[])ipOut.getPixels();
                //Paint the spots in the pixel array
                for(int j = 0; j < h; j++){
                        for(int i = 0; i < w; i++){
                                for(int k = 0; k < K; k++){
                                        float r2 = (i - X[k])*(i - X[k]) + (j - Y[k])*(j - Y[k]);
                                        float R2 = 2*u[k]*u[k];
                                        float ratio = r2/R2;
                                        if(ratio < 10){
                                                pixels[i + w*j] += (float)(I[k]*Math.exp(-ratio));
                                        }
                                }
                        }
                }
                //Set the brightness and contrat limits to the full range of values in the ImageProcessor
                ipOut.setMinAndMax(0,0);
                //Create an ImagePlus and show it.
                ImagePlus impOut = new ImagePlus("Spots",ipOut);
                impOut.show();
                //Chage the lookup table
                IJ.run("Fire");
        }
}

On Mar 8, 2012, at 1:17 PM, Christophe Leterrier wrote:

> Hi Bob,
>
> Thanks a lot for this suggestion. Actually it helped me understand the
> FloatProcessor class (I'm really learning the IJ API step by step !).
> So it allows to manipulate images with no pixel sampling, and then
> sample them by converting them to pixellated images? Maybe it could be
> also possible to use the GaussianBlur.blurFloat() method to generate
> individual gaussian sub-images and adding them to create the "big"
> Floatprocessor where the reconstruction is done ?
>
> Christophe
>
> On Thu, Mar 8, 2012 at 17:41, Robert Dougherty <[hidden email]> wrote:
>> Christophe,
>>
>> It may be slow, but you could use the explicit formula for the Gaussian.  Suppose your fluorophores are at (X[k], Y[k]), k = 0,...K-1, where X and Y are, say, 32-bit numbers scaled to the range of 0..w-1 and 0..h-1, where w and h are the width and height of the image in pixels. Suppose further that the uncertainty in the location is u[k] and the intensity is I[k].
>>
>> ImageProcessor ipOut = new FloatProcessor(w,h);
>> float[] pixels = ipOut.getPixels();
>> for(j = 0, j < h; j++){
>>        for(i = 0; i < w; i++){
>>                for(k = 0; k < K; k++){
>>                        float r2 = (i - X[k])*(i - X[k]) + (j - Y[k])*(j - Y[k]);
>>                        float R2 = 2*u[k]*u[k];
>>                        float ratio = r2/R2;
>>                        if(ratio < 10){
>>                                pixels[i + w*j] += (float)(I[k]*Math.exp(-ratio));
>>                        }
>>                }
>>        }
>> }
>>

Robert Dougherty, Ph.D.
President, OptiNav, Inc.
1414 127th Place NE #106
Bellevue, WA 98005
Tel. [phone]
FAX [phone]
www.optinav.com
[hidden email]

[lechristophe]

Hi Robert and IJers,

Thanks to Robert's help, I got the gaussian-localizations image
reconstruction script working.

Here is the relevant code part. The generateImage2 function takes as
arguments the square image size s, and an array Loc containing  an X
(x coordinates of localizations), Y (y coordinates of localizations)
and U (sigma of localization uncertainty) arrays, then returns the
generated float image ImageProcessor

function generateImage2(s, Loc) {
       var ipOut = new FloatProcessor(s,s);
       var X = Loc[0];
       var Y = Loc[1];
       var U = Loc[2];
       var K = X.length;
       pixels = ipOut.getPixels();
       for (var k = 0; k < K; k++){
               jmin = Math.max(0, Math.round(Y[k] - 4 * U[k]));
               jmax = Math.min(s-1, Math.round(Y[k] + 4 * U[k]));
               imin = Math.max(0, Math.round(X[k] - 4 * U[k]));
               imax = Math.min(s-1, Math.round(X[k] + 4 * U[k]));
               for(var j = jmin; j <= jmax; j++){
                      for(var i = imin; i <= imax; i++){
                              r2 = (i + 0.5 - X[k]) * (i + 0.5 - X[k]) + (j
+ 0.5 - Y[k]) * (j + 0.5 - Y[k]);
                              R2 = 2 * U[k] * U[k];
                              ratio = r2 / R2;
                              if(ratio < 10){
                                      pixels[i + s*j] += 1000 *
Math.exp(-ratio);
                              }
                       }
               }
       IJ.showProgress(k, K);
       }
       ipOut.setPixels(pixels);
       return ipOut;
}


Thanks to Robert Dougherty for his help. If you see any possible
improvement to the code, just let me know!

Christophe

On Thu, Mar 8, 2012 at 17:41, Robert Dougherty <[hidden email]> wrote:

> Christophe,
>
> It may be slow, but you could use the explicit formula for the Gaussian.  Suppose your fluorophores are at (X[k], Y[k]), k = 0,...K-1, where X and Y are, say, 32-bit numbers scaled to the range of 0..w-1 and 0..h-1, where w and h are the width and height of the image in pixels. Suppose further that the uncertainty in the location is u[k] and the intensity is I[k].
>
> ImageProcessor ipOut = new FloatProcessor(w,h);
> float[] pixels = ipOut.getPixels();
> for(j = 0, j < h; j++){
>        for(i = 0; i < w; i++){
>                for(k = 0; k < K; k++){
>                        float r2 = (i - X[k])*(i - X[k]) + (j - Y[k])*(j - Y[k]);
>                        float R2 = 2*u[k]*u[k];
>                        float ratio = r2/R2;
>                        if(ratio < 10){
>                                pixels[i + w*j] += (float)(I[k]*Math.exp(-ratio));
>                        }
>                }
>        }
> }
>
> Then you can convert ipOut to some other type, if you so desire, and display it.
>
> If you also have fluorophore intensities, I[k], then you would multiply the exponential by I[k] before or after converting it to float.
>
> Bob
>
>
> On Mar 8, 2012, at 2:05 AM, Christophe Leterrier wrote:
>
>> Hi,
>>
>> Note : I crosspost this on both ImageJ mailing list and Fiji-users
>> list, as I think both could contribute to useful advices.
>>
>> I'd like to implement a script (or maybe a plugin ?) for
>> reconstructing images from "pointillist" super-resolution microscopy
>> experiments.
>>
>> Basically, I already use a plugin that for detection and fitting of
>> the individual fluorophore localizations, and I end up with a table
>> for each localization, the XY coordinates (with a 0.01 nm numerical
>> precision) + a localization uncertainty around these coordinates (
>> typically ~10 nm).
>>
>> From this localizations table, what I would like to do is to
>> reconstruct an image where each localization would be plotted as a
>> gaussian centered on its precise coordinates, with a width
>> corresponding to the localization uncertainty. Of course this would be
>> sampled to an arbitrary pixel size (let's say 5 nm), so each gaussian
>> would end up as a bunch of square pixels centered on the subpixel
>> center of the gaussian.
>>
>> Anyone knows how I could do that in ImageJ? Notably, is there a simple
>> way to generate a sub-pixel centered gaussian directly on the
>> arbitrary-sized pixel array ?
>>
>> Thanks for your help,
>>
>> Christophe
>>
>>
>> --
>> Christophe Leterrier
>> Researcher
>> Axonal Domains Architecture Team
>> CRN2M CNRS UMR 7286
>> Aix Marseille University, France
>
> Robert Dougherty, Ph.D.
> President, OptiNav, Inc.
> 1414 127th Place NE #106
> Bellevue, WA 98005
> [phone]
> FAX [phone]
> www.optinav.com
> [hidden email]