Metabolic dyshomeostasis induced by SARS-CoV-2 structural proteins reveals immunological insights into viral olfactory interactions

One of the most common symptoms in COVID-19 is a sudden loss of smell. SARS-CoV-2 has been detected in the olfactory bulb (OB) from animal models and sporadically in COVID-19 patients. To decipher the specific role over the SARS-CoV-2 proteome at olfactory level, we characterized the in-depth molecular imbalance induced by the expression of GFP-tagged SARS-CoV-2 structural proteins (M, N, E, S) on mouse OB cells. Transcriptomic and proteomic trajectories uncovered a widespread metabolic remodeling commonly converging in extracellular matrix organization, lipid metabolism and signaling by receptor tyrosine kinases. The molecular singularities and specific interactome expression modules were also characterized for each viral structural factor. The intracellular molecular imbalance induced by each SARS-CoV-2 structural protein was accompanied by differential activation dynamics in survival and immunological routes in parallel with a differentiated secretion profile of chemokines in OB cells. Machine learning through a proteotranscriptomic data integration uncovered TGF-beta signaling as a confluent activation node by the SARS-CoV-2 structural proteome. Taken together, these data provide important avenues for understanding the multifunctional immunomodulatory properties of SARS-CoV-2 M, N, S and E proteins beyond their intrinsic role in virion formation, deciphering mechanistic clues to the olfactory inflammation observed in COVID-19 patients.


Introduction
The coronavirus disease 2019 (COVID- 19) is an ongoing viral pandemic caused by the rapid spread of severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) (1). It is well known that infection is transmitted by respiratory droplets with an incubation period of approximately 5 days (2), however, there are still many gaps regarding the molecular mechanisms disrupted during the infection and the immune system response triggered against this pathogen. The approximately 30,000nucleotide SARS-CoV-2 RNA genome encodes for 14 major open-reading frames (ORFs). The 4 structural proteins (S (spike), membrane (M), envelope (E), and nucleocapside (N)) are expressed from subgenomic RNAs transcribed from the RNA genome. In addition, SARS-CoV-2 expresses another 14 non-structural proteins and 8 accessory proteins which participate in RNA replication-transcription, and pathogenesis (3). The clinical presentation of COVID-19 varies from asymptomatic, mild to moderate disease in the majority of cases (4). However, fatal consequence can occur in patients with comorbidities (5,6). SARS-CoV-2 infection in human patients include fever and respiratory symptoms, such as fever, dry cough, dyspnea, and pneumonia (4). In addition, accumulating evidence has demonstrated the appearance of other symptoms of clinical importance, such as neurological affectation and olfactory and gustatory dysfunctions, among others (7)(8)(9)(10)(11)(12). Importantly from a diagnostic point of view, these changes usually occur earlier than other symptoms, suggesting a potential specific predictive value (13,14).
Although the molecular mechanisms of olfactory and neurological disorders in SARS-CoV-2-infected patients are not fully understood, there is growing evidence that the olfactory loss observed in COVID-19 occurs through the infection of supporting cells in the olfactory system (15). The key pathological mechanism for SARS-CoV-2 to invade the host cell involves the S protein attachment to the human angiotensinconverting enzyme 2 (ACE2) cell surface protein. Then, the virus uptake is facilitated by the action of the TMPRSS2 protease (16). Importantly, sustentacular cells highly express both ACE2 and TMPRSS2 proteins, while little expression is found in olfactory neurons (15). Therefore, these sustentacular cells might be the main targets for SARS-CoV-2 infection, leaving the neurons vulnerable and depleted of nutrients. However, other studies which focused on the olfactory epithelium from SARS-CoV-2infected patients with acute loss of smell, point out that olfactory sensory neurons, supporting cells, and immune cells are major sites of SARS-CoV-2 infection (17). Certain coronaviruses have the ability to pass from the OE through the cribriform plate to infect the OB, a second stage of the olfactory pathway (18,19). In fact, several reports have detected evidences of SARS-CoV-2 at the OB (20). However, the neuroinvasive potential of SARS-CoV-2 and future neurological vulnerabilities remain controversial (21,22).
Omics technologies have allowed to increase our knowledge about the molecular imbalance triggered by the SARS-CoV-2 in multiple human organs as well as in-vitro systems (23-26). However, there are still major gaps in the pathways leading to loss of smell caused by SARS-CoV-2, in addition to the molecular consequences from SARS-CoV-2 structural proteins at the olfactory level. Bearing in mind that the olfactory loss is considered an early predictor of COVID-19, we consider that the methodologies simulating SARS-CoV-2 infection in olfactory cells in combination with -omics strategies is a straightforward approach to elucidate the host and viral mechanisms involved in SARS-CoV-2-induced anosmia. These approaches can also identify potential therapeutic targets for non-invasive intranasal treatments. In this study, we have characterized the molecular dyshomeostasis induced by each SARS-CoV-2 structural protein (S, N, M, E) in murine OBC1 cells as an olfactory in-vitro system. The effects induced by the expression of GFP-tagged SARS-CoV-2 structural proteins uncovered a tangled crosstalk between multiple biofunctions and pathways as well as an immunological response differentially modulated by each SARS-CoV-2 structural protein. Our results also identified novel viral-host protein interactions partially shared by the structural components of the SARS-CoV-2 proteome.

Materials
The following reagents and materials were used.

Generation of M, E, N, S and GFP fusion proteins
Coding sequences for M, E, N and S proteins were retrieved from the SARS-CoV-2 reference genome (https://www.ncbi. nlm.nih.gov/nuccore/1798174254). Coding sequences were further codon-optimized for their expression in human cells, and recombinant synthetic genes fused to the GFP-coding sequence at the 3-termini were ordered from GeneArt (ThermoFisher). Genes were further cloned into pDUAL-PuroR lentivectors (27,28) under the transcriptional control of the SFFV promoter. pDUAL-PuroR lentivectors express puromycin resistance under the transcriptional control of the human ubiquitin promoter.

Cell culture and generation of transfected cell lines
Immortalized murine olfactory bulb OBC1 cells (ABM. T0235) were cultured in DMEM (Gibco) supplemented with 10% FBS (Merck Millipore) and 1% penicillin/streptomycin (ABM) and grown in a 5% CO 2 humidified atmosphere at 37°C . For transfection, OBC1 cells were seeded in 6-well plates at a confluency of 80%. Twenty microliters of Fugene HD (Promega) were incubated separately in 50 ul of DMEM without FBS and then mixed with three micrograms of M-GFP, E-GFP, N-GFP, S-GFP and control GFP plasmid DNA, followed by 15 min incubation at room temperature. The cells were transfected in a 3ml total volume for 72 hours. After three days, the culture medium was replaced with DMEM supplemented with 1ug/ml puromycin (Thermo). Transfection efficiency was assessed using flow cytometry and to achieve a more specific population, cell sorting based on GFP expression was performed among all experimental conditions.

RNA sequencing (RNA-seq) and data analysis
Briefly, total RNA was extracted and purified using a RNeasy Mini Kit (Qiagen, Hilden, Germany) following manufacturer's instructions. Sequencing libraries were prepared by following the Illumina Stranded Total RNA Prep with Ribo-Zero Plus (Illumina Inc., San Diego, CA) from 100 ng of total RNA, that has been depleted by following the instructions. All libraries were run in a HiSeq1500 PE100 lane in Rapid mode, pooled in equimolar amounts to a 10nm final concentration. The library concentration was measured by Qubit 3.0 (Invitrogen) and library size ensured by capilar electrophoresis in Fragment Analyzer (AATI). The quality of the RNAseq results was initially assessed using FastQC v0.11.9 (http://www.bioinformatics. babraham.ac.uk/projects/fastqc/) and MultiQC v1.9 (http:// multiqc.info/). The raw reads were trimmed, filtered for those with a Phred quality score of at least 25 and all adapters were removed with TrimGalore v0.5.0 (https://www.bioinformatics. brabraham.ac.uk/projects/trim_galore/). Trimmed reads were analyzed with SortMeRNA v2.1 software (https://bioinfo.lifl.fr/ RNA/sortmerna/) (29) to delete the 18S and 28S rRNA to eliminate the rRNA residues that could remain undepleted by the chemical treatment in the library preparation. Clean reads were aligned versus the Mus Musculus reference genome (release GRCm38.p6/GCA_000001635.8, ftp://ftp.ensembl.org) using HISAT2 v2.2.1 (https://daehwankimlab.github.io/hisat2/) (30) with default parameters. Resulting alignment files were quality assessed with Qualimap2 (http://qualimap.bioinfo.cipf.es) (31) and sorted and indexed with Samtools software (32). After taking a read count on gene features with the FeatureCounts tool (http:// subread.sourceforge.net) (33), quantitative differential expression analysis between conditions was performed by DESeq2 (34), implemented as R Bioconductor package, performing read-count normalization by following a negative binomial distribution model. In order to automate this process and facilitate all group combination analysis, the SARTools pipeline (35) was used. All resultant data was obtained as HTML files and CSV tables, including density count distribution analysis, pairwise scatter plots, cluster dendrograms, Principal Component Analysis (PCoA) plots, size factor estimations, dispersion plots and MA and Volcano plots. The resulting CSV file, including raw counts, normalized counts, Fold-Change estimation and dispersion data was annotated with additional data from the Biomart database (https://www.ensembl.org/ biomart/martview/346d6d487e88676fd509a1b9a642edb2). In order to control the False Discovery Rate (FDR), the p-values were amended by Benjamini-Hochberg (BH) multiple testing corrections. Those features showing corrected p-values below the 0.05 threshold and FoldChange values >1.5 or <0.5 were considered up-or down-regulated genes, respectively.

Proteomics and data analysis 2.5.1 Protein extraction
Culture medium was removed, and cells were washed with 1X cold PBS. Then, pellet was resuspended in lysis buffer containing 7M urea, 2M tiourea and 50Mm DTT and incubated in ice for 30 min, vortexing each 10 min. After a sonication step, the lysate was centrifuged for 20 minutes at 20000xg at 15°C. The supernantant was then added to a new Eppendorf and protein concentration was assessed using the Bradford assay (Biorad).

SWATH-mass spectrometry proteomics: MS/MS library generation
As an input for generating the SWATH-MS assay library, a pool of 15 samples (2.5µg/sample) derived each cell replicate was used. Protein extracts (30 µg) were diluted in Laemmli sample buffer and loaded into a 1.5 mm thick polyacrylamide gel with a 4% stacking gel casted over a 12.5% resolving gel. Total gel was stained with Coomassie Brilliant Blue and 12 equals slides from the pooled sample was excised from the gel and transferred into 1.5 mL Eppendorf tubes. Protein enzymatic cleavage was carried out with trypsin (Promega; 1:20, w/w) at 37°C for 16 h. Peptide mixture was dried in a speed vacuum for 20 min. Purification and concentration of peptides was performed using C18 Zip Tip Solid Phase Extraction (Millipore). Peptides recovered from ingel digestion processing were reconstituted into a final concentration of 0.5µg/µL of 2% ACN, 0.5% FA, 97.5% MilliQ-water prior to mass spectrometric analysis. MS/MS datasets for spectral library generation were acquired on a Triple TOF 5600+ mass spectrometer (Sciex, Canada) interfaced to an Eksigent nanoLC ultra 2D pump system (SCIEX, Canada) fitted with a 75 mm ID column (Thermo Scientific 0.075 × 250mm, particle size 3 mm and pore size 100 Å). Prior to separation, the peptides were concentrated on a C18 precolumn (Thermo Scientific 0.1 × 50mm, particle size 5 mm and pore size 100 Å). Mobile phases were 100% water 0.1% formic acid (FA) (buffer A) and 100% Acetonitrile 0.1% FA (buffer B). Column gradient was developed in a gradient from 2% B to 40% B in 120 min. Column was equilibrated in 95% B for 10 min and 2% B for 10 min. During all process, precolumn was in line with column and flow maintained all along the gradient at 300 nl/min. Output of the separation column was directly coupled to nano-electrospray source. MS1 spectra was collected in the range of 350-1250 m/z for 250 ms. The 35 most intense precursors with charge states of 2 to 5 that exceeded 150 counts per second were selected for fragmentation, rolling collision energy was used for fragmentation and MS2 spectra were collected in the range of 230-1500 m/z for 100 ms. The precursor ions were dynamically excluded from reselection for 15 s. MS/MS data acquisition was performed using AnalystTF 1.7 (Sciex) and spectra files were processed through ProteinPilot v5.0 search engine (Sciex) using Paragon ™ Algo-rithm (v.4.0.0.0) for database search. To avoid using the same spectral evidence in more than one protein, the identified proteins were grouped based on MS/MS spectra by the Progroup ™ algorithm, regardless of the peptide sequence assigned. The protein within each group that could explain more spectral data with confidence was depicted as the primary protein of the group. FDR was performed using a non-lineal fitting method (36) and displayed results were those reporting a 1% Global FDR or better.

SWATH-mass spectrometry proteomics: quantitative analysis
Protein extracts (20 µg) from each sample were reduced by addition of DTT to a final concentration of 10mM and incubation at room temperature for 30 minutes. Subsequent alkylation by 30 mM iodoacetamide was performed for 30minutes in the dark. An additional reduction step was performed by 30mM DTT, allowing the reaction to stand at room temperature for 30 min. The mixture was diluted to 0.6M urea using MilliQ-water, and after trypsin addition (Promega) (enzyme:protein, 1:50 w/w), the sample was incubated at 37°C for 16h. Digestion was quenched by acidification with acetic acid. The digestion mixture was dried in a SpeedVac.
Purification and concentration of peptides was performed using C18 Zip Tip Solid Phase Extraction (Millipore). The peptides recovered were reconstituted into a final concentration of 0.5µg/µL of 2% ACN, 0.5% FA, 97.5% MilliQ-water prior to mass spectrometric analysis. For SWATH-MS-based experiments the instrument (Sciex TripleTOF 5600+) was configured as described by Gillet et al. (37). Briefly, the mass spectrometer was operated in a looped product ion mode. In this mode, the instrument was specifically tuned to allow a quadrupole resolution of Da/mass selection. The stability of the mass selection was maintained by the operation of the Radio Frequency (RF) and Direct Current (DC) voltages on the isolation quadrupole in an independent manner. Using an isolation width of 16 Da (15 Da of optimal ion transmission efficiency and 1 Da for the window overlap), a set of 37 overlapping windows were constructed covering the mass range 450-1000 Da. In this way, 1 mL of each sample was loaded onto a trap column (Thermo Scientific 0.1 × 50mm, particle size 5 mm and pore size 100 Å) and desalted with 0.1% TFA at 3 mL/ min during 10 min. The peptides were loaded onto an analytical column (Thermo Scientific 0.075 × 250mm, particle size 3 mm and pore size 100 Å) equilibrated in 2% acetonitrile 0.1% FA. Peptide elution was carried out with a linear gradient of 2 to 40% B in 120 min (mobile phases A:100% water 0.1% formic acid (FA) and B: 100% Acetonitrile 0.1% FA) at a flow rate of 300 nL/min. Eluted peptides were infused in the mass spectrometer. The Triple-TOF was operated in swath mode, in which a 0.050 s TOF MS scan from 350 to 1250 m/z was performed, followed by 0.080 s product ion scans from 230 to 1800 m/z on the 37 defined windows (3.05 s/cycle). Collision energy was set to optimum energy for a 2 + ion at the center of each SWATH block with a 15 eV collision energy spread. The mass spectrometer was always operated in high sensitivity mode. The resulting ProteinPilot group file from library generation was loaded into PeakView ® (v2.1, Sciex) and peaks from SWATH runs were extracted with a peptide confidence threshold of 99% confidence (Unused Score ≥ 1.3) and a FDR lower than 1%. For this, the MS/MS spectra of the assigned peptides was extracted by ProteinPilot, and only the proteins that fulfilled the following criteria were validated: (1) peptide mass tolerance lower than 10 ppm, (2) 99% of confidence level in peptide identification, and (3) complete b/y ions series found in the MS/MS spectrum. Only proteins quantified with at least two unique peptides were considered. The quantitative data obtained by PeakView ® were analyzed using Perseus software 1.6.14 version (38) for statistical analysis and data visualization. MS data and search results files were deposited in the Proteome Xchange Consortium via the JPOST partner repository (https://repository.jpostdb.org) (39) with the identifier PXD027645 for ProteomeXchange and JPST001274 for jPOST (for reviewers: https://repository. jpostdb.org/preview/6417285016102b4b16aaa0; Access key: 3355). Interactome and pathway analysis were performed using Metascape (40) and machine learning-based bioinformatic QIAGEN's Ingenuity Pathway Analysis (IPA, QIAGEN Redwood City, www.qiagen.com/ingenuity)

Protein arrays
For the secretome analysis, a dot-blot protein array was used for cytokine profiling (Abcam). Briefly, membranes with 62 cytokine antibodies were blocked with the manufacturer's blocking bufier at room temperature (RT) for 30 min, and incubated o/n with 1 ml of undiluted cell-cultured media from OBC1 transfected cells (n = 3/condition). After washing, a biotinylated anti-cytokine antibody mixture was added to the membranes followed by incubation with HRP-conjugated streptavidin and then exposed to the manufacturer's peroxidase substrate.

Western-blotting
Equal amounts of protein (10 µg) were resolved in 4-15% stain free SDS-PAGE gels (BioRad). Protein extracts derived from OBC1 cells were electrophoretically transferred onto nitrocellulose membranes using a Trans-blot Turbo transfer system (up to 25V, 7min) (BioRad). Membranes were probed with primary antibodies at 1:1000 dilution in 5% nonfat milk or BSA according to manufacturer instructions. After incubation with the appropriate horseradish peroxidase-conjugated secondary antibody (1:5000), the immunoreactivity was visualized by enhanced chemiluminiscence (Perkin Elmer) and detected by a Chemidoc MP Imaging System (Bio-Rad). Equal loading of the gels was assessed using Stain-free imaging technology. Thus, protein normalization was performed by measuring total protein directly on the gels used for western blotting. After densitometric analyses (Image Lab Software Version 5.2; Bio-Rad), optical density values were expressed as arbitrary units and normalized to total protein levels.

Results and discussion
Multiple studies have pointed out that SARS-CoV-2 infection affects the chemosensory processing in animals and humans, especially the olfaction. SARS-CoV-2 is able to reach the OE and the OB in mice and golden Syrian hamsters (17,41). In rhesus monkeys, SARS-CoV-2 also invades the CNS primarily via the OB (42). Although it has been proposed that the olfactory nerve is not a likely route to brain infection in COVID-19 patients (43), several reports have shown the presence of SARS-CoV-2 (mRNA/protein levels or viral particles) in multiple brain areas including the OB derived from COVID-19 patients (20,(44)(45)(46)(47)(48)(49)(50). However, the involvement of OB cells remains controversial (51) In particular, the S protein has the capacity to cross the blood-brain barrier and reach the OB in mice (52). In view of these preclinical and neuropathological data, an indepth olfactory molecular characterization is necessary to unveil the missing links in the mechanisms of the sudden smell impairment induced by SARS-CoV-2 infection. To address this gap in knowledge, we have evaluated in-depth the OB cell response upon the intracellular expression of the SARS-CoV-2 structural proteins. The S glycoprotein is the main determinant of viral entry in target cells (52). The M glycoprotein is key for viral particle assembly (53,54), whereas the E protein is a multifunctional protein, acting as a viroporin in its oligomeric form, participating on viral assembly and virion release (55,56). The N protein packages the positive strand viral genome RNA into a helical ribonucleocapsid (RNP), being involved during virion assembly through its interactions with the membrane protein M and the viral genome (53,57). In general, SARS-CoV-2 structural protein factors are essential in the viral genome replication, receptor attachment and virion formation, promoting the viral spreading and pathogenesis.

Molecular derangements induced by SARS-CoV-2 structural proteins in OB cells
To examine the molecular alterations due to the presence of SARS-CoV-2 structural proteins in an olfactory cell system, GFP-tagged SARS-CoV-2 structural proteins S, M, N and E were expressed in OBC1 cells ( Figure 1A). It is important to note that these cells derived from the OB, structure of the forebrain that is located just above the nasal cavity, corresponding to the second stage of the olfactory pathway (23). After GFP-positive cell sorting, the olfactory metabolic imbalance induced by each SARS-CoV-2 protein was monitored across all experimental conditions by RNA-seq and shotgun proteomics. At transcriptome level, 1506, 1440, 359 and 1237 differential expressed genes were found between control-GFP and GFP-E, GFP-M, GFP-N and GFP-S SARS-CoV-2 proteins respectively (p-value adjusted <0.01; fold-change: 50%) (Supplementary Table 1). In contrast, around 15-23% of the quantified proteome was modulated by GFP-tagged SARS-CoV-2 structural proteins. Specifically, 435, 387, 594 and 584 proteins differential expressed proteins (DEPs) were found between Control-GFP and GFP-S, GFP-N, GFP-M and GFP-E respectively (p-value <0.05; fold-change 30%) (Supplementary Table 2). GFP-tagged SARS-CoV-2 proteins were detected by Western-blotting (Supplementary Figure 1A) and massspectrometry (Supplementary Table 2) except the GFP-E protein, due to its small size and the high content of aliphatic hydrophic aminoacids in its sequence (36% of Leu and Val) (58). Multiple commonalities and differences were observed when transcriptomic and proteomic datasets derived from the expression of GFP-tagged SARS-CoV-2 structural proteins  were compared (Figures 1B, C; Supplementary Figures 1B, 2A). Focusing on the most affected proteins, we identified a downregulation in Ak4 (Adenylate kinase 4), a protein that controls ATP levels through the regulatory activity of AMPK which plays a protective role against oxidative stress (59, 60). A detailed functional clustering revealed that SARS-CoV-2 structural proteins independently impact on multiple pathways at transcriptomic and proteomic levels ( Supplementary  Figures 1C, 2B; Supplementary Tables 3, 4). Comparing the functionality derived from both omic approaches, all SARS-CoV-2 structural proteins commonly interferes with multiple biofunctions such as extracellular matrix organization, deregulation of lipid metabolism, hemostasis, non-integrin membrane interactions and signaling by receptor tyrosine kinases ( Figure 1D; Supplementary Table 5). However, SARS-CoV-2 structural proteins commonly induced specific alterations with different impact on transcriptomic homeostasis and proteostasis. Considering the most statistically significant enriched pathways, SARS-CoV-2 proteins interfere with transcriptional activities associated to chondroitin/dermatan sulphate metabolism, platelet homeostasis, PPAR-dependent gene expression, opioid signaling, glycosylation and cytokine, IGF and PDGF signaling ( Figure 1D). At proteostatic level, SARS-CoV-2 structural proteins modulate aminoacid metabolism, biological oxidations, interleukin-12 signaling and processing of pre-mRNA and rRNA, between others ( Figure 1D). To decipher specificities and singularities associated to the expression of each SARS-CoV-2 structural protein at olfactory level, functional analyses were performed excluding the common deregulated proteome and taking into account only the mapping of unique DEPs associated to each SARS-CoV-2 protein (Supplementary Table 6). Interestingly, densely connected complexes with regulatory roles in GTP hydrolysis, calcium uptake into mitochondria (processing SMDT1), nuclear pore complex formation, vesicle-mediated transport, antigen presentation and peroxisomal protein import were commonly altered by structural SARS-CoV-2 proteins (Figure 2, upper). However, protein components of these complexes were differentially targeted by each structural viral protein (Figure 2, lower). Inferring more biologically interpretable results, GFP-M significantly modulated protein mediators located in organelle membranes that are mainly involved in trans-Golgi network vesicle budding, ferroptosis and cellular response to chemical stress ( Figures 3A, B). GFP-S independently impacted on nuclear proteins regulating cholesterol synthesis, DNA elongation and mRNA splicing between others (Figures 3A, B). GFP-E-dysregulated proteins presented a widespread localization, affecting a plethora of non-related functions such as mRNA processing, RAB geranylgeranylation and integrin interactions, whereas GFP-N specifically altered lipid, glutathione and glyoxilate metabolism (Figures 3A, B). In contrast to recent results with full viral particles, these data generated at the organelle-complex-pathway axis demonstrate the enriched multifactorial functions associated to each SARS-CoV-2 protein.
Interestingly, 14-18% of the perturbations generated by SARS-CoV-2 M, N and E proteins (at the level of protein-coding gene) in OBC1 cells were common with data derived from pulmonary cells expressing these structural viral proteins (26) ( Figure 3C). However, only 3% of the proteome modulated by the S protein was shared across studies, indicating that a great proportion of proteostatic events intercepted by structural SARS-CoV-2 proteins are cell-dependent. In addition, part of the altered molecular routes identified in our in-vitro system have been recently detected in the OB from COVID-19 subjects (62) (Supplementary Figure 3), pointing out specific pathways directly targeted by SARS-CoV-2 structural proteins beyond the CNS affectation and/or the systemic immune response to SARS-CoV-2 infection.

SARS-CoV-2 structural proteins differentially modulate the activation profile of survival routes and mitochondrial homeostasis in OB cells
To enhance the analytical outcome of proteostatic alterations, proteome-scale interaction networks were constructed merging the differential proteomes induced by SARS-CoV-2 structural proteins. As shown in Figure 4, both p38 MAPK and Akt appeared as principal nodes in protein interactome maps. Although much effort has been spent on characterizing the pleiotropic effects of SARS-CoV-2 infection in multiple biological contexts, there is no detailed information about the impact of each structural SARS-CoV-2 protein on the activation dynamics of the cell survival. Even though changes in p38 MAPK and Akt expression were not detected in our proteomic analysis, the alteration of some targets might be indicative of a dysfunctional state when SARS-CoV-2 structural proteins are continuously expressed. Subsequent experiments were performed to analyze the activation state of p38 MAPK and Akt signaling pathways across GFP-SARS-CoV-2 structural protein expressing-OBC1 cells ( Figure 5). As shown in Figure 5A, the presence of SARS-CoV-2 structural proteins induced a p38 MAPK and Akt inactivation. Specifically, the independent expression of SARS-CoV-2 S and N proteins inactivated both kinases. In contrast, the M protein specifically inactivated Akt and the E protein interfered with p38 MAPK. To complement our signaling mapping, other stress-responsive kinases relevant to OB homeostasis (63) were checked. N and S proteins also induced the inactivation of GSK3 and the SEK1 stress-activated protein kinase (SAPK) axis ( Figure 5B). However, S protein expression induced the inactivation of PKA whereas the expression of the N protein enhanced the active phosphorylation of the PKA catalytic subunit ( Figure 5B). No appreciable changes were observed in the activation profile of  Figure 4). These data demonstrate that SARS-CoV-2 structural proteins differentially impact on the survival potential of OB cells. Many viral infections trigger the activation of the p38 MAPK and Akt signaling pathways for an efficient replication, including SARS-CoV-2 (64)(65)(66). p38 MAPK and Akt inactivation induced by the expression of SARS-CoV-2 structural proteins could be part of an antiviral mechanism triggered in olfactory cells. In fact, pharmacological inhibition of p38 MAPK has antiviral efficacy (67,68) and Akt inhibitors have been proposed as candidate drugs (69). Therefore, our results are in agreement with evidence pointing to both kinases as potential  (73). However, the mechanistic connection with SARS-CoV-2 is not known. We observed an enhanced PKA activity when SARS-CoV-2 N protein is overexpressed in OBC1 cells. Interestingly, the coe x p r e s s i o n o f P K A a n d N p r o t e i n g e n e r a t e s a poliphosphorylated N dimers able to sequestrate cellular 14-3-3 isoforms (74). It has been recently demonstrated the existence of physical and functional interactions between SARS-CoV-2 and host  Frontiers in Immunology frontiersin.org mitochondria, serving this organelle as a subcellular platform for anti-SARS-CoV-2 immunity (75). Although mitochondrial proteostasis was modified by all SARS-CoV-2 structural proteins ( Figure 3B), our network analysis uncovered a generalized down-regulation of mitochondrial complex I subunits specifically induced by the expression of SARS-CoV-2 E and M proteins (Figures 6A, B). Mitochondrial prohibitin complex (Phb1/Phb2 subunits) modulates respiratory complex assembly and mitochondrial homeostasis, triggering anti-oxidant effects (76, 77). According to Biogrid database (https:// thebiogrid.org/), Phb complex interacts with SARS-CoV-2 proteome at multiple levels. Specifically, Phb1 interacts with ORF9B/C, NSP2, NSP3, NSP6 and NSP8 whereas NSP2, NSP3, NSP6, ORF6, ORF9B, E and M proteins are Phb2 interactors. As shown in Figure 6C, a depletion in Phb1 levels was observed in E protein expressing OBC1 cells. In contrast, a significant increment in Phb2 levels was evidenced when S and N proteins were produced. In general, Phb1 repression triggers a Frontiers in Immunology frontiersin.org concomitant reduction of its partner Phb2 and viceversa (78), indicating that SARS-CoV-2 structural proteins has the capacity to impact on mitochondrial homeostasis differentially interfering with the functional interdependency of Phb subunits.

SARS-CoV-2 structural proteins induce a dissimilar immunological effectome in OB cells
The vast amount of information about the relationship between immune system and COVID-19 points out that multiple aberrations of innate and acquired immunity are induced by SARS-CoV-2 infection. Specifically, astrogliosis, microgliosis and infiltration by cytotoxic T lymphocytes in the OB have been characterized in COVID-19 subjects (45). This olfactory neuroinflammation is different from the observed in other brain areas (79). In addition, it has been proposed that olfactory ensheating cells may produce glia transit tubules through which cytokines and chemokines might migrate across the OE-OB axis (22). However, little is known about the immunological mechanisms triggered by each SARS-CoV-2 structural protein at olfactory level. In agreement with previous reports (26), specific gene/protein subsets modulated by each SARS-CoV-2 structural protein profile a pro-inflammatory signature (Supplemental Table 7, Supplementary Figure 5), inducing specific changes in the OBC1 immunological effectome. As shown in Figure 7, S protein expression impacted on regulators of IL-6, IL-10 and IL-12 and in the differentiation of mature B cells. Alterations in protein-coding genes relevant in inflammasomes, lymphocyte migration, response to IL-7 and B cell activation were induced by the expression of GFP-N protein expression in OBC1 cells (Figure 7). On the other hand, the M protein had the potential to interfere with lymphocyte activation, T cell differentiation, macrophage differentiation as well as in the IL-10, IL-27, and IL-35 signaling, whereas SARS-CoV-2 E protein specifically impacts on adaptive immune response, IL-17 pro-duction, B cell receptor signaling and natural killer cell differentiation between other processes (Figure 7).
We wanted to partially map the consequences over the intracellular immunological effectome dependent on SARS-CoV-2 structural proteins. Cytokines and growth factors secreted by each stable OB cell line constitutively expressing a B A FIGURE 5 Levels and residue-specific phosphorylation of (A) p38 MAPK, Akt and (B) PKA, SEK1 and GSK3 in GFP-SARS-CoV-2 structural protein expressing-OBC1 cells. Equal loading of the gels was assessed by stain-free digitalization. Panels show histograms of band densities. Data are presented as mean ± SEM from four independent biological replicates. *P < 0.05 vs. GFP; **P < 0.01 vs. GFP; *** P < 0.001 vs. GFP (a.u: arbitrary units).
GFP-tagged SARS-CoV-2 structural protein may provide novel insights into the inflammatory imbalance produced by SARS-CoV-2. As shown in Figure 8, GFP-tagged M and E protein expression in OBC1 cells induced more extracellular changes. Part of the 62 secreted cell-cell signaling molecules analyzed were commonly altered between the structural proteins ( Figure 8). Specifically, the increment of RANTES (CCL5) and MIP3a (CCL20) levels in the OBC1 secretome was independently induced by all structural proteins. Both chemokines are elevated in sera from COVID-19 patients (80). IFN gamma and MCP5 (CCL12) extracellular levels were also increased in SARS-CoV-2 S, M and E protein-expressing OBC1 cells (Figure 8). Although most of the differential soluble mediators have been previously related to the immune cell recruitment associated to COVID-19 cytokine storm (81), these data contribute to the better understanding of the multifunctional immunomodulatory properties of SARS-CoV-2 M, N, S and E proteins beyond their intrinsic role in virion formation.

Proteotranscriptomic data integration by machine learning unravels TGF-beta signaling route as a confluent activation node by SARS-CoV-2 structural proteins
It has been previously demonstrated that sophisticated machine learning (ML) approaches are useful workflows to integrate multi-omics data with the aim to discover complex interconnections between different type of entities (82,83). To uncover regulatory complex composed by molecular effectors potentially modulated by SARS-CoV-2 structural proteins at olfactory level, transcriptomic and proteomic datasets were integrated and mined by a machine-learning approach in which multiple entities (genes, proteins, upstream regulators, pathways, diseases) were interconnected in the form of knowledge graphs. Considering the IPA Z-score as a statistical parameter based on the match between expected relationship direction and experimental gene/protein expression, functional graphs were constructed. As shown in Figure 9, common hubs B C A FIGURE 6 Functional protein interactome maps for differentially expressed proteins in GFP-SARS-CoV-2 E and M protein-expressing OBC1 cells, revealing a mitochondrial impairment (A, B). OB expression of Phb complex by Western blot (C). Equal loading of the gels was assessed by stain-free digitalization. Panels show histograms of band densities. Data are presented as mean ± SEM from four independent biological replicates. *P < 0.05 vs. GFP; **P < 0.01 vs. GFP; *** P < 0.001 vs. GFP (a.u: arbitrary units).
Laché n-Montes et al. 10.3389/fimmu.2022.866564 Frontiers in Immunology frontiersin.org were predicted for each SARS-CoV-2 structural protein such as TGFB1/B3, SMAD and TEAD proteins, EGF, SPP1 and EDN1. It has been previously observed that SARS-CoV-2 N protein interacts with SMAD proteins enhancing TGF-beta signaling (84) whereas SARS-CoV-2 S protein triggers a transcriptional response associated to TGF-beta signaling (85). Interestingly, a therapeutic strategy focused on the inactivation of TGF-beta using integrin inhibitors has been proposed to mitigate COVID-19 severity (86). EDN1 (Endothelin-1) has been considered a neuroprotective factor that participates in the olfactory response modulation through the uncoupling of gap junctions (87,88). Interestingly, elevated EDN1 levels together with a sustained inflammation has been observed 3 months after recovery from acute COVID-19 symptoms (89). SPP1 (osteopontin) participates in the OB synaptic plasticity (90) and acts as a molecular brake regulating neuroinflammatory response to chronic viral infections (91). Involved in the enhancing production of interferon-gamma and interleukin-12, SPP1 is also overproduced in COVID-19 patients (92). Although our study has uncovered the existence of commonalities and differences in the SARS-CoV-2 structural protein functionality, potential limitations exist that warrant discussion. It is important to note that SARS-CoV-2 structural and non-structural protein interactions are absent in our study (93,94). We are aware that using our multi-omic approach, the molecular dimension generated by posttranslational The expression of SARS-CoV-2 structural proteins differentially modulate the immunological capacity of OBC1 cells. Processes related to immune system were directly extracted from the functional analysis of proteotranscriptomic datasets (see Supplementary Table 6 for more details). Extracellular cytokine profiling of GFP-SARS-CoV-2 structural protein-expressing OBC1 cells. The analysis of 62 cytokines/growth factors was performed in the cell media derived from all biological conditions using a dot-blot protein array method. Four independent experiments were performed. Data are presented as mean ± SEM. *p < 0.05, **p < 0.01, and ***p < 0.001 vs GFP. modifications (i.e phosphorylation) not only in host proteins but also in the SARS-CoV-2 structural proteins (26) has not been considered in this study, hampering the characterization of potential substrates modulated by the virus-hijacked kinase activation profiles. Due to our experimental design, the overexpression of exogenous proteins may generate drawbacks concerning protein misfolding, localization and regulation as well as intrinsic limitations associated to GFP expression systems (95). However, our three-layered approach (transcriptomicsproteomics-secretomics) coupled to Systems Biology has Machine learning-based analysis proposes a potential activation of TGFbeta/SMAD axis commonly induced by all SARS-CoV-2 structural proteins. The functional networks use a precomputed table containing inferred relationships between molecules, functions, diseases, and pathways obtained and scored by a machine learning algorithm operating entirely on prior knowledge. The heuristic graph algorithm present in IPA software was optimized to create a manageable network that brings together the most significantly activated (in orange; positive Z-score) or inhibited (in blue; negative z-score) upstream regulators, diseases, functions, and pathways from the differential proteotranscriptomic alterations demonstrated the capacity to detect therapeutic targets potentially useful against COVID-19. Hubs detected by our network analysis such as p38 MAPK, Akt and GSK-3 further experimentally validated in our olfactory in vitro system, have been proposed as therapeutic targets. In particular, pharmacological inhibition of p38 MAPK using drugs and compounds such as ARRY-797 (phase 2/3), MAPK13-IN-1 (preclinical), SB203580 (preclinical) and ralimetinib (phase 2) has antiviral activity (Bouhaddou, M. et al., 2020). The Akt inhibitor MK-2206 as well as specific GSK-3 inhibitors (COB-187, LY2090314, and AZD-1080) have been also proposed as candidate drugs against SARS-CoV-2 (Fagone, P et al., 2020; Ghazanfari D et al., 2022). Based on the olfactory cellular system used, additional experiments are needed to verify the specific role of the SARS-CoV-2 structural proteome in different human olfactory cellular contexts (23, 96). As shown in previous reports performed at cellular and tissular levels (24,26), the application of proteomics in different olfactory cell layers as well as in olfactory areas directly derived from COVID-19 individuals, would increase our understanding of not only the early smell impairment associated to SARS-CoV-2 infection but also the olfactory recovery potential in COVID-19 patients.

Data availability statement
Transcriptomic and proteomics data and results files were deposited in GEO accession GSE182849 and ProteomeXchange Consortium via the JPOST partner with the identifier PXD027645 for ProteomeXchange and JPST001274 for jPOST.  Comparative analysis between deregulated pathways observed in OBC1 cells expressing SARS-CoV-2 structural proteins and altered pathways observed at the level of the OB from COVID-19 subjects.

SUPPLEMENTARY FIGURE 5
Differential immunological effectome (genes/proteins) induced by the expression of SARS-CoV-2 structural proteins in OBC1 cells.

SUPPLEMENTARY TABLE 4
Gene-Ontology (GO) functional analysis derived from proteomics.

SUPPLEMENTARY TABLE 5
Top-20 most significant disrupted pathways. Enrichment and heatmap SUPPLEMENTARY TABLE 6 Functional analysis of DEPs exclusively modulated by a SARS-CoV-2 structural protein.

SUPPLEMENTARY TABLE 7
Immunological functionality derived from the interlocking of transcriptome and proteome datasets.