Nicotinamide adenine dinucleotide (NAD) is a crucial metabolite for all living cells, existing in the oxidized (NAD⁺) and/or reduced (NADH) forms. It functions as a coenzyme in various cellular redox reactions and serves as a substrate for various NAD⁺-consuming enzymes (Belenky et al., 2007; Verdin, 2015). NAD has been shown to play pivotal roles in numerous cellular processes and functions, including energy metabolism, mitochondrial function, cellular signaling, calcium homeostasis, aging, longevity, and cell death (Ziegler, 2000; Belenky et al., 2007; Verdin, 2015; Cuenoud et al., 2020). Therefore, the in vivo measurement and quantification of intracellular NAD levels in the human brain have attracted significant interest in numerous scientific investigations on brain aging (Braidy et al., 2019; Lautrup et al., 2019; Griffiths et al., 2020), cognitive functions (Zhao et al., 2021; Campbell, 2022), and neurodegenerative diseases such as Alzheimer’s (Ruan et al., 2018; Hou et al., 2021; Wang et al., 2021), Parkinson’s (Pérez et al., 2021; Kriebs, 2022; Mischley et al., 2023), and Huntington diseases (Lautrup et al., 2019; Lundt and Ding, 2021).

Despite the profound importance of intracellular NAD metabolism in human health and disease, assessing brain intracellular NAD levels in-situ has been challenging. Recent advancements in in vivo magnetic resonance (MR) spectroscopy/spectroscopic imaging (MRS/MRSI)-based techniques have enabled non-invasive measurement of NAD levels in the human brain (Lu et al., 2014b; Zhu et al., 2015; Lu et al., 2016; de Graaf et al., 2017; Bagga et al., 2020; Skupienski et al., 2020; Dziadosz et al., 2022). However, most of the current studies were restricted to the examination of a small brain region or using averaged spectra for a few brain regions, with limited exploration of NAD distributions across various brain regions over the entire human brain (Lu et al., 2014b; Zhu et al., 2015; Das et al., 2021; Ren et al., 2023). These limitations stem from several long-standing technical hurdles. The intracellular NAD concentrations in the brain reside within the sub-millimolar range, yielding an exceptionally low signal-to-noise ratio (SNR; Lu et al., 2014b; Zhu et al., 2015). Furthermore, as a high-dimensional imaging method, MRSI typically demands a large number of measurements to encode both spatial and spectral information. Therefore, using existing MRSI methods, especially phosphorus-31 MRSI (³¹P-MRSI) to obtain high-quality maps of NAD concentrations throughout the entire brain could take prohibitively long scan times on the order of hours (Santos-Díaz and Noseworthy, 2020).

In recent years, a range of techniques has been developed to tackle these challenges and advance ³¹P-MRSI applications for mapping brain metabolites. The emergence of ultrahigh field MR systems has yielded significantly enhanced SNR and spectral resolution (Lei et al., 2003; Moser et al., 2012; Barker, 2014; Ruhm et al., 2021). For instance, a 2.8-fold sensitivity increase has been shown in ³¹P-MRS at 7 T compared with 3 T (Rodgers et al., 2014; Lu et al., 2014a). Specialized acquisition methods such as proton decoupling and nuclear overhauser effect (NOE) have been utilized to further improve spectral quality and SNR (Luyten et al., 1989; Lei et al., 2003); and fast scanning trajectories like echo-planar or spiral trajectories have been used to enhance data acquisition speed (Valkovič et al., 2016; Korzowski and Bachert, 2018). Advanced post-processing methods like low-rank denoising have further contributed to noise reduction (Ma et al., 2017; Clifford et al., 2020; Korzowski et al., 2020). Nevertheless, quantitative imaging of the intracellular NAD level throughout the entire human brain still remains challenging given its exceptionally low concentration.

The purpose of this study is to demonstrate the feasibility of whole-brain NAD mapping by synergistically integrating ultrahigh field ³¹P-MRSI with an advanced denoising method based on the probabilistic subspace model (Chen et al., 2020; Li et al., 2021a; Zhang et al., 2023). This method leverages ultrahigh field, spectral priors, spatial constraints, and statistical priors to significantly enhance the sensitivity of ³¹P-MRSI for NAD imaging. Its performance in terms of SNR, accuracy, and reproducibility of the resulting NAD measurements were evaluated. The potential of mapping tissue metabolic heterogeneity and revealing age dependence of NAD levels were also preliminarily assessed. The results indicate that this method holds a significant promise for mapping the NAD content in the human brain.

Figure 1 shows a comparison of resulting ³¹P-MRSI spectra, SNR maps, and metabolite maps (on one set of “2.3-cc data”) using different processing methods to illustrate the denoising efficacy of the presented method. On the representative spectra from a selected voxel, we can see many metabolite peaks other than PCr and ATPs were buried under the noise level in the raw ³¹P-MRSI data; the basic low-rank denoising method performed well in noise reduction but the NAD peak was still difficult to be distinguished from noise; while the results from the proposed method achieved a better noise reduction and the NAD peak became clearly above the noise level. The SNR maps also affirmed the significant improvement in SNR compared with original data and basic low-rank denoising. The mean SNRs over the brain were 12.65 ± 1.91 dB, 17.05 ± 1.67 dB, and 22.49 ± 1.48 dB for raw MRSI data, low-rank denoising, and the proposed method, respectively. Compared with the original noisy data, the presented method almost provided an around 10-fold sensitivity enhancement, which might make it possible for reliable detection of NAD signals in the voxel base. As shown in Figure 1, given the relatively high concentration of PCr, its spatial intensity maps were of good quality and consistent using different methods. But the spatial maps of low-concentration metabolites like GPC and NAD were heavily contaminated by noise without denoising. The NAD maps using basic low-rank denoising showed significant noise reduction but there were still noticeable spatial variations induced by noise. In contrast, the results produced using the proposed method showed further reduced noise fluctuations and a relatively homogenous distribution of NAD within the human brain. The signal intensities of NAD over the brain were 0.52 ± 0.18, 0.49 ± 0.11, and 0.42 ± 0.08 (in institutional unit) for these three methods, respectively. As expected, the proposed method produced the smallest spatial variations.

One complete set of whole-brain metabolite intensity maps (4 representative ³¹P-MRSI slices from another set of “2.3-cc data”) is shown in Figure 2, including PCr, αATP, γATP, βATP, Pi, GPE, GPC, PE, PC, and NAD. Most of the intensity maps, including NAD, were of high quality with very minor effects of noise, except PC, which had the lowest signal intensity in the acquired ³¹P-MRSI data (as shown in the spectra in Figure 1). The NAD map suggested that intracellular NAD concentrations appeared to have a relatively homogeneous distribution over the brain, with only slight enhancements in midbrain or thalamic regions. It is worth noting that PCr had high signals in some subcutaneous regions, but other metabolites did not, which was suspected due to the scalp muscles with a higher PCr concentration.

Figure 3 shows results of the computational simulation to illustrate accuracy of the proposed method. Representative spectra and NAD estimates using different processing methods were compared. Similar to the observations in Figure 1, raw measurements without denoising had very noisy estimates of NAD intensity, the low-rank denoising provided a reduced noise level, while the proposed method offered the best SNR. With a significant noise reduction, the NAD intensity map and localized spectrum produced using the proposed method were very close to the ground truth. Quantitatively, rRMSEs of NAD estimates referring to the ground truth were 48.7, 20.5, and 12.4% for the raw data, low-rank denoising, and the proposed method, respectively. In this realization of noise simulation, SNR enhancement of the proposed method led to an improved accuracy on NAD estimation by around six-fold, and the resulting relative errors were around 12%.

Figure 4A displays the quantified NAD concentration maps (in the unit as mM) from all seven subjects in the “2.3-cc data.” These NAD concentrations are displayed in tripolar views and show similar spatial distributions in the brain. A summary of the NAD concentrations across these subjects is shown in Figure 4B. The mean NAD concentration of each subject was around 0.4 mM, and the standard deviation was around 0.1 mM. These statistics were reasonably consistent among these seven subjects studied to date. Figure 4 also includes the reproducibility results from two repeated scans on the same subject. One pair of spectra in the same spatial voxel from these two measurements were displayed and compared in Figure 4C. A good agreement was observed between these two spectra, including the small NAD peaks, as indicated in the zoom-in regions between −5 and −10 ppm. The R2-plot in Figure 4D includes quantitative NAD concentrations of all voxels within the brain. Compared with the identical line, there was no significant bias between these two NAD measurements, and the NAD concentrations were mostly distributed from 0.30 mM to 0.50 mM. Additionally, the Pearson’s correlation coefficient between these two measurements was 0.72 and the averaged coefficient of variation was 6.1%, further suggesting an excellent reproducibility between these two scans, even though the NAD concentration is very low.

Figure 5 presents the SNR analysis on one high-resolution “1.0-cc data” using different processing methods. Compared with the “2.3-cc data,” this set of data was acquired with higher resolution and a shorter scan time, thus intrinsically having a lower SNR. The SNR differences with the “2.3-cc data” could be noticed from the representative spectra and SNR maps compared with Figure 1. But similarly, the proposed method still outperformed the raw MRSI data and basic low-rank denoising method and produced spectra of a largely improved quality, with the NAD resonance visible above the noise level. The calculated SNRs over the brain were 9.65 ± 3.06 dB and 18.22 ± 2.61 dB before and after denoising, approximately one half and one third of the “2.3-cc data,” respectively.

Figure 6 shows the high-resolution metabolite maps obtained from one set of “1.0-cc data,” including intensity maps of PCr, αATP, GPC, and NAD. With the degraded SNR compared with the “2.3-cc data,” the proposed method still produced reasonable although relatively noisier NAD estimates. Nevertheless, benefitting from the higher resolution, the heterogeneity of metabolite distributions between brain tissues became more visible. More specifically, PCr exhibited a clear higher concentration in gray matter than white matter, while the ATPs and NAD showed relatively uniform distributions between brain tissues. The regression analysis results of this data are shown in Figure 7. From the plots and regression curves, we found that higher gray matter fractions were associated with higher PCr levels, lower GPC levels, uniform ATP distributions, and similar NAD levels, which matched the observations in literature reports (Ruhm et al., 2021).

The age dependence of NAD level was revealed using the “1.0-cc data.” The NAD levels, which were signal intensities of NAD over the brain, of the 14 subjects and their ages were displayed in Figure 8. As we can see, cluster of the data with ages below 30 showed the highest overall NAD levels, and the NAD levels declined with ages above 30. The linear regression curve also confirmed this trend of NAD decline as age increases. The whole brain averaged spectra before and after denoising of these subjects were displayed in Supplementary Figure S5. This NAD decline with aging matched the biological expectation (Zhu et al., 2015), and could serve as strong indirect evidence to further support the feasibility of our ³¹P-MRSI method in mapping brain NAD content.

Using ³¹P-MRSI techniques, the feasibility for mapping phosphorus metabolite distributions over the human brain with reasonable sensitivity and resolution has been demonstrated previously (Korzowski and Bachert, 2018; Korzowski et al., 2020, 2021; Rowland et al., 2020; Ruhm et al., 2021). However, for low-concentration metabolite signals such as NAD, it is still very challenging to measure and quantify due to the limited SNR in individual voxels. In this study, we applied a probabilistic subspace-based approach to process the ³¹P-MRSI data acquired at ultrahigh field and achieved significant noise reduction, which made it possible to reliably map and quantify the low-level NAD content across the entire human brain at 7 Tesla. Computational simulation and in vivo experiments were carried out to evaluate the SNR performance, accuracy, and reproducibility of the imaging method in NAD measurement. From the in vivo results, brain NAD concentrations around 0.4 mM were consistently found among multiple subjects. This averaged NAD level over the brain is in good agreement with the reported values obtained from human occipital and other brain regions (Zhu et al., 2015; Ruhm et al., 2021).

This feasibility study showed a relatively homogenous distribution of NAD levels over the brain. On the one hand, in the “2.3-cc data,” the spatial resolution is limited. Even though no hamming window was applied in processing to preserve spatial resolutions, partial volume effects could still partly contribute to the homogenous distribution we observed. On the other hand, the contrast analysis on the “1.0-cc data” also implied an insignificant NAD contrast between gray matter and white matter, especially compared with PCr and GPC. Even though only a small cohort was included in this analysis and further studies will be needed to confirm, the experimental evidence on the current data we have up to now consistently suggested such a homogenous distribution of brain NAD concentrations. This homogenous distribution of NAD might indicate a potentially uniform need for maintaining NAD homeostasis across the brain as ATP does, and might suggest that NAD plays a consistent role in supporting cellular functions throughout the entire brain.

This study showed a noticeable decline of brain NAD levels with aging (see Figure 8). This trend aligned with previous ex vivo studies and in vivo reports focusing on occipital and frontal lobes (Gomes et al., 2013; Zhu et al., 2015; Kim et al., 2017), suggesting a decline in mitochondrial metabolism efficiency and functionality during aging process. Successful reveal on the age dependence of NAD supported the feasibility of our method in measuring NAD changes under different biological conditions, which may have many potential applications in neurodegenerative diseases and aging-related disorders. But one limitation of this set of data is the limited size of cohort. Only 14 data were included for now and only three of them have an age above 30. In the future, including more subjects, especially more elderly volunteers, will strengthen the observations and evidence. In addition, more data in a wider age range may make it possible to understand better the age dependency of NAD concentration, even its variations in different age stages. It should also be noted that only the signal intensities after spectral fitting instead of quantified NAD concentrations were analyzed for the “1.0-cc data.” This is because NOE was employed in this dataset, which made the absolute quantification of metabolites challenging. In order to achieve absolute quantification, additional calibration scans and systematic studies will be needed, which is beyond the scope of this paper.

The partial separability model-based method has been widely used for noise reduction in not only ³¹P-MRSI, but also MRSI of other nuclei including ¹H, ¹³C, ²H, etc. (Liang, 2007; Nguyen et al., 2013; Lam and Liang, 2014; Li et al., 2017; Guo et al., 2019). Compared with the basic low-rank approximation method used in previous ³¹P-MRSI studies (Korzowski et al., 2020; Ruhm et al., 2021), the probabilistic subspace method used in this work has several key features for enhanced noise reduction. First, the presented method pre-determined the temporal basis functions from the group data, in contrast to the low-rank approximation method which only depends on a single noisy data. One set of representative spectra of these basis functions was displayed in Supplementary Figure S4. This strategy leverages signal correlations across all the collected data to improve the SNR for subspace estimation, thus improving accuracy. Second, the presented method imposed two additional constraints on the spatial coefficients, including a total variation regularization and a statistical distribution constraint. The spatial regularization imposes anatomically weighted spatial smoothness on the coefficients, which has been widely used in various imaging reconstruction applications (Haldar et al., 2008; Lam and Liang, 2014; Zhao et al., 2015; Liu et al., 2018; Wang et al., 2020). The statistical distributions absorb prior information from the group data, which provide soft boundary constraints on the spatial coefficients to prevent them from going outside of the feasible distributions due to noise. These statistical constraints have been applied and proved effective for reliable estimations of myelin water components and B₁ (radiofrequency magnetic field) mapping (Li et al., 2021a; Zhang et al., 2023).

The noise reduction provided by the probabilistic subspace model can be utilized in many applications of ³¹P-MRSI techniques. On the one hand, the SNR enhancement can benefit the imaging resolution and/or speed, especially for those metabolites with relatively high concentrations (e.g., >1 mM). The “1.0-cc data” in this study has shown the potential in increasing spatial resolution, even for low-concentration metabolites like NAD. For applications where high-concentration metabolites are of interest or high resolutions are not required, the scan time can be significantly reduced, possibly to the level close to typical MR imaging scans. On the other hand, the SNR benefits can be applied to measure metabolites with even lower concentrations. For example, a spectrum averaging 9 adjacent voxels in the “2.3-cc data” is displayed in Figure 9. With the improved SNR, the doublet of NAD⁺ (which should be a doublet of doublet, but it is only visible as doublet in in vivo data due to wide linewidth) plus the singlet of NADH result in an apparent triplet, interestingly, similar with the ³¹P MRS acquired in the human occipital lobe with ¹H-decoupling at 4 Tesla (Lu et al., 2014b). This result indicates great promise to separate these two metabolites on the voxel base, thus, to estimate and imaging the NAD redox ratio (=[NAD⁺]/[NADH]) in human brain (Lu et al., 2014b). Moreover, the Uridine Diphosphate Glucose (UDPG) resonance peaks at around −9.8 ppm were also above the noise level, whose concentrations were reported as around 0.1 mM (Ren et al., 2015). These observations demonstrated the potential of mapping the [NAD⁺]/[NADH] ratio and/or [UDPG] in the entire brain, possibly with a relatively lower resolution. This is one important direction of our future work.

The computation simulation study in Figure 3 used a group averaged ³¹P-MRSI data as the ground truth. This provided a good reference with realistic spectral lineshapes and allowed for quantitative evaluation on the performance of different processing methods. But using group-averaged spectra could potentially average out subtle variations of individual data. Therefore, we performed an additional comparison using the individual data. Specifically, one of the raw “2.3-cc data” was used, and every 24 voxels of the data were averaged to generate a large-voxel data with a relatively good SNR in each voxel. The proposed denoising method was applied on the raw small-voxel data and the same voxel-averaging was performed to generate the large-voxel data for comparison. The NAD estimates and localized spectra from these two data are shown in Figure 10. The results were reasonably consistent between the large-voxel raw data and denoised data, and the rRMSE on NAD estimates was 5.7%. This also supports the feasibility of the denoising method.

There are several potential improvements that can be made to achieve a more accurate quantification of NAD concentration. First, the excitation RF profile has not been corrected in the NAD quantification due to its small impact. For the “2.3 cc” data, the excitation hard pulse was centered in the middle of the αATP and γATP peaks. Therefore, the correcting factors of excitation profile for γATP and NAD are 1.038 and 1.071, respectively. Correcting this factor could slightly improve accuracy of NAD quantification. Second, the current spectral quantification model assumed a singlet resonance for the NAD peak, while the NAD peak should include NAD⁺, NADH, and several UDPG components (Ren et al., 2020; Ren and Sherry, 2021). Separating these components requires exceptional SNR and more complicated signal modeling. Being able to separate these components will further improve the accuracy of NAD quantification.

The technical advancements in MR hardware and systems could further improve the sensitivity of such ³¹P-MRSI methods. The ultrahigh field systems used in this study were 7 T scanners, while there are several MR human systems with higher field strength, such as 9.4 and 10.5 T. Higher field strength could further boost the SNR, thus enhancing the imaging capability of ³¹P-MRSI techniques. We used volume ³¹P RF coils in this study. It is well known that ³¹P RF array coil with multiple receive channels could provide better SNR than the volume coil (van Uden et al., 2019; Rowland et al., 2020; Peeters et al., 2021). Combined with a higher field strength, better RF coils, and the advanced processing methods used in this work, the ³¹P-MRSI technology has the potential to provide sufficient resolution and speed for practical applications in research and clinical settings.

The presented capability of non-invasively mapping NAD levels in entire human brain is expected to be highly valuable and have broad implications in many biomedical and clinical research fields. It will also allow evaluation of the alteration of NAD levels in diseased brain regions due to neurological disorders, neurodegenerative diseases, or age-related cognitive decline. Mapping NAD contents in these conditions may uncover specific patterns of regional NAD-related metabolic dysregulation associated with disease progression, potentially opening avenues for targeted interventions and therapeutic strategies. In addition to its impact on disease research, this technique could serve as a powerful imaging tool for neurodevelopmental studies, which allows for the exploration of how brain NAD dynamics contribute to the maturation and developmental progression of the brain.

In conclusion, by combining the ultrahigh field 3D ³¹P-MRSI with advanced probabilistic subspace-model based post-processing method, we demonstrated for the first time the feasibility of whole-brain NAD mapping in healthy human subjects at 7 T. This approach would allow non-invasively monitoring cerebral NAD and its changes in brain regions under various brain conditions.

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

The studies involving humans were approved by the Institutional Review Board of the University of Minnesota; Institutional Review Board at the University of Pittsburgh. The studies were conducted in accordance with the local legislation and institutional requirements. The participants provided their written informed consent to participate in this study.

RG: Formal analysis, Investigation, Methodology, Writing – original draft, Writing – review & editing. SY: Data curation, Formal analysis, Writing – review & editing. HW: Data curation, Investigation, Writing – review & editing. YL: Conceptualization, Formal analysis, Methodology, Writing – review & editing. YZ: Data curation, Methodology, Writing – review & editing. Z-PL: Conceptualization, Investigation, Resources, Supervision, Writing – review & editing. WC: Conceptualization, Formal analysis, Project administration, Resources, Supervision, Validation, Writing – review & editing. X-HZ: Conceptualization, Data curation, Funding acquisition, Resources, Supervision, Validation, Writing – review & editing.