# Patent

## Title
Methods and systems for image data processing

## Abstract

Embodiments of the computer-implemented methods, storage mediums, and systems may be configured to determine locations of particles within a first image of the particles. The particles may have fluorescence-material associated therewith. The embodiments may include calculating a transform parameter, and the transform parameter may define an estimated movement in the locations of the particles between the first image of the particles and a second image of the particles. The embodiments may further including applying the transform parameter to the locations of the particles within the first image to determine movement locations of the particles within the second image.

## Background
PRIORITY CLAIM
This application claims the benefit of priority to U.S. Provisional Patent Application 61/627,851 entitled “Methods and Systems for Image Data Processing” filed Oct. 18, 2011, which is hereby incorporated by reference in its entirety.


BACKGROUND OF THE INVENTION
1. Field of the Invention
This invention generally relates to methods and systems for image data processing. Certain embodiments relate to methods and systems for performing one or more steps for processing multiple images of particles to account for movement of the particles between the images.
2. Description of the Related Art
The following descriptions and examples are not admitted to be prior art by virtue of their inclusion within this section.
Imaging using detectors such as charged coupled device (CCD) detectors is employed in several currently available instruments in biotechnology applications. Such applications may require taking multiple images of particles. In these multiple images of the particles, the particles may appear to move. In images of particles taken close together in time (or perhaps at the same time), the particles may appear to shift or move. Accordingly, it would be desirable to develop methods and systems for data processing of images of particles to account for movement of particles between images.
SUMMARY OF THE INVENTION
The problem outlined above may be in large part addressed by computer-implemented methods, storage mediums, and systems for performing one or more steps associated with data image processing of particles. The following are mere exemplary embodiments of the computer-implemented methods, storage mediums, and systems and are not to be construed in any way to limit the subject matter of the claims.
Embodiments of the computer-implemented methods, storage mediums, and systems may be configured to determine locations of particles within a first image of the particles, wherein the particles have fluorescence-material associated therewith; calculate a transform parameter, wherein the transform parameter defines an estimated movement in the locations of the particles between the first image of the particles and a second image of the particles; and apply the transform parameter to the locations of the particles within the first image to determine movement locations of the particles within the second image.
In some embodiments of the methods, storage mediums, and systems, the transform parameter includes a radial component and a constant component. The radial component may be proportional to a distance between the locations of particles within the first image and a center of the first image.
In some embodiments of the methods, storage mediums, and systems, calculating the transform parameter may include estimating estimated locations of particles within the second image and calculating potential transform parameters based on the locations of the particles within the first image and the estimated locations of the particles in the second image. In some embodiments, calculating the transform parameter may further include determining an optimal transform parameter based on the potential transform parameters.
In some embodiments of the methods, storage mediums, and systems, estimating the estimated locations of the particles within the second image may include determining maximal integral locations based on the second image and the locations of the particles within the first image. Moreover, in some embodiments, determining the optimal transform parameter comprises using a Hough transform.
Some embodiments of the methods, storage mediums, and systems, may further include calculating an error component based on a force between the particles.
The terms “a” and “an” are defined as one or more unless this disclosure explicitly requires otherwise.
The term “substantially” and its variations are defined as being largely but not necessarily wholly what is specified as understood by one of ordinary skill in the art, and in one non-limiting embodiment “substantially” refers to ranges within 10%, preferably within 5%, more preferably within 1%, and most preferably within 0.5% of what is specified.
The terms “comprise” (and any form of comprise, such as “comprises” and “comprising”), “have” (and any form of have, such as “has” and “having”), “include” (and any form of include, such as “includes” and “including”) and “contain” (and any form of contain, such as “contains” and “containing”) are open-ended linking verbs. As a result, a method or device that “comprises,” “has,” “includes” or “contains” one or more steps or elements possesses those one or more steps or elements, but is not limited to possessing only those one or more elements. Likewise, a step of a method or an element of a device that “comprises,” “has,” “includes” or “contains” one or more features possesses those one or more features, but is not limited to possessing only those one or more features. Furthermore, a device or structure that is configured in a certain way is configured in at least that way, but may also be configured in ways that are not listed.
Other features and associated advantages will become apparent with reference to the following detailed description of specific embodiments in connection with the accompanying drawings.



BRIEF DESCRIPTION OF THE DRAWINGS
The following drawings form part of the present specification and are included to further demonstrate certain aspects of the present invention. The invention may be better understood by reference to one or more of these drawings in combination with the detailed description of specific embodiments presented herein.
FIG. 1 is a schematic diagram illustrating a cross-sectional view of one embodiment of a system configured to acquire and process images of particles;
FIG. 2 is a flowchart outlining a method for processing images of particles;
FIGS. 3A-C illustrate potential grid distortion between a first image of particles and a second image of particles;
FIG. 4 illustrates a free body diagram of a particle;
FIGS. 5A-5B illustrate a specific embodiment of determining a maximal integral location;
FIGS. 6A-6D illustrate a vote space used in determining an optimal transform parameter; and
FIG. 7 illustrates an additional free body diagram of a particle.



DETAILED DESCRIPTION
Although embodiments are described herein with respect to particles, it is to be understood that the systems and methods described herein may also be used with microspheres, polystyrene beads, microparticles, gold nanoparticles, quantum dots, nanodots, nanoparticles, nanoshells, beads, microbeads, latex particles, latex beads, fluorescent beads, fluorescent particles, colored particles, colored beads, tissue, cells, micro-organisms, organic matter, non-organic matter, or any other discrete substances known in the art. The particles may serve as vehicles for molecular reactions. Examples of appropriate particles are illustrated and described in U.S. Pat. No. 5,736,330 to Fulton, U.S. Pat. No. 5,981,180 to Chandler et al., U.S. Pat. No. 6,057,107 to Fulton, U.S. Pat. No. 6,268,222 to Chandler et al., U.S. Pat. No. 6,449,562 to Chandler et al., U.S. Pat. No. 6,514,295 to Chandler et al., U.S. Pat. No. 6,524,793 to Chandler et al., and U.S. Pat. No. 6,528,165 to Chandler, which are incorporated by reference as if fully set forth herein. The systems and methods described herein may be used with any of the particles described in these patents. In addition, particles for use in method and system embodiments described herein may be obtained from manufacturers such as Luminex Corporation of Austin, Tex. The terms “particles” and “microspheres” and “beads” are used interchangeably herein.
In addition, the types of particles that are compatible with the systems and methods described herein include particles with fluorescent materials attached to, or associated with, the surface of the particles. These types of particles, in which fluorescent dyes or fluorescent particles are coupled directly to the surface of the particles in order to provide the classification fluorescence (i.e., fluorescence emission measured and used for determining an identity of a particle or the subset to which a particle belongs), are illustrated and described in U.S. Pat. No. 6,268,222 to Chandler et al. and U.S. Pat. No. 6,649,414 to Chandler et al., which are incorporated by reference as if fully set forth herein. The types of particles that can be used in the methods and systems described herein also include particles having one or more fluorochromes or fluorescent dyes incorporated into the core of the particles.
Particles that can be used in the methods and systems described herein further include particles that in of themselves will exhibit one or more fluorescent signals upon exposure to one or more appropriate light sources. Furthermore, particles may be manufactured such that upon excitation the particles exhibit multiple fluorescent signals, each of which may be used separately or in combination to determine an identity of the particles. As described below, image data processing may include classification of the particles, particularly for a multi-analyte fluid, as well as a determination of the amount of analyte bound to the particles. Since a reporter signal, which may represent the amount of analyte bound to the particle, is typically unknown during operations, specially dyed particles, which not only emit fluorescence in the classification wavelength(s) or wavelength band(s) but also in the reporter wavelength or wavelength band, may be used in conjunction with the systems described herein.
The methods described herein generally include analyzing images of particles and processing data measured from the images to determine the location of the particles within the images. Subsequent processing of the one or more characteristics of the particles can be performed according to the methods described in U.S. Pat. No. 5,736,330 to Fulton, U.S. Pat. No. 5,981,180 to Chandler et al., U.S. Pat. No. 6,449,562 to Chandler et al., U.S. Pat. No. 6,524,793 to Chandler et al., U.S. Pat. No. 6,592,822 to Chandler, and U.S. Pat. No. 6,939,720 to Chandler et al. as well as U.S. patent application Ser. No. 11/534,166 to Roth et al. which are incorporated by reference as if fully set forth herein.
Turning now to the drawings, it is noted that FIG. 1 is not drawn to scale. In particular, the scale of some of the elements of the figure is greatly exaggerated to emphasize characteristics of the elements. Some elements of the system have not been included in the figures for the sake of clarity.
One embodiment of a system configured to generate, acquire, or supply images of particles and to process the images according to embodiments of methods described herein is shown in FIG. 1. The system shown in FIG. 1 may be used in applications such as multi-analyte measurement of particles. The system includes an imaging subsystem that includes light source 10. Light source 10 may include one or more light sources such as light emitting diodes (LED), lasers, arc lamps, incandescent lamps, or any other suitable light sources known in the art. In addition, or alternatively, the light source may include more than one light source (not shown), each of which is configured to generate light a different wavelength or a different wavelength band. One example of an appropriate combination of light sources for use in the system shown in FIG. 1 includes, but is not limited to, two or more LEDs. Light from more than one light source may be combined into a common illumination path by a beam splitter (not shown) or any other suitable optical element known in the art such that light from the light sources may be directed to the particles simultaneously. Alternatively, the imaging subsystem may include an optical element (not shown) such as a reflecting mirror and a device (not shown) configured to move the optical element into and out of the illumination path depending on which light source is used to illuminate the particles. In this manner, the light sources may be used to sequentially illuminate the particles with different wavelengths or wavelength bands of light. The light source(s) may also illuminate the substrate from above, rather than below the substrate (not shown).
The light source(s) may be selected to provide light at wavelength(s) or wavelength band(s) that will cause the particles or material coupled thereto to emit fluorescence. For instance, the wavelength(s) or wavelength band(s) may be selected to excite fluorescent dyes or other fluorescent materials incorporated into the particles and/or coupled to a surface of the particles. In this manner, the wavelength(s) or wavelength band(s) may be selected such that the particles emit fluorescence that is used for classification of the particles. In addition, the wavelength(s) or wavelength band(s) may be selected to excite fluorescent dyes or other fluorescent materials coupled to the particles via a reagent on the surface of the particles. As such, the wavelength(s) or wavelength band(s) may be selected such that the particles emit fluorescence that is used to detect and/or quantify reaction(s) that have taken place on the surface of the particles.
As shown in FIG. 1, the imaging subsystem may include optical element 12 that is configured to direct light from light source 10 to substrate 14 on which particles 16 are immobilized. In one example, optical element 12 may be a collimating lens. However, optical element 12 may include any other appropriate optical element that can be used to image light from light source 10 onto substrate 14. In addition, although the optical element is shown in FIG. 1 as a single optical element, it is to be understood that optical element 12 may include more than one refractive element. Furthermore, although optical element 12 is shown in FIG. 1 as a refractive optical element, it is to be understood that one or more reflective optical elements may be used (possibly in combination with one or more refractive optical elements) to image light from light source 10 onto substrate 14.
Particles 16 may include any of the particles described above. Substrate 14 may include any appropriate substrate known in the art. The particles immobilized on substrate 14 may be disposed in an imaging chamber (not shown) or any other device for maintaining a position of substrate 14 and particles 16 immobilized thereon with respect to the imaging subsystem. The device for maintaining a position of substrate 14 may also be configured to alter a position of the substrate (e.g., to focus the imaging subsystem onto the substrate) prior to imaging. Immobilization of the particles on the substrate may be performed using magnetic attraction, a vacuum filter plate, or any other appropriate method known in the art. Examples of methods and systems for positioning microspheres for imaging are illustrated in U.S. patent application Ser. No. 11/270,786 to Pempsell filed Nov. 9, 2005, which is incorporated by reference as if fully set forth herein. The particle immobilization method itself is not particularly important to the method and systems described herein. However, the particles are preferably immobilized such that the particles do not move perceptibly during the detector integration period, which may be multiple seconds long.
As shown in FIG. 1, the imaging subsystem may include optical element 18 and beam splitter 20. Optical element 18 is configured to focus light from substrate 14 and particles 16 immobilized thereon to beam splitter 20. Optical element 18 may be further configured as described above with respect to optical element 12. Beam splitter 20 may include any appropriate beam splitter known in the art. Beam splitter 20 may be configured to direct light from optical element 18 to different detectors based on the wavelength of the light. For example, light having a first wavelength or wavelength band may be transmitted by beam splitter 20, and light having a second wavelength or wavelength band different than the first may be reflected by beam splitter 20. The imaging subsystem may also include optical element 22 and detector 24. Light transmitted by beam splitter 20 may be directed to optical element 22. Optical element 22 is configured to focus the light transmitted by the beam splitter onto detector 24. The imaging subsystem may further include optical element 26 and detector 28. Light reflected by beam splitter 20 may be directed to optical element 26. Optical element 26 is configured to focus the light reflected by the beam splitter onto detector 28. Optical elements 22 and 26 may be configured as described above with respect to optical element 12.
Detectors 24 and 28 may include, for example, charge coupled device (CCD) detectors or any other suitable imaging detectors known in the art such as CMOS detectors, two-dimensional arrays of photosensitive elements, time delay integration (TDI) detectors, etc. In some embodiments, a detector such as a two-dimensional CCD imaging array may be used to acquire an image of substantially an entire substrate or of all particles immobilized on a substrate simultaneously. In this manner, all photons from the illuminated area of the substrate may be collected simultaneously thereby eliminating error due to a sampling aperture used in other currently available systems that include a photomultiplier tube (PMT) and scanning device. In addition, the number of detectors included in the system may be equal to the number of wavelengths or wavelength bands of interest such that each detector is used to generate images at one of the wavelengths or wavelength bands.
Each of the images generated by the detectors may be spectrally filtered using an optical bandpass element (not shown) or any other suitable optical element known in the art, which is disposed in the light path from the beam splitter to the detectors. A different filter “band” may be used for each captured image. The detection wavelength center and width for each wavelength or wavelength band at which an image is acquired may be matched to the fluorescent emission of interest, whether it is used for particle classification or the reporter signal. In this manner, the imaging subsystem of the system shown in FIG. 1 may be configured to generate multiple images at different wavelengths or wavelength bands simultaneously. Although the system shown in FIG. 1 includes two detectors, it is to be understood that the system may include more than two detectors (e.g., three detectors, four detectors, etc.). As described above, each of the detectors may be configured to generate images at different wavelengths or wavelength bands simultaneously and/or sequentially by including one or more optical elements for directing light at different wavelengths or wavelength bands to the different detectors simultaneously and/or sequentially.
In addition, although the system is shown in FIG. 1 to include multiple detectors, it is to be understood that the system may include a single detector. The single detector may be used to generate multiple images at multiple wavelengths or wavelength bands sequentially. For example, light of different wavelengths or wavelength bands may be directed to the substrate sequentially, and different images may be generated during illumination of the substrate with each of the different wavelengths or wavelength bands. In another example, different filters for selecting the wavelength or wavelength bands of light directed to the single detector may be altered (e.g., by moving the different filters into and out of the imaging path) to generate images at different wavelengths or wavelength bands sequentially.
The imaging subsystem shown in FIG. 1, therefore, is configured to generate a plurality or series of images representing the fluorescent emission of particles 16 at several wavelengths of interest. In addition, the system may be configured to supply a plurality or series of digital images representing the fluorescence emission of the particles to a processor (i.e., a processing engine). In one such example, the system may include processor 30. Processor 30 may be configured to acquire (e.g., receive) image data from detectors 24 and 28. For example, processor 30 may be coupled to detectors 24 and 28 in any suitable manner known in the art (e.g., via transmission media (not shown), each coupling one of the detectors to the processor, via one or more electronic components (not shown) such as analog-to-digital converters, each coupled between one of the detectors and the processor, etc.). Preferably, processor 30 is at least configured to process and analyze these images to determine one or more characteristics of particles 16 such as a classification of the particles and information about a reaction taken place on the surface of the particles. The one or more characteristics may be output by the processor in any suitable format such as a data array with an entry for fluorescent magnitude for each particle for each wavelength. Specifically, the processor may also (e.g. alternatively or additionally) be configured to perform one or more steps of the method embodiments described herein to process and analyze the images.
Processor 30 may be a processor such as those commonly included in a typical personal computer, mainframe computer system, workstation, etc. In general, the term “computer system” may be broadly defined to encompass any device having one or more processors, which executes instructions from a memory medium. The processor may be implemented using any other appropriate functional hardware. For example, the processor may include a digital signal processor (DSP) with a fixed program in firmware, a field programmable gate array (FPGA), or other programmable logic device (PLD) employing sequential logic “written” in a high level programming language such as very high speed integrated circuits (VHSIC) hardware description language (VHDL). In another example, program instructions (not shown) executable on processor 30 to perform one or more steps of the computer-implemented methods described herein may be coded in a high level language such as C#, with sections in C++ as appropriate, ActiveX controls, JavaBeans, Microsoft Foundation Classes (“MFC”), or other technologies or methodologies, as desired. The program instructions may be implemented in any of various ways, including procedure-based techniques, component-based techniques, and/or object-oriented techniques, among others.
Program instructions implementing methods such as those described herein may be transmitted over or stored on a storage medium. The storage medium may include but is not limited to a read-only memory, a random access memory, a magnetic or optical disk, or a magnetic tape. For each image, all located particles and the values and/or statistics determined for each identified particle may be stored in a memory medium within the storage medium. The image processing methods described herein may be performed using one or more algorithms. As described in more detail below, the algorithms may be complex and, therefore, may be best implemented through a computer (e.g., processing device). As such, the methods described herein and particularly in reference to FIG. 2 may be referred to as “computer-implemented methods” and, thus, the terms “method” and “computer-implements method” may be used interchangeably herein. It is noted that the computer-implemented methods and program instructions of the systems described herein may, in some cases, be configured to perform processes other than those associated with methods described herein and, therefore, the computer-implemented methods and program instructions of systems described herein are not necessarily limited to the depiction in the figures.
According to one embodiment, a computer-implemented method for image data processing includes one or more of the following steps (i.e., high level operations): determining locations of particles within a first image, calculating a transform parameter, and applying the transform parameter to the locations of the particles within a first image to determine movement locations of the particles within the second image. In some embodiments, these steps may be performed sequentially in the order listed above.
FIG. 2 illustrates a method 200 illustrating an exemplary sequence of steps for image processing. As shown, the method 200 may include acquiring 202 data for multiple images of particles having fluorescence-material associated therewith, wherein each of the multiple images corresponds to a different wavelength band. In some cases, the data may be acquired directly from an imaging system, but in other cases, the data may be acquired from a storage medium. In either case, the data may be representative of multiple images taken at different wavelengths as noted above. Exemplary wavelengths that may be used may correspond to different color channels, such as but not limited to red for classification channel 1, red (e.g., a same or different wavelength of red) for classification channel 2, green for the reporter channel. As further noted above, in order to accommodate each color channel, the particles used for the method described herein may be specially dyed to emit at all wavelengths or in all wavelength bands of interest. In particular, in order to measure both classification and reporter signals within the multiple images, the methods described herein may be performed using specially dyed particles, which not only emit fluorescence in the classification wavelength(s) or wavelength band(s), but also in the reporter wavelength or wavelength band.
In some embodiments, an image taken a first bandwidth may be referred to as a “first image,” and a subsequent/simultaneous image taken at the same or a different bandwidth may be referred to as a “second image.” In preferred embodiments, the first image may be related to a classification channel (e.g., CL1 or CL2), and the second image may be related to a reporter channel (RP). As described above, the first image and the second image may be taken successively (in any order) or simultaneously.
In embodiments of the subsystem described with respect to FIG. 1, the particles, when viewed across one or more images, may appear to move. The movement of the particles between images may at least be caused by lens distortion and/or chromatic aberration. That is even though particles may not actually shift or move between multiple images they may appear to move based on lens distortion and or chromatic aberration. With respect to FIG. 1, the lens distortion and/or chromatic aberration may be a result of the filters (e.g., a filter wheel), light source 10, and/or one or more of optical elements 12, 18, 22 and 26. FIGS. 3A, 3B, and 3C demonstrate a movement that may be caused by grid distortion. FIG. 3A illustrates the grid distortion in a classification channel image (e.g., a first image) due to the lens distortion. As shown in the figure, the maximum distortion in the image occurs at the corners. FIG. 3B illustrates the grid distortion in a reporter channel image (e.g., a second image). As shown in this figure, the maximum distortion in the image also occurs at the corners. FIG. 3C illustrates an overlay of one corner of the reporter channel image and the classification channel image. As shown, the FIG. 3C illustrates a apparent movement towards the center of the images from the classification channel image to the reporter channel image. Thus, as a result of the lens distortion and chromatic aberration, a particle may appear to move between the classification channel image and the reporter channel image. As shown in FIG. 3, the lens distortion is a contributor to the radial movement. Light passing through any glass (e.g., optical elements, the chamber, and the like) may refract different wavelengths differently—like a prism. Variations on the chamber (e.g., the top plate) or other optical elements may cause chromatic aberration to vary as well.
After acquiring 202 data for multiple (e.g., at least two) images, the method 200 may proceed by determining the locations of particles within a first image of the particles. As discussed throughout, the particles may have florescence-material associated therewith. Moreover, in embodiments, the “first image” may refer specifically to a classification channel image. One of skill in the art will recognize a variety of image processing techniques to determine locations of particles within a classification channel image including peak searching and similar methods. For example, a variety of methods are discussed in U.S. patent application Ser. No. 11/534,166 to Roth et al.
In some embodiments, determining the location of particles within a classification channel image may be easier than determining the location of particles within a reporter channel image. As described above, the classification channel image may be configured to illuminate the particles themselves while the reporter channel image may be configured to illuminate the substances (e.g., analyte) bound to the particle. As such, a peak search (or similar algorithm) in a classification channel image may closely reveal the location of a particle. At the same time, a peak search (or similar algorithm) may reveal the location of the analyte—which may or may not correlate to the location of the particle. Rather, in some instances, such an algorithm may reveal the location of the analyte on the edge of the particle or even a different particle.
In some embodiments of the method 200, the method may include calculating 206 a transform parameter. A transform parameter defines an estimated movement in the locations of particles between the first image of the particles and the second image of the particles. Using the transform parameter would allow one to determine the location of particles in the second image as a function of the location of the particles in the first image.
In certain embodiments, the transform parameter may include a radial component and a constant component. For example, FIG. 4 illustrates one embodiment of a radial movement and a constant movement that may be used to define the transform parameter. The constant movement may be induced by the optical effects of the filter wheel. As shown in FIG. 4, p0 illustrates the location of a particle in the first image, and similarly p1 illustrates the location of the particle in the second image. The movement can be defined as the combination of two components: (1) s the vector of the constant movement, and (2) a the radial movement. The radial component may be proportional to the distance between the location of particles within the first image and a center O of the image. As such, the location of p1 may be determined using Equation 1:

p1=T(p0)=p0+s+a(O−p0)  (1)

Using Equation 1 and the transform parameter (e.g., as defined by a and s), the location of particle in the first image may be determined. As defined above, Equation 1 (using components a and s) reflects a specific embodiment of the disclosed method. A transform parameter may include only a single component or even multiple components. The relationship between p0 and p1 need not be linear, and may even be non-linear.
In some embodiments, calculating 206 a transform parameter may include estimating 208 estimated locations of particles within a second image. As described above, an algorithm such a peak search algorithm may or may not find the location of particle in the second image. Such an algorithm may be used to estimate the locations of particles within a second image. Another method, referred to as the “maximal integral location” is discussed in more detail below.
Based on one or more of these estimated locations, potential transform parameters may be calculated 210. For example, a single pair of points (e.g., p0 and an estimated p1) may be used to define one or more potential transform parameters. In an embodiment of a transform parameter comprising more than one component, however, may require more than one pairs of points to determine a potential transform parameter. More than one pairs of points may be used to define a set of potential transform parameters—where each pair may define a single potential transform parameter. In some embodiments, by analyzing the set of potential transform parameters, an optimal transform parameter may be determined 212. In a simple embodiment, determining an optimal transform parameter may include taking the average, mean, mode, or the like of the set of potential transform parameters. Another method, using a Hough transform is described in more detail below.
In some embodiments, estimating 208 the locations of particles within the second image may include determining the maximal integral location based on the second image and the locations of the particles within the first image. FIG. 5A illustrates one embodiment of finding the maximal integral location (e.g., an estimated/p1). Finding the maximal integral location includes estimating the location of p1 in the second image based the location of p0 in the first image. In some embodiments, finding the maximal integral location includes finding the location of the point p0 in the second image. The original location of the particle p0 is illustrated in FIG. 5A. Finding the maximal integral location may further include analyzing the image a certain distance around p0. This certain distance around p0 may define an area to be analyzed. In FIG. 5A, the part of the image contained within the square may be analyzed. More specifically, in the figure, the square is defined by the three pixels in each direction from p0. In various embodiments this distance may be any number of pixels or other metric. In some embodiments, the certain distance around p0 may be defined by as a circle (e.g., based on a radius) around p0 instead.
Finding the maximal integral location may further include for one or more points q in the area to compute the integrated intensity centered at q. As shown in the figure, the set of pixels in the areas to be analyzed (e.g., the square box) may define a set of points q. For each point g, the integrated intensity is computed. In the specific embodiment of the figure, the area to be integrated is defined by the circle (e.g., with a radius of 2.5 pixels). In other embodiments, the area to be integrated may be defined by a square (e.g., with a half side length of 2.5 pixels). In some embodiments, the value of q that maximizes the integral intensity is estimated to be the location of the particle at in the second image (e.g., p1).
As shown in FIG. 5A, in this specific embodiment, the maximal integral location method correctly identifies the location of the particle in the second image (e.g., p1). In contrast, in FIG. 5B, the maximal integral location method does not correctly identify the location of the particle in the second image. In this figure, there is both a “dim” particle and a “bright” particle in the second image. Given the location p0 of the dim particle in the first image, the actual location of the particle in the second image is identified “correct bead center” in the figure. However, due to interference from the adjacent particle, the point p1 is found as the maximal integral location. The estimated location of the particle shown is FIG. 5(b) is an outlier. Whereas the pair of points (p0 and p1) from FIG. 5A may be used to calculate a “correct” transform parameter, the pair of points from FIG. 5B may result in an “incorrect” transform parameter.
The method steps discussed with respect to FIG. 5 are explained in a specific embodiment with respect to Equations 2 and 3. Equation 2 defines f(p) the sum of pixels in the second image (e.g., RP) about point p with integration radius r. Equation 3 defines K(p0,m) as the set of integrated intensities in the second image (e.g., RP) at all points p1 less than m distance from p0.

f(p)=Σ{q1|[q−p]|<r}RP(q)  (2)

K(pa,m)={f(p1):∀p1;∥<m}  (3)

Given the center of a particle p0 in the first image, the maximal integral location of this particle in the second image may be defined (in this specific embodiment) as the location p1 where f(p1) is the maximum of K(p0,m). Here, m is the maximum detectable movement of a particle from the first image to the second image.
After estimating 208 the locations of particles within a second image, a set of potential transform parameters may be calculated, and this set of transform parameters may be used to determine an optimal transform parameter. As discussed with respect to FIG. 5, some of the pairs of points (e.g., p0 and the estimated p1) may be “correct.” That is the estimated p1 corresponds to the actual location of the particle in the second image. Likewise, some of the pairs of points will be “incorrect.” Since many of the pairs of points will be right (e.g., in the set of all pairs it is more likely than not that p1 was estimated correctly), an analysis of the set can reveal an optimal transform parameter.
In a specific embodiment, potential transform parameters may be calculated based on each pair of points. Further, a transform may be used (e.g., a Hough transform) where each pair of points “votes” for potential transform parameters. Thus, the optimal transform parameter would be the potential transform parameter that gets the most overall votes. Alternative algorithms that may be used include 1) a random sample consensus (RANSAC) of p0, p1 pairs and 2) an elimination of variables by first specifically choosing p0,p1 close to O to eliminate the effect of a, solving for s, and then using both s and the remaining p0,p1 pairs to solve for a.
A specific embodiment of 210 calculating potential transform parameters and determining 212 an optimal transform parameter is described below with respect to Equations 4 and 5. This embodiment uses a transform parameters with a radial and constant component as defined by Equation 1. Equation 4 solves Equation 1 for s, and Equation 5 solves Equation 1 for a.

s=p2−p0−a(O−p0)  (4)





  
    
      
        a
        =
        
          
            
              
                (
                
                  o
                  -
                  
                    p
                    0
                  
                
                )
              
              T
            
            ⁢
            
              (
              
                
                  p
                  1
                
                -
                
                  p
                  0
                
                -
                s
              
              )
            
          
          
            
              
              
                
                
                  o
                  -
                  
                    p
                    0
                  
                
                
              
              
            
            2
          
        
      
    
    
      
        (
        5
        )
      
    
  




The pseudo code below introduces a parameter V (e.g., a Hough space or a vote space).















1) Initialize the array V to 0



2) For each value a in A do









For each pair p0, p1 do









Compute s(p0, p1, a)



Coerce s to closest quantized value s′ in S



Increment V[s′][a]









3) For each value s in S do









For each pair p0, p1 do









Compute a(p0, p1, s)



Coerce a to closest quantized value a′ in A



Increment V[s][a′]










In a specific embodiment, vectors may vary in the set of S. For example, the set S may be defined to vary from <−3, −3> to <3,3>, and S may quantized at every 1⁄5 pixel. Such a set of S will yield 30 possible variations. One of skill in the art will recognize that the set S may be defined to be larger or smaller either by increasing/decreasing the range of S or changing the quantization factor. The scalar a may vary in the set of A. For example, the set A may be defined to vary from −0.0004 to 0.0004 and quantized at every 0.0001 steps. Such a set of A has a size of 80. In embodiments where an image has 2048×2048 pixels, every 0.001 change in a corresponds to roughly 1 pixel of radial movement at the image boundary (since the distance from O to the edge of the image is 1024 pixels). A range of −0.004 to 0.004 could enable the detection of approximately 4 pixels of radial movement. One of skill in the art will recognize that the set A may be changed by changing the range of quantization. Increasing the range of either S or A could result in detecting larger components of radial and scalar movement. Moreover, using a finer quantization could result in a finer determinations of an optimal transform parameter.
FIG. 6 illustrates a particular example of portions of array V. Since s is a two-dimensional vector, it is represented on the x and y axis. As shown, S varies from −3 to 3 along both axes. The component a is represented on the z axis (vertically), and A varies from −0.004 to 0.004. Each of the FIGS. 6A, 6B, 6C, and 6D represent a slice of array V. As shown, FIG. 6C represents the highest peak in the vote space. This peak occurs at quantized location a=0.0009, s=<−0.4, 0.4>. The quantized location, however, may not be the “optimal transform parameter.” In some embodiments, computing the center of mass inside a volume of vote space about the peak location yields the optimal transform parameter: a*=9.07×10−4 and s*=<−0.3684, 0.3800>.
In some embodiments, the method 200 may further include applying 214 the transform parameter to the second image. In such an embodiment, the calculated transform parameter (e.g., the optimal transform parameter) may be used to determine the location of the particles in the second image based on the location of the particles in the first image. Using Equation 1 for example, the optimal transform parameter values a and s may be used to determine the movement of each of the particles.
In some embodiments, prior to applying 214 the transform parameter to the second image, an error vector may be calculated. Thus, the error vector may be account for at least one source of error in the calculation of the transform parameter. Specifically, the error vector may take into account the affect of neighboring particles. Moreover, neighboring particles may apply a force on a given particle causing them to move with respect to time. The error vector is defined in Equation 6.

ε=p1−T(p0;a*;s*)=p1−p0−s*−a*(O−p0)  (6)

As illustrated in FIG. 7, each particle q may exert a force on p0. More specifically, this force may include a magnetic force between the particles. The force exerted by a particle q may have a direction as defined by Equation 7, and the force exerted by a particle q may have a magnitude as defined by Equation 8. The magnitude of the force exerted by a particle q is inversely proportional to the square of the distance from q to p0. In addition to calculating the square of the distance between q and p0, Equation 8 introduces a variable g. The total error vector combining the forces exerted on p0 by all neighboring beads q is summed together in Equation 9.

Direction=±(q−p0)  (7)





  
    
      
        Magnitude
        =
        
          g
          
            
              
              
                
                
                  q
                  -
                  
                    p
                    0
                  
                
                
              
              
            
            2
          
        
      
    
    
      
        (
        8
        )
      
    
  
  
    
      
        ɛ
        =
        
          g
          ⁢
          
            
              ∑
              
                q
                ∈
                Q
              
            
            ⁢
            
              
                (
                
                  q
                  -
                  
                    p
                    0
                  
                
                )
              
              
                
                  
                  
                    
                    
                      q
                      -
                      
                        p
                        0
                      
                    
                    
                  
                  
                
                3
              
            
          
        
      
    
    
      
        (
        9
        )
      
    
  




Once an optimal value of g has been found, Equation 1 can be modified to account for the error vector. Equation 10 below accounts for the radial component, the constant component, and the error component:




  
    
      
        
          p
          1
        
        =
        
          
            T
            ⁡
            
              (
              
                p
                0
              
              )
            
          
          =
          
            
              p
              0
            
            +
            s
            +
            
              a
              ⁡
              
                (
                
                  O
                  -
                  
                    p
                    0
                  
                
                )
              
            
            +
            
              g
              ⁢
              
                
                  ∑
                  
                    q
                    ∈
                    Q
                  
                
                ⁢
                
                  
                    (
                    
                      q
                      -
                      
                        p
                        0
                      
                    
                    )
                  
                  
                    
                      
                      
                        
                        
                          q
                          -
                          
                            p
                            0
                          
                        
                        
                      
                      
                    
                    3
                  
                
              
            
          
        
      
    
    
      
        (
        10
        )
      
    
  




Similar to the calculation of the optimal transform parameter, the error component can be determined by estimating a set of potential values for the error component and calculating the optimal value of the error component. Equation 11 illustrates how to calculate g based on a given point p surrounded by particles q within a given radius to form the set Q




  
    
      
        
          
            Letting
            ⁢
            
              
            
            ⁢
            W
          
          =
          
            
              ∑
              
                q
                ∈
                Q
              
            
            ⁢
            
              
                
                  (
                  
                    q
                    -
                    
                      p
                      0
                    
                  
                  )
                
                
                  
                    
                    
                      
                      
                        q
                        -
                        
                          p
                          0
                        
                      
                      
                    
                    
                  
                  3
                
              
              ⁢
              
                :
              
            
          
        
        ⁢
        
          

        
        ⁢
        
          g
          =
          
            
              w
              s
              T
            
            
              
                
                
                  
                  w
                  
                
                
              
              2
            
          
        
      
    
    
      
        (
        11
        )
      
    
  




An optimal value of g can be calculated using the following pseudo code:















1) For each p0,p1 pair do









Compute ε



Define Q to be neighboring particles “close” to p0



Compute W



Compute g



If g is within a specified bounds (e.g., −20 to 0), record in









G









Select the optimal value g* from G










As discussed above, the optimal values of a, s, and g may then be used to determine the movement of particles between a first image and a second image.
It will be appreciated to those skilled in the art having the benefit of this disclosure that this invention is believed to provide computer-implemented methods, storage mediums, and systems for image data processing. Further modifications and alternative embodiments of various aspects of the invention will be apparent to those skilled in the art in view of this description. Accordingly, this description is to be construed as illustrative only and is for the purpose of teaching those skilled in the art the general manner of carrying out the invention. It is to be understood that the forms of the invention shown and described herein are to be taken as the presently preferred embodiments. Elements and materials may be substituted for those illustrated and described herein, parts and processes may be reversed, and certain features of the invention may be utilized independently, all as would be apparent to one skilled in the art after having the benefit of this description of the invention. Changes may be made in the elements described herein without departing from the spirit and scope of the invention as described in the following claims.

## Claims

1. A computer-implemented method, comprising:
determining locations of particles within a first image of the particles, wherein the particles have fluorescence-material associated therewith;
calculating a transform parameter, wherein the transform parameter defines an estimated shift in the locations of the particles between the first image of the particles and a second image of the particles, and wherein the transform parameter comprises a radial component and a constant component, the radial component proportional to a distance between the locations of the particles within the first image and a center of the first image;
applying the transform parameter to the locations of the particles within the first image to determine movement locations of the particles within the second image.

2. The computer-implemented method of claim 1, wherein calculating the transform parameter comprises:
estimating estimated locations of particles within the second image;
calculating potential transform parameters based on the locations of the particles within the first image and the estimated locations of the particles in the second image.

3. The computer-implemented method of claim 2, wherein calculating the transform parameter further comprises determining an optimal transform parameter based on the potential transform parameters.

4. The computer-implemented method of claim 2, wherein estimating the estimated locations of the particles within the second image comprises determining maximal integral locations based on the second image and the locations of the particles within the first image.

5. The computer-implemented method of claim 3, wherein determining the optimal transform parameter comprises calculating a Hough transform.

6. The computer-implemented method of claim 1, further comprising calculating an error component based on a force between the particles.

7. A non-transitory storage medium comprising program instructions which are executable by a processor for:
determining locations of particles within a first image of the particles, wherein the particles have fluorescence-material associated therewith;
calculating a transform parameter, wherein the transform parameter defines an estimated movement in the locations of the particles between the first image of the particles and a second image of the particles, and wherein the transform parameter comprises a radial component and a constant component, the radial component proportional to a distance between the locations of particles within the first image and a center of the first image;
applying the transform parameter to the locations of the particles within the first image to determine movement locations of the particles within the second image.

8. The non-transitory storage medium of claim 7, wherein calculating the transform parameter comprises:
estimating estimated locations of particles within the second image;
calculating potential transform parameters based on the locations of the particles within the first image and the estimated locations of the particles in the second image.

9. The non-transitory storage medium of claim 8, wherein calculating the transform parameter further comprises determining an optimal transform parameter based on the potential transform parameters.

10. The non-transitory storage medium of claim 8, wherein estimating the estimated locations of the particles within the second image comprises determining maximal integral locations based on the second image and the locations of the particles within the first image.

11. The non-transitory storage medium of claim 9, wherein determining the optimal transform parameter comprises calculating a Hough transform.

12. The non-transitory storage medium of claim 7, further comprising calculating an error component based on a force between the particles.

13. A system configured to acquire and process image data, comprising:
an imaging subsystem configured to image, at different wavelength bands, particles disposed within the imaging subsystem; and
a non-transitory storage medium comprising program instructions which are executable by a processor for:
storing data acquired for multiple images of the particles, wherein each of the multiple images corresponds to a different wavelength band;
determining locations of particles within a first image of the particles, wherein the particles have fluorescence-material associated therewith;
calculating a transform parameter, wherein the transform parameter defines an estimated movement in the locations of the particles between the first image of the particles and a second image of the particles, and wherein the transform parameter comprises a radial component and a constant component, the radial component proportional to a distance between the locations of particles within the first image and a center of the first image;
applying the transform parameter to the locations of the particles within the first image to determine movement locations of the particles within the second image.

14. The system of claim 13, wherein calculating the transform parameter comprises:
estimating estimated locations of particles within the second image;
calculating potential transform parameters based on the locations of the particles within the first image and the estimated locations of the particles in the second image.

15. The system of claim 14, wherein calculating the transform parameter further comprises determining an optimal transform parameter based on the potential transform parameters.

16. The system of claim 14, wherein estimating the estimated locations of the particles within the second image comprises determining maximal integral locations based on the second image and the locations of the particles within the first image.

17. The system of claim 15, wherein determining the optimal transform parameter comprises calculating a Hough transform.