Automatic segmentation and measurement methods of living stomata of plants based on the CV model

Background The stomata of plants mainly regulate gas exchange and water dispersion between the interior and external environments of plants and play a major role in the plants’ health. The existing methods of stomata segmentation and measurement are mostly for specialized plants. The purpose of this research is to develop a generic method for the fully automated segmentation and measurement of the living stomata of different plants. The proposed method utilizes level set theory and image processing technology and can outperform the existing stomata segmentation and measurement methods based on threshold and skeleton in terms of its versatility. Results The single stomata images of different plants were the input of the method and a level set based on the Chan-Vese model was used for stomatal segmentation. This allowed the morphological features of the stomata to be measured. Contrary to existing methods, the proposed segmentation method does not need any prior information about the stomata and is independent of the plant types. The segmentation results of 692 living stomata of black poplars show that the average measurement accuracies of the major and minor axes, area, eccentricity and opening degree are 95.68%, 95.53%, 93.04%, 99.46% and 94.32%, respectively. A segmentation test on dayflower (Commelina benghalensis) stomata data available in the literature was completed. The results show that the proposed method can effectively segment the stomata images (181 stomata) of dayflowers using bright-field microscopy. The fitted slope of the manually and automatically measured aperture is 0.993, and the R2 value is 0.9828, which slightly outperforms the segmentation results that are given in the literature. Conclusions The proposed automated segmentation and measurement method for living stomata is superior to the existing methods based on the threshold and skeletonization in terms of versatility. The method does not need any prior information about the stomata. It is an unconstrained segmentation method, which can accurately segment and measure the stomata for different types of plants (woody or herbs). The method can automatically discriminate whether the pore region is independent or not and perform pore region extraction. In addition, the segmentation accuracy of the method is positively correlated with the stomata’s opening degree.


Background
Stomata are structures that are found in the aboveground parts of all terrestrial plants and account for approximately 95% of gas exchange [1]. The stomata are composed of pores formed by a pair of guard cells. Stomata regulate the exchange of water vapour and CO 2 between the plant and the atmosphere through changes in the aperture of the pores. They play a pivotal role in controlling the balance between the water loss and carbon gain [2][3][4][5][6]. Therefore, research on stomatal behaviour has been a hot topic in the field of botany. To study the behaviour of stomata, one must be able to calculate the morphology of the pores and quantitatively describe the behaviour of the pores.
There are many methods to measure the morphology of the pores. Omasa and Onoe [7] first conducted stomata measurement studies on sunflowers (Helianthus annuus L cv. Russian Mammoth) using Fourier transform filtering to remove the scanning lines and an unsharp masking and thresholding method for the pore region extraction. The parameters of the pores (major and minor axes and area) are obtained by translating and rotating the principal components of the extraction region. This method performs the transformation in the frequency domain of the image, and the calculation costs are large. In addition, the threshold adjustment needs to be performed manually according to the different stomata images, and it is not an automatic segmentation method. Sanyal et al. [8] used the watershed method to perform the stomatal segmentation of SEM images, and used their morphology to remove noise and the Sobel operator to extract the stomatal edges. This is an edge-based method that performs poorly when the edges of the pores are discontinuous or there is more noise. Laga et al. [9] used template matching technology to identify wheat stomata, but the calculation of the stomatal morphological parameters is complex, and the template parameters are only valid for wheat stomata. New templates need to be constructed when handling the stomata of other species. Therefore, the algorithm is less versatile. Liu [10] uses the MSERs (maximal stable external regions) to detect ellipses in grape stomatal images and identify stomata, but this method requires manual interaction when selecting the ellipses, which is a semi-automatic method that assumes that the stomata shape is elliptical. Jayakody [11] proposed a method for automatic pore detection and measurement for grapevines. It combined the threshold with the morphological skeleton to extract the stomatal edges and calculate the parameters for ellipse fitting such as the major and minor axis of the pores and the eccentricity. The algorithm assumes that the stomata are close to the centre of the region of interest (ROI) and that the length of the skeleton remnant needs to be specified. If the pore size (different species) changes, the algorithm will fail to identify the pores. Therefore, the algorithm is only applicable to grape species. Toda et al. [12] used facial recognition technology to perform automatic pore measurement of the stomata for dayflowers (Commelina benghalensis). This method uses adaptive threshold technology to calculate the pore parameters. However, the algorithm needs to manually define the ranges of the parameters (such as the area, solidity, major axis length, and centroid coordinates). When the sizes or shapes of the stomata change, the algorithm will fail to identify the pores. Therefore, the above methods are only applicable to a specific plant, such as sunflowers, grapes, wheat or dayflowers, and the data collection process requires picking the leaves and tearing off the epidermis to make a micro-slide and acquire the images using bright-field microscopy. Once the pore size and the type or the conditions of image acquisition are substantially changed, the method mentioned above will not successfully identify the pores.
In this paper, we aim to develop a general method for the automatic segmentation and measurement of plant stomata, which is a level set method based on the Chan-Vese (CV) model [13]. Pore region extraction is accomplished by evolving the energy function, and the pores' morphological parameters are obtained using the ellipse fitting technique. The method does not rely on image gradient information and can overcome the edge leakage and poor anti-noise ability based on the edge information model. It is suitable for segmentation and measurement of any plant stomata and does not require prior information about the stomata. In addition, this is the first time that research has been conducted on the behaviour of living stomata of plants in situ.
The remainder of this paper is organized as follows. In the "Methods" section, detailed steps and examples of stomata segmentation and pore measurement are given. The experimental results of the study and comparisons with existing methods are presented in the "Results" section. The discussion and conclusion are given in the last section.

Methods
The purpose of this study was to develop a general automatic stomata segmentation and pore measurement method. It requires the input of a single stomata image (manual or automatic detection) and outputs the morphology parameters of the pores. The method comprises 5 steps: (1) detect and crop a single stomata image as the input, (2) convert the image to greyscale, (3) conduct level set segmentation based on the CV model (with the ROI image centre as the initial position), (4) conduct region shape analysis, and (5) conduct post-processing and ellipse fitting.
The overall flow of the method is shown in Fig. 1. Each 1000× microscope image contains approximately 10-20 stomata, and each 500× magnification image contains approximately 30-40 stomata. First, the stomata are detected using the Faster-RCNN [14], and then the detected stomata are cropped to generate a single stomata image. Second, the single stomata image is used as an input for the stomata segmentation algorithm, and the level set segmentation based on the CV model is implemented (the ROI centre is taken as the initial position).The output is the stomata pore region. Next, shape analysis is performed on the binary image of the pore region to determine whether the segmentation result is an isolated region or not. If it is not an isolated region, morphological post-processing will be executed in order to disconnect the non-isolated regions and obtain an isolated pore region. Conversely, if it is an isolated region, then one can proceed directly to the next step: ellipse fitting. Finally, boundary extraction (red line) and ellipse fitting (yellow line) are performed on the isolated stomata pore region, and the morphological parameters (the major and minor axes, area, opening degree, etc.) of the pore are output. The integral algorithm is automatically operated in series without the prior information assumption for the stomata, and it is independent of the stomata pore dimensions, the opening degree and the data collection conditions. It is a general segmentation method. The details of each step are described below. It should be noted that the acquisition of the clear image of the full field of view is a key preparation step, and the depth composition technique is used in this paper. The image acquisition experiment parameters are detailed in the image acquisition section.

Architecture of the detection model
The living stomata were detected using the Faster regionbased convolutional neural networks (Faster-RCNN), which is a deep learning method for the detection of objects in natural images [15,16]. The network architecture is shown in Fig. 2 [14].
The Faster R-CNN network consists of the following four modules. (1) Feature extraction. The Faster R-CNN uses VGG16 to extract feature maps of the input images, and these feature maps will be used for the subsequent region proposal network (RPN) layers and fully connected layers. (2) RPN. The RPN network is used to generate the region proposals. First, several of anchor boxes are generated. Then, softmax is used to determine whether these anchors belong to the foreground or background (whether it is the object or not) after they are cropped and filtered. In addition, another branch named the bounding box regression is used to correct the anchor boxes and form more precise proposals. (3) ROI pooling. The proposals that are generated by the RPN and the feature map that are obtained in the last layer of the VGG16 are combined in order to generate a fixed size proposal feature map, and are prepared for the next step, which is the fully connected operation of target recognition and location. (4) Classification and regression. The fully connected operation of a fixed-size feature map is performed, and then the object classification is completed using softmax. Furthermore, the smooth L 1 loss is used to complete the bounding box regression operation and get the exact location of the object. The open source codes [17] of Faster R-CNN detection algorithm is implemented under the TensorFlow framework. The data set consists of 1000 stomata images with 500 images for each of two magnifications (500× and 1000×). The resolution of the images is resized to 800 × 600 in order to train the Faster R-CNN. The ratio of the number of images in the training set, verification set and test set is set to 400/400/200 prior to the model's training, and the training was performed with the following parameters: the learning rate that was used in the first 50,000 iterations was 0.001, the rate in the last 20,000 iterations was 0.0001, there were 200 epochs, and the batch size was 64. The loss function is where L cls is the Softmax Loss, L reg is Smooth L 1 Loss, = 10 , and the details of other parameters can be obtained from literature [14].

Validation of the detection model
1000 images of stomata were used in the datasets, Which include Chinese Necklace Poplar and Black poplar. We use 400 images as the training set, and another 400 images as the verification set and the remained 200 images as the test set.
The loss curve for training the model is shown in Fig. 3. We detected 1290 stomata from 1316 (Ground truth of the 200 images) stomata using this model successfully, and the recall is 98.02%. The precision is 100%. There is no false detection for the test set.
In order to evaluate the performance of the Faster-RCNN for different data, the data of dayflower provided in the literature [12] was tested using the model. The r-p curve of the model are shown in Fig. 4. The visualization results of the dayflower by the detection model are given in Fig. 5. The recall and the precision is 86.31% and 83.59%, respectively.

Stomata segmentation
The Chan-Vese model is a geometric active contour segmentation model based on the level set method, which is an improved method of the classical level set [13,18]. The model uses the average value of the greyscale inside and outside the region as the energy function in order to where C 1 = average (I(x, y)), outside(C); and C 2 = average (I(x, y)), inside (C).
The length (C) is the length of the closed contour C.
The area (C) is the inner area of C. λ 1 and λ 2 are the weight coefficients of the respective energy terms. The first two, μ and ν, are called "smooth terms", and they control that the curve maintains a certain smoothness during the evolutionary process. The latter two are called "fitting terms", and they mainly move the segmentation curve to the edge of the image in order to minimize the fitting error. The position of the final segmentation contour C and the unknown constants C 1 and C 2 are obtained by minimizing the energy function.
The steps of the CV model algorithm are as follows.
(2) Calculate the average grey value of the foreground C 1 and background C 2 . The stomata segmentation algorithm that is proposed assumes that the stomata pores are close to the centre of the ROI, and the centre of the single stoma image is used as the initial point of the evolution of the energy function in the CV model. The whole segmentation process does not need to manually specify the initial position, and does not require prior information about the stomata, which is a fully automatic method for stomata segmentation. The flow chart of the segmentation algorithm is shown in Fig. 2.
A visual representation of the level set segmentation algorithm for living stomata is shown in Fig. 3.
Our aim is to find the stomata pore regions. After segmentation by the CV model, the stomata pore region is obtained and the algorithm outputs the binary image of the stomata pore. Next, the geometric parameters of the pores are calculated.

Stomata pore measurement
The stomata pore regions that are segmented by the CV model maybe non-independent regions, which have branches or connections to other regions(see Fig. 5b). Therefore, a morphology post-processing technique is used to disconnect the non-independent regions before the pore measurement, and the algorithm must automatically determine whether the segmented region is independent or non-independent. The stomata pore measurement algorithm is composed of region shape analysis, morphological post-processing, region filling, boundary extraction, ellipse fitting, and parameter calculations. The outputs are the length of the major and minor axes, the opening degree and the pore area. The flow chart is shown in Fig. 4.
The region shape analysis is a key step for this algorithm. It can determine whether the segmented region is independent or not. If it is a non-independent region, a post-processing step is added to disconnect the connection. In this paper, morphology erosion and dilation are used to disconnect the non-isolated regions. The Solidity attribute is used to automatically determine whether the segmented region is an independent region or not. (For a non-independent region, Solidity < 0.85, which represents the concavity and convexity of the region.) A visual representation of this algorithm is shown in Fig. 5.
The Solidity is calculated as follows: where Area is the area of region, and Convex Area is the area of the smallest convex polygon. It indicates the degree of Solidity (concavity and convexity) of the region.
For most stomata images in our data, the stomata regions that are segmented by the CV model are independent and do not need the post-processing step, as shown in Fig. 6. In this situation, the Solidity of the region is greater than 0.85, and the ellipse fitting and parameter calculations can be directly accomplished for the segmented stomata pore region.

Data acquisition
In the experiment, a large depth-of-field microscope observation system, a VHX-2000 [19], from the Keyence Corporation, was used to capture the fully focused images of the living stomata of poplar leaves. The species are one-year-old black poplars (Populus simonii × P. nigra) and Chinese white poplars (Populus lasiocarpa Oliv). Using the depth composition technology of the VHX-2000 system, the images of stomata in the full field of view are taken. The depth composition parameters are set as follows. The depth-of-field interval is 2-4 μm Solidity = Area/Convex Area and the number of images included is 20-25 frames, which depend on the flatness of the leaves. The resolution of the captured images is 1600 × 1200 pixels, and we have 4.8 pixels/μm(1000×). A total of 51 microscope images were collected at two magnifications (1000× and 500×), and the total number of stomata is 708. The stoma image acquisition is shown in Fig. 7. Figure 7b is the fully focused microscope images for the living stomata using dark-field microscopy. In addition, a leaf clamp must be used to fix the plant leaves.

Pores measurement
To evaluate the segmentation performance of the method, the stomata pores were manually segmented by ImageJ under the supervision of botany experts, and the ellipse fitting and corresponding parameters of the pores were measured by ImageJ as the Ground-Truth. The results of this annotation were compared with the segmentation and measurement results of our method.
Ten of the 708 stomata failed to be identified in the test, and thus, the experimental results show that our method can be used to segment and measure the stomata pores at 1000× and 500× magnifications with high accuracy. The segmentation results of our method are shown in Fig. 8.
According to the geometric knowledge, the semi-major axis and semi-minor axis of a and b for the fitted ellipses are first obtained, and then other parameters are deduced using these two parameters. The calculation formulas are as follows: Area: A = πab . Eccentricity: The measurement errors of the major and minor axes of the stomata pores, the area, and the opening degree are obtained. The results are listed in Table 1.
In addition, this paper analyses the relationship between the major and minor axes measurement errors and the pore opening degree, as shown in Table 2. The total error is equal to the square root of the major and minor axes errors.
It can be seen from Table 2 that the segmentation accuracy of the method is positively correlated with the stomata's opening degree. The greater the stomata opening degree is, the higher the accuracy of the method for measuring the lengths of the major and minor axes. The  Fig. 6 Flowchart of the pore region segmentation for living stomata opening dimension of the stomata has some influence on the algorithm's performance. The larger the opening degree is, the higher the measurement accuracy of the method. We also give the experimental results that characterize the consistency of the algorithm. The relationship between the minor axis length's measurement accuracy of our method and that of the manual measurement method is shown in Fig. 9 (for 698 stomata). The regression line's slope is 1.015 and the R 2 is 0.9848; therefore, the algorithm has very good consistency.

Algorithm comparison
We refer to the stomata as the foreground and the nonstomatal area as the background of the images. The microscope images of the living stomata are feature-rich and include the epidermis, veins, trichomes, stains, etc.   Therefore, the stomata that are segmented by the threshold method are not often independent regions and are connected with other areas. In severe cases, even with post-processing techniques such as skeleton extraction and morphological techniques, it is difficult to isolate the stomata regions. The adaptive threshold method in the literature [12] and skeleton extraction combined with ellipse fitting [11] cannot segment the living stomata in the complex image background, especially for low contrast images (such as the living stomata data in this paper), as shown in b and c of Fig. 10. Although the threshold method is simple, it requires good-quality images with clear backgrounds. Often, a complicated post-processing step is required in order to extract the pore regions, and the prior information of the stomata is needed (such as the area range and the stomata perimeter range). Our method based on the CV model proposed in this paper can successfully segment the stomata pores in the complex background, as shown in d of Fig. 10.
The proposed method can also achieve good segmentation accuracy of the stomata with clear backgrounds (data publicly available in literature [12]). In this paper, the stomata microscope image data that are disclosed in the literature [12] are tested using our method. Due to the significant non-uniform illumination in the bright field environment, the stomata segmentation error is large if the original image is directly segmented using the level set method. Therefore, non-uniform illumination correction is necessary before the segmentation step. Since a reflection region in an image is an area with a higher pixel value, the areas with larger pixel values are replaced with some smaller values to remove the reflections. In this paper, the average grey value of the image is used to replace the pixel value of the reflective area in order to eliminate the reflection, which is followed by the CV model segmentation, and then the stomata pore region can be obtained. The segmentation results with the reflection removal algorithm are shown in Fig. 11.
In Fig. 11, c is the segmentation result by the CV model without non-uniform illumination correction. It can be seen from c that there is a notch in the segmentation results and that the segmentation error is large. Therefore, it is necessary to remove the non-uniform illumination before the pores segmentation step. d is the greyscale image of the stomata after the reflection is removed. When the CV model is run on the images in d, the stomata pore regions that are segmented are relatively intact, as shown in e. The last row f compares the To verify the algorithm's segmentation performance, our method was tested using 188 stomata of dayflower from 16 microscope images from the literature [12]. Among these stomata, 7 stomata failed to be segmented, and 181 stomata were successfully segmented by our method. The fitting relationship between the aperture value (minor axis length) and the manual measurement results of the 181 stomata apertures is shown in Fig. 12.
It can be seen from Fig. 11 that the fitting slope of our method and the manual segmentation is 0.993, and the R 2 value is 0.9828, which slightly outperforms the segmentation results of literature [12] (slope = 1.0485 and R 2 = 0.98215). More importantly, the CV-based segmentation method that is proposed in this paper does not require a priori information about the stomata, such as the ranges of the stomata pore area, the major axis length and the pore perimeter. It is an unconstrained, fully automatic stomata segmentation and measurement method. The method has good versatility and is not limited to a specific plant. It is able to segment the stomata of any plant, including non-elliptical shaped stomata.

Discussion
It can be seen from the experimental results that our method can obtain high segmentation accuracy for most of the stomata images.
Non-uniform illumination correction is not needed for the living stomata in the dark-field microscope images collected in this paper. However, it is necessary for the bright-field images. The reason may be that reflections often result from the stomata pores in the bright-field images, but there are fewer reflections that are generated from the pores in the dark-field images.  The stomata that are discarded in our method are shown in Fig. 13. These include blurred stomata, stomata with too small of opening degrees, and obscured stomata. Figure 13a shows the stomata images collected under bright and dark fields. The guard cells are also regarded as stomata pore regions, which causes the stomata to be over-segmented. The main reason for the stomata oversegmentation is that the stomata opening degree is too small, and the initial position of the CV does not fall on the pore. Therefore, segmentation results in a stomata apparatus together with the guard cells. b is the three obscured stomata in the dark lighting field. The stomata regions obtained by our method are discontinuous. Post-processing is required to obtain the complete pore region. c is the blurred stomata case. The algorithm will fail to segment these images (Figs. 14, 15, 16, 17).
In addition, the method in this paper only assumes that the pores are located near the centre of the ROI. This assumption is normally true for a detection algorithm. However, for pores with a very small opening (minor axis length < 0.83 μm and approximately 4 pixels for our living stomata images), the initial position of the evolution of the CV model (the centre of the image) usually falls outside the pores, and the segmentation result is oversegmented. Therefore, the stomata with small degrees of opening and deviations from the centre of the ROI cannot be effectively segmented by our method. Fortunately, when the stomata pore opening is too small, the CO 2 gas exchange is negligible for the stomata, and the stomata can be considered closed.
The proposed method is a region-based level set method, which finds the boundary of the segmentation by calculating the maximum difference between the average grey values inside and outside the region. It does not depend on the gradient information of the image, overcomes the edge leakage of the gradient method and is susceptible to noise. It is a perfect region-based level set segmentation method. However, since the CV-based segmentation method is an

Conclusions
In this paper, a general method for the automatic segmentation and measurement of plants' living stomata based on the CV model is proposed. The method consists of five parts: reflective removal pre-processing(if necessary), CV model-based segmentation, non-independent pore region discrimination, morphological post-processing and ellipse fitting. In this paper, the measurement accuracy of the major-axis length, minoraxis length, area, eccentricity and opening degree of the living stomata were 95.68%, 95.53%, 93.04%, 99.46%, and 94.32%, respectively, after segmenting and measuring the 698 stomata of poplar leaves. The versatility of the algorithm is better than those of existing methods of stomata segmentation. In addition, the algorithm's consistency is very good for both the bright-field and dark-field stomata images of different datasets, and the R 2 is greater than 0.98. This research will be extended to test living stomata of other plant species in the future.