Uncovering protein glycosylation dynamics and heterogeneity using deep quantitative glycoprofiling (DQGlyco)

Abstract

Protein glycosylation regulates essential cellular processes such as signaling, adhesion and cell-cell interactions; however, dysregulated glycosylation is associated with diseases such as cancer. Here we introduce deep quantitative glycoprofiling (DQGlyco), a robust method that integrates high-throughput sample preparation, highly sensitive detection and precise multiplexed quantification to investigate protein glycosylation dynamics at an unprecedented depth. Using DQGlyco, we profiled the mouse brain glycoproteome, identifying 177,198 unique N-glycopeptides-25 times more than previous studies. We quantified glycopeptide changes in human cells treated with a fucosylation inhibitor and characterized surface-exposed glycoforms. Furthermore, we analyzed tissue-specific glycosylation patterns in mice and demonstrated that a defined gut microbiota substantially remodels the mouse brain glycoproteome, shedding light on the link between the gut microbiome and brain protein functions. Additionally, we developed a novel strategy to evaluate glycoform solubility, offering new insights into their biophysical properties. Overall, the in-depth profiling offered by DQGlyco uncovered extensive complexity in glycosylation regulation.

Fig. 1. Deep profiling of protein glycosylation in human cell lines and mouse brain.

a, DQGlyco workflow. b, Average numbers of unique glycopeptides, glycosites and glycoproteins identified from LC–MS/MS triplicates of unfractionated samples from human cell lines (HelaK and HEK293T cells) and mouse brain. Data are presented as the mean values ± s.d. c, Glycopeptide enrichment specificity, calculated as the ratio of glycopeptide spectrum matches to all peptide spectrum matches (PSM) for unfractionated samples from human cell lines (HeLaK and HEK293T; n = 3 technical replicates per cell line) and mouse brain (n = 3 technical replicates), as well as from PGC-fractionated mouse brain samples. Data are presented as the mean values ± s.d. d, Total number of identified N-glycopeptides (unique sequence and glycan composition), N-glycosites and N-glycoproteins from PGC-fractionated mouse brain samples. e, UniProt and PhosphoSitePlus annotation of N-glycosites identified in mouse brain samples. f, Total number of unique glycopeptides of the PGC-fractionated mouse brain sample per binned m/z range compared to recent large-scale glycoproteomics studies (Riley et al. and Liu et al.). g, Total number of O-glycopeptides and O-glycoproteins identified in the PGC-fractionated mouse brain sample. h, GO enrichment analysis (STRINGdb) for proteins found to be only N-glycosylated, only O-glycosylated or both N-glycosated and O-glycosylated for two selected terms (intracellular organelle and plasma membrane).

Fig. 2. Protein glycosylation composition and microheterogeneity in the mouse brain.

a, Frequency of the number of unique glycan compositions per site in the PGC-fractionated mouse brain N-glycosylation data (on average, 18 glycan compositions were identified per site). b, Total number of unique glycopeptides, glycosites and glycoproteins per glycan class for the PGC-fractionated mouse brain sample. c, Frequencies of glycan classes identified in this study and in a recent glycoproteomics study using lectin-based enrichment. d, Distinct glycan classes are significantly overrepresented or underrepresented on specific functional protein domains from the InterPro database (log₂ odds ratio > 0 or log₂ odds ratio < 0, Fisher’s exact test, adjusted P value < 1 × 10⁻³; examples shown). e, Glycosylation profiles of sites belonging to the same protein are significantly more correlated compared to when glycan compositions are shuffled across the whole glycoproteome (Kendall rank correlation coefficient for n = 29,932 site pairs per distribution, two-sided Wilcoxon rank-sum test, P < 2.2 × 10⁻¹⁶). Box plots indicate the median and the first and third quartiles. Whiskers extend from the hinges to the largest value no further than 1.5 × the interquartile range. Data points beyond the end of the whiskers are plotted individually. f, Results of the GO enrichment analysis (STRINGdb) of proteins with sites displaying low, medium or high microheterogeneity (n = 1–2, 3–11 and 11+ unique glycan compositions per site, respectively) for selected significant terms (two-sided Fisher’s exact test, FDR < 0.05). g, Frequency of glycosites or glycopeptides with a close known phosphosite (±5 residues) per glycan class for O-glycosylation and N-glycosylation data.

Fig. 3. N-glycosylation in a structural context.

a, Glycosylated asparagines highlighted as spheres in the AlphaFold-predicted structure of the coxsackievirus and adenovirus receptor (P78310). Colors indicate the model predicted local distance difference test (pLDDT). b, Total number of identified glycans and relative glycan type composition per site in the same protein. Topological domains and IDRs are also highlighted (IDR definition in Methods). c, Validation AUROC for the predictive model of lowly and highly glycosylated asparagines that uses structural features as input (details in Methods; low-N = 2,203 and high-N = 747 in mouse; low-N = 1,196 and high-N = 488 in human). d, Feature importance plot. The x axis indicates the F score for every feature on the y axis. e, pPSE (24 Å, 180°) for highly and lowly glycosylated asparagines (low-N = 2,203 and high-N = 747 in mouse; low-N = 1,196 and high-N = 488 in human). f, Log-transformed number of glycoforms on asparagines within extracellular (N = 2,813) or cytoplasmic (N = 335) topological domains. g, Average AlphaMissense score for lowly (N = 1,196) and highly (N = 488) glycosylated asparagines on data from human cell lines. All box plots indicate the median and the first and third quartiles. Whiskers extend from the hinges to the largest value no further than 1.5 × the interquartile range. Data points beyond the end of the whiskers are plotted individually. All P values indicate the comparison of distributions using a two-sided Wilcoxon rank-sum test.

Fig. 4. Time course of fucosylation inhibition.

a, Experimental design. b, PCA of the normalized TMT glycopeptide intensities per treatment condition, time point and replicate. PC, principal component. c, Glycopeptide modulation over time. Glycopeptides significantly regulated for at least one time point were grouped in four clusters depending on their kinetics of regulation using neural gas clustering. d, Proportion of glycan classes present in the four clusters. e, Fucosylated peptides being significantly downregulated for at least one time point were grouped into three clusters using neural gas clustering. f, Glycoproteins having fucosylated glycoforms belonging to the three clusters. Here, the mean log₂ fold changes of the fucosylated glycopeptides belonging to each cluster are displayed. FC, fold change. g,h, The fucosylated glycoforms of the C-type mannose receptor 2 (MRC2) and the neuronal cell adhesion molecule (NRCAM) exhibit different rates of downregulation.

Fig. 5. Systematic characterization of surface-exposed glycoforms.

a, Intact living HEK293 cells are treated with either PNGase or proteinase K. Surface-exposed glycoforms are affected while intracellular glycoforms remain protected. b, Number of glycoforms on affected sites (at least two glycoforms changing) significantly changing in abundance upon enzymatic treatments, per subcellular compartment annotation of proteins. c,d, Characterization of surface-exposed glycoforms belonging to the zinc ion channel ZIP6 and to the cell adhesion protein DSC2. Each data point represents a unique glycoform, with the color representing the glycan class. The shape of the data point reflects significance, while the dotted line represents the effect size cutoff (log₂ 0.6 = 1.5 fold change). e, Fold change of glycoforms, separated per glycan class, on sites for which at least two glycoforms were affected by the enzymatic treatments. High-mannose glycans were significantly less exposed than the other glycan classes (two-sided t-test). The horizontal lines represent the median. f, Comparison of fold changes of glycopeptide abundance upon proteinase K (y axis) and PNGase treatments (x axis), for sites affected by both enzymes (sites with at least two glycopeptides significantly changing in abundance upon treatment) and glycopeptides identified in both experiments. g, Number of times a given glycan composition was identified on a glycopeptide significantly changing in abundance upon proteinase K (y axis) and PNGase (x axis) treatments. The size of the dot represents the frequency with which this composition was identified as surface exposed.

Fig. 6. Quantitative profiling of the glycoproteome across tissues.

a, Experimental design. Glycoproteomes of the brain, liver and kidney tissues belonging to two male mice colonized (each in technical triplicates) with a human gut microbiome were profiled. b, PCA of the glycopeptide intensities after protein abundance normalization. c, Spearman correlation of glycopeptide intensities among different samples after protein abundance normalization. d, Pairwise Pearson correlation of the relative intensities of all glycopeptides on a given glycosite between tissues (for glycosites with at least five glycopeptides quantified) identifying glycosites with low tissue specificity (high correlation) and sites with high tissue specificity (low correlation). e, Density plot of site correlation values across tissues. The mean (0.669) is represented as a dotted line and chosen as a cutoff value for glycosite tissue specificity. f, Density plot of site correlation values for the same tissue across two mice. g, Fraction of tissue-specific sites (correlation < 0.669) within a given glycoprotein (for proteins with at least three glycosites with a correlation value). h, Correlation values across tissues for the three sites quantified on the insulin receptor protein. i, Functional domain enrichment depending on glycosite tissue specificity (InterPro database, log₂ odds ratio > 0 or log₂ odds ratio < 0, two-sided Fisher’s exact test, adjusted P value < 0.001; examples shown). j, Glycosite correlation values across tissues when only one type of glycan class is considered (two-sided, Wilcoxon rank-sum test). The horizontal lines represent the median. k, Glycosite correlation values across tissues per subcellular compartments. The horizontal lines represent the median.

Fig. 7. Systematic characterization of glycoform biophysical properties.

a, Solubility proteome profiling experimental workflow, enabling the proteome-wide characterization of the in vivo biophysical properties of glycoforms by assessing their solubility in a nondenaturing detergent (NP40) and in a SDS. b, Glycoform solubility compared to the solubility of the general protein population per glycan class (two-sided t-test). The horizontal lines represent the median. c–e, Solubility profiling of glycoforms of proteins encoded by Cdh13 (c), Thy1 (d) and Cd47 (e). Each dot represents a glycoform while the solid lines represent the solubility of the general protein population, based on nonmodified peptides.

Fig. 8. Impact of the gut microbiome on the mouse brain glycoproteome.

a, Three groups of six adult germ-free mice (three males and three females) were colonized by B. uniformis, by a defined eight-member community of human gut commensal or remained germ-free for 2 weeks. b, Spearman correlation of glycopeptide intensities between biological replicates after protein abundance normalization. c, PCA of normalized reporter intensities of the 18 brain samples. d, N-glycopeptide regulation in the mouse brain depending on gut microbiome composition. e, Overlap of significantly regulated proteins at the glycoproteome and proteome level. f, Mitochondrial proteins involved in cellular respiration and G protein γ-subunits regulated in protein abundance. g, Cytoskeleton proteins change in thermal stability upon gut microbiome colonization (n, number of proteins in each family; actins: Actb, Acta2 and Actg1; dyneins: Dnm1, Dnm2 and Dnm3; myosins: Myh10, Myh11, Myh9 and Myl12b; α-tubulins: Tuba1b, Tuba4a, Tuba8 and Tubal3; β-tubulins: Tubb1, Tubb2, Tubb2b, Tubb3, Tubb4a, Tubb4b, Tubb5 and Tubb6; γ-tubulins: Tubg1 and Tubg2). h, ΔAUC of thermal stability of proteins in brains of germ-free and eight-member colonized mice. Glycoproteins with at least one glycopeptide changing significantly in abundance upon colonization and proteins with GO term ‘structural molecule activity’ had higher thermal stability compared to other glycoproteins (two-sided t-test) The horizontal lines represent the median. i, N-glycosylation modulation on proteins involved in neurotransmission. j, Cntnap1 changes in thermal stability upon gut microbiome colonization. k, Site-specific regulation of N-glycosylation on the Grin2a glutamate receptor. Each vertical line represents one glycopeptide and each dot represents one replicate (3 conditions × 6 biological replicates). l, Eaat2 encodes a membrane-bound transporter, glycosylated at the N205 and N215 sites on an extracellular loop (blue). Human structure combined with AlphaFold predictions for the extracellular region is shown. The light-gray structure represents the homotrimer (Protein Data Bank 7VR8). Yellow planes depict the membrane orientation; dark-blue and dark-gray regions depict the predicted structures aligned to the template structure. Spheres indicate the N205 and N215 glycosites. m, The glycoforms on the two glycosylation sites of the Eaat2 protein are differentially modulated upon gut microbiota colonization. n, Solubility of the Eaat2 glycoforms from the solubility experiment. Gut-microbiome-modulated glycoforms are highlighted. The horizontal lines represent the median.