Probabilistic data integration identifies reliable gametocyte-specific proteins and transcripts in malaria parasites

Plasmodium gametocytes are the sexual forms of the malaria parasite essential for transmission to mosquitoes. To better understand how gametocytes differ from asexual blood-stage parasites, we performed a systematic analysis of available ‘omics data for P. falciparum and other Plasmodium species. 18 transcriptomic and proteomic data sets were evaluated for the presence of curated “gold standards” of 41 gametocyte-specific versus 46 non-gametocyte genes and integrated using Bayesian probabilities, resulting in gametocyte-specificity scores for all P. falciparum genes. To illustrate the utility of the gametocyte score, we explored newly predicted gametocyte-specific genes as potential biomarkers of gametocyte carriage and exposure. We analyzed the humoral immune response in field samples against 30 novel gametocyte-specific antigens and found five antigens to be differentially recognized by gametocyte carriers as compared to malaria-infected individuals without detectable gametocytes. We also validated the gametocyte-specificity of 15 identified gametocyte transcripts on culture material and samples from naturally infected individuals, resulting in eight transcripts that were >1000-fold higher expressed in gametocytes compared to asexual parasites and whose transcript abundance allowed gametocyte detection in naturally infected individuals. Our integrated genome-wide gametocyte-specificity scores provide a comprehensive resource to identify targets and monitor P. falciparum gametocytemia.

may require interventions that specifically aim to prevent malaria transmission. Malaria transmission depends on male and female gametocytes, the sexually reproducing forms of the Plasmodium parasite that are ingested by blood-feeding Anopheles mosquitoes. In the mosquito gut, gametocytes may complete the parasite's reproductive cycle and, following sporogonic development, render the mosquito infectious. Factors that govern gametocyte production and infectivity remain poorly understood. Whilst recent studies have shed light on the processes controlling gametocyte commitment 7,8 , commitment and maturation of gametocytes may differ between infections and over the course of infections, under influence of environmental and host factors 9,10 . A better understanding of gametocyte dynamics during infections, as well as the development of tools to monitor or target gametocytes, may be informed by high-throughput protein and transcriptome studies 11,12 . In the past 15 years, a number of large-scale studies on Plasmodium gametocytes have been reported: the proteome of Pf and the rodent malaria parasite Plasmodium berghei (Pb) have been examined by mass spectrometry [13][14][15][16][17][18][19][20][21] , and the transcriptome of both species by microarray and RNA sequencing 15,20,[22][23][24][25][26][27][28] . These studies differed in their focus and resolution in examining (sexual) developmental stages and each faced challenges in detecting low abundance proteins 29 and by the purity of parasite populations 13,14 . The use of fluorescent parasites and fluorescence-assisted sorting of staged parasites have recently permitted a better discrimination of proteins in either male or female gametocytes 16,17,20 and have allowed more detailed comparisons of Plasmodium life-stages. However, individual studies are still vulnerable to imperfect sample purity, and other sources of uncertainty such as correct gene identification for accurate peptide assignment. These technical and methodological challenges lead to discrepancies between individual studies and hamper firm conclusions about gametocyte-specificity of proteins and transcripts.
We utilized the numerous published proteomics and transcriptomics Plasmodium data sets in a comprehensive data integration framework to obtain a consensus of gametocyte-specific transcripts and proteins. Our data integration approach is an adaptation of the naïve Bayesian classifiers that have previously been applied in the prediction of protein interactions and components of cellular systems 30,31 . The framework calculates probabilities that any given transcript or protein is gametocyte-specific given the evidence presented across the total of transcriptomics and proteomics data. A key aspect of the methodology is that it takes into account the predictive power of each contributing data set: it assigns weights to data sets based on their ability to distinguish gold standard lists of gametocyte and asexual proteins. These we have constructed using existing literature where life-stage specificity was confirmed using classical "non-omics" approaches (e.g. protein detection in immunofluorescence-assays, functional/genetic studies), followed by expert curation. The most informative data (from datasets with the highest discriminative power against the gold-standard lists) will thus contribute most to the predictions, while less informative data are down-weighted. This allows for (i) the resolution of conflicting evidence without disregarding data, and (ii) the construction of a transparent scoring system in which the relative contribution of each data set is directly visible.
Using this approach, we propose a robust gametocyte-specificity score for all Pf genes that allows a consensus list of gametocyte-specific genes at protein and transcript level. We illustrate the utility of our findings by examining naturally acquired responses to newly identified gametocyte-specific proteins in gametocyte-carriers and non-carriers by protein microarray. In addition, we confirmed gametocyte-specificity for a selection of gametocyte-specific transcripts using culture material from geographically distinct Pf strains and samples from naturally infected malaria patients.

Results
Weighted integration of proteomics and gene expression data. Using Bayesian statistics, we integrated Plasmodium mass spectrometry and transcript datasets from 18 different studies (Table 1) on P. falciparum (Pf; n = 14), P. berghei (Pb; n = 3) and P. vivax (Pv; n = 1). Since gametocyte biology differs between Plasmodium species, scores were calculated for the total set of Plasmodium studies and for Pf only. Unsupervised clustering of genes based on peptide counts or mRNA expression resulted in grouping according to data acquisition method rather than parasite stage ( Fig. 1), illustrating the necessity of a supervised approach to discriminate between gametocyte-specific and non-gametocyte genes. To objectively assess the value of individual data sets and allow their assembly into a gametocyte-specificity score, we created a gold standard that served as a benchmark for every sample. This gold standard was collected from literature review and comprises two lists; one of asexually expressed proteins, mainly blood stage but also sporozoite and liver stage, and one of known gametocyte proteins (Supplementary Table S1). A gametocyte-specificity score was then derived for each gene by comparing its expression in all studies to the relative expression of the gametocyte and asexual gold standards in those samples ( Supplementary Fig S1-2). Proteins or RNAs detected in a study with high discriminative power for gametocyte and asexual gold standard genes (ratio of gametocyte to asexual gold standard genes) received higher gametocyte-specificity scores than those detected in a study with lower discriminative power. The individual log-transformed scores per gene were combined for proteomics and transcriptomics data separately. Scores for Pf-only studies (Supplementary Table S2) and all combined data sets were highly correlated ( Fig. 2D; Pearson's r = 0.9867 and r = 0.9514 for proteomics and transcriptomics, respectively) and the latter were used in the remainder of the manuscript.
The distribution of scores for proteins and transcripts are presented for the two sets of genes of the gold standard as well as for all other genes ( Fig. 2A). As expected, the gametocyte and asexual gold standard set of genes are perfectly separated by their respective proteomics-derived scores and show only little overlap in their transcriptomics-derived scores (Fig. 2B). The shift in the density peak of the proteomics compared to transcriptomics is due to an inherent property of the method that gives a negative score for gametocyte-specificity to all proteins that were not detected in the proteomics studies (n = 1583). The highest scoring 100 genes for proteomics and transcriptomics contained 26 (63.4%) and 15 (36.6%) of the 41 gold standard gametocyte genes, respectively (Fig. 2C), indicating both the respective discriminating power of the gold standard and that many other genes are as specific for gametocytes as the highest gold standard representatives. Translationally repressed genes are common in late stage female gametocytes 32,33 and are detectable by high transcriptomics and low proteomics score. In our analysis 461 genes have this profile (Supplementary Table S3), including genes that are known to be translationally repressed like Pfs28 32,34 and 186 genes with a previously reported bias towards expression in female gametocytes 20 .
Cross-validation illustrates the improved predictive power of the integrated data. Ten-fold cross-validation was performed using random subsamples of the gold standard lists to predict the ranks of left-out genes. The resulting proteomics ranking shows near perfect sensitivity, with all but two gold standard gametocyte genes ranking higher than the gold standard asexual genes (Fig. 3). The added value of our integrated approach is illustrated by the receiver operating characteristic curve where the integration of data sets gave higher sensitivity and areas under the curve for both proteomics and transcriptomics than any individual study (Fig. 3). Using the Bayesian integration based on the complete gold standards, we ranked all Pf proteins by giving them a gametocyte-specificity score (Supplementary Table S2). All proteins with a score >5 (n = 602) were considered gametocyte-specific. Most of these have not consistently been described as "specific" or "enriched" in gametocytes in the original data sets ( Fig. 4A and Supplementary Table S4). Previous studies defined 315 13 to 1725 20 proteins as gametocyte-specific for Pf. Not only did our integrated approach lead to a better recovery of gold standard listed known gametocyte proteins, we also identified 178 genes with undescribed function as gametocyte-specific (Supplementary Table S2). We further identify a number of proteins as gametocyte-specific even though they had been reported as asexual by previous studies (Fig. 4B).
A recent proteomics study of male and female Pf gametocytes 35 , not included in our original analysis, was used to test the robustness of our scores. When we included this data set in our final Bayesian proteomics scores, both gametocyte scores and gene ranks before and after addition of this data set were highly correlated (Pearson's r = 0.997 and Spearman's rho = 0.995, respectively). Furthermore, the top 100 gametocyte proteins did not change and the top 602 proteins were 96% identical (578 of 602). Taken together, cross-validation and independent data suggest that the integrated gametocyte-specificity score is robust and contains potential novel gametocyte markers.
Predicted gametocyte-specific proteins are recognized by gametocyte-carriers. As an illustration of the utility of gametocyte-specific proteins as markers of gametocyte exposure, we utilized protein microarray data from a study that aimed to characterize the immune profile associated with transmission-reducing immunity in naturally infected gametocyte carriers (Stone, Campo et al. accepted manuscript 36 ). For the current study, we compared responses to our gold standard gametocyte genes (n = 40) and novel gametocyte genes from our 100 highest scoring proteins that were on the array (n = 30). Antibody prevalences for these genes were compared between Gambian gametocyte carriers and Gambians who carried asexual parasites but not gametocytes as determined by microscopy. Antibody responses to the predicted gametocyte-specific proteins were significantly higher in gametocyte carriers (p = 0.005), while for the gold standard antigens this difference was less significant (p = 0.058) (Fig. 5A, Mann-Whitney U test). When antigens were analysed individually, a significantly higher antibody prevalence in gametocyte carriers was detected for five novel gametocyte antigens (Fig. 5B, p < 0.05 in Fisher's exact, corrected for multiple testing, Supplementary Table S5, Supplementary Fig S3). Only two of these five have an assigned function -a DNA ligase, and Gamete egress and sporozoite-traversal protein (GEST). For two of the three remaining Plasmodium proteins, we were able to predict a function based on homology, using the sensitive homology detection tool HHpred 37 . PF3D7_1251000 is homologous to the co-chaperone HSP20 heat shock protein and PF3D7_1439600 is homologous to the MLRQ subunit of complex IV of the oxidative phosphorylation, underlining the enrichment in mitochondrial proteins as discussed below (Supplementary Table S6 includes homology predictions for all conserved, highly gametocyte-specific Plasmodium proteins).

Gametocyte-specific RNA transcripts detect (submicroscopic) gametocyte carriage. Of the 100
highest-scoring transcripts, 15 non-gold standard candidates were selected for qRT-PCR validation based on their gametocyte scores in a preliminary analysis ( Table 2). Mature gametocytes of four Pf strains from different geographical origins were compared to asexual blood stage parasites. The minimum transcript abundance difference (Ct Asexuals − Ct Gametocytes ) ranged from 4.76 to 14.95 ( Fig. 6A and Supplementary Table S7 qPCR & primers), reflecting 27.1 to 31,500-fold higher transcript numbers in gametocytes compared to asexual parasites and confirming pronounced upregulation of all selected targets in gametocytes. With a very conservative threshold of 1,000-fold enrichment in three of the four strains tested (Fig. 6A), eight of the 15 tested transcripts were highly specific to gametocytes. Transcript abundance in ring-stage parasites was assessed and compared to Pfs25 mRNA, an established and highly abundant yet intron-less female gametocyte specific transcript 11,38 . Five out of eight gametocyte specific transcripts were undetectable in asexual ring stages at ≤10 5 parasites/mL, similar in specificity to Pfs25 (Fig. 6B); the five most sensitive gametocyte markers detected gametocytes across the range of 10 2 -10 6 gametocytes/mL (Fig. 6C). In RNA samples from a previously reported clinical trial conducted in Kenya 39 , all eight gametocyte markers detected gametocytes at densities below 10 3 /mL (Fig. 6D).
Gametocyte-specific proteins are enriched for cytoskeletal movement and metabolism functions. To uncover novel characteristics underlying gametocyte function, we analyzed over-represented gene ontology (GO) terms in our integrated consensus gametocyte proteins. The 100 highest ranked proteins were examined for enrichment of GO terms that reflect specific biological functions. Microtubule based movement, metabolism of carboxylic acids and metabolism of nucleic acids were highly enriched among gametocyte Figure 1. Clustered data sets used in this study with genes ranked according to their protein or transcript expression. Level of expression as detected in the respective samples with unique peptide counts for MS data and percentiles for transcriptomics. The studies are clustered using complete linkage according to their overall gene expression similarities (Euclidean distance). See Table 1 for study keys. Distribution of asexual(a)/ gametocyte (g) samples (red/blue) is shown in top bar, proteomics (P) and transcriptomics (T) (dark/light grey) in lower bar.
proteins ( Supplementary Fig. S4). Four out of the six putative Pf dynein heavy chain proteins are found back among the 100 most gametocyte-specific proteins, alongside a tubulin gamma chain and a tubulin chaperone. The importance of DNA elongation and ligation processes is reflected in GO term associations as well as in antibody response to a DNA ligase (Fig. 5B). The "classic" GO term enrichment calculation was complemented by a rank-based gene set enrichment analysis (GSEA). GSEA uses all Pf genes and their respective (proteomics derived) gametocyte-specificity score and thus does not include an arbitrary cutoff of the proteins that are or are not gametocyte-specific. It confirmed the above-mentioned results, and in addition to those the terms   Table 1 for study keys. "mitochondrial protein complex" (GO:0098798) and "TCA" (GO:0006099) were enriched, stressing both mitochondrial location and processes. Although male gametocytes carry pre-synthesized proteins to rapidly form eight motile gametes upon activation in the mosquito midgut, we did not observe flagellum associated terms in the GSEA. The reason for this is that current GO term annotation for Pf has only one gene with the "cilium" GO term for cellular component (GO:0005929; PF3D7_1025500), one for "axoneme" (GO:0005930; PF3D7_0828700) and none with the biological-process terms "cilium or flagellum-dependent cell motility" (GO:0001539) or "axoneme assembly" (GO:0035082). We supplemented the GO annotation with a list of 28 Pf cilium genes (Methods). The newly assembled "cilium" GO term now acquired the highest enrichment score in the GSEA (Supplementary Fig. S5). This may reflect the formation of the flagella of the microgamete but may also (partially) reflect intracellular trafficking and or be associated with genome replication as the term has overlap with the genes annotated for microtubule processes (Supplementary Table S5 GO terms).

Discussion
Combining proteomic and transcriptomic data from 18 sources, we present an integrated consensus score for gametocyte-specific proteins and transcripts. We predict 602 gametocyte-enriched proteins of which 186 are currently without ascribed function. We illustrate the potential utility of our gametocyte score by providing evidence for differential recognition of gametocyte proteins by naturally infected gametocyte carriers and the sensitive detection of mRNA of novel gametocyte transcripts in field samples.
The gametocyte proteome of P. falciparum (Pf ) has been assessed repeatedly. Individual lists of gametocyte-specific proteins 13,14,16,17,19,20 have unavoidable limitations related to comparator (asexual) parasite stages, sample purity, assay sensitivity and arbitrary cut-offs used to define gametocyte-specificity and show only partial agreement. To acquire a more robust gametocyte-specificity score, we integrated data from these individual studies, along with studies of purified asexual parasites and related Plasmodium species. Including gene expression data from multiple species generally increases the likelihood that the combined gene expression data reflect underlying biology, as observed in the Apicomplexa 40,41 . We applied a Bayesian classifier first applied to 'omics data by Jansen and colleagues 30 and adapted by Van der Lee and colleagues to identify genes involved in anti-viral immune responses 31 . The probabilistic approach combines the evidence from all studies in an unbiased way, without giving a priori preference of one study over another. Instead, the measurements of all studies were weighted inherently during the scoring process by assessing the retrieval of a gold standard set of genes. As these gametocyte and asexual gold standard sets are of central importance to the study, they have undergone expert curation (see Methods and Acknowledgments). The power of this integrative approach lies not only in weighting data sets by the retrieval of gold standard genes but also in the opportunity to exclude proteins from the gametocyte-specific list by appreciating their presence in (other) asexual samples. A further strength of the approach is that it allows the ranking of gametocyte proteins that have only been reported in a subset of studies. Our integration of data sets reveals that 602 proteins are likely to be specific to gametocytes although very few  gametocyte-specific proteins were detected in every underlying dataset and seven proteins had never before been reported as gametocyte-specific. A general limitation of all mass spectrometry (MS) studies is their bias toward highly abundant proteins. Proteins with low-level expression may be missed in a bulk proteome analysis. After integration of the MS studies listed in Table 1, 1583 Pf proteins were never detected, representing approximately 28% of all proteins encoded by Pf. Some of these might be of too low abundance or expressed during sporozoite or gamete, ookinete and liver stage, which are underrepresented or not included in our data, respectively. New advances in MS that include the sensitive detection of peptides from currently understudied Plasmodium life stages may shed light on these currently uncharacterized genes. In addition, approaches that focus specifically on post translational modifications like phosphorylation of proteins as has been done for asexual parasites 42-45 may add new lines of evidence towards gametocyte-specific functions of proteins. Our approach suggests that the currently available MS data is sufficiently comprehensive to identify stage-specific proteins when analysed in an integrative approach. We examined this directly by incorporating a new Pf gametocyte MS study 35 in our scoring. The authors reported 44 new gametocyte-specific proteins that were not reported by earlier studies. We compared this data set to our integrated data set and found 24 of the 44 had been detected in one or more erythrocytic stages or sporozoites 14,18,42,46,47 while 11 others had been identified in a (single) gametocyte sample before (Supplementary Fig. S6). Importantly, the scores and top 100 gametocyte genes remained unaltered by integrating this new dataset.
The ranking of gametocyte-specificity that we provide here can i) aid in understanding the biology of this life stage and ii) improve diagnostics related to gametocyte exposure and carriage. Regarding gametocyte biology, our high-ranking gametocyte-specific genes are enriched for mitochondrial, metabolism and microtubule processes and DNA replication, supporting the quality of the data integration. The enrichment of mitochondrial localization and process is consistent with what we know about the enlarged mitochondrion of gametocyte stages 48 and increased activity of the citric acid cycle 49 . DNA replication terms are highly enriched which is consistent with what happens in the subsequent life stage in which the (micro)gamete rapidly duplicates its genomic DNA three times. Regarding the use of the gametocyte score to inform gametocyte diagnostics, diagnostics can directly detect nucleic acids specific to gametocytes 11 or detect antibody responses reflecting past/recent exposure as is increasingly used for asexual P. falciparum and P. vivax parasites 50,51 . We use our integrated gametocyte list to explore its utility for both approaches. We validated 15 transcript targets in four different Pf strains, comparing transcript abundance in gametocytes and asexual parasites. All tested targets were enriched in gametocytes. Five targets were tested for their sensitivity and can recognize 100 gametocytes/mL, while the signal is undetectable when fewer than 10 5 -10 6 ring-stage parasites/mL are present. In practical terms, these markers may be used to reliably detect gametocytes at densities well below the microscopic threshold of detection in samples without high-densities of asexual parasites, similar to the gametocyte marker that is currently most widely used, the female gametocyte-specific Pfs25 38 .
As an alternative approach to the detection of gametocyte carriage in populations, we utilized a gametocyteenriched protein microarray (Stone, Campo et al. accepted manuscript 36 ) to determine antibody responses to genes that we here describe as highly gametocyte-specific. The bacterial expression system used for the array has known limitations with the expression of conformational proteins 52 and should thus be considered a 'rule in' rather than 'rule out' approach to immune recognition. Moreover, the array was constructed with the aim of detecting surface proteins or exported proteins whilst our list does not require these characteristics. Only 30 of our top 100 novel gametocyte antigens were thus printed on this array. Antibody responses to five gametocyte proteins were significantly more prevalent in gametocyte carriers than in carriers of the asexual blood stage only. This is the first evidence that antibody responses may be indicative of current gametocyte carriage. Importantly, the dichotomization of gametocyte-exposed and non-exposed individuals was based on a single time-point screening for gametocytes by microscopy. Microscopy has a low sensitivity for detecting gametocytes that commonly circulate at low densities 53 and several of the asexual parasite carriers are likely to have had preceding or concurrent low densities of circulating gametocytes. Antibody prevalence in the group classified as gametocyte-negative by microscopy may thus be associated with concurrent low-density gametoctytemia and/or long-lived antibody responses acquired following previous gametocyte exposure. The presently analysed samples thus do not allow any conclusions on a possible role of submicroscopic gametocyte densities in boosting or maintaining antibody responses to gametocyte antigens. Refined studies with longitudinal sampling and gametocyte detection by sensitive qRT-PCR methodologies are needed to formally assess antibody kinetics in relation to gametocyte exposure and determine whether recent markers of exposure to blood stage antigens 50 can be complemented by a set of markers for recent or long-term gametocyte exposure.
We described the assembly of a curated gold standard set of gametocyte and asexual proteins and used this new resource to rank the likelihood of all Pf proteins and transcripts being specific to the gametocyte stage. Data from 18 publicly available studies were integrated to resolve partially conflicting evidence. The resulting consensus lists can be used for guidance of future investigations as we have shown the value of our predictions by in vitro validation.

Materials and Methods
Assembly of a gold standard for gametocyte and asexual proteins to weigh whole proteome/ transcriptome data sets. To build a gold standard against which the performance of individual data sets could be assessed, we identified proteins that are known to be expressed in either asexual parasites (mostly blood stage, also including sporozoites and liver stage) or gametocytes. This list was initially informed by literature review (Supplementary Table S1) for expression in the respective stages as detected by immunofluorescence assays and/or western blot, supplemented with P. falciparum blood stage or transmission blocking vaccine candidates. This initial list was then communicated with experts (including the authors DAB, PA, FS, CJJ, SMK, TWAK, MM, CD, RS and TB) and edited. If additional proteins were suggested for inclusion in the list, published evidence was requested and examined prior to inclusion of the proteinThe final asexual gold standard list contains 46 proteins; the final gametocyte list contains 41 proteins. These gold standard lists (Supplementary Table 1) represent the balance between very strict inclusion criteria and sufficient set size to evaluate the quality of all data sets integrated. We tested for the detection of these proteins or transcripts in the respective samples, using a Bayesian statistics approach that we have successfully applied previously for genes involved in anti-viral immune defense 31 .
Data selection and integration. Data sets that measured protein and transcript abundance in Pf gametocytes were balanced with data from other life stages and supplemented with studies of the rodent malaria parasite P. berghei. One MS study on P. vivax was included as it is based on of ex vivo blood material as opposed to all other studies that used in vitro cultivated parasites. Unique peptide counts were retrieved from plasmoDB (version 28) in which sequenced peptides from the published studies are always mapped to the most recent genome annotation, or supplementary material of the respective studies. Many aspects determine how well a protein is represented in proteomics data that are obtained via MS, like its length or posttranslational modifications. Remapping original MS data to newly annotated genomes improves the quality of the predicted proteins 19 . We were however not able to retrieve those data from the studies 13,14,16 and therefore decided to take those proteins at face value. Notice that also these early studies contribute significantly to our integrated lists.
Expression percentiles were retrieved from plasmoDB (version 28) or calculated from raw data in the respective supplementary material. Gametocyte samples were summarized if applicable (using the maximum peptide count/expression percentile of different stages or male and female gametocytes) as were asexual samples, only considering the highest expression in any sample or time point.
For MS and transcriptomics data sets, separate scores for gametocyte-specificity of any Pf gene have been calculated. In brief, protein or transcript expression has been categorized from absent to high expression levels as given by number of unique peptides or expression percentiles, respectively. For each of the respective bins, a score was calculated depending on the relative retrieval of gametocyte and asexual gold standard genes. The log ratio of these retrieved genes defined the score for all other genes within the same bin. The final gametocyte score calculates as the prior probability of a gene being gametocyte-specific that is updated using the contributions of the data sets: where gametocyteGS and asexualGS are the fractions of retrieved gametocyte and asexual gold standard genes in sample i, respectively. We used a pseudocount of 1 if necessary to prevent division by zero if none of the gold standard genes was retrieved in this specific sample and bin. It was assumed that the likelihood of a gene to be either gametocyte or asexual specific is equally high, thus the (log-transformed) prior equals 0 and the final score depends solely on the integrated data. In the selection of a set of proteins that we assigned to be gametocyte-specific we chose a cutoff score of 5.0 (proteomics-derived). The cutoff score of 5.0 can be interpreted as: a gene has to be 2 5 = 32 times more likely to be gametocyte-specific than asexual specific. The score of 5.0 was based on the behavior of the gold standard genes. Out of the 41 gametocyte gold standard genes, 37 have a score higher than 5.0, while none of the asexual gold standard genes do. When applicable, genes from Pb and Pv were treated as their respective Pf orthologs as retrieved from plas-moDB 54 , to be able to integrate all data sets. When no ortholog is known, the respective non-Pf data sets did not contribute to the score of this particular gene. Scores using Pf data exclusively were also calculated (Supplementary Table S2 includes all scores and rankings with expression information from all integrated studies).

Cross-validation of the scoring method.
We performed a ten-fold cross-validation to assess the predictive performance of the integrated gametocyte-specificity score (i.e. its ability to discriminate known gametocyte vs. asexual genes). For that, we subsampled both gold standard gene sets ten times (folds), without replacement (i.e. each gene is selected exactly once). Then for each fold we re-weighed and integrated the data sets based on nine-tenth of gold standard genes, and collected the ranks of the one-tenth of genes that were left out in that particular fold. A ROC curve was constructed based on those ranks. Using the same strategy, ROC curves for individual data sets that comprised both gametocyte and asexual samples were constructed for comparison.
Protein microarray to measure humoral immune responses. A protein microarray that was enriched for gametocyte proteins was produced and probed as described earlier for a study aiming to unravel the immune signature of naturally acquired transmission-reducing immune responses in gametocyte carriers (Stone, Campo et al. accepted manuscript 36 ). As a control group, Gambian asexual parasite carriers without gametocytes detectable by microscopy were included in the probing. For the current study array data from this control group (n = 63) and Gambian gametocyte carriers were used (n = 164) [55][56][57][58][59][60] . All of these 227 individuals were sampled during a period of intense malaria transmission intensity in The Gambia and likely had (multiple) previous malaria infections 55,61 . For these populations, responses to 30 newly defined highly gametocyte-specific antigens (from 24 genes) were compared between gametocyte carriers and non-carriers (Mann-Whitney U test). Seropositivity for each of the antigens was determined using a mixture model-based cutoff and related to gametocyte carriage using Fisher's exact test, corrected for multiple testing (Benjamini-Hochberg) of a total of 70 antigens (including 40 antigens from the gametocyte gold standard).
Transcript abundance in different life stages and strains. The abundance of 15 predicted gametocytespecific targets was measured in asexual parasites and gametocytes of four different Pf strains from in vitro culture. The targets were selected from the 100 highest scoring transcripts to account for uncertainties about the absolute scoring of transcriptomics data with a protein-based gold standard. We do not assume a clear hierarchy between these top 100 scoring transcripts and consider any of these genes highly gametocyte-specific. The 15 highest-ranking non-gold standard genes were selected based on a preliminary analysis of the data, and contain genes that are currently not annotated as well as genes with known protein function in gametocytes (PUF1 62 and Ccp4 63 ). In the final generation of the gametocyte-scores, all validation genes were retained in the top 100 scoring genes. The Pf strains used are of West African (NF54, NF166, NF175) and Southeast Asian origin (NF135). All strains were cultured and synchronized as described previously 20 . Using established standard curves, the same concentrations of parasites were compared for Ct values in qRT-PCR (for primers, see Supplementary Table S7). Extracted nucleic acids were DNase-treated before reverse transcription when introns were absent from the targets. Initial comparison was between mixed asexual blood stage parasites (considering the lowest Ct measured in any strain and replicate) and stage V gametocyte (highest Ct measured per strain and replicate). Promising targets with a high Ct Asexuals − Ct Gametocytes were further examined in serial dilutions of stage V gametocytes and synchronized asexual material of the strain NF54 (10, 20, 30, 40 hours post invasion, resembling early rings, late rings, trophozoites and schizonts, respectively). All qRT-PCR reactions were analyzed in technical triplicates, from biological triplicates (NF54) or duplicates (remaining strains).
RNA samples were used from a clinical malaria trial conducted in Western Kenya 39 . Samples from days 3 and 7 after treatment were selected to ensure a range of (low-density) gametocyte carriage to test qRT-PCR sensitivity.

Go term enrichment in top 100 proteins or rank-based enrichment. Current GO term annotation for
Pf was retrieved from plasmoDB (release 30) and analyzed using the topGO R package 64 for the enrichment of terms in the 100 highest scoring proteins versus all Pf proteins. Semantic clustering of the significant GO terms in Biological Process ontology was done with the Revigo webtool for Pf 65  Pf proteins and their ranks and scores was performed using the software available at http://software.broadinstitute.org/ gsea/downloads.jsp. Based on the cilium genes reported by the Syscilia consortium 66 , we assembled a "cilium" GO term of mixed ontology (GO:9999999) for Pf with 28 predicted orthologs (Supplementary Table S6). Data availability. All data generated or analysed during this study are included in this published article (and its Supplementary Information files).
• The gold standard lists (Supplementary Table S1) • Bayesian gametocyte scoring for proteomics and transcriptomics data. Includes Pf only-scores and expression values for any gene and individual data set (Supplementary Table S2) • Potential translationally repressed genes with high transcriptomics score (>7) and low proteomics score (<−10, Supplementary Table S3) • Overview over previously reported gametocyte-specificity per study (Supplementary Table S4) • Seroprevalence in The Gambia in gametocyte carriers and non-carriers (Supplementary Table S5) • Function predictions for highly gametocyte-specific proteins with lacking annotation (Supplementary  Table S6) • Transcript validation in 15 targets, including primer sequences for qRT-PCR (Supplementary Table S7) • GO term analyses with cilium genes (Supplementary Table S8).