Sequential deconstruction of composite drug transport in metastatic breast cancer

See allHide authors and affiliations

Science Advances  24 Jun 2020:
Vol. 6, no. 26, eaba4498
DOI: 10.1126/sciadv.aba4498


It is challenging to design effective drug delivery systems (DDS) that target metastatic breast cancers (MBC) because of lack of competent imaging and image analysis protocols that suitably capture the interactions between DDS and metastatic lesions. Here, we integrate high temporal resolution of in vivo whole-body PET-CT, ex vivo whole-organ optical imaging, high spatial resolution of confocal microscopy, and mathematical modeling, to systematically deconstruct the trafficking of injectable nanoparticle generators encapsulated with polymeric doxorubicin (iNPG-pDox) in pulmonary MBC. iNPG-pDox accumulated substantially in metastatic lungs, compared to healthy lungs. Intratumoral distribution and retention of iNPG-pDox varied with lesion size, possibly induced by locally remodeled microenvironment. We further used multiscale imaging and mathematical simulations to provide improved drug delivery strategies for MBC. Our work presents a multidisciplinary translational toolbox to evaluate transport and interactions of DDS within metastases. This knowledge can be recursively applied to rationally design advanced therapies for metastatic cancers.


Metastatic breast cancer (MBC) is refractory to traditional anticancer therapies and is presently incurable (1, 2), necessitating innovation in drug discovery and delivery of chemotherapeutics. Owing to their small size, heterogeneity, and dispersed nature, MBCs are difficult to detect and treat, accounting for >90% of breast cancer–related deaths (3, 4). Despite notable advances in engineering and development of complex macromolecular and nano- and micrometer-sized drug delivery systems (DDS) for primary solid tumors, only limited studies have evaluated the delivery, accumulation, and interactions of such systems in metastatic settings. A series of biological barriers, such as entrapment in reticuloendothelial system organs, hemorheological and hemodynamic limitations, endothelial extravasation, and tumor cell membrane, hinders the accumulation of drugs in tumors and reduces the efficacy of conventional small molecule– and nanoparticle (NP)–based therapeutics (5). Multicomponent and multifunctional NPs have been developed to improve therapeutic efficacy while reducing nonproductive systemic distribution by introducing active targeting moieties or pH-activatable drug release (6, 7). However, these systems have shown limited efficacy, and none have been translated to the clinic. We have demonstrated sequential and successful negotiation of most, if not all, biological barriers by iNPG-pDox composed of porous silicon microparticle-based injectable nanoparticle generator (iNPG) that was packaged with poly(lactic-co-glycolic acid) polymer-doxorubicin (pDox) conjugates. pDox molecules self-assemble into NPs upon release from iNPG, and iNPG-pDox treatment significantly inhibited tumor metastasis, including functional cure in 40 to 50% of the treated mice with pulmonary MBC (8). While such tremendous outcomes are unprecedented in MBC, it is imperative that their physiological voyage be traced and a deeper mechanistic understanding of their transport processes be gained to aid and expedite their preclinical testing and clinical translation.

Molecular imaging modalities can elucidate a wealth of disease-relevant information across length and time scales (9). Although relatively under-explored, multiscale molecular imaging, when applied rationally and integrated across subcellular to whole-body levels, is particularly well suited to dissect the pharmacokinetics (PK) of administered therapies and hence predict and monitor their therapeutic responses and efficacy (10). As such, in conjunction with mathematical modeling of mass and momentum transport through the body, molecular imaging modalities form an essential pillar of the transport oncophysics approach to cancer (11, 12).

In this work, we demonstrate the use of positron emission tomography–computed tomography (PET-CT), combined with optical imaging (OI), high-resolution microscopy, and mathematical modeling, in deconstructing the biodistribution kinetics of composite drug, iNPG-pDox, in a quantitative, stepwise manner from whole body to organ and finally to tissue and cellular levels in a murine model of pulmonary MBC (Fig. 1). We adopted a modular, orthogonal surface engineering approach that is independent of carrier type, cargo, formulation method, or delivery route to modify iNPG-pDox for subsequent PET-CT and OI/microscopy, respectively. The early events in the systemic distribution of iNPG-pDox were monitored up to 3 hours postinjection (p.i.) with high spatiotemporal resolution via longitudinal PET-CT. iNPG-pDox demonstrated passive tumor tropism, accumulating rapidly and substantially in metastatic lungs. Further, we integrated in vivo imaging with mathematical modeling to mechanistically describe the PK of iNPG-pDox. The semimechanistic model, based on the conservation of mass and law of mass action (13, 14), captured the essential transport processes and biophysical interactions that govern the in vivo disposition of the particles and their differential kinetics in the lungs of naïve and diseased animals. Subsequently, multiplexed, volumetric OI of intact tumor-bearing and naïve lungs at 24 hours p.i. revealed an apparent divergence in trafficking patterns of iNPG and pDox that varied with the size of metastatic lesions. Spatially heterogeneous distribution of the two components, afforded by the release of pDox NP from the iNPGs, was further explored with high-resolution confocal microscopy. While the iNPGs tended to marginate in the tumor-associated pulmonary blood vessels, pDox NPs were found to infiltrate into the metastatic lesions.

Fig. 1 Schematic of the experimental workflow combining multiscale imaging modalities to deconstruct the transport kinetics of composite DDS, iNPG-pDox.

Whole-body longitudinal PET-CT allowed visualization of spatiotemporal distribution of [64Cu]NOTA-iNPG-pDox in murine models of 4T1 pulmonary metastasis. Nondestructive whole-organ multiplexed OI on excised lungs revealed intra-organ trafficking pattern of dual-fluorescent AF647-iNPG-pDox. High-magnification, high-resolution microscopy delineated the spatial heterogeneity in iNPG and pDox distribution at the cellular level. Representative images in the bottom panel illustrate the data collected sequentially during the experimental workflow.

We believe that this study provides an evolving toolbox of experimental and theoretical approaches that link clinically relevant imaging modalities such as PET-CT with orthogonal experimental approaches such as whole-organ OI and high-resolution microscopy and complementary mathematical modeling, to aid the development, characterization, and assessment of high-performance experimental and investigational new therapeutics for metastatic cancers in preclinical and translational studies. We contend that such interdisciplinary retrospective analyses not only enable multiscale assessment of a therapeutic’s behavior in vivo but also make reverse engineering possible to guide future rational design of novel and effective anticancer drugs. Being agnostic to therapeutic class or morphology, administration route, animal species, or disease type, our approach can be easily customized for a wide range of studies spanning preclinical to clinical spheres. The strategy can be harnessed for rapid screening and comparison among potential therapeutic formulations, thus promising immense implications in preclinical and clinical settings.


Engineering multistage delivery vectors for multimodal imaging

iNPG-pDox were surface-engineered in a stepwise manner to render them amenable to multimodal imaging spanning multiple length scales. In preparation for longitudinal whole-body PET-CT, amine-modified, discoidal porous silicon particles (called iNPG), with dimensions 2.6 μm diameter by 700 nm thickness and 40- to 80-nm-sized pores (15, 16), were first conjugated to p-SCN-Bn-NOTA (1,4,7-triazacyclononane-N,N′,N″-triacetic acid), a bifunctional divalent macrocyclic metal chelator via amine-isothiocyanate chemistry (Fig. 2A) (17). The conjugation strategy was iteratively refined to induce minimal alteration in vector morphology and PK behavior while ensuring sufficient signal-to-noise ratio. The probes were labeled with a PET radionuclide, copper-64 (64Cu, half-life; t1/2 = 12.7 hours) (17), with ~58% radiolabeling yield and ~99% radiochemical purity (fig. S1B). Intermediate probes NOTA-iNPG and [64Cu]NOTA-iNPG depicted similar physicochemical characteristics (size, shape, and surface charge) to native iNPGs (Fig. 2, B and C, and fig. S1A). Before injections, [64Cu]NOTA-iNPG was loaded with pDox monomers, with ~25% encapsulation efficiency [matching that of native iNPG-pDox (8)], and the final constructs, [64Cu]NOTA-iNPG-pDox (Fig. 2B), depicted no major changes in morphology (fig. S1A) or surface charge (Fig. 2C) when compared to iNPG-pDox. Radiostability of [64Cu]NOTA-iNPG-pDox was investigated in vitro before systemic studies. 64Cu isotope remained stably chelated to the probe, with only minor detachment (~93.5% intact 64Cu) after 6 hours of incubation in mouse serum (Fig. 2D and fig. S1C). Broadening of the analyte peak in radio-chromatogram 24 hours after incubation (fig. S1C) was attributed to gradual degradation of iNPG-pDox (fig. S1E), consistent with previous observations with porous silicon vectors (15).

Fig. 2 Surface engineering and characterization of iNPG-pDox.

(A) Stepwise illustration of [64Cu]NOTA-iNPG-pDox synthesis. Partially oxidized and amine-modified silicon microparticles were conjugated to a radiometal chelator (NOTA-SCN), followed by radiolabeling with 37 MBq 64CuCl2 (outset in yellow depicts the conjugation chemistry; a similar protocol was adopted to synthesize AF647-iNPG). [64Cu]NOTA-iNPG was loaded with polymeric doxorubicin (pDox) cargo. (B) Representative scanning electron micrographs depict disc-like, porous morphology and pDox-filled channels in [64Cu]NOTA-iNPG-pDox. Scale bars, 2 μm. (C) Surface charge of the composite DDS at different stages of surface modification. Radiolabeling or NIR dye (AF647) conjugation does not significantly alter the zeta potential of naïve iNPG or iNPG-pDox. (D) Serum stability of [64Cu]NOTA-iNPG-pDox over the 24-hour incubation period obtained by radio–thin-layer chromatography (radio-TLC). (E) Absorbance spectra of dual-fluorescent AF647-iNPG-pDox indicate peaks typical of Dox (exmax: 490 nm) and NIR dye AlexaFluor647 (exmax: 650 nm) absorbance. a.u., arbitrary units.

We used a similar surface engineering approach to introduce near-infrared (NIR) dye AlexaFluor647 (AF647; excitation: 640 nm, emission: 680 nm) on iNPG microcarriers to investigate the intra-organ and intratissue distribution of iNPG-pDox. AF647-conjugated iNPG could be loaded with pDox monomers with 26.8 ± 5.9% encapsulation efficiency, similar to that of native iNPG (8). The final fluorescent probe, AF647-iNPG-pDox, depicted the optical characteristics of both Dox and NIR dye (Fig. 2E and fig. S2, A and B), indicating multiplexing capability for ex vivo OI and confocal microscopy. The conjugation chemistry was intricately tailored to avoid any changes in the physicochemical characteristics of the composite drug (Fig. 2C and fig. S1A), previously optimized to ensure greatest accumulation in metastatic organs, including liver and lungs (8). pDox NPs released from AF647-iNPG-pDox retained the optical characteristics of free Dox drug (fig. S2C). When incubated in serum, AF647-iNPG-pDox demonstrated an initial burst release of ~44% pDox NPs within 3 hours, increasing gradually to ~55 and 75% cumulative release by days 1 and 5, respectively, concomitant with the gradual degradation of iNPG (fig. S2D).

Longitudinal whole-body PET-CT depicts spatiotemporal distribution of [64Cu]NOTA-iNPG-pDox

A series of multiscale in vivo and ex vivo experiments was designed to comprehensively explore the distribution kinetics of complex multicomponent iNPG-pDox. At the whole-body or systemic level, we first investigated the early events in the spatiotemporal distribution of [64Cu]NOTA-iNPG-pDox in murine models of 4T1 MBC. Age- and species-matched unchallenged mice were used as control. Harnessing the exquisite temporal resolution of dynamic PET-CT, a 30-min scan was acquired simultaneously with an intravenous bolus injection of 60 μCi [64Cu]NOTA-iNPG-pDox (corresponding to a dose of 6 mg kg−1 Dox). The image sequence was reconstructed into six frames: 6 × 300 s to visualize the spatial distribution of the particles in real time. Combined with a terminal static PET-CT scan at 3 hours p.i., our PET-CT protocol could adequately map the changes in the trafficking pattern of [64Cu]NOTA-iNPG-pDox from blood pool to metastatic lungs and mononuclear phagocytic system (MPS) organs (Fig. 3).

Fig. 3 Whole-body PET-CT imaging reveals spatiotemporal kinetics of [64Cu]NOTA-iNPG-pDox.

(A and B) Representative coronal images reconstructed from longitudinal whole-body PET-CT scans of 4T1 metastasis-bearing (A) and naïve (B) mice. Imaging protocol included an initial 30-min dynamic scan (5 frames × 300 s) and a terminal static scan at 3 hours p.i. (C and D) 3D-rendered maximum intensity projection (MIP) of the terminal PET-CT scan (blue outline) depicting the spatial distribution of [64Cu]NOTA-iNPG-pDox in 4T1 tumor-bearing (A) and naïve (B) mice. (E) Bioluminescence imaging (BLI) performed before PET-CT imaging depicting strong signal from 4T1-Luc in lungs to confirm the presence of pulmonary metastasis. (F and G) Axial CT image of the thoracic cavity [green dashed outline in (A) and (B)] depicting the difference between 4T1 tumor-laden and tumor-free lungs (yellow outline). (H and I) Whole-tissue microscopy of H&E-stained sections from 4T1 metastatic and healthy lungs, harvested after terminal PET-CT scans.

Representative coronal PET-CT images (Fig. 3, A and B) depict the evolution of the whole-body distribution profile of [64Cu]NOTA-iNPG-pDox in tumor-bearing and naïve mice. A representative three-dimensional (3D)–rendered maximum intensity projection (MIP) is also illustrated (Fig. 3, C and D). [64Cu]NOTA-iNPG-pDox accumulated visibly and predominantly into 4T1 tumor-laden lungs [presence of tumors was confirmed by bioluminescence imaging (BLI) before PET-CT scans (Fig. 3E) as well as coregistered CT (Fig. 3F)] as early as 5 min p.i. (demarcated by green box in Fig. 3A). Axial PET-CT images showed a strong signal from the thoracic cavity at 3 hours, indicating reasonable retention of the radiolabeled particles in the diseased lungs (fig. S3A and movie S1). In contrast, radiolabeled pDox NPs injected directly in the absence of carrier iNPGs depicted >10-fold lower accumulation and faster clearance in metastatic lungs by 3 hours p.i. when compared to multistage iNPG-pDox DDS (fig. S3, C to E). Healthy tumor-free lungs (Fig. 3G) showed greatly attenuated accumulation of [64Cu]NOTA-iNPG-pDox (Fig. 3B) even at early time points, which diminished rapidly by 3 hours, possibly due to the washout of particles over time with blood flow (Fig. 3D and movie S2). [64Cu]NOTA-iNPG-pDox cleared rapidly from the blood pool, evident from the weak cardiac activity at 5 min p.i. and negligible signal at later time points (fig. S3A).

Tumor tropic transport of [64Cu]NOTA-iNPG-pDox was validated microscopically via immunohistochemistry (IHC). Lung tissues were harvested after the terminal PET-CT scans and subsequently processed for cross-sectional, high-magnification whole-tissue imaging (Fig. 3, H and I) and higher-resolution microscopy (fig. S4). iNPGs (golden brown discoidal particles, marked by yellow arrows) localized in the lung microcapillaries. Consistent with the in vivo PET-CT observations, significantly higher number of iNPGs could be observed in the tumor-bearing lungs (fig. S4A) compared to healthy control (fig. S4B). The results further validate that the systemic transport observed noninvasively with PET-CT reflected the actual trafficking pattern of the composite DDS and was not influenced by potential radionuclide dislocation.

Besides the lungs, organs of the MPS, liver and spleen, were the major sites of [64Cu]NOTA-iNPG-pDox sequestration in both groups, consistent with the observations for intravenously injected nano- and microparticles (movies S1 and S2) (5, 18). The in vivo observations were confirmed by ex vivo gamma-counting analysis and hematoxylin and eosin (H&E) staining (fig. S5, A and B). Absence of iNPG-pDox from the heart and no visible structural disruption of the cardiac tissues are important, in view of the widely reported instances of severe cardiomyopathy in patients administered with free or liposomal formulations of doxorubicin (19, 20).

Quantification of longitudinal trafficking of [64Cu]NOTA-iNPG-pDox

Sensitive, accurate, reliable, and robust quantification of the radionuclide signal adds a strong value to PET-CT–based noninvasive assessment of NP transport kinetics in vivo for preclinical and translational applications (10). Regions of interest (ROIs) were defined for each major anatomical site, lungs, heart (surrogate for blood pool activity), liver, spleen, kidney, and flank muscle (baseline control), and integrated across anterior to posterior planes of view. With 50.07 ± 2.92% injected dose per gram of tissue (%ID g−1) at 5 min p.i., 4T1 metastatic lungs were the largest initial recipient organ for intravenously injected [64Cu]NOTA-iNPG-pDox (Fig. 4, A and C, and fig. S5, C and E). In contrast, the accumulation of [64Cu]NOTA-iNPG-pDox in naïve lungs was significantly lower (26.2 ± 3.3 and 7.5 ± 2.5 %ID g−1, at 5 min and 3 hours p.i., respectively; Fig. 4, B and D, and fig. S5, C and F). A ~3-fold high accumulation of the composite drug in tumor-laden lungs (21.0 ± 5.3 %ID g−1 at 3 hours p.i.), compared to naïve lungs (fig. S5D), and the consequent enhancement in exposure [area under the curve (AUC)0-3h = 6066.6 ± 598.6 and 2282.0 ± 141.1 %ID g−1·hour in tumor-bearing and naïve lungs, respectively; Fig. 4E] suggested intrinsic tumor tropism of iNPG-pDox. Plasma time-activity curves (TACs) depicted a biphasic decline with an initial distribution phase and a terminal elimination phase. Thus, a two-compartment pharmacokinetic model best described the plasma concentration kinetics of iNPG-pDox particles in both healthy and diseased animals (fig. S5G). We observed that the particles have a short elimination half-life [t1/2, β≈ 75 min (healthy), 19 min (diseased)] and an even shorter distribution half-life (t1/2, α≈ 1 to 1.5 min) (table S1). The relatively shorter t1/2, β, larger clearance, Cl, and smaller plasma bioavailability, AUCP0 (table S1), observed in diseased animals compared to the healthy ones can be attributed to the increased pulmonary retention of the particles in the former (Fig. 4, C to E). Retention of [64Cu]NOTA-iNPG-pDox in lungs was computed as the ratio of lung-to-plasma activity, indicating a >5-fold enhancement in the tumor-bearing group when compared to naïve, by the 3-hour time point (fig. S5H).

Fig. 4 Quantitative analysis of systemic transport of [64Cu]NOTA-iNPG-pDox.

(A and B) Time-activity curves (TACs) representing temporal evolution of [64Cu]NOTA-iNPG-pDox disposition in major organs, liver, spleen, kidney, muscle, and lungs [4T1 tumor-bearing (A) or healthy control (B)]. (C and D) Lung accumulation of [64Cu]NOTA-iNPG-pDox in individual animals in tumor-bearing (C) or control (D) groups. (n = 3; solid line, mean; shaded region, SD) (E) Area under the curve (AUC) obtained from TACs (0 to 3 hours) in (A) and (B) for major organs in the two groups (n = 3, means ± SD; two-tailed t test; **P < 0.01, ***P < 0.001; ns, not significant). (F) Schematic of the multicompartment model developed to understand the in vivo disposition kinetics of [64Cu]NOTA-iNPG-pDox. Red arrows, plasma flow; green arrow, phagocytic uptake; yellow arrow, hepatobiliary excretion; and black arrows, association and dissociation of the particles with capillary endothelium (in lungs) and macrophages (in MPS). The MPS has an additional extravascular subcompartment (blue box). (G and H) Model fits obtained for major organs; Pearson correlation coefficient R > 0.99. Data represent means ± SD. (I) Schematic depicting the margination of iNPG in pulmonary microvessels resulting in enhanced transient interactions with the endothelium, characterized by phenomenological parameters kon,L and koff,L. (J) Table displaying average values of kon,L and koff,L derived from least squares fitting of the model to PET-CT measurements of [64Cu]NOTA-iNPG-pDox whole-body distribution. (K) High-magnification micrograph depicting the near-wall localization of a single [64Cu]NOTA-iNPG (red arrow) in a vWF-stained (brown) tumor-associated vessel in 4T1 metastasis-bearing lungs excised after 3 hours of PET-CT scan. Scale bar, 10 μm.

To understand these observations mechanistically, we developed a multicompartment, semimechanistic mathematical model to study the in vivo behavior of the microparticulate iNPG-pDox DDS (Fig. 4F). As shown in Fig. 4 (G and H) and fig. S6 (A and B), the model satisfactorily fits the mean experimental data (Pearson correlation coefficient R > 0.99) in both the naïve and diseased groups (fig. S10, C and D), providing reliable estimates of the unknown model parameters (table S2). First-order phenomenological rate constants (kon, L and koff, L) were introduced in the lung compartment (Fig. 4F) as lumped parameters that account for physicochemical (particle-related), hemorheological (blood composition–related), and hemodynamic (blood flow–related) factors responsible for governing the association and dissociation of the particles to and from the capillary wall (Fig. 4, I and J) to explain the differential retention of the particles in diseased and naïve lungs. As determined from the model, the ratio of kon, L to koff, L, which indicates the tendency of the particles to associate with the pulmonary microvasculature, is ≪1 in the naïve mice, while conversely, this ratio is >1 in mice with tumor metastases in the lungs (table S2). Inferences drawn from the mathematical model were verified by IHC on von Wildebrand Factor (vWF)–stained lung tissues. A representative high-magnification, high-resolution micrograph in Fig. 4K displays the margination of a plateloid iNPG microparticle (indicated by a red arrow) with the vascular endothelial cells of a tumor-associated microcapillary, possibly aided by local hematocrit changes or vascular remodeling.

TACs generated from the tomographic PET-CT images of liver, spleen, kidney, and muscle demonstrated similar trafficking trends for nontarget organs in the two groups (Fig. 4, A and B, and fig. S5, E and F). Within the first 5 min, 33.86 ± 12.36 and 36.26 ± 5.04 %ID g−1 is seen in the liver, while spleen accumulated 14.5 ± 4.7 and 14.8 ± 5.10 %ID g−1 in diseased and naïve animals, respectively, accounting for >70 %ID being sequestered by the MPS compartment. Further, the estimated value (~103 hour−1) of the kon, MPS parameter (table S2), which characterizes the attachment process between particles and macrophages in the MPS, indicates that the interaction occurs instantaneously, and the ratio of kon, MPS to koff, MPS parameter being >1 suggests that the attachment of particles to the macrophages is strongly favored over their detachment. Together, the rapid macrophage-assisted hepatic and splenic clearance of the particles from the plasma (21) contributes to the extremely short distribution half-life of [64Cu]NOTA-iNPG-pDox seen above. Further, the estimated phagocytic rate constant (kp ~1.26 hour−1) suggests a time scale of ~48 min for phagocytosis, which is in the same order of magnitude as the elimination half-lives of the particles, indicating that phagocytosis-driven hepatobiliary excretion is the primary excretion mechanism of iNPG-pDox particles.

To investigate the relative importance of model parameters on governing the delivery of iNPG-pDox to lungs, we performed a local sensitivity analysis with area under the pulmonary concentration kinetics curve (AUCL) as the model output of choice (fig. S6, E and F). Variation in the kon, L and koff, L parameters had relatively less influence on the kinetics of [64Cu]NOTA-iNPG-pDox in the healthy lungs, indicating that the plasma flow parameters and MPS-related parameters (kon, MPS, koff, MPS, and kp) could sufficiently explain the observed kinetics. In contrast, concentration kinetics in tumor metastasis in the lungs were found to be highly dependent on kon, L and koff, L, in addition to the other parameters. This analysis further corroborates our hypothesis that metastasis-induced local alterations in hemorheological and hemodynamic conditions significantly affect the kon, L and koff, L parameters, which may be manifested in the form of increased interactions between the particles and the vessel wall, resulting in larger exposure to the metastatic lungs.

Multiplexed OI depicts mesoscopic spatial heterogeneity in iNPG and pDox distribution

Spatial distribution of iNPG and pDox is expected to diverge temporally, as self-assembled pDox NPs (8) are released from iNPG microparticles arrested in the capillaries of tumor-laden lungs and infiltrate the metastatic tumors. Our in vitro studies indicated that >50% of loaded pDox was released from iNPGs within 24 hours of incubation in serum (fig. S2). Accordingly, we chose this time point to assess the differences in intra-organ trafficking patterns of the carrier and payload (Fig. 5). Dual-fluorescent AF647-iNPG-pDox were developed and intravenously injected in 4T1-GFP-Luc (green fluorescent protein) tumor-bearing and naïve mice (dose: 6 mg kg−1 of Dox). Animals were euthanized 24 hours p.i., and lungs were excised for nondestructive whole-organ multiplexed OI. Images were acquired over 18 filter sets encompassing the fluorescence spectra of GFP [exmax (maximum excitation): 488 nm, emmax (maximum emission): 510 nm], Dox (exmax: 460 nm, emmax: 550 nm), and AF647 (exmax: 640 nm, emmax: 680 nm). Figure 5 (A and B) shows spectrally unmixed single channel and merged multiplexed images of excised lungs from healthy and diseased groups, respectively. Obvious GFP signal was observed from tumor lesions in the metastatic lungs (mono green channel images in Fig. 5B). Image analysis indicated that 45 ± 15% of the lung surface area was covered by 4T1 nodules. AF647-iNPG depicted significantly higher retention in metastatic lungs compared to the naïve group at 24 hours p.i., indicated by higher average radiant efficiency (P < 0.05; Fig. 5C), dose-normalized signal intensity (P < 0.05; Fig. 5D), and surface area coverage (P < 0.01; Fig. 5E). Consequently, the delivery, and hence the fluorescence intensity of pDox in the tumor-bearing lungs, was significantly enhanced (P < 0.05; red bars in Fig. 5C). Dose-normalized pDox fluorescence from the control lungs was found to be >2.5-fold lower than the metastatic lungs (P < 0.01; Fig. 5D), while surface area coverage was >2-fold less (P < 0.01; Fig. 5E). Released pDox NPs demonstrated greater spatial overlap with GFP-expressing 4T1 tumors than AF647-iNPG (Fig. 5G and fig. S7).

Fig. 5 Multiplexed whole-organ imaging reveals intra-organ spatial trafficking of AF647-iNPG-pDox.

(A and B) Spectrally unmixed single channel and merged images of intact lungs, excised from (A) healthy and (B) 4T1 tumor-bearing mice, 24 hours p.i. of AF647-iNPG-pDox, depicting fluorescence from 4T1-GFP (green), AF647-iNPG (blue), and pDox NP (red). Scale bar, 5 mm. (C to E) Semiquantitative ROI analysis of fluorescence signals from AF647-iNPG and pDox in diseased and healthy lungs; (C) average radiant efficiency [(p s−1 cm−2 sr−1)/(μW cm−2)], (D) dose-normalized signal intensity, and (E) area fraction (%). (F and G) Pixel-by-pixel line profile analysis of spectrally unmixed signals from healthy (F) and tumor-bearing (G) lungs indicating GFP (green), iNPG (blue), and pDox (red) intensities along the white dashed arrows in (A) and (B), respectively. (H) Line profiles across individual tumor lesions (white circles) depicting diminishing signals from colocalized iNPG and pDox with increasing size of lesion (I to V). Lesion dimensions presented as means ± SD (n = 3 measurements). (I) Heat map profile quantifying the GFP (4T1 tumor), AF647-iNPG, and pDox signals across tumor lesions of increasing size (I to V).

A tumor size–dependent iNPG-pDox transport heterogeneity was discerned from the multiplexed fluorescence images (Fig. 5B, merged). Thus, we sought to perform a semiquantitative, pixel-by-pixel line profile analysis to measure the spatial overlap between the GFP, AF647-iNPG, and pDox signals across individual metastatic lesions of varying sizes (Fig. 5B, white dashed circles). Small tumors showed substantial colocalization of GFP, iNPG, and pDox fluorescence (Fig. 5H, I). However, increase in tumor size (green curves, Fig. 5H, I to V) was accompanied by a rapid decrease in overlapping AF647-iNPG intensity (blue curves, Fig. 5H, I to V). Similarly, overlap between pDox (red line) and GFP was more pronounced in small- to medium-sized metastatic lesions, visibly higher than that of AF647-iNPG and GFP (Fig. 5H, I to III). A gradual reduction in pDox signal was observed in larger tumors >1 mm in diameter (Fig. 5, B and H, IV and V), highlighted through heat map profiling of average fluorescent intensities across the diametric line profiles (Fig. 5I). Notably, AF647-iNPG signal peaked in the peritumoral regions, flanking the GFP signal peak, while the peak of the pDox signal coincided with that of GFP, suggesting that the two components of our multistage delivery vector occupied distinct regions within the metastatic organ.

High-resolution microscopy delineates intratissue spatial heterogeneity in iNPG-pDox distribution

We next dissected the trafficking patterns of AF647-iNPG vectors and pDox NPs at tissue and cellular scale with respect to tumor lesion size and associated microvasculature (indicated by endothelial biomarker CD31, green). High-resolution confocal microscopy on cryo-sectioned lung tissues confirmed the mesoscopic observations of the preceding experiment, whereby greater accumulation of AF647-iNPG, and consequently of pDox NPs, could be observed in the vicinity of tumor cell clusters and micrometastasis (Fig. 6, A and B, and fig. S8A) than in macrometastases (Fig. 6C and fig. S8B). Higher-magnification imaging showed dispersed pDox NPs (red) released locally from clusters of AF647-iNPG microparticles (pink) arrested in the tumor-associated blood vessels (zoomed-in image in Fig. 6B and fig. S8A) over a distance of few micrometers, as depicted in the normalized line profile analysis (Fig. 6D, II). Although the spatial shift in pDox signal is clearly visible, a partial overlap between pDox and iNPG signals can also be discerned (Fig. 6D, I), indicating that a fraction of pDox was released by 24 hours (corroborating our in vitro drug release and whole-organ OI observations).

Fig. 6 High-resolution microscopy reveals tumor size–dependent transport heterogeneity of AF647-iNPG-pDox.

Confocal images depict greater accumulation of AF647-iNPG (pink) and released pDox NPs (red) in (A) 4T1 tumor cell clusters and (B) micrometastases. The tissues are stained for CD31 endothelial marker (green) to demarcate the blood vessels and counterstained with DAPI (blue). Scale bar, 50 μm. Zoomed-in image of area boxed in yellow is shown on the right depicting the release of pDox NPs from iNPG adhering to a nearby capillary. Scale bar, 5 μm. (C) A single macrometastatic nodule depicting locally remodeled vascular architecture inducing peritumoral accumulation of iNPG (pink), while released pDox NPs (red) diffuse into the core. Scale bars, 100 and 50 μm, as indicated. (D) Line profiles obtained along the white dashed arrows shown in (A) to (C) indicate spatial distribution of CD31 (green), iNPG (pink), and pDox (red). (E) Numerical solution of pDox transport model (diffusion equation) showing the spatial distribution of pDox NPs along the tumor diameter at different time points p.i. Note that 0 on the x axis denotes the vascular end of the tumor periphery. The mass of pDox in the tumor was fixed at 0 %ID at 0 hour. (F) Distribution of AF647-iNPG-pDox in healthy lungs. Zoomed-in image of area boxed in yellow depicts release of pDox from a solitary iNPG near pulmonary blood vessels. Scale bar, 5 μm. (G to J) Image analysis indicates significant difference in area coverage (G and H) and fluorescence intensity (I and J) of iNPG and pDox between diseased and healthy tissues. Bars represent means ± SEM; two-tailed t test; **P < 0.01, ***P < 0.001.

Macrometastases (≥500 μm in diameter) demonstrated starkly different distribution of the two components. Metastatic growth was accompanied by a gradual but significant remodeling of pulmonary vasculature (Fig. 6C and fig. S8B), when compared to healthy, tumor-free lungs (Fig. 6F), although no significant difference could be observed in area fraction or signal intensity of CD31 stain (fig. S8C). Peritumoral areas were markedly distinguished by rearranged air space pattern and deviation from the mesh-like architecture of vessels. Notably, though co-option of alveolar capillaries and larger vessels occurred at higher frequency and could be clearly demarcated throughout tumor nodules of varying sizes (Fig. 6, A to C, and fig. S8, A and B), metastatic growth beyond 1 to 2 mm in diameter additionally presented dilated, tortuous, and tubular vessels both in the core and in periphery (fig. S8B), possibly due to vessel remodeling and angiogenesis (22). AF647-iNPG was mainly found confined to the rearranged peritumoral capillaries, as the cores remained conspicuously devoid of larger microparticles (Fig. 6C). pDox NPs released from AF647-iNPG were found both in the periphery and diffused into the tumor cores, displaying a gradual gradient in intensity (indicated by high-magnification line profile analysis from periphery to core; Fig. 6D, III). Notably, tumors ~1 mm in diameter depicted almost complete absence of iNPGs, while accumulation of pDox NPs was also greatly diminished (fig. S8B), corroborating our observations from whole-organ OI (Fig. 5H, IV and V).

In accordance with systemic PET-CT and multiplexed OI, accumulation of AF647-iNPG-pDox in naïve lungs was significantly lower (Fig. 6F), in terms of both area fraction coverage (Fig. 6, G and H) and fluorescence signal intensity of the two components (Fig. 6, I and J). Notably, released pDox NPs occupied a noticeably larger area of metastatic lungs than iNPGs (~2-fold increase in area fraction; Fig. 6, G and H).

To understand the intratumoral transport of pDox released from the vasculature-bound iNPG-pDox, we solved the pDox transport model numerically from 0 to 24 hours p.i. (fig. S9A). As shown in Fig. 6E, the model simulates the spatiotemporal distribution of released pDox along the tumor diameter. Simulation results show that by the 24-hour time point, the maxima of the pDox diffusion curve lies around 0.35 mm, corroborating with our observations (Figs. 5H and 6, A and B). In addition, distance penetrated by pDox into the tumor interstitium increased with the passage of time; however, this was accompanied by a reduction in the maxima of the pDox curve. This indicates that as time progressed, pDox diffusion into the tumor interstitium moved the maxima of the curve forward, but the magnitude of the maxima declined, arguably due to reduction in the influx of pDox at the boundary. This can be attributed to the time-dependent decrease in the bound iNPG-pDox at the vascular end, as evident from our findings of PET-CT imaging (Fig. 3) and multicompartment pharmacokinetic model (Fig. 4G). The above observations are qualitatively consistent with the line profiles of pDox distribution obtained through experimental OI measurements (Figs. 5H and 6D).

We further used the model to investigate strategies to improve distribution of pDox in the metastases. For this, we analyzed the effect of parameters that govern the availability and transport of pDox on the spatiotemporal distribution of pDox along the tumor diameter. As shown in fig. S9B, the smaller size of released pDox NPs appears to improve the distribution across the tumor, which is an effect of improved diffusivity of pDox through the interstitium. Further, the release rate of pDox from iNPG-pDox seems to affect the maxima of the pDox curve, such that the larger the release rate constant, the greater the magnitude of the maxima (fig. S9C). Similarly, a larger value of the ratio of kon,L to koff,L of iNPG-pDox interaction with lung vasculature increases the maxima of the pDox curve (fig. S9D). This indicates that increasing the binding tendency of iNPG-pDox to the lung vasculature (possibly through active targeting) can enhance the availability of pDox for improved tumor delivery. Thus, on the basis of the above findings, by optimizing the characteristics of iNPG-pDox and/or pDox payload, we can potentially improve the intratumoral transport and delivery of pDox to pulmonary MBC.


Drug delivery approaches and systems have evolved tremendously in the past decades, owing to advances in physical sciences and engineering as applied to oncology (23). Many recent breakthroughs in experimental therapeutics have a multicomponent framework, formulated to act synergistically toward improving the delivery and efficacy of the active drug ingredient. In contrast, traditional chemotherapy that focuses on nonspecifically modulating the dosage of small molecule drugs is often met with only marginal improvements in treatment response and significant side effects (6). However, the empirical nature of drug development and preclinical evaluation paradigms has greatly impeded the rational design and translation of effective multicomponent therapeutics. Multiplexed, multiscale imaging strategies, as presented in this study, can provide a means to precisely monitor the transport of complex delivery systems within the body and evaluate, screen, and compare different formulations at preclinical, translational, and clinical levels. Further, this framework can be extended to reveal spatiotemporal determinants of efficacy of innovative interventions such as therapeutic vaccines (24) and immunomodulatory drugs (25) that require different sites of action for different components.

In this retrospective study, we harnessed two naturally orthogonal complements, PET-CT and OI, to thoroughly characterize the multiscale negotiation of physiochemical and biological barriers by our composite drug, iNPG-pDox, in an effort to unravel the reason behind its efficacious treatment of pulmonary MBC, reported earlier (8). Despite significant progress, metastatic cancer that accounts for >90% of cancer-related deaths (3) presents a formidable and as-yet incurable challenge for even the most advanced therapeutic technologies. Therefore, visualization and reverse engineering of the spatiotemporal kinetics of nanomedicines that effectively target metastases can allow directed development and rational improvement of other advanced therapies in the future.

We used PET-CT imaging to visualize the early events in the whole-body trafficking of iNPG-pDox. Because of its widespread application in diagnosis, all major hospitals and medical centers are equipped with PET-CT scanners, streamlining the logistical operations during translational and clinical studies. Excellent spatiotemporal resolution afforded by clinical PET-CT scanners, whole-body imaging capability in humans and nonhuman primates, and flexible and modular nature of PET-CT make it highly suitable for translational and clinical evaluation of experimental drugs and DDS. We used 64Cu for its short half-life t1/2 ≈ 12.7 hours that is more suitable for short circulation iNPG-pDox DDS, while longer-circulating pDox NPs were radiolabeled with zirconium-89, t1/2 ≈ 78.4 hours, using the same facile radiochemistry, which attests the versatile and inherently modular nature of our presented framework (26). Further, we show that experimental techniques like ex vivo OI and microscopy that offer exquisite spatial resolution and multiplexed tracking of different fluorescent labels, but have been traditionally limited by their qualitative nature, can serve as powerful accomplices to quantitative PET-CT and mathematical modeling to comprehensively interrogate drug kinetics and predict viable strategies for improvement of drug delivery and efficacy. Alongside the macroscopic information yielded by systemic imaging modalities, knowledge of the local microdistribution patterns of nanomaterials within an organ or tissue or at subcellular levels can serve as a direct determinant and predictor of therapeutic efficiency.

While noncompartment and compartmental PK analyses are routinely used to characterize the in vivo disposition of nano- and microparticles, a multicompartment, whole-body scale model is more appropriate to obtain a comprehensive understanding of their biodistribution and clearance, including kinetics at the target site. The model developed in this study is a reduction of a traditional physiologically based pharmacokinetic model (27) for small molecules, as guided by the in vivo data and an understanding of the physical limits to transport processes imposed by particle size. Despite its simplistic nature, the model is fairly universal for similarly sized particles, provides useful insights into the in vivo kinetics of iNPG-pDox, and can potentially be used for interspecies scaling to make prospective predictions in humans or other species.

Given the extensive surface area of the pulmonary vasculature (28) and the high margination propensity of plateloid silicon microparticles (2931), it is expected that iNPG-pDox will have a tendency to interact with the endothelial cells of capillaries in the lungs. This can lead to the formation of reversible nonspecific bonds at the particle-cell interface, thereby slowing down or arresting the particles in the lung microvasculature. As mentioned previously, the kon,L and koff,L parameters were incorporated to encompass the physicochemical, hemorheological, and hemodynamic factors affecting these interactions. Since the physicochemical particle characteristics remain unchanged, our results suggest that the presence of metastases may alter the local hemorheology and hemodynamic conditions in diseased mice, thus increasing iNPG-pDox margination (30, 32, 33) and adhesion to the microvasculature in diseased lungs. Further, local changes in blood viscosity, hematocrit, and vessel architecture can increase the flow resistance of blood in tumors (34), thereby lowering the shear rate that governs the dislodging of particles from the vessel walls and increasing the tendency of the particles to remain attached to the capillary wall in metastatic sites. High-resolution, high-magnification microscopic analysis validated this hypothesis, where distinct remodeling of the pulmonary vascular and bronchial network (fig. S8), as well as deregulated deposition of extracellular stromal components such as collagen and smooth muscle actin (fig. S10), could be observed, possibly induced by the evolving and invasive metastatic growth.

One concern about engineering of complex systems to make them amenable to multiplexed imaging can be the apparent influence on their therapeutic efficacy. However, as demonstrated, our imaging technology necessitated minimal surface engineering of the composite drug, producing negligible change in their physicochemical properties. The chemistry is inherently modular and allows ready adaptation to any imaging modality of choice, as well as test system including NPs, polymers, biologics, or even cell-based delivery vectors. Furthermore, the proposed technology requires no additional sophisticated resources for development and characterization of engineered therapeutics, making it widely applicable to most basic and translational research settings. Combined with research efforts across multiple disciplines, our presented framework can serve as a powerful tool to design rational DDS and evaluate the nature and efficacy of downstream therapeutic responses. Overall, we believe that the present multiplexed and multiscale monitoring approach can be an excellent complement to current translational workflows for rapid and accurate evaluation of advanced complex therapies.


Surface modification of iNPG-pDox

Discoidal porous silicon particles (iNPG) with dimensions 2.6 μm by 700 nm were fabricated as described previously (16). Particles were then partially oxidized and modified by 2% (v/v) 3-aminopropyltriethoxysilane (APTES) in isopropanol for 48 hours at 55°C to confer primary amines on the particle surface. Polymeric doxorubicin (pDox) was separately synthesized by conjugating the drug Dox to side chains of poly(l-glutamic acid) via a pH-sensitive hydrazone linker following previous protocols (8).

For radiolabeling, iNPG was first conjugated to NOTA (Macrocyclics, Dallas, TX) (17). Briefly, 3 billion APTES-modified (-NH2-bearing) iNPG particles were dissolved in 0.5 ml of N,N′-dimethylformamide, followed by addition of 0.8 mg of p-SCN-Bn-NOTA in 0.5 ml of dimethyl sulfoxide (DMSO). The mixture was sonicated until a homogeneous suspension was observed, followed by the addition of 20 μl of diisopropylethylamine. The reaction was allowed to proceed at 25°C for 4 hours under vigorous shaking at 700 rpm. The final conjugate was washed by centrifugation successively, twice with DMSO, followed by isopropanol and water. As-synthesized NOTA-iNPG was suspended in water and freeze-dried for further usage. For radioactive labeling, NOTA-iNPG was suspended in 0.1 M Na2CO3 solution, and the pH was adjusted to 5.5. Copper-64 [64CuCl2; 37 megabecquerel (MBq); Washington University, St. Louis, MO] was added, and the mixture was incubated at 37°C for 2 hours. Excess non-chelated 64Cu was removed by centrifugation at 4000 rpm. To obtain the final constructs, polymeric pDox monomers (3 mg) were loaded into [64Cu]NOTA-iNPG constructs, in DMSO solution (70 mg ml−1), following vigorous vortexing for 4 hours. Unloaded pDox was removed by centrifugation. Radio–thin-layer chromatography (radio-TLC; Bioscan AR-2000 Radio-TLC Imaging Scanner, Eckert & Ziegler, Valencia, CA) using 50 mM EDTA as the mobile phase was used to calculate radiochemical yields and assess time-dependent radiostability in 1× mouse serum.

A similar procedure was adopted to develop fluorescent iNPG by replacing chelator NOTA with NIR dye, AlexaFluor AF647 (excitation: 640 nm; emission: 680 nm, Sigma, St. Louis, MO) to form AF647-iNPG. AF647-iNPG microparticles were filled with pDox monomers as described above to produce AF647-iNPG-pDox for multiplexed OI.

To perform 89Zr radiolabeling of pDox NPs, pDox was first conjugated to p-SCN-Bn-Deferoxamine (DFO; Macrocyclics, Dallas, TX) via amine-isothiocyanate conjugation chemistry as described above (26). Final conjugates were obtained after filtration through PD-10 desalting columns (GE Healthcare, USA). As-synthesized pDox-DFO NPs were suspended in 0.1 M HEPES buffer, and the pH was adjusted to 7 to 7.5. Zirconium-89 (89Zr oxalate; 37 MBq; Washington University, St. Louis, MO) was added, and the mixture was incubated at 37°C for 2 hours. Excess non-chelated 89Zr was removed by PD-10 columns to obtain the final radiolabeled product [89Zr]DFO-pDox.


Morphology of iNPG, iNPG-pDox, NOTA-iNPG, [64Cu]NOTA-iNPG-pDox, and AF647-iNPG-pDox was visualized on a scanning electron microscope (Nova NanoSEM 230, Thermo Fisher Scientific, USA) after each step of chemical conjugation and drug loading. To prepare samples for scanning electron microscopy, 20 μl of sample (105 ml−1 in Milli-Q water) was applied on aluminum stubs (Ted Pella Inc., USA) and allowed to dry at room temperature, followed by thermal pretreatment at 90°C to remove any traces of moisture. Conditions for scanning electron microscopy operation were kept constant across all samples: high vacuum, 3 × 106 Torr; acceleration voltage, 5 kV; spot size, 3; working distance, 5 mm. Surface charge was measured with Zetasizer Nano ZS (Malvern Instruments Ltd., Worcestershire, UK) using Milli-Q water as the dispersant, following software instructions. Absorbance and emission spectra of fluorescent probes were obtained in aqueous solution with an ultraviolet-visible-NIR spectrophotometer (Synergy H4 Hybrid, BioTek, USA). All measurements were performed in triplicate.

In vitro degradation study of [64Cu]NOTA-iNPG-pDox

To test whether NOTA conjugation or radiolabeling will influence degradability of iNPG, [64Cu]NOTA-iNPG-pDox were incubated in 1× mouse serum at 37°C for 1 day. Aliquots were drawn at different time points and centrifuged at 4000 rpm and washed with phosphate-buffered saline (PBS) once. Supernatant was discarded, and dried pellets were stored as is, until complete radioactive decay. “Cold” samples were then visualized on a scanning electron microscope to assess the degree of degradation.

In vitro pDox release

One billion AF647-iNPG-pDox vectors were incubated in 1× mouse serum and gently shaken at 37°C over a period of 5 days. Aliquots were collected at different time points and centrifuged at 4000 rpm. Supernatant (released pDox NPs) was collected and measured by UV absorbance.

Animal models

All animal studies were performed in accordance with the guidelines from the Institutional Animal Care and Use Committee at the Houston Methodist Research Institute. Murine breast cancer cell line 4T1, transfected with plasmid carrying both luciferase and GFP genes (4T1-GFP-Luc), was cultured in RPMI 1640 medium, supplemented with 10% fetal bovine serum and 1% penicillin and streptomycin antibiotics at 37°C and 5% CO2. An experimental pulmonary metastasis model was generated by injecting 1 × 105 4T1-GFP-Luc cells into the tail vein of syngeneic 5- to 6-week-old female Balb/c mice (Envigo, Indianapolis, IN). Tumor growth was regularly monitored using BLI (IVIS Spectrum, PerkinElmer). Tumor-bearing animals were used roughly 10 days later for all imaging studies. Age-matched, tumor-free naïve mice were used as control.

Whole-body static and dynamic PET-CT imaging and image analysis

Noninvasive longitudinal PET-CT imaging was performed on a micro PET-CT Inveon rodent model scanner (Siemens Medical Solutions USA, Inc.). For dynamic PET-CT imaging, 60 to 80 μCi [64Cu]NOTA-iNPG-pDox was injected per mouse using a homemade catheter system installed in the tail vein, concurrently with the initiation of the scan. A 30-min dynamic sequence was acquired and framed into six frames: 6 × 300 s. Following this, static scans were performed at 3 hours p.i. After the terminal scans, all animals were euthanized, and major organs were harvested and wet-weighed, and tissue radioactivity was measured on a gamma counter (PerkinElmer). Immediately afterward, lungs and selected tissues (liver, spleen, heart, kidney, and flank muscle) were fixed in 4% paraformaldehyde overnight for immunohistochemical analysis. For longitudinal static PET-CT imaging with [89Zr]DFO-pDox, 60 to 80 μCi [89Zr]DFO-pDox (dose: 6 mg/kg of Dox) was injected via the tail vein, and static PET-CT scans were acquired at 0.5, 3, 24, and 48 hours p.i.

Image reconstruction and ROI analysis were performed on decay-corrected, DICOM-formatted images using vendor software (Inveon Research Workplace). ROIs were kept constant during the analysis of all frames on a single dynamic scan. Volume of ROIs did not vary more than 5% between different scans. All data are presented as percent injected dose per gram of tissue (%ID g−1).

Pharmacokinetic analysis

A two-compartment PK model was used to estimate the systemic half-lives (distribution, t1/2, α; elimination, t1/2, β) of iNPG-pDox particles by least squares fitting of the closed-form solution of the model to the plasma concentration kinetics. The following biexponential equation is the analytical solution of a two-compartment PK modelCP=C1ek1t+C2ek2t(1)where CP is the plasma concentration of the particles; C1 and C2 are the intercepts of the back-extrapolated initial distribution phase and terminal excretion phase, respectively, with slopes k1 and k2 on a semilogarithmic scale used to calculate distribution (t1/2, α = ln (2)/k1) and elimination (t1/2, β = ln (2)/k2) half-lives, respectively. Clearance (Cl), the volume of plasma cleared of particles per unit time, was estimated asCl=DAUCP0(2)where D represents the injected dose of particles (100% ID) and AUCP0 is the area under the curve of the plasma concentration kinetics curve from 0 to ∞ time. For the calculation of Cl, plasma density was assumed to be 1 g ml−1.

Mathematical model

We developed a semimechanistic mathematical model, which is a reduced physiologically based pharmacokinetic model, based on conservation of mass and law of mass action in a physiologically relevant network of compartments representing the dominant distribution sites of iNPG-pDox in vivo. The model was formulated as a system of ordinary differential equations and contains both mechanistic and phenomenological parameters that describe the key processes responsible for iNPG-pDox biodistribution and clearance, which were either obtained directly from published literature or estimated via nonlinear regression (see table S3). Briefly, the model consists of five compartments [lungs, plasma, MPS (comprising liver and spleen), muscle, and kidney], which represent the total volumes (Vi) of the given organ or tissue and are connected in an anatomical fashion via plasma (Qi) flow rates. MPS in the model is subcompartmentalized into vascular (VMPS,v) and extravascular (VMPS,e) volumes, consistent with the MPS being recognized as a sink-like tissue (35), where particles are typically phagocytized and translocated from the vascular space into the extravascular tissue. Phagocytosis was empirically modeled as a two-step process: (i) reversible attachment of iNPG-pDox with the macrophages, characterized by the first-order rate constants kon,MPS and koff,MPS for association and dissociation, respectively, and (ii) engulfment and translocation of the particles into the extravascular space characterized by the phenomenological first-order rate constant kp. These two steps can vary significantly in their time scales, with the former being much faster than the latter (36). From the extravascular subcompartment of the MPS, the particles are excreted via hepatobiliary excretion governed by the bile flow rate B.

Interaction of iNPG-pDox with pulmonary vasculature in the lung compartment was modeled by introducing first-order phenomenological rate constants, kon, L and koff, L, to describe the association and dissociation of the particles to and from the capillary wall, respectively.

Model equations

The model, based on conservation of mass and law of mass action, consists of the following system of equations and initial conditions.

Plasma. VPdCPdt=CP(QL+QMPS+QM+QK)+CLfreeQL+CMPS,vfreeQMPS+CMQM+CKQKCP(0)=C0(M1)where CP is the concentration of particles in the plasma compartment, which has a volume VP; Qi represents plasma flow rate of compartment i (where i = L, MPS, M, and K for lungs, MPS, muscle, and kidneys, respectively); CLfree is the free particle concentration in the lungs, i.e., particles not associated with the endothelial cells; CMPS,vfree is the free particle concentration in the MPS vasculature, i.e., not associated with the macrophages; CM and CK represent particle concentration in muscle and kidney, respectively; and C0 is the concentration of particles in plasma at time t = 0. Note that QMPS is composed of the blood flow from hepatic artery and portal vein (spleen, pancreas, small intestine, and large intestine).

Free-flowing particles (not associated with endothelial cells)

VLdCLfreedt=CPQLCLfreeQLCLfreeVLkon,L+CLboundVLkoff,LCLfree(0)=0(M2)where VL is the volume of the lungs; kon, L and koff, L are the first-order rate constants that characterize the nonspecific association and dissociation of the particles with and from the lung microvasculature, respectively; and CLbound represents the concentration of the particles bound to the lung vasculature.

Bound particles (associated with endothelial cells)


Mononuclear phagocytic system.
Vascular free (free-flowing particles in vascular volume)

VMPS,vdCMPS,vfreedt=CPQMPSCMPS,vfreeQMPSCMPS,vfreeVMPS,vkon,MPS+CMPS,vboundVMPS,vkoff,MPSCMPS,vfree(0)=0(M4)where VMPS,v is the volume of the vascular space in the MPS; kon,MPS and koff,MPS are the first-order rate constants that characterize the association and dissociation of the particles with and from the macrophages, respectively; and CMPS,vbound is the macrophage-bound particle concentration in the MPS vasculature.

Vascular bound (particles bound to macrophages in vascular volume)

VMPS,vdCMPS,vbounddt=CMPS,vfreeVMPS,vkon,MPSCMPS,vboundVMPS,vkoff,MPSCMPS,vboundVMPS,vkpCMPS,vbound(0)=0(M5)where kp is the first-order rate constant that characterizes phagocytosis and translocation of particles into the extravascular subcompartment.

Extravascular (particles in extravascular volume)

VMPS,edCMPS,edt=CMPS,vboundVMPS,vkpCMPS,eB, CMPS,e(0)=0(M6)where CMPS,e is the particle concentration in the extravascular subcompartments of the MPS; VMPS,e is the volume of the extravascular space in the MPS; and B is the bile flow rate responsible for hepatobiliary excretion of the particles.

Muscle. VMdCMdt=CPQMCMQM, CM(0)=0(M7)where VM is the volume of the muscle.

Kidneys. VKdCKdt=CPQKCKQK, CM(0)=0(M8)where VK is the volume of the kidneys. Note that renal excretion is negligible for particles >10 nm (37); therefore, the kidney compartment in our model is not an excretory compartment, and hepatobiliary excretion forms the sole basis of iNPG-pDox excretion.

Sensitivity analysis

Local sensitivity analysis was performed to evaluate the dependence of iNPG-pDox exposure to the lungs (i.e., AUC for lung concentration kinetics, AUCL) on various model parameters, excluding the ones that were physiological constants derived from literature. Analysis was performed by determining the sensitivity coefficient (SC) for each tested parameterSC=AUCpertAUCoptimAUCoptimPpertPoptimPoptim(3)where Poptim is the optimized value of the parameter obtained through model fitting, Ppert is the perturbed value of the parameter (between −50 and +50% of Poptim), AUCoptim is the value of AUCL obtained from the optimized parameters, and AUCpert is the value of AUCL obtained using perturbed parameters (parameters were perturbed one at a time). Trapezoidal numerical integration (via the built-in MATLAB function trapz) was used to estimate the AUCs.

Ex vivo whole-organ multiplexed OI

A separate cohort of tumor-bearing and tumor-free mice was used to assess the macroscopic intra-organ distribution of AF647-iNPG-pDox vectors. AF647-iNPG-pDox (dose: 6 mg kg−1 of Dox) was intravenously injected, and the animals were euthanized 24 hours p.i. Lungs were harvested and immediately imaged ex vivo on IVIS Spectrum (PerkinElmer). Images were acquired for GFP, Dox, and AF647 fluorescence with 18 filter sets, predefined in vendor software (Living Image, PerkinElmer) using the same acquisition conditions for all animals. Spectral unmixing and data analysis were performed in automatic mode on Living Image software, using four-component principal components analysis. Contribution from autofluorescence was subtracted during the analysis. Volume- and area-matched ROIs were drawn on spectrally unmixed images for each channel and presented as average radiant efficiency [units: (p s−1 cm−2 sr−1)/(μW cm−2)]. Image analysis, including pixel-by-pixel line profile analysis and area fraction quantification, was performed on FIJI [National Institutes of Health (NIH)] software where appropriate. Automated thresholding was applied after segmenting the images into individual color channels.

IHC and immunofluorescence microscopy

Tissues harvested after terminal PET-CT scans at 3 hours p.i. of [64Cu]NOTA-iNPG-pDox were fixed in 4% paraformaldehyde overnight, followed by washing in ethanol and PBS, paraffin embedding, and sectioning into 6-μm-thick sections. H&E, vWF (vessel), Col-I (collagen), and α-SMA (α–smooth muscle actin) staining was performed on separate sections to analyze metastasis morphology and particle distribution. Wide-field whole-tissue images were obtained on EVOS Auto FL System (Thermo Fisher Scientific, USA). High-magnification, high-resolution images were obtained on a Nikon Eclipse Ti microscope with constant acquisition settings.

Tissues harvested for ex vivo OI, 24 hours p.i. of AF647-iNPG-pDox, were rinsed in PBS and preserved in TissueTek OCT compound at −80°C. Frozen tissue blocks were cryosectioned into 6-μm-thick sections for immunofluorescence microscopy and image analysis. Briefly, tissue sections were incubated in cold acetone for 10 min, followed by blocking with 1:1 horse-goat serum, supplemented with 0.1% Triton X-100, for 1 hour. Primary antibodies (see table S4) were added to the tissue for overnight incubation at 4°C. Tissues were washed in PBS (three times, 5 min each time), and secondary antibodies were added for 1 hour. After washing three times in PBS, 4′,6-diamidino-2-phenylindole (DAPI; Life Technologies) was added for 10 min, followed by further washing and mounting with Gold anti-fade mounting medium. Stained sections were observed under Nikon Ti confocal microscope. Acquisition settings were kept constant across groups. All analyses were performed using Nikon Elements software and FIJI (NIH) image analysis software. Gaussian blur correction was applied across all images. For quantification, segmented images were automatically thresholded, and fluorescence intensities and fraction area coverage were quantified.

Intratumoral transport modeling of pDox NPs

A diffusion-based transport model of pDox transport in the tumors was developed with the premise that iNPG-pDox particles that are bound to the pulmonary vasculature release pDox into the tumor interstitium, which then diffuses from the tumor periphery toward the tumor core (fig. S9A). Advective intratumoral transport of pDox was ignored in this problem under the assumption that high interstitial fluid pressure can hamper convection of interstitial fluid. As a simplification, this process is modeled as the 1D diffusion of pDox along the tumor diameter, starting from the tumor periphery. For this, we solved a 1D diffusion equation (Eq. 4) with mixed boundary conditions. The boundary condition on the vascular end of the tumor domain is defined by the temporal kinetics of bound iNPG-pDox (obtained from the multicompartment pharmacokinetic model) and the release kinetics of pDox (obtained from the data given in fig. S2D). According to the data in fig. S2D, first-order kinetics governs the fraction of pDox, contained in the bound iNPG-pDox, released into the interstitium. No-flux boundary condition is applied on the opposite end of the tumor domain; i.e., pDox is assumed to not exit the tumor domainNpDOXt=D2NpDOXx2, for t in (0,T)and x in (0,L),NpDOX(0,t)=(1ekrelt)VLCLbound(t)(Dirichlet boundary condition)NpDOXx (L,t)=0 (Neumann boundary condition)NpDOX(x,0) (Initial condition)(4)where NpDOX represents the mass of pDox (units, %ID), D is the diffusion coefficient of pDox, krel is the first-order release rate constant of pDox from iNPG-pDox, VL is the volume of lungs, and CLbound(t) is the concentration of iNPG-pDox bound to the pulmonary vasculature.

Parameter values are given in table S3, and the model was solved numerically in MATLAB through the finite difference method using an explicit scheme (forward Euler).

Statistical analysis

Statistical analyses were performed using GraphPad Prism 8.0 (GraphPad Software Inc., CA, USA), MATLAB R2018a, and Microsoft Excel. Animals were randomized before imaging studies. Two-tailed Student’s t test (α = 0.05) was used to evaluate statistical differences between the two groups. Multiple comparison tests were performed with Holm-Sidak adjustment using GraphPad Prism, where required. Data are presented as means ± SD or SEM as noted. P values less than 0.05 were considered statistically significant; *P < 0.05, **P < 0.01, ***P < 0.001, ****P < 0.0001. The model system of equations was solved through numerical integration using the nonstiff ODE (ordinary differential equation) solver ode45 in MATLAB. Least squares fitting of the models (semimechanistic and two-compartment) to experimental data was performed using the “Levenberg-Marquardt” algorithm in MATLAB.


Supplementary material for this article is available at

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.


Acknowledgments: We would like to thank J. Gu of the Electron Microscopy Core at the Houston Methodist Research Institute for help with scanning electron microscopy imaging. All animal experiments were performed according to the IACUC guidelines of the Houston Methodist Research Institute. Funding: This work was partially supported by NIH grants U54CA210181, R01CA193880, R01CA222959, U01CA196403, U01CA213759, R01CA226537, and R01CA222007; U.S. Department of Defense grant W81XWH-17-1-0389; National Science Foundation grants DMS-1716737 and DMS-1930583; and the Ernest Cockrell Jr. Presidential Distinguished Chair. Author contributions: S.G., M.F., and H.S. conceived and designed the experiments. S.G. performed surface modification of iNPG-pDox, material characterization, 64Cu and 89Zr labeling, radiochemistry analysis, PET-CT imaging, necropsies, gamma counting, PET-CT data analysis, tissue extraction and processing, ex vivo IVIS imaging and analysis, tissue staining, confocal imaging, image processing and analysis, and manuscript writing. G.Z. and Z.H. synthesized and characterized pDox. P.D., V.C., and Z.W. developed and executed mathematical modeling and wrote and edited the manuscript. S.N. performed animal injections and necropsy. Z.L. provided assistance with PET-CT experiments. X.L. provided silicon particles. H.S. and M.F. edited the manuscript. Competing interests: M.F. and X.L. are inventors of a U.S. patent from The University of Texas-Houston (patent no. US 8,920,625 B2, granted 30 December 2014). M.F., H.S., and G.Z. are inventors of a nonprovisional patent application from the Houston Methodist Hospital (PCT/US 2014/0010879 A1). Both inventions are related to the scientific research described herein. More recently, M.F. has formed a company (BrYet LLC) and has secured certain commercial rights on said patents. Though not an executive officer, he remains a majority shareholder, director, and scientific advisor of the company. All other authors declare that they have no competing interests. Data and materials availability: All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. Additional data related to this paper may be requested from the authors.

Stay Connected to Science Advances

Navigate This Article