CN107844683B - Method for calculating concentration of digital PCR (polymerase chain reaction) liquid drops - Google Patents
Method for calculating concentration of digital PCR (polymerase chain reaction) liquid drops Download PDFInfo
- Publication number
- CN107844683B CN107844683B CN201710966704.8A CN201710966704A CN107844683B CN 107844683 B CN107844683 B CN 107844683B CN 201710966704 A CN201710966704 A CN 201710966704A CN 107844683 B CN107844683 B CN 107844683B
- Authority
- CN
- China
- Prior art keywords
- image
- candidate
- droplet
- liquid drop
- calculating
- 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.)
- Active
Links
- 239000007788 liquid Substances 0.000 title claims abstract description 106
- 238000000034 method Methods 0.000 title claims abstract description 22
- 238000007847 digital PCR Methods 0.000 title claims abstract description 12
- 238000003752 polymerase chain reaction Methods 0.000 title description 6
- 108090000623 proteins and genes Proteins 0.000 claims abstract description 24
- 238000005516 engineering process Methods 0.000 claims abstract description 7
- 238000007781 pre-processing Methods 0.000 claims abstract description 6
- 238000004364 calculation method Methods 0.000 claims description 13
- 238000001914 filtration Methods 0.000 claims description 11
- 238000004891 communication Methods 0.000 claims description 10
- 238000003708 edge detection Methods 0.000 claims description 10
- 230000002146 bilateral effect Effects 0.000 claims description 9
- 230000000877 morphologic effect Effects 0.000 claims description 9
- 238000006243 chemical reaction Methods 0.000 claims description 5
- 238000004458 analytical method Methods 0.000 claims description 4
- 238000000265 homogenisation Methods 0.000 claims description 3
- 238000003709 image segmentation Methods 0.000 claims description 3
- 230000009466 transformation Effects 0.000 claims description 3
- 238000001514 detection method Methods 0.000 abstract description 7
- 238000010586 diagram Methods 0.000 description 9
- 239000000523 sample Substances 0.000 description 4
- 238000000605 extraction Methods 0.000 description 2
- 230000001965 increasing effect Effects 0.000 description 2
- 238000003909 pattern recognition Methods 0.000 description 2
- 238000005070 sampling Methods 0.000 description 2
- 230000004544 DNA amplification Effects 0.000 description 1
- 238000012408 PCR amplification Methods 0.000 description 1
- 238000005260 corrosion Methods 0.000 description 1
- 230000007797 corrosion Effects 0.000 description 1
- 230000002708 enhancing effect Effects 0.000 description 1
- 239000007850 fluorescent dye Substances 0.000 description 1
- 239000012634 fragment Substances 0.000 description 1
- 238000010606 normalization Methods 0.000 description 1
- 238000005191 phase separation Methods 0.000 description 1
- 239000000126 substance Substances 0.000 description 1
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 1
Images
Classifications
-
- G06T5/70—
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/13—Edge detection
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/155—Segmentation; Edge detection involving morphological operators
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
- G06V40/00—Recognition of biometric, human-related or animal-related patterns in image or video data
- G06V40/10—Human or animal bodies, e.g. vehicle occupants or pedestrians; Body parts, e.g. hands
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
- G16B20/00—ICT specially adapted for functional genomics or proteomics, e.g. genotype-phenotype associations
Abstract
The invention discloses a method for calculating the concentration of a digital PCR liquid drop, which comprises the following steps: s1, acquiring a digital image; s2, image preprocessing: image quality is improved by utilizing image processing technologies such as image enhancement and the like; s3, droplet positioning: determining the position of each droplet in the image; s4, calculating the volume of the liquid drop: calculating the volume of each liquid drop according to the area of each liquid drop, and further obtaining the total volume V; s5, calculating the concentration of the liquid drop: substituting the total volume V into a formula to calculate the target gene concentration m:
wherein q is the negative rate of the target gene. The invention can effectively detect the quantity of the liquid drops and calculate the volume of the liquid drops, and realizes the detection of the concentration of the target gene, thereby successfully realizing the digital PCR gene detection.
Description
Technical Field
The invention relates to the technical field of image processing, in particular to a method for identifying liquid drops in a liquid drop picture by gene detection and calculating target gene concentration.
Background
The PCR technology is used for DNA amplification in gene detection, and the invention is directed to a liquid drop digital PCR system, which is developed to the third generation PCR technology-digital PCR. The PCR reaction system is obtained by oil-water two-phase separation and takes liquid drops as a unit, each liquid drop is an independent PCR reaction unit, and after PCR amplification is carried out, the liquid drops can be marked by a fluorescent probe for marking a specific gene fragment. Droplets with a fluorescent signal are positive droplets and droplets without a fluorescent signal are referred to as negative droplets. And (3) photographing by using an industrial high-definition camera, detecting the number and total number of the positive droplets and the volume of each droplet by using an image pattern recognition technology, and finally calculating the target gene concentration of the sample. In the droplet chip, the droplet density is large, the droplet size is different, the background and the background of the negative droplet are not greatly distinguished, and the accurate positioning of the droplet number and the calculation of the droplet volume are important conditions for the successful detection of the digital PCR gene.
Disclosure of Invention
The invention provides a method for calculating the concentration of digital PCR liquid drops, which can effectively detect the quantity of the liquid drops and calculate the volume of the liquid drops and the concentration of a target gene.
A method of calculating a digital PCR droplet concentration, comprising:
s1, acquiring a digital image;
s2, image preprocessing: image quality is improved by utilizing image processing technologies such as image enhancement and the like;
s3, droplet positioning: determining the position of each droplet in the image;
s4, calculating the volume of the liquid drop: calculating the volume of each liquid drop according to the area of each liquid drop, and further obtaining the total volume V;
s5, calculating the concentration of the liquid drop: substituting the total volume V into a formula to calculate the target gene concentration m:
wherein q is the negative rate of the target gene.
The following are preferred technical schemes of the invention:
the step S2 includes:
s21, reading in images of various formats in a file reading mode;
s22, performing graying processing on the read image by using a corresponding conversion formula, and converting the read image into a grayscale image;
s23, carrying out median filtering and denoising processing on the image;
s24, carrying out image homogenization treatment of top hat transformation on the image;
and S25, performing bilateral filtering fuzzy operation and edge enhancement operation on the image to finish the preprocessing of the image.
The step S3 includes:
s31, carrying out image segmentation on the image to obtain a droplet generation area image;
s32, carrying out canny edge detection on the image of the liquid drop generation area;
s33, performing multiple morphological operations on the image after edge detection;
s34, performing connected domain analysis on the image after each morphological operation to obtain all possible candidate droplet positions;
and S35, for each possible candidate droplet position, eliminating false candidate droplet positions of which the candidate droplet center is located in the background region, the edge region and the plurality of droplet centers are distributed on the same droplet by utilizing the surrounding neighborhood region, and obtaining the final droplet position.
The step S35 includes:
s351, expanding the center position of each possible candidate droplet to the left, the right, the upper and the lower respectively to obtain a partial image with the size of W x H;
s352, removing candidate liquid drops of which the center positions are located in the background area;
s353, removing candidate liquid drops of which the center positions are located in the edge areas of the liquid drops;
and S354, removing the candidate liquid drop with the center position of the liquid drop on the same liquid drop.
The step S352 includes:
s3521, carrying out OTSU threshold operation on the local image to obtain a threshold ThrdT;
s3522, when the gray value of the droplet center position is smaller than the threshold ThrdT, the candidate droplet center is considered to be located in the background area, and the candidate droplet center is removed.
The step S353 includes:
s3531, calculating pixel values of 8 neighborhood pixels of the residual candidate liquid drop center positions by utilizing 4 neighborhoods;
s3532, when one of the pixel values of 8 neighborhood pixels is smaller than a threshold value ThrdT, the center of the candidate liquid drop is considered to be positioned in the edge area, and the candidate liquid drop is removed;
the step S354 includes:
s3541, searching a connected area where the background is located for the candidate liquid drops which are not located in the background area and not located in the edge area;
s3542, finding the maximum gray value of the image;
s3543, searching all critical threshold values with disconnected foreground and background;
s3544, judging whether the processed candidate liquid drops are communicated with other candidate liquid drops in the binarization process under the condition of a critical threshold value;
s3545, removing the background communication area without participating in subsequent calculation;
s3546, removing a small connected region without participating in subsequent calculation;
s3547, filling a small communication area to avoid interference on subsequent calculation;
s3548, cutting off the connection part of the candidate liquid drop to be processed and other liquid drops;
s3549, extracting the boundaries of the candidate droplets to be processed by using morphological operations;
s3550, radius calculation in 360 directions is carried out on the extracted boundary;
s3551, for every two candidate droplets, calculating the distance between the two candidate droplets, and if the distance between the two candidate droplets is smaller than the sum of the radii of the two candidate droplets, considering that the two candidate droplets are positioned on the same droplet, and removing one of the candidate droplets.
The step S3550 includes:
s35501, scattering the boundary image in the 360-degree direction from the central position pt 1;
s35502, for each direction, recording a first point pt2 that changes from 0 to 255;
s35503, using formula (7), calculates the euclidean distance between pt1 and pt2 as the radius r in that direction:
wherein (x)1,y1) Is the coordinate of point pt1, (x)2,y2) And (4) obtaining the radius value of each droplet in the 360-degree direction by taking the coordinates of pt2 as r-D.
The step S3551 includes:
s35511, for every two candidate droplets, calculating the direction dir of the two candidate droplets;
s35512, calculating the Euclidean distance D1 between the two candidate droplets;
s35513, obtaining radii r1, r2, r1 of the two droplets in the dir direction respectively, wherein r2 is a radius value of the center position of one candidate droplet in the dir direction, and r2 is a radius value of the center position of the other candidate droplet in the dir direction;
s35514, comparing the size between D1 and r1+ r2, if D1< r1+ r2, then the two candidate drops are positioned on the same drop, and removing the candidate drop position with smaller gray value;
and S35515, removing all overlapped candidate droplet positions to obtain final positioning, and completing positioning of droplets.
The step S4 includes:
s41, performing up-down and left-right expansion on each droplet position to obtain a local image with the size of W x H;
s42, carrying out binarization operation by using OTSU;
s43, carrying out OTSU binarization operation again on the background area to obtain a threshold value ThrdV;
s44, performing binarization operation on the local image by using a threshold value ThrdV;
s45, setting the pixel value of the adhesion position of the binary image to 0 by using the information of adhesion of the liquid drop and other liquid drops obtained in the above way;
s46, extracting the boundary of the binary image;
s47, calculating the number of the points of the connected domain where the liquid drop is positioned, and adding the number of the points on the boundary to be used as the area A of the liquid drop;
s48, calculating the volume v of the droplet by using the formula (8)i：
vi＝A*H (8)
Wherein H is height;
s49, calculating the total volume of the target genes by using the formula (9):
wherein n is the total number of the liquid drops;
thereby obtaining the total volume V.
Compared with the prior art, the invention has the following advantages:
and (3) photographing by using an industrial high-definition camera, detecting the number and total number of the positive droplets and the volume of each droplet by using an image pattern recognition technology, and finally calculating the target gene concentration of the sample. In the liquid drop chip, the liquid drop density is large, the liquid drop size is different, and the background of the liquid drop are not greatly distinguished, the invention can accurately position the liquid drop quantity and calculate the liquid drop volume through image processing, and can effectively detect the liquid drop quantity and calculate the liquid drop volume and the target gene concentration m, thereby successfully realizing the digital PCR gene detection.
Drawings
Fig. 1 is an image of the graying processing of step S22;
FIG. 2 is the bilateral filtered image of step S25;
FIG. 3 is the image after the edge enhancement in step S252;
fig. 4 is an image after the edge detection of canny in step S32;
fig. 5 is an image after the first expansion in step S332;
fig. 6 is the image after the second expansion in step S332;
fig. 7 is an image after the last expansion in step S332;
FIG. 8 shows all possible candidate drop positions of step S34;
fig. 9 is a partial image of the center of one droplet expanded from left to right and up and down in step S351;
fig. 10 is a partial image binarized at step S35411;
FIG. 11 is a diagram illustrating the step S3522 after the background region is removed from the center of the droplet candidate;
FIG. 12 is a schematic diagram of step S3532 after the center of the droplet candidate at the edge region is removed;
fig. 13 is a histogram of the partial image of step S35421;
FIG. 14 is a schematic diagram of the local image when the foreground and the background are not connected in step S35435;
fig. 15 is a partial image after the small connected region is filled in step S35471;
fig. 16 is a partial image in which the droplet joint is removed in step S35485;
fig. 17 is a partial image after the boundary extraction in step S35502;
FIG. 18 is a schematic diagram of step S35514 in which the centers of the candidate droplets on the same droplet are removed;
fig. 19 is a schematic diagram of calculating 8 neighborhood pixels for the remaining candidate droplet center positions using 4 neighborhoods in step S3531.
Detailed Description
The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the drawings in the embodiments of the present invention, and it is obvious that the described embodiments are only a part of the embodiments of the present invention, and not all of the embodiments. All other embodiments, which can be derived by a person skilled in the art from the embodiments given herein without making any creative effort, shall fall within the protection scope of the present invention.
A method of calculating a digital PCR droplet concentration, comprising:
s1, acquisition of digital images: acquiring a digital image of a target area by using a photosensitive element such as a CCD (charge coupled device);
s2, image preprocessing: image quality is improved by utilizing image processing technologies such as image enhancement and the like;
s3, droplet positioning: determining the position of each droplet in the image;
s4, calculating the volume of the liquid drop: calculating the volume of each liquid drop according to the area of each liquid drop, and further obtaining the total volume;
s5, calculating the concentration of the liquid drop: and calculating the concentration of the target gene.
The step S2 includes:
s21, reading in images of various formats in a file reading mode;
s22, converting the read image into a gray image by performing a graying process using a corresponding conversion formula, wherein the step S22 includes:
the read-in image is grayed according to a corresponding formula, such as: the RGB color image can be grayed out as follows:
GrayV＝0.2989*R+0.5870*G+0.1140*B(1)
wherein, R, G, B are R, G, B components of the color image respectively, GrayV is the converted gray value;
the grayed image is shown in fig. 1;
s23, performing median filtering denoising on the gray-scale image, wherein the step S23 includes:
s231, performing sliding operation on each pixel point in the image according to a sliding frame of n × m, wherein the values of n and m are 7 for example;
s232, sorting all pixel values in the sliding frame (sorting from large to small or sorting from small to large, and taking a middle value);
s233, replacing the gray value of the position where the center of the sliding frame is located in the image with the sorted median value;
s24, performing top hat transformation (image homogenization processing) on the image, wherein the step S24 includes:
s241, performing opening operation of firstly corroding and then expanding the image by using the structural element n x n, wherein the value of n is 50 for example;
s242, subtracting the image after the opening operation from the original image to obtain a top-hat converted image;
s25, carrying out bilateral filtering fuzzy operation and edge enhancement operation on the image, wherein the step S25 comprises the following steps:
s251, performing bilateral filtering operation on each liquid drop generation area image, and keeping boundary information while removing noise, wherein the range of a space domain Gaussian kernel is [5,15], and the range of a color domain Gaussian kernel is [10,20 ];
the bilateral filtered image is shown in fig. 2;
s252, performing a sharpening operation on the image after the bilateral filtering operation, and enhancing edge information, where the step S252 includes:
s2521, high utilizationThe gaussian filtering performs gaussian blurring on the bilateral filtered image, and the template used in the gaussian blurring in step S2521 is a pass function
When δ is 1, sampling is performed with respect to the center, and the resulting filter template is as follows:
in the formula, x and y are x and y coordinate values of the image respectively, and sigma refers to variance.
S2522, subtracting the image after Gaussian blur from the image after bilateral filtering to obtain a difference image, and recording the difference image as a template, namely a template image;
s2523, adding the template image to the original image to obtain a sharpened image, wherein step S2523 includes:
the sharpened image is calculated using the following formula:
g(x，y)＝f(x，y)+k*gmask(x，y) (3)
wherein g (x, y) is the sharpened image, f (x, y) is the original image, k (k > 0) is the weighting factor, high boost filtering is performed when k > 1, and the contribution of unsharp template is de-emphasized when k < 1, where k is in the range [5,20 ]],gmask(x, y) is a template image, which is obtained by the following formula:
wherein the content of the first and second substances,
the image is a Gaussian blurred image;
the edge enhanced image is shown in fig. 3.
The step S3 includes:
s31, carrying out image segmentation on the image to obtain a region of interest, wherein the step S31 comprises the following steps:
s311, positioning the sample inlet hole and the oil storage hole to obtain a rough processing area;
s312, removing the sampling hole and the oil storage hole area;
s313, finding the position of the flow channel by projection;
s314, dividing the region to be processed into two liquid drop generating regions through the flow channel, namely the region of interest;
s32, carrying out edge detection on the preprocessed image, wherein the step S32 comprises the following steps:
for the preprocessed image, edge detection is carried out by using a canny edge detection operator, and two threshold value ranges in the detection operator are respectively [30,80], [90,240 ];
the image after canny edge detection is shown in fig. 4;
s33, performing multiple morphological operations on the edge-detected image, wherein the step S33 includes:
s331, performing multiple cycle expansion operations on the image after edge detection, and terminating the cycle when a new connected region is not generated any more;
s332, an expansion operation, using a rectangular structure of 3 × 3, the structure being:
the image after the first expansion is shown in fig. 5, the image after the second expansion is shown in fig. 6, and the image after the last expansion is shown in fig. 7;
s34, performing connected domain analysis on the image after each morphological operation to obtain all possible candidate droplet positions, wherein the step S34 includes:
s341, analyzing the connected domain of the image after the expansion operation, recording the information of each connected region, and obtaining all possible candidate droplet positions, wherein 4-way communication is used for operation;
all possible candidate drop positions are shown in fig. 8;
s35, for each possible candidate droplet position, using its surrounding neighborhood region to eliminate the false candidate droplet positions where the candidate droplet center is located in the background region, the edge region and the multiple droplet centers are distributed on the same droplet, and obtain the final droplet position, where the step S35 includes:
s351, expanding the center position of each possible candidate droplet to the left, the right, the upper and the lower respectively to obtain a partial image with the size of W x H (the range of W and H is 100-200, and is determined according to specific situations);
the partial image of one of the liquid drops after the center is expanded from left to right and up and down is shown in fig. 9;
s352, removing the candidate liquid drop with the center position of the liquid drop positioned in the background area, wherein the step S352 comprises the following steps:
s3521, performing OTSU threshold operation on the local image (OTSU refers to a binary method commonly used in image processing and named as Otsu method, and is the prior art) to obtain a threshold ThrdT (ThrdT is a threshold obtained according to OTSU);
s3522, when the gray value of the center position of the droplet is smaller than ThrdT, the candidate droplet center is considered to be located in the background area, and the candidate droplet center is removed;
the schematic diagram of the candidate droplet center after the background region is removed is shown in fig. 11;
s353, removing the candidate liquid drop with the liquid drop center position located in the liquid drop edge area, wherein the step S353 comprises the following steps:
s3531, for the remaining candidate droplet center positions, the pixel values of the 8 neighboring pixels are calculated using the 4 neighbors, and the structure is shown in fig. 19:
s3532, when one of the pixel values of 8 neighborhood pixels is smaller than the threshold thrdT, the center of the candidate liquid drop is considered to be positioned in the edge area, and the candidate liquid drop is removed;
the schematic diagram of the candidate droplet with the center at the edge area removed is shown in fig. 12;
s354, removing the candidate droplet whose droplet center position is located on the same droplet, where the step S354 includes:
s3541, for candidate droplets that are not located in the background region and not in the edge region, finding a connected region where the background is located, where the step S3541 includes:
s35411, carrying out binarization operation on the local image by using a threshold value ThrdT;
the binarized partial image is shown in fig. 10;
s35412, connected domain analysis is conducted on the binary image;
s35413, storing the number of pixels of each connected domain and the position information of the pixels to obtain background information for subsequent use;
s3542, finding the maximum gray value of the image, wherein the step S3542 comprises the following steps:
s35421, calculating a histogram of the local image;
the histogram of the partial image is shown in fig. 13;
s35422, carrying out normalization on the histogram to a range [0,100 ];
s35423, from 255 to a threshold ThrdT, obtaining the gray value position with the gray value number larger than 5 and the number of the following continuous 3 gray values also larger than 5, recording the gray value, and taking the gray value as the maximum gray value maxGrayV of the local image;
s3543, finding a critical threshold value at which all foreground and background are not communicated, wherein the step S3543 includes:
s35431, sequentially carrying out binarization operation of the threshold value from threshold value ThrdT +1 to maxGrayV on the local image, and increasing the threshold value by 1 each time;
s35432, recording the information of the connected domains in each binarization process, and storing the information in an array neighbor InfoRecord, wherein the array neighbor InfoRecord comprises the information of the total number of the connected domains, the number of points of each connected domain, the central position, the label and the like for subsequent use;
s35433, comparing the connected domain information after each binarization with the background connected domain information, and if the number of the connected domains of the two is equal or fluctuates in a small range (the difference of the number of pixels is about 10) and the labels are matched, determining that the background region is separated from the foreground region and is not fused;
s35434, counting the number of the foreground areas and the background areas that are not fused, if the total number is equal to the number of the background areas, then, we consider that the foreground areas and the background areas are not fused under the threshold, find a critical threshold T, and record the position start _ pos of the threshold in the neighborInfoRecord (described in S35431 is a loop operation, each time the threshold is increased by 1, when the condition "the total number is equal to the number of the background areas", we consider that the current threshold is the critical threshold, which is denoted as T, and the start _ pos is a subscript in the neighborInfoRecord of the group), and exit the loop operation;
s35435, recording the residual connected domain information, enabling the threshold to change from a critical threshold T to maxGrayV, recording the connected domain information of the change process, and recording the total number of neighborInfo records as num;
a schematic diagram of a local image when the foreground is not communicated with the background is shown in fig. 14, and the center position is located in the third droplet in the first row;
s3544, determining whether the processed candidate droplet is connected to other candidate droplets in the binarization process under the condition of the critical threshold, wherein the step S3544 includes:
s35441, calculating a label of a communication area where the center of the candidate liquid drop to be processed is located;
s35442, if the tag is non-zero, i.e., it is not in the background region under the threshold T, then determining whether the droplet is communicated with other droplets, wherein the step S35442 includes:
s354421, storing the positions of all the pixels with labels of label;
s354422, traversing the neighborInfo record from start _ pos +1 to num, recording a new connected domain label set of all points with labels being label in the change process, and counting the number of points of each label in the label set and the total number Sum;
s354423, if Sum is larger than 54, more than two labels exist, and the number of the labels with the number of the communication domain points larger than 0.3 Sum is larger than 2, the candidate liquid drops to be processed are considered to be communicated with other liquid drops, and circulating in such a way, and recording the position pos where new communication liquid drops do not appear any more;
s354424, if there is no candidate droplet in communication with other droplets, the final cut threshold final _ T is equal to T, otherwise final _ T is equal to pos-start _ pos + T (i.e., pos minus start _ pos plus T);
s35443, if the label is zero, finding a neighborhood point of which the label is non-zero and which does not fall into the background area with the threshold value of thrdT by using the four neighborhood areas, calculating label of the communication area where the label is located, and judging whether the liquid drop is communicated with other liquid drops by using the neighborhood point;
s3545, removing the background connected region without participating in subsequent calculation, wherein the step S3545 comprises the following steps:
s35451, when the threshold value is final _ T, comparing the neighborInfo record with the background area when the threshold value is thrdT;
s35452, if the number of the connected domain points of the two connected domain points is equal or fluctuates in a small range (the difference of the number of the pixels is about 10), and a label with a background point is matched with a label in the neighborInfo record, the connected domain is regarded as a background domain, and the connected domain is removed;
s3546, removing small connected regions and not participating in subsequent calculation, wherein the step S3546 comprises the following steps:
s35461, if the number of the points of the connected area is less than 27, namely the minimum drop radius is determined to be 3, regarding the remaining connected area as an invalid area, and removing the invalid area;
s3547, filling a small connected area to avoid interference on subsequent calculation, wherein the step S3547 comprises the following steps:
s35471, if the number of the points in a certain connected region is less than 27, setting the pixel values of all the pixel points in the region to be 255;
the partial image after the small connected region is filled is shown in fig. 15, and the center position is located in the third droplet in the first row;
s3548, cutting off the connection part of the candidate liquid drop to be processed and other liquid drops, wherein the step S3548 comprises the following steps:
s35481, scattering the liquid drops to be processed from the center of the liquid drops in the direction of 360 degrees;
s35482, for each direction, find the first point pt1 where 0 changes to 255;
s35483, continue in this direction, find the first point pt2 that changes from 255 to 0;
s35484, if pt2 is located in the background area, in the direction, the drop to be treated is not adhered with other drops;
s35485, if pt2 is located in a non-background area, in the direction, the liquid drops to be processed are adhered to other liquid drops, and the gray value of the adhered pixel points is set to be 0;
the partial image with the droplet junction removed is shown in fig. 16, and the center position is located in the third droplet in the first row;
s3549, extracting the boundary of the candidate liquid drop to be processed by using morphological operation, wherein the step S3549 comprises the following steps:
s35491, traversing the processed binary image, and setting the gray value of the pixel point as 1 if the gray value of the pixel point is equal to 255;
s35492, the utilization structure is:
carrying out corrosion operation on the image;
s35493, subtracting the corroded image from the original image, setting the pixel value of the corresponding pixel point to be 255 if the pixel value of the corresponding pixel point is 1, and setting the other pixel values to be 0 to obtain a boundary image;
s3550, calculating a radius of the extracted boundary in 360 directions, wherein the step S3550 includes:
s35501, scattering the boundary image in the 360-degree direction from the central position pt 1;
s35502, for each direction, recording a first point pt2 that changes from 0 to 255;
the partial image after the boundary extraction is shown in fig. 17, and the center position is located in the third droplet in the first row;
s35503, calculating the euclidean distance between pt1 and pt2 as the radius r in the direction using the following formula:
the formula listed is a general Euclidean distance calculation formula, where (x)1,y1) Is the coordinate of point pt1, (x)2,y2) And (4) obtaining the radius value of each droplet in the 360-degree direction by taking the coordinates of pt2 as r-D.
S3551, calculating the distance between every two candidate droplets, and if the distance between them is smaller than the sum of their radii, then it is considered that they are located on the same droplet, and one of them is removed, and the step S3551 includes:
s35511, for every two candidate droplets, calculating the direction dir of the two candidate droplets;
s35512, calculating the Euclidean distance D1 between the two candidate droplets;
s35513, obtaining radii r1, r2 of the two droplets in the dir direction respectively (the formula in S35503, it is calculated that 360 radius values are supplied in the 360 degree direction, and r1, r2 are radius values in the dir direction, that is, r1 is a radius value of the center position of one candidate droplet in the dir direction, and r2 is a radius value of the center position of the other candidate droplet in the dir direction);
s35514, comparing the size between D1 and r1+ r2, if D1< r1+ r2, then the two candidate drops are positioned on the same drop, and removing the candidate drop position with smaller gray value;
a schematic diagram of the candidate droplet after the center of the candidate droplet is removed on the same droplet is shown in fig. 18;
s35515, removing all overlapped candidate droplet positions to obtain final positioning (completing positioning of droplets);
the step S4 includes:
s41, performing up-down and left-right expansion on each droplet position to obtain a local image with the size of W x H (the range of W and H is 100-200, and the local image is determined according to specific situations);
s42, carrying out binarization operation by using OTSU;
s43, carrying out OTSU binarization operation again on the background area to obtain a threshold ThrdV (ThrdV is a threshold obtained according to OTSU);
s44, performing binarization operation on the local image by using a threshold value ThrdV;
s45, setting the pixel value of the adhesion position of the binary image to 0 by using the information of adhesion of the liquid drop and other liquid drops obtained in the above way;
s46, extracting the boundary of the binary image;
s47, calculating the number of the points of the connected domain where the liquid drop is positioned, and adding the number of the points on the boundary to be used as the area A of the liquid drop;
s48, calculating the volume v of the droplet by using the formula (8)i：
vi＝A*H (8)
Wherein H is height; the H value is currently estimated to be 3/4 × r, r being the average radius of the droplet;
s49, the total volume of the target genes is as follows:
wherein n is the total number of droplets.
The step S5 includes:
s51, calculating the target gene concentration m according to the Poisson distribution of the large sample by using the following formula:
wherein q is the negative rate of the target gene. The value of q is calculated from the images, which are different from one image to the next, ranging from 0 to 1, and V is the aforementioned volume.
Claims (6)
1. A method for calculating the concentration of a digital PCR droplet, comprising:
s1, acquiring a digital image;
s2, image preprocessing: image quality is improved by utilizing image processing technologies such as image enhancement and the like;
s3, droplet positioning: determining the position of each droplet in the image;
s4, calculating the volume of the liquid drop: calculating the volume of each liquid drop according to the area of each liquid drop, and further obtaining the total volume V;
the step S3 includes:
s31, carrying out image segmentation on the image to obtain a droplet generation area image;
s32, carrying out canny edge detection on the image of the liquid drop generation area;
s33, performing multiple morphological operations on the image after edge detection;
s34, performing connected domain analysis on the image after each morphological operation to obtain all possible candidate droplet positions;
s35, for each possible candidate droplet position, eliminating false candidate droplet positions of which the candidate droplet center is located in a background region, an edge region and a plurality of droplet centers are distributed on the same droplet by utilizing a surrounding neighborhood region of the candidate droplet position to obtain a final droplet position;
the step S35 includes:
s351, expanding the center position of each possible candidate liquid drop to the left, the right, the upper and the lower respectively to obtain a local image;
s352, removing candidate liquid drops of which the center positions are located in the background area;
s353, removing candidate liquid drops of which the center positions are located in the edge areas of the liquid drops;
s354, removing candidate liquid drops of which the center positions are located on the same liquid drop;
the step S352 includes:
s3521, carrying out OTSU threshold operation on the local image to obtain a threshold ThrdT;
s3522, when the gray value of the center position of the liquid drop is smaller than a threshold ThrdT, the candidate liquid drop center is considered to be located in the background area, and the candidate liquid drop center is removed;
s5, calculating the concentration of the liquid drop: substituting the total volume V into a formula to calculate the target gene concentration m:
wherein q is the negative rate of the target gene.
2. The method according to claim 1, wherein the step S2 includes:
s21, reading in images of various formats in a file reading mode;
s22, performing graying processing on the read image by using a corresponding conversion formula, and converting the read image into a grayscale image;
s23, carrying out median filtering and denoising processing on the image;
s24, carrying out image homogenization treatment of top hat transformation on the image;
and S25, performing bilateral filtering fuzzy operation and edge enhancement operation on the image to finish the preprocessing of the image.
3. The method according to claim 1, wherein the step S354 comprises:
s3541, searching a connected area where the background is located for the candidate liquid drops which are not located in the background area and not located in the edge area;
s3542, finding the maximum gray value of the image;
s3543, searching all critical threshold values with disconnected foreground and background;
s3544, judging whether the processed candidate liquid drops are communicated with other candidate liquid drops in the binarization process under the condition of a critical threshold value;
s3545, removing the background communication area without participating in subsequent calculation;
s3546, removing a small connected region without participating in subsequent calculation;
s3547, filling a small communication area to avoid interference on subsequent calculation;
s3548, cutting off the connection part of the candidate liquid drop to be processed and other liquid drops;
s3549, extracting the boundaries of the candidate droplets to be processed by using morphological operations;
s3550, radius calculation in 360 directions is carried out on the extracted boundary;
s3551, for every two candidate droplets, calculating the distance between the two candidate droplets, and if the distance between the two candidate droplets is smaller than the sum of the radii of the two candidate droplets, considering that the two candidate droplets are positioned on the same droplet, and removing one of the candidate droplets.
4. The method according to claim 3, wherein the step S3550 comprises:
s35501, scattering the boundary image in the 360-degree direction from the central position pt 1;
s35502, for each direction, recording a first point pt2 that changes from 0 to 255;
s35503, using formula (7), calculates the euclidean distance between pt1 and pt2 as the radius r in that direction:
wherein (x)1,y1) Is the coordinate of point pt1, (x)2,y2) And (4) obtaining the radius value of each droplet in the 360-degree direction by taking the coordinates of pt2 as r-D.
5. The method according to claim 4, wherein the step S3551 comprises:
s35511, for every two candidate droplets, calculating the direction dir of the two candidate droplets;
s35512, calculating the Euclidean distance D1 between the two candidate droplets;
s35513, obtaining radii r1, r2, r1 of the two droplets in the dir direction respectively, wherein r2 is a radius value of the center position of one candidate droplet in the dir direction, and r2 is a radius value of the center position of the other candidate droplet in the dir direction;
s35514, comparing the size between D1 and r1+ r2, if D1< r1+ r2, then the two candidate drops are positioned on the same drop, and removing the candidate drop position with smaller gray value;
and S35515, removing all overlapped candidate droplet positions to obtain final positioning, and completing positioning of droplets.
6. The method according to claim 1, wherein the step S4 includes:
s41, performing up-down and left-right expansion on each liquid drop position to obtain a local image;
s42, carrying out binarization operation by using OTSU;
s43, carrying out OTSU binarization operation again on the background area to obtain a threshold value ThrdV;
s44, performing binarization operation on the local image by using a threshold value ThrdV;
s45, setting the pixel value of the adhesion position of the binary image to 0 by using the information of adhesion of the liquid drop and other liquid drops obtained in the above way;
s46, extracting the boundary of the binary image;
s47, calculating the number of the points of the connected domain where the liquid drop is positioned, and adding the number of the points on the boundary to be used as the area A of the liquid drop;
s48, calculating the volume v of the droplet by using the formula (8)i：
vi＝A*H (8)
Wherein H is height;
s49, calculating the total volume of the target genes by using the formula (9):
wherein n is the total number of the liquid drops;
thereby obtaining the total volume V.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710966704.8A CN107844683B (en) | 2017-10-17 | 2017-10-17 | Method for calculating concentration of digital PCR (polymerase chain reaction) liquid drops |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710966704.8A CN107844683B (en) | 2017-10-17 | 2017-10-17 | Method for calculating concentration of digital PCR (polymerase chain reaction) liquid drops |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN107844683A CN107844683A (en) | 2018-03-27 |
| CN107844683B true CN107844683B (en) | 2020-05-12 |
Family
ID=61661626
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201710966704.8A Active CN107844683B (en) | 2017-10-17 | 2017-10-17 | Method for calculating concentration of digital PCR (polymerase chain reaction) liquid drops |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN107844683B (en) |
Families Citing this family (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110689926A (en) * | 2018-06-19 | 2020-01-14 | 上海交通大学 | Accurate detection method for high-throughput digital PCR image droplets |
| CN108956400A (en) * | 2018-07-20 | 2018-12-07 | 润电能源科学技术有限公司 | A kind of Measure Method of Pulverized Coal and device based on image recognition |
| CN109358026B (en) * | 2018-09-13 | 2021-04-06 | 中国科学院苏州生物医学工程技术研究所 | Fluorescent liquid drop detection method and device and server |
| US11305284B2 (en) * | 2018-11-26 | 2022-04-19 | Tokitae, LLC | Determining a bulk concentration of a target in a sample using a digital assay with compartments having nonuniform volumes |
| CN109903296B (en) * | 2019-02-15 | 2021-06-01 | 领航基因科技（杭州）有限公司 | LBP-Adaboost algorithm-based digital PCR (polymerase chain reaction) droplet detection method |
| CN111951294A (en) * | 2019-12-04 | 2020-11-17 | 上海孚清生物科技有限公司 | Micropore traversing method and device of PCR chip and storage medium |
| CN112766075B (en) * | 2020-12-31 | 2021-12-28 | 中国冶金地质总局矿产资源研究院 | Hyperspectral remote sensing black and odorous water body grading method based on semi-supervised learning strategy |
| CN114262733A (en) * | 2022-01-10 | 2022-04-01 | 深圳麦科田生物医疗技术股份有限公司 | Micro-drop digital PCR (polymerase chain reaction) fluorescent signal processing method |
| CN114540469A (en) * | 2022-01-11 | 2022-05-27 | 深圳大学 | Digital nucleic acid quantification method based on non-uniform volume liquid drops and image processing |
| CN114486841B (en) * | 2022-04-07 | 2022-06-24 | 南京长健生物科技有限公司 | Micro fluorescence immunoassay analyzer and image processing method thereof |
| CN114638832B (en) * | 2022-05-19 | 2022-09-23 | 深圳市中科先见医疗科技有限公司 | DPCR liquid drop fluorescence detection method based on watershed algorithm |
Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101887066A (en) * | 2009-05-15 | 2010-11-17 | 中国科学院空间科学与应用研究中心 | System and method for measuring droplet evaporation rate |
-
2017
- 2017-10-17 CN CN201710966704.8A patent/CN107844683B/en active Active
Patent Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101887066A (en) * | 2009-05-15 | 2010-11-17 | 中国科学院空间科学与应用研究中心 | System and method for measuring droplet evaporation rate |
Non-Patent Citations (2)
| Title |
|---|
| 基于图像处理的黑液浓度测量;闫友房;《中华纸业》;20111231;全文 * |
| 自吸离散式数字PCR芯片检测装置的关键技术研究与系统开发;袁茂恺;《中国优秀硕士学位论文全文数据库信息科技辑》;20160815;第39页第1段-第51页第2段，图3.5-3.16 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN107844683A (en) | 2018-03-27 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN107844683B (en) | Method for calculating concentration of digital PCR (polymerase chain reaction) liquid drops | |
| CN107545239B (en) | Fake plate detection method based on license plate recognition and vehicle characteristic matching | |
| CN106651872B (en) | Pavement crack identification method and system based on Prewitt operator | |
| Yousif et al. | Toward an optimized neutrosophic K-means with genetic algorithm for automatic vehicle license plate recognition (ONKM-AVLPR) | |
| Lee et al. | Automatic extraction of characters in complex scene images | |
| US6587586B1 (en) | Extracting textual information from a video sequence | |
| Yoon et al. | Blob extraction based character segmentation method for automatic license plate recognition system | |
| CN107665327B (en) | Lane line detection method and device | |
| CN109460735B (en) | Document binarization processing method, system and device based on graph semi-supervised learning | |
| CN111489337B (en) | Automatic optical detection pseudo defect removal method and system | |
| Bedruz et al. | Real-time vehicle detection and tracking using a mean-shift based blob analysis and tracking approach | |
| Ten Kate et al. | Mid-range and distant vehicle detection with a mobile camera | |
| CN110598566A (en) | Image processing method, device, terminal and computer readable storage medium | |
| CN110648330B (en) | Defect detection method for camera glass | |
| CN108509950B (en) | Railway contact net support number plate detection and identification method based on probability feature weighted fusion | |
| Azad et al. | New method for optimization of license plate recognition system with use of edge detection and connected component | |
| Azad et al. | A novel and robust method for automatic license plate recognition system based on pattern recognition | |
| CN111898627A (en) | SVM cloud microparticle optimization classification recognition method based on PCA | |
| US20080310685A1 (en) | Methods and Systems for Refining Text Segmentation Results | |
| CN112580447A (en) | Power line detection method based on edge second-order statistics and fusion | |
| US9332154B2 (en) | Image binarization using dynamic sub-image division | |
| Aung et al. | Automatic license plate detection system for myanmar vehicle license plates | |
| CN112289377A (en) | Method, apparatus and computer program product for detecting bright spots on an image | |
| CN110084117B (en) | Document table line detection method and system based on binary image segmentation projection | |
| Yoon et al. | Blob detection and filtering for character segmentation of license plates |
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 | ||
| CB02 | Change of applicant information |
Address after: Room 1114, No. 341, Jin Jun Road, Jianggan District, Jianggan District, Hangzhou, Zhejiang Applicant after: PILOT GENE TECHNOLOGIES(HANGZHOU) Co.,Ltd. Address before: Room 1114, No. 341, Jin Jun Road, Jianggan District, Jianggan District, Hangzhou, Zhejiang Applicant before: HUADONG MEDICINE (HANGZHOU) GENE TECHNOLOGY Co.,Ltd. |
|
| CB02 | Change of applicant information | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant | ||
| TR01 | Transfer of patent right |
Effective date of registration: 20211216 Address after: 518000 g1316, Lianxing building, building B, Yihua new village, district 46, Haifu community, Xin'an street, Bao'an District, Shenzhen, Guangdong Province Patentee after: Pilot medical technology (Shenzhen) Co.,Ltd. Address before: 310000 room 1114, Jin Jun Road, 341 Shui Xiang Road, Jianggan District, Hangzhou, Zhejiang. Patentee before: PILOT GENE TECHNOLOGIES (HANGZHOU) Co.,Ltd. |
|
| TR01 | Transfer of patent right |