Comet Assay Analysis

The single-cell gel electrophoresis (or comet) assay is a rapid and sensitive fluorescence microscopic method for the detection of primary DNA damage at the individual cell level. During the last 10 years, it has gained widespread acceptance for genotoxicity testing. Compared to other genotoxicity tests, the comet assay is an inexpensive and simple method that is well suited to analysis and measurement by means of computer-assisted techniques [Christoph Helma and Maria Uhl, Introduction in “A Public Domain Image Analysis Program for the Single Cell Gel Electrophoresis (Comet) Assay,” Mutation Research 466 (2000) 9-15].

See http://cometassay.com for additional information and references.

This loads the package.

<< ImageProcessing`

Here an example image is loaded (data the courtesy of Blake Whitaker, University of Southern Maine).

img = ImageRead["comet.tif"];

This defines the regions of interest as user-selected rectangular areas in the graphical coordinate system of the Mathematica front end. The coordinate pairs represent the lower-left and upper-right vertices of the bounding rectangles for each of six regions in the image.

\(\(v = 
        Round[{{{9.71295, 48.9237}, {161.548, 118.5148}}, 
{{98.1746, 
                120.8957}, {228.921, 189.16}}, 
{{191.072, 
                37.37961}, {336.58, 102.9163}}}];\)\)

In order to draw the regions of interest on the original image, define a utility function that takes two vertices and returns a Line graphics primitive depicting the bounding rectangle.

toLine[{{ax_, ay_}, {bx_, by_}}] := 
      Line[{{ax, ay}, {bx, ay}, {bx, by}, {ax, by}, {ax, ay}}]

Here is a the original image data showing the regions of interest.

Show[{Graphics[img], 
          Graphics[{RGBColor[1, 0, 0], toLine /@ v}]}];

Proceed with extracting the comets from the original image.

comet = ImageTake[img, 
              Sequence @@ #] & /@ (Reverse[Transpose[#]] & /@ 
              v);

This shows the individual comets.

Show[GraphicsArray[Graphics /@ comet]];

The original images are in three-channel RGB format. Clearly the green channel signal is strongly dominant, which is easily verified by displaying the planar format of an example comet.

Show[Graphics[PlanarImageData[comet[[1]]]]];

A closer examination reveals the presence of weak signals in the red and blue channels, especially in the tail section of the comet (interestingly, this effect is not discussed in the scientific literature). Since all computer-aided quantification results are traditionally presented in the context of gray-scale images, proceed with a color transformation.

cometG = Round[ToGrayLevel /@ comet];

This displays the histograms of the example comet images.

<< Graphics`Graphics`

DisplayTogetherArray[ShowImageHistogram[ImageHistogram[#], 
              Ticks -> {{50, 100, 150}, None}] & /@ 
          cometG];

Next is a typical calculation. For example, the interest here is to obtain the area of the head of the comet (“AH”), the total density of the head (“DNAH”), and the mean intensity of the head (“MH”). Usually, the same set of measurements is obtained for the tail, and the two results are contrasted. In the literature on comet assays, the head is identified by the pixels with the highest intensity. As is commonly used, the threshold value is 85 percent of the maximum intensity in the image.

0.85 Max /@ cometG

The images are binarized, adding a morphological filtering operation to smooth the borders of the head regions.

bin = BinaryOpen[#, {{1, 1, 1}, {1, 1, 1}, {1, 1, 1}}] & /@ 
          Thread[Threshold[cometG, 0.85 Max /@ cometG]];

Show[GraphicsArray[Graphics /@ %]];

This continues with identifying the regions and measuring the required properties, additionally including calculations of the centroids and enclosing circles.

hd = Where[#, 1] & /@ bin;

AH = RegionArea /@ hd

DNAH = 
      MapThread[(Plus @@ Extract[#1[[1]], #2] &), {cometG, 
          hd}]

MH = DNAH/(Length /@ hd) // N

This displays the results in a tabular fashion.

\({{"image #", "area", "density", "mean intensity"}, {Range[3], AH, DNAH, MH}} // TableForm\)

The graphical elements used to annotate the comet images are calculated.

gr = GraphicsRotate /@ 
            N[{Point[RegionCenter[#]], Circle @@ RegionCircle[#]}] & /@ 
        hd

MapThread[
        Show[{Graphics[#1], 
              Graphics[{RGBColor[1, 0, 0], PointSize[0.03], #2}]}, 
            DisplayFunction  -> Identity] &, {comet, gr}];
Show[GraphicsArray[%]];

This is an alternative visualization of a comet image.

ListPlot3D[cometG[[3, 1]], Mesh -> False, 
        Boxed -> False, PlotRange -> {0, 160}, 
        ViewPoint -> {1.552,  -2.866,  0.910}];

Finally, another property of a comet, the so-called comet profile, is calculated and displayed. It is a plot of the cross-sectional density (density measured perpendicularly to the lengthwise dimension of the comet).

Tr /@ Transpose[#[[1]]] & /@ cometG;
Show[
        GraphicsArray[ListPlot[#, PlotJoined -> True, 
                DisplayFunction -> Identity] & /@ %]];