Identification of pollen taxa by different microscopy techniques

Melissopalynology is an important analytical method to identify botanical origin of honey. Pollen grain recognition is commonly performed by visual inspection by a trained person. An alternative method for visual inspection is automated pollen analysis based on the image analysis technique. Image analysis transfers visual information to mathematical descriptions. In this work, the suitability of three microscopic techniques for automatic analysis of pollen grains was studied. 2D and 3D morphological characteristics, textural and colour features, and extended depth of focus characteristics were used for the pollen discrimination. In this study, 7 botanical taxa and a total of 2482 pollen grains were evaluated. The highest correct classification rate of 93.05% was achieved using the phase contrast microscopy, followed by the dark field microscopy reaching 91.02%, and finally by the light field microscopy reaching 88.88%. The most significant discriminant characteristics were morphological (2D and 3D) and colour characteristics. Our results confirm the potential of using automatic pollen analysis to discriminate pollen taxa in honey. This work provides the basis for further research where the taxa dataset will be increased, and new descriptors will be studied.


Introduction
Pollen is a natural component of honey, but due to its specific shape characteristics, it is successfully used as a method to identify the geographical and botanical origin of honey, or to determine a single-species honey. This issue is studied by the scientific branch called mellisopalynology. Honey analysis is based on visual inspection according to international standards of the International Honey Committee and this technique was described in detail in the work by Von Der Ohe et al. [1]. However, this method of evaluation is prone to be affected by subjective errors in analyses [1,2].
A number of different microscopic techniques are used for pollen identification. Morphological characteristics [3], especially the shape characteristics (length, width, circularity, and shape factor or length/width ratio) [4,5], are used most commonly. In light microscopy (LM) the morphological (pollen shape) and surface structures (exina) are used for the identification as well [6]. Other microscopic techniques include scanning electron microscopy (SEM) which allows a more accurate determination of the type of pollen grain based on differences in surface structures [7,8]. Transmission electron microscopy can provide information on the morphology of internal structures, but with regard to the complexity of processing, it is used mainly for descriptive studies in the field of palynology [9]. On the other hand, light microscopy is able to recognize far more taxa on species level than scanning electron microscopy [8]. Highlighting of surface structures is possible using other microscopic techniques of light microscopy, such as phase contrast (Ph), Nomarski phase, and polarization [8].
All microscopic techniques have some advantages and disadvantages. In combination they can provide complementary information about the observed objects. The BF is preferably for observations of transparent light-absorbing specimens showing a medium optical density. High-density or opaque specimens are not passed through by the transmitted light, thus fine details situated inside them or localized at the specimen's surface cannot be well perceived. DF is preferably for observation of low-density specimens. Phase boundaries are highly contrasted by light diffraction. Lateral resolution of DF is maximized, because the incoming circular light runs to the specimen in an oblique direction [10]. DF is also able to visualize very small structures that are beyond the resolution limit, where visible diffraction patterns are generated [11]. Ph is used for visualisation of transparent specimens in light microscopy while the peripheral zones of such specimens appear in lower brightness and the central parts are highlighted up to the brightness of the surrounding medium [12,13].
Identification of the pollen taxon is based on pollen grain morphology, and the accuracy of identification may vary according to the microscopic technique used. Ullah et al. [14] used light microscopy and SEM to study pollen micro-morphology of 18 species belonging to seven genera of Alsinoideae subfamily. And their results showed that pollen grains of this subfamily have the common feature of pantoporate ornamentation. Another study applied scanning electron microscopy to observe and describe the ultrastructure of the exine ornamentation and aperture structure of the Lamioideae taxa [15]. Also the study [16] confirmed palyno-anatomical features in monocot taxa, where the correct identification of the 24 monocot taxa used the light and scanning electron microscopy.
With the increase of computer power, it is possible to use image analysis process for pollen identification. Image analysis transfers image content into feature description, such as colour, shape, and texture [17]. The colour of the pollen grain is represented by the mean values of the RGB components (Red, Green, Blue) or the hue value which, in contrast to RGB, is independent of brightness. The shape is described by geometric parameters such as shape factor, circularity, and dispersion which, unlike the area, are independent of the size of the object. The texture is generally calculated as spatial variation of the brightness intensity of the pixels [18]. But also contour-inner segmentation can be used as textural description [19]. Some authors developed different texture description. The use of image analysis for the identification and classification of pollen grains has been validated by a number of authors in various modifications [6,[20][21][22][23]. Another study confirmed slight differences in the number of pollen grains in the case of human visual inspection or an automatic analysis [24]. However, only a few studies have been devoted to the identification of pollen grains in honey.
In these studies, bright field microscopy (BF) was used to recognize pollen grains. This microscopic technique is also used by the harmonized method [1]. The dark field microscopy (DF) [23], fluorescence technique [25], was validated for pollen analysis in palynology. In the study by Jacinto-Pimienta et al. [26], melissopalynology was performed using the phase contrast microscopic technique (Ph). No studies have yet been performed on the analysis of pollen in honey by other microscopic techniques using image analysis.
For the comparison of microscopic techniques, a quantitative melissopalynological technique was newly applied, this technique uses the filtration of pollen grains on a microcellulose filter, which is one of the harmonized analyses of honey according to the International Honey Committee [1]. The aim of this study is to compare three microscopic techniques (BF, DF and Ph) and their applicability for the classification of pollen grains by image analysis and application of new EDF (Extended depth of focus) and Volume descriptors for pollen classification.

Materials and methods
Four types of honey originating from the Czech Republic were obtained in 2019 from various honey collections. These included spring blossom honey (SBH), acacia honey (AH), lime tree honey (LTH), and honeydew honey (HH). Dominant taxa in SBH were Brassica napus (Brassica), Salix and Pyrus/prunus. Dominant taxa in AH were Brassica, Salix, Pyrus/prunus, Phacelia campanularia (Phacelia), Robinia pseudoacacia (Robinia). In LTH, dominant taxa were Brassica, Salix Pyrus/prunus, Phacelia, Tilia. HH contained Brassica, Salix, and Asteraceae pollen taxa. The identification of pollen taxa was performed using pollen atlases, scanning electron microscopy and confirmed by an independent examination (Intertek Group plc, GER).
Honey samples were processed for the quantitative pollen analysis in compliance with the harmonized method [1]. A 25 mm/3 μm Millipore 3 μm filter (Merck KGaA, GER) was used in a vacuum filtration assembly (ThermoFisher, USA). After filtration, the filter was dried at 40˚C for 12 hours and then mounted with solacryl (Merci, CZE).
The samples were analysed under the Eclipse Ci-L microscope (Nikon, JPN) with motorized stage of Proscan III (Prior, USA). Images were captured by the DFK 23U274 camera (Imaging Source, GER). Three microscopic techniques were used with 400x magnification. The used method was bright field light microscopy, dark field microscopy, and phase contrast microscopic technique. For dark field dry condenser with a Planchromat 40x lens was used for DF; C-C phase Contrast condenser and CFI Plan Achromat DL 40X lens, N.A. 0.65, W.D. 0.56 mm, Ph2 was used for Ph. Since it is not a fully automated system, a USB4000-UV-VIS-ES spectrophotometer (Ocean Optics Inc., USA) was used to set the same lighting conditions. An identical profile of light passing through the slide in the range of 420-750 nm in a place with no pollen grains or other impurities was ensured for each measurement. The pollen grains were captured in 5 different focal planes with distances of 8 μm from each other and an extended depth of focus image (EDF) was created to be used for image analysis.
The complete procedure for image acquisition comprised several steps. After solacryl mounting, microscopic slides were scanned in random positions in order to capture 1000 pollen grains for each microscopic technique. The selection of monitored botanical taxa followed from previous studies and represented taxa typically contained in Czech honey [27,28].
For taxa identity confirmation, the samples were examined by SEM. After filtration, the sample was dried in a desiccator for 24 hours and applied to an aluminium double-sided adhesive tape and followed by nanocoating with a 10 nm layer of gold using the sputter coater device Q150R ES (Quorum Technologies, UK). Scanning was performed using MIRA3 (Tescan, a.s., CZE). Pollen grains was scanned in different magnifications from 4 kx to 40 kx in high vacuum. Acceleration voltage 5.0 keV and secondary electron detector were used for picture acquisition. The procedure for processing the images and the evaluation of the data obtained by image analysis is described in Fig 1. Image analysis was performed using NIS-Elements AR 5.20 (Laboratory Imaging, CZE). Binary mask was created semi-automatically. The binarization required manual outline corrections. But a higher level of robustness, which is expected in dual segmentation, is better for our comparative hypothesis [19].
Thresholding method to the RGB colour space was applied, the thresholding range is indicated in Table 1. Object smoothing filter, and object cleaning filter, opening algorithms were used as postprocessing operation. Object cleaning filter removed small objects from the image.
Opening algorithms are based on erosion followed by dilation to smooth object. Object smoothing filter smooths the binary image contours. Similar thresholding methods were used by other authors [19,29].
The measured morphological descriptors included 2D descriptors of Area, Equatorial Diameter (EqDiameter), Perimeter Contour, Mean Chord, Length, Width, Min Feret, Max Feret 90, Circularity, Elongation, Shape Factor, Convexity, Roughness, and 3D descriptors of Volume Equatorial Sphere (VolumeEqSphere), Volume Equatorial Cylinder (VolumeEqCylinder). Morphology parameters are created by mathematical operations. They simplified the image data for the subsequent operations. The exact formulas are presented in detail (S1 Table). All parameters were calculated as standard functions. The 3D descriptors are specific for the software used. VolumeEqSphere and VolumeEqCylinder are given by Eqs 1 and 2 respectively.
Roughness ¼ Convex Hull perimeter Perimeter Eq 3 2D descriptor corresponded with surface texture was calculated by Eq 3. The used colour descriptors included Mean Intensity, Intensity Variation, Mean Red, Mean Green, Mean Blue, Hue Typical, Hue Variation, Mean Saturation, Mean Brightness, Bright Variation, Mean Density, and Density Variation. The colour features are arithmetic mean values on RGB. Hue is converted from RGB to the HSV colour-space. Descriptors were calculated as standard features, for more details see the S1 (S2 Table). And finally, the used extended depth of focus (EDF) descriptors included EDF surface, EDF Roughness, and EDF-Z distance. The automated and semiautomated microscope systems provide EDF image generation function and corresponding 3D reconstruction function by EDF technology. In our study, NIS-Elements EDF reconstruction system (Laboratory Imaging, CZE) was used. The EDF descriptors were simplified according to equations presented in detail. The EDF parameters corresponding to surface texture are EDF Surface and EDF Roughness. They were calculated by Eqs 4 and 5, respectively. EDF-Z distance corresponds to height of the reconstructed object (Eq 6).
EdfRoughness ¼ EdfSurface Area Eq 5 1000 pollen grains were selected in each sample to be classified, and the above mentioned descriptors were measured for them. Following the semi-automated thresholding, pollen grains of dominant taxa in honey samples were classified.
The data were processed statistically using the 2014.5.03 XLSTAT software (Addinsoft, USA). The normality test confirmed the normal distribution of the data. The ANOVA Tukey HSD test was used to compare the pollen characteristics. The multidimensional analysis was provided according to paradigm for object recognition where factor analysis and discrimination analysis are used for object recognition [30].
Factor analysis was used for multi-dimensional description. Factor analysis uses the correlation structure of the measured variables to create "new" variables (factors). The aim is to use the "strong" correlation of the measured variables to reduce the dimensionality of the data while maintaining as much information as possible as well as high level interpretability of the newly created factors. The basis of this method are eigenvalues and eigenvectors of the correlation matrix of the measured variables. The eigenvalues determine the amount of information retained from all the original variables in a given factor, and the corresponding eigenvectors indicate the transformation relationship of the original data to factors. In addition, the product of the square root of the eigenvalue by the corresponding eigenvector gives the correlation of the factors with the original variables. Thus, this procedure produces the same number of factors as the original variables, where each factor is a linear combination of all measured variables. The reduction of the dimension then consists of ignoring the factors with a sufficiently small eigenvalue (usually omitting the factors with an eigenvalue less than 1). Another requirement for a new factor representation is to preserve all factors containing the maximum square of the correlation with some of the original variables.
Discriminant analysis was used to compare the microscopic methods used except for SEM. The aim of the discriminant analysis is to assign the obtained observations based on their measured characteristics to predetermined groups. The idea is to find a data representation that minimizes variability within groups and maximizes variability between groups and then use it for assignment. To perform the discriminant analysis, the decomposition of the covariance matrix into a part corresponding to the intragroup variability and a part corresponding to the intergroup variability is used. Eigenvalues and eigenvectors are again used for the calculation, but this time from the matrix product of the inversion of the intragroup covariance matrix and the intergroup covariance matrix. Again, new variables are created as a combination of all the original variables. The meaning of the eigenvectors is the same as in the factor analysis, but the eigenvalues here express the "discriminatory power" of the resulting variable. The assignment itself then takes place according to the distance from the centres of the predefined groups. The measure of the success of discrimination is then the confusion matrix, where the correct groups are in the rows, the groups assigned using the analysis are in the columns, and inside there is the frequency of each category.

One-dimensional description
Identification of pollen grains by image analysis can be performed using various descriptors. Basic shape characteristics are commonly used [4,6], nevertheless, other characteristics resulting from the textural properties of pollen grains can be utilized as well [19]. For applications other than mellisopalynology, the benefits of EDF for improving discrimination of the objects studied were also confirmed [31]. Another parameter obtainable from the image data is the colour characteristic [32]. A total of 30 descriptors were measured, of which 13 were morphological descriptors, 2 were volume descriptors, 12 were colour descriptors including inner surface texture information, and 3 were EDF descriptors. Shape characteristics are used to describe the geometric parameters of an object and are determined by a simplified mathematical model. For pollen grain discrimination, they represent an important parameter, both, for human evaluation [33] as well as for automatic systems [5,19].
Basic morphological descriptors are often used for pollen identification. Pollen grains are 3D structures which are described by conventional microscopic methods as a plane structure.
The plane description factors are presented in Tables 2 and 5. The important parameter is length which was defined as mean length of the major axis [34], similar to polar axis [35]. The differences in length for evaluated pollen taxa is described in Table 2. The BF shows differences between 4 taxa (p<0.05). The DF shows differences between 5 taxa and Ph shows differences between 5 taxa (p<0.05). The reason is the halo effect around a pollen grain, which enables better automatic detection of object edges. The benefits of the DF technique for pollen identification was also confirmed by several authors [23,36]. Similar number of separated groups was in Ph where object plasticity improved edge detection. The microscopic methods have not been compared to each other in one dimensional description, because different individual pollen grains are analysed. The comparison of our results with two common public atlases ( Table 2) showed also differences between atlases and microscopic methods. Statistical analysis was not provided because only mean and SD value are available, and besides, different number of pollen grains was used for these descriptions according to the author's methodology.
Both private databases of pollen grains and publicly available pollen atlases describing characteristic shape properties accompanied by photographs of pollen [1,39] serve as a support for manual evaluation most commonly. So far, private collections [19,32,40], which are then used to discriminate pollen taxa, are used for image analysis. In comparison with other authors, for Ph, the Brassica pollen length was 24.48 μm, Prunus pollen length was 36.94 μm, Asteraceae represents a wide group the pollens sized within the range of 33. .08 [41]. For BF, the size of Brassica was measured at 21-30 μm, for Salix 16-27 μm, for Asteraceae (single species Helianthus) 27-31 μm [19], for the Prunus 55.15 μm [5]. On the other hand, size and morphological similarities were likely related to the confusion in taxa classification for DF [24]. Size similarity has also been confirmed in our study ( Table 2, Fig 2). The morphological characteristics of pollen can differ between botanical species (Fig 2) as well as within botanical species. Morphological differences within a species have been described by several authors and include  [42]. The length is the common descriptor used for pollen classification based on image processing and manual classification. The length of the pollen can be affected by sample preparation which causes differences between values published https://doi.org/10.1371/journal.pone.0256808.g002 [37,38], but also between the results of this study. Another factor that affects pollen grain size is moisture level of the cytoplasm [43]. For the size changing, and important factor is also the type of the apertures which also provide routes for transfer of water and other substances [35]. The aperture can also be used for precise identification of the pollen taxa [44]. This feature was not evaluated in this study, primarily due to the need of a specific algorithm for identification. But apertures have impact on pixel colour, and in our study they are included in colour descriptors. Pollen grain variability ranges considerably and a clear classification according to pollen grain size is not possible in botanical taxa of similar size. Differences in the length measurement of the pollen grain (Table 2). In addition to a similar length of pollen grain in Tilia and Robinia, both of these botanical taxa are morphologically similar (Fig 2). Therefore, we also assume their similar variability in the spatial arrangement. Likewise, the similarity can be explained for Brassica and Tilia, where some projections have the same length in planar view. The quantitative technique of filtered grains was used in our study [1], with regard to light scattering, information on the real colour of pollen grains is preserved and the halo effect around pollen grain remains (Fig 2A-2D, DF). As documented in Fig 2, the image information about the pollen grains on the filters differs from conventional microscopic techniques. The structural as well as colour characteristics of the pollen grains stay, but a lower contrast between the pollen grain and the background is achieved primarily in DF compared to techniques using pollen grain centrifugation. In addition to easier quantification, the speed of this analysis is also an advantage of this filtration technique. However, with regard to the need to make the filter transparent using immersion oil or solacryl, it is necessary to take the morphological changes caused by dehydration, both during drying before mounting and the effect of the organic solvent in the mounting medium. into account. With regard to the stabilizing structures of the pollen wall, these changes are not significant for all taxa.
Other morphological descriptors, some of which have lower impact on the pollen size variability, can be easily calculated also by the image processing. Using the simple geometrical measurement, the shape and ornamentation of pollen can be analysed [32]. The length descriptor as an important parameter was discussed separately ( Table 2). The other descriptors used in our study include Area, EqDiameter, Perimeter Contour, Mean Chord, Width, Min Feret, Max Feret 90, Circularity, Elongation, Shape Factor, Convexity, Roughness.
Our results show that circularity, elongation, shape factor, and convexity are not proper for pollen discrimination (Table 3). They divided the pollen taxa into three different groups for Ph, DF, and BF, except for Circularity and Elongation where 2 separate groups were formed for BF, and Elongation with 2 separate groups for DF and Shape Factor with 4 separate groups (Table 3). However, the shape descriptors are not reflected properly in the planar view, namely with regard to the Z-view, which can eliminate hight variability for some botanical taxa [19].
The shape descriptors of Circularity, Elongation, Shape Factor, and Convexity describe specific characteristics regarding the geometry of a pollen grain. These shape descriptors are invariant with respect to translation, rotation, and scaling. These measurements are clearly defined, relatively easy to use and there is a fast and easily understood procedure for their computing [45]. Therefore, these parameters specify the discrimination of individual pollen grains and are also used in descriptive palynological studies. Common pollen shapes are circular, elliptic, lobate, triangular, and polygonal [46], which are described in our study by the parameter of Circularity. The Elongation is more representative for common description of the pollen terminology know as P/E (polar axis/equatorial axis). P/E is widely used system for a shape classes definition. For the spheroidal pollen in range of 0.88-1.14, for suboblate 0.75-0.88 and for subprolate 1.14-1.33 [35]. The Shape Factor described object roughness. It is based on morphometric measurement, this characteristic represented pollen ornamentation. In the morphometric principle it is able to detect spine only in distances according to the resolution of the microsystems used. Other morphological descriptors show higher ability to divide taxa into separate groups (Table 4). These descriptors are size dependent and they provided basic information of the pollen shape [32].
For BF and Ph, pollen was classified into five groups for Area, EqDiameter, and Perimeter Contour. In contrast, for DF all taxa were classified into 6 different groups (p <0.05). The Mean Chord shows the variable differences of taxa groups in microscopic methods comparison. DF divided taxa to 6 groups, BF to 5 groups, and Ph to 4 groups. The Equatorial Diameter is a parameter used also for pollen characterisation by human evaluation. In comparison, EqDiamether is slightly different because it is calculated as a diameter of a circle with the same area as the measured object (S1 Table), not as a line, lying in the equatorial plane [35]. The results are in conformity with other authors who confirmed DF as a suitable method for pollen discrimination [23,47]. The chord length describes object shape which is represented by a set of n-vertex polygons [48]. In our study, it is characterised as a mean value measured at 0, 45, 90 and 135 degrees directions (S1 Table). The chord length is highly relevant to the properties related to heterogeneous materials [49,50], in conformity with these, our results confirmed that it is also usable for pollen classification in BF and DF.
The complementary morphological descriptors were able to classify the pollen taxa into separate groups. Six different groups were formed by Width, Min Feret, and Max Feret 90 descriptors for BF and DF (p <0.05). Ph reached a lower number of separate groups for Width and Min Feret, nevertheless, MaxFeret90 separated six different groups similarly to BF and DF (Table 5). Our results confirmed that these descriptors can be used for classification and they show slight differences between individual microscopic methods. Feret diameter is also confirmed in literature as a relevant descriptor for pollen classification by DF [51] and by BF [5]. P/E ratio is commonly used for pollen description, where P stands for polar diameter and E stands for equatorial diameter. In our morphometric characterization, it is represented by Elongation (Max Feret 90/Min Feret) but spatial orientation of the microspore is not reflected. The Width is a complementary descriptor to the Length. It is appropriated for elongated pollen taxa and for dried pollen. The useful length descriptors for particle size classification are also Min Feret and Max Feret 90 and they can be also use for description of object elongation [52]. The Max Feret 90 is a complementary parameter with equatorial diameter used for pollen characterisation in pollen terminology [35].
The 3D descriptor was able to divide the pollen taxa to separate groups depending on the microscopic technique used (Table 6). For DF, the taxa were divided into 6 groups, for BF and Ph, they were classified into 5 significant different groups (p <0.05). 3D discrimination for pollen taxa was confirmed for BF [40,53,54] where various characteristic were used. Our results confirmed that 3D reconstruction can be used for discrimination of pollen taxa also for DF and Ph, where differences between taxa were found (p <0.05).
As already mentioned before, pollen grains are 3D structures. Various methods for 3D acquisition and visualization of pollen was used [40,53,54]. 3D descriptors in our study are calculated from several focal planes. The descriptor used was VolumeEqSphere, and VolumeEq-Cylinder and they represented a simplified 3D model of irregular pollen shape (Table 6). This model is used for size characterization of a particle in 3D reconstructed objects [55].
Important factor for pollen determination by image processing is the surface of the pollen grain [6,21] formed by the pollen wall. In our study, this descriptor is based on colour and shape characteristics. A comparison of the differences within the observed taxa is summarised  (Table 7). Colour characteristic represented textural information calculated as the spatial organisation of RGB levels of an image (Mean Brightnes, Intensity Variation, and Bright Variation). The shape characteristics represented textural information on Roughness, EDF Surface and EDF Roughness of pollen grains which correspond to the exine layer. In particular, the outer exine layer-sexine-creates various arrangements that form the typical pollen surface (columellae, tectum, supratectal elements) [56]. In 2D, the exine sculpture is most characterized by the morphometric parameter of Roughness and in 3D by EDF Surface and EDF Roughness.
Exact correspondence with the exine sculpture, however, was not demonstrated for either parameter. This is because the accurate exine visualization can only be achieved at higher resolutions, typically using SEM microscopic techniques [57]. In light microscopy, the differences cannot be detected with high accuracy. Sometimes, there is even a different classification of

PLOS ONE
Pollen identification by microscopy techniques exine sculpture depending on the microscopic technique used. Typical exine sculptures are shown in Fig 3. The suitability of surface texture for pollen identification was also confirmed [6]. EDF derivate descriptor without texture correspondence is EDF-Z, these represented mean distance in reconstructed Z-view. This distance is complementary to the length of the pollen grain but in the direction towards the optical view. The results of pixel brightness, hue, and pixel distribution are summarised in Table 8. The Mean Brightness, Intensity Variation, and Bright Variation divided the taxa into 4, 3, and 5 groups (p <0.05) for BF, DF, and Ph, respectively. Hue Typical represented the colour of the pollen and it divided the taxa into 4, 5, 4 separate groups (p <0.05) for BF, DF, and Ph, respectively. Surface texture of a pollen grain affects the distribution of pixel brightness values of the observed object. This distribution can be characterized by various parameters (Mean Value, Coefficient of Variation, Askewness and Kurtosis). The Mean Value and Coefficient of Variation were used in our study. The coefficient of Variation including Intensity and Bright Variation. These descriptors represented characteristic elements on pollen surface [58]. Notably, Lacuna, aperture, annulus, arcus, atrium, sculpture elements, tectum, and columella are used for precise identification of pollen taxa by human evaluators. In our study, variation value characterises all surface elements based on differences between pixel intensity and brightness.
The use of colour characteristics to identify and discriminate microscopic objects has been validated in several studies [32,40,53,59,60]. Although colour characteristics are important for pollen grain discrimination, their transfer to other systems is difficult. The reason is the influence of the light source, microscopic glass, microscope optics, RGB sensor, and the influence of post-processor processing on the colour of the digital image [61]. For proper transferability, the optical system of the microscope as well as camera parameters must allow accurate standardization of the device [62]. In our study, only one optical system and one camera were used, and for each measurement, the microscopic system was set to the reference spectrum using a spectrophotometer as a standard.
In our study, RGB colour descriptors confirm significant differences between taxa (p <0.05) for microscopic methods (Table 9). Mean value of RGB and Mean Density were used as colour descriptors. The highest variability between taxa was confirmed for Ph. Different groups for Mean Red and Mean Green reached 3, for Mean Blue reached 5 (p <0.05). The BF shows 3 different groups for Mean Red and Mean Blue, 4 for Mean Green (p <0.05). The DF shows 3 different groups for Mean Blue and 2 groups for Mean Red and Mean Blue (p <0.05). The derivate function for RGB is Mean Density. Mean Density shows the highest variability for Ph (5 groups) and similar for BF and DF (3 groups). In contrast to [63], where BF was used and only Mean Red was relevant, our results show that Mean Blue and Mean Green colour descriptors can be relevant for pollen classification depending on the microscopic method used. This reflects human perception, where the RGB colour system is used. By red-greenblue combination, every possible colour is defined [59].
Transferred RGB colour space in to the HSB (hue, saturation, brightness) colour space is shown in Table 10. Our results show that Hue Variation had a lower variability between taxa in BF. For DF and Ph, the taxa were separated into 4 different groups (p <0.05) ( Table 10). Mean Saturation was more variable between the taxa. For BF, the taxa were separated to 5 different groups, for DF and Ph into 4 different groups (p <0.05). The brightness represented also textural features as described above (Table 8). The mean intensity classified taxa for Ph, into 5 groups, for BF into 4 groups, and for DF into 3 different groups (p<0.05). Mean Intensity was used for pollen classification using autofluorescence image analysis [64]. HSB colour space was successfully used for pollen discrimination in BF [63]. Our results show that other microscopic methods can be used as well (Table 10).

PLOS ONE
The sequence of the measured descriptors for morphological and colour characteristics of pollen has been described above. However, none of the descriptors allows clear identification of the studied taxa in any microscopic methods used. In the case of multiple descriptors, use of discriminant analyses is the most appropriate, which has also been successfully validated for pollen discrimination by a number of authors [4,19,21,24,60].
It can also be assumed that the more taxa studied, the more characteristics will be necessary for the recognition. In our study, honey from the Czech Republic was used, which is typically affected by large-scale farming [28]. In comparison, honey can be collected also from country with fragmented agricultural activity [65], from national parks [66,67] or other region with high floral diversity, where typically one family is represented by multiple species, wild cultivars and endemism. These types of honey are also typical by high content of pollen with similar shape, polar and equatorial length, P/E, and pollen ornamentation. In these types of honey, the precise identification should be provided. In general, identification of the type of surface, lacune, apertures, thickness of apertures, ornamentation is provided often by means of a combination of different microscopic techniques.

Multi-dimensional description
Due to the number of characteristics measured and the dependence of some of them, it is possible to reduce the size of the input data before. This was done by factor analysis. Already the first two BF factors retain 56.26% of the information contained in the original thirty characteristics (68.35% for DF and 60.67% for Ph). In all cases (BD, DF, and Ph), all factors corresponding to eigenvalues greater than 1 and all factors containing the maximum of the squared correlation with one of the original variables were used for the discriminant analysis (S3 Table. Factor variable correlation matrix). Due to the size of the eigenvalues relevant to other factors, the information from them can be declared insignificant and can be considered "as noise". The first 6 factors retaining 91.1% of the original information were used for discriminant analysis of the data obtained by DF; and 7 factors retaining 89.23% and 89.09% information for DF and BF, respectively. The results of the factor analysis are illustrated by the scree plots in Fig 4. Factor 1 for BF consists mainly of morphometric parameters (S3 Table). The measured morphometric parameters can be divided into 2D characteristics (Area, EqDiameter,  Table). In addition to colour characteristics identical to BF, (Hue Typical, Hue Variation, Mean Saturation, Bright Variation, Density Variation) contribute to Factor 2 as well. In contrast, Mean Red has a lower impact on this factor in DF. Recalculated morphological parameters (Convexity, Circularity, Roughness) also have a significant share in Factor 2. For Ph, the characteristics of Factor 1 (S3 Table) are identical to BF except for the EDF characteristic (EDF Surface) which has less weight. In addition to the colour characteristics identical to BF, the calculated colour characteristics (Intensity Variation, Hue Variation, Mean Saturation, Bright Variation) also contribute to Factor 2.
The significance of most of the measured morphological characteristics in the first factor confirms the suitability of determining pollen grains using size and shape characteristics. These characteristics are used both, in the manual evaluation of pollen grains [68,69] as well as in the evaluation of pollen grains by image processing.
Although morphological characteristics significantly contribute to Factor 1, they are not sufficient to discriminate the assessed taxa.
The colour characteristics of pollen grains play a significant role in F2 and also contribute to F4 and F5 in the case of BF. For DF and Ph, colour characteristics contribute significantly to F2 as well as to other factors to a lesser extent (S3 Table).
However, the first two factors are not sufficient for the accurate classification of pollen grains. In addition to the already mentioned colour characteristics, shape characteristics also contribute to other factors in BF.

Multi-dimensional discriminant analysis
The microscopic methods were contrasted based on the correct classification rate (CCR) of the discriminant analysis (DA).
The overall CCR of pollen grain discrimination based on the created model differed between the microscopic techniques used; the Ph reached 93.05%, DF 91.02%, and BF 88.88% (Figs 4-6).
For the Ph microscopic technique, the lowest CCR was achieved for the Salix taxon (60%) (Fig 5). We assume that the reason is a similar size and similar reticulate exine ( Table 2, Fig 3) to Brassica and thus the classification into this taxon occurs incorrectly. Likewise, for the same reasons obviously, Tilia is incorrectly classified into the Brassica taxon with CCR 74.07%. Incorrect classification also occurs in Pyrus/Prunus 78.46% and Robinia 90%, where we assume that the reason is mainly the similar tricolpate shape of the pollen grain, see Fig 2. The use of Ph imaging technique thus proves to be the most suitable method for discriminating the studied pollen taxa (Fig 5). Ph is used in melissopalynological analysis as one of the microscopic methods [8,41] that allows 3D imaging of unstained structures. It can also be successfully used, for example, to identify honeydew elements [70]. In our study, its applicability for the semiautomated image processing method was also confirmed.
For DF microscopy, the least accurate discrimination was for the Salix taxon (56.91%), similar to Ph (Fig 6). Most misclassified pollens were classified as the Brassica taxon. Lower CCR was also for Tilia (74.70%) when most incorrectly classified pollen grains were assigned to the Brassica taxon. The pollen surface of these taxa forms a reticulate exine and discrimination using textural parameters is therefore difficult. CCR for Pyrus/prunus was 78.24% due to incorrect classification of part of the pollen into the Brassica and Robinia taxa. We attribute this

PLOS ONE
Pollen identification by microscopy techniques incorrect discrimination to the tricolpate shape similar to Robinia and to the size similar to Pyrus/prunus and Brasica sp. Another factor of incorrect identification is also attributed to the shape similar to Brassica in some projections that are given by the polar or equatorial view. That is typical for heteropolar pollens.
The DF technique of automatic as well as semi-automatic analysis of pollen grains has been validated by several authors for pollen discrimination. In the study by Sevillano et al. [36], an automatic system for acquisition and identifying pollen grains was researched and the CCR of the DF method reached 97.86%. Lagerstrom et al. [24] reached 81.2% CCR when using DF. The difference in CCR values are caused not only by the different classification systems used but also by the different number of identified images and the number of classified pollen taxa.
For the BF technique, the lowest CCR was reached for Tilia (44.83%) and pollen grains that were most often incorrectly classified as Brassica. Even in this case, we attribute this to the heteropolarity of the Tilia pollen grain. Similar to Ph and DF, Salix was partially incorrectly classified into Brassica (Fig 7). With this technique, the CCR was below 50% (47.24%). BF is the most commonly used technique for evaluating pollen by both, visual inspection [1,39] and automatic analysis [18,24,[71][72][73]. Differences in CCR vary between authors and range from 64% to 100% [6,63].
The lowest CCR for all microscopic methods was for Salix and Tillia taxa. Given the diversity of pollen taxa in honey, it is likely that even if automated analyses of pollen grains are established, precise identification of some taxa will be necessary. Either imaging-based methods such as SEM [21,74] or transmission electron microscopy [75] may be considered most appropriate. However, not all taxa can be distinguished from each other by melissopalynology methods [76] or by chemical markers [77]. However, specific identification of pollen taxa is also possible by methods based on other principles. Methods based on DNA detection can be considered as promising. When the DNA sequence is compared instead of the pollen grain morphology. For example, the chloroplast trnL barcoding fragment [78,79], rbcL DNA [80,81]. However, DNA methods also have their disadvantages which lead to a reduction in https://doi.org/10.1371/journal.pone.0256808.g007 discriminatory power [82,83]. DNA methods allow for the examination of a larger number of pollen grains [80], on the other hand, the establishment of automated analysis will allow for a comparable number of grains to be examined and therefore a more accurate result.
Our study confirmed that discrimination of 7 botanical taxa in all microscopic methods using a membrane filter is possible. In contrast to paleoclimatology, for which this application is problematic [20], this quantitative method can be also used for automatic discrimination of pollen grains in melissopalynology. The discriminant model was created from 2243, 2419, and 2482 images of pollen grains for Ph, DF, and BF, respectively. The difference in the number of pollen grains is due to the different ability of the system to automatically focus and recognize the pollen grain. Conventional microscopic image acquisition systems are optimized for BF or fluorescence. Thus, most pollen grains were detected in BF; on the contrary, detection of some pollen grains was problematic for the Ph system and the number of detected and subsequently analysed pollen grains was therefore the lowest. In the case of routine capturing, the detection and focusing mechanism can be adjusted for Ph conditions. The BF microscopic method achieved the lowest CCR, however, it is the most commonly used method in the literature, for which there is a number of accompanying materials that facilitate the identification of pollen grains, especially for the human evaluator. Improvement of CCR in all compared techniques can be expected with the introduction of other forms of classification techniques, such as artificial neural networks [4] or convolutional neural networks [36].

Conclusions
Various pollen descriptors were introduced in this study. The differences between pollen taxa were confirmed for microscopic methods used. Morphometric, volumetric, and colorimetric descriptors were mostly relevant for the discrimination with diverse variability for microscopic methods. 3D descriptors classified pollen taxa mostly in dark field, colour descriptors mostly in phase contrast. Likewise, EDF descriptors were helpful for the classification mostly in the dark field microscopy. The EDF Surface classified the most pollens into separate groups, compared to EDF Roughness and EDF-Z. None of one-dimensional descriptors can discriminate all pollen taxa into separate groups. For this reason, the multidimensional analysis was used. Discrimination analyses show different CCR. All microscopic methods can be used to classify pollen in honey into the appropriate taxa. The best technique for classification was phase contrast microscopy, where the CCR values achieved 93.05%. For dark field microscopy CCR value was 91.02%, and for bright field microscopy CCR reached 88.88%. However, although the method achieved a high CCR for the selected taxa, a large number of new taxa can be expected in the honey that will require human evaluation of pollen morphology for their correct identification. Our results provide background for future research where the taxa dataset will be increased, and new representative descriptors will be studied for pollen classification in honey.
Supporting information S1