Research Articles

Download PDF (2.72 MB)

TOMOGRAPHY, June 2020, Volume 6, Issue 2:129-138
DOI: 10.18383/j.tom.2020.00005

Computational Modeling of Interstitial Fluid Pressure and Velocity in Head and Neck Cancer Based on Dynamic Contrast-Enhanced Magnetic Resonance Imaging: Feasibility Analysis

Eve LoCastro1, Ramesh Paudyal1, Yousef Mazaheri1, Vaios Hatzoglou2, Jung Hun Oh1, Yonggang Lu3, Amaresha Shridhar Konar1, Kira vom Eigen1, Alan Ho4, James R. Ewing5, Nancy Lee7, Joseph O. Deasy1, Amita Shukla-Dave1

Departments of 1Medical Physics and2Radiology, Memorial Sloan Kettering Cancer Center, New York, NY;3Department of Radiology, Medical College of Wisconsin, Milwaukee, WI;4Department of Medicine, Memorial Sloan Kettering Cancer Center, New York, NY;Departments of 5Neurology and6Neurosurgery, Henry Ford Hospital, Detroit, MI; and7Department of Radiation Oncology, Memorial Sloan Kettering Cancer Center, New York, NY


We developed and tested the feasibility of computational fluid modeling (CFM) based on dynamic contrast-enhanced magnetic resonance imaging (DCE-MRI) for quantitative estimation of interstitial fluid pressure (IFP) and velocity (IFV) in patients with head and neck (HN) cancer with locoregional lymph node metastases. Twenty-two patients with HN cancer, with 38 lymph nodes, underwent pretreatment standard MRI, including DCE-MRI, on a 3-Tesla scanner. CFM simulation was performed with the finite element method in COMSOL Multiphysics software. The model consisted of a partial differential equation (PDE) module to generate 3D parametric IFP and IFV maps, using the Darcy equation and Ktrans values (min−1, estimated from the extended Tofts model) to reflect fluid influx into tissue from the capillary microvasculature. The Spearman correlation (ρ) was calculated between total tumor volumes and CFM estimates of mean tumor IFP and IFV. CFM-estimated tumor IFP and IFV mean ± standard deviation for the neck nodal metastases were 1.73 ± 0.39 (kPa) and 1.82 ± 0.9 × (10−7 m/s), respectively. High IFP estimates corresponds to very low IFV throughout the tumor core, but IFV rises rapidly near the tumor boundary where the drop in IFP is precipitous. A significant correlation was found between pretreatment total tumor volume and CFM estimates of mean tumor IFP (ρ = 0.50, P = 0.004). Future studies can validate these initial findings in larger patients with HN cancer cohorts using CFM of the tumor in concert with DCE characterization, which holds promise in radiation oncology and drug-therapy clinical trials.


In solid tumors, interstitial fluid pressure (IFP) is elevated owing to increased permeability of abnormally formed tumor blood vessels and a lack of functional lymphatic removal pathways (1, 2). Elevated IFP serves to nullify the hydrostatic pressure differential between vasculature and the interstitium, hampering extravasation of drug into the tumor (3). Previous studies have shown that elevated IFP is known to have serious implications for the effective delivery of anticancer drugs (46).

Elevated tumor IFP has been invasively assessed in previous clinical and preclinical studies and is not routinely measured in clinical practice (710). The commonly accepted gold standard procedure for acute measurement of IFP is the wick-in-needle (WIN) method, which was established by Fadnes et al. (11) in rats. A previous report (12) attests its performance under different levels of hydration in dogs. In clinical settings, Gutmann et al. (8) observed high IFP values in head and neck (HN) tumors using the WIN technique (0.53–4.39 kPa [4–33 mmHg]). Use of these direct methods is limited to measuring IFP in locations where the tumor can be accessed easily. The development of a noninvasive method to estimate IFP as a surrogate biomarker in HN tumors, based on dynamic contrast-enhanced (DCE) MRI, can be a significant step forward if it can be incorporated into the early assessment of tumor response to therapy.

Multimodal therapy in various combinations of surgery, radiation, or chemotherapy is used in the management of HN cancers, depending on the primary site and TNM stage (1315). For nasopharyngeal carcinoma (NPC), treatment options include neoadjuvant chemotherapy (NAC), as it has high incidence of distant metastatic disease (16). Concurrent chemoradiation is the main approach for nonsurgical treatment of locoregionally advanced HN cancers (1719). RECIST v1.1 is the commonly used guideline tool to assess patient response to treatment and is limited by usage of bidimensional lesion measurements on anatomical imaging (20). It may be insightful to include positron emission tomography (PET) metabolic criteria (PERCIST) (21), which requires an additional fluorodeoxyglucose (18F-FDG) PET scan. MRI can provide both anatomical and functional results, which can be used to model tumor IFP estimates (22, 23).

DCE-MRI uses low-molecular-weight gadolinium-based contrast agents (CAs) to measure the influx of CA across the capillary wall (24). Pharmacokinetic modeling of the DCE data provides physiologically relevant biomarkers of perfusion and permeability for clinical endpoints such as prognosis and prediction of response in HN cancers (25, 26). For estimates of IFP in clinical HN settings, there is a need to develop computational fluid modeling (CFM) methods that use noninvasive techniques, such as DCE-MRI, and may help to understand the mechanisms of drug or radiation dose delivery in relation to tumor IFP (27), of particular relevance in regions not readily accessible by WIN.

CFM studies based on DCE-MRI, although in early stages of development, have introduced a noninvasive approach toward assessing IFP in preclinical settings (23, 28, 29) and in clinical studies on brain tumors (22, 30). The method implemented in the aforementioned studies all use Ktran, influx of CA, as an input in the Starling equation of fluid exchange accounting for both osmotic and hydrostatic pressure, and capillary hydraulic conductivity (23, 28, 31). The CFM assumes a porous medium to simulate interstitial fluid transport through tumor tissue extracellular space (pores). To the best of our knowledge, application of CFM to explore IFP in clinical HN cancers is yet unexplored, where the aforementioned modeling method has the potential to help treating physicians determine optimal radiation treatment dose for tumor control and/or assess early treatment response of novel targeted therapies, which are unmet challenges to date.

The purpose of the present clinical study is to develop and test the feasibility of quantitative CFM estimation of IFP and IFV from DCE-MRI in patient with HN cancer with locoregional lymph node (LN) metastases.


Patient Studies

The institutional review board approved and granted a waiver of informed consent for this retrospective clinical study, which was compliant with the Health Insurance Portability and Accountability Act. A total of n = 22 patients (male, 21; female, 1) with HN squamous cell carcinoma were enrolled between June 2014 and October 2015. Twenty patients were human papillomavirus–positive (HPV+) and 2 were HPV-negative (HPV−); all had neck nodal metastases with histologically proven squamous cell carcinoma. Patient characteristics are given in Table 1. Figure 1 shows the workflow pipeline steps for the CFM presented in this study.

Table 1.

Patient Characteristics

Characteristic Value (%)
Number of patients 22
Number of LN 38
Mean age (years) 59
Range (years) 44−69
    Male 21 (95)
    Female 1 (5)
Location of Primary Tumor
    Base of tongue 15 (68)
    Oropharynx 5 (23)
    Unknown 2 (9)
HPV Status
    Positive 21 (95)
    Negative 1 (5)
Figure 1.

Workflow for magnetic resonance imaging (MRI)-based computational fluid modeling (CFM) simulations. The patient undergoes a magnetic resonance (MR) examination with dynamic contrast-enhanced imaging. The images are contoured by a neuroradiologist on all slices to properly demarcate the 3-dimensional tumor structure. The images are then processed through extended Tofts model (ETM) using MRI-QAMPER software. A representative patient's individualized arterial input function (AIF) and tissue signals are plotted here. The ETM Ktran map is generated and incorporated into the simulation, along with the 3-dimensional mesh of the tumor region of interest (ROI). Finally, the CFM solves the dynamical equation to generate estimates of interstitial fluid pressure (IFP) in the tumor mesh domain.


DCE-MRI Data Acquisition

MRI studies were performed using a neurovascular phased-array coil on a 3-Tesla (T) scanner (Ingenia; Philips Healthcare; The Netherlands). The standard magnetic resonance (MR) acquisition parameters were as follows: multiplanar (axial, coronal, and sagittal) T2-weighted (T2w), fat-suppressed, fast spin-echo images (repetition time [TR] = 4000 milliseconds; echo time [TE] = 80 milliseconds; number of averages [NA] = 2; matrix = 256 × 256; slice thickness = 5 mm; field of view [FOV] = 20–24 cm), and multiplanar T1-weighted (T1w) images (TR = 600 milliseconds; TE = 8 milliseconds; NA = 2; slice thickness = 5 mm; matrix = 256 × 256; FOV = 20–24 cm). Standard T1w and T2w imaging were followed by multiple flip angle T1w imaging (to measure precontrast T1 [ie, T10]), and dynamic T1w imaging.

Precontrast T1w images were acquired using a fast spoiled gradient recalled (SPGR) echo sequence with acquisition parameters as follows: TR = 7 milliseconds; TE = 2.7 milliseconds; flip angles θ = 5°, 15°, 30°, and single-excitation; NA = 1, pixel bandwidth = 250 Hz/pixel, FOV = 20–24 cm2, 256 × 128 matrix zero-filled to 256 × 256 during image reconstruction.

Dynamic acquisition was performed with parameters identical to precontrast T1w imaging using a flip angle of 15° as mentioned above. Antecubital vein catheters delivered a bolus of 0.1 mmol/kg gadolinium-based clinically approved contrast at 2 mL/s, followed by saline flush using an MR-compatible programmable power injector (Spectris; Medrad, Indianola, PA) after acquiring the first 5–6 precontrast images. The entire node was covered contiguously with 5-mm-thick slices, zero gap, yielding 8–10 slices and 40 or 60 phases were acquired with ≤8-second temporal resolution. Total acquisition time was ≤8 minutes.

DCE-MRI Pharmacokinetic Analysis

The signal intensity in a voxel, measured from a SPGR T1w acquisition, is given by the following equation (32):

where S(t) is the voxel signal intensity at time t, M0 is the equilibrium magnetization of the protons, θ is the flip angle, TR is repetition time and TE is echo time. R1(t) (R1 = 1/T1) and R2*(t) (R2* = 1/T*2) are the time courses of longitudinal and transverse relaxation rates, respectively.

For the R1(t) calculation, TE ≪ T2*, and so approximating e-TER2*(t)1, then Equation [1] becomes the following equation (33):


The solution for the relaxation rate R1 (t) follows from Equation [2] and is given by:


In the limit of fast water exchange, the change in water proton relaxation rate (ie, ΔR1 =(R1-R10) [s−1]) owing to CA relaxivity, r1 = 4.0 [(mM)−1−1], is linearly related to the tissue CA concentration, Ct (mM) as follows:

where R10 is the precontrast longitudinal relaxation rate.

The extended Tofts model (ETM) is based on a 2-compartment model (vascular space and extravascular extracelluar space [EES]). The ETM expression for modeling Ct(t) is given by the following equation (24):

where, Ktran (min−1) is the volume transfer constant of CA, Cp(t) is the time-course of plasma CA concentration, which is called arterial input function (AIF). The rate constant describing CA back-flux into the vascular space from the EES is defined as kepKtran/ve, where ve and vp are the volume fractions of the EES and blood plasma, respectively.

Individual AIF for each patient was determined from an algorithm that uses cross correlation to rate arterial voxel time series in comparison with an ideal AIF signal. Signal from arterial voxels with strongest similarity to ideal AIF was averaged to be used for each patient AIF (34). Regions of interest (ROIs) on tumors were contoured manually by an experienced neuroradiologist on a late phase of the T1w DCE images using ITK-SNAP (35). The extent of necrosis was evaluated by the same neuroradiologist on T2w and postcontrast T1w MRI, reflecting tumor heterogeneity. The total tumor volume was calculated on T2w images for the 38 neck LN metastases of the 22 patients.

T10 values were calculated from multiple flip angle precontrast T1w images (36). The time course of tissue concentration data, Ct(t), was fitted (Equation [5]) using a nonlinear curve fitting technique that minimizes the sum of squared errors (SSEs) between model fit and data. The fitting procedure produces parametric maps of Ktran, ve, and vp, and parameter estimation values were bounded in the fitting routine as follows: Ktran ∈ [0, 5] (min−1), ve, and vp ∈ [0, 1]. Automation of AIF findings, quantitative pharmacokinetic analysis of DCE-MRI data, and generation of parametric maps were performed using in-house MRI-QAMPER software (Quantitative Analysis Multi-Parametric Evaluation Routines) approved by NCI Quantitative Imaging Network for Technical Benchmark (December 2019) written in MATLAB (The MathWorks, Inc., Natick, MA), as detailed in Paudyal et al. (37).

CFM Mathematical Model

The fluid mechanics of the system are given by the Navier–Stokes hydrodynamic equation. The mass-balance describing the continuity of flow (1/s) in and out of interstitial space with source and sink conditions gives a generalized description of fluid velocity from forces acting on, within, and between mass elements (38).


Where u is the fluid velocity vector (m/s), pi is the (interstitial) fluid pressure (Pa), K is the permeability of the porous membrane (m2, describing the packing of the extracellular matrix), ρ is the interstitial fluid density (kg/m3), μ is interstitial fluid viscosity (Pa*s), additional forces acting on the fluid volume element may be represented as F (N), and I and T are the notation for the identity matrix and the transpose of the operator, respectively.

Fluid movement through space is approximated with low Reynolds number flow (39) and modeled under assumption of steady-state velocity. Neglecting both friction within fluid and exchange of momentum between fluid and solid phases, the case of an incompressible fluid (·u=0) in the Navier–Stokes equation (Equation [6]) simplifies to Darcy's law where interstitial fluid velocity (IFV) is related to the gradient in IFP (pi):


This Darcy velocity u describes bulk fluid movement within the interstitial space and expresses velocity as the gradient of instantaneous pressure, proportional to the local hydraulic conductivity KHKμ:


Fluid source-sink assumptions are as follows: fluid can enter the interstitium via the vascular compartment (ϕv); fluid is drained from the system via lymphatic pathways (ϕL) (1).


Blood flow transport ϕv across capillary walls was regulated according to the Starling law (1, 31):


Where LP is the hydraulic conductivity of the capillary wall, or vessel permeability, S/V is microvascular surface area per unit volume,pV is the blood pressure in the microvessel, pi is IFP, πV is osmotic pressure in microvasculature, πi is osmotic pressure in interstitial space, and σT is the osmotic reflection coefficient.

Lymphatic drainage of fluid is assumed to be available in normal tissue outside of the metastatic LN (40, 41), dependent on the difference between local interstitial pressure,pi,and pressure of the lymphatic system, pL and the lymphatic clearance rate, LpLSLV:


In the subsequent sensitivity analysis, we observe the change in IFP profile if we create a simulated sink, or removal, term in the tumor.

Combining Equations [9], [10], and [11], the continuity equation (becomes):


Numerical values for physical parameters were based on the default selections from prior literature, listed with references in Table 2.

Table 2.

Tissue and Vascular Parameters Used in Simulations

Parameter Unit Description (# References) Value
Lp m Pa−1 s−1 Vessel permeability (1, 23) 2 × 10−11 (tumor)
3 × 10−12 (normal)
LpLSL/V Pa−1 s−1 Lymphatic filtration coefficient (23) 1 × 10−7
K m2 Pa−1 s−1 Hydraulic conductivity (23) 1.9 × 10−12 (tumor)
3.8 × 10−13 (normal)
S/V m−1 Microvascular surface area per unit volume (1, 47) 2 × 104 (tumor)
7 × 103 (normal)
pV Pa Microvascular pressure (59) 2,300
πi Pa Osmotic pressure in interstitial space (1, 46, 47) 3,230 (tumor)
1,330 (normal)
πV Pa Osmotic pressure in microvasculature (46) 2,670
σT Unitless Average osmotic reflection coefficient for plasma (1, 50) 0.82 (tumor)
0.91 (normal)

Computational Fluid Modeling

As mentioned above, the influx of fluid through the capillary wall to the interstitial space is described by the continuity equation, consisting of the source and sink terms. The ETM-derived Ktran is incorporated into the Starling equation (Equation [11]), and scaled by the mean tumor value, <Ktran> to account for the heterogeneous bulk vascular influx (23, 28). The final expression for the continuity equation (Equation [12]) in terms of the dependent variable interstitial pressure, pi:


The continuity partial differential equation (PDE) (Equation [13]) was implemented using the COMSOL CFM simulation PDE module. Solving (Equation [13]) provides the basis for estimation of pi and 3D parametric maps of IFP and IFV.

The physiological 3D mesh model was generated from each patient's T1w DCE tumor images. The contoured tumor ROIs in the model were converted to binary masks and extended by 10 pixels to represent normal tissue around tumor. The ROIs for tumor with normal surrounding tissue were resliced to be 1 mm3 isotropic in MATLAB using the NIfTI toolbox (42), and converted to stereolithography (STL) file format. STL files were imported into the simulation software and interpreted as boundary meshes for the model. In the present study, a simplified geometry was used to represent a complex LN tumor structure. On imaging, intranodal tumor necrosis is observed, which is a mix of tumor necrosis, keratin pooling, fibrous tissue, edema, viable tumor cells, and possibly hemorrhage (43). Usually, a localized group of nodes is present in an expected nodal draining area for a specific primary tumor (43). Once extracapsular extension has occurred, the LN tumor can extend to invade surrounding tissue (43, 44).

ETM-estimated Ktran maps were resliced to matching space to ensure accurate spatial interpolation (45), and converted from units of min−1 to s−1 to carry out calculations in standard SI units; the Ktran values were interpreted in COMSOL as a scalar field over the mesh domain. The Ktran field values were normalized by the mean Ktran value in the tumor ROI (23). Numerical values for physical constants in normal and tumor tissue were respectively assigned to the appropriate regions of the 3D STL domain mesh, as given in Table 2. These physiological values are taken from those previously cited in literature (1, 23, 46, 47). To generate IFP maps, a stationary solution of the PDE in Equation [13] was computed on the 3D extended domain ROI, and a no-flux condition was implemented at mesh boundary edges.

Simulation was conducted using the general coefficient form PDE module in a commercial multiphysics software package (COMSOL Inc., Stockholm, Sweden) which uses finite element method to solve PDE equations.

IFP and IFV Sensitivity Analysis

A sensitivity analysis was performed by parameter variation in the CFM models to test the influence on IFP and IFV estimates as detailed in a preclinical setting (23). The parameter values in the CFM were sequentially varied in separate simulation runs: vessel permeability (in normal and tumor tissue, Lp), lymphatic clearance (LpLSL/V in normal tissue and tumor), and ratio of hydraulic conductivity (allowing KH,tumor to vary relative to KH,normal). The range of high and low values was chosen to scale 2× and 0.5× relative to the baseline value (Table 2). Values were assigned ([LpLSL/V]t/[LpLSL/V]n,baseline = 0.125 and 0.5) to study the effect of fluid clearance from the tumor region (either via lymphatics or vascular reabsorption).

Statistical Analysis

Spearman correlation coefficient (ρ) was calculated between the total tumor volume and CFM estimates of IFP and IFV. The results were reported as mean ± standard deviation. In addition, heterogeneity measures from the histogram were calculated as skewness and kurtosis. A P-value < .05 was considered statistically significant.


All 22 patients with HN cancer underwent MRI pretreatment (pre-TX), and 38 neck LN metastases were analyzed. Out of the 22 patients, 10 patients had ≥2 metastatic nodes. Figure 2 shows the acquired MR images of a representative patient (patient #1) with a right lateral LN (total tumor volume = 14.26 cm3).

Figure 2.

T1-weighted postcontrast image of patient #1 with right lateral neck nodal metastasis (yellow outline) (A); T1-weighted (T1w) postcontrast image with overlaid ROI for analysis (green denotes normal tissue, red denotes viable tumor tissue, blue denotes necrotic core) (B); Ktran map of tumor ROI (C); preview of the imported geometry mesh of the tumor (blue) for patient #1 for CFM (D).


The CFM estimation maps of IFP and IFV are shown in Figure 3 for representative patients #1, 2, 3 with and 4 with HN cancer. The presence of cystic necrosis for the nodes in patients #1 and 4 was evidenced by hypointense signal in postcontrast T1w images and hyperintensity in the T2w images. In contrast, the nodes of patients #2 and 3 did not exhibit necrosis on imaging. The subtle differences in IFP and IFV profiles are evident for the necrotic tumors.

Figure 3.

Visualization of IFP and interstitial fluid volume (IFV) maps as estimated by simulation in patient #1, 2, 3, and 4. The profiles were calculated along the lines drawn on the corresponding maps. Row 1: T1w postcontrast image (inset: T2w view of tumor) (V1 = 14.26 cm3, V2 = 25.20 cm3, V3 = 31.27 cm3, V4 = 28.80 cm3); Row 2: Estimated IFP map of tumor and surrounding normal tissue (mean tumor pressure P¯1= 1.73 kPa, P¯2 = 2.05 kPa, P¯3 = 2.09 kPa,P¯4 = 1.77 kPa); Row 3: IFP profile along a vertical bisector line (tumor boundary denoted by dashed line); Row 4: Estimated IFV map of tumor and surrounding normal tissue (mean tumor velocity u¯1 = 1.60 × 10−7m/s, u¯2 = 1.80 × 10−7m/s, u¯3 = 1.67 × 10−7m/s, u¯4 = 1.51 × 10−7m/s); Row 5: IFV profile along vertical bisector line, with tumor boundary (dash).


In Figure 3, mean IFP (pmean) for the cystic node in patient #1 was 1.73 kPa with maximum IFP at tumor core pmax = 2.23 kPa; the IFP profile demonstrated slowly changing pressure with little curvature throughout the tumor core. The resulting IFV profile for patient #1 (IFV gradient of IFP, by the Darcy law) results in negligible velocity throughout the tumor core. The estimated pmean and pmax for the metastatic node in patient #2 was 2.05 kPa, and 2.64 kPa, respectively. Values of pmean and pmax for patients #3 and 4 were 2.09, 1.65 kPa (pmean) and 2.67, 2.04 kPa (pmax), respectively. In the IFP profiles of patients #2 and 3 with non-necrotic tumors, a moderate curvature is seen, and the resulting IFV is more variable and elevated in the tumor core.

Table 3 summarizes the mean(± standard deviation) values for total tumor volume (cm3), Ktran (min−1), IFP (kPa), and IFV (m/s) obtained from the LN metastases. The mean total tumor volume for the 38 nodes was V = 15.83 ±10.80 cm3 ranging from 1.36 to 40.02 cm3.

Table 3.

Summary of DCE-MRI Ktrans- and CFM-estimated IFP and IFV from 22 Patients with Head and Neck Squamous Cell Carcinoma

Parameter (Unit) Value (Mean ± SD)
Total Tumor Volume (cm3) 15.83 ± 10.80
Ktrans (min−1) 0.02 ± 0.02
Interstitial Fluid Pressure
(IFP [kPa])
1.73 ± 0.39
IFP Skewness −0.48
IFP Kurtosis 3.10
Interstitial Fluid Velocity
(IFV [×10−7 m/s])
1.82 ± 0.89
IFV Skewness 0.79
IFV Kurtosis 3.84

Figure 4 exhibits the scatter plot of pre-TX total tumor volume versus mean CFM-estimated IFP. Significant correlation (Spearman) was found between pre-TX total tumor volume and the mean IFP (ρ = 0.5, P = .004). No significant correlation was found between mean tumor IFV and total tumor volume (ρ = −0.1, P = .5).

Figure 4.

Plot of total tumor volume, calculated from T2-weighted image, and mean tumor IFP as estimated by COMSOL simulation, for all pretreatment neck lymph node (LN) metastases. A significant correlation between tumor volume and the mean intratumor pressure is found (ρ = 0.5, P = .004).


Figure 5 shows the of IFP and IFV profiles from the sensitivity analysis as carried out for representative patients #1 and 2. The parameter with the greatest effect on amplitude of IFP values was the tumor vessel permeability (Lp). For the necrotic tumor in patient #1, changes to parameters affected the magnitude of IFP but the shape and slope of the IFP field did not vary greatly and hence had little influence on the internal IFV. The non-necrotic node in patient #2 showed greater IFP sensitivity to all changes in parameters, with more variability in the magnitude of predicted IFV, and maintaining IFV appreciably higher.

Figure 5.

Sensitivity analyses of IFP and IFV profiles for patients #1 (A, B) and #2 (C, D) as a result of varying model parameters; in each case, the high value is 2× the default parameter (as listed in Table 2), and the low value is 0.5× the default. Subscript t, n denote values in tumor and normal, respectively.



This is the first feasibility study to apply CFM using DCE data to estimate IFP/IFV in HN LN metastases. A 3D mesh tumor ROI, derived from the contours on DCE imaging, and the continuity partial differential equation (to) describe the physiological processes (Starling law and Darcy law), were combined and solved with the finite element method (COMSOL Multiphysics computing environment) using the 3D Ktran values. The Ktranvoxel values in the CFM were normalized by the mean Ktran values in each tumor ROI. CFM-generated IFP and IFV maps depict heterogeneity within the tumor, stemming from incorporation of Ktran to inform fluid influx rate variability from the spatial CA extravasation, similar to Pishko et al. (23).

The CFM-estimated mean tumor IFP showed a moderate but significant correlation with the total tumor volume at pre-TX. This finding is in good agreement with previous clinical studies in HN, cervical, and breast cancer that directly measured IFP by using an invasive WIN method (8, 48, 49). Across different tumor types, the mean IFP value may vary but the relationship between IFP and tumor volume still holds (8, 48, 49). IFP is elevated within the tumor and decays sharply at its boundary, and in accordance with the Darcy law, IFV attains a maximum value at the tumor boundary and movement is directed outward from the tumor (1, 50, 51). Most of the high-velocity fluid movement exists near the periphery of the tumor owing to a precipitous fall-off of pressure (52). Tumors with mild-to-severe necrotic fractions in their interior, as in patients #1 (severe) and 4 (mild), tend to quickly reach maximum IFP, starting off a gradual plateau and varying less inside the tumor compared with patients #2 and 3, which attained a central local maximum of pressure within the core (Figure 3). The physical impact of cystic necrosis on the IFV is shown in these contrasting profiles. The variability of IFP in the non-necrotic tumor in contrast to necrotic tumor may lead to a significant change in IFV of tumor core.

Sensitivity analyses illustrated the influence that tissue-descriptive parameters, such as hydraulic conductivity, KH, or the effectiveness of lymphatic vessels, have on estimates of IFP and IFV in the CFM. The analysis also assesses the stability of the results, in spite of estimated values from literature. In the core of the necrotic node, IFV was not changed from the baseline despite changes to the model parameters. In contrast, the IFV in the node with a large viable fraction showed greater tendency for variation depending on parameter values. Differing from the conditions of the main analysis, we considered lymphatic functionality for the tumor core. Tumor drainage introduces a distinct change in curvature of both patients' IFP profiles and is the only factor that introduces variability to IFV values in the core of the necrotic tumor. In the patient with necrotic tumor, increasing hydraulic conductivity ratio KH,t/KH,n has no noticeable effect from baseline value, but in the non-necrotic case, a lower-ratio KH,t/KH,n results in lower overall IFP with the greatest IFV compared with all other parameter changes. In contrast, for the non-necrotic tumor, a lower lymphatic clearance in normal tissue has no impact on IFP/IFV from baseline, but it effectively increases the IFP in the necrotic tumor. In general, any fluid delivered into a tumor will be mostly reliant on diffusive transport or be trapped owing to the approximately zero internal IFV. Intravenous drug delivery effectiveness may be affected by the presence of necrotic core and clearance mechanisms therein.

The effect of mesh resolution or partitioning size on resulting finite element method simulations was shown to be minimal. CFM studies in Pishko et al. (23) and Magdoom et al. (53) reanalyzed the same preclinical DCE data set using different quantization schemes for the mesh; their result shows the resulting pressure maps were comparable. Although the absolute IFP magnitude was affected by differences in spatial discretization due to changes in voxel size, the overall profile trends in IFP maps were not changed.

Known limitations of this study include literature-based parameter values in tumor and normal tissue (23), a shortcoming of all similarly constructed IFP simulations (1, 5457). Physiologically accurate approaches in the future would aim to include additional imaging-based biomarkers to characterize the model parameters. A B1 nonuniformity acquisition is required for high-field DCE-MRI. The accuracy of ETM-derived Ktran values would be improved by inclusion of a B1 correction (58). Individual AIF measurements for large patient cohort studies can be a challenge, which can be addressed by using an average AIF (34). Lastly, the repeatability of the metric and comparison to a gold standard WIN IFP measurement needs to be assessed before the IFP derived from CFM method can be established as a prognostic/predictive biomarker.

In conclusion, our initial findings establish that DCE-MRI-based CFM can noninvasively provide estimates of tumor IFP and IFV in neck nodal metastases. After validation, IFP estimates from CFM may help treating physicians better understand and plan optimal radiation treatment dose for tumor control, and/or assess early treatment response of novel targeted therapies.


[1] Abbreviations:


arterial input function


contrast agent


computational fluid modeling


dynamic contrast-enhanced magnetic resonance imaging


extravascular extracellular space


extended Tofts model


field of view


head and neck


interstitial fluid pressure


interstitial fluid volume


lymph node


partial differential equation


regions of interest


spoiled gradient recalled echo


sum of squared errors




echo time


repetition time


number of averages




Equal contribution: E.L. and R.P. contributed equally to this work.

This work was supported by U01 CA211205 and in part through the NIH/NCI Cancer Center Support Grant P30 CA008748. We would like to thank the MRI technologists for their great efforts in helping to perform the MRI examinations, and Mr. Christian Czmielewski for his helpful contribution to patient enrollment and data management. This work was supported by the MSKCC internal IMRAS grant, U01 CA211205, and in part through the NIH/NCI Cancer Center Support Grant: P30 CA008748.


    Baxter LT, Jain RK. Transport of fluid and macromolecules in tumors. I. Role of interstitial pressure and convection. Microvasc Res. 1989;37:77–104.
    Leu AJ, Berk DA, Lymboussaki A, Alitalo K, Jain RK. Absence of functional lymphatics within a murine sarcoma: a molecular and functional evaluation. Cancer Res. 2000;60:4324–4327.
    Jain RK. Delivery of molecular medicine to solid tumors. Science. 1996;271:1079–1080.
    Jain RK. Barriers to drug delivery in solid tumors. Sci Am. 1994;271:58–65.
    Heldin CH, Rubin K, Pietras K, Ostman A. High interstitial fluid pressure - an obstacle in cancer therapy. Nat Rev Cancer. 2004;4:806–813.
    Lunt SJ, Fyles A, Hill RP, Milosevic M. Interstitial fluid pressure in tumors: therapeutic barrier and biomarker of angiogenesis. Future Oncol. 2008;4:793–802.
    Boucher Y, Kirkwood JM, Opacic D, Desantis M, Jain RK. Interstitial hypertension in superficial metastatic melanomas in humans. Cancer Res. 1991;51:6691–6694.
    Gutmann R, Leunig M, Feyh J, Goetz AE, Messmer K, Kastenbauer E, Jain RK. Interstitial hypertension in head and neck tumors in patients: correlation with tumor size. Cancer Res. 1992;52:1993–1995.
    Boucher Y, Salehi H, Witwer B, Harsh GR, Jain RK. Interstitial fluid pressure in intracranial tumours in patients and in rodents. Br J Cancer. 1997;75:829–836.
    Taghian AG, Abi-Raad R, Assaad SI, Casty A, Ancukiewicz M, Yeh E, Molokhia P, Attia K, Sullivan T, Kuter I, Boucher Y, Powell SN. Paclitaxel decreases the interstitial fluid pressure and improves oxygenation in breast cancers in patients treated with neoadjuvant chemotherapy: clinical implications. J Clin Oncol. 2005;23:1951–1961.
    Fadnes HO, Reed RK, Aukland K. Interstitial fluid pressure in rats measured with a modified wick technique. Microvasc Res. 1977;14:27–36.
    Wiig H, Reed RK, Aukland K. Measurement of interstitial fluid pressure in dogs: evaluation of methods. Am J Physiol. 1987;253(2 Pt 2):H283–90.
    Argiris A, Li S, Savvides P, Ohr JP, Gilbert J, Levine MA, Chakravarti A, Haigentz Jr M, Saba NF, Ikpeazu CV, Schneider CJ, Pinto HA, Forastiere AA, Burtness B. Phase III randomized trial of chemotherapy with or without bevacizumab in patients with recurrent or metastatic head and neck cancer. J Clin Oncol. 2019;37:3266–3274.
    Beckham TH, Barney C, Healy E, Wolfe AR, Branstetter A, Yaney A, Riaz N, McBride SM, Tsai CJ, Kang J, Yu Y, Chen L, Sherman E, Dunn L, Pfister DG, Tan J, Rupert R, Bonomi M, Zhang Z, Lobaugh SM, Grecula JC, Mitchell DL, Wobb JL, Miller ED, Blakaj DM, Diavolitsis VM, Lee N, Bhatt AD. Platinum-based regimens versus cetuximab in definitive chemoradiation for human papillomavirus-unrelated head and neck cancer. Int J Cancer. 2019.
    Rosenthal DI, Istenmaa DA, Glatstein E. A review of neoadjuvant chemotherapy for head and neck cancer: partially shrunken tumors may be both leaner and meaner. Int J Radiat Oncol Biol Phys. 1994;28:315–320.
    Han JE, Yi SK, Wang S, Erman A, Bearelly S, Sindhu S, Robbins JR, Bauman J, Hsu CC. Neoadjuvant chemotherapy improves survival compared with concurrent chemoradiation alone in nasopharyngeal carcinoma patients with N3 disease. Head Neck. 2019;41:4076–4087.
    Ang KK, Harris J, Garden AS, Trotti A, Jones CU, Carrascosa L, Cheng JD, Spencer SS, Forastiere A, Weber RS. Concomitant boost radiation plus concurrent cisplatin for advanced head and neck carcinomas: radiation therapy oncology group phase II trial 99-14. J Clin Oncol. 2005;23:3008–3015.
    Dilling TJ, Bae K, Paulus R, Watkins-Bruner D, Garden AS, Forastiere A, Kian Ang K, Movsas B. Impact of gender, partner status, and race on locoregional failure and overall survival in head and neck cancer patients in three radiation therapy oncology group trials. Int J Radiat Oncol Biol Phys. 2011;81:e101–9.
    Setton J, Caria N, Romanyshyn J, Koutcher L, Wolden SL, Zelefsky MJ, Rowan N, Sherman EJ, Fury MG, Pfister DG, Wong RJ, Shah JP, Kraus DH, Shi W, Zhang Z, Schupak KD, Gelblum DY, Rao SD, Lee NY. Intensity-modulated radiotherapy in the treatment of oropharyngeal cancer: an update of the Memorial Sloan-Kettering Cancer Center experience. Int J Radiat Oncol Biol Phys. 2012;82:291–298.
    Eisenhauer EA, Therasse P, Bogaerts J, Schwartz LH, Sargent D, Ford R, Dancey J, Arbuck S, Gwyther S, Mooney M, Rubinstein L, Shankar L, Dodd L, Kaplan R, Lacombe D, Verweij J. New response evaluation criteria in solid tumours: revised RECIST guideline (version 1.1). Eur J Cancer. 2009;45:228–247.
    Wahl RL, Jacene H, Kasamon Y, Lodge MA. From RECIST to PERCIST: evolving considerations for PET response criteria in solid tumors. J Nucl Med. 2009;50:122S–150S.
    Bhandari A, Bansal A, Singh A, Sinha N. Perfusion kinetics in human brain tumor with DCE-MRI derived model and CFD analysis. J Biomech. 2017;59:80–89.
    Pishko GL, Astary GW, Mareci TH, Sarntinoranont M. Sensitivity analysis of an image-based solid tumor computational model with heterogeneous vasculature and porosity. Ann Biomed Eng. 2011;39:2360–2373.
    Tofts PS, Brix G, Buckley DL, Evelhoch JL, Henderson E, Knopp MV, Larsson HBW, Lee T-Y, Mayr NA, Parker GJM, Port RE, Taylor J, Weisskoff RM. Estimating kinetic parameters from dynamic contrast-enhanced T(1)-weighted MRI of a diffusable tracer: standardized quantities and symbols. J Magn Reson Imaging. 1999;10:223–232.
    Shukla-Dave A, Lee NY, Jansen JFA, Thaler HT, Stambuk HE, Fury MG, Patel SG, Moreira AL, Sherman E, Karimi S, Wang Y, Kraus D, Shah JP, Pfister DG, Koutcher JA. Dynamic contrast-enhanced magnetic resonance imaging as a predictor of outcome in head-and-neck squamous cell carcinoma patients with nodal metastases. Int J Radiat Oncol Biol Phys. 2012;82:1837–1844.
    Kim S, Loevner LA, Quon H, Kilger A, Sherman E, Weinstein G, Chalian A, Poptani H. Prediction of response to chemoradiation therapy in squamous cell carcinomas of the head and neck using dynamic contrast-enhanced MR imaging. AJNR Am J Neuroradiol. 2010;31:262–268.
    Lee N, Schoder H, Beattie B, Lanning R, Riaz N, McBride S, Katabi N, Li D, Yarusi B, Chan S, Mitrani L, Zhang Z, Pfister DG, Sherman E, Baxi S, Boyle J, Morris LGT, Ganly I, Wong R, Humm J. Strategy of using intratreatment hypoxia imaging to selectively and safely guide radiation dose de-escalation concurrent with chemotherapy for locoregionally advanced human papillomavirus-related oropharyngeal carcinoma. Int J Radiat Oncol Biol Phys. 2016;96:9–17.
    Zhao J, Salmon H, Sarntinoranont M. Effect of heterogeneous vasculature on interstitial transport within a solid tumor. Microvasc Res. 2007;73:224–236.
    Steuperaert M, Debbaut C, Carlier C, De Wever O, Descamps B, Vanhove C, Ceelen W, Segers P. A 3D CFD model of the interstitial fluid pressure and drug distribution in heterogeneous tumor nodules during intraperitoneal chemotherapy. Drug Deliv. 2019;26:404–415.
    Bhandari A, Bansal A, Singh A, Gupta RK, Sinha N. Comparison of transport of chemotherapeutic drugs in voxelized heterogeneous model of human brain tumor. Microvasc Res. 2019;124:76–90.
    Starling EH. On the absorption of fluids from the connective tissue spaces. J Physiol. 1896;19:312–326. sp000596.
    Bagher-Ebadian H, Jain R, Nejad-Davarani SP, Mikkelsen T, Lu M, Jiang Q, Scarpace L, Arbab AS, Narang J, Soltanian-Zadeh H, Paudyal R, Ewing JR. Model selection for DCE-T1 studies in glioblastoma. Magn Reson Med. 2012;68:241–251.
    Aryal MP, Nagaraja TN, Brown SL, Lu M, Bagher-Ebadian H, Ding G, Panda S, Keenan K, Cabral G, Mikkelsen T, Ewing JR. Intratumor distribution and test-retest comparisons of physiological parameters quantified by dynamic contrast-enhanced MRI in rat U251 glioma. NMR Biomed. 2014;27:1230–1238.
    Shukla-Dave A, Lee N, Stambuk H, Wang Y, Huang W, Thaler HT, Patel SG, Shah JP, Koutcher JA. Average arterial input function for quantitative dynamic contrast enhanced magnetic resonance imaging of neck nodal metastases. BMC Med Phys. 2009;9:4.
    Yushkevich PA, Piven J, Hazlett HC, Smith RG, Ho S, Gee JC, Gerig G. User-guided 3D active contour segmentation of anatomical structures: significantly improved efficiency and reliability. Neuroimage. 2006;31:1116–1128.
    Deoni SC, Josseau MJ, Rutt BK, Peters TM. Visualization of thalamic nuclei on high resolution, multi-averaged T1 and T2 maps acquired at 1.5 T. Hum Brain Mapp. 2005;25:353–359.
    Paudyal R, Lu Y, Hatzoglou V, Moreira A, Stambuk HE, Oh JH, Cunanan KM, Nunez DA, Mazaheri Y, Gonen M, HO A, Fagin JA, Wong RJ, Shaha A, Tuttle RM, Shukla-Dave A. Dynamic contrast-enhanced MRI model selection for predicting tumor aggressiveness in papillary thyroid cancers. NMR Biomed. 2020;33(1):e4166.
    Dullien FAL. Porous media: fluid transport and pore structure. NY: Academic Press; 1979.
    Yao W, Li Y, Ding G. Interstitial fluid flow: the mechanical environment of cells and foundation of meridians. Evid Based Complement Alternat Med. 2012;2012:853516.
    Baxter LT, Jain RK. Vascular permeability and interstitial diffusion in superfused tissues: a two-dimensional model. Microvasc Res. 1988;36:108–115.
    Padera TP, Kadambi A, di Tomaso E, Carreira CM, Brown EB, Boucher Y, Choi NC, Mathisen D, Wain J, Mark EJ, Munn LL, Jain RK. Lymphatic metastasis in the absence of functional intratumor lymphatics. Science. 2002;296:1883–1886.
    Shen J. Tools for NIfTI and ANALYZE image 2014. Available from:
    Sakai O, Curtin HD, Romo LV, Som PM. Lymph node pathology. Benign proliferative, lymphoma, and metastatic disease. Radiol Clin North Am. 2000;38:979–998.
    Som PM. Detection of metastasis in cervical lymph nodes: CT and MR criteria and differential diagnosis. AJR Am J Roentgenol. 1992;158:961–969.
    Magdoom KN, Pishko GL, Kim JH, Sarntinoranont M. Evaluation of a voxelized model based on DCE-MRI for tracer transport in tumor. J Biomech Eng. 2012;134:091004.
    Stohrer M, Boucher Y, Stangassinger M, Jain RK. Oncotic pressure in solid tumors is elevated. Cancer Res. 2000;60:4251–4255.
    Tan WH, Wang F, Lee T, Wang CH. Computer simulation of the delivery of etanidazole to brain tumor from PLGA wafers: comparison between linear and double burst release systems. Biotechnol Bioeng. 2003;82:278–288.
    Haider MA, Milosevic M, Fyles A, Sitartchouk I, Yeung I, Henderson E, Lockwood G, Lee TY, Roberts TPL. Assessment of the tumor microenvironment in cervix cancer using dynamic contrast enhanced CT, interstitial fluid pressure and oxygen measurements. Int J Radiat Oncol Biol Phys. 2005;62:1100–1107.
    Nathanson D, Nelson L. Interstitial fluid pressure in breast cancer, benign breast conditions, and breast parenchyma. Ann Surg Oncol. 1994;1:333–338.
    Baxter LT, Jain RK. Transport of fluid and macromolecules in tumors. II. Role of heterogeneous perfusion and lymphatics. Microvasc Res. 1990;40:246–263.
    Simonsen TG, Lund KV, Hompland T, Kristensen GB, Rofstad EK. DCE-MRI-derived measures of tumor hypoxia and interstitial fluid pressure predict outcomes in cervical carcinoma. Int J Radiat Oncol Biol Phys. 2018;102:1193–1201.
    Ewing JR, Nagaraja TN, Aryal MP, Keenan KA, Elmghirbi R, Bagher-Ebadian H, Panda S, Lu M, Mikkelsen T, Cabral G, Brown SL. Peritumoral tissue compression is predictive of exudate flux in a rat model of cerebral tumor: an MRI study in an embedded tumor. NMR Biomed. 2015;28:1557–1569.
    Magdoom KN, Pishko GL, Rice L, Pampo C, Siemann DW, Sarntinoranont M. MRI-based computational model of heterogeneous tracer transport following local infusion into a mouse hind limb tumor. PLoS One. 2014;9 e89594.
    Sarntinoranont M, Rooney F, Ferrari M. Interstitial stress and fluid pressure within a growing tumor. Ann Biomed Eng. 2003;31:327–335.
    Pozrikidis C, Farrow DA. A model of fluid flow in solid tumors. Ann Biomed Eng. 2003;31:181–194.
    Smith JH, Humphrey JA. Interstitial transport and transvascular fluid exchange during infusion into brain and tumor tissue. Microvasc Res. 2007;73:58–73.
    Liu LJ, Brown SL, Ewing JR, Schlesinger M. Phenomenological model of interstitial fluid pressure in a solid tumor. Phys Rev E Stat Nonlin Soft Matter Phys. 2011;84:021919(2 Pt 1).
    Sengupta A, Gupta RK, Singh A. Evaluation of B1 inhomogeneity effect on DCE-MRI data analysis of brain tumor patients at 3T. J Transl Med. 2017;15:242.
    Boucher Y, Jain RK. Microvascular pressure is the principal driving force for interstitial hypertension in solid tumors: implications for vascular collapse. Cancer Res. 1992;52:5110–5114.


Download the article PDF (2.72 MB)

Download the full issue PDF (12.51 MB)

Mobile-ready Flipbook

View the full issue as a flipbook (Desktop and Mobile-ready)