CN112927182A - Digital PCR microarray image analysis method - Google Patents

Digital PCR microarray image analysis method Download PDF

Info

Publication number
CN112927182A
CN112927182A CN202011638314.6A CN202011638314A CN112927182A CN 112927182 A CN112927182 A CN 112927182A CN 202011638314 A CN202011638314 A CN 202011638314A CN 112927182 A CN112927182 A CN 112927182A
Authority
CN
China
Prior art keywords
image
channel
region
units
area
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Withdrawn
Application number
CN202011638314.6A
Other languages
Chinese (zh)
Inventor
唐骏
余跃
肖旻
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Xiamen University of Technology
Original Assignee
Xiamen University of Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Xiamen University of Technology filed Critical Xiamen University of Technology
Priority to CN202011638314.6A priority Critical patent/CN112927182A/en
Publication of CN112927182A publication Critical patent/CN112927182A/en
Withdrawn legal-status Critical Current

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/0002Inspection of images, e.g. flaw detection
    • G06T7/0012Biomedical image inspection
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T3/00Geometric image transformations in the plane of the image
    • G06T3/40Scaling of whole images or parts thereof, e.g. expanding or contracting
    • G06T3/4007Scaling of whole images or parts thereof, e.g. expanding or contracting based on interpolation, e.g. bilinear interpolation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • G06T7/11Region-based segmentation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • G06T7/136Segmentation; Edge detection involving thresholding
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • G06T7/187Segmentation; Edge detection involving region growing; involving region merging; involving connected component labelling
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/60Analysis of geometric attributes
    • G06T7/62Analysis of geometric attributes of area, perimeter, diameter or volume
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10056Microscopic image
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10064Fluorescence image

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Health & Medical Sciences (AREA)
  • Geometry (AREA)
  • General Health & Medical Sciences (AREA)
  • Medical Informatics (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • Radiology & Medical Imaging (AREA)
  • Quality & Reliability (AREA)
  • Investigating, Analyzing Materials By Fluorescence Or Luminescence (AREA)

Abstract

The invention provides a digital PCR microarray image analysis method, which comprises the following steps: inputting a first channel fluorescence image, and preprocessing the first channel fluorescence image; carrying out image segmentation on the preprocessed fluorescence image; identifying effective units of the segmented image and extracting mass center coordinates of the effective units; positioning and extracting effective units in the second and third channel images according to the coordinates, and counting the fluorescence signal intensity of the effective units in the second and third channel images according to the effective units; and drawing a scatter diagram and a histogram of the fluorescence intensity of the reaction units in the second and third channel images according to the signal intensity, counting the proportion of the positive units, and calculating the concentration of the target nucleic acid molecules in the original reaction liquid by combining Poisson distribution. According to the invention, the mask is established by the first channel image so as to indirectly analyze the second channel image and the third channel image, so that the multichannel digital PCR microarray image under non-uniform illumination can be effectively analyzed, the accuracy of image analysis is improved, and absolute quantitative analysis of target nucleic acid molecules is realized.

Description

Digital PCR microarray image analysis method
Technical Field
The invention relates to the field of image processing, in particular to a digital PCR microarray image analysis method.
Background
The digital PCR is a nucleic acid detection technology with high sensitivity, high specificity and absolute quantification, and utilizes the concept of 'divide-and-conquer' to subdivide a PCR reaction system into tens of thousands of micro-reaction systems, so that all the reaction systems independently carry out PCR circulating amplification reaction, after the reaction is finished, the proportion of positive units is counted, and the concentration of target nucleic acid molecules in the original reaction liquid can be calculated by utilizing Poisson distribution. In recent years, with the progress of chip manufacturing technology, chips containing tens of thousands of micro reaction cells are manufactured by micromachining technology. The chip in the Quantastudio 3D chip-based digital PCR system introduced by Life Techmologices, Inc. in 2013 already contains twenty thousand microcell structures. Compared with qPCR, the method can realize absolute quantitative analysis without depending on a standard curve as a reference, and has been widely applied to the fields of molecular diagnosis, virus detection, food safety and the like. However, the statistical accuracy of the positive units after PCR reaction is not high, and at present, image processing technology is mainly used to identify the positive units in the digital PCR microarray image, and with the increasing update of image analysis technology, the analysis of fluorescence images mainly includes: image preprocessing, image segmentation, image recognition, post-data analysis and the like, wherein an image segmentation algorithm plays a key role. Currently, the gridding analysis method, the automatic cell positioning method and the deep learning technology are mainly used for segmentation, and the method has the defect that the segmentation can only be carried out on orthogonal and vertically distributed microarray images, and the method is not suitable for non-orthogonal images. And the segmentation threshold value is calculated by an Otsu method or the segmentation is carried out by calculating double threshold values by the Otsu method under the condition that the illumination of the image is uniform. And the image segmentation accuracy for uneven illumination and blurred images is low.
Therefore, the development of a high-accuracy digital PCR microarray image analysis method is of great significance in solving the problems of inaccurate segmentation, low recognition rate, weak robustness and the like of the existing digital PCR image analysis technology.
Disclosure of Invention
The present invention is directed to a digital PCR microarray image analysis method to solve the above-mentioned problems.
In order to achieve the above object, an embodiment of the present invention provides a digital PCR microarray image analysis method, including
Inputting a first channel fluorescence image, and preprocessing the first channel fluorescence image;
performing image segmentation on the preprocessed first channel fluorescence image to obtain a segmented image;
identifying effective units of the segmented image and extracting centroid position coordinates of the effective units;
positioning and extracting effective units in the images of the second channel and the third channel according to the coordinates of the mass center of the effective units, and counting the fluorescence signal intensity of the effective units in the images of the second channel and the third channel according to the effective units in the images of the second channel and the third channel;
drawing a scatter diagram and a cylindrical diagram of the fluorescence intensity distribution of the reaction units in the second channel image and the third channel image according to the fluorescence signal intensity of the effective units in the second channel image and the third channel image, and counting the proportion of the positive units;
and calculating the concentration of the target nucleic acid molecules in the original reaction solution based on the proportion of the positive units and the Poisson distribution.
Further, inputting the first channel fluorescence image, and preprocessing the first channel fluorescence image specifically comprises:
carrying out gray scale processing on the first channel fluorescence image;
carrying out bilinear interpolation processing on the image subjected to the gray processing to obtain a bilinear interpolation result;
performing image enhancement processing on the first channel fluorescence image based on the bilinear interpolation result;
and sequentially carrying out median filtering and Gaussian filtering on the image after the image enhancement so as to finish preprocessing.
Further, the image enhancement of the first channel fluorescence image based on the bilinear interpolation result specifically includes:
performing adaptive histogram equalization on the first channel fluorescence image to obtain an equalized image;
carrying out top-hat transformation on the equalized image to obtain a top-hat transformation result;
performing bottom-cap transformation on the equalized image to obtain a bottom-cap transformation result;
and summing the top-hat transformation result and the bilinear interpolation result, and then subtracting the result from the bottom-hat transformation result to complete image enhancement.
Further, the image segmentation is performed on the preprocessed first channel fluorescence image to obtain a segmented image specifically as follows:
segmenting the preprocessed first channel fluorescence image by using an Otsu method to obtain a first segmentation image;
expanding the first segmentation image to obtain a maximum connected region image;
performing logical AND operation on the maximum connected region image and the first segmentation image to obtain a region-of-interest segmentation image;
acquiring a gray-scale image of an under-segmentation area in the segmentation image of the area of interest, and preprocessing the gray-scale image;
according to the image information entropy improvement Otsu method, segmenting the preprocessed under-segmentation region again to obtain a second segmentation image;
and merging the first segmentation image and the second segmentation image to finish image segmentation, wherein merging is finished by using logical operation.
Further, the acquiring a gray map of an under-segmented region in the segmentation map of the region of interest and the preprocessing thereof specifically include:
setting an upper limit of an area threshold value to screen out a large-area from the under-segmentation area;
expanding the large-area to obtain an expanded connected area gray scale image;
and preprocessing the expanded connected region gray level image line.
Further, the above-mentioned method of improving Otsu according to image information entropy is to perform segmentation again on the preprocessed under-segmented region to obtain a second segmented image specifically:
the maximum gray value and the minimum gray value in the under-segmentation region after the pre-processing are respectively set as Zmax and Zmin, and the average value of the maximum gray value and the minimum gray value is calculated as an initial segmentation threshold T0:
Figure BDA0002877411380000051
calculating the information entropy of the global image:
Figure BDA0002877411380000052
wherein E represents the information entropy of the global image, h (i) represents the number of pixels with the gray value of i, and N represents the total number of pixel points of the whole image;
reducing gray values and pixels at two ends of the image gray level at the speed of step 1, and setting the minimum gray value and the maximum gray value after reduction as alpha and beta respectively;
calculating the information entropy of the (alpha, T0-1) and (T0, beta) intervals, which are respectively:
Figure BDA0002877411380000053
Figure BDA0002877411380000054
calculating the proportion of the number of pixels in the (alpha, beta) interval and the proportion of the number of pixels in the (alpha, beta) interval in the (alpha, T0-1) interval and the (T0, beta) interval in the total number of pixels in the (alpha, beta) interval, wherein the proportion is respectively as follows:
Figure BDA0002877411380000055
Figure BDA0002877411380000056
calculating (alpha, beta) interval average information entropy:
E2=r0E0+r1E1
if any number greater than 0 exists, when the minimum gray value is alpha and the maximum gray value is beta, the iteration condition is satisfied:
ΔE=|E-E2|<ε;
calculating a threshold value when the variance between the classes of the interested gray scale interval is maximum according to the interested gray scale interval (alpha, beta);
and segmenting the under-segmentation region according to the threshold value.
Further, the identifying the effective unit of the segmented image and extracting the centroid position coordinates of the effective unit specifically include:
according to the upper and lower limits of the area threshold value, t is respectively1And t2Dividing the region of interest segmentation graph into a small area region, a normal area region and a large area region, wherein the area is smaller than t1Is divided into a small area A1 with an area t1And t2Is divided into a normal area A2 with an area greater than t2Is divided into a large area A3;
setting a circularity threshold, and merging regions which do not meet the circularity threshold in the region A2 into the region A1;
acquiring a connected region around the region A3, merging the connected region into the region a1 to form a new region a 12;
excluding the connected region around region A3 from region a2, resulting in region a 22;
carrying out edge detection on the preprocessed first channel fluorescence image by using a canny operator, and combining the region A12 with an edge detection result to obtain a region A13;
expanding the area A12, and performing logical AND operation with the area A13 to form an area A14;
filling and opening the area A14 to obtain an area A15;
setting an eccentricity threshold, and removing a connected region which does not meet the eccentricity threshold range in the region A15 to obtain a region A16;
merging the region a16 with the region a22 to obtain an image a 4;
calculating the average gray value of all the areas in the area A22; setting the upper and lower limits of a gray threshold value to be k respectively according to the average gray value1、k2The gray value in all connected domains in the image A4 is less than k1And is greater than k2The gray value of the pixel position of (a) is set to 0, and an image A41 is obtained;
calculating the sizes of holes of all connected domains in the image A41; setting a hole threshold, and removing a connected domain which does not meet the hole threshold in the image A41 to obtain an image A42;
calculating the circularity sizes of all connected regions in the image A42, and removing the connected regions which are not in the threshold range according to a circularity threshold to obtain an image A43;
removing the connected domains with too large and too small areas in the image A43 to complete the identification of the effective units;
and marking the centroid of the finally identified connected domain on the gray scale image of the first channel fluorescence image.
Further, the positioning and extracting the effective units in the second channel image and the third channel image according to the coordinates of the effective unit centroid positions, and counting the fluorescence signal intensities of the effective units in the second channel image and the third channel image according to the effective units in the second channel image and the third channel image specifically include:
and counting the number and the centroid positions of the effective units in the first channel image, and reading the fluorescence signal intensities of the effective units in the second channel image and the third channel image.
Furthermore, according to the fluorescence signal intensities of the effective units in the images of the second channel and the third channel, drawing a scatter diagram and a histogram of the fluorescence intensity distribution of the effective units in the images of the second channel and the third channel, and counting the proportion of the positive units:
drawing a histogram according to the fluorescence signal intensities of all effective units of the second channel image and the third channel image;
searching two extreme point coordinates with maximum peak span in the histogram, taking the mean value of the two extreme point coordinates as a positive judgment threshold, if the fluorescence intensity is greater than the threshold, judging the two extreme point coordinates as a positive unit, otherwise, judging the two extreme point coordinates as a negative unit;
and counting the number of positive units and the total number of effective units according to the judgment threshold and the fluorescence signal intensity of the effective units, and drawing a scatter diagram, a histogram and a two-parameter scatter diagram of the fluorescence intensity of the images of the second channel and the third channel.
Further, the concentration of the target nucleic acid molecule in the original reaction solution is calculated based on the proportion of the positive units and in combination with poisson distribution, and the calculation formula is as follows:
Figure BDA0002877411380000081
wherein p represents the ratio of the number of positive units to the total number of valid units, NtThe total number of wells of the chip is shown, and V is the total volume of the reaction solution.
The invention provides a multichannel digital PCR microarray image analysis method, which comprises the steps of inputting a first channel fluorescence image and preprocessing the first channel fluorescence image; performing image segmentation on the preprocessed fluorescence image to obtain a segmented image; identifying effective units of the segmented image and extracting centroid position coordinates of the effective units; positioning and extracting effective units in the images of the second channel and the third channel according to the coordinates of the mass center of the effective units, and counting the fluorescence signal intensity of the effective units in the images of the second channel and the third channel according to the effective units in the images of the second channel and the third channel; drawing a scatter diagram and a histogram of the fluorescence intensity distribution of the reaction units in the second channel image and the third channel image according to the fluorescence signal intensity of the effective units in the second channel image and the third channel image, and counting the proportion of the positive units; and calculating the concentration of the target nucleic acid molecules in the original reaction solution based on the proportion of the positive units and the Poisson distribution. The invention establishes a mask by the first channel image so as to indirectly analyze the second channel image and the third channel image, can effectively analyze the multichannel digital PCR microarray image under the non-uniform illumination, improves the accuracy of image analysis and realizes the absolute quantitative analysis of target nucleic acid molecules.
Drawings
In order to more clearly illustrate the technical solutions of the embodiments of the present invention, the drawings that are required to be used in the embodiments will be briefly described below, it should be understood that the following drawings only illustrate some embodiments of the present invention and therefore should not be considered as limiting the scope, and for those skilled in the art, other related drawings can be obtained according to the drawings without inventive efforts.
FIG. 1 is a flow chart of a method for analyzing a digital PCR microarray image according to an embodiment of the present invention.
FIG. 2 is a first channel fluorescence image of the present invention;
FIG. 3 is a graph of the effect of the first channel fluorescence image after adaptive histogram equalization according to the present invention;
FIG. 4 is a diagram illustrating enhanced fluorescence images of a first channel according to the present invention;
FIG. 5 is a graph of the effect of the invention after median filtering in a first channel fluorescence image;
FIG. 6 is a graph of the effect of the first channel fluorescence image after Gaussian filtering according to the present invention;
FIG. 7 is an Otsu segmented image of the present invention;
FIG. 8 is a diagram of the present invention after removing peripheral impurities;
FIG. 9 is a gray scale diagram of an Otsu method segmentation under-segmentation region according to the present invention;
FIG. 10 is an under-segmented region segmentation image of the present invention;
FIG. 11 is a final segmentation of a first channel image of the present invention;
FIG. 12 is an image of a large area segmented by Otsu method of the present invention;
FIG. 13 is an image of a large area and its surrounding connected domain of the present invention;
FIG. 14 is an image of a normal area region excluding a large area of the surrounding connected component region according to the present invention;
FIG. 15 is a small area region of the present invention;
FIG. 16 is a large area surrounding region passing edge detection image of the present invention;
FIG. 17 is a comparison of eccentricity screens of the present invention;
FIG. 18 is an image of the present invention before hole, circularity, and area screening;
FIG. 19 is a comparison of the present invention after hole, circularity, and area screening;
FIG. 20 is a label diagram of the finally identified valid cell of the present invention;
FIG. 21 is a fluorescence image of a second channel of the present invention;
FIG. 22 is a fluorescence image of the invention and of a third channel;
FIG. 23 is a gray scale diagram of the effective cells of the second channel image extracted from the first channel image according to the present invention;
FIG. 24 is a gray scale diagram of the effective cells of a third channel image extracted from a first channel image according to the present invention;
FIG. 25 is a digital PCR image analysis APP interface of the present invention.
Detailed Description
In order to make the objects, technical solutions and advantages of the embodiments of the present invention more apparent, the embodiments of the present invention will be described clearly and completely with reference to the accompanying drawings, and it is to be understood that the described embodiments are some, but not all embodiments of the present invention. All other embodiments, which can be obtained by a person skilled in the art without any inventive step based on the embodiments of the present invention, are within the scope of the present invention. Thus, the following detailed description of the embodiments of the present invention, presented in the figures, is not intended to limit the scope of the invention, as claimed, but is merely representative of selected embodiments of the invention.
Referring to fig. 1, an embodiment of the present invention provides a digital PCR microarray image analysis method, including:
and S11, inputting the first channel fluorescence image and preprocessing the first channel fluorescence image.
In this embodiment, the first channel fluorescence image is first subjected to gray scale processing; and then carrying out bilinear interpolation processing on the image subjected to the gray processing to obtain a bilinear interpolation result.
And carrying out image enhancement processing on the first channel fluorescence image based on the bilinear interpolation result. The image enhancement processing process specifically comprises the following steps:
performing adaptive histogram equalization on the first channel fluorescence image to obtain an equalized image, referring to fig. 3;
carrying out top-hat transformation on the equalized image to obtain a top-hat transformation result;
performing bottom-cap transformation on the equalized image to obtain a bottom-cap transformation result;
and summing the top-hat transformation result and the bilinear interpolation result, and then subtracting the result from the bottom-hat transformation result, and referring to the figure 4, finishing the image enhancement.
The image after image enhancement is subjected to median filtering and gaussian filtering in sequence, with reference to fig. 5 and 6, to complete preprocessing.
And S12, performing image segmentation on the preprocessed first-channel fluorescence image to obtain a segmented image.
In the present embodiment, the preprocessed first channel fluorescence image is first segmented using the Otsu method to obtain a first segmented image. Specifically, a threshold value Thre1 which maximizes the variance between image classes is calculated by an Otsu method, and a global image is segmented by the threshold value Thre1, wherein the formula is as follows;
Figure BDA0002877411380000121
wherein f is1(x, y) is the image to be segmented, f* 1(x, y) is the segmented image, L is the image maximum gray level; otsu method segmentation of image parametersRefer to fig. 7.
In this embodiment, the first segmented image is then dilated to obtain a maximum connected region image. The specific process is as follows:
and solving the areas of all connected domains, sequencing the areas in an increasing mode, and taking the average value of the middle 80% of the areas as the average area of the normal unit.
And (3) solving an equal diameter R corresponding to the area average value of the normal unit:
Figure BDA0002877411380000122
in the formula, S is an area.
The segmentation image is expanded by taking the equal diameter R as the size of an expansion structure element.
And (4) counting the area of all the connected domains after expansion, reserving the connected domain with the largest area, and removing other connected domains.
In the present embodiment, logical and operation is performed on the maximum connected region image and the first divided image described above to obtain a region-of-interest divided map, referring to fig. 8.
In this embodiment, a grayscale image of an under-segmented region in the segmentation map of the region of interest is obtained, refer to fig. 9, and is preprocessed; the pretreatment process comprises the following steps:
setting an upper limit of an area threshold value to screen out a large area region; for example, 1.5 times the normal area may be set as the upper area threshold t1, and 0.5 times may be set as the lower area threshold t 2; and screening out the under-divided areas according to the upper limit of the area threshold. It should be noted that the upper and lower limits of the area threshold may be other parameters, and these schemes are all within the scope of the present invention.
In this embodiment, the large-area region is expanded to obtain a gray scale map of the expanded connected region; for example, the size of the inflatable structure may be 0.8 times the isodiametric size of the normal cells.
In this embodiment, the post-expansion connected component gray-scale map line is preprocessed. The pretreatment specifically comprises:
and carrying out bilinear interpolation processing on the expanded connected region gray-scale image to obtain a bilinear interpolation result.
And carrying out adaptive histogram equalization on the expanded connected region gray level image.
And performing top-hat transformation on the image subjected to the self-adaptive histogram equalization to obtain a top-hat transformation result.
And performing bottom-cap transformation on the image subjected to the self-adaptive equalization to obtain a bottom-cap transformation result.
And summing the top-hat transformation result and the bilinear interpolation result, and then carrying out difference on the result and the bottom-hat transformation result.
And sequentially carrying out median filtering and Gaussian filtering on the summed and differentiated image so as to finish preprocessing the expanded connected region gray level image line.
In the present embodiment, the preprocessed under-divided region is again divided according to the image information entropy improvement Otsu method to obtain a second divided image. The specific segmentation process is as follows:
the maximum gray value and the minimum gray value in the under-segmentation region after the pre-processing are respectively set as Zmax and Zmin, and the average value of the maximum gray value and the minimum gray value is calculated as an initial segmentation threshold T0
Figure BDA0002877411380000131
Calculating the information entropy of the global image:
Figure BDA0002877411380000141
wherein E represents the information entropy of the global image, h (i) represents the number of pixels with the gray value of i, and N represents the total number of pixel points of the whole image;
and reducing the gray values at two ends of the image gray level and pixels thereof at the speed of step 1, and setting the minimum gray value and the maximum gray value after reduction as alpha and beta respectively.
Calculating the information entropy of the (alpha, T0-1) and (T0, beta) intervals, which are respectively:
Figure BDA0002877411380000142
Figure BDA0002877411380000143
calculating the proportion of the number of pixels in the (alpha, beta) interval and the proportion of the number of pixels in the (alpha, beta) interval in the (alpha, T0-1) interval and the (T0, beta) interval in the total number of pixels in the (alpha, beta) interval, wherein the proportion is respectively as follows:
Figure BDA0002877411380000144
Figure BDA0002877411380000145
calculating (alpha, beta) interval average information entropy:
E2=r0E0+r1E1
if any number greater than 0 exists, when the minimum gray value is alpha and the maximum gray value is beta, the iteration condition is satisfied:
ΔE=|E-E2|<ε
where ε may be taken to be 0.3E.
According to the interested gray level interval (alpha, beta), calculating the threshold value Thte when the variance between the interested gray level interval classes is maximum2(ii) a According to the threshold value Thte2Segmenting the under-segmented region
Figure BDA0002877411380000151
In the formula f2(x, y) is an image to be segmented of the under-segmented region,
Figure BDA0002877411380000152
for the image after the under-divided region division, L is the maximum gray level of the image, and the division result refers to fig. 10.
And merging the first segmentation image and the second segmentation image to remove the phenomenon of disconnection between the reaction units due to expansion, and completing image segmentation, referring to fig. 11, wherein the merging is completed by using logical operation.
And S13, identifying effective units of the segmented image and extracting the centroid position coordinates of the effective units.
In the embodiment, firstly, the region of interest is divided into three types, namely a small area region, a normal area region and a large area region according to an upper area threshold t1 and a lower area threshold t 2; wherein the area is smaller than t1Is divided into a small area A1 at t1And t2Is divided into normal area A2, which is larger than t2Is divided into a large area a3, see fig. 12.
Setting a circularity threshold, and merging regions which do not meet the circularity threshold in the region A2 into the region A1; specifically, the circularities of all the cells in the area a2 are obtained, and the calculation formula of the circularities is:
Figure BDA0002877411380000153
where S and C are the area and perimeter of each cell, respectively. When F is 1, the region is circular, and when F is greater than 1, the region is other shapes. By setting the circularity threshold to 0.86, regions that are not within the circularity range are merged into the region a 1.
The connected regions around region A3 were obtained and merged into region a1 as shown in fig. 13, forming a new region a 12.
The connected region around region A3 was excluded from region a2 to give region a22, region a22, as shown in fig. 14.
Edge detection was performed on the preprocessed first channel fluorescence image using the canny operator to obtain region a 13.
The region a12 (see fig. 15) is combined with the edge detection result to obtain a region a 13.
The region a12 is expanded so that the expanded structure size may be 0.2 times the normal dot level size, and the expanded image is logically anded with the image of the region a13 to form the region a14, with the result being referred to fig. 16.
The region a14 is filled and opened to obtain a region a 15.
Setting the eccentricity threshold to 0.6, screening out connected regions that do not meet the eccentricity threshold, and removing connected regions that do not meet the eccentricity threshold range in region a15, resulting in region a16, as shown in fig. 17.
Region a16 was combined with region a22 to obtain region a4 after the preliminary screening, see fig. 18.
Calculating the average gray value of all the cells in the area A22, and recording as g; setting the upper and lower limits of a gray threshold value to be k respectively according to the average gray value1、k2E.g. the mean gray value is 1g, then k1May be 0.5g, k2It may be 1.5g, and setting the gray values of all connected domains in image a4 to be less than 0.5g and greater than 1.5g to be 0, resulting in image a 41.
And (3) calculating the sizes of the holes of all the connected regions in the image A41, wherein the size of the holes refers to the number of pixels with the gray value of 0 in the connected regions, setting the hole threshold value to be 12, and removing the connected regions with the sizes of the holes larger than the threshold value in the connected regions to obtain an image A42.
And calculating the circularity sizes of all connected regions in the image A42, and removing the connected regions which do not meet the circularity threshold value according to the circularity threshold value to obtain an image A43.
Connected components of image a43 that are too large and too small in area are removed as shown in fig. 19.
The centroid of the finally identified connected domain is marked on the grayscale map of the first channel as shown in fig. 20.
And S14, positioning and extracting effective units in the images of the second channel and the third channel according to the coordinates of the mass center positions of the effective units, and counting the fluorescence signal intensity of the effective units in the images of the second channel and the third channel according to the effective units in the images of the second channel and the third channel.
In the present embodiment, the first channel image is used as a reference image, the second channel image, as shown in fig. 21, and the third channel image, as shown in fig. 22, are used as floating images, and image registration is performed by using mutual information as a metric.
According to the recognized centroid positions and the number of the effective units of the first channel image, the effective unit positions in the second channel image and the third channel image are located, and the effective units in the second channel image and the third channel image are extracted, as shown in fig. 23 and fig. 24;
and (4) counting the fluorescence signal intensity of all effective units according to the effective units in the second channel image and the third channel image.
And S15, drawing a scatter diagram and a histogram of the fluorescence intensity distribution of the reaction units in the second channel image and the third channel image according to the fluorescence signal intensity of the effective units in the second channel image and the third channel image, and counting the proportion of the positive units.
In this embodiment, a histogram is drawn according to the fluorescence signal intensities of all effective sampling points of the second channel image and the third channel image, two maximum points with the largest peak span in the histogram are found, the mean value of the maximum points is taken as a positive determination threshold, if the fluorescence intensity is greater than the threshold, the determination is made as a positive unit, otherwise, the determination is made as a negative unit.
And counting the number of positive units and the total number of effective units according to the judgment threshold and the fluorescence value of the reaction units. And drawing a scatter diagram, a histogram, a two-parameter scatter diagram and a chip effect diagram of fluorescence values of the images of the second channel and the third channel, and displaying the images by manufacturing an APP interface, as shown in a drawing area in FIG. 25.
S16: and calculating the concentration of the target nucleic acid molecules in the original reaction solution based on the proportion of the positive units and the Poisson distribution.
In this example, the concentration of the target nucleic acid molecule in the original reaction solution is calculated by the following formula:
Figure BDA0002877411380000181
wherein p represents the ratio of the number of positive units to the total number of valid units, NtThe total number of wells of the chip is shown, and V is the total volume of the reaction solution.
Calculating the confidence interval of the concentration according to the statistical principle:
Figure BDA0002877411380000182
where lambda is the mean of the poisson distribution,
Figure BDA0002877411380000183
is the standard deviation of the poisson distribution and n is the number of significant elements in the image. The image parameters of the second channel and the third channel are shown in a parameter area below an interface of FIG. 25, and meanwhile, the threshold value can be adjusted by manually clicking the interface in the interface, so that the positive threshold value can be more accurately judged, and the absolute quantitative precision of the digital PCR is improved.
The method for analyzing the multichannel digital PCR microarray image provided by the embodiment comprises the following steps: inputting a first channel fluorescence image, and preprocessing the first channel fluorescence image; performing image segmentation on the preprocessed fluorescence image to obtain a segmented image; identifying effective units of the segmented image and extracting centroid position coordinates of the effective units; positioning and extracting effective units in the second channel image and the third channel image according to the effective unit mass center position coordinate, and counting the fluorescence signal intensity of the effective units in the second channel image and the third channel image according to the effective units in the second channel image and the third channel image; drawing a scatter diagram and a histogram of the fluorescence intensity distribution of the reaction units in the second channel image and the third channel image according to the fluorescence signal intensity of the effective units in the second channel image and the third channel image, and counting the proportion of the positive units; and calculating the concentration of the target nucleic acid molecules in the original reaction solution based on the positive unit proportion and the Poisson distribution. The mask is established on the first channel image so as to indirectly analyze the second channel image and the third channel image, the multichannel digital PCR microarray image under non-uniform illumination can be effectively analyzed, the accuracy of image analysis is improved, and the absolute quantitative analysis of target nucleic acid molecules is realized.
While the invention has been particularly shown and described with reference to a preferred embodiment, it will be understood by those skilled in the art that various changes in form and detail may be made therein without departing from the spirit and scope of the invention as defined by the appended claims.

Claims (10)

1. A digital PCR microarray image analysis method, comprising:
inputting a first channel fluorescence image, and preprocessing the first channel fluorescence image;
performing image segmentation on the preprocessed fluorescence image to obtain a segmented image;
identifying effective units of the segmented image and extracting centroid position coordinates of the effective units;
positioning and extracting effective units in the second channel image and the third channel image according to the effective unit mass center position coordinates, and counting the fluorescence signal intensity of the effective units in the second channel image and the third channel image according to the effective units in the second channel image and the third channel image;
drawing a scatter diagram and a histogram of the fluorescence intensity distribution of the reaction units in the second channel image and the third channel image according to the fluorescence signal intensity of the effective units in the second channel image and the third channel image, and counting the proportion of the positive units;
and calculating the concentration of the target nucleic acid molecules in the original reaction solution based on the proportion of the positive units and the Poisson distribution.
2. The method for analyzing digital PCR microarray images according to claim 1, wherein the inputting of the first channel fluorescence image and the pre-processing thereof are specifically:
carrying out gray scale processing on the first channel fluorescence image;
carrying out bilinear interpolation processing on the image subjected to the gray processing to obtain a bilinear interpolation result;
performing image enhancement processing on the first channel fluorescence image based on the bilinear interpolation result;
and sequentially carrying out median filtering and Gaussian filtering on the image after image enhancement to finish preprocessing.
3. The radix digital PCR microarray image analysis method of claim 2, wherein the image enhancement of the first channel fluorescence image based on the bilinear interpolation result is specifically:
performing adaptive histogram equalization on the first channel fluorescence image to obtain an equalized image;
carrying out top-hat transformation on the equalized image to obtain a top-hat transformation result;
performing bottom-cap transformation on the equalized image to obtain a bottom-cap transformation result;
and summing the top-hat transformation result and the bilinear interpolation result, and then subtracting the result from the bottom-hat transformation result to complete image enhancement.
4. The method for analyzing digital PCR microarray image according to claim 1, wherein the image segmentation is performed on the preprocessed first channel fluorescence image to obtain a segmented image, and the segmented image is specifically:
segmenting the preprocessed first channel fluorescence image by using an Otsu method to obtain a first segmentation image;
expanding the first segmentation image to obtain a maximum connected region image;
performing logical AND operation on the maximum connected region image and the first segmentation image to obtain an interesting region segmentation image;
acquiring a gray-scale image of an under-segmentation area in the segmentation image of the area of interest, and preprocessing the gray-scale image;
according to the image information entropy improvement Otsu method, segmenting the preprocessed under-segmentation region again to obtain a second segmentation image;
and merging the first segmentation image and the second segmentation image to finish image segmentation, wherein merging is finished by using logical operation.
5. The method for analyzing the digital PCR microarray image according to claim 4, wherein the obtaining of the grayscale map of the under-segmented region in the segmentation map of the region of interest and the pre-processing thereof are specifically:
setting an upper limit of an area threshold value to screen out a large-area from the under-segmentation area;
expanding the large-area to obtain an expanded connected area gray scale image;
and preprocessing the expanded connected region gray level image line.
6. The method for analyzing the digital PCR microarray image according to claim 4, wherein the pre-processed under-segmented region is segmented again according to the Otsu method for improving the entropy of the image information to obtain the second segmented image specifically comprises:
the maximum gray value and the minimum gray value in the under-segmentation region after the pre-processing are respectively set as Zmax and Zmin, and the average value of the maximum gray value and the minimum gray value is calculated as an initial segmentation threshold T0:
Figure FDA0002877411370000031
calculating the information entropy of the global image:
Figure FDA0002877411370000032
wherein E represents the information entropy of the global image, h (i) represents the number of pixels with the gray value of i, and N represents the total number of pixel points of the whole image;
reducing gray values and pixels at two ends of the image gray level at the speed of step 1, and setting the minimum and maximum gray values after reduction as alpha and beta respectively;
calculating the information entropy of the (alpha, T0-1) and (T0, beta) intervals, which are respectively:
Figure FDA0002877411370000041
Figure FDA0002877411370000042
calculating the proportion of the pixel number of the (alpha, T0-1) interval and the pixel number of the (T0, beta) interval to the total pixel number of the (alpha, beta) interval, wherein the proportion is respectively as follows:
Figure FDA0002877411370000043
Figure FDA0002877411370000044
calculating (alpha, beta) interval average information entropy:
E2=r0E0+r1E1
if any number greater than 0 exists, when the minimum gray value is alpha and the maximum gray value is beta, the iteration condition is satisfied:
ΔE=|E-E2|<ε;
calculating a threshold value when the variance between the classes of the interested gray scale interval is maximum according to the interested gray scale interval (alpha, beta);
and segmenting the under-segmentation region according to the threshold value.
7. The method for analyzing digital PCR microarray image according to claim 4, wherein the identifying the effective units of the segmented image and extracting the coordinates of the centroid positions of the effective units are specifically:
according to the upper and lower limits of the area threshold value, t is respectively1And t2Dividing the region-of-interest segmentation map into a small area region, a normal area region and a large area region, wherein the area is smaller than t1Is divided into a small area A1 with an area t1And t2Is drawn as normal area region a2,area greater than t2Is divided into a large area A3;
setting a circularity threshold, and merging regions which do not meet the circularity threshold in the region A2 into the region A1;
acquiring a connected region around the region A3, merging the connected region into the region a1 to form a new region a 12;
excluding the connected region around region A3 from region a2, resulting in region a 22;
carrying out edge detection on the preprocessed first channel fluorescence image by using a canny operator, and combining the region A12 with an edge detection result to obtain a region A13;
expanding the area A12, and performing logical AND operation with the area A13 to form an area A14;
filling and opening the area A14 to obtain an area A15;
setting an eccentricity threshold, and removing a connected region which does not meet the eccentricity threshold range in the region A15 to obtain a region A16;
merging the region a16 with the region a22 to obtain an image a 4;
calculating the average gray value of all the areas in the area A22; setting the upper and lower limits of a gray threshold value to be k respectively according to the average gray value1、k2The gray-level values in all connected domains in the image A4 are less than k1And is greater than k2The gray value of the pixel position of (a) is set to 0, and an image A41 is obtained;
calculating the sizes of holes of all connected domains in the image A41; setting a hole threshold, and removing a connected domain which does not meet the hole threshold in the image A41 to obtain an image A42;
calculating the circularity sizes of all connected regions in the image A42, and removing the connected regions which are not in the threshold range according to a circularity threshold to obtain an image A43;
removing the connected domains with overlarge and undersize areas in the image A43 to complete the identification of the effective units;
and marking the centroid of the finally identified connected domain on the gray scale image of the first channel fluorescence image.
8. The digital PCR microarray image analysis method of claim 1, wherein the effective units in the second and third channel images are located and extracted according to the effective unit centroid position coordinates, and the statistics of the fluorescence signal intensities of the effective units in the second and third channel images according to the effective units in the second and third channel images specifically comprises:
and counting the number and the centroid positions of the effective units in the first channel image, and reading the fluorescence signal intensities of the effective units in the second channel image and the third channel image.
9. The method for analyzing digital PCR microarray image of claim 1, wherein the scatter plot and the histogram of the fluorescence intensity distribution of the effective units in the second channel image and the third channel image are plotted according to the fluorescence signal intensity of the effective units in the second channel image and the third channel image, and the ratio of the positive units is counted as follows:
drawing a histogram according to the fluorescence signal intensities of all effective units of the second channel image and the third channel image;
searching two extreme point coordinates with maximum peak span in the histogram, taking the mean value of the two extreme point coordinates as a positive judgment threshold, judging as a positive unit if the fluorescence intensity is greater than the threshold, and otherwise, judging as a negative unit;
and counting the number of positive units and the total number of effective units according to the judgment threshold and the fluorescence signal intensity of the effective units, and drawing a scatter diagram, a histogram and a two-parameter scatter diagram of the fluorescence intensity of the images of the second channel and the third channel.
10. The method of claim 1, wherein the concentration of the target nucleic acid molecule in the original reaction solution is calculated based on the proportion of the positive units and the poisson distribution, and the calculation formula is as follows:
Figure FDA0002877411370000071
wherein p represents the ratio of the number of positive units to the total number of valid units, NtThe total number of wells in the chip was indicated, and V was the total volume of the reaction solution.
CN202011638314.6A 2020-12-31 2020-12-31 Digital PCR microarray image analysis method Withdrawn CN112927182A (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011638314.6A CN112927182A (en) 2020-12-31 2020-12-31 Digital PCR microarray image analysis method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011638314.6A CN112927182A (en) 2020-12-31 2020-12-31 Digital PCR microarray image analysis method

Publications (1)

Publication Number Publication Date
CN112927182A true CN112927182A (en) 2021-06-08

Family

ID=76163277

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011638314.6A Withdrawn CN112927182A (en) 2020-12-31 2020-12-31 Digital PCR microarray image analysis method

Country Status (1)

Country Link
CN (1) CN112927182A (en)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116051556A (en) * 2023-03-31 2023-05-02 至美时代生物智能科技(北京)有限公司 Micro-fluidic chip reaction hole image recognition method and system based on relative coordinates
CN117851957A (en) * 2024-03-07 2024-04-09 徐州市检验检测中心 PCR (polymerase chain reaction) detection fluorescence intensity data processing system

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116051556A (en) * 2023-03-31 2023-05-02 至美时代生物智能科技(北京)有限公司 Micro-fluidic chip reaction hole image recognition method and system based on relative coordinates
CN117851957A (en) * 2024-03-07 2024-04-09 徐州市检验检测中心 PCR (polymerase chain reaction) detection fluorescence intensity data processing system
CN117851957B (en) * 2024-03-07 2024-05-17 徐州市检验检测中心 PCR (polymerase chain reaction) detection fluorescence intensity data processing system

Similar Documents

Publication Publication Date Title
CN107316077B (en) Automatic adipose cell counting method based on image segmentation and edge detection
Wang et al. Automated blob detection using iterative Laplacian of Gaussian filtering and unilateral second-order Gaussian kernels
CN103048329B (en) A kind of road surface crack detection method based on active contour model
CN103543277B (en) A kind of blood group result recognizer based on gray analysis and category identification
CN108961230B (en) Identification and extraction method for structural surface crack characteristics
US11538261B2 (en) Systems and methods for automated cell segmentation and labeling in immunofluorescence microscopy
CN111027443B (en) Bill text detection method based on multitask deep learning
CN112927182A (en) Digital PCR microarray image analysis method
CN115861320B (en) Intelligent detection method for automobile part machining information
CN114820625B (en) Automobile top block defect detection method
CN110598030A (en) Oracle bone rubbing classification method based on local CNN framework
CN114549528B (en) Micro-droplet digital PCR (polymerase chain reaction) droplet detection method and system
Bariamis et al. Unsupervised SVM-based gridding for DNA microarray images
CN110517273B (en) Cytology image segmentation method based on dynamic gradient threshold
CN110148126B (en) Blood leukocyte segmentation method based on color component combination and contour fitting
CN113012124B (en) Shoe print hole and embedded object feature detection and description method
CN109147932B (en) Cancer cell HER2 gene amplification analysis method and system
CN115294377A (en) System and method for identifying road cracks
CN116978543A (en) Artificial intelligent auxiliary marrow tumor pathological diagnosis device
CN113850792A (en) Cell classification counting method and system based on computer vision
CN115115939B (en) Remote sensing image target fine-grained identification method based on characteristic attention mechanism
CN116402822A (en) Concrete structure image detection method and device, electronic equipment and storage medium
CN115082379A (en) Activated sludge phase contrast microscopic image floc and filamentous bacterium segmentation method
CN115170507A (en) Grouting pipe surface defect detection method and system based on image data
CN111815591B (en) Lung nodule detection method based on CT image

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
WW01 Invention patent application withdrawn after publication

Application publication date: 20210608

WW01 Invention patent application withdrawn after publication