The N-glycan profile at diagnosis was evaluated to explore the prognostic potential of total plasma N-glycome as a biomarker for COVID-19 disease severity. Patients were classified at diagnosis as mild (n=56), severe (n=105) and critical (n=35) (24). Total plasma N-glycan composition was determined by LC-MS/MS analysis of RFMS-labelled glycans, confirming structures with MS/MS data. A total of 36 structures were identified in all patients ( Supplementary Table 1 ) and 13 additional structures were detected but could not be accurately annotated. The most abundant structures in all the samples were biantennary and were specifically Neu5Ac2Hex5HexNAc4 followed by Neu5Ac1Hex5HexNAc4 and Fuc1Hex3HexNAc4. A total of 22 N-glycans were found to be statistically significant between groups at diagnosis ( Figure 1 ). These 22 N-glycans were further tested as potential biomarkers for COVID-19 for their ability to predict the disease course at diagnosis.

A random forest analysis was used to evaluate the effectiveness of the statistically significant total plasma N-glycans as predictive biomarkers to indicate disease severity in COVID-19 at diagnosis. Results indicated that the top 3 N-glycans with the highest discriminatory power between groups were Fuc1Hex5HexNAc5 and Hex6HexNAc3, which decrease with increasing severity, and Neu5Ac1Hex4HexNAc3, which increases with severity ( Figure 2 ).

Subsequently, prediction models were built to evaluate the predictive efficacy of some of the N-glycans with the highest classification accuracy. Receiver operating characteristic (ROC) curves were generated through a binary logistic regression analysis comparing the mild vs. critical groups, severe vs. critical groups and mild vs. severe groups ( Table 1 ; Figure 3 ). The ROC curve analysis was performed for an individual N-glycan and for the ratio of two N-glycans. The ratio was computed as it was able to provide more information than the two corresponding N-glycans alone. These ratios were computed for all possible N-glycan pairs and the top-ranked ratio (based on p-values) was selected for further biomarker analysis.

The ratio of Fuc1Hex5HexNAc4/Hex10HexNAc2 was the top-ranked ratio (p<0.001), showing a high capacity to differentiate between mild and critical patients, with an AUC of 0.881 and a specificity and a sensitivity of 90.9 and 73.1%, respectively. Similarly, Fuc1Hex5HexNAc5, individually, was also optimal to distinguish between mild from critical patients (AUC = 0.88; p<0.001) ( Table 1 ; Figure 3A ). In the random forest analysis, Fuc1Hex5HexNAc5 was considered to be the N-glycan with the highest classification accuracy.

Regarding the other two pairwise comparisons, severe vs critical and mild vs severe, more biomarkers should be considered to rapidly identify the clinical progression of patients at an early stage, as AUCs ranged between 0.666 and 0.793. Next, a regression model was performed in a randomly selected set of patients to validate the accuracy of these aforementioned N-glycans to properly predict mild and critical outcomes at diagnosis. Results indicated that the ratio Fuc1Hex5HexNAc4/Hex10HexNAc2 (AUC=0.862, p<0.001, specificity=78.8% and sensitivity=88.0%) and Fuc1Hex5HexNAc5 (AUC=0.85, p<0.001, specificity=90.9% and sensitivity=72.0%) could be suitable biomarkers to distinguish mild from critical outcomes ( Figure 3B ).

The N-glycan profile at 4 weeks postdiagnosis was evaluated to explore the potential of total plasma N-glycome to act as a biomarker for COVID-19 disease course prediction. A total of 18 N-glycans were found to be statistically significant between groups at 4 weeks postdiagnosis, 15 of which were also significantly altered at diagnosis. The progression of these 15 N-glycans at diagnosis and at 4 weeks for the different groups of severity is shown in Figure 4 . The aforementioned N-glycans (Fuc1Hex5HexNAc5 and Hex10HexNAc2) that were tested as potential biomarkers using a regression model were altered at diagnosis and at 4 weeks postdiagnosis. Fuc1Hex5HexNAc5 showed the same tendency at both time points, decreasing with increasing severity. Differently, Hex10HexNAc2 increased with increasing severity at diagnosis and at 4 weeks postdiagnosis. Interestingly, for both N-glycans, the significance of the alterations was higher at diagnosis than at 4 weeks. The remaining 13 N-glycans altered at both time points showed similar tendencies at diagnosis and at 4 weeks postdiagnosis ( Figure 4 ). More specifically, all the statistically significant fucosylated structures decreased with increasing severity at diagnosis and at 4 weeks, except for Neu5Ac3Fuc1Hex6HexNAc5 which increased with increasing severity.

Three additional N-glycans were altered only at 4 weeks between severity groups, specifically Hex4HexNAc5 and NeuAc1Hex4HexNAc5, which increased with severity and Hex7HexNAc2, which decreased with severity. To obtain a global picture of the observed significant alterations at diagnosis and at 4 weeks postdiagnosis, glycan structures were grouped by their traits in terms of bisecting GlcNAc (N-acetylglucosamine), oligomannose-type, fucosylation, galactosylation and sialylation ( Figure 5 ). Total plasma bisecting GlcNAc structures significantly decreased in critical patients when compared with mild and severe ones at diagnosis, whereas the trend shifted at 4 weeks postdiagnosis, increasing in severe patients compared to mild ones. A similar pattern was observed for fucosylation and galactosylation, both decreasing with increasing severity at diagnosis and at 4 weeks postdiagnosis. Notably, a greater significance between groups of severity was observed at diagnosis than at 4 weeks. Differently, oligomannose structures increased with severity at diagnosis and showed no significant alterations between groups at 4 weeks postdiagnosis. No significant changes were observed in sialylation.

A total of 196 patients with SARS-CoV-2 infection who had positive polymerase chain reaction (PCR) confirmation within the first 21 days of infection make up the COVID-19 patient cohort. Patients were categorised into 3 groups of severity (mild (n=56), severe (n=105) and critical (n=35)), according to the inclusion criteria described in “Diagnosis and Treatment Protocol for COVID-19 Patients (version 8 trial)” (24). Additionally, disease progression was evaluated at 4 weeks postdiagnosis for 122 of the aforementioned patients (mild (n=22), severe (n=82) and critical (n=18)). The severity classification changed for 49 patients during the 4 weeks postdiagnosis, with 16 patients changing from mild to severe, 32 from severe to critical and one from critical to severe. Data from patients were stored in a database with information about the hospitalization, including the symptoms that were present at the time of admission, radiological findings, the severity of pneumonia, the need for oxygen therapy, the medical treatment received, as well as demographic information and previously diagnosed diseases of interest (9). Patients’ demographic and relevant clinical data are summarised in Supplementary Table 2 . Briefly, the cohort median age increased with COVID-19 severity and the female sex predominated in the mild and critical groups. Hypertension was the most common comorbidity in all groups and more predominant in severe and critical patients. Drug administration and the need for oxygen and aggressive treatments were consistent with the severity of the disease. A total of 25 patients from the cohort died from COVID-19. The serum biochemical composition was also characterised at the time of admission in the entire cohort, including both routine and inflammatory parameters. The most significant differences in the blood pattern were observed between mild and critical patients (9).

The sampling protocol consisted of a clinical evaluation, blood cell count, and standard biochemical parameters at inclusion (baseline). Serum samples were stored at -80°C at BioBank - Institut d’Investigació Sanitària Pere Virgili (IISPV) facilities until analysed (9).

Protocols were carried out in accordance with the recommendations of the Ethical and Scientific Committees from each participating institution and were approved by the Committee for Ethical Clinical Research following the rules of Good Clinical Practice from the IISPV (079/2020, CEIm IISPV). The CEIm IISPV is an independent committee, which oversees the correct adherence to the ethical standards governing clinical trials and research projects that are carried out in our environment, specifically in terms of its methodology, ethics and laws. Its members include both health and non-health professionals. All subjects or their relatives gave written informed consent in accordance with the Declaration of Helsinki (9).

Ammonium bicarbonate (ABC), ammonium formate and formic acid (LC-MS grade) were purchased from Sigma-Aldrich (St. Louis, MO, USA). Acetonitrile (LC-MS grade) was purchased from Merck (Darmstadt, Germany). The water used throughout the study was purified with a Milli-Q system from Millipore (Burlington, MA, USA).

The release and labelling of N-glycans were done using a method previously described (32). Briefly, total protein was quantified with the Bradford assay. Sample denaturation, de-N-glycosylation, labelling with Rapifluor-MS (RFMS), and purification were performed in accordance with the Waters Corporation “GlycoWorks RapiFluor-MS N-Glycan Kit Care and Use Manual” (p/n 715004793). Briefly, 15 µg of protein were denatured at 90°C for 3 min in the presence of 5% (w/v) RapiGest. De-N-glycosylation was then conducted at 50°C for 5 min, by adding 1.2 μL of Rapid PNGase F. Subsequently, the digested samples were directly labelled with 6 μL of RFMS reagent without the need for purification. The reaction proceeded at room temperature for 5 min before using a GlycoWorks μElution Plate for the SPE Clean-up procedure. Glycans were eluted with 200 mM ammonium acetate in 5% acetonitrile.

The LC-MS/MS analysis was carried out using a previously described analytical method (32). Derivatised samples were analysed using an Agilent UHPLC 1290 Infinity Series coupled to an Agilent qTOF/MS 6550 Series (Agilent Technologies, Santa Clara, CA). N-glycans were separated on a Waters ACQUITY UPLC BEH amide column (2.1 mm × 150 mm i.d., particle size 1.7 μm), using 50 mM ammonium formate solution (mobile phase A) and 100% acetonitrile (mobile phase B). The analysis was performed at a flow rate of 0.4 mL/min, with an injection volume of 20 µL, and the column temperature set at 60°C. The chromatographic gradient began by ramping mobile phase A from 25% to 46% over a period of 35 minutes. From 35 to 36.5 min, the gradient ramped from 46 to 100% solvent A and the flow rate was lowered to 0.2 mL/min. 100% solvent A was held constant from 36.5 to 39.5 min, after which the percentage of solvent A decreased to 25%, from 39.5 min to 43.1 min. The flow rate was then increased back to 0.4 mL/min from 43.1 to 47.6 min and solvent A was held constant at 25%. Lastly, the parameters were held constant from 47.6 to 55.0 min. The qTOF operated in positive electrospray ionisation mode (ESI+), and mass spectra were recorded between m/z 300–1700 at 1.5 spectra/s. The source conditions were set as follows: nebuliser gas at 25 psi, gas temperature at 200°C, gas flow at 12 L/min, sheath gas temperature at 250°C, sheath gas flow at 12 L/min, capillary voltage at 3500 V, and nozzle voltage at 500 V. Additionally, tandem mass experiments (MS/MS) using data dependence acquisition at a collision energy of 30 eV from the 10 most intense ions were used for identification purposes.

The MS data were first processed using Agilent MassHunter Qualitative Analysis B.07 software. For the identification of N-glycans, total ion chromatograms (TIC) containing MS/MS fragmentation data were deconvoluted using the “find molecular feature” algorithm, which detected chromatographic peaks considered to be N-glycans. The resulting list of entities containing the MS/MS data was exported and loaded to Simglycan software for molecular and structural elucidation. Simglycan is a high-throughput structural identification tool that uses a built-in database with theoretical fragmentation profiles to predict the structure of glycans (33). Additionally, the Simglycan results were matched against the GlycoStore database (https://www.glycostore.org) to refine the identification results. GlycoStore provided elution property information for over 850 unique structures including standardised retention times, expressed as glucose units (GU), for the RFMS-labelled glycans (34). Following the identification process, 36 N-glycans were identified in human plasma and 13 additional structures were detected but could not be annotated. Subsequently, the exact masses [M+H]⁺, [M+2H]²⁺ or [M+3H]³⁺ of these structures were extracted on all samples using Agilent Mass Hunter Quantitative software (B.07) to generate a refined matrix of quantitative data for statistical purposes. Finally, the obtained data matrix, which included the peak area for each identified N-glycan, was normalised by dividing the peak area of each N-glycan by the sum of the peak areas of all N-glycans within each respective sample.

Non-parametric Kruskal-Wallis test was used to assess significant differences between groups of severity, employing Mass Profiler Professional (MPP) software v.15.1 (Agilent Technologies, Massachusetts, USA). Graphical representations were generated with GraphPad Prism software (version 9.0, GraphPad Inc., San Diego, CA) and Metaboanalyst 5.0. Both SPSS (version 21.0, SPSS Inc., Chicago, IL) and Metaboanalyst were employed to perform random forest analyses, and create binary logistic regression models and ROC curves to evaluate the potential accuracy of the selected biomarkers for predicting COVID-19 severity. Kruskal-Wallis test and χ2 test were used to determine the demographic and clinical features differences between groups of severity. The Spearman’s correlation coefficient, used to measure the association between the study cohort variables and groups of severity, was performed with SPSS. Results were considered statistically significant at p < 0.05.