Automatic Image Classication Using Neural Networks Increases Accuracy for Allergenic Pollen Monitoring

Monitoring of airborne pollen concentrations provides an important source of information for the globally increasing number of hay fever patients. Airborne pollen are traditionally counted under the microscope, but with the latest developments in image recognition methods, automating this process has become feasible. A challenge that persists, however, is that many pollen grains cannot be distinguished beyond the genus or family level using a microscope. Here, we assess the use of a Convolutional Neural Network (CNN) to increase taxonomic accuracy for airborne pollen. As a case study we use the nettle family (Urticaceae), which contains two main genera (Urtica and Parietaria) common in European landscapes which pollen cannot be separated by trained specialists. While pollen from Urtica species have very low allergenic relevance, those from several species of Parietaria are severely allergenic. We collect pollen from both fresh as well as from herbarium specimens and use these to train the CNN model VGG16. The model shows that Urticaceae pollen can be distinguished with 98.3% accuracy. We then apply our model on Urticaceae pollen collected from aerobiological samples and show that the genera can be condently distinguished, despite the more challenging input images that are often overlain by debris. Our method can also be applied to other pollen families in the future and will thus help to make allergenic pollen monitoring more specic.


Introduction
Pollen allergies are on the rise globally, with worldwide approximately 10 to 30% of adults and 40% of children affected 1,2 . For patients the symptoms include a runny nose, sneezing and itchy eyes, mouth or skin. Control measures and medication are readily available, but to alleviate the symptoms most e ciently, exposure to allergens should be kept to a minimum 3 . Therefore, for more and more people, fast and accurate monitoring of airborne pollen provides an essential early warning system 4,5 . Pollen concentrations in the air are monitored using samplers that collect airborne pollen on sticky tape, e.g.
Hirst type samplers 6 . These tapes are microscopically inspected for their pollen content, a process that requires highly trained specialists. Moreover, although the allergenic pollen from some plants can be monitored at the species level (e.g. species of plantain, Plantago L. 7 ), many other pollen grains cannot be accurately identi ed to this level. In many taxa, only a genus-or family-level identi cation is possible using current microscopic methods 8 . This is problematic since different species and even genera within the same family can possess very different allergenic pro les. An extra challenging factor in airborne pollen identi cation from Hirst samples is that they are collected directly from the air. In contrast to acetolyzed pollen 9 , these still contain all organic material, and de ning features are less apparent 10 .
This identi cation challenge is exempli ed in the case of the nettle family (Urticaceae). Pollen grains produced by all species from the genus Urtica L. (stinging nettles) have a low allergenic pro le 11 , while pollen from several species of Parietaria L. (pellitory) are a major cause of hay fever and asthma, in particular P. judaica L. and P. o cinalis L. 12,13 . These pellitory species are native to the Mediterranean, but throughout the second half of the twentieth century, a range expansion occurred through north-eastern Europe, the Americas and Australia as a result of anthropogenic distribution and climate change 14,15 . Parietaria sensitization is highly different per geographic area, but has been reported to reach 80% in southern Italy while a value of 13% was found in the United Kingdom 16 . Species of Parietaria ower throughout the year but their main owering peaks occur in May-June and August-October, which overlaps with the owering season of Urtica species (June -October) 17 . Cross-reactivity is present between species of Parietaria, but is absent between the genera Urtica and Parietaria 11,18,19 . Parietaria pollen is microscopically indistinguishable from that of Urtica and their contribution to the total airborne pollen load is currently not assessed in either native or expanded range 20 .
Pollen grains from Urtica and Parietaria species have a simple morphology: they are small (~11-20µm), rounded to slightly ellipsoidal tri-, tetra-or zonoporate with a psilate to scabrate surface ornament and small pores. Most species have an annulus around the pore, i.e. a thickening of the otherwise very thin exine and a germination area called the oncus (lens-shaped body located in the apertural region) 7 . The only species of Urticaceae that can be distinguished in aerobiological samples is Urtica membranacea due to its small size (~10-12μm) and high number of pores (usually more than six 21 . The main difference between the pollen of Urtica and Parietaria are the slightly smaller size and coarser surface ornamentation of Parietaria, and a more angular outline and more pronounced annulus of Urtica 22 . Neural networks have been used successfully to manage both the tasks of differentiating pollen from non-pollen debris as well as correctly identifying different taxa (for an overview please refer to 23 ).
Automatic image recognition can however also be used to improve identi cation of pollen taxa that are di cult to distinguish using traditional methods. Subtle variations in morphology that are not readily apparent through microscopic investigation may be consistently detected by neural networks. This has for example been shown for the highly similar pollen of black spruce (Picea mariana (Mill.) Britton, Sterns & Poggenb.) and white spruce (Picea glauca (Moench) Voss) using machine learning 24 and for pollen of 10 species of the thistle genus Onopordum L. using an arti cial neural network 25 . However, neural networks have so far not been tested for improvement of taxonomic resolution in unacetolyzed pollen in aerobiological samples.
Here we use the Convolutional Neural Network (CNN) model VGG16 to distinguish morphologically similar, unacetolyzed pollen from the nettle family. We collect pollen from all species of Urticaceae present in the Netherlands (Urtica dioica, U. membranacea, U. urens, Parietaria judaica and P. o cinalis), freshly collected as well as from herbaria, and create a pollen image reference dataset. We rst optimize the performance of the CNN model using data augmentation and 10-fold cross validation. The model is then applied to Urticaceae pollen from three aerobiological samples with high Urticaceae pollen count. We use one sample from the Leiden University Medical Centre (LUMC), Leiden, the Netherlands as well as one sample each from Lleida and Vielha, Catalonia, Spain (ICTA-UAB). In the Netherlands, stinging nettles (Urtica) are highly abundant and therefore it is expected that most Urticaceae pollen will be from this genus. Urtica is also expected to be dominant in Vielha, while in the direct surroundings of Lleida, Parietaria is very abundant.
The main objectives of this study are (1) to see whether a CNN model can distinguish morphologically similar unacetolyzed pollen of two common genera in the Urticaceae family that have highly differing allergenic pro les; (2) to test whether the trained model can be successfully applied on aerobiological samples containing more complex and for the model before unseen input images.

Model Performance
The trained VGG16 model accurately identi ed pollen to the genus level for 97.5% of the test images for Urtica and 98.5% for Parietaria ( Figure 1A). For Parietaria four images were misclassi ed (all to Urtica), while seven were misclassi ed for Urtica (all to Parietaria). The species Urtica membranacea was con dently distinguished from all other Urticaceae species (100%), but distinction at the species-level was not possible for any of the other Urtica and Parietaria species. This is because the distinguishing features of pollen from these species (e.g. exine ornamentation) could not be resolved in the used image projections. The average accuracy for the model using the two genera (Urtica and Parietaria), and the species Urtica membranacea is 98.3%. We nd that using a hard voting approach in our model improves the accuracy (Correct Classi cation Rate, CCR) by 0.58% compared with calculating the average CCR of 10 models ( Figure 1B; see Materials and Methods for full explanation). Including data augmentation into the model further improved the accuracy of the model by 0.95% on the test data set.
For all species, pollen were collected from a minimum of four different plants. Looking at the raw pollen images from the different plants, we identi ed intra-speci c differences that result from natural variability within each species. To test whether VGG16 learned the pollen-speci c distinguishing features rather than sample-speci c details, we produced feature maps ( Figure 2). Despite the highly variable input images of unacetolyzed pollen from different plants, the model consistently learned features such as edges in the rst convolutional layers, while ner features such as pores and annuli were learned in deeper layers. Table 1 shows the results of the CNN on Urticaceae pollen from an aerobiological sample from Leiden, the Netherlands, as well as from Lleida and Vielha, Catalonia, Spain. We set the identi cation threshold at a value of 60% as derived from the model test images, and therefore the CNN also returned unknown images (see Supplementary Table S2 for the full results). For the sample from Leiden, 85.7% of the Urticaceae pollen were identi ed as Urtica, with only minor presence of Parietaria (4.5%). The sample from Lleida shows dominance of Parietaria pollen grains (81.0%) while 14,3% of the Urticaceae pollen grains were classi ed as Urtica. Finally, for Vielha we nd a mixture of ~70% Urtica and ~20% Parietaria. No Urtica membranacea pollen grains were identi ed in any of the samples. On average, unknown images account for 8.7% of the total images when using 60% identity threshold. When using a stricter identity threshold (e.g. 70%, see Table 1), the unknown image category increases to an average value of 13.5%.

Discussion
This study demonstrates incorporating neural networks to increase the taxonomic resolution of pollen grains in aerobiological samples. The feature maps in Figure 2 show that the trained deep learning model VGG16 looks at the traditionally used morphological features to distinguish Urtica from Parietaria pollen. The characteristic thickening of the exine around the pores of Urtica shows the highest activation in the deeper convolutional layers. The distinct thickening is missing in Parietaria pollen, and the model instead focuses on the pollen outline. As expected, the only species to be distinguished by our model is Urtica membranacea which shows a slightly angular outline due to the larger numbers of pores ( Figure 2D). For the other species used in this study, no distinction was possible even though it has been shown that pollen from species of Urtica (U. dioica and U. urens) and Parietaria (P.judaica and P. o cinalis) can be separated based on differences in their exine ornamentation 22 . These differences can, however, only be imaged using specialized microscopy methods such as SEM or phase-contrast imaging, and are very hard to visualize using bright eld microscopy. Furthermore, these features are obscured when pollen grains are not acetolyzed. For our purposes, this species level distinction is not relevant as no known differences in allergenicity are known between either the species of Urtica 11 or Parietaria 18 . This is the rst time a deep learning model has been used to increase the taxonomic accuracy of unacetolyzed pollen. The model represents a signi cant improvement of earlier attempts in distinguishing Urticaceae pollen using automatic image classi cation. In a previous study using shape and texture features, pollen from three Urticaceae species could be distinguished from another with a 89% accuracy 26 , though only a small image dataset was used to train the model (i.e. 100 images per species). Similar results were obtained by 27 where shape features were used with a minimum distance classi er to obtain a 86% accuracy between three species of Urticaceae. Because not all species of Urticaceae were included and a low amount of training images was used, these studies have limited applicability to the highly diverse pollen encountered in aerobiological slides.
Deep learning models have shown similar accuracy rates to ours on larger and more varied pollen datasets as well, but these either focussed on the family level 28 or on insect-collected pollen for honey analysis 29,30 . Increasing the taxonomic resolution of pollen grains has been achieved by incorporating an extensively trained deep learning model with super-resolution microscopy on a case study of fossil pollen 31 . Similarly, incorporating SEM images has been found to allow for highly accurate distinction of pollen types 32 . These microscopy methods, however, are often much more expensive than using Light Microscopy and require extensive sample preparation. Moreover, nearly all of these studies work with acetolyzed pollen that allow easier recognition of distinguishing features, and used pollen collected from a single location. To validate our model, we tested it on Urticaceae pollen from aerobiological samples collected from different locations in Spain and the Netherlands. Most of the pollen grains from the sample from Leiden, the Netherlands were identi ed by the deep learning model as Urtica, with only a low number of images identi ed as Parietaria. While Parietaria plants are relatively abundant around the sampling location in Leiden and were owering on the chosen date, its pollen is most likely simply outnumbered by the much larger number of nettles in the area. For Lleida (Catalonia), where pellitory plants are abundantly present, Parietaria pollen grains dominated the assemblage, while the sample from Vielha showed a mixed assemblage. The number of unknown images was the highest for the sample from Vielha (11.5%), which is most likely the result of the presence of more debris on the pollen grains making a certain identi cation impossible. In all aerobiological slides, debris on top of or below the pollen grains was observed in different focal plains. Nevertheless, the model still successfully classi ed most of the pollen grains, and in most cases with high con dence (Supplementary Table S2). This shows the potential broad application of this method and opens up opportunities to study both seasonal as well as long-term yearly dynamics of Parietaria versus Urtica abundance of airborne pollen, as well as using this method to distinguish other morphologically similar species of allergenic importance from different families (e.g.

Betulaceae, Amaranthaceae, Oleaceae).
A limitation of our method is that currently pollen from aerobiological slides have to be located manually. In other systems like the commercially available Classifynder system, pollen are automatically located and imaged using dark eld imaging after which a simple neural network classi es the pollen 33 . This is also the case for the BAA500 system used by e.g. Oteros et al. 34 , that was particularly developed for recognizing and classifying unacetolyzed airborne pollen for hay fever predictions. While both systems achieve automated and accelerated pollen counting, our method instead particularly increases the accuracy of information useful for allergy prevention by making it more speci c.

Conclusions
In conclusion, using a combination of an image-processing work ow and a su ciently trained deep learning model, we were able to differentiate unacetolyzed pollen grains from two genera and one species in the nettle family. These are genera that are indistinguishable with current microscopic methods but possess different allergenic pro les, and thus the ability to differentiate them is of medical signi cance. Our method can be more broadly applied to distinguish pollen from similarly challenging allergenic plant families and can thus help in producing more accurate pollen spectra to improve the forecasts for allergy sufferers.

Collection of Pollen
Pollen were collected from all ve species of Urticaceae found in the Netherlands. In the genus Urtica, the native species U. dioica L. (common nettle) and U. urens L. (small nettle) are ubiquitous in nitrogen rich moist areas, ditches, woodlands, disturbed sites and roadsides. The exotic Mediterranean species U. membranacea is rarely encountered, though is included in this study since its range is expected to increase due to the effects of global warming. The genus Parietaria is represented in the Netherlands by the species P. judaica L. (pellitory of the wall) and P. o cinalis L. (upright pellitory) that both occupy rocky substrates, mainly in the urban environment 15 . Moreover, P. judaica has shown a big increase in abundance over the past decades, e.g. in the Netherlands (Supplementary Figure S1), but also in many other parts of the world.
Pollen from all Urticaceae species was either freshly obtained or collected from herbarium specimens (Naturalis Biodiversity Center; see Supplementary Table S1) To produce palynological reference slides, thecae of open owers were carefully opened on a microscopic slide using tweezers. A stereo microscope was mounted in a fume hood to avoid inhalation of the severely allergenic pollen of Parietaria species. Non-pollen material was manually removed to obtain a clean slide. The pollen were mounted using a glycerin:water:gelatin (7:6:1) solution with 2% phenol and stained with Safranin (0.002% w/v). These represent the same conditions also used in airborne pollen collected by a Hirst type sampler. Cover slips were secured with para n.

Pollen Image Capture
A total of 6,472 individual pollen grains were scanned, using a minimum of 1000 individual pollen grains for each species, though numbers varied (Supplementary Table S1). The system used for the imaging was a Zeiss Observer Z1 (inverted microscope) linked to a Hamamatsu EM-CCD Digital Camera (C9100), located at the Institute of Biology Leiden (IBL). Grayscale images were used, since the pollen were stained to increase contrast and not for species recognition.
The imaging procedure was as follows: on each microscope reference slide containing only pollen of one species of Urticaceae, an area rich in pollen was identi ed by eye and this area was automatically scanned using multidimensional acquisition with the Zeiss software Zen BLUE. For areas that were very rich in pollen, a user-de ned mosaic was created consisting of all the tiles to be scanned (e.g. 20x20 tiles), while a list of XY positions was used for microscopic slides less rich in pollen. Because pollen grains are 3-D shapes, catching all important features can only be achieved using different focal levels, so-called 'Z-stacks'. A total of 20 Z-stacks were used in this study with a step size of 1.8 µm. The settings used for scanning were a Plan Apochromat 100x (oil) objective and numerical aperture 0.55 with a bright eld contrast manager. To maintain similar conditions in the image collection process, the condenser was always set to 3.3V with an exposure time of 28ms.

Reference Pollen Image Library
All images were post-processed in ImageJ (Fiji) 36 using the script Pollen_Projector (https://github.com/pollingmarcel/Pollen_Projector). The input for this script is a folder containing all raw pollen images (including all Z-stacks), and the output is a set of projections for each individual pollen grain that are subsequently used as input for the VGG16 deep learning model.
Pollen_Projector identi es all complete, non-overlapping pollen grains and extracts them as stacks from the raw Z-stack. This is achieved using binarization on the raw images to detect only those rounded objects with a circularity >0.3 and a size larger than 5µm. Out-of-focus images within each group of 20 Zstack slices were removed using a threshold for minimum and maximum pixel values. The conventional input of a convolutional neural network is a three-channel image. In colour images RGB channels are commonly used, but since we use grayscale images, three different Z-stack projections were chosen to represent the three different channels. The projections used are Standard Deviation, Minimum Intensity and Extended Focus. Standard Deviation creates an image containing the standard deviation of the pixel intensities through the stack, where positions with large differences appear brighter in the nal projection. Minimum intensity takes the minimum pixel value through the stack and uses that for the projection. Finally, the Extended Focus projection was created using the 'Extended_Depth_of_Field' ImageJ macro of Richard Wheeler (www.richardwheeler.net) 37 . This macro takes a stack of images with a range of focal depths and builds a 2D image from it using only in focus regions of the images. A schematic overview of the processes behind the Pollen_Projector script is shown in Supplementary Figure S2.

Convolutional Neural Network
Convolutional Neural Networks (CNN) are widely used in the eld of computer vision for image classi cation, object detection, facial recognition, autonomous driving, etc. For our work we used the pretrained VGG16 network 38 in Keras 39 . Compared with traditional neural networks and shallow convolutional neural networks, VGG16 has deeper layers that extract more representative features from images. A feature extractor and classi er are two key structural parts of the CNN that perform the classi cation task. In order to train a CNN model, a large number of tagged data sets are fed into a model, to train the model to learn more features to be able to distinguish the images.
The VGG16 network contains 13 convolutional layers that form ve blocks, which generate features from images in the feature extraction phase. During the training process, parameters of convolutional layers were derived from the pre-trained network on the ImageNet dataset. Subsequently, three fully connected layers were built and added to the convolutional layers to classify the different classes (Supplementary Figure S3). To improve the effectiveness, robustness and generalization ability of the VGG16 model, as well as to prevent over tting, 10-fold cross-validation was applied in the training process. The complete pollen image dataset was split into a training data set (90%) and a test data set (10%). The training data set was split into ten subsets and used for 10-fold cross-validation. For each fold, the number of epochs was set to 30. The accuracy of the model converged at this point and the model is therefore found not to be over tting. For each fold, a model was trained using nine of the folds as training data and validated on the remaining part of the data. After the 10-fold cross-validation, 10 models were obtained and tested on the test dataset. In order to try to improve the accuracy of the model further, we further compared hard voting (or majority voting) to the standard average measure for the ten models. Hard voting is summing the votes for class labels from each model and predicting the class with the most votes.
In order to quantify model accuracy, several commonly used performance measures were used: where TP refers to true positives, TN to true negatives, FP to false positives and FN to false negatives. Recall is the number of true positives divided by the total number of elements that belong to the correct class, which is the sum of the true positives and false negatives. The F1-score is the weighted average of the precision and recall. The correct classi cation rate (CCR) re ects the accuracy of the model.

Data augmentation
A large number of images for each class is required to train a deep learning model, as the performance will increase when more variation is fed to the model. Due to the nature of the images investigated in this study, the model was sensitive to small changes, since the differences between the genera are very subtle. Therefore, data augmentation was tested to try to increase the variety of pollen images used as input. We selected the augmentation options brightness and ip. Brightness range was set from 0.1 to 2, with <1 corresponding to a darker image and > 1 to a brighter image. Horizontal-and vertical ip were also applied randomly (Supplementary Figure S4). The results were compared to running the model without data augmentation.   territory, city or area or of its authorities, or concerning the delimitation of its frontiers or boundaries. This map has been provided by the authors.

Supplementary Files
This is a list of supplementary les associated with this preprint. Click to download.