Imaging-based characterization of convective tissue properties

Abstract

Convective transport is an important phenomenon for nanomedicine delivery. We present an imaging-based approach to recover tissue properties that are significant in the accumulation of nanoparticles delivered via systemic methods. The classical pharmacokinetic analysis develops governing equations for the particle transport from a first principle mass balance. Fundamentally, the governing equations for compartmental mass balance represent a spatially invariant mass transport between compartments and do not capture spatially variant convection phenomena. Further, the parameters recovered from this approach do not necessarily have direct meaning with respect to the governing equations for convective transport. In our approach, a framework is presented for directly measuring permeability in the sense of Darcy flow through porous tissue. Measurements from our approach are compared to an extended Tofts model as a control. We demonstrate that a pixel-wise iterative clustering algorithm may be applied to reduce the parameter space of the measurements. We show that measurements obtained from our approach are correlated with measurements obtained from the extended Tofts model control. These correlations demonstrate that the proposed approach contains similar information to an established compartmental model and may be useful in providing an alternative theoretical framework for parameterizing mathematical models for treatment planning and diagnostic studies involving nanomedicine where convection dominated effects are important.

Keywords:

• Modeling (i.e., heat transfer, ultrasound, EM, integrated, treatment planning); mass transport; mixture theory; Darcy law

1. Introduction

The unique physiochemical and optical properties of nanoparticles make them ideal candidates for numerous biomedical applications including drug delivery, bioimaging, tissue engineering, and biosensors. Factors that influence the tumoral accumulation of nanoparticles remain an active area of interest [1,2]. Many delivery strategies for particulate systems have relied on the enhanced permeability and retention (EPR) effect for particles to passively cross the tumor endothelial barrier [3]. A common approach to improve the delivery of nanoparticles into the target site is by prolonging the retention time in the blood [4,5]. Relatively smaller size nanoparticles (<10 nm) are cleared by renal excretion directly and very large particles (>100 nm) are absorbed by liver cells [6]. Nanoparticles between 10 and 100 nm can pass through the endothelial fenestrae in the liver and spleen and may reach the target organ via the blood vasculature. The relationship between nanoparticle circulation time and retention in tumors is complex [7]. Systemic delivery of nanoparticle accumulation in a target tumor results from both passive transport within the convection dominated effects of blood flow within the large vasculature as well as the active biological selection and intra-endothelial transport once the particles are near the target. Recent studies have demonstrated that passive transport through inter-endothelial gaps is not responsible for the transport of nanoparticles into solid tumors [8]. Macrophages [9], as well as dendritic cells, neutrophils, and monocytes [1], have a significant role in the nanoparticle accumulation in tumors. This manuscript presents a pharmacokinetic model of nanoparticle transport within the scope of passive transport.

In particular, the pharmacokinetic model presented views the nanoparticle delivery problem as convective transport through a porous media. Imaging provides anatomical information that drives the clinical workflow and provides a natural setting for model developments and practical application. Here, model developments are with respect to computed tomography (CT). CT imaging is a widely accepted tool for use in a variety of medical imaging applications. Compared to other common imaging modalities such as magnetic resonance imaging (MRI), ultrasound (US), and positron emission tomography (PET), CT is regarded as being accessible, inexpensive, and well-tolerated by imaged subjects. The high spatial and temporal resolution makes CT a common first-line imaging modality for diagnosis and staging of primary and metastatic tumors, anticancer therapy prediction and delivery, and monitoring recurrence following therapy [10,11]. The wide availability of standardized imaging protocols enables reproducible results between scanners and institutions. In CT imaging, ionizing radiation is used to penetrate tissue with inherent contrast arising from differences in attenuation specific to each tissue type. This attenuation, expressed in Hounsfield units (HU), can be artificially increased with the addition of high-density contrast material. Measured HU values are directly proportional to the concentration of contrast in the defined area and provide the basis for quantitative CT applications [12]. While CT imaging can be performed with or without contrast, the use of contrast enables quantitative dynamic and functional imaging when used with techniques such as CT perfusion imaging. Tracer kinetic models of imaging contrast motivate our modeling approach. Contrast enhancement observed during CT perfusion imaging is used as a surrogate for tracking nanoparticle transport within this manuscript. Imaging measurements provide the data to calibrate our mathematical formulation.

High fidelity mathematical modeling of nanoparticle transport kinetics involves multiphysics fluid flow within large vessels and porous living tissue. Common mathematical formulations [13–17] consider the uptake and exchange of nanoparticles among major anatomical compartments including the liver, blood plasma, and tumor. For example, the plasma compartment represents particles in circulation and available to bind to the tumor. The liver compartment represents particles trapped by the liver to be excreted. Analysis of nanoparticle kinetics is commonly developed from first principle mass transport theory and results in a system of nonlinear ordinary differential equations describing the time-varying kinetics of the nanoparticle biodistribution. Such models incorporate several biological parameters that necessitate systematically-designed experimental data to calibrate these model parameters.

The widely applied Tofts model is a compartmental model that was designed for tissues with negligible blood volume [18]. Tofts models assume equilibrium of contrast media between the blood plasma and the extravascular–extracellular space (EES) as well as isodirectional permeability [19]. The notion of permeability in the Tofts’ approach is interpreted as the rate transfer constant between a blood vessel and the EES. The feeding vessels within the tissue are assumed to provide a spatially homogeneous arterial input function source term to the governing ordinary differential equation. Implicitly, this assumes that the time scale for the transport between imaging visible vessels and vessels not visible on imaging is less than the sequential time between a dynamic imaging acquisition. The extended Tofts models build upon the Tofts model and include additional parameters for intravascular signal contributions. However, permeability is still interpreted as the rate transfer constant between a blood vessel and the EES. Fundamentally, the governing equations for compartmental mass balance represent a spatially invariant mass transport between compartments and do not capture spatially variant convection phenomena.

In this manuscript, we view the transport problem as convective flow through porous tissue. While Darcy flow models have been applied to model mass transport in porous biological tissues [20–27], the key idea of this manuscript is to provide a framework for directly measuring permeability in the sense of Darcy flow through porous tissue. In the case where the porous media does not provide an imaging signal, NMR measurements may be applied to directly measure permeability and porosity [28] from an MR visible fluid. Here, the porous biological tissue provides an imaging signal. Obtaining permeability in the sense of Darcy flow is important for applications of nanotechnology in which the governing equations are dominated by convective transport [29,30]. Permeability measurements within the presented framework are obtained from the same convection dominated equations. Dimensional analysis of the parameters obtained are directly related to units of flow rates, pressure gradients, permeability, and, thus, have direct meaning for parameterizing mathematical models for treatment planning and diagnostic studies. Pixels parameter recovery results in high dimensional inverse problems that are computationally expensive [31]. An iterative clustering algorithm is applied pixel-wise and evaluated for data reduction [32].

2. Methods

2.1. Signal model

A mathematical formulation for mass transport within vascularized porous living tissue shares much in common with the mathematics, aerospace, biomedical, geosciences, and petroleum engineering mixture theory literature [33–36]. Mixture theory provides a natural framework to consider nanoparticle transport through porous tissue. At each point in space, $x\in \Omega ,$ we consider the mixture of porous tissue, ${\phi }_{\text{tissue}},$ and blood that does, ${\phi }_{\text{nano}},$ and does not, ${\phi }_{\text{blood}},$ contain the nanoparticles: (1) $${\phi }_{\text{tissue}}+{\phi }_{\text{nano}}+{\phi }_{\text{blood}}=1\forall x\in \Omega .$$(1)

Mass is conserved for each component in the mixture with and without nanoparticles flowing with velocity, $v\left[\frac{m}{s}\right]:$ $$\frac{\partial {\rho }_{\iota }{\phi }_{\iota }}{\partial t}+\nabla ·\left({\rho }_{\iota }{\phi }_{\iota }v\right)=0\mathrm{\hspace{1em}}\mathrm{\hspace{1em}}\iota \in \left\{\text{nano},\mathrm{ }\text{blood}\right\}.$$

Mass density of each component is denoted ${\rho }_{\iota }\left[\frac{kg}{m^3}\right].$ We rewrite the mass balance in terms of the saturation: (2) $$\frac{\partial s_{\iota }}{\partial t}+\nabla ·\left(s_{\iota }v\right)=0\mathrm{\hspace{1em}}\mathrm{\hspace{1em}}\iota \in \left\{\text{nano},\mathrm{ }\text{blood}\right\},$$(2) (3) $$s_{\text{nano}}+s_{\text{blood}}=1\mathrm{\hspace{1em}}\mathrm{\hspace{1em}}{\phi }_{\text{tissue}}=1-\varphi \mathrm{\hspace{1em}}\mathrm{\hspace{1em}}{\phi }_{\text{nano}}=\varphi s_{\text{nano}}\mathrm{\hspace{1em}}\mathrm{\hspace{1em}}{\phi }_{\text{blood}}=\varphi s_{\text{blood}},$$(3)

where $\varphi$ represents the tissue porosity or void space. Saturation of a given constituent, $s_{\iota },\iota \in \left\{\text{blood},\mathrm{ }\text{nano}\right\},$ represents the ratio of the void volume filled with the constituent to the total of the void volume in the porous medium. Summing the mass balance (2) for each constituent, $\iota \in \left\{\text{nano},\mathrm{ }\text{blood}\right\},$ and applying the constraint that the saturations sum to one,(3), the flow is seen to be incompressible: (4) $$\nabla ·v=0.$$(4)

Characteristic lines of convective transport are assumed in the direction of the normal to imaging visible vessels, n, that are the source of the flow. The speed of the flow along the characteristic lines is denoted a. The governing equation for tracking nanoparticle transport is a first order hyperbolic equation: $$\begin{array}{c}\frac{\partial s_{\text{nano}}}{\partial t}+\frac{\partial s_{\text{nano}}}{\partial n}a=0\mathrm{\hspace{1em}}\mathrm{\hspace{1em}}a=n·v=\left|\right|v\left|\right|\\ s_{\text{nano}}(0,x)=s_0(x)\mathrm{\hspace{1em}}\mathrm{\hspace{1em}}{s}_{\text{nano}}(t,0)=s_{\text{AIF}}(t)\end{array}.$$

Here, the initial condition s₀ may be obtained directly from the imaging data and is the only source of nanoparticle transport. The inlet boundary condition from the arterial input function, $s_{\text{AIF}},$ is obtained by placing a region of interest (ROI) on the contrast enhancement observed in the aorta. The vector pointing from the vessel centerline implicitly defines the normal vector, n, from the centerline point of the closest feeding vessel, x₀: $$n=\frac{x-x_0}{\left|\right|x-x_0\left|\right|}.$$

A signed distance map [37] is used to estimate the closest perpendicular distance to the vessel input. This perpendicular distance approximation implicitly assumes that both imaging visible and non-visible vessels contribute to the flow and that vessel branches exist that are perpendicular to the imaging visible vessels. The feeding vessels within the tissue are assumed to supply the nanoparticle system input as initial conditions and boundary conditions for our system. Under these assumptions, the solution is analytic and of the form $s(x,t)=f(x-at).$ The function, f, is arbitrary and is determined by the initial conditions and boundary conditions: (5) $$s_{\text{nano}}(t,x)=s_0\left(x-at\right)=s_{\text{AIF}}\left(-x/a+t\right).$$(5)

The flow speed, $a[m/s],$ is obtained as the distance value, $\left|\right|x-x_0\left|\right|,$ divided by the contrast bolus arrival time (BAT). The bolus arrival time is computed from the peak-gradient method available in the 3 D Slicer PkModeling extension [38]. Flow speed is assumed piecewise constant to allow for tissue heterogeneity: (6) $$a(x)=\sum _i{\alpha }_i{\psi }_i(x).$$(6)

A superpixel model [32] is used to obtain piecewise constant imaging regions. A linear iterative clustering algorithm with a grid size of 20 voxels is used to generate the superpixel model.

2.2. Flow model

Following Darcy’s law, the flow speed a is the assumed result of a pressure gradient along the flow direction: (7) $$v·n=-\frac{\kappa }{\mu }\nabla p·n=a.$$(7)

The flow velocity, v, is proportional to the pressure gradient, p, across the tissue and in the direction of the normal, n, from the source vessel. Tissue permeability is denoted $\kappa [m^2].$ The viscosity of blood is known, $\mu =3.5·{10}^{-3} [Pa·s].$ Poiseuille’s law [39] provides a relationship between the flow rates and pressure gradient in the vessels: (8) $$Q=\frac{\pi R^4}{8\mu }\nabla p·n.$$(8)

Flow rates, Q, in rabbit liver have been measured as $177[ml/\mathit{min}]$ [40]. The average vessel radius, R, was extracted from image Hessian based filters applied to the image subtraction of the pre-contrast image from the maximum intensity projection image [39]. The pressure gradient is assumed spatially homogeneous across the tissue and vessels. Substituting Equation (8)(8) $$Q=\frac{\pi R^4}{8\mu }\nabla p·n.$$(8) into Equation (7)(7) $$v·n=-\frac{\kappa }{\mu }\nabla p·n=a.$$(7) provides an estimate of the tissue permeability: (9) $$\kappa =\frac{a\pi R^4}{8Q}.$$(9)

Within this manuscript, Equation (9)(9) $$\kappa =\frac{a\pi R^4}{8Q}.$$(9) estimates the permeability as proportional to the flow speed. However, for further comparison with existing literature, it is worth noting that alternative measurements of the flow rate, Q, may be obtained from indicator dilution theory [41] and blood volume measurements [42]. Flow rates are derived from the conservation of flow across the tissue and vessel interface. The flow rate input is equal to the flow rate out of the tissue. A constant velocity over each interface, $\partial \Omega ,$ is assumed. The flow rate, Q is derived from the integral form of mass conservation applied to constituent ι at the tissue interface: (10) $$\begin{array}{c}{m}_{\iota }(T)-m_{\iota }(0)={\int }_{\Omega }\rho \phi (x,T)dV-{\int }_{\Omega }\rho \phi (x,0)dV={\int }_0^T{\int }_{\partial \Omega }\rho \phi v·n\mathrm{\hspace{1em}}\mathrm{\hspace{1em}}\iota \in \left\{\text{nano},\mathrm{ }\text{blood}\right\}\hfill \end{array}.$$(10)

Under the assumption that the time interval t = 0 represents the time prior to the nanoparticle arrival and t = T is less than the transit time out of the tissue Ω, outflow is not considered. Within the time interval $[0,T],$ the flow rate into the tissue may be obtained: (11) $$m_{\text{nano}}(T)={\int }_0^T\rho {\phi }_{\text{nano}}dt{\int }_{\partial \Omega }{v}_{\text{in}}·\mathit{ndA}\mathrm{\hspace{1em}}\mathrm{\hspace{1em}}\Rightarrow \mathrm{\hspace{1em}}\mathrm{\hspace{1em}}Q={\int }_{\partial \Omega }{v}_{\text{in}}·\mathit{ndA}=\frac{m(T)}{{\int }_0^T\rho {\phi }_{\text{nano}}dt}.$$(11)

Here, steady state input velocity, $v_{\text{in}},$ is assumed so that the velocity may be factored out of our the integral in time. The tracer distribution is assumed homogeneous across the feeding vessel such that the concentration factors out of the spatial integrals. At this step, the various analysis methods depart in assumptions for calculation of the flow rate, Q. Indicator dilution theory [41] assumes that the mass in the tissue is known from the injected concentration bolus, $m(T)=m_0.$ Calculations of blood volume (BV) measurements use the area under the curve for the concentration in the tissue [42]. Dimensional analysis of the units of each term are used to guide the algebraic relationships seen in the literature. The average flow rate, $Q,$ may be related to BV as (12) $$Q\equiv \frac{1}{\Delta t}{\int }_0^{T}Q=\frac{1}{\Delta t}{\int }_{\Omega }BV(x)dV\mathrm{\hspace{1em}}\mathrm{\hspace{1em}}BV(x)\equiv \frac{{\int }_0^T{\phi }_{\text{nano}}(x,t)dt}{{\int }_0^T\rho {\phi }_{\text{nano}}dt}.$$(12)

2.3. Compartmental model

In the assessment of tumor biology and angiogenesis, compartmental analysis of CT perfusion imaging provides a noninvasive volumetric assessment of tumor vasculature [11,43–46]. Correlation between parametric maps obtained from compartmental analysis and microvessel structure has been used as a predictor of disease severity and treatment outcomes [47–50].

We apply the PkModeling extension for 3 D Slicer [51] to compute the extended Tofts compartmental model parameters as a control for our speed recovery in Equation (5)(5) $$s_{\text{nano}}(t,x)=s_0\left(x-at\right)=s_{\text{AIF}}\left(-x/a+t\right).$$(5) : (13) $$s_{\text{nano}}(t)=K^{\text{trans}}{\int }_0^t{s}_{\text{aif}}(u)e^{-\frac{K^{\text{trans}}}{v_e}(t-u)}du+\mathit{fpv}{s}_{\text{aif}}(t).$$(13)

The solution to the compartmental model was obtained from a nonlinear least square curve fit to the extended Tofts model analytical solution parametrized by $K_{\text{trans}},$ v_e, and fpv. Here, $K_{\text{trans}}[1/s]$ represents the influx forward volume transfer constant (into EES from plasma), v_e is the fractional volume of EES per unit volume of tissue, fpv represents the fractional plasma volume of the arterial input present at each voxel [19]. Contrast enhancement is linearly proportional to the concentration of contrast agent in the tissue [42]: (14) $${\phi }_{\text{nano}}\propto I(x,t)-I_b(x),$$(14)

where$I_b(x)$ represents a baseline value of the tissue before contrast enhancement. The correlation between the Tofts model parameters and speed a estimates are measured.

2.4. Experimental methods

Animal procedures were performed following a protocol approved by the MD Anderson Cancer Center Institutional Animal Care and Use Committee per recommendations in the National Institutes of Health Guide for Care and Use of Laboratory Animals. Four male New Zealand White rabbits (Charles River Laboratories, Massachusetts, USA) weighing 2.4-3 kg underwent VX2 tumor inoculation in the left lateral lobe of the liver. The four animals are identified as id1, id2, id3, and id4 in our analysis. Tumor implantation was performed via the administration of freshly harvested and prepared VX2 tumor fragments from an in-house VX2 donor line. VX2 fragments 3-4 mm in diameter were implanted in the liver at a single site through an 18 GA needle attached to a 1 ml syringe. Tumors were allowed to grow for 10-14 days until the tumors reached approximately 1 cm diameter. Following sufficient tumor growth, all animals underwent a sham procedure involving the insertion and removal of a catheter (20 GA x 1.88 in). CT imaging occurred at the time of the sham procedure and 24 h post-procedure.

Before imaging, rabbits were anesthetized with isoflurane (1-5%)/oxygen (1.5 L/min) administered via endotracheal tube and laryngeal mask. CT imaging was performed using a multi-slice detector CT (SOMATOM Definition Edge, Siemens Healthineers, Erlangen, Germany) to determine contrast transport through the liver and tumor. A vac fix immobilization device was used to maintain supine positioning on the scanner table. Perfusion imaging was acquired with 4, 6, or 8 ml of 320 mg/ml iodixanol contrast (Visipaque, GE Healthcare, Cork, Ireland) injected at a rate of at 2 ml/sec into the marginal ear vein using a power injector (Medrad Envision, CT, Connecticut, USA). Scan parameters were 80 kVp, 264 mA, 570 ms rotation time, and pitch of 0.5. Scans were acquired over the entire liver before contrast injection began, and at delays of 2, 4, and 6 s following the injection start; the total scanning time was 87 s. A 15 cm field-of-view was reconstructed from the raw projection data with a slice thickness of 1.5 mm and an interval of 1.0 mm with a B20f reconstruction kernel.

3. Results

Transport of iodinated contrast agent is used as a surrogate for nanoparticle transport in the tissue. Datum used in the in vivo evaluation of our model is shown in Figure 1. Figure 1(a) illustrates the time history of the contrast uptake at the aorta used for as the arterial input, $s_{\text{AIF}}.$ Regions of interest used in the data analysis are outlined in Figure 1(b). Imaging visible vessels are outlined in red. Pixel-wise distance is measured from the perpendicular distance of the closest vessel outlined. The blue region outlines the VX2 tumor growth. Green outlines liver parenchyma tissue.

Figure 1. In vivo data. (a) The time history of the contrast in the aorta illustrates the bolus of contrast arrival time. The different lines illustrate three representative pixels. (b) The labeled regions illustrate the imaging visible vessels (red), the tumor (blue), and the liver parenchyma tissue (green). Pixelwise distance from the nearest vessel boundary provides the spatial length scale for our governing hyperbolic equation for convective transport.

Linear iterative clustering segmentation of the rabbit liver resulted in 174, 176, 111, and 156 super-pixel regions within the livers of animal id1, id2, id3, and id4, respectively. Two iterations of iterative clustering were performed in each data set for convergence. The average volume of each superpixel was 960 ± 26 mm³. A representative example of the pixelation of the liver is shown in Figure 2. Super-pixel regions form the piecewise constant basis for inversion of the tissue heterogeneity.

Figure 2. Domain decomposition of liver anatomy with simple linear iterative clustering (SLIC) segmentation. An (a) axial view and (b) 3 D visualization of the domain decomposition is shown.

Figure 3 provides a visualization the data input to our algorithm as well as a visualization of the extended Tofts model results. Figure 3(a) illustrates spatial variations in the fractional plasma volume, fpv. Fractional plasma volume is a dimensionless quantity $[1].$ Negative values of the fractional plasma volume are seen in the curve fit of the extended Tofts model (13) to the data. Figure 3(b) illustrates the pixelwise bolus arrival time in seconds. Voxel neighborhoods of pixels are seen to have the same arrival time. Figure 3(c) illustrates surface of the imaging visible vessels. Flow speed is estimated as the distance to the nearest vessel surface divided by the arrival time. Correlations between the extend Tofts model parameters and our measurement flow speed are summarized in Table 1. The fractional plasma volume, fpv, extended Tofts parameter shows a positive correlation with our flow speed measurements for all animals.

Figure 3. An overview of our algorithmic approach is shown for animal id1. (a) The extended Tofts model provides a control for the proposed approach. Spatial variations in the tissue fractional plasma volume parameter fpv within the parenchyma and vessels is shown. Negative values are the result of the curve fit of the analytical model in Equation (13)(13) $$s_{\text{nano}}(t)=K^{\text{trans}}{\int }_0^t{s}_{\text{aif}}(u)e^{-\frac{K^{\text{trans}}}{v_e}(t-u)}du+\mathit{fpv}{s}_{\text{aif}}(t).$$(13) to the time-intensity data at each pixel. (b) Spatial variations in the bolus arrival time is measured in seconds. (c) Distance information is computed with respected to imaging visible vessels. An example of the 3 D vessel segmentation used to calculate the signed distance transform $x-x_0$ is shown. (d) The bolus arrival time combined with the distance information provide estimates of the flow speed. Correlation between the super pixel recovered speed and the average of the pixelwise is shown. Within these pixels the average speed represented by the super pixel is expected to provide dimension reduction of the parameter optimization space while maintaining accuracy of the solution field.

Figure 4. Tofts model, fpv, versus flow speed estimates are provided. The fractional plasma volume is a dimensionless $[1]$ parameter. Animal id1, id2, id3, and id4 are shown in (a), (b), (c), and (d); respectively. The correlation, r, between the fractional plasma volume and flow speed is provided within each dataset.

Table 1. Summary of extended Tofts analysis.

	id1	id2	id3	id4
$r(a,K_{\text{trans}})$	–0.27	–0.32	0.35	–0.05
$r(a,v_e)$	–0.16	–0.30	0.23	–0.38
r(a, fpv)	0.43	0.42	0.73	0.41

The observed correlation, r, between the measured flow speed (6) and the $K_{\text{trans}},$ v_e, fpv, extended Tofts model parameters (13) are shown for each animal: id1, id2, id3, and id4.

Figure 3(d) plots the average pixel flow speed against the flows speed recovered directly from superpixel domain decomposition for animal id1. The x-axis represents the pixelwise speed averaged across a given super pixel illustrated in Figure 2. The y-axis represents the speed recovered from our model fits applied to the raw data averaged across a superpixel. A Pearson correlation of 0.95 (p <.05) was observed and indicates that the superpixel domain decomposition is a good approximation of the tissue heterogeneity. Table 2 provides a global summary of the pixelwise average permeability estimate with each animal liver, id1, id2, id3, and id4. Permeability was calculated according to Equation (9)(9) $$\kappa =\frac{a\pi R^4}{8Q}.$$(9) . The average vessel radius extracted from the surface of the imaging visible vessels shown in Figure 3(c) was 1.09 mm, 1.07 mm, 1.10 mm, and 1.09 mm for animal id1, id2, id3, and id4; respectively. Figure 4 provides a comprehensive quantitative summary of the correlations measured between the pixelwise average speed and the extended Tofts fpv parameter

Table 2. Summary of permeability measurements.

	id1	id2	id3	id4
$\kappa [m^2]$	1.74e–07	1.72e–07	2.83e–07	2.79e–07

The pixelwise average permeability estimates with the liver parenchyma tissue is shown for each animal: id1, id2, id3, and id4.

4. Discussion

Families of tracer kinetic models may be grouped as model-independent based on Fick’s Law, deconvolution analysis, compartmental modeling, and modeling that accounts for convective transport [42]. Each has assumptions appropriate for domain-specific applications. Here we assume that convective transport from the vessels directly supplies the flow within the tissue. Importantly, our approach allows us to recover permeability properties with respect to porous media flow. Table 2 presents initial steps and demonstrates the feasibility in the recovery of convection related properties from imaging data. As a reference, permeability measurements on the order of $1e^{-7}{m}^2$ are similar to gravel or fractured rock within the oil and gas communities [52]. Tofts model proposes a spatially invariant compartmental formulation. The resulting measurements of the Tofts formulation have different units, [1/s] for $K_{\text{trans}}$ and dimensionless for v_e and fpv, with respect to transport equations within porous tissue. Calibrating and validation of the same governing equations have the benefit of direct interpretability. Although our approach applied the model to contrast-enhanced imaging, imaging visible measurements of nanoparticle concentrations may be used to directly recover permeability measurements of interest for transport within porous tissue.

This manuscript presents a pharmacokinetic model of nanoparticle transport within the scope of passive transport. Permeability estimates (Table 2) are with respect to contrast flow through porous tissue. Antibody targeting and interactions with the immune system [1] are not considered in this approach. However, the key ideas of the framework may also be coupled with mechanisms for active transport of nanoparticles into the target tumor. This may be modeled as a sink term with appropriate rate constants within our governing equations. For example, within a mixture theory framework, the concentration of immune cells may be included as an additional constituent within Equation (1)(1) $${\phi }_{\text{tissue}}+{\phi }_{\text{nano}}+{\phi }_{\text{blood}}=1\forall x\in \Omega .$$(1) . Rate constants for the uptake of the nanoparticles by the immune system formulates the interaction as a first-order relationship between the presence of nanoparticles and uptake by the immune system. Additional model complexity may incorporate nanoparticle interactions with blood cells, plasma proteins and complement proteins. Serum proteins present in the bloodstream bind to the nanoparticles and results in their rapid recognition and uptake of nanoparticles by the mononuclear phagocyte system. Premature clearance of nanoparticles by the mononuclear phagocyte system due to their size and nonspecific binding limit the delivery efficacy of nanoparticles into a specific target site [4]. Additional constituents couples additional equations and parameters to our governing equations. Solution of the governing equations coupled with immune cells share much in common within the field of tumor growth modeling [53]. These models are based on the conservation laws of continuum physics complemented by source and flux characterizations that trigger and control cancer development and decline through mathematical representations of the events listed in the Hallmarks of Cancer [54].

The presented methodology is applicable to flow rates characteristic of systemic delivery as well as convection-enhanced delivery. Convection enhanced delivery applied to nanomedicine utilizes direct injections of the nanotechnology into the interstitial space to achieve the desired therapeutic or diagnostic effect. Planning software for convection-enhanced delivery is likely to be useful in guiding in deciding catheter placement [26]. Our parameter recovery approach provides a methodology to recover patient-specific transport parameters with direct meaning with respect to convection-dominated governing equations. Either a literature value for pressure may be used to recover the permeability, or the framework may be applied with patient-specific blood pressure measurements. The direct correlation between flow speed and permeability is expected to provide a reliable reproducible methodology for Darcy’s Law based predictions.

Limitations in these measurements arise from the tradeoffs in computational efficiency and model fidelity. The analytical solution in Equation (5)(5) $$s_{\text{nano}}(t,x)=s_0\left(x-at\right)=s_{\text{AIF}}\left(-x/a+t\right).$$(5) facilitates an efficient numerical solution but is inherently one-dimensional from the closest perpendicular distance to the nearest vessel source. Sufficiently near a vessel source, this approximation is likely to be appropriate. Future efforts will evaluate tradeoff of increase model fidelity. Finite element methods are needed to curve fit the 4 D enhancement data to the partial differential equations for mass balance, Equation (2)(2) $$\frac{\partial s_{\iota }}{\partial t}+\nabla ·\left(s_{\iota }v\right)=0\mathrm{\hspace{1em}}\mathrm{\hspace{1em}}\iota \in \left\{\text{nano},\mathrm{ }\text{blood}\right\},$$(2) , coupled with incompressible flow, Equation (4)(4) $$\nabla ·v=0.$$(4) . The finite element method further requires a computation mesh that accurately conforms the tissue anatomy. Adjoint methods [31] are needed to provide analytic gradients with respect to the domain decomposition for efficient model calibration to the data. The presented results are a theoretically sound simplification of the 4 D partial differential equation constrained optimization using a finite element approach. Here, we hypothesize that any solution bias induced by our simplified approach may be reconciled with appropriate calibration and cross-validation to imaging data. Thus, the overall correlation between the recovered permeability and permeability of a high fidelity finite element approach is anticipated to be high but with much less cost.

Future work will also incorporate numerical simulations of the flow through the vessels into the transport problem. The problem of particle flow through the vasculature has been considered for varying particle diameters. Navier–Stokes computational fluid dynamic simulations are considered for flow through the vasculature. For large particle size relative to the vessels, perturbations of the fluid flow from the presence of the particles should be considered. However, for nanomedicine, the particles are small relative to capillary diameters ($\approx$ 10 µm). The fluid force on the particle should be considered, however, the inverse particle momentum exchange is not considered. This one-way coupling approach has been shown to save significant computational resources while providing accurate results in this small nanoparticle to vessel diameter setting [55]. For systemic delivery, tracking of flow streamlines from the fluid flow simulations provides a mechanism for downstream targeting of the desired delivery [56]. Runge-Kutta methods are used to solve the differential equations with a varying time step. Ultimately, for a given application and data quality, multiple models may be considered feasible to describe the data. Future efforts will apply Bayesian model selection approaches to a family of models to balance model complexity with the goodness of fit [18].

Disclosure statement

No potential conflict of interest was reported by the author(s).

Additional information

Funding

This work was supported in part by the National Institutes of Health [grant number R01CA201127]. AC Crouch supported by a training fellowship from The University of Texas Health Science Center at Houston Center for Clinical and Translational Sciences TL1 Program [grant number TL1 TR003169]. EA Thompson was supported by the Society of Interventional Radiology Foundation Allied Scientist Training Grant and CPRIT Research Training Grant [grant number RP170067]. The support of the High Performance Computing for research facility at the University of Texas MD Anderson Cancer Center for provided computational resources that have contributed to the research results reported in this paper.