dataesr/scientific-paragraphs-categorization huggingface dataset:
k-NN based classifier Following  [Weber et al. , 2021] , we built a k-Nearest Neighbours (k-NN) based classifier to distinguish TCRs recognizing a specific antigen from generic ones. It works by computing the Levenshtein distance of the CDR3β under analysis with respect to a set of positive examples (training set of antigen-specific CDR3β) and to a set of negative examples (training set of bulk repertoire CDR3β). In the computation of the Levenshtein distance, the V and J segments are used as well (different V or J segments increase the distance by 1). Then the average distances from the k nearest neighbours are computed for both the positive and negative examples, and a score is computed as the difference between the two average distances. The only parameter of this model, k, has been fixed for the largest dataset of positive examples (TCRs reactive to the NLVPMVATV peptide) by cross-validation.
We split the positive and negative data in 50 independent training and test sets (with respectively 80% and 20% of the sequences); from each training set, a portion of 20% of its data is held out as validation set. Next k, for each partition, is fixed by maximizing the model performance on the validation set, and the model with the best k is evaluated on the test set. For the other datasets of positive examples, which are smaller, we used k = 26, obtained as the k for which the average AUC on the validation set is maximal (Figure  S18 ). The performance of the k-NN based classifier was evaluated as described in section 5. 5.
NetTCR-2. 0 NetTCR-2. 0  [Montemurro et al. , 2021]  is a convolutional neural network method trained from sets of pairs of peptides and cognate TCR sequences (β chain only or α + β chains). To re-train and evaluate it with our datasets in the β chain only version, we replaced the sets of CDR3β specific to the 4 epitopes of interest (M1 58 , pp65 495 , BMLF1 280 , S 269 ) with our VDJdb-derived datasets. We kept all the other NetTCR training data (the 3 peptides different from our 4 peptides of interest with their specific CDR3β as well as their control data for these additional peptides). Re-training NetTCR-2. 0 with only the 4 epitope-specific sets of CDR3β we considered gave AUCs in the same range of values as for the case with 7 peptides. In all cases our training of NetTCR-2. 0 uses the same default settings as  [Montemurro et al.
Classification performance Having found the optimal diffRBM architecture (section 5. 3) for the immunogenicity model, for each of the 3 HLA types considered, we trained 50 HLA-specific models on the original training sets (consisting of the 80% of the full datasets available) and we assessed their average performance over the corresponding 50 choices of the test set. In particular, we tested the ability of the HLA-specific immunogenicity models to identify new immunogenic peptides by the Receiver Operating Characteristic curve (ROC). For each of the 50 repetitions, we assigned scores of immunogenicity predicted by a given HLA-specific model (given by the diffRBM units' score L d of equation 11) to the sequences of the test set of positives (immunogenic peptides with the HLA type under consideration) and of negatives (a test set of non-immunogenic peptides presented by the same HLA).
Varying the threshold score value to discriminate positives from negatives, we obtained the ROC curve describing the fraction of immunogenic peptides predicted by the models' scores, against the fraction of predicted non-immunogenic ones. We took the Area Under the Curve (AUC) as a metric of the models' ability to discriminate immunogenic from non-immunogenic peptides. We performed the same validation for all the RBM-based approaches (section 5. 1, Figure  6 ), using their corresponding output scores (equation 9 for the full RBM, 10 for background RBM). The performance of the full RBM obtained through the score (9) is by far and large equivalent to the one of an RBM with the same hyperparameters entirely trained, in one step only, from the selected dataset (Figure  S10 ), confirming that there is a gain in performance with the differential learning strategy only when we focus on the differential units and their parameters.
, 2017]  or  [Dean et al. , 2015] ), randomly drawn at each repetition with the same size as the positive test set and with no overlap with the 10 6 sequences of the training set of the background model. All the approaches (diffRBM and alternatives versions, SONIA, k-NN, NetTCR-2. 0, TCRex) are trained and tested on 50 independent random partitions of both positives and negatives into training and test sets, and the performance shown in Figures 6, S16, S17, S19 is the average AUC over these 50 partitions. As negative set to train the supervised approaches (k-NN, NetTCR-2. 0) we took, for each training repetition, a subset of the bulk-repertoire dataset from  [Emerson et al. , 2017]  with the same size as the positive training set.
PRIME tool PRIME  [Schmidt et al. , 2021]  was downloaded from https://github. com/GfellerLab/PRIME. Being it not possible to re-train PRIME on our datasets for a fair comparison of performance, we simply evaluated it on the set of immunogenic and non-immunogenic peptides we collected and we obtained a discrimination performance with AUC=0. 53 for peptides specific to HLA-A*02:01, AUC=0. 45 for HLA-B*07:02 and AUC=0. 56 for HLA-B*35:01. Since we set out to predict immunogenicity conditioned on binding to a given HLA allele, PRIME was run in the mono-allelic mode (e. g. with option -a A0201 for the case of HLA-A*02:01); in general, however, different results are obtained by adding more alleles, where then the best presenting allele according to the predictor is taken.
4. 2) and hence represent position-specific amino acid frequencies. In the case of the immunogenicity model, the prediction by (18) simply describes the log ratio between the frequency of amino acids in immunogenic peptides relative to the amino acid frequency of all presented peptides from the background dataset (we refer to it as 'AA frequency ratio (to all presented)', see Figure  2C-D ); in the case of the TCR specificity model  (18)  gives the log ratio between the frequency of amino acids in peptide-specific CDR3β repertoires relative to the amino acid frequency in bulk repertoire (we refer to it as 'AA frequency ratio (to bulk)', see Figure  5B-C ). In both cases the prediction indicated as 'Amino acid frequency' is obtained by estimating and ranking the single-site factors: T i (σ i ) = log P i (σ i ) (19) In the case of the immunogenicity model, we have also learnt an independent-site model from nonimmunogenic data (that we call P NI (σ i )) and we have looked at the single-site factors: T i (σ i ) = log P i (σ i ) P NI i (σ i ) (20) which represents a log ratio of amino acid frequency in immunogenic data relative to non-immunogenic data (this prediction is labelled 'AA frequency ratio (to non-immunogenic)' in Figure  2C -D and 3B).
e. , whereby TCR response could not be recovered even at the highest concentrations). For a mutation in sequence σ at position i changing σ i to σ i , we estimated the model prediction of the mutation cost as: Predicted mutation cost = T i (σ i ) -T i (σ i ) (25) where we took T i as the single-site factors of (  17 ) for lethal mutations (Figure  3E ) and of background RBM (21) for non-lethal mutations (Figure  3G ). To assess whether the distribution towards positive values for lethal mutations is significantly higher than the expectation for generic, non immunogenicity-impacting mutations, we estimated the mutation cost distribution of a 'control' case (see Figure  3E ) where we drew at random 3000 HLA-A*02:01presented peptides from the background dataset, and we calculated the costs of all possible amino acid substitutions at each position.
85). This result is confirmed by scoring the peptides via NetMHCpan4. 1 (AUC=0. 85). Excluding this case, the average diffRBM units AUC=0. 68, which is comparable to the AUC obtained via standard validation (average AUC=0. 66, Figure  4B ). Upper row: the AUC of classification of immunogenic vs non-immunogenic peptides given by the full RBM (trained in part on the background dataset and in part on the selected data and indicated as 'RBM (differential)') is essentially the same as the one of an RBM with the same hyperparameters trained only on selected data (denoted as 'RBM (non-differential)'). The comparison is made separately for each HLA. In these scatter plots every point gives the AUC of discrimination of the two models compared for one of the 50 random partitions into training and test sets of the corresponding HLA-specific dataset (Methods).
Bottom row: we consider the AUC of immunogenic vs non-immunogenic classification obtained through the difference of scores between the diffRBM units trained on immunogenic peptides and the ones trained on nonimmunogenic peptides. The scatter plots show that such performance (on the y-axis), as it should, is equivalent to the one attained by the difference of the full RBM log-likelihoods (x-axis). The hyperparametric search for the classifier is based on the optimal performance on the validation dataset at discriminating immunogenic vs non-immunogenic peptides as measured by the AUC. We show the mean AUC over 50 random partitions into training set (A) and validation set (B). The labels on the x-axis denote different classifier architectures, which are described in Table  S2. Performance at discriminating immunogenic vs non-immunogenic peptides (like in Figure  4 ) where a reweighting scheme based on sequence similarity is applied during training (Methods).
Given that NetTCR2. 0  [Montemurro et al. , 2021]  does not account for V and J type in its input, we report the performance of the diffRBM units, full RBM and background RBM for the 4 epitope-specific models under consideration when only the CDR3β amino acid sequence has been used to train and test the models. AUCs are averages over 50 training-test partitions and error bars give the corresponding standard deviation (Methods). A: Entropy S (expressed in nats, see Methods) of the space of HLA-specific presented antigens (evaluated by background RBM) and of HLA-specific immunogenic antigens (evaluated by the full RBM) for the 3 HLAs. Error bars represent the uncertainty on the estimated entropy coming from limited sampling and was calculated as in  [Marchi et al. , 2019]. B: The entropy of the background dataset (CDR3β bulk repertoire) obtained from background RBM is compared to the entropy of epitope-specific CDR3β obtained from the full RBM models of TCR specificity to GILGFVFTL, NLVPMVATV, GLCTLVAML and YLQPRTFLL.
However, the uncertainty from sampling is substantial (especially for the smallest datasets HLA-B*07:02, HLA-B*35:01 and YLQPRTFLL), and prevents us from estimating the reduction in diversity with high accuracy. The first 3 rows refer to the HLA-specific models of immunogenicity (entropy values plotted in A) while the last 4 rows refer to the epitope-specific models of TCR specificity (entropy values in B). 17.
Figure1: Cartoon of the differential RBM (diffRBM) learning approach. A: The parameters of background RBM (gray) are learnt from the 'background' sequence dataset. B: The diffRBM units (gold) are learnt from a small subset of 'selected' sequences. C: We consider the application of diffRBM to modeling peptide immunogenicity or T-cell receptor (TCR) antigen specificity, whereby the background dataset consists, respectively, of all antigens presented by a given Human Leukocyte Antigen class I complex (HLA) or of generic TCRs from the bulk repertoire. D: The selected sequences correspond to HLA-specific antigens validated to be immunogenic or to TCRs that are antigen-specific responders. The inferred parameters associated to the diffRBM units allow one to identify putative contact positions in the peptide-HLA-TCR structure (E) and more generally to assign scores that distinguish the selected from the background sequences (F).
Figure3: DiffRBM units encode molecular-level features of immunogenicity. A: Enrichment in amino acid usage in immunogenic peptides relative to non-immunogenic ones calculated as a log frequency ratio across the central positions (4-8) for each HLA type. The color code used in the x-axis labels indicates amino acid properties: red = negatively charged (E,D), blue = positively charged (H, K, R), purple = non-charged polar hydrophilic (N, T, S, Q), golden = aromatic (F, W, Y), black = aliphatic hydrophobic (I, L, M, V), green = cysteine (C), grey = tiny (A, G, P) or other. B: Distribution of the diffRBM single-site factors evaluated across the HLA-A*02:01-specific immunogenic sequences with W at position 5. For all amino acid/position combinations that were found to be key to immunogenicity, the diffRBM single-site factors assume predominantly positive values (see also FigureS6A).
Figure4: Performance at discriminating immunogenic vs non-immunogenic peptides. A: For a given HLA (e. g. HLA-A*02:01), we assign diffRBM units' scores to held-out sets of immunogenic and non-immunogenic antigens and we measure the discrimination performance via the Area Under the Curve (AUC), see Methods. B: Performance of diffRBM units, full RBM, background RBM and other methods for the HLA-A*02:01-specific dataset as measured by the AUC. We have divided the methods into unsupervised and supervised ones (color code). Unsupervised methods are meant as the methods trained on immunogenic (or presented) only peptides to learn an unsupervised model. Supervised methods are trained on both immunogenic and non-immunogenic to optimally discriminate them. The approach 'diffRBM units (difference)' is intermediate: it exploits the annotation of peptides as immunogenic and non-immunogenic but it is not optimized for the discrimination task.
(resulting in the choice of HLA-A*02:01, HLA-B*07:02 and HLA-B*35:01). The size of the final datasets of immunogenic peptides is: |D s | = 1682 for HLA-A*02:01, |D s | = 258 for HLA-B*07:02, |D s | = 215 for HLA-B*35:01. Sets of non-immunogenic peptides consist of 2301 sequences (HLA-A*02:01), 807 (HLA-B*07:02), 166 (HLA-B*35:01).
spectrometry techniques. The resulting training dataset sizes for the background models are |D b | = 4265 for HLA-A*02:01, |D b | = 1006 for HLA-B*07:02, |D b | = 1211 for HLA-B*35:01. For consistency with the type of datasets used to train the presentation model by RBM-MHC, scores of presentation from the algorithm NetMHCpan4. 1.
FigureS11: Hyperparametric search for the classifier of immunogenicity. The hyperparametric search for the classifier is based on the optimal performance on the validation dataset at discriminating immunogenic vs non-immunogenic peptides as measured by the AUC. We show the mean AUC over 50 random partitions into training set (A) and validation set (B). The labels on the x-axis denote different classifier architectures, which are described in TableS2.
FigureS20: Model-based entropy estimation. A: Entropy S (expressed in nats, see Methods) of the space of HLA-specific presented antigens (evaluated by background RBM) and of HLA-specific immunogenic antigens (evaluated by the full RBM) for the 3 HLAs. Error bars represent the uncertainty on the estimated entropy coming from limited sampling and was calculated as in [Marchi et al. , 2019]. B: The entropy of the background dataset (CDR3β bulk repertoire) obtained from background RBM is compared to the entropy of epitope-specific CDR3β obtained from the full RBM models of TCR specificity to GILGFVFTL, NLVPMVATV, GLCTLVAML and YLQPRTFLL. C: The entropies calculated from background RBM and the full RBM plotted in A and B is lower than the one estimated from independent-site models of the same data (Entropy PWM), because RBM models can account for correlations between sequence sites, hence for additional constraints on sequence diversity.
Colors are the same as in A, B. D: Reduction in diversity in selected data. For each type of data, given the entropy S, exp(S) gives an estimate of the effective size of the corresponding sequence space, hence the ratio exp(S RBM )/ exp(S Background ) gives an estimate of the effective fraction of sequences retained in selected data compared to background ones. The table shows that the reduction in diversity quantified by this fraction is typically less than one order of magnitude, with the exception of the YLQPRTFLL-specific CDR3β repertoire. However, the uncertainty from sampling is substantial (especially for the smallest datasets HLA-B*07:02, HLA-B*35:01 and YLQPRTFLL), and prevents us from estimating the reduction in diversity with high accuracy. The first 3 rows refer to the HLA-specific models of immunogenicity (entropy values plotted in A) while the last 4 rows refer to the epitope-specific models of TCR specificity (entropy values in B).
Acknowledgements We thank  Julien Racle ,  Nathanael Spisak ,  Francesco Camaglia ,  Marta Łuksza ,  Morten Nielsen  and  Vadim Karnaukhov  for discussions. We are particularly grateful to  Giulio Isacchini  for help with the data and software. A. D. G. is supported by the  European Union 's  Horizon 2020 research and innovation programme  under the  Marie Sklodowska-Curie  grant agreement No  101026293. JFdCD acknowledges funding from  ANR-19  Decrypted  CE30-0021-01  and  ANR-17 RBMPro CE30-0021-01. This research was also supported by the  European Research Council   COG 724208  and  ANR-19-CE45-0018  ' RESP-REP ' from the  Agence Nationale de la Recherche.
Simple Summary: Wing venation traits are used to identify honey bee subspecies. While several wing-based tools are available, they suffer from weaknesses that were addressed by the recently developed software DeepWings©. This software allows fully automated identification of wing images in a friendly, free, and rapid manner. Here, we sought to test DeepWings© on 14,816 wing images representing 2601 colonies sampled in the native areas of three widespread subspecies in Europe: the Iberian honey bee (Apis mellifera iberiensis), the dark honey bee (Apis mellifera mellifera), both belonging to the M lineage, and the Carniolan honey bee (Apis mellifera carnica), belonging to the C lineage. DeepWings© classification of these colonies largely matched the endemic M and C lineages, with proportions of 71.
One of the most intensive assessments of wing shape variation is provided by the Discriminant Analysis with Numerical Output (DAWINO) method. This method is based on 30 characters extracted from vein lengths, their ratios, and vein angles. On the other side of the spectrum are the methods that only require estimations of the Cubital Index (CI), Hantel Index (HI), and/or Discoidal Shift Angle (DSA). These are popular amongst many beekeepers involved in the conservation of A. m. mellifera  [34]  as they are simple and can be implemented by semi-automatic software, such as BeeMorph (http://www. hockerley. plus. com/ accessed on 1 July 2022 ) or CBeeWings (http://www. cybis. se/cbeewing/ accessed on 1 July 2022). The problem is that these software packages provide a less reliable identification than that obtained from character-intensive methods.
In addition, they require manual annotation of the vein junctions from which the indexes are calculated, which is time-consuming and prone to error. Another way of assessing wing shape variation is through wing geometric morphometrics (WGM). This method is recognized as robust and reliable in insect taxonomy and is widely used in honey bee subspecies identification for varying purposes, including conservation  [36] [37] [38] [39] [40] [41] [42]. WGM uses the coordinates defined by 19 landmarks located in the vein junctions to capture variation in wing shape  [43]  and can be implemented by the software IdentiFly  [36]. However, IdentiFly is a semi-automatic tool that requires several steps before the wings can be fully identified, which makes its use difficult for the layman. The most recent advance in WGM uses deep learning to automatically extract the 19 landmarks from the right forewing of honey bee workers  [44].
The approach is implemented by the software DeepWings©, which allows fully automated identification of honey bees. DeepWings© is very friendly, requiring only that the users drag the wing images into a file drop zone. Then, for each analyzed wing, it automatically retrieves the classification probabilities for the top three subspecies and the estimation of CI, HI, and DSA along with the landmark coordinates. DeepWings© is offered as a free web service at https://deepwings. ipb. pt. (accessed on 1 November 2022) Herein, we employed DeepWings© to identify 2601 colonies from the analysis of 14,816 wings. These colonies were located in 15 countries, covering the native ranges of A. m. iberiensis, A. m. mellifera and A. m. carnica. In addition, we compared the wing shape data with molecular data for a subset of A.
m. ligustica, or A. m. caucasia and corresponding lineage. The software retrieves the three highest prediction probabilities for each wing individually and for an average wing estimated from the wing batch dragged into the file drop zone. DeepWings© constructs the average wing by averaging the coordinates of each of the 19 landmarks across all the wings processed in one batch  [44]. When all the wings from a colony are simultaneously uploaded, as is carried out here, the estimated average wing represents the colony and the classification at the colony level can be retrieved for the average wing. We used the wings output at both individual and colony levels. The wings were assigned to one of the five subspecies based on the highest prediction probability, even if the probability was low.
Association between Wing Data and Molecular Data The association between DeepWings© classification and that obtained from molecular markers (microsatellites and SNPs) was assessed for 1214 colonies sampled from the native ranges of A. m. iberiensis (Portugal, Spain) and A. m. mellifera (France, Ireland, Wales, Russia). For each colony, the highest prediction probability of belonging to M lineage, as inferred from the 19 landmark coordinates using DeepWings©, was compared against the corresponding M-lineage membership proportion, as inferred from different sets of microsatellites or SNPs (Table  S1 ) using the software Structure  [45]  or Admixture  [46] , respectively. The molecular dataset was generated in previous works (see Table  S1 ).
Statistical Analysis The probability data obtained from the average wing for the subspecies identified by DeepWings© did not follow a normal distribution, as per the Kolmogorov-Smirnov test. Accordingly, the summary statistics were presented for each location as medians and interquartile ranges (median, interquartile range [IQR]; Tables  S2  and S3 ). The distributions of the probability data points were compared among the identified subspecies in each dataset using the Mann-Whitney U test or the Kruskal-Wallis test, followed by the Dunn's multiple comparison test with statistical significance levels (p) adjusted by Bonferroni. The association between the probability of belonging to the M lineage, as inferred by DeepWings© from wing shape data, and the membership proportion in the M lineage, as inferred from microsatellite or SNP data, was assessed using the Spearman's rank-order correlation coefficient (r).
Classification of the Total Wings Dataset A total of 14,816 worker forewings, representing 2601 colonies, were processed by DeepWings©. From these, 856 (5. 8%) were rejected by the software as the 19 landmarks could not be annotated due to different image problems (Figure  S1 ). The rejected subset included 106 wings from Poland and 7 from Hungary. Since the colonies of these two countries were mostly represented by single wings, the total number of classified colonies decreased to 2488. In 82. 8% (709) of the discarded wings, rejection was mainly due to the very low resolution of the images and noisy background, leading DeepWings© to aggregate close landmarks. The remaining 147 (17. 2%) images displayed some kind of corruption, including missing landmarks (9, 1. 1%), folded or twisted wings (9, 1.