Cerebrovascular recovery drives restoration of neurometabolite levels after mild COVID-19

Abstract

SARS-CoV-2 infection can trigger broad acute-stage endothelial dysfunction that is often followed by neurocognitive post-acute sequelae (PASC) chronically. It is thus of interest to evaluate long-term post-infection cerebrovascular function dynamics and their effects of neuronal health. Using functional magnetic resonance imaging (fMRI), we examined cerebrovascular reactivity (CVR) magnitude and delay across a broad post-infection timeline (3–59 months) in 69 participants who previously had mild cases of COVID-19. We also assessed the relationships between CVR and neurometabolite markers of neuronal health using magnetic resonance spectroscopy (MRS) in the thalamic region and the corona radiata. Increasing time since infection (TSI) was associated with shorter CVR delay in global gray (GM) and white matter (WM), with no effects of TSI on CVR magnitude. Parcellation of the GM revealed TSI-dependent decreases in CVR delays in nine of the 10 examined GM regions. CVR delay was inversely related to total choline (tCho), creatine (tCr), and N-acetylaspartate (tNAA) levels in the thalamic region, but not in the corona radiata. The results suggest slower cerebrovascular reactivity follows mild COVID-19 and eventually resolves spontaneously, albeit on a protracted timeline. Further, improved levels of tCho, tCr, and tNAA in the GM are associated with this functional cerebrovascular recovery.

Introduction

COVID-19 is caused by the SARS-CoV-2 virus and is characterized by widespread and potentially life-threatening acute-stage pulmonary effects. Significant long-term consequences of this virus in other organs have also emerged, and can manifest themselves as a chronic multi-system condition known as long COVID, or post-acute sequelae of COVID-19 (PASC).^1–3 Mounting evidence implicates SARS-CoV-2 in vascular damage affecting multiple organs, inflicted both directly and indirectly, as having a role in lasting post-infection symptomology.^4–6 This is consistent with the elevated rates of acute-stage myocardial infarctions and strokes as well as chronic vascular complications that can follow COVID-19.^7,8 The role of post-infection damage to the brain’s vasculature in particular has recently been implicated in the neurocognitive symptoms of long PASC.^9–11 These cerebrovascular impairments, while detrimental on their own, can be further compounded by endothelial and epithelial damage elsewhere in the body, such as the intestines, with harmful gut-derived compounds proliferating systemically and potentially gaining access to the brain.^12,13 In moving toward elucidating long-term vascular sequelae of COVID-19, here we describe the cerebrovascular reactivity (CVR) dynamics across the post-infection timeline in a group of individuals who previously had mild, non-hospitalized COVID-19. In addition, we evaluate the effects that these vascular processes may have on neuronal health, as measured via levels of major neurometabolites in the brain.

The SARS-CoV-2 virus binds to angiotensin-converting enzyme 2 (ACE2) receptors, which are widely expressed by the lungs as well as other organs.^14,15 Extensive ACE2 expression by endothelial cells facilitates direct viral interactions with the vasculature, causing receptor downregulation and damage to the blood vessel walls.^16,17 Endothelial integrity can be further compromised indirectly by the effects of systemic cytokines and oxidative stress in the context of hyperinflammatory cascades triggered by the infection.^18,19 Impairment of the endothelium’s critical role in brain arteriolar function can be exacerbated by capillary contractile dysfunction that relies on pericytes, which also interact with SARS-CoV-2 and have been implicated in cognitive decline.^20–22 While endothelial injury is well-characterized in severe acute-stage COVID-19, it remains unclear whether mild COVID-19 produces cerebrovascular or neuronal changes that may persist for months or years after infection. In moving toward improved characterization of long COVID-19 pathogenesis that often follows unremarkable acute-stage infections, we focused on examining cerebrovascular function and its neuronal implications after mild, non-hospitalized COVID-19 cases.

We used functional magnetic resonance imaging (fMRI) of CVR to measure the ability of brain blood vessels to respond to hypercapnia, which was induced by respiring of a gas mixture of CO₂—a rapid and potent vasodilator. ²³ Breath-hold paradigms have previously been used to induce hypercapnia during neuroimaging, yielding reliable and clinically relevant indicators of cerebrovascular health. ²⁴ However, such designs can introduce motion artifacts and variance due to individual differences in pulmonary function and task compliance, with fixed respired CO₂ and end-tidal CO₂ targeted stimulus paradigms emerging as the preferred CVR administration methods. ²⁵ Our approach to use highly controlled, end-tidal CO₂ targeted delivery of CO₂ gas mixture via a mask during normal respiration addresses these limitations, while also allowing for additional procedures to monitor participant safety. With hypercapnia induced in an intermittent manner over several blocks of alternating administration of CO₂ gas mixture and normal room air, the paradigm allows for the measurements of two main response variables: CVR magnitude and delay. Specifically, CVR magnitude refers to the increase of the blood oxygenation level dependent (BOLD) signal in response to a rise in arterial partial pressure of CO₂, and CVR delay is the temporal lag of this response relative to the onset of CO₂ administration. Collectively, these two metrics evaluate the ability of brain vasculature to respond to exogenous vasoactive stimuli, which is highly reflective of general cerebrovascular endothelial health.^26,27 As such, CVR offers unique information about brain vascular function, and can effectively evaluate the trajectories of cerebrovascular health following mild SARS-CoV-2 infection.

Our overarching hypothesis is that SARS-CoV-2 inflicts endothelial damage that subsequently results in lasting cerebrovascular dysfunction. Testing this necessitates accurate identification and recruitment of healthy control participants who have never had COVID-19—an increasingly challenging task given that the majority of the population has been infected with SARS-CoV-2 at some point. As such, in this analysis we examined the post-infection trajectories of CVR dynamics in a cohort consisting entirely of participants who previously had mild COVID-19, with time since infection (TSI) being the main predictor of interest. Our specific predictions were that CVR magnitude and delay would be diminished and long, respectively, in the early post-infection timeframe, and that long-term normalization would be observed as increases in CVR magnitude and decreases in its delay.

In addition to evaluating cerebrovascular dynamics following COVID-19, we also examined the associations of post-infection CVR metrics with neuronal tissue health. Specifically, we quantified the effects of CVR magnitude and delay on total concentrations of choline-containing compounds (tCho), creatine (tCr), and N-acetylaspartate (tNAA) in the brain. These neurometabolites are highly abundant in brain tissue and serve as reliable markers of neuronal health.^28–30 Choline-containing compounds maintain cellular lipid membranes, support synthesis of neurotransmitter acetylcholine, and are glial markers.^31,32 While choline is recycled via membrane turnover, its de novo brain supply is entirely dependent on transport from blood, highlighting the importance of healthy cerebrovascular function for cholinergic homeostasis.^33,34 Creatine and phosphocreatine play important roles in supporting neuronal metabolism by buffering adenosine triphosphate (ATP), and while ATP is synthesized locally via astrocyte–neuron coupling, around 50% is transported into the brain from the circulation, thus being partially dependent on healthy brain vascular function. ³⁵ Neuronal mitochondria synthesize NAA, which is involved in myelination and lipid synthesis and serves as a reliable indicator of overall brain tissue health. While all NAA is produced locally, its levels are indirectly influenced by cerebrovascular malnourishment, resulting in mitochondrial dysfunction and myelin lipid synthesis. ³⁶ Collectively, tCho, tCr, and tNAA levels serve as excellent markers of neuronal health, and are expected to be affected by any potential dysregulation of cerebrovascular function. We hypothesized that levels of major metabolite concentrations would be associated with any post-infection temporal changes in CVR metrics, as these cerebrovascular dynamics are highly reflective of oxygen delivery that supports ATP synthesis and mitochondrial function in brain tissue.

Methods

Study participants

The study protocol was approved by the University of New Mexico Health Sciences Center Human Research Protections Office in accordance with the Belmont Report, and all individuals consented to study participation. Sixty-nine participants (age = 43.7 ± 14.4) with histories of non-hospitalized SARS-CoV-2 infections, as established via self-reported positive tests, were recruited from the Albuquerque community. If recalled, exact dates of infection were recorded, with the first of the month used if only the month and year were provided. To minimize acute infection effects, inclusion was limited to individuals who had COVID-19 at least 3 months prior to enrollment. If vaccinated for SARS-CoV-2, participants were enrolled at least 14 days after the last dose to avoid transient vaccine effects. The study group TSI range was 3.13–59.17 months (M = 18.19, SD = 13.82), with a 94% vaccination rate and 97% of participants having had only one infection.

Exclusion criteria included: serious psychiatric, developmental, neurological, gastrointestinal, metabolic, or vascular conditions, moderate-to-severe traumatic brain injury, tobacco use, and history of or ongoing substance use disorder. ³⁷ These data are part of an ongoing study on gut–brain dynamics, so exclusion criteria related to gastrointestinal function also applied: use of motility agents <2 days before the study, laxatives <1 week, probiotics <2 weeks, antibiotics <1 month, and bowel purgatives <2 weeks. Participants followed identical dietary restrictions and avoided exercise before study visits to stabilize their gut microbiomes and homogenize any vascular effects of gut-derived compounds, adhering to the North American Consensus for the Hydrogen and Methane-Based Breath Testing in Gastrointestinal Disorders. ³⁸

Anatomical MRI acquisition

Neuroimaging sequences were collected during the same session on a Siemens Prisma 3T located at the Mind Research Network in Albuquerque, NM, USA. A high-resolution structural T1 volume was obtained using echo time (TE) = 1.61, 3.47, 5.33, 7.19, and 9.05 ms; repetition time (TR) = 2530 ms; inversion time (TI) = 1200 ms; flip angle = 7°; field of view (FOV) = 256 mm; and voxel size = 1.0 mm³. The structural scans were used for MRS voxel placements in the thalamic and corona radiata regions and for calculating CVR values within these voxel masks. Structural volumes were processed using FreeSurfer to isolate high-fidelity subcortical GM structure for MRS metabolite quantification. ³⁹ Anatomical scans were also used for fMRI normalization, and were segmented using statistical parametric mapping (SPM) to separately evaluate GM and WM CVR metrics, as GM and WM respond to hypercapnia differently.^40,41

CVR data acquisition

The hypercapnia paradigm was adopted from Lu et al., utilizing a RespirAct system for precise gas delivery. ⁴² Participants were fitted with a mask connected to the gas delivery device that blended the CO₂ mixture (5% CO₂, 21% O₂, remainder N₂). Participants were instructed to breathe normally throughout the CVR experiment. The RespirAct system monitored each participant’s breath-by-breath inspiration to adjust CO₂ gas delivery in real-time, targeting exhaled end-tidal CO₂ (ETCO₂) partial pressure of 10 mmHg relative to baseline. The hypercapnia paradigm was administered via a block design of four 35-s “on” and “off” conditions, alternating between CO₂ gas and room air. Continuous ETCO₂ data were recorded and used as a regressor in modeling the BOLD response to hypercapnia. The CVR fMRI parameters were: TE = 37 ms, TR = 800 ms, 52° flip angle, FOV = 208 mm, and voxel size = 2.0 × 2.0 × 2.0 mm, 72 slices. Field mapping scans were acquired with reverse-phase encoding directions (A-P and P-A).

CVR data processing

Neuroimaging data were processed using open-source toolboxes, including AFNI, advanced normalization tools (ANTs), FMRIB software library (FSL), SPM, and custom MATLAB scripts.^43–45 The fMRI data were de-spiked and corrected for motion, slice-timing, and susceptibility-induced distortions using the field mapping scans. Motion regressors, their derivatives, and mean framewise displacements were used to assess motion, with no scans having had mean frame-to-frame displacement exceeding 3 mm. The functional data were aligned to the T1, normalized to MNI space, and smoothed using a 6 mm kernel. ⁴⁶ The ETCO₂ data were processed using peak-to-peak interpolation, and values were sampled form the resulting time course at time points corresponding to the fMRI volume acquisitions. The resulting ETCO₂ data were expanded into a series of time-lagged regressors for modeling vascular reactivity. For each voxel, the regressor with the highest Pearson correlation with the BOLD signal was selected, and the associated correlation coefficient and its corresponding lag time reflected the CVR fit and delay, respectively. Voxel-wise CVR magnitude was then calculated using the fit measure as %BOLD signal change per mmHg CO₂, and delay was defined as the time of maximal BOLD response relative to the ETCO₂ trace.

The CVR data were examined using both global and region of interest-based (ROI) approaches. For global analyses, GM and WM were assessed separately, as the WM exhibits lower CVR magnitude and longer delay. ⁴⁰ The GM and WM were segmented using SPM. Following normalization, MNI atlas-defined ROIs to assess whether separate brain regions differ in their post-infection CVR trajectories. ⁴⁷ For each participant, voxel-wise CVR magnitude and delay values were averaged within each area of interest to produce summary metrics for the following 12 regions: global GM, global WM, left/right (l/r) thalamus, l/r/ frontal, l/r parietal, l/r/ occipital, and l/r temporal lobes.

MRS data acquisition

Two MRS voxels were scanned: one placed over the left thalamus and another over the left corona radiata (CR), consistent with similar experiments published previously.^48,49 The thalamic region was chosen to differentially target the GM and minimize WM contributions. The ovoid shape of the thalamus maximizes the GM:WM ratio within the voxel, as doing so is challenging in the thin cortical GM. Further, the thalamus is one of the most CVR-responsive regions in the brain, being vascularized via the posterior cerebral artery that is 43% and 23% more reactive to hypercapnia than the anterior and middle cerebral arteries, respectively. ⁵⁰ Importantly, thalamic neurometabolites have been shown to be susceptible to various neurological conditions, making it a suitable candidate for examining any alterations due to cerebrovascular deficiencies.^51–55 The CR voxel placement was chosen for analogous anatomical reasons as the thalamic voxel, maximizing WM volumes within the imaged area. Collectively, the thalamic and CR placements facilitated our goal to make differential inferences about vascular effects on neuronal health in the GM and WM.

A semi-Localization by Adiabatic SElective Refocusing (semi-LASER) sequence was used to obtain MRS data, which were part of a broader diffusion MRS experiment that included acquisition of diffusion-weighted spectra. ⁵⁶ Given our interest in conventional MRS metabolite concentrations, only the gradient-free b0 (0 s/mm²) spectra were used for this analysis. The following parameters were used: TR = 5000 ms, TE = 120 ms, spectral width = 2500 Hz, number of complex points = 2048, number of excitations (NEX = 32). Additional scans without VAriable Power and Optimized Relaxations delays (VAPOR) water suppression were collected at each voxel location (NEX = 4). ⁵⁷

MRS data processing

The MRS data were processed using MRSpa, LCModel, and MATLAB.^58,59 The MRSpa toolbox was used to apply corrections for phase jumps in the free induction decay (FID) signal and eddy currents. ⁶⁰ Frequency and phase instabilities were adjusted by minimizing the cross-correlations of frequency and phase differences between the 32 spectral averages. The averaged spectra were assessed for quality via linewidths (15 Hz threshold) of unsuppressed water peaks and NAA peak signal-to-noise ratios (SNR = 5 threshold), calculated using the FID Appliance toolbox.^61,62 LCModel was used to estimate the peak areas for choline compounds (tCho; choline, phosphocholine, and glycerophosphocholine), creatine (tCr; creatine and phosphocreatine), and N-acetylaspartate (tNAA, N-acetylaspartyl-glutamate and acetyl moiety of NAA). Cramér–Rao lower bound (CRLB) metabolite concentration variance estimates (%SD) were used as data reliability indicators, and datasets with tCho, tCr, or tNAA %SD > 20 were excluded. ⁶³

Subject-specific CSF, GM, and WM volume fractions were calculated within each MRS voxel to account for metabolite and water T1 and T2 relaxation time biases across different tissues (Figure 1). An MRSpa component based on SPM was used to segment the brain tissue in the CR voxel. Since deep GM segmentation was often pixelated, FreeSurfer was used to isolate deep GM regions that may fall within the thalamic voxel: thalamus, globus pallidus, hypothalamus, and brainstem. These structures were used to calculate the GM and WM composition of the thalamic voxels in each participant. The resulting tissue type fractions within each MRS voxel were used to adjust tCho, tCr, and tNAA concentrations using metabolite- and region-specific T1 and T2 relaxation times. The relaxation times were based on existing work using Point RESolved Spectroscopy (PRESS) MRS, and were multiplied by factors accounting for tissue- and metabolite-specific relaxation time differences between the PRESS and semi-LASER sequences.^64–66 A similar metabolite quantification approach described elsewhere lists the specific relaxation values and the formula used to adjust metabolite concentrations. ⁴⁹

Figure 1.

Th and CR MRS voxel placements in a representative participant, with tissue segmentations within each voxel showing that they were mostly occupied by their target tissue types (Th = GM, CR = WM). Specific tissue volumes obtained from these segmentations were also used to adjust the MRS concentrations for tissue-specific T1 and T2 relaxation times of water and individual metabolites (tCho, tCr, and tNAA). The MR spectra acquired from the Th and CR voxels are presented from the same participant. Vertical axes are not shown due to signal amplitudes being in arbitrary units.

Statistical analyses

Multivariate linear regression analyses were conducted in SPSS to assess the effects of TSI on CVR metrics, with TSI and age used as the predictor variables and either GM magnitude, GM delay, WM magnitude, or WM delay as the dependent variables. TSI and age variables were log-transformed to address positive skews. Age was included due to improved model qualities, as per Akaike Information Criteria (AIC). To evaluate the effects of CVR metrics on metabolite levels, a similar approach utilized a series of linear regression models. The age variable was omitted due to non-significance and lack of importance as per the AIC analyses. The thalamic voxel MRS data were used as a proxy for global GM metabolite profile, while the CR data as an analogous proxy for global WM. As such, the following models were used to examine the relationships between GM CVR metrics and thalamic metabolite levels: tCho_Thal ~ TSI + GM magnitude, tCr_Thal ~ TSI + GM magnitude, tNAA_Thal ~ TSI + GM magnitude, tCho_Thal ~ TSI + GM delay, tCr_Thal ~ TSI + GM delay, and tNAA_Thal ~ TSI + GM delay. Likewise, the relationship between WM CVR metrics and CR metabolite levels was evaluated using: tCho_CR ~ TSI + WM magnitude, tCr_CR ~ TSI + WM magnitude, tNAA_CR ~ TSI + WM magnitude, tCho_CR ~ TSI + WM delay, tCr_CR ~ TSI + WM delay, and tNAA_CR ~ TSI + WM delay.

Results

The group CVR averages were: GM magnitude = 0.211 ± 0.051 %BOLD/mmHg, GM delay = 13.699 s, WM magnitude = 0.107 ± 0.025 %BOLD/mmHg, WM delay = 14.115 ± 2.180 s. The data were normally distributed and contained no outliers (3 SD cutoff; GM magnitude lower bound/upper bound (LB/UB) = 0.056/0.366 %BOLD/mmHg; GM delay LB/UB = 6.315/21.217 s; WM magnitude LB/UB = 0.031/0.183 %BOLD/mmHg; WM delay LB/UB = 7.456/20.883 s).

Regression analyses with age and TSI as predictors and either CVR delay or magnitude as the outcome variables revealed differential results. The model examining delay was significant in both global gray and white matter (Figure 2; GM F (2, 68) = 4.916, p = 0.010, R²_Adj = 0.103; WM F (2, 68) = 5.772, p = 0.005, R²_Adj = 0.123). In both the GM and WM models, TSI was the only term that significantly predicted delay (TSI × GM delay β = −0.326, t = −2.830, p = 0.006; age × GM delay β = 0.182, t = 1.575, p = 0.120; TSI × WM delay β = −0.338, t = −2.971, p = 0.004; age × WM delay β = 0.215, t = 1.883, p = 0.064). Increasing TSI significantly predicted shorter CVR delay in 9/10 examined brain regions after a false discovery rate multiple comparisons correction (p_FDR < 0.05). ⁶⁷ In predicting CVR magnitude, the overall effect of age and TSI was not observed in either tissue type (GM F (2, 68) = 1.577, p = 0.214, R²_Adj = 0.017; WM F (2, 68) = 0.672, p = 0.514, R²_Adj = −0.010). Neither TSI nor age were predictive of CVR magnitude (TSI × GM magnitude β = 0.034, t = 0.280, p = 0.781; age × GM magnitude β = −0.214, t = −1.771, p = 0.081; TSI × WM magnitude β = 0.082, t = 0.674, p = 0.503; age × WM magnitude β = −0.122, t = −0.994, p = 0.324). TSI was not associated with CVR magnitude in any of the examined GM regions.

Figure 2.

The effects of TSI on CVR magnitude and delay in global GM and WM. The scatter plots in the middle panel show that the delay of the CVR response to CO₂ became significantly shorter with increasing TSI in both global GM and WM. However, as seen in the scatter plots in the top panel, CVR magnitude showed no significant associations with TSI in either global GM or WM. The bottom panel shows effects of TSI on CVR magnitude and delay in individual brain regions. None of the regions displayed significant changes in CVR magnitude with increasing TSI. All regions except for the right frontal lobe displayed significant decreases in CVR delay with increasing TSI.

Out of the 69 CVR sessions, usable MRS was collected for 44 thalamic and 43 CR voxel placements (both voxels N = 35, thalamic only N = 9, CR only = 8, no MRS N = 17). Water spectral linewidths and NAA SNR were: LW_Thal = 10.52 ± 1.61 Hz, LW_CR = 8.19 ± 2.28 Hz, SNR_Thal = 26.68 ± 5.45, SNR_CR = 58.46 ± 17.42. The MRS data were normally distributed and contained no outliers based on the 3 SD cutoff criterion (tCho_Thal LB/UB = 0.500/2.847 µMol/g; tCr_Thal LB/UB = 1.887/10.616 µMol/g; tNAA_Thal LB/UB = 4.595/17.927 µMol/g; tCho_CR LB/UB = 0.671/2.260 µMol/g; tCr_CR LB/UB = 2.073/7.456 µMol/g; tNAA_CR LB/UB = 5.291/12.936 µMol/g). The observed mean concentrations and CRLBs for each metabolite and voxel location were as follows: tCho_Thal = 1.67 ± 0.39 µMol/g, CRLB = 6.05% ± 1.45%; tCr_Thal = 6.25 ± 1.46 µMol/g, CRLB = 7.36% ± 2.08%; tNAA_Thal = 11.26 ± 2.22 µMol/g, CRLB = 3.73% ± 1.19%; tCho_CR = 1.47 ± 0.26 µMol/g, CRLB = 3.93% ± 0.80%; tCr_CR = 4.81 ± 0.89 µMol/g, CRLB = 5.05% ± 0.95%; tNAA_CR = 9.14 ± 1.24 µMol/g, CRLB = 2.47% ± 0.94%.

Given the TSI effects on CVR delay, we examined whether CVR delay measures were predictive of tCho, tCr, and tNAA concentrations within the two thalmic and CR MRS voxels. Regression analyses were performed using tCho, tCr, and tNAA levels as the dependent variables, with predictors being TSI and MRS region-specific CVR delays (delay values averaged across all voxels within either the thalamic or CR MRS region). The TSI variable was included due to the specific interest of whether metabolite levels are restored as a factor of time. However, Akaike Information Criteria (AIC) analyses rejected the age variable due to lack of importance in all models. Age showed no effect on any metabolites in either examined region (tCho_Thal: β = 0.116, t = 0.755, p = 0.455; tCr_Thal: β = −0.041, t = −0.266, p = 0.791; tNAA_Thal: β = −0.064, t = −0.417, p = 0.679; tCho_CR: β = 0.068, t = 0.434, p = 0.667; tCr_CR: β = 0.144, t = 0.931, p = 0.357; tNAA_CR: β = −0.048, t = −0.304, p = 0.762). Within the thalamic MRS voxel, regression models were significant for all three metabolites, with shorter CVR delay being predictive of higher tCho, tCr, and tNAA concentrations. In this region, the TSI term was not predictive of any metabolites’ levels. Neither CVR delay nor TSI were predictive of any metabolite concentrations within the CR voxel. Further, CVR magnitude was not associated with any metabolite levels in either MRS voxel. These results are summarized in Table 1 and Figure 3.

Table 1.

Summary of the regression analyses examining the effects of TSI and CVR magnitude and delay on tCho, tCr, and tNAA concentrations within the thalamic and CR MRS voxels. The overall regression models were significant only with CVR delay and TSI as predictor terms and thalamic neurometabolite levels as the dependent variables. In these models, the TSI had no significant associations with any metabolite concentrations, while CVR delay was significantly predictive of tCho, tCr, and tNAA concentrations.

Region	Metabolite	Multivariate			Univariate
		R 2 Adj	F (df)	p	Term	β	t	p
Delay of CVR
Thalamus	tCho	0.130	4.201 (2, 41)	0.022	TSI	0.165	1.151	0.256
					CVR delay	−0.355	−2.471	0.018
	tCr	0.152	4.867 (2, 41)	0.013	TSI	0.052	0.366	0.716
					CVR delay	−0.428	−3.015	0.004
	tNAA	0.123	4.018 (2, 41)	0.026	TSI	0.090	0.621	0.538
					CVR delay	−0.382	−2.650	0.011
Corona radiata	tCho	−0.050	0.002 (2,40)	0.998	TSI	0.008	0.049	0.961
					CVR delay	−0.001	−0.007	0.995
	tCr	0.015	1.328 (2, 40)	0.276	TSI	−0.121	−0.738	0.465
					CVR delay	−0.265	−1.620	0.113
	tNAA	−0.047	0.050 (2, 40)	0.951	TSI	−0.049	−0.291	0.773
					CVR delay	−0.037	−0.219	0.828
Magnitude of CVR
Thalamus	tCho	0.014	1.295 (2,41)	0.285	TSI	0.175	1.091	0.282
					CVR magnitude	0.121	0.750	0.458
	tCr	0.062	2.421 (2,41)	0.101	TSI	0.004	0.027	0.979
					CVR magnitude	0.324	2.064	0.045
	tNAA	−0.010	0.793 (2, 41)	0.459	TSI	0.098	0.603	0.550
					CVR magnitude	0.137	0.841	0.405
Corona radiata	tCho	0.008	1.176 (2, 40)	0.319	TSI	−0.014	−0.092	0.927
					CVR magnitude	−0.237	−1.532	0.133
	tCr	−0.048	0.039 (2, 40)	0.961	TSI	−0.024	−0.150	0.881
					CVR magnitude	0.035	0.222	0.826
	tNAA	−0.021	0.569 (2, 40)	0.571	TSI	−0.052	−0.331	0.743
					CVR magnitude	−0.163	−1.041	0.304

CVR: cerebrovascular reactivity; TSI: time since infection.

Figure 3.

Shorter CVR delay was significantly related to concentrations of tCho, tCr, and tNAA in the thalamic MRS voxel. However, no significant effects of CVR delay on any neurometabolite concentrations were observed within the CR MRS voxel.

Discussion

Our results demonstrated a significant decrease in CVR delay as a factor of time since SARS-CoV-2 infection, suggesting that brain vasculature may be functionally impaired in the early post-infection phase, recovering spontaneously over the following months and years. Importantly, this was observed in individuals who had mild COVID-19, highlighting that cerebrovascular dysregulation may not only follow serious acute-stage infections. These findings bolster the evidence for the involvement of cerebrovascular impairments in neurocognitive symptoms of long COVID.^9,68–70 It remains of interest whether cerebrovascular dysfunction is driving long COVID symptomology alone, or if it emerges due to a combination of cerebrovascular deficits with SARS-CoV-2-induced damage in other organs. Future studies should evaluate whether a post-infection cooccurrence of gut epithelial and brain endothelial damage facilitates unchecked proliferation of gut-derived compounds systemically, ultimately accessing the brain.

CVR analyses are often framed in the context of vascular dysfunction having effects on neuronal health, yet few studies have evaluated this link directly.^71–75 We showed that GM CVR delay was associated with tCho, tCr, and tNAA levels in the thalamic region, which served as a proxy for global GM. Unlike CVR delay itself, the concentrations of these metabolites did not exhibit spontaneous normalization as a function of post-infection time, but rather appeared to depend on vascular recovery. While all examined metabolites are markers of neuronal health, they differ in the mechanisms through which they rely on cerebrovascular function. In the context of COVID-19, the relationship between longer CVR delay and lower NAA may reflect diminished mitochondrial function.^76–79 Constrained ATP synthesis is coupled with oxygen-dependent dysregulation of its buffering by Creatine, undermining its distribution to meet neuronal energy demands.^80,81 Vascular choline undersupply jeopardizes healthy cell membrane turnover, which is exacerbated by the aforementioned under-production of ATP and distribution. ⁸² Of note are existing reports of increased post-COVID-19 choline concentrations—findings attributed to broad dysregulation of cholinergic brain metabolism.^83,84 In the context of chronic post-COVID conditions, the observed neurometabolic pattern underscores the role of vasculature in driving neuronal tissue recovery and potential amelioration of neurocognitive symptoms.

Intermittent hypercapnia has been demonstrated as an effective method for enhancing glymphatic action—the mechanism that clears unwanted waste products from the brain.^85,86 Similar to the CVR on/off block paradigm described herein, the resulting oscillatory cycles of cerebrovascular dilation and constriction appear to drive influx of CSF into and egress out of the brain, carrying potentially neurotoxic metabolic byproducts out of the interstitial space. Glymphatic dysfunction has been discussed extensively in the context of neurodegenerative conditions that involve aggregation of unwanted proteins and peptides in the brain.^87,88 Given the emerging evidence that cerebral vasculature may exhibit functional impairments following COVID-19, it is of interest to leverage our design to evaluate whether the glymphatic system is also affected. Importantly, such endeavors may point to better mechanistic understanding of the recent links made between long COVID and higher risk of dementia.^89–91 Specifically, diminished post-COVID glymphatic clearance may exacerbate the accumulation of dementia-related proteins and peptides in at-risk groups. Relatedly, glymphatic dysfunction would also result in insufficient clearance of neurotoxic compounds that originate outside of the brain. ⁹² For example, neurotoxic substances produced by the gut microbiome may proliferate systemically and access the brain via impaired intestinal and blood–brain barriers, accumulating and driving long COVID pathology.

While it remains of interest to contrast post-infection vascular dynamics between long COVID patients and those who experienced no long-term symptoms after infection, the main objective of the described analysis is to characterize chronic physiological CVR trajectories in a sufficiently powered sample without granular parcellation of post-infection symptomology. Given the downward trend in CVR delay across time and its associated improvements in neurometabolite levels, the intuitive hypothesis would be that long COVID patients experience longer CVR delays to exogenous vascular stimuli than healthy individuals. Indeed, the CVR literature overwhelmingly points to longer CVR delays being associated with a multitude of negative outcome measures. ⁹³ However, we do not dismiss the possibility of the opposite scenario in long COVD patients, whereas it is possible for them to exhibit shorter CVR delays reflective of hyper-reactive vasculature. Patients with PASC commonly report symptoms consistent with postural tachycardia syndrome (POTS)—a condition with subtypes that can involve reduced peripheral vascular resistance, excessive vasodilation, and reduced stress response.^94–96 While several studies have reported CBF autoregulation impairments in POTS based solely on metrics that are similar to CVR magnitude, the timing of vascular responses has not been explicitly evaluated.^97–99 As such, while the described data offers insight into the long-lasting cerebrovascular sequelae of SARS-CoV-2, their extrapolation to the long COVID patient population warrants caution, specifically with regard to faster CVR delay being necessarily indicative of healthy vasculature.

Given that the majority of the population has had COVID-19 at some point, study recruitment of individuals who have never been infected is challenging. ¹⁰⁰ In the absence of such a control group, using a known infection date as the anchor relative to which the CVR metrics were examined is the next best alternative for drawing inferences about any potential COVID-induced cerebrovascular damage in the early post-infection stages. As such, these challenges limit our ability to empirically categorize early post-infection CVR delays as aberrant or impaired. However, the highly significant downward trends in CVR delay in both GM and WM across time offers compelling evidence that this may be the case, and that brain vascular function exhibits spontaneous recovery over a protracted timeline. The observed trends were also quite consistent, with no significant heteroskedasticity observed in the TSI × CVR delay relationship in either the GM or the WM, as revealed by Breusch Pagan tests (GM: F (2, 66) = 0.862, p = 0.427; WM: F (2, 66) = 1.427, p = 0.247). Despite this, we note that the described cross-sectional analysis is inherently vulnerable to confounds by pre-existing conditions that may have been overlooked by the stringent participant exclusion criteria implemented in our study, such as individual differences in premorbid vascular health, inflammatory burden, or other lifestyle factors.

Other studies utilizing the RespirAct system in similar fMRI paradigms offer some insight into normal CVR response values; however, most were conducted after 2020 and did not characterize prior SARS-CoV-2 infections in their participants. One 2019 study using RespirAct reports a mean time-to-peak GM CVR delay of 6–9 s in a group of 19 healthy young adults (20 ± 2 years). ¹⁰¹ A 2015 transfer function analysis of the BOLD response to CO₂ utilized a RespirAct-controlled CO₂ administration in healthy participants, and while CVR delays in absolute unit(s) were not presented, the negative phase angles of ~−0.3 to −0.6 radians, as inferred from figure colormaps, yield a CO₂–BOLD delay of ~5–10 s at a frequency of 0.01 Hz. ¹⁰² These values correspond to the shortest CVR delays observed in our data. Collectively, given the variability in equipment setups, study samples, and analysis methods, it is difficult to derive an accurate normative range for healthy CVR delay from the literature.

Single-voxel spectroscopy necessitates the selection of a brain region to be examined, limiting the inferences that can be made about the whole brain. In this case, we conducted MRS in two spatially separated voxels with the goal of them being proxies for global GM (thalamic voxel) and global WM (CR voxel).The thalamic MRS voxel was dominated by GM (GM = 73.5% ± 4.2%, WM = 26.3% ± 4.3%, CSF = 0.2 ± 0.4), while the CR voxel consisted mostly of WM (GM = 16.4% ± 7.5%, WM = 82.7 ± 8.0, CSF = 0.9% ± 0.9%). Therefore, our selection of the thalamic and CR MRS voxel placements to differentially capture GM- and WM-specific processes was largely successful, considering the practical difficulties of positioning a cubic MRS voxel within the brain to capture any one specific tissue type. To evaluate the extent to which each voxel reflected its target tissue of interest (thalamic = GM, CR = WM), we examined whether global GM CVR delay is differentially predictive of tCho, tCr, and tNAA levels within the Thalamic MRS voxel relative to the CR MRS region. To counterbalance this, we conducted analogous analyses using global WM CVR delay as the predictor variable. Indeed, TSI-corrected global GM delay was significantly predictive of metabolite concentrations in the thalamic MRS voxel, but not in the CR (tCho_Thal: β = −0.308, t = −2.090, p = 0.043; tCr_Thal: β = −0.435, t = −3.053, p = 0.004; tNAA_Thal: β = −0.411, t = −2.869, p = 0.006; tCho_CR: β = 0.055, t = 0.329, p = 0.744; tCr_CR: β = −0.155, t = −0.930, p = 0.358; tNAA_CR: β = 0.020, t = 0.120, p = 0.905). We found that global WM CVR delay was also significantly predictive of tCr and tNAA concentrations in the thalamic voxel, but was not predictive of any metabolite levels in the CR voxel (tCho_Thal: β = −0.278, t = −1.857, p = 0.070; tCr_Thal: β = −0.405, t = −2.780, p = 0.008; tNAA_Thal: β = −0.376, t = −2.566, p = 0.014; tCho_CR: β = 0.019, t = 0.110, p = 0.913; tCr_CR: β = −0.188, t = −1.125, p = 0.267; tNAA_CR: β = −0.055, t = −0.325, p = 0.747). Collectively, our interpretation of these findings is that the thalamic region may serve as a reasonable proxy for assessing metabolic effects of global cerebrovascular dynamics in both white and gray matter. Nevertheless, concrete generalizability of deep-GM thalamic neurometabolite findings to global cortical GM was not assessed, and should be interpreted with caution. The CR region appeared to offer little utility for capturing neurometabolite changes in the context of global cerebrovascular recovery. This discrepancy may reflect diminished vascularization and BOLD SNR in the WM relative to the GM, with downstream effects of WM vascular recovery subsequently being slow and limited in their detectability.

The variability in disease severities caused by different SARS-CoV-2 strains is well established, with Omicron producing notably milder acute-stage infections than its predecessors. ¹⁰³ While the described study did not directly characterize the specific viral strains that its participants were infected with, we used infection dates to estimate whether each participant has been infected with a pre- or post-Omicron SARS-CoV-2 variant. We conducted several post-hoc analyses to evaluate the potential effects of viral variant severities on the observed CVR and MRS measures. With the emergence of Omicron being a pivotal point in the pandemic that was followed by progressively diminishing overall acute-stage COVID-19 severities, we used the estimated date of January 1, 2022 as the cutoff point after which Omicron became the dominant strain in the United States. ¹⁰⁴ Using this date, participants were classified as either having pre- or post-Omicron infections. As expected, the majority of the participants were in the latter group (pre-Omicron N = 19, post-Omicron N = 50). Nevertheless, to explore whether this broad grouping of strain severities had any potential to confound the results of our study, we used a series of ANCOVA models to examine TSI- and age-corrected effects of pre- versus post-Omicron grouping on the primary outcome measures of CVR (magnitude, delay) and MRS (thalamic tCho, tCr, tNAA; CR tCho, tCr, tNAA). The results indicated no significant differences between the pre- and post-Omicron groups on any measures of CVR (GM magnitude: F (1, 65) = 0.043, p = 0.837; GM delay F = 0.112, p = 0.739; WM magnitude: F = 0.603, p = 0.440; WM delay F = 0.196, p = 0.660) or MRS (tCho_Thal: F = 3.374, p = 0.131; tCr_Thal: F = 1.694, p = 0.201; tNAA_Thal: F = 1.785, p = 0.189; tCho_CR: F = 0.149, p = 0.702; tCr_CR: F = 0.207, p = 0.652; tNAA_CR: F = 0.275, p = 0.603).

The described study does not include neurocognitive sequelae of the suggested cerebrovascular recovery as a function of post-infection time. Quantification of such relationships is key to characterizing the role of cerebrovascular dysfunction in long COVID pathogenesis. Our group is actively working on the disentanglement of complex multi-faceted symptomology profile of PASC in the context of the described cerebrovascular patterns, yet these associations remain beyond the scope of this report. Specifically, our goal here was to determine whether a foundational long-term post-infection pattern exists in the dynamics of cerebrovascular function. To this end, the study described a compelling time-dependent trend of spontaneous, albeit protracted, recovery of brain vasculature following COVID-19. These findings offer novel and useful information to this field of work, highlighting that CVR dynamics may be altered for months or years after infection, while also offering a degree of reprieve to PASC patients, as the described cerebrovascular aberrations ultimately do appear to resolve on their own with time. In summary, future work on linking these findings with specific symptoms as well as interventions that can accelerate the process of vascular recovery are worthwhile pursuits for alleviating PASC symptoms.

We demonstrated in individuals with histories of mild COVID-19 that the ability of brain vasculature to respond to exogenous vasoactive stimuli is delayed in the early relative to late post-infection timeframes. Diminishing CVR delay over time suggests spontaneous cerebrovascular recovery and that endothelial damage inflicted by SARS-CoV-2 may be reversible. While it remains unknown if this is the case in PASC, our findings potentially implicate endothelial dysfunction in this condition. Importantly, we demonstrated that recovery of brain vasculature corresponds to the restoration of neuronal health, bolstering the feasibility of using vascular interventions to alleviate neurocognitive symptoms in PASC patients.