Multiparametric Renal Magnetic Resonance Imaging: Validation, Interventions, and Alterations in Chronic Kidney Disease

Background: This paper outlines a multiparametric renal MRI acquisition and analysis protocol to allow non-invasive assessment of hemodynamics (renal artery blood flow and perfusion), oxygenation (BOLD T2*), and microstructure (diffusion, T1 mapping). Methods: We use our multiparametric renal MRI protocol to provide (1) a comprehensive set of MRI parameters [renal artery and vein blood flow, perfusion, T1, T2*, diffusion (ADC, D, D*, fp), and total kidney volume] in a large cohort of healthy participants (127 participants with mean age of 41 ± 19 years) and show the MR field strength (1.5 T vs. 3 T) dependence of T1 and T2* relaxation times; (2) the repeatability of multiparametric MRI measures in 11 healthy participants; (3) changes in MRI measures in response to hypercapnic and hyperoxic modulations in six healthy participants; and (4) pilot data showing the application of the multiparametric protocol in 11 patients with Chronic Kidney Disease (CKD). Results: Baseline measures were in-line with literature values, and as expected, T1-values were longer at 3 T compared with 1.5 T, with increased T1 corticomedullary differentiation at 3 T. Conversely, T2* was longer at 1.5 T. Inter-scan coefficients of variation (CoVs) of T1 mapping and ADC were very good at <2.9%. Intra class correlations (ICCs) were high for cortex perfusion (0.801), cortex and medulla T1 (0.848 and 0.997 using SE-EPI), and renal artery flow (0.844). In response to hypercapnia, a decrease in cortex T2* was observed, whilst no significant effect of hyperoxia on T2* was found. In CKD patients, renal artery and vein blood flow, and renal perfusion was lower than for healthy participants. Renal cortex and medulla T1 was significantly higher in CKD patients compared to healthy participants, with corticomedullary T1 differentiation reduced in CKD patients compared to healthy participants. No significant difference was found in renal T2*. Conclusions: Multiparametric MRI is a powerful technique for the assessment of changes in structure, hemodynamics, and oxygenation in a single scan session. This protocol provides the potential to assess the pathophysiological mechanisms in various etiologies of renal disease, and to assess the efficacy of drug treatments.


INTRODUCTION
Magnetic Resonance Imaging (MRI) offers the possibility to non-invasively assess the structure of the kidney as well as renal function in a single scan session. This article outlines the development of a quantitative functional multiparametric renal MRI protocol to probe hemodynamics (total and regional blood flow, perfusion), oxygenation [Blood Oxygen Level Dependent (BOLD) T 2 * imaging], and microstructure (diffusion weighted imaging, longitudinal relaxation time T 1 mapping) and describes associated analysis methods. This multiparametric MRI protocol is applied in healthy participants, to assess both reproducibility and the field strength dependence of MRI parameters between 1.5 and 3 Tesla (T). In addition, studies are performed in healthy participants to evaluate the possibility of using hypercapnia and hyperoxia to monitor changes in renal BOLD and T 1 reactivity. Finally, pilot data demonstrating the feasibility of this multiparametric protocol in Chronic Kidney Disease (CKD) is shown.
The kidney is an intricate organ which regulates electrolytes, acid-base balance, and blood pressure and filters blood to remove water soluble waste products (Skorecki et al., 2016). Regulation of renal tissue oxygenation is complex, because renal blood flow is not only needed to prevent hypoxic injury but is inextricably linked to the requirement for glomerular filtration. The kidney's response to hypoxia cannot simply be an increase in renal blood flow, as this would also increase oxygen demand; a number of hemodynamic mechanisms are required to regulate the fine balance between oxygen delivery and consumption, and these can be measured using multiparametric MRI.
Oxygen delivery is determined by arterial blood flow, which is regulated by arterial blood pressure and intrarenal vascular resistance, local tissue perfusion and blood oxygen content (Evans et al., 2008). Arterial blood is supplied to the kidney via the renal artery, and blood flow can be estimated from phase contrast MR. The renal artery sequentially divides into segmental, interlobar, arcuate, and interlobular arteries before finally reaching the afferent arterioles that supply the glomeruli. The renal microcirculation varies depending on the location of each nephron within the cortex. In the outer cortex, glomerular efferent arterioles give rise to a capillary network that surrounds the tubules, important for reabsorption of water and electrolytes. In contrast, the efferent arterioles supply the medulla and give rise to the vasa recta, the long unbranched capillary loops that run into the inner medulla associated with the loop of Henle, as well as capillaries in the outer medulla. This facilitates concentration of urine in the medulla, but also has consequences for oxygenation.
The majority of arterial blood delivered to the kidney is directed toward the renal cortex, which primarily is responsible for filtration and tubular reabsorption; 5-15% is delivered to the medulla, whose purpose is concentration of the urine by maintaining a hypertonic environment. Despite the much lower proportional blood flow, the absolute blood flow to the medulla is still significant, reflecting large total renal blood flow (∼20% of cardiac output), and a number of physiological and pathological conditions can produce significant redistribution of renal blood flow. However, a significant cortico-medullary oxygen gradient exists, with the inner medulla having a tissue oxygenation (pO 2 ) as low as 10 mmHg, compared to 50 mmHg in the cortex. MR potentially provides a non-invasive method to assess this change in tissue oxygenation.
Tubular epithelial transport allows the kidney to regulate volume and composition of urine, but has significant energy and oxygen requirements. Reabsorption of sodium is the main determinant of oxygen consumption (Blantz et al., 2007;Thomson and Blantz, 2008), so that oxygen consumption is related to renal function, with reductions in glomerular filtration rate (GFR) resulting in a lower filtered load and lower requirement for sodium reabsorption. Renal perfusion is driven primarily by the need to maintain glomerular filtration rather than oxygenation, and therefore arteriovenous shunting of oxygen occurs to prevent hyperoxic tissue injury (via periglomerular shunts and between arterial and venous limbs of the vasa recta).
Although the kidney can reduce oxygen consumption in response to hypoxia, the lower pO 2 in the medulla increases its propensity to ischemic damage, which is considered a key pathogenic event in acute kidney injury (AKI) and CKD (Venkatachalam et al., 2010). Since multiple interacting mechanisms operate in concert to provide tight regulation of intrarenal oxygenation, dysfunction of these mechanisms may contribute to the pathogenesis of kidney disease. For example, vascular morphologic changes may occur such as, capillary rarefaction as well as factors that affect regional blood flow and oxygen diffusion e.g., upregulation of the renin-angiotensin or sympathetic nervous systems (Adler et al., 2004). Alternatively, other changes may impact regional oxygen utilization such as alterations in global and single nephron GFR, drugs interfering with glomerular hemodynamics or tubular transport and hydration status. In addition to the complexity of kidney function, renal diseases such as CKD are diverse in terms of pathophysiological processes, etiology and outcomes, highlighting the need for multiparametric MRI measures. Renal perfusion and tissue oxygenation appear central integrating factors in kidney disease, highlighting the need to perform a combined assessment of these parameters, and regardless of the nature of initial insult, fibrosis is the final common pathway.
The potential use of complementary MRI techniques to noninvasively assess multiple parameters to provide a wealth of information on renal blood flow and regional perfusion, tissue oxygenation, and degree of fibrosis, as well as behavior in low or high oxygen or carbon dioxide, will undoubtedly aid the understanding of kidney disease. Prior to the use of MRI in kidney disease, the reproducibility of MRI measures and their dependence on different factors must be understood in normal tissue.
Here, we assess the inter-subject variability, repeatability, and field strength dependence of multiparametric MRI measures in healthy participants. Physiological modulations such as, hyperoxia and hypercapnia are performed, and pilot data are shown to illustrate the feasibility of detecting changes in MR measures in CKD.

Study Design
Studies were carried out according to the principles of the Declaration of Helsinki. Healthy participant studies were approved by the Local Ethics Committee and patient studies were approved by the East Midlands Research Ethics Committee. Written informed consent was obtained from all participants.
Imaging was performed on either a 1.5 or 3 T Philips whole body MR scanner. Data is presented from studies that use the multiparametric renal MRI protocol, comprising quantification of renal blood flow and perfusion, renal oxygenation, and markers of renal microstructural change due to fibrosis/inflammation. All data was collected with subjects fasted for at least 2 h prior to their scan.

Variability, Repeatability, and Field Strength Dependence in Healthy Participants
Here, we evaluate the variation in MRI measures within normal tissue of a healthy participant cohort, specifically we assess renal artery and renal vein blood flow [as measured with phase contrast (PC)-MRI], kidney perfusion [as measured with arterial spin labeling (ASL)], T 1 measures [and a comparison of readout schemes: spin echo-echo planar imaging (SE-EPI) and balanced fast field echo (bFFE)], tissue oxygenation (from BOLD T 2 * ), diffusion weighted imaging (DWI), and total kidney volume. Data collated across a number of studies are first shown, giving a cohort of 127 participants (88 male) with mean age of 41 ± 19 years. This data is then divided into two groups comprising healthy participants <40 years and >40 years [see Section Application in Chronic Kidney Disease (CKD)]. In addition the field strength dependence of MR relaxation times is assessed. Since clinical MR scanners at both 1.5 and 3 T are now widely available, the field dependence of MR relaxation time measures of T 1 (using both SE-EPI and bFFE) and T 2 * for renal cortex and renal medulla was assessed.
A subset of 11 participants (age 20-28 years, body mass index 20-26 kg/m 2 ) had two/three repeat 3 T scans at the same time of day and after an overnight fast to limit diurnal and dietary variability. To determine the between session repeatability of MRI measures, the intra-subject Coefficient of Variation (CoV; defined as the standard deviation/mean) and intra class correlation (ICC, average measures, two-way random, absolute agreement) were assessed.

Physiological Modulation in Healthy Participants
Physiological modulations, such as gas enrichment by hypercapnia, hyperoxia, or carbogen (hypercapnic-hyperoxia; Milman et al., 2013) may provide a more sensitive marker to assess changes in renal oxygenation and microcirculation reactivity and functionality, and changes in these parameters associated with pathology. Here, we assess the change in MRI parameters in healthy participants in response to hypercapnia and hyperoxia. We induced hypercapnia and hyperoxia using a sequential gas delivery breathing circuit and a prospective, feedforward gas delivery system (Respiract TM , Thornhill Research Inc., Toronto, Canada) to control and monitor end-tidal oxygen (P ET O 2 ) and carbon dioxide (P ET CO 2 ) partial pressures. Hypercapnia was targeted at P ET CO 2 ∼6 mmHg above the subjects' baseline value whilst keeping P ET O 2 constant at the subjects' resting value, the paradigm comprised 5 min normoxia and 5 min of hypercapnia. Hyperoxia was targeted at P ET O 2 ∼500 mmHg with P ET CO 2 targeted to remain constant at the subjects' resting value. The paradigm comprised 5 min normoxia and 5 min of hyperoxia, P ET O 2 was increased/decreased over a 1 min transition period.
BOLD T 2 * was measured at 3 T using a multi-gradient echo Fast Field Echo (mFFE) sequence in six healthy participants (3 male, mean age 25 years, range 22-28 years) during the hyperoxia and hypercapnia challenge. T 1 was measured during a hyperoxic challenge in five healthy participants (3 male, mean age 26 years, range 22-31 years) using an inversion recovery sequence with modified respiratory triggering and a bFFE readout at 3 T.

Application in Chronic Kidney Disease (CKD)
To demonstrate the feasibility of use of the multiparametric MR protocol in patients, 11 patients with CKD Stage 3 or 4 were scanned (inclusion criteria: estimated GFR 15-66 ml/min/1.73 m 2 , age 18-85 years). Baseline blood pressure and estimated GFR of the patients was recorded. The complete multiparametric protocol was performed comprising of localizer scans, PC-MRI, ASL, T 1 , T 2 * , and DWI data as described below. All scans were acquired in approximately 45 min. Figure 1 outlines the key MRI parameters within the multiparametric protocol, these measures can all be performed within a 45 min scan session. All mapping data are collected with matched geometry with slices in a coronal-oblique plane through the long axis of the kidneys, allowing automated interrogation of the resulting multiparametric maps. All data is acquired using respiratory triggering or an end-expiration breath hold to ensure data is acquired at the same point in the respiratory cycle. Each of the parameters within this protocol are outlined below.

Localizers and Kidney Volume Assessment
Balanced turbo field echo (bTFE) scans are acquired in three orthogonal planes (30 slices of 1.75 × 1.75 × 7 mm 3 resolution, data collected in single breath hold per orientation). These scans provide a localizer to allow accurate planning of subsequent images, and segmentation of these images yields total kidney volume.

Phase Contrast (PC)-MRI to Assess Renal Artery and Vein Blood Flow
Prior to the PC-MRI acquisition, an angiogram is acquired to plan the placement of the PC-MRI renal artery slice to ensure that it is positioned prior to any bifurcations of the artery.
PC-MRI is then used for the measurement of blood flow in the renal arteries and veins (Debatin et al., 1994;Schoenberg et al., 1997;Bax et al., 2005;Park et al., 2005;Dambreville et al., 2010). PC-MRI is performed using a single slice TFE image placed perpendicular to the vessel of interest. Multiple phases are collected across the cardiac cycle when imaging the renal artery (20 phases) and renal vein (15 phases). Imaging parameters use a flip angle of 25 • , reconstructed resolution 1.2 × 1.2 × 6 mm 3 , and velocity encoding of 100/50 cm/s for renal artery and vein, respectively. Each measurement is acquired during a single 15-20 s breath hold, dependent on the subjects' heart rate.
To measure renal cortex perfusion, we have implemented a respiratory-triggered FAIR (Flow-sensitive Alternating Inversion Recovery) ASL scheme. For imaging, we use either a SE-EPI or bFFE readout. Typical imaging parameters at 3 T are a post label delay (PLD) time of 1,800 ms (depending on choice of readout scheme and field strength Buchanan et al., 2015), 40 label/control pairs, 288 × 288 mm field of view, 3 × 3 × 5 mm 3 voxel resolution. A SE-EPI readout provides good spatial coverage, allowing multiple slices to be acquired in a short acquisition time over the ASL signal curve (five slices in ∼ 300 ms at 3 T). A bFFE readout provides the benefit of high spatial resolution, typically 1.5 mm in-plane spatial resolution and 5 mm slice thickness, however this can limit slice coverage due to the increased acquisition time per slice (∼ 250 ms slice spacing at 3 T). Since ASL is a subtraction technique, we use respiratory triggering to minimize the effects of respiratory motion leading to misalignment or blurring. It is important to take into account the arrival time of the blood to the tissue when quantifying perfusion, particularly in disease where the arrival time can be increased, resulting in an apparent reduction in perfusion. A separate scan to assess the arrival time of the blood to the tissue is acquired, by collecting ∼4 label/control pairs at shorter PLD times (500, 700, 900, 1,100 ms). A base equilibrium M 0 scan and T 1 map are also required for accurate perfusion quantification. Depending on respiratory rate, scan time for 40 label/control pairs of ASL data is approximately 6 min, with a further 2 min for assessment of arrival time and collection of a base M 0 scan.

Longitudinal Relaxation Time T 1 Mapping
The assessment of the longitudinal relaxation time T 1 of tissue is essential for the quantification of ASL perfusion. Recently, T 1 mapping alone has been shown to provide an important parameter by which to evaluate fibrosis (due to the association of collagen with supersaturated hydrogel) or inflammation (interstitial edema, cellular swelling). T 1 has been shown to correlate well with fibrosis and edema in the myocardium (Iles et al., 2008;Jellis and Kwon, 2014), liver (Hoad et al., 2015;Tunnicliffe et al., 2017), and more recently in the kidney (Friedli et al., 2016).
Here, an inversion recovery sequence with a modified respiratory triggering scheme (Figure 2) has been developed to minimize respiratory-induced abdominal motion between images of differing contrast collected across the range of inversion times required to compute a T 1 map. The respiratory trigger is applied at the peak of inspiration in the respiratory cycle and the image is then acquired at a constant time following this trigger, during the flat end-expiration period of the respiratory FIGURE 2 | Modified respiratory triggered inversion recovery sequence shown for (A) short (TI 1 ) and (B) long (TI 2 ) inversion time. By altering the variable delay, Tv, each image acquisition is collected at a constant time (Tv + TI) following the respiratory trigger. First arrow indicates the respiratory trigger ("Respiratory Trigger"), second the inversion pulse ("Inversion"), and the shaded block indicates the image acquisition readout ("Acquisition"). cycle. A variable delay, Tv, is introduced between the respiratory trigger and the inversion pulse which is followed by the inversion time, TI, between the inversion pulse and image acquisition. By holding the total time period Tv + TI constant, this results in all image readouts for all inversion times being collected at a constant time of (Tv + TI) following the respiratory trigger and as such all images are aligned across the inversion times. The time (Tv + TI) is chosen to be at the end-expiration period of the respiratory cycle to minimize any potential motion artifacts.
Here, we use either a SE-EPI or bFFE readout scheme for T 1 mapping. In general, the same readout scheme as is used for the ASL acquisition is chosen. Importantly, the chosen image readout scheme has an impact on the measured T 1 value. A SE-EPI readout scheme provides a "true" T 1 value, whereas a bFFE readout scheme results in an "apparent" T 1 , shorter than the "true" T 1 due to the influence of transverse relaxation rates (T 2 /T 2 * ; Schmitt et al., 2004). At 3 T, we typically collect 13 inversion times of 200, 300, 400, 500, 600, 700, 800, 900, 1,000, 1,100, 1,200, 1,300, and 1,500 ms in a total scan time of <3 min. For a multi-slice bFFE readout, the temporal slice spacing is longer (∼250 ms at 3 T) than for a SE-EPI (∼60 ms at 3 T) readout, and therefore the dynamic range of TIs can be increased by acquiring the scans ascend, descend and interleaved slice order.
An alternative scheme for T 1 mapping is to use a modified look-locker inversion recovery (MOLLI) sequence originally developed for cardiac T 1 mapping. This typically involves the acquisition of a cardiac-gated single-shot MOLLI sequence using a bFFE readout [23] with a 3(3)3(3)5 sampling pattern collected in a breath hold. However, this acquisition scheme is not best suited to the kidney, since it is cardiac triggered, requires a number of breath holds for complete coverage of the kidneys, and does not match the ASL acquisition readout scheme, and so it is not implemented in our multiparametric protocol.

Diffusion Weighted Imaging (DWI)
DWI assesses the thermally induced Brownian motion of water within tissues, which can be quantified from the Apparent Diffusion Coefficient (ADC). ADC may also be affected by factors such as tubular flow and capillary perfusion, which can be better distinguished using the IntraVoxel Incoherent Motion (IVIM) model to quantify pure diffusion (D; Le Bihan et al., 1988). In DWI, at least two single-shot echo-planar images are acquired without and with diffusion weighting gradients (b-values) from which molecular diffusion can be quantified and spatially mapped. It is important to note that the quantification of the ADC is affected by the b-values acquired. In this multiparametric protocol, DWI data is acquired with a SE-EPI readout at multiple b-values (for example, 11 b-values of 0, 5, 10, 20, 30, 50, 100, 200, 300, 400, 500 s/mm 2 ). The highest b-value is chosen such that the echo time (TE) does not become so long as to limit the signalto-noise ratio (SNR) of the image. Typically a 288 × 288 mm field of view is used with 3 × 3 × 5 mm 3 voxel resolution which has a minimum TE of 56 ms. This sequence is acquired with respiratory triggering such that the image readouts are collected at the end-expiration period. For 11 b-values, the acquisition time is approximately 8 min.

Blood Oxygenation Level Dependent (BOLD) Imaging to Assess Tissue Oxygenation
BOLD MRI exploits the paramagnetic properties of deoxygenated blood, which acts to shorten the transverse relaxation time constant (T 2 * )-alternatively expressed as the relaxation rate R 2 * (1/T 2 * )-a measure which provides an indirect non-invasive assessment of oxygen content. Higher R 2 * (or lower T 2 * ) is an indicator of lower tissue pO 2 . BOLD MRI is more sensitive at detecting changes in medullary compared to cortical pO 2 due to their relative positions on the oxygen dissociation curve-cortical pO 2 lies near the plateau of the hemoglobin oxygenation curve and medullary pO 2 lies on the linear part of the curve, thus a large change in local pO 2 is needed to cause a similar change in R 2 * for the cortex compared to the medulla. The use of BOLD MRI to measure renal oxygenation has been extensively studied. However, it should be highlighted that a number of other factors, such as, hydration status, dietary sodium intake, and susceptibility effects also alter BOLD R 2 * (Pruijm et al., 2017), this can make it difficult to draw definite conclusions from its independent use. For a review of this technique, see Pruijm et al. (2017). In this multiparametric protocol, BOLD T 2 * data is acquired using a mFFE sequence with multiple slices. Typical imaging parameters are 1.5 mm in-plane resolution, 5 mm slice thickness, initial TE 5 ms, TE spacing 3 ms, 12 echoes, flip angle 30 • . Each measurement is acquired in a single ∼17 s breath hold.

Kidney Volume Assessment
Analyze9 software (AnalyzeDirect, Overland Park, KS) is used to define a region of interest around the kidneys on each bTFE localizer image slice. Total kidney volume can then be calculated by summing across all slices, typically the coronal slices are used for organ volume measures. Analysis time is approximately 10 min.

PC-MRI Renal Blood Flow Assessment
A region of interest is placed over the vessel using Q-flow software (Philips Medical Systems, Best, NL). Mean flow velocity (cm/s), mean cross-sectional area of the lumen (mm 2 ), and hence mean bulk renal blood flow (ml/s) over the cardiac cycle, are calculated for each vessel. Total perfusion of each kidney can then be calculated by correcting the renal blood flow to kidney volume. Analysis time is approximately 2 min per vessel.

Multiparametric Interpretation
Combining multiparametric MRI maps adds considerable insight into the underlying physiology. We have developed a multiparametric image analysis program (MATLAB, The Mathworks Inc., Natick, MA) that generates and combines the parametric ASL perfusion, T 1 , diffusion, and BOLD T 2 * maps in the same data space. The multiparametric maps can then be used to perform multivariate analysis of structural and hemodynamic measures in automated regions of interest in the cortex and medulla.

Mapping perfusion from ASL data
Individual perfusion weighted difference images (controllabel) are calculated, inspected for motion (exclude >1 voxel movement) or realigned, and averaged to create a single perfusion-weighted ( M) map. M, T 1 maps (see below), and M 0 maps are then used in a kinetic model (Equation 1; Buxton et al., 1998) to calculate tissue perfusion (f ) maps (in ml/100 g tissue/min). T 1,blood is assumed to be 1.55 s at 3 T and 1.36 s at 1.5 T (Dobre et al., 2007), whilst λ, the blood-tissue partition coefficient, is assumed to be 0.8 ml/g for kidney. Analysis time is approximately 10 min.

Longitudinal relaxation time (T 1 ) mapping
Inversion recovery data is fit on a voxel-by-voxel basis to Equation (2) to generate a "true" T 1 map for the SE-EPI readout, for the bFFE readout an "apparent" T 1 map is obtained. Analysis time is approximately 3 min of user intervention, and up to 1 h processing time on a standard pc.
Mapping ADC, D, D * , and f p from DWI data DWI data are fit to form ADC maps (in mm 2 /s) by taking the log of the exponential signal decay (Equation 3). In addition, since the DWI data is collected at a number of b-values, it is possible to model the bi-exponential IVIM model (Equation 4). In the IVIM model, D (in mm 2 /s) is the pure tissue molecular diffusion coefficient representing the diffusion coefficient of slow or non-perfusion-based molecular diffusion, D * (in mm 2 /s) is the pseudodiffusion coefficient which is the fast or perfusion-based molecular diffusion representing intravoxel microcirculation or perfusion, and f p is the perfusion fraction (%) of the voxel (Le Bihan et al., 1988;Koh et al., 2011). To fit data to the IVIM model, D was first fit to Equation (3) for b-values of >200 s/mm 2 , this assumes that the pseudodiffusion component D * can be neglected above this value. Second, f p was determined from the zero intercept of this fit. Finally, D * was obtained from the monoexponential fit using the precalculated values of D and f p (Suo et al., 2015). Analysis time is approximately 5 min of user intervention, and up to 5 min processing time on a standard pc.
BOLD MRI to map T 2 * /R 2 * mFFE data are fit voxelwise using a weighted echo time (TE) fit to form T 2 * /R 2 * maps from the log of the exponential signal decay (Equation 5). Analysis time is approximately 5 min.

Interpretation of multiparametric maps
Binary whole kidney masks are formed from the manual segmentation of the base equilibrium M 0 scan or T 1 map. To distinguish renal cortex and medulla, a histogram of T 1 values across both kidneys is formed (with a bin size of 20 ms). Two peaks in the histogram, originating from the renal cortex and medulla, can be identified from which to form separate renal cortex and renal medulla masks. This segmentation procedure is illustrated for both a healthy participant and CKD patient in Figure 3. It should be noted that T 1 values are elevated in CKD [see Section Application in Chronic Kidney Disease (CKD)], but sufficient corticomedullary differentiation remains to segment the cortex from medulla. These binary cortex and medulla masks can then be applied to each parametric map (perfusion, T 1 , ADC, D, D * , and f p ) to interrogate identical regions of interest in which to assess mean values of each parameter. Importantly, to assess heterogeneity of measures and remove bias, a Gaussian curve fit can be applied to the histogram to determine both the mode and full-width-at-half-maximum (FWHM) of renal cortex and medulla parameter values across one or both kidneys (Rossi et al., 2012). The assessment of corticomedullary differentiation (medulla-cortex) in MRI parameters also provides important information. Analysis time is approximately 10 min.

RESULTS
All results given are the mean and standard deviation across participants.
Variability, Repeatability, and Field Strength Dependence in Healthy Participants Figure 4 shows example multiparametric MRI maps for a single healthy participant collected at 3 T, illustrating that the maps can be combined in the same data space and allow assessment of heterogeneity across the kidney. Table 1 provides the mean and associated standard deviation for MRI parameters collected across the cohort of healthy participants, with the number of subjects included in each analysis provided, and a comparison to literature values. Table 2 shows the field strength dependence of longitudinal (T 1 ) and transverse (T 2 * ) relaxation times. As expected, T 1values are longer at 3 T compared with 1.5 T for both the SE-EPI and bFFE readout schemes. It should be noted that the "apparent" T 1 measured using a bFFE readout scheme is shorter than the "true" T 1 measured using a SE-EPI readout. The corticomedullary differentiation of T 1 can be seen to be greater at 3 T compared with 1.5 T. Conversely, the transverse relaxation time (T 2 * ) is longer at 1.5 T. Table 3 provides the CoV and ICCs for the repeatability study at 3 T. The CoV is low for T 1 (< 2.9%), ADC (2.9%), T 2 * (4.1%), and kidney volume (4.2%). The ICCs were high for cortex perfusion (0.801), cortex and medulla T 1 (0.848 and 0.997 using SE-EPI), renal artery flow (0.844) and total kidney volume (0.985). Figure 5A shows the T 2 * mode and FWHM in renal cortex and medulla at normoxia and during hypercapnia or during hyperoxia. During hypercapnia, there was a trend for a decrease in the T 2 * mode in the renal cortex (P = 0.098, paired t-test), but the T 2 * FWHM did not change. The T 2 * mode in the renal medulla did not change, but the T 2 * FWHM was found to increase (P = 0.02, paired t-test). During hyperoxia, there was no change in T 2 * mode or FWHM in either the renal cortex or renal medulla. Figure 5B shows the "apparent" T 1 mode and FWHM for renal cortex and medulla at normoxia and during hyperoxia. During hyperoxia, there was no significant change in the "apparent" T 1 mode of renal cortex or medulla, but the FWHM increased in the renal cortex (P = 0.009, paired t-test) and medulla (P = 0.092, paired t-test).

Application in Chronic Kidney Disease (CKD)
All 11 CKD patients had glomerular kidney disease, Table 4 provides demographic details of the patients and divides the healthy participants into young (<40 years) and older (>40 years) age groups for comparison. Table 5 provides the MRI results for each of these groups.
Renal artery blood flow was significantly reduced in the older healthy participants compared to the young healthy participants, though no difference is seen between the older participants and CKD patients. In CKD, renal cortex perfusion and renal vein blood flow were lower than in older healthy participants. T 1 SE-EPI in the renal cortex was significantly higher in CKD patients compared to older healthy participants, and corticomedullary T 1 differentiation was reduced in CKD patients compared to older healthy participants. T 2 * measured in the renal cortex and medulla was not significantly different in CKD patients compared with healthy participants. In this patient cohort, renal cortex ADC, D and total kidney volume in CKD patients were also not significantly different to healthy participants.

DISCUSSION
This article has demonstrated acquisition and analysis methods to perform multiparametric assessment of the kidneys in healthy participants and CKD patients.

Variability, Repeatability, and Field Strength Dependence in Healthy Participants
We provide a comprehensive summary of MRI parameter values for healthy participants, results are in agreement with values reported across separate studies in the literature ( Table 1). When comparing T 1 measures for the renal cortex and medulla to literature values, it is important to consider the MR field strength and readout scheme used for the image acquisition. Here, we show the expected T 1 increase with field strength ( Table 2). Further, the computed T 1 value is dependent on the image readout scheme, with a shorter "apparent" T 1 measured for a bFFE readout compared to a SE-EPI readout, due to the influence of transverse relaxation on the bFFE readout. The T 1 of the medulla was higher than that of the cortex, resulting in clearly visualized corticomedullary differentiation.
The CoV of T 1 measures is very low, <3% for cortex and medulla (Table 3). Cutajar et al. reported CoVs of between 0.3 and 11.5% for repeatability of renal cortex T 1 measures on the same day (Cutajar et al., 2012). Gillis et al. showed no significant difference between visits for repeated measures of renal cortex T 1 using a MOLLI method (Gillis et al., 2014). However, MOLLI has some compromises, it is a cardiac gated scheme which provides poor sampling of the inversion recovery curve, requires a breath hold per slice, and since it uses a bFFE readout, its "apparent" T 1 value is also affected by tissue fat content at 3 T (Mozes et al., 2016).
PC-MRI measures of renal artery and vein blood flow have a reasonably high CoV, as previously described (Bax et al., 2005;Khatir et al., 2014). This is likely a result of placement of the imaging slice. In contrast, ASL renal cortex perfusion is a voxel-wise measure and this is shown to have a low CoV and high ICC. Cutajar et al. reported CoVs of between 1.8 and 12.1% for repeatability of renal cortex perfusion measures on the same day (Cutajar et al., 2012) and Chowdhury et al. reported a within session CoV of 3.3% (Chowdhury et al., 2012), but to our knowledge there have been no CoVs reported for measures collected between visits. Gillis et al. showed no significant differences between visits for repeated measures of renal cortex perfusion (Gillis et al., 2014). Thoeny et al. showed that the measured value of ADC is affected by the choice of b-values (Thoeny et al., 2005). Using only low b-values (0-100 s/mm 2 ) will result in a high calculated ADC, whilst high b-values (500-1,000 s/mm 2 ) will result in a low calculated ADC. Using a wide range of b-values provides the least variation in ADC between healthy participants. Here, we use b-values of between 0 and 500 s/mm 2 and show comparable results to Thoeny et al. (2005). Cutajar et al. found no significant difference in ADC between sessions, their ADC values were higher than we report, likely due to their acquisition using only two b-values (Cutajar et al., 2011). The value of both ADC and D had a low CoV, whilst D * and f p had poor repeatability.
A wide range of renal cortex and medulla T 2 * (R 2 * ) values are reported in the literature. T 2 * decreases with increasing field strength and is longer in the renal cortex compared to the renal medulla, indicating the hypoxic state of the medulla. The T 2 * /R 2 * values we present are in agreement with several studies (Li et al., 2004a;Ding et al., 2013;Khatir et al., 2014;Piskunowicz et al., 2015;van der Bel et al., 2016), whilst others give lower (Simon-Zoula et al., 2006;Park et al., 2012) or higher (Li et al., 2004b) R 2 * -values. Khatir et al. (2014) measured similar between session CoVs to those we present here.
There is discrepancy in the literature of the effect of breathing 100% oxygen on T 2 * . Some studies have shown no change in T 2 * in the renal cortex (Jones et al., 2002;O'Connor et al., 2009;Khatir et al., 2014;Niendorf et al., 2015) or medulla (Jones et al., 2002), whilst a small number of studies show a small increase in T 2 * in the renal cortex (Winter et al., 2011;Ganesh et al., 2016) and medulla (Donati et al., 2012;Khatir et al., 2014). At normal levels of inspired oxygen, the body maintains hemoglobin levels in arterial blood near to saturation level. During hyperoxia, a higher fraction of HbO 2 /Hb and a reduction in blood volume could both be expected to contribute to a small increase in T 2 * . As an alternative, hypercapnic-hyperoxia has been shown to cause a marked 50% increase in renal T 2 * -weighted signal intensity, suggesting this method provides enhanced sensitivity (Milman et al., 2013).
For T 1 , previous studies have shown a decrease in the renal cortex on breathing 100% oxygen, which is equivalent to ∼600 mmHg (Jones et al., 2002;O'Connor et al., 2007O'Connor et al., , 2009Ganesh et al., 2016). Here, we used our modified respiratory triggered scheme to measure T 1 and independently controlled end-tidal concentrations of oxygen and carbon dioxide (constant to ∼0.1 mmHg). No significant difference in the mode of T 1 was found between hypercapnia and normoxia. The controlled gas delivery was equivalent to breathing ∼80% oxygen, and this may explain the smaller T 1 change seen in our data. It should be noted that breathing 100% oxygen can lead to hypocapnia (Becker et al., 1996) resulting in a reduction in flow.
To our knowledge, no studies of human kidneys have used hypercapnia. Winter et al. showed no change in T 2 * at 1.5 T in the rabbit renal cortex when inspiring 10% carbon dioxide (balance air; Winter et al., 2011), whilst Ganesh et al. show a decrease in T 2 * at 3 T when inspiring 10-30% carbon dioxide (21% oxygen, balance nitrogen; Ganesh et al., 2016). Milman et al. showed that hypercapnia induced by 5% CO 2 inhalation caused a marked decline in hemodynamic response imaging maps, based on changes in the signal intensity of a T 2 * -weighted image, resembling results of studies in the liver (Milman et al., 2013). The level of inspired carbon dioxide in this work is significantly lower than 10%, this may explain why our T 2 * decrease did not reach significance.  Alternative mechanisms of physiological modulation to assess renal oxygenation and microcirculation reactivity and functionality include water loading, sodium loading, or drug administration (e.g., angiotensin, furosemide, saline). Studies have shown that water loading results in an increase in BOLD T 2 * in the medulla (Prasad et al., 1996;Prasad and Epstein, 1999;Tumkur et al., 2006a;Vivier et al., 2013;Ding et al., 2015), this is thought to be due to the production of endogenous prostaglandin PGE2 in the medulla which decreases deoxyhemoglobin (Hb) levels, but it is not possible to distinguish between changes in oxygen supply and oxygen consumption. Similar more pronounced results have been shown following administration of furosemide (a sodium pump inhibitor; Prasad et al., 1996;Li et al., 2004a;Tumkur et al., 2006b;Vivier et al., 2013), coupled with a larger increase in urinary output (Vivier et al., 2013). Interestingly, T 2 * is not altered in older subjects after water loading (Prasad and Epstein, 1999) or furosemide administration (Epstein and Prasad, 2000).

Chronic Kidney Disease
The standard clinical assessment of renal function is the estimated GFR (eGFR) calculated from serum creatinine concentration. However, this is a late marker of renal dysfunction, is often discordant with tissue damage, is subject to hemodynamic fluctuation, and cannot be used to assess individual kidney function. Kidney biopsy has sampling error associated with the small specimen size, and comes with associated risks of an invasive procedure. This pilot study has assessed the use of multiparametric MRI in CKD patients, potentially providing a number of techniques by which to assess kidney structure and function. Renal blood flow and renal cortex perfusion was lower in CKD patients compared with healthy participants. T 1 values were increased in both renal cortex and medulla compared to healthy participants, though primarily in cortex, resulting in a loss of corticomedullary differentiation.
There have been a number of previous studies assessing changes in individual MRI parameters related to hemodynamics and structure in CKD patients (Inoue et al., 2011;Michaely et al., 2012;Xin-Long et al., 2012;Khatir et al., 2014Khatir et al., , 2015Milani et al., 2016). Studies have compared perfusion in CKD patients with healthy participants and found perfusion to be lower in CKD patients (Rossi et al., 2012;Tan et al., 2014). Gillis et al. showed that the T 1 relaxation time was longer in CKD patients compared to healthy participants (Gillis et al., 2016). Further, ADC values have been shown to be reduced in CKD compared to healthy participants (Goyal et al., 2012). A recent study using DWI and T 1 mapping has demonstrated changes in both kidney ADC and T 1 in animal models and humans with CKD (Friedli et al., 2016). Prior studies have shown conflicting changes in measures of oxygenation in CKD, with some groups reporting a reduction in oxygenation in CKD, whilst others report no differences in cortical or medullary R 2 * (Pruijm et al., 2014). Khatir et al. showed similar cortical and medulla R 2 * values at baseline between patients and controls. But on inspiring 100% oxygen, R 2 * significantly decreased in the renal cortex of CKD patients with no change in R 2 * was observed in healthy participants. Medullary R 2 * increased in both patients and controls on inspiring 100% oxygen (Khatir et al., 2015). Pruijm et al. (2014) assessed patients with CKD and arterial hypertension (Pruijm et al., 2014), no difference in R 2 * was seen between the patient group and healthy participants at baseline. However, following administration of furosemide, a blunted R 2 * decrease was seen in patients compared with healthy participants. Xin-Long et al.
(2012) measured the corticomedullary differentiation in R 2 * in healthy participants and CKD patients and found an increased differentiation in CKD patients compared to healthy participants (Xin-Long et al., 2012).

Limitations
It is important to consider the different factors which can impact on reported MR measures. Inconsistent BOLD results have been widely documented between published studies, whilst Michaely et al. showed that in a study of 280 subjects, R 2 * correlated poorly with eGFR (Michaely et al., 2012). This is likely since R 2 * is only an estimator of oxygenation, and is confounded by many other factors, with it being suggested that changes in the blood volume fraction considerably influences renal T 2 * (Niendorf et al., 2015). Estimates of total renal blood flow need to consider kidney volume to also compute total perfusion, and in CKD patients the shrinkage of the kidney should be considered, which can mean that blood flow per kidney is preserved. However, this correction does not take into account that the cortex and medulla may not lose volume at the same rate.
ICC's are high for some MRI parameters presented-T 1 , perfusion, renal artery flow and ADC-but other values are relatively low, presently hampering the introduction of these methods in clinical practice. Currently, MRI is expensive and multiple breath hold methods cannot be used in older, frail patients. Here, our multiparametric protocol includes a limited number of breath holds, with ASL, T 1 , and DWI data collected using respiratory triggered acquisitions. In this study, interobserver variability was not assessed since the post-processing is automated, including ROI placement. Further automation of this pipeline could be included and with the introduction of greater processing power, maps could be computed online at the scanner. In future, functional sodium technology to provide information on renal concentrating capacity will provide a further additional measure for multiparametric protocols (Maril et al., 2006). At this point, studies showing that MRI parameters can predict hard outcomes, such as, end stage renal disease, death or rapid decline of kidney function are necessary. For ultimate use in the clinic, MRI protocols need to be time efficient, and so it will be important to define key MRI parameters of high ICC which can be used for clinical assessment.

CONCLUSIONS
This paper has outlined a multiparametric MRI acquisition and analysis protocol for assessment of renal structure, hemodynamics and oxygenation. No other modality can combine non-invasive techniques to provide such a comprehensive evaluation of renal function as MRI. Studies showing that MRI has added value to simply monitoring serum creatinine and proteinuria in kidney disease, and that MRI can provide similar information as a kidney biopsy are now eagerly awaited. The ability of early identification of patients at risk of progressing to end-stage kidney disease and protocols to assess the efficacy of treatments would improve clinical outcome, be of cost benefit for society and improve life quality for the patients.

AUTHOR CONTRIBUTIONS
EC: study design, acquisition, analysis and interpretation of data, statistical analysis, drafting of manuscript, final approval of manuscript, accountable for all aspects of the work; CEB: study design, acquisition, analysis and interpretation of data, critical revision of manuscript, final approval of manuscript, accountable for all aspects of the work; CRB: acquisition, analysis and interpretation of data, critical revision of manuscript, final approval of manuscript, accountable for all aspects of the work; BP: acquisition, analysis and interpretation of data, final approval of manuscript, accountable for all aspects of the work; HM: study design, interpretation of data, final approval of manuscript, accountable for all aspects of the work; MT: study design, interpretation of data, final approval of manuscript, accountable for all aspects of the work; NS: study design, interpretation of data, critical revision of manuscript, final approval of manuscript, accountable for all aspects of the work; SF: study design, acquisition, analysis and interpretation of data, drafting of manuscript, final approval of manuscript, accountable for all aspects of the work.

FUNDING
The authors acknowledge the financial support from the National Institute for Health Research Nottingham Digestive Diseases Biomedical Research Unit, Nottingham University Hospitals NHS Trust and University of Nottingham, the Medical Research Council Confidence in Concept Award and the Dr. Hadwen Trust. The Dr. Hadwen Trust (DHT) is the UK's leading nonanimal biomedical research charity that exclusively funds and promotes human-relevant research that replaces the use of animals whilst supporting the progress of medicine.