ORIGINAL RESEARCH

published: 24 June 2021

doi: 10.3389/fonc.2021.655781

Enhanced Drug Delivery to Solid

Tumors via Drug-Loaded
Nanocarriers: An Image-Based
Computational Framework
Farshad Moradi Kashkooli 1, M. Soltani 1,2,3,4*, Mohammad Masoud Momeni 1
and Arman Rahmim 5,6
1 Department of Mechanical Engineering, K. N. Toosi University of Technology, Tehran, Iran, 2 Department of Electrical and

Computer Engineering, Faculty of Engineering, School of Optometry and Vision Science, Faculty of Science, University of
Waterloo, Waterloo, ON, Canada, 3 Advanced Bioengineering Initiative Center, Multidisciplinary International Complex, K. N.
Toosi University of Technology, Tehran, Iran, 4 Centre for Biotechnology and Bioengineering (CBB), University of Waterloo,
Waterloo, ON, Canada, 5 Departments of Radiology and Physics, University of British Columbia, Vancouver, BC, Canada,
6 Department of Integrative Oncology, BC Cancer Research Institute, Vancouver, BC, Canada

Objective: Nano-sized drug delivery systems (NSDDSs) offer a promising therapeutic

technology with sufficient biocompatibility, stability, and drug-loading rates towards
efficient drug delivery to solid tumors. We aim to apply a multi-scale computational
Edited by: model for evaluating drug delivery to predict treatment efficacy.
Yuxia Tang,
Nanjing University, China
Methodology: Three strategies for drug delivery, namely conventional chemotherapy
Reviewed by:
(one-stage), as well as chemotherapy through two- and three-stage NSDDSs, were
Anju Gupta, simulated and compared. A geometric model of the tumor and the capillary network was
University of Toledo, United States obtained by processing a real image. Subsequently, equations related to intravascular and
Francesco Giansanti,
University of L’Aquila, Italy interstitial flows as well as drug transport in tissue were solved by considering real
*Correspondence: conditions as well as details such as drug binding to cells and cellular uptake. Finally, the
M. Soltani role of periodic treatments was investigated considering tumor recurrence between
msoltani@[Link]
treatments. The impact of different parameters, nanoparticle (NP) size, binding affinity of
Specialty section:
drug, and the kinetics of release rate, were additionally investigated to determine their
This article was submitted to therapeutic efficacy.
Pharmacology of Anti-Cancer Drugs,
a section of the journal Results: Using NPs considerably increases the fraction of killed cells (FKCs) inside the
Frontiers in Oncology tumor compared to conventional chemotherapy. Tumoral FKCs for two-stage DDS with
Received: 19 January 2021 smaller NP size (20nm) is higher than that of larger NPs (100nm), in all investigate release
Accepted: 26 May 2021
Published: 24 June 2021
rates. Slower and continuous release of the chemotherapeutic agents from NPs have
Citation:
better treatment outcomes in comparison with faster release rate. In three-stage DDS, for
Moradi Kashkooli F, Soltani M, intermediate and higher binding affinities, it is desirable for the secondary particle to be
Momeni MM and Rahmim A (2021)
released at a faster rate, and the drug with slower rate. In lower binding affinities, high
Enhanced Drug Delivery to Solid
Tumors via Drug-Loaded release rates have better performance. Results also demonstrate that after 5 treatments
Nanocarriers: An Image-Based with three-stage DDS, 99.6% of tumor cells (TCs) are killed, while two-stage DDS and
Computational Framework.
Front. Oncol. 11:655781.
conventional chemotherapy kill 95.6% and 88.5% of tumor cells in the same
doi: 10.3389/fonc.2021.655781 period, respectively.

Frontiers in Oncology | [Link] 1 June 2021 | Volume 11 | Article 655781

Moradi Kashkooli et al. Delivery of Drug-Loaded Nanocarriers to Tumors

Conclusion: The presented framework has the potential to enable decision making for new
drugs via computational modeling of treatment responses and has the potential to aid
oncologists with personalized treatment plans towards more optimal treatment outcomes.
Keywords: solid tumors, drug delivery, nanomedicine, drug-loaded nanocarriers, tumor penetration, image-based
model, treatment efficacy

INTRODUCTION small NPs experience higher transvascular and interstitial

transports (10). Investigation of the various functions of NPs
Encapsulation of drugs in nanoparticles (NPs) can enhance the with different physiological features in combating biological
delivery of therapeutic and diagnostic agents to tumors, and at barriers has demonstrated that NPs require to be dynamically
the same time, reduce their accumulation in healthy tissue (1). In adapted to these obstacles, resulting in the emergence and
treatments, this issue results in improved local drug development of multi-stage DDSs (11). In such a system, the
concentrations at disease sites while decreasing systemic NPs can change their size in response to stimuli at different
toxicity. On the other hand, conventional small molecule stages (e.g., initial size of 100 nm and secondary size of 10 nm)
compounds usually distribute randomly among all tissues (2). (12). Wong et al. (13) suggested a multi-stage system comprising
However, NPs using advantages such as enhanced permeability the primary particle containing smaller secondary particles,
and retention effect (EPR) have provided moderate which subsequently carry the therapeutic agent. Eventually,
enhancements over small molecule therapeutics regarding therapeutic agents contained in the secondary particles must
patient survival and care (3). A qualitative analysis of the be released to reach cancer cells. If an additional stage was
pharmacokinetics of NP-mediated doxorubicin (DOX) in included in the conventional NPs, it could improve drug
patients with cancer was carried out by Gabizon et al. (4) distribution into the tumor and also tumor penetration, as well
comparing with an equal DOX dosage delivered in free form. as increase treatment efficacy. The secondary particles release
The findings distinctly indicated that NP-mediated delivery is their cargo within the tumor, triggered by exposure to external
capable of decreasing plasma clearance and enhancing DOX (magnetic and electric fields, acoustic, etc.) or internal stimuli
concentration inside tumors (4). In vivo tests on human prostate (TME properties such as enzymes, pH, etc.).
carcinoma showed that NP-mediated DOX may improve Since the introduction of the first model for chemotherapy,
treatment effectiveness because of decreased systemic the use of mathematical and computational models has become
elimination, enhanced penetration in tumor, and prolonged widespread in examining different DDSs (14–19). These models
NP existence with a slow release rate of drug (5). The clinical can provide guidance on required composition and preparation
advantages of NP-encapsulated DOX is also evaluated in patients methods for administration of DDSs. Various computational
with metastatic breast cancer (6). The results suggested that NP- models have been employed for simulation of NSDDSs to
mediated DOX was capable of significantly reducing examine efficacy, understand biological phenomena, and select
cardiotoxicity and improving therapeutic efficacy. Nevertheless, optimal treatment plans. Based on the investigated spatial and
despite multiple promising ideas and experimental outcomes, temporal scales, these models are classified to discrete,
NPs do not commonly reach tumors efficiently in clinical trials continuous, and hybrid models (8). A summary is presented in
(7). Successful drug delivery depends on the ability of NPs to Table 1 of important studies conducted on employing
deep penetration into tissue, achieve an efficient spatial mathematical modeling of drug-containing NPs for drug
distribution, ensure its proper binding affinity, and adequately delivery to solid tumors (20–45).
release the drugs (8). However, progress in NPs will need a As seen in Table 1, there exist a number of gaps in the
number of technical and physiological barriers to be realized and literature regarding usage of drug-containing NPs for delivery of
overcome, such as opsonization and nonspecific protein drug to solid tumors. With this motivation, in the present study,
adsorption, non-specific uptake by cells and organs the main contribution is to apply a computational model of three
encompassing the immune system, targeting and penetrating DDSs, namely (i) a conventional chemotherapy system (one-
the tumor microenvironment (TME), and gaining access to stage DDS) and (ii) a conventional two-stage DDS containing
cancer cells for intracellular drug delivery (9). Physical NP and chemotherapy (two-stage DDS), and (iii) a three-stage
circumstances influencing the function of NPs facing these DDS containing a primary NP, a secondary NP, and the
obstacles are still understood imperfectly, and the chemotherapy; on a real image of vascularized tumor. To this
nanomedicine lacks a detailed explanation of physical end, intravascular and interstitial fluid flows as well as drug
principles to guide logical NP designs that can overcome transport equations are solved by considering actual conditions
biological barriers in TME (1). in tumor. Then, various parameters, NP size, binding affinity of
For most desirable efficiency, sufficient amounts of the drug ligands to the receptors of cells, and the kinetics of release
therapeutic agents must reach tumors in order to eliminate rate, are studied to determine their effects on treatment efficacy.
cancer cells, and at the same time, they must not induce NPs employed in this study could for instance be drug-loaded
considerable side effects on normal tissues. Generally, relatively nanocarriers, liposomes, or magnetic NPs. Additionally, the

Frontiers in Oncology | [Link] 2 June 2021 | Volume 11 | Article 655781

Frontiers in Oncology | [Link]

Moradi Kashkooli et al.

TABLE 1 | A summary of important studies conducted on employing mathematical modeling of drug-containing NPs for drug delivery to solid tumors.

Reference / Subject Geometry & Findings Study gap

Year simulation method

El-Kareh Developing a mathematical model and A PK/PD model Authors recommended that shorter injection duration might enhance treatment -Spatial distribution is not considered.
and Secomb applying it to compare the efficacy of performance, and explored it by computational studies. The treatment outcome was
(20)/2000 different administration modes of both predicated according to peak intracellular concentration over the whole treatment period.
free DOX and thermo-sensitive A comparison between bolus injection and continuous infusion demonstrated that
liposome (TSL) encapsulated DOX. duration of infusion had a great effect on the treatment outcome. Optimal duration is
dependent on cellular pharmacokinetics. Drug release rate from non-TSLs is an effective
parameter so that if this rate is optimized, the efficacy of non-TSLs is slightly fewer than
continuous infusion.
Zhang et al. Developing a mathematical model A 2D-0D* model & Compared to liposomes, diffusion of free drug plays a greater role in drug transport to -Avascular model;
(21)/2009 coupling heat and mass transfer to finite element tumor, as the free drug diffusivity is higher than that of liposomes. Hyperthermia alone -Normal tissue is not considered;
investigate spatiotemporal distributions method (FEM) only increases drug accumulation in the tumor periphery, and the TCs in the central area -Real image of tumor is not considered.
of drug that are released from the are barely damaged because of weak diffusion. Necrosis or apoptosis of the TCs can
liposome. importantly affect the penetration of drug and must be taken into account in modeling of
drug diffusion to precisely simulate the treatment effect. Combination of radio-frequency
ablation and liposomal DOX delivery demonstrates more efficacious therapeutic result,
particularly for larger tumors.
Hendricks et al. Presenting a multi-scale mathematical A PK model Authors illustrated that, for varying tumor transport features, there exist a regimen where -Spatial distribution is not considered.
(22)/2012 model of Liposomal DOX delivery for liposomal and conventional DOX deliver identical amounts of dox to tumor cell nuclei.
quantifying the role of parameters They also showed that liposome PKs and tumor deposition (which reflects vascular
related to tumor and drug in drug permeability) are highly variable.
delivery to solid tumors
Chauhan et al. Investigating the effect of normalizing A 2D-1D model & Decreasing the vessel-wall pore size via normalization reduces the IFP in tumors, allowing -Real image of tumor is not considered;
3

(23)/2012 blood vessels of tumor for enhancing FEM small NPs to enter them more quickly. However, enhanced steric and hydrodynamic -Concentration equations are solved
nanomedicine delivery in a size- hindrances, also associated with smaller pores, make it more difficult for large NPs to without taking into account the binding
dependent method enter tumors. It was suggested that smaller (∼12 nm) NPs are ideal for treating cancer affinity and cellular internalization.
because of their better penetration into the tumor.
Gasselhuber et al. Proposing a mathematical model for A PK model While stealth-DOX led to high concentrations in tumor in comparison with free-DOX, just a -Spatial distribution is not considered.
(24)/2012 comparison of Conventional minor fraction was bioavailable, resulting in little cellular uptake. Optimum time constants
chemotherapy, TSLs, and stealth of release for maximum cellular uptake for stealth-DOX and TSLs are obtained.
liposomes
Zhan and Xu Employing a mathematical modeling for A 2D-0D model & The model was applied to idealized geometry of tumor, and comparisons have been -Real image of tumor is not considered;
(25)/2013 TSL delivery of DOX to solid tumor FEM performed between continuous infusion of DOX and TSL-mediated delivery. Authors -Avascular model;
illustrated that TSL-mediated delivery performs better in reducing concentration of drug -Lack of capillary network and
in healthy tissues. Compared with direct infusion, TSL delivery results a much higher assumption of uniform vascular density

Delivery of Drug-Loaded Nanocarriers to Tumors

peak intracellular concentration of DOX, which may enhance fraction of killed cells in in the tissue.
tumor thereby improving the treatment impact of the drug.
June 2021 | Volume 11 | Article 655781

Stylianopoulos et al. Developing a mathematical platform for A 2D-1D model & The model simulations offer that electrostatic repulsion has a small effect on the -Real image of tumor is not considered
(26)/2013 NP delivery to solid tumors considering FEM transcapillary transport of NPs. Conversely, electrostatic attraction generated even by and instead, a mathematical model is
electrostatic interactions between the small cationic charges can result in a two-fold enhancement in the transvascular flux of used for angiogenesis;
NPs and the negatively-charged vessel- NPs into the tumor interstitium. For each size of NP, there exist an amount of charge -Concentration equations are solved
wall pores. density above which a sharp enhancement in transcapillary transport is simulated. ignoring binding affinity and cellular
internalization.
Kim et al. Overviewing different mathematical ─ (Review paper) Authors overviewed the state of mathematical modeling approaches that address ─
(27)/2013 frameworks of anti-cancer drug phenomena regarding drug delivery. They described how different types of models were
penetration into solid tumor employed to predict spatial-temporal drug distributions in solid tumor, to simulate various
approaches to overcome obstacles to drug delivery, or to optimize treatment programs.

(Continued)
Frontiers in Oncology | [Link]

Moradi Kashkooli et al.

TABLE 1 | Continued

Reference / Subject Geometry & Findings Study gap

Year simulation method

They also discussed how integration of in silico modeling with in vivo or clinical data
can provide better tools to understand the drug transport process.
Stylianopoulos et al. Employing mathematical modeling to A 2D-1D tumor Adjusting the release kinetics and binding affinities of drug results in enhanced drug -Real image of tumor is not considered
(28)/2015 examine the effect of drug features on model & delivery. Smaller NPs have better treatment efficacy than bigger ones. and instead, a 1D network is used for
the distribution and efficacy of NPs and FEM angiogenesis;
also investigating two multi-stage NP
delivery systems.
Stylianopoulos Design considerations for nano- ─ (Review paper) Authors evaluated different design parameters that can be regulated to optimize DDS, ─
and Jain therapeutics in oncology. suggested specific design approaches that should optimize delivery to most tumors,
(29)/2015 and discussed under which circumstances active targeting would be advantageous.
Chou et al. Developing a mathematical model of A 2D-0D tumor The efficacy of anti-cancer drug delivery was determined by the interplay of the -Lack of capillary network;
(30)/2017 tumor according to interstitial fluid flow model & microvascular density and NP size. All NPs and chemotherapeutic drugs have a limited -Real image of tumor is not considered;
and particle transport to study the drug FEM concentration in the necrotic zone of tumor, where transport of drug is only through -Concentration equations are solved
transport and cumulative diffusion. Using NPs as anti-cancer drug carriers is generally a better option compared to ignoring binding affinity and cellular
concentrations in a tumor. molecular chemotherapeutic agent due to its higher therapeutic efficacy on tumor and internalization.
lower damage to healthy tissue.
Zhan and Wang Investigating the convection-enhanced A 3D-0D model & Liposomes are able to increase the accumulation and penetration of drug in the -Avascular model;
(31)/2018 delivery of liposome containing DOX FVM convection enhanced delivery treatment. Transport of liposome is affected by convection -Lack of capillary network and
under different circumstances in an rather than diffusion. The effective delivery volume has nonlinear relation with the release- assumption of uniform vascular density
MRI-based brain tumor model. rate of drug. in the tissue.
Shamsi et al. Proposing a computational model for A 2D-0D tumor A great enhancement in the intratumoral concentration of magnetic NPs compared to -Avascular model;
(32)/2018 magnetically-assisted drug delivery model & free drugs. The success of magnetic drug targeting in larger tumors (10–20 mm in size) -Lack of capillary network and
4

approach to assess the penetration of FEM is found to be significantly due to the strength of magnetic field and tumor-magnet assumption of uniform vascular density
drug into peritoneal tumors nodules and distance while these two parameters are less important in small tumors. in the tissue.
improve intraperitoneal (IP) -Concentration equations are solved
chemotherapy. ignoring binding affinity and cellular
internalization.
Stylianopoulos et al Reengineering the TME to enhance the ─ (Review paper) Authors discussed the mechanics of both solid and fluid components of tumor, focusing ─
(33)./2018 efficacy of drug delivery from on how they prevent the delivery of drug and create an abnormal TME that promotes
computational modeling to bench to tumor growth and resistance to treatment. They also provide strategies to re-engineer
bedside. the TME by normalizing the vessels of tumor and the ECM to enhance the treatment of
cancer. Eventually, they summarized different mathematical approaches that have
provided insights into the physical obstacles against efficient cancer treatment and
suggested novel methods to overcome these impediments.

Delivery of Drug-Loaded Nanocarriers to Tumors

Huang et al. Presenting a mathematical modeling for A PK/PD model Authors demonstrated that complicated relationships between the related factors (various -Spatial distribution is not considered;
(34)/2019 spatial–temporal distribution of chemotherapy drugs, release rate constants, and heating duration) and the predicted
June 2021 | Volume 11 | Article 655781

chemotherapy drug in TSL-mediated treatment result, making it difficult to identify the best parameter set. a model-based
DDSs optimization approach is presented to overcome this challenge. Optimization showed
that the best result would be obtained with a low drug release rate at physiological
temperature, combined with a moderate to high release rate at mild hyperthermia and
1 h heating post- injection.
Rezaeian et al. IP injection of TSL DOX with the A 2D-0D tumor Using TSL-DOX delivery is efficacious than conventional chemotherapy. Adjusting the -Avascular model;
(35)/2019 triggered release by mild hyperthermia model & TSL size must be carried out according to the vessel wall permeability. Smaller TSLs -Lack of capillary network and
caused by high intensity focused FEM have better treatment efficacy. TSL-DOX delivery system in smaller tumors is less assumption of uniform vascular density
ultrasound. beneficial compared to larger ones. in the tissue.
Shamsi et al. A review of computational modeling of ─ (Review paper) Authors investigated different theoretical modeling approaches as influential tools to ─
(36)/2019 nano-engineered DDSs. furnish future design and development of DDSs.

(Continued)
Frontiers in Oncology | [Link]

Moradi Kashkooli et al.

TABLE 1 | Continued

Reference / Subject Geometry & Findings Study gap

Year simulation method

He et al. Developing a mathematical modeling to A PK model Authors quantified the effect of influencing parameters on the efficacy of tumor delivery, -Real image of tumor is not considered;
(37)/2019 analyze nanomedicine distributions in the magnitude of heterogeneous distribution, and the EPR effect. They also compared the -Spatial distribution is not considered;
solid tumors spatial distributions of the NPs and the free drugs within tumors. The model predicted -Lack of capillary network and
high degrees of distributional heterogeneity for both NPs and free drugs. They found that assumption of uniform vascular density
diffusion coefficient of NPs was the most efficient factor in decreasing the NPs in the tissue.
distributional heterogeneity but it has moderate impact on the free drugs.
Dogra et al. Overviewing different mathematical ─ (Review paper) Authors provided an overview on mathematical modeling works that have been applied ─
(38)/2019 modeling about application of towards a better insights of nano-bio interactions for enhancing the efficacy of drug
nanomedicine in cancer treatment. delivery to tumor.
Tehrani et al. Conducting numerical simulation to A 2D-0D tumor Injection process has an essential impact on distribution of MNPs. Sufficient diffusion time -Real image of tumor is not considered;
(39)/2020 investigate the impacts of diffusion of model & can enhance the ablation zone after thermal therapy. Balance between diffusion time and -Avascular model;
MNPs on microwave ablation FEM size of MNPs can enhance the efficacy of therapy. -Drug transport equations are solved
treatment. ignoring binding affinity and cellular
internalization.
Wirthl et al. Presenting a multi-phase tumor growth A 2D-0D tumor This study allows investigation of the properties and of the limitations of NP delivery to -Real image of tumor is not considered;
(40)/2020 model to examine NP delivery to solid model & solid tumors, which currently complicate the translation of NP therapy to a clinical trials -Drug transport equations are not
tumors FEM investigated.
Wijeratne and Proposing a mathematical framework of A 3D-1D tumor This model allows for drug features (e.g., size and binding affinity) to be explicitly defined, -Real image of tumor is not considered.
Vavourakis dynamic growth of solid tumor, drug model & thus facilitating investigation into the interaction between the changing TME and cytotoxic
(41)/2020 delivery, and angiogenesis. FEM and NP drugs. They predict a heterogeneous distribution of NPs after delivery; that NPs
need a ECM with high porosity to cause tumor regression; and that transcapillary fluid
velocity is dependent on porosity of ECM, and implicitly on the drug size.
5

Dogra et al. Conducting sensitivity analysis to Physiologically Degradation rate of NPs, size of NPs, blood viscosity of tumor, vascular fraction of tumor, -Real image of tumor is not considered;
(42)/2020 characterize the effective parameters on based PK model and tumor vascular porosity of tumor are effective factors in governing kinetics of NPs -Simulation is performed without taking
low delivery of NP to tumor and high within the interstitial space of tumor. into account the binding affinity and
off-target accumulation of NPs by cellular uptake.
whole-body NP pharmacokinetics
Shojaee et al. Effect of NP size, magnetic intensity, A 2D-2D tumor Magnetic field and size changes has a moderate impact on the drug penetration to the -Real image of tumor is not considered;
(43)/2020 and tumor distance on the distribution model & tumor. The dense ECM, elevated IFP, and the availability of capillary network have -Simulation is performed without
of the MNPs in a TME FEM negative influences on the MNP distribution. The size and the magnetic field are the two considering the binding affinity and
most promised factors for enhancing the convection term in the tumor area. cellular uptake.
Stillman et al. Presenting in silico modeling of ─ (Review paper) Authors investigated latest outcomes in multi-scale modeling of NP transport obstacles, ─
(44)/2020 nanomedicine for cancer, across scales as well as existing software packages, with the goal of focusing the wider research
and transport obstacles. community in building a common computational platform that able to overcome some

Delivery of Drug-Loaded Nanocarriers to Tumors

of the current barriers facing effective design of NPs.
Moradi A review of different mathematical ─ (Review paper) Investigation of various issues regarding the use of NPs as vehicles of anticancer drug ─
June 2021 | Volume 11 | Article 655781

Kashkooli et al. modeling approaches for NSDDSs. delivery: specifically, administration into the circulation system, transvascular transport,
(45)/2021 distribution in the extracellular matrix, cellular internalization, and release of drug from NPs.

*A 2D-0D means that the geometry of tumor and microvascular network are considered 2-dimensional and 0-dimentional (avascular), respectively.
Moradi Kashkooli et al. Delivery of Drug-Loaded Nanocarriers to Tumors

impact of successive treatment cycles is numerically examined into the blood, and subsequently, NP characteristics cause them
considering tumor recurrence between two consecutive to be released in accordance with a pre-designed controllable
treatments for three investigated DDSs. procedure. NP formulations have benefits over commonly used
chemotherapies since they can combine many diagnostic and
therapeutic factors (48). They are also associated with
significantly lower side effects, owing to their capacity for
MATERIAL AND METHODS optional accumulation in tumorous tissue. To eradicate TCs,
NPs have to first reach the tumor via the vascular system, then
In this section, a detailed description of mathematical models
extravasate from the relatively leaky regions of the microvessels
and their assessment are investigated. First, the mathematical
into the TME. Subsequently, NPs release their cargo in a
equations are presented, including equations governing the
controlled manner, while staying put at the ECM. Drug
interstitial fluid flow, solute transport and cellular uptake of
molecules diffuse into the tumor slightly more than into
chemotherapeutic agents, and NP delivery system. Then,
normal tissue and can bind to cancer cells and/or TME
parameters of model, relevant physical and transport
components eventually, becoming internalized by cells (28). A
properties, and NP-related calculations are also provided.
schematic of drug delivery mechanisms considered in the current
Finally, numerical methods to assess the mathematical models,
study is shown in Figure 1.
model parameters, boundary conditions (BCs), input images as
geometry, computational domain, as well as solution strategy
are presented. Governing Equations
The mathematical models for delivery of drugs to solid tumors
consist of equations as:
Delivery Mechanisms
Two well-known delivery mechanisms, namely chemotherapeutic i. Mass and momentum conservation for interstitial fluid flow;
delivery and delivery through drug loaded NP, are considered in ii. Mass transport for the free, bound, and internalized drug;
the present study. In the following, their mechanisms are iii. Mass transport for nano-encapsulated drug; and
described in detail. iv. Treatment efficacy for intracellular drug concentration.
After intravenous (IV) injection, the chemotherapeutic drugs
(here, DOX) are transported to the tumor site through the blood Detailed descriptions of each aspect of mathematical
vessels. Subsequently, they pass via the tumor vessel wall and modeling are provided next.
travel the remaining distance from the vessel wall to the cancer The dynamics process in drug delivery includes binding and
cells. In tumor ECM, free molecules of drug agents can bind to unbinding of drug agent ligands with receptors of the cells at the
receptors of the cells, unbind or get internalized (28, 46). rates of KON and KOFF, respectively; exchange of drug between
In general, only a small fraction of an injected drug reaches capillary network and interstitium, and influx/efflux of drugs
tumor tissue, the remaining being cleared from the body. One from interstitium to TCs, internalization of drugs to cellular
possible way to overcome this problem is to target drug delivery space, and finally cell-killing by drug agents. The cell-killing rate
by encapsulating the anticancer drug in a nanocarrier (47). In is governed by a model according to the predicted intracellular
such a system, drugs which encapsulated inside NPs are injected concentration of anticancer drugs. In the case of two-stage and

A B C

FIGURE 1 | Schematic of drug delivery mechanisms considered in the current study. (A) one-stage DDS or conventional chemotherapeutic delivery, (B) two-stage DDS
(i.e., NP delivery), and (C) three-stage DDS.

Frontiers in Oncology | [Link] 6 June 2021 | Volume 11 | Article 655781

Moradi Kashkooli et al. Delivery of Drug-Loaded Nanocarriers to Tumors

three-stage DDSs, additional equations describing the transport ∂C

of NPs in ECM are required. Mathematical modeling of fluid = ∇ · ½Deff ∇ C − ∇ · ½vi C + (FB − FL ) (6)
∂t
flow and solute transport in interstitium, convection-diffusion- in which ФB is the drug transport rate through the microvessels
reaction (CDR) modeling of drug transport in the extracellular into the interstitial space, and ФL is the drug transport rate from
space (for conventional chemotherapy, two-stage, and three- interstitial space into the lymph system. ФB is defined according
stage DDSs), as well as FKCs and tumor-cell survival equations to Patlak’s model, as (18, 49):
are presented in the following sections, respectively.
PS P
FB = fB (1 − sf )Cp + (Cp − Cf ) Pe e (7)
Fluid Flow in Interstitium V e −1
First, the interstitial fluid flow equations are solved to provide the
basic biomechanical environment for transport of drug. Darcy fB (1 − sf )
equation, which demonstrates the relationship between Pe = (8)
P VS
interstitial fluid velocity (IFV) and interstitial fluid pressure
(IFP), is employed to describe interstitial fluid flow in a porous Pe demonstrates the Peclet number, sf is the coefficient of
environment as follows (16): filtration reflection, P represents the permeability coefficient of
microvessels, and Cp is the injected drug concentration.
vi = −k ∇ Pi (1) The drug transport rate through lymphatic microvessels has
been considered only in healthy tissue as (18, 49):
where vi and k are the IFV and the hydraulic conductivity.
Considering the presence of source/sink terms in biological FL = fL C (9)
tissues, the continuity equation is modified as (7):
A bolus injection of chemotherapeutic agents, representing
∇ ·vi = fB − fL (2) the initial vascular concentration of the drug, is modeled as (50):
in which fB is the rate of fluid flow from the microvessels to the Cp = C0 exp ( − t=kd ) (10)
extracellular matrix (ECM) and and fL is the rate of fluid flow
from ECM to lymph system, defined as (7): where C0 and kd are the initial concentration and blood

 circulation decay, respectively.
S
fB = Lp ðPB − Pi − ss ðpB − pi ÞÞ (3) Convection-Diffusion-Reaction Modeling of Drug
V
Transport in the Interstitium

 The drugs exist in different forms as: NPs in the interstitial space
S
fL = LPL ðP − PL Þ (4) (CN), free drugs in the interstitial space (CF), bound drugs (CB),
V L i
and intracellular drugs (CINT) (20, 49). In tissues, with regards to
in which S/V demonstrates the surface area per unit volume of the existence of mass flow source/sink, a system of equations,
microvessels. LPL ( VS )L and PL are the lymphatic filtration which is called CDR equations, is utilized to represent the
coefficient and lymphatic pressure, respectively. process of drug delivery. The general block diagram of the
Combining Eq. (2) with Eq. (1), we arrive at: model of current study considering NPs and chemotherapeutic
drugs for two-stage DDS is shown in Figure 2.
−k ∇2 Pi = fB − fL
Conventional Chemotherapy


S The CDR equations for conventional chemotherapy are as
= Lp ðPB − Pi − ss ðpB − pi ÞÞ follows (17, 41):
V

 ∂ CF
S ∂t = −v ∇ CF + Deff ∇2 CF − j1 KON Crec CF + KOFF CB + (FB − FL )
− LPL ðP − PL Þ (5)
V L i Free drug
∂ CB
∂t = 1
j KON Crec CF − KOFF CB − KINT CB
(11)
Transport of Drug in the Interstitium Bound drug
The comprehensive model for drug transport in the interstitial ∂ CINT
= KINT CB
∂t
space includes:
Drug internalized into the cell
- Transport in ECM by diffusion and convection mechanisms, in which CF demonstrates free drug concentration, CB bound
- Transport across microvessels by diffusion and convection drug concentration, CINT internalized drug concentration, and
mechanisms, and Crec is the concentration of cell-surface receptor. In this
- Binding to cells and internalization. equations, v is the IFV, D is the diffusion coefficient, KON is
the association rate, KOFF is disassociation rate, KINT, represents
The drug transport equation for biological tissue can be the cellular internalization; and j demonstrates the tumor
written as (18): volume fraction available for the drugs.

Frontiers in Oncology | [Link] 7 June 2021 | Volume 11 | Article 655781

Moradi Kashkooli et al. Delivery of Drug-Loaded Nanocarriers to Tumors

FIGURE 2 | Block diagram of the current study for computational modeling of drug transport of two-stage DDS.

Two-Stage Drug Delivery System in which CN1 and CN2 are the primary and secondary NP
Systemically administered NPs, as demonstrated in Figures 1 concentrations, respectively. Krel1 is the rate constant for the
and 2, are transported to tumor sites through the circulation release of the secondary NP from the primary one and Kel2 the
system, undergo transvascular extravasation followed by rate constant for the drug release from the secondary NP. a and
distribution in the interstitial space, and are finally delivered to b are the number of secondary NPs released by the primary and
cancer cells. For transport of an NP containing drug particles released by the secondary NP, respectively.
chemotherapeutic agents (CN), the system of equations is
adjusted thus (51): Fraction of Killed Cells
Despite its cardiotoxicity, the drug DOX is frequently employed
∂ CN
∂t = −vi ∇ CN + DN ∇2 CN − Krel CN + (FV − FL ) in tumor treatment. DOX is a standard-of-care, DNA-damaging
Nano − carrier agent employed in the treatment of multiple tumors (e.g.,
∂ CF
bladder, breast, and lung cancers). Using the internalized drug
∂t = aKrel CN − vi ∇ CF + D ∇2 CF − j1 KON Crec CF + KOFF CB concentration, the efficacy of drugs was calculated according to
Free drug the empirical equation obtained for DOX (52). The FKCs
∂ CB
(12) parameter is defined, as (53):
∂t = 1
j KON Crec CF − KOFF CB − KINT CB
Bound drug FKC = 1 − SF = 1 − exp ( − w · CINT ) (14)
∂ CINT
∂t = KINT CB in which SF is the fraction of cells remaining after treatment and
Internalized drug w is a fitting parameter specified for DOX based on the results of
experiments (54).
in which Krel, DN, and a are respectively the drug release rate
from the carrier, its diffusion coefficient, and the number of Tumor-Cell Survival
chemotherapy molecules contained in the nanocarrier. The number of TCs after a period of time (ni), which is obtained
through Gompertz’s model, depends on intervals between
Three-Stage Drug Delivery System chemotherapy sessions (t) (55). Gompertz equation is a
Three-stage NP drug delivery system is an efficient option to function of three parameters (as demonstrated in Eq.)15(): the
overcome the barriers of two-stage DDS in ECM. The system of number of TCs surviving after the ith treatment (Ni), the number
equations for three-stage DDS, as demonstrated in Figure 1C, is of saturated cells after a very long period (N∞), and eventually the
as following (32): rate of tumor progression (b).
∂ CN1

 
∂t = −vi ∇ CN1 + DN1 ∇2 CN1 − Krel1 CN1 + (FV − FL ) N∞
ni (t) = Ni exp Ln ½1 − expð−bt Þ (15)
Nano − carrier 1 Ni
∂ CN2
∂t = aKrel1 CN1 − vi ∇ CN2 + DN2 ∇2 CN2 − Krel2 CN2 The number of surviving TCs after each therapeutic phase
Nano − carrier 2 will be examined as a criterion of treatment efficacy evaluation.
∂ CF This criterion is obtained by using the FKCs according to Eq.
∂t = bKrel2 CN2 − vi ∇ CF + D ∇2 CF − j1 KON Crec CF + KOFF CB
(13) (14). The number of remaining TCs as a result of the difference
Free drug between the number of cells after (Ni) and before (N0) treatment
∂ CB
∂t = 1
j KON Crec CF − KOFF CB − KINT CB will also be examined as a criterion for evaluating treatment
Bound drug
efficacy. In the present model, SF is also defined, as (30):
∂ CINT
∂t = KINT CB Ni
SF = (16)
Internalized drug N0

Frontiers in Oncology | [Link] 8 June 2021 | Volume 11 | Article 655781

Moradi Kashkooli et al. Delivery of Drug-Loaded Nanocarriers to Tumors

The initial numbers of TCs here, adopted from York et al. TABLE 2 | Parameters of interstitial transport used in numerical simulations.
(55), are N0 = 5 × 109, N∞ = 3.1 × 1012, and b = 0.0283 month-1. Parameter Unit Description Value Ref.
The number of cells for healthy tissue was considered to be N1 =
4.64 × 1012 (55). The cell number is assumed to depend on tumor pB [mmHg] Oncotic pressure of 20 (Normal) (57)
microvessels 20 (Tumor)
volume; i.e., when the cell number decreases, the tumor shrinks
pi [mmHg] Oncotic pressure of interstitial 10 (Normal) (57)
(R reduces). The ratio of the density of healthy tissue to the fluid 15 (Tumor)
density of TC is considered to be 0.2 (56), and it is also assumed ss – Coefficient of average osmotic 0.91 (57)
that the microvascular density distribution in the computational reflection (Normal)
field does not change after each treatment. In addition, the 0.82 (Tumor)
Lp [cm/ Hydraulic conductivity of the 0.36×10-7 (57)
regrowth of healthy tissue cells was examined by using Eq. (15)
((mmHg)*s)] microvessel wall (Normal)
with the assumption that the growth rate (b) of healthy tissue is 2.8×10-7
half that of the tumor (30). (Tumor)
LpLSL/V [1/ Coefficient of Lymph filtration 1.33×10-5 (58)
(mmHg*s)] (Normal)
0 (Tumor)
Model Parameters k [cm2/ Hydraulic conductivity of 8.53×10-9 (58)
Interstitial and Drug Transport Parameters (mmHg*s)] interstitium (Normal)
Tables 2 and 3 demonstrate both the interstitial and DOX drug 4.13×10-8
transport parameters, respectively, including both tumor and (Tumor)
healthy tissues. Table 4 represents the baseline state parameters PL [Pa] Hydrostatic pressure of lymph 0 (58)
vessels
for NSDDS for 20 nm particle size and 200 nm vessel-wall pore
(VWP) size.
where F is the partition coefficient and is defined as (53, 60):
NP-Related Calculations
F = (1 − l)2 (23)
The hydraulic conductivity, vascular permeability, and reflection
coefficient are calculated for the NPs by using the theory of in which l is the ratio of drug particle size to the VWP size, as
particle transport through cylindrical pores (28, 59): follows (53, 59):
rs
g ro2 l= (24)
Lp = (17) ro
8mL
Kt and Ks coefficients in Eqs. (21) and (22) are given by (53,
g HDo 59):
P= (18)
L

sf = 1 − W (19)
in which g is the fraction of the surface area of a porous vessel- TABLE 3 | Parameters for chemotherapy drug (DOX) applied to computational
wall, ro is the pore radius, h is the viscosity of water at 310 K, and modeling.
L is the vessel-wall thickness. Do represents a particle diffusion
Parameter Unit Description Value Ref.
coefficient in a free solution at 310 K, given by the Stokes-
2 -10
Einstein relationship as in (53, 59): D [m /s] Coefficient of diffusion 1.58×10 (58)
(Normal)
Kb T 3.40×10-10
Do = (20) (Tumor)
6phrs P [m/s] Microvessel permeability 3.75×10-7 (58)
in which Kb, T, and rs are the Boltzmann constant, the coefficient (Normal)
3.00×10-6
temperature, and the diffusing particle radius, respectively.
(Tumor)
H and W are diffusive and convective hindrance factors, sf – Filtration reflection coefficient 0.35 (52)
respectively, and related to hydrodynamic and electrostatic KON [m3/ Binding rate constant 15 (32)
interactions. Neglecting electrostatic interactions (E=0), H and (mole s)]
W are reduced to (53, 59): KOFF [1/s] Unbinding rate constant 8×10-3 (32)
KINT [1/s] Internalization rate constant 5×10-5 (32)
j – Volume fraction of tumor 0.4 (32)
6pF
H= (21) available to drugs
Kt Crec [M] Cell-surface receptors 1×10-5 (32)
concentration
Kd [Min] Half-life of drug in plasma 6 (32)
F(2 − F)Ks
W= (22) w [m3/ Survival constant of cancer cells 0.6603 (51)
2Kt mole]

Frontiers in Oncology | [Link] 9 June 2021 | Volume 11 | Article 655781

Moradi Kashkooli et al. Delivery of Drug-Loaded Nanocarriers to Tumors

TABLE 4 | Parameters of baseline state for NP drug delivery for 20 nm particles and 200 nm VWP size.

Parameters Unit Description Value Ref.

2 -12
D [m /s] Diffusion 7×10 (30)
j - Volume fraction of tumor available to drugs 0.05 (32)
KON [m3/(mole s)] Binding rate constant 15 (32)
KOFF [s-1] Unbinding rate constant 8×10-3 (32)
KINT [s-1] Cellular uptake rate constant 5×10-5 (32)
Kel [s-1] Release rate constant 2.1×10-6 (32)
Kd [min] Blood circulation decay constant 1320 (32)
a – Number of particles in the NP carrier 20 (32)
Crec [M] Concentration of cell-surface receptors 1×10-5 (32)

! " ! #
Kt 9 2 pffiffiffi 2 an Two boundaries are assumed in the computational field: the
− 52
= p 2(1 − l) 1 + o (1 − l)n
inner boundary (between the tumor and the normal tissue), and
Ks 4 n=1 bn
the outer boundary (at the outer edges of computational field). For
! the inner boundary, the continuity BC is considered for IFV, IFP, as
4 an+3
+o l2 (25) well as concentration and its flux. For the outer boundary, in which
n=0 bn+3 the IFP is constant, the Dirichlet BC is used for fluid flow, and the
open boundary is employed for concentration (50). The input and
Parameter values used in the Eqs. (17) to (25) are
output pressure amounts for the parent vessels are selected as:
demonstrated in Table 5.
PInlet,1 = 25 mmHg, PInlet,2 = 25 mmHg, and POutlet =10 mmHg,
Model Geometry and Boundary Conditions according to realistic physiological conditions reported in literature
In this study, a tumor model (as shown in Figure 3A) is (50, 61).
employed as an input geometry based on real image of a
tumor with a capillary network surrounded by healthy tissue, Solution Strategy
extracted from Roudnicky et al. (60). About the circumstance of There exist two distinct phases for solving the present problem:
this input image, it should be mentioned that this type of tumor steady-state and transient. Calculations related to blood flow in a
is inoculated in a mice by injecting human A431 squamous cell model with a discretized capillary network in the computational
carcinoma (SCC) cells (60). This image was taken about 12 days field provides a system of non-linear equations. Thus, an iterative
after tumor inoculation. In fact, 2 days after inoculation, approach was applied to solve the fluid flow equations. The blood
thrombospondin-2 (TSP2), an anti-angiogenic matricellular flow and interstitial fluid flow were solved concurrently, where
protein that inhibits tumor growth and angiogenesis, was IBP and IFP were coupled via Starling’s equation (Eq. (3) in
injected for 10 days and subsequently this image was taken. supporting file). Obtained values for IBP, IFP, and IFV were
After image-processing, a computational field is considered with employed for solving the transient equations of CDR to achieve
the existence of a tumor in the middle of the domain as well as different drug concentrations (CN1, CN2, CF, CB, and CINT) as well
the parent vessels (Figure 3B). as FKCs.

TABLE 5 | Parameter values used for NP-related calculations.

Parameter Description Value Ref.

−6
L Vessel-wall thickness 5×10 m (31)
h Water viscosity at 310K 7×10−4 Pa∙s (31)
g Fraction of surface area of vessel-wall occupied by pores 1×10−4 [-] (51)
a1 1st coefficient for Kt -73/60 [-] (59)
a2 2nd coefficient for Kt 77.293/50.400 [-] (59)
a3 3rd coefficient for Kt -22.5083 [-] (59)
a4 4th coefficient for Kt -5.617 [-] (59)
a5 5th coefficient for Kt -0.3363 [-] (59)
a6 6th coefficient for Kt -1.216 [-] (59)
a7 7th coefficient for Kt 1.647 [-] (59)
b1 1st coefficient for Ks 7/60 [-] (59)
b2 2nd coefficient for Ks -2.227/50.400 [-] (59)
b3 3rd coefficient for Ks 4.0180 [-] (59)
b4 4th coefficient for Ks -3.9788 [-] (59)
b5 5th coefficient for Ks -1.9215 [-] (59)
b6 6th coefficient for Ks 4.392 [-] (59)
b7 7th coefficient for Ks 5.006 [-] (59)

Frontiers in Oncology | [Link] 10 June 2021 | Volume 11 | Article 655781

Moradi Kashkooli et al. Delivery of Drug-Loaded Nanocarriers to Tumors

FIGURE 3 | (A) Real image of tumor, and (B) computational field considered in numerical simulation which is obtained by image-processing of realistic image.

After image-processing, geometries were meshed and different parameters are investigated in detail. Subsequently, results
analyzed utilizing the COMSOL Multiphysics software-version of computational modeling for delivery of NPs are presented.
5.5a. The coupled nonlinear set of the above-mentioned
governing equations and also the BCs were assessed through Validation of the Results
FEM. A segregated approach is applied to solve the equations In this study, validation is performed with the same governing
with the time-step of 0.1 [s] and relative tolerance of 0.001. A six- equations as the literature (28) for the FKCs over time (Figure 4).
fold drop of residuals is chosen as the criterion for convergence. As is clear, there is a correspondence between the results of the
For solving drug delivery equations in vascularized tumors, a present study and those in the literature so that by considering
Core (TM) i24 CPU @ 3 GHz with 32 GB RAM system is used. the real geometry and physics, the FKCs has a similar trend.
However, considering similar conditions, its value has about 4%
difference, due to differences in tumor geometry, computational
domain, and structure of capillary network. One should note the
RESULTS AND DISCUSSION 2D geometry of the tumor and capillary network investigated in
this study, while Stylianopoulos et al. (28) utilized 1D for the
Characteristics of tumors are determinant factors in transport of capillary network.
drug and final efficacy of treatment. Since chemotherapy drugs
are carried by the circulatory system, the tumor vasculature Baseline Model Analysis
features play an important act in delivery of drug. In the present The spatial-temporal distribution of the non-dimensional
study, a tumor model incorporating details of capillary network concentrations of drug that is taken up by TCs in the tumor and
distribution is employed to evaluate the impact of heterogeneous surrounding normal tissue for three investigated DDSs ─single-
distribution of capillary network on delivery of drug. This is followed stage (conventional chemotherapy), two-stage, and three-stage ─
by studies of transport of drug in tumors with microvascular density are shown in Figures 5–7. These non-dimensional concentrations
(MVD). First, the results of chemotherapeutic agent delivery for a (Cei ) were obtained at a given time by dividing that concentration at
case study, extracted from real image of tumor, are proposed, and any point by the maximum concentration in the entire field. Total

Frontiers in Oncology | [Link] 11 June 2021 | Volume 11 | Article 655781

Moradi Kashkooli et al. Delivery of Drug-Loaded Nanocarriers to Tumors

FIGURE 4 | Comparison of the results with previously published study (28) using FKCs over time.

concentration is defined as the sum of various drug concentrations. Analysis of NP Drug Delivery System
In single-state DDS, free drug concentration reduces over time and Parameters of the TME that inhibit the delivery of NPs into the
the drug gradually turns into a bound drug. Then, the bound drug tumor include the size-dependency of the transport both across
slowly enters the intracellular space and is consumed there. It is also the tumor microvessel wall and then via the interstitial space of
obvious that drug concentration in the tumor is greater than that in tumor. Transport across the tumor vessel wall is determined by
the normal tissue. The reason for this is the greater extravasation the relative size of the particle compared to the VWP size. On the
rate from microvasculature in the tumor area and additionally the other hand, not only the size but also other parameters of drug,
highly-dense MVD in this region. The highest value for total such as the drug release kinetics might play a crucial role in the
concentration takes place in the tumor zone and this value is outcome of the therapy. Therefore, one of the main goals of the
several times higher than the concentration in the normal tissue. present study is to determine under what circumstances two-
In two-stage DDS, the loads of NPs are released in the interstitium stage and three-stage DDSs can be beneficial, relative to one
in a free drug form and the remaining process is similar to single- another or to conventional chemotherapy. In the following, the
state system. In three-stage DDS, the primary NPs release the results are presented for two- and three-stage DDSs. Then, a
second NPs, and the remaining process is similar to previous parameter study is carried out to examine the effect of three
two-stage DDS.C eINT represents the internalization of drugs to the important parameters ─size of NPs, binding affinity, release rate
cellular space, determining the succeed of drug delivery in killing of drug─ on the real geometry of tumor.
cancer cells (i.e., leading to higher FKCs). Another important factor
for efficient drug delivery is uniform distribution of drug in tumor, Two-Stage Drug Delivery System
expressing the penetration depth of drugs extravasated from The effect of release rate of drug for two sizes of NPs are
microvascular network. Comparison of the distribution of demonstrated in Figure 8. Overall, FKCs for NPs with the size
different drug concentrations in three investigated systems of 20nm is higher than that of 100nm, in all investigate release
demonstrates that three-stage DDS provides much more uniform rates of drug. The main reason is the greater blood half-life of
drug distribution than the other two ones. After that, two-stage smaller NPs compared to larger ones. Therefore, a 20nm NP is
system has more efficient drug distribution compared to single- expected to circulate in the blood for a longer time in comparison
stage DDS. From Figures 5–7, it is also clear that in conventional with a 100nm particle, further improving efficacy of delivery
chemotherapy, there exist a small concentration of drug in healthy systems. Moreover, it is demonstrated that slower and
tissue, leading to side effects; whereas, in both NSDDSs, the side continuous release of the chemotherapeutic agents from NPs
effects are negligible. have better treatment results compared to faster release rate.

Frontiers in Oncology | [Link] 12 June 2021 | Volume 11 | Article 655781

Moradi Kashkooli et al. Delivery of Drug-Loaded Nanocarriers to Tumors

FIGURE 5 | Spatiotemporal distributions of concentrations of chemotherapy drug in tumor and its surrounding normal tissue with increasing time. These non-dimensional
e
concentrations C INT were calculated at a given time by dividing that concentration at any point of geometry by the maximum concentration in the whole domain.

FIGURE 6 | Spatiotemporal distributions of concentrations of drug-loaded NPs in tumor and its surrounding normal tissue with increasing time.

NP delivery systems may have other functions in addition to and anticancer activity of a DDS. Depending on various factors
acting as carriers of drug. One of these functions is to control the such as the NP formulation, fabrication method, surrounding
release rate of therapeutic agents from the NPs (2). Release rate environment, etc., release rate can vary across a wide range. For
from drug-loaded nanocarriers directly determines the toxicity instance, stealth liposomes can provide sustainable release over

Frontiers in Oncology | [Link] 13 June 2021 | Volume 11 | Article 655781

Moradi Kashkooli et al. Delivery of Drug-Loaded Nanocarriers to Tumors

FIGURE 7 | Spatiotemporal distributions of concentrations of three-stage NPs in tumor and its surrounding normal tissue with increasing time.

weeks, while TSLs are designed to release their payloads in a reduces the need for higher doses of the drug, which results in
short time. Hence, to cover all these wide ranges, the release rate reducing side effects.
is changed from 2.1×10-6 to 1×10-3 [1/s] in this study. A very
rapid drug release before the NPs have penetrated deep into the Three-Stage NP Drug Delivery System
tumor may result in non-uniform drug distribution, while a very In addition to enhancing the tumoral penetration depth, multi-
slow release may cause most NPs to be cleared out before stage NP delivery systems can provide further tunability in the
reaching the tumor as well as it may cause multi-drug spatial delivery control to solid tumors (8, 28). Due to the
resistance (29). Moreover, release rate of drug should be physiological obstacles that a chemotherapeutic agent must
regulated based on maximum efficacy in TC kill and minimum encounter, a multi-stage DDS can improve treatment efficacy
side effects. Sustained release with the goal of drug delivery over a by altering its physical properties (including shape, size,
long period of time is important for drugs that are rapidly flexibility, charge, and/or surface coating) to suit the transport
metabolized and excreted. Sustained release can stabilize across each obstacle (8). As interest in increasingly complicated
plasma concentrations of the drug at a constant level, thus DDSs grows, we face a corresponding challenge to set the model

A B
1 1

0.8 0.8

0.6 0.6
FKCs

FKCs

0.4 0.4

0.2 0.2

0 0
2.1E-6 1.18E-4 1.0E-3 2.1E-6 1.18E-4 1.0E-3

Release rate (1/s) Release rate (1/s)

FIGURE 8 | Comparison of treatment efficacy of two NP sizes for different release rates. KON was set to 15[m3/[(mol.s)] in all cases. (A) 100 nm, (B) 20 nm.

Frontiers in Oncology | [Link] 14 June 2021 | Volume 11 | Article 655781

Moradi Kashkooli et al. Delivery of Drug-Loaded Nanocarriers to Tumors

parameters. In this study, we had presented two common be affected (29). Indeed, there is a competition between diffusion
scenarios for multi-stage DDS: (i) a 100-nm particle, which in the ECM (which depends on NP size) and binding affinity
released secondary 10-nm particles; and (ii) a 20-nm particle, and/or the release rate of drugs. However, the one of the main
which released secondary 5-nm particles. conclusions is that all these three parameters should be
Based on Figure 9 for the three-stage DDSs, we assumed considered together and studying their impact alone does not
multiple different possible states for release rate constants (i.e., provide the right policy on optimal design of NSDDSs.
Kel1 and Kel2) of the secondary particle and the drug. The results
enable a set of observations as to how release rate constants affect Treatment Evaluation
the efficacy of drugs. The four most predominant are: (i) overall, In this section, the treatment outcomes of three investigated drug
the second scenario (NP1 = 20nm and NP2 = 5nm) has higher delivery systems are calculated by considering 1-month drug-free
treatment efficacy, almost in all investigated states, compared to breaks between treatments. As demonstrated in Figure 10, three-
the first scenario (NP1 = 100nm and NP2 = 10nm); (ii) in low stage system has better treatment outcome than two-stage and
binding rates, the high release rates have better performance; one-stage (i.e., conventional chemotherapy) systems, implying
(iii) in moderate and high binding rates, the NP release must superiority based on higher rate of killing TCs and shorter
have high release rates and the drug release must have the lower treatment time, simultaneously. For the conventional
release rates; (iv) the least treatment efficacy occurs when both chemotherapy, the respective efficacies of the first, second, and
release rates are slow. third stages of treatment in reducing the TCs are about 36.78%,
A high binding affinity of drug agents to receptors of cancer- 29.67%, and 24.41%. This rate is 20.24% for the fourth stage,
cell leads to a higher drug concentration in the intracellular space 14.9% for the fifth stage, 5.9% for the sixth stage, 4.65% for the
of tumor, and it simultaneously decreases the drug concentration seventh stage, 3.42% for the sixth stage, 2.73% for the ninth stage,
in tumor tissue. In other words, with a higher binding affinity, and 2.11% for the tenth stage. Tumor regrowth between the first
more drugs are taken up by cancer cells, reducing the drug level and second stages of treatment, between the second and third
in the interstitium. Binding affinity plays different roles in the stages, and up to the end of the process in the tenth stage are
delivery of smaller or larger particles because in smaller particles calculated to be 13.28%, 10.57%, 8.06%, 5.59%, 3.4%, 2.85%,
the diffusion mechanism is more dominant than convection, 2.4%, 2.17%, and 2.03%, respectively. At the last stage of
while in larger particles the convection via vessel-walls is more treatment, it is obvious that the regrowth of TCs is about
dominant in the transport across blood vessels. It should be 2.03%, while the treatment efficacy is about 2.11%, implying
mentioned that a very high binding affinity of the NPs leads to that for this size of tumor using conventional chemotherapy
aggregation nearby the vessels that NPs are extravasated from. alone is not the best choice and need more cycles of treatments,
On the other hand, for transport from microvessels to tissue, NP so an adjuvant therapy would be suggested at this size of tumor.
efficacy depends on the kinetics of drug release from particles In general, after seven treatments, 5.54% of the initial tumor cells
(28). Eventually, the released drug may also rapidly bind to the remain for conventional chemotherapy.
cells, causing heterogeneous and incomplete distribution within Results indicate that after 5 treatments with three-stage
tumor (62). After simultaneous investigation of the impacts of system, 99.6% of TCs are killed, while two-stage and one-stage
NP size, binding affinity, and the release rate of drugs, it is system respectively kill 95.6% and 88.5% of TCs in the same
predicted that NP penetration from tumor microvessels would period. It should be mentioned that the rate of killing TCs for the

A B
Kel1 = Kel2 =1E-3 Kel1 = Kel2 =1E-3
Kel1 = Kel2 = 2.1E-6 Kel1 = Kel2 = 2.1E-6
1.0 Kel1 = 2.1E-6, Kel2 = 1E-3 1.0 Kel1 = 2.1E-6, Kel2 = 1E-3
Kel1 = 1E-3, Kel2 = 2.1E-6 Kel1 = 1E-3, Kel2 = 2.1E-6

0.8 0.8

0.6 0.6
FKCs
FKCs

0.4 0.4

0.2 0.2

0.0 0.0
0.15 15 1500 0.15 15 1500
3
K_ON (m /(mol.s)) K_ON (m3/(mol.s))

FIGURE 9 | Comparison of treatment efficacy of two multi-stage delivery scenarios for different values of KON, Kel1, and Kel2. (A) NP1=100nm and NP2=10nm,
(B) NP1=20nm and NP2=5nm.

Frontiers in Oncology | [Link] 15 June 2021 | Volume 11 | Article 655781

Moradi Kashkooli et al. Delivery of Drug-Loaded Nanocarriers to Tumors

FIGURE 10 | Survival rate of TCs for one (conventional), two and three-stage drug delivery systems.

two-stage system reaches 98.9% after 7 cycles of treatment. This investigate release rates. Slower and continuous release of the
rate is also 94.46% for one-stage system after 10 treatment cycles. chemotherapeutic agents from NPs have better treatment
As a consequence, the treatment efficacies for both the NSDDSs outcomes in comparison with faster release rate.
have shown significant improvements compared to conventional • For three-stage DDSs, in intermediate and higher binding
chemotherapy. Results of the case study show that the survival affinities, it is desirable for the secondary particle to be
rates of TCs after several treatment cycles offer significant released with faster rate, and the drug with slower rate. In
prognostic insight about the survival rate and the cell regrowth lower binding affinities, the high release rates have better
percentage. Being able to evaluate the efficacy of a treatment performance
scenario using multi-scale computational modeling is very • Three-stage system has better treatment results relative to
helpful in clinical situations. two-stage and one-stage systems, reaching 99.6% effectiveness
in killing TCs after 5 treatments; while two-stage and one-
stage system respectively kill 95.6% and 88.5% of TCs in the
CONCLUSION same period.

In the present work, a mathematical model of drug transport is Overall, a mathematical framework has been developed for
proposed to predict delivery of chemotherapeutic drugs, either in drug delivery in both free and encapsulated forms to provide
their free form or encapsulated in NPs. Parameters of the model qualitative and mechanistic understanding of transport of drug
describing physiological and biological characteristics of tissue in solid tumors. The effect of treatment outcomes considering
and drug are extracted from experimental data in the literature. tumor recurrence between the two presented models is highly
The model has been applied to realistic tumor reconstructed complementary to PK/PD models as it considers both the spatial
from actual image of vascularized tumor to: (i) understand the and temporal scales at the same time, while PK/PD models
transport steps of non-encapsulated and encapsulated drugs, and merely consider the temporal scales. Moreover, different
(ii) elucidate the impact of drug properties on drug delivery and equation parameters in the present study have physiological
treatment. The main findings are as follows: and bio-chemical implications (e.g., intravascular pressure,
microvascular density, microvascular diameter, permeability,
• Spatiotemporal drug distribution illustrates that the diffusion coefficient to name a few), whereas PK/PD methods
complexity of the capillary network is the main factor for might involve a parameter to express several parameters.
non-uniform drug distribution in the tumor. It should be mentioned that anti-angiogenic agent (here,
• The FKCs of tumor for two-stage DDS with smaller size of TSP2) has effects on real tumor geometry and changes of
NPs (20nm) is higher than that of larger ones (100nm), in all permeability, but we have no discussion about TSP2 in the

Frontiers in Oncology | [Link] 16 June 2021 | Volume 11 | Article 655781

Moradi Kashkooli et al. Delivery of Drug-Loaded Nanocarriers to Tumors

present study. On the other hand, if we have three images for DATA AVAILABILITY STATEMENT
three different states (negative control, treatment with TSP2 and
without treatment), the mathematical modeling can include the The original contributions presented in the study are included in the
new geometry of microvascular networks. This is an interesting article. Further inquiries can be directed to the corresponding author.
suggestion for future studies. Additionally, this method can be
applied on each images of tumor because the geometry of tumor
can be obtained after image-processing; however, in this study,
we have chosen an actual image from the literature due to the AUTHOR CONTRIBUTIONS
lack of experimental set-up in our group.
Protection systems for NPs such as PEGylation as well as Data curation, FM. Investigation, FM, MS, and MM.
adjustable parameters of NP design, all try to enhance drug half- Methodology, FM. Project administration, FM and MS.
life in the circulation system, which can increase drug delivery Resources, FM. Software: FM and MM. Supervision, MS and
into the tumor and reduce side effects to healthy tissue. AR. Validation: FM. Visualization, FM and MM. Writing –
Furthermore, NPs can be manipulated to release their loaded original draft, FM. Writing – review and editing, MS, MM, and
drugs if exposed to a specific external (e.g., ultrasound or AR. Funding acquisition, AR. All authors contributed to the
magnetic field) or internal (e.g., pH or enzyme) stimulus. Multi article and approved the submitted version.
stimuli-responsive DDSs ─combination of internal and external
stimuli─ have not only succeeded in targeted drug delivery but
also in the multi-modal cancer diagnosis and treatment. In
addition to enhancing efficacy of encapsulation, these systems FUNDING
may enhance the drugs half-life. Our presented approach has the
potential of modeling these targeting systems. For example, by The authors wish to acknowledge funding from the Natural
simultaneously solving bio-heat and wave equations with the Sciences and Engineering Research Council of Canada (NSERC)
presented equations, the drug delivery through thermo-sensitive Discovery Grant RGPIN-2019-06467, as well as the BC
NPs can be modeled. Cancer Foundation.

REFERENCES 10. Popović Z, Liu W, Chauhan VP, Lee J, Wong C, Greytak AB, et al. A
Nanoparticle Size Series for In Vivo Fluorescence Imaging. Angewandte
1. Ng TSC, Garlin MA, Weissleder R, Miller MA. Improving Nanotherapy Chemie (2010) 122:8831–4. doi: 10.1002/ange.201003142
Delivery and Action Through Image-Guided Systems Pharmacology. 11. Stylianopoulos T, Wong C, Bawendi MG, Jain RK, Fukumura D. Multistage
Theranostics (2020) 10(3):968–97. doi: 10.7150/thno.37215 Nanoparticles for Improved Delivery Into Tumor Tissue. Methods enzymol
2. Moradi Kashkooli F, Soltani M, Souri M. Controlled Anti-Cancer Drug Elsevier (2012) 508:109–30. doi: 10.1016/B978-0-12-391860-4.00006-9
Release Through Advanced Nano-Drug Delivery Systems: Static and 12. Niu Y, Zhu J, Li Y, Shi H, Gong Y, Li R, et al. Size Shrinkable Drug Delivery
Dynamic Targeting Strategies. J Controlled Release (2020) 327:316–49. doi: Nanosystems and Priming the Tumor Microenvironment for Deep
10.1016/[Link].2020.08.012 Intratumoral Penetration of Nanoparticles. J Controlled Release (2018)
3. Decuzzi P. Facilitating the Clinical Integration of Nanomedicines: The Roles 277:35–47. doi: 10.1016/[Link].2018.03.012
of Theoretical and Computational Scientists. ACS nano (2016) 10:8133–8. 13. Wong C, Stylianopoulos T, Cui J, Martin J, Chauhan VP, Jiang W, et al.
doi: 10.1021/acsnano.6b05536 Multistage Nanoparticle Delivery System for Deep Penetration Into Tumor
4. Gabizon A, Catane R, Uziely B, Kaufman B, Safra T, Cohen R, et al. Prolonged Tissue. Proc Natl Acad Sci (2011) 108:2426–31. doi: 10.1073/pnas.1018382108
Circulation Time and Enhanced Accumulation in Malignant Exudates of 14. Baxter LT, Jain RK. Transport of Fluid and Macromolecules in Tumors. I.
Doxorubicin Encapsulated in Polyethylene-Glycol Coated Liposomes. Cancer Role of Interstitial Pressure and Convection. Microvasc Res (1989) 37:77–104.
Res (1994) 54:987–92. doi: 10.1016/0026-2862(89)90074-5
5. Vaage J, Barberá -Guillem E, Abra R, Huang A, Working P. Tissue 15. Soltani M, Chen P. Numerical Modeling of Fluid Flow in Solid Tumors. PloS
Distribution and Therapeutic Effect of Intravenous Free or Encapsulated One (2011) 6:e20344. doi: 10.1371/[Link].0020344
Liposomal Doxorubicin on Human Prostate Carcinoma Xenografts. Cancer 16. Sefidgar M, Soltani M, Raahemifar K, Sadeghi M, Bazmara H, Bazargan M,
(1994) 73:1478–84. doi: 10.1002/1097-0142(19940301)73:5<1478::AID- et al. Numerical Modeling of Drug Delivery in a Dynamic Solid Tumor
CNCR2820730526>[Link];2-1 Microvasculature. Microvasc Res (2015) 99:43–56. doi: 10.1016/
6. Batist G, Ramakrishnan G, Rao CS, Chandrasekharan A, Gutheil J, Guthrie T, [Link].2015.02.007
et al. Reduced Cardiotoxicity and Preserved Antitumor Efficacy of Liposome- 17. Moradi Kashkooli F, Soltani M, Rezaeian M, Taatizadeh E, Hamedi M-H.
Encapsulated Doxorubicin and Cyclophosphamide Compared With Image-Based Spatio-Temporal Model of Drug Delivery in a Heterogeneous
Conventional Doxorubicin and Cyclophosphamide in a Randomized, Vasculature of a Solid Tumor—Computational Approach. Microvasc Res
Multicenter Trial of Metastatic Breast Cancer. J Clin Oncol (2001) 19:1444– (2019) 123:111–24. doi: 10.1016/[Link].2019.01.005
54. doi: 10.1200/JCO.2001.19.5.1444 18. Moradi Kashkooli F, Soltani M, Hamedi MH. Drug Delivery to Solid Tumors
7. Sun D, Zhou S, Gao W. What Went Wrong With Anticancer Nanomedicine With Heterogeneous Microvascular Networks: Novel Insights From Image-
Design and How to Make it Right. ACS Nano (2020) 14(10):12281–90. doi: Based Numerical Modeling. Eur J Pharm Sci (2020) 151:105399. doi: 10.1016/
10.1021/acsnano.9b09713 [Link].2020.105399
8. Jain RK, Stylianopoulos T. Delivering Nanomedicine to Solid Tumors. Nat 19. Baxter LT, Jain RK. Transport of Fluid and Macromolecules in Tumors. (II)
Rev Clin Oncol (2010) 7(11):653–64. doi: 10.1038/nrclinonc.2010.139 Role of Heterogeneous Perfusion and Lymphatics. Microvasc Res (1990)
9. Lane LA. Physics in Nanomedicine: Phenomena Governing the In Vivo 40:246–63. doi: 10.1016/0026-2862(90)90023-K
Performance of Nanoparticles. Appl Phys Rev (2020) 7:011316. doi: 20. El-Kareh AW, Secomb TW. A Mathematical Model for Comparison of Bolus
10.1063/1.5052455 Injection, Continuous Infusion, and Liposomal Delivery of Doxorubicin to

Frontiers in Oncology | [Link] 17 June 2021 | Volume 11 | Article 655781

Moradi Kashkooli et al. Delivery of Drug-Loaded Nanocarriers to Tumors

Tumor Cells. Neoplasia (New York NY) (2000) 2:325. doi: 10.1038/ 39. Soltani M, Tehrani MHH, Moradi Kashkooli F, Rezaeian M. Effects of
[Link].7900096 Magnetic Nanoparticle Diffusion on Microwave Ablation Treatment: A
21. Zhang A, Mi X, Yang G, Xu LX. Numerical Study of Thermally Targeted Numerical Approach. J Magn Magn Mater (2020) 514:167196. doi: 10.1016/
Liposomal Drug Delivery in Tumor. J heat transfer (2009) 131(4):043209. doi: [Link].2020.167196
10.1115/1.3072952 40. Wirthl B, Kremheller J, Schrefler BA, Wall WA. Extension of a Multiphase
22. Hendriks BS, Reynolds JG, Klinz SG, Geretti E, Lee H, Leonard SC, et al. Tumour Growth Model to Study Nanoparticle Delivery to Solid Tumours.
Multiscale Kinetic Modeling of Liposomal Doxorubicin Delivery Quantifies PloS One (2020) 15(2):e0228443. doi: 10.1371/[Link].0228443
the Role of Tumor and Drug-Specific Parameters in Local Delivery to Tumors. 41. Wijeratne PA, Vavourakis V. A Quantitative in Silico Platform for Simulating
CPT Pharmacometrics Syst Pharmacol (2012) 1(11):e15. doi: 10.1038/ Cytotoxic and Nanoparticle Drug Delivery to Solid Tumours. Interface Focus
psp.2012.16 (2019) 9: (3):20180063. doi: 10.1098/rsfs.2018.0063
23. Chauhan VP, Stylianopoulos T, Martin JD, Popović Z, Chen O, Kamoun WS, 42. Dogra P, Butner JD, Ramı́rez JR, Chuang Y-L, Noureddine A, Brinker CJ, et al.
et al. Normalization of Tumour Blood Vessels Improves the Delivery of A Mathematical Model for Nanomedicine Pharmacokinetics and Tumor
Nanomedicines in a Size-Dependent Manner. Nat Nanotechnol (2012) 7:383– Delivery. Comput Struct Biotechnol J (2020) 18:518–31. doi: 10.1016/
8. doi: 10.1038/nnano.2012.45 [Link].2020.02.014
24. Gasselhuber A, Dreher MR, Rattay F, Wood BJ, Haemmerich D. Comparison 43. Shojaee P, Niroomand-Oscuii H, Sefidgar M, Alinezhad L. Effect of
of Conventional Chemotherapy, Stealth Liposomes and Temperature- Nanoparticle Size, Magnetic Intensity, and Tumor Distance on the
Sensitive Liposomes in a Mathematical Model. PloS One (2012) 7(10): Distribution of the Magnetic Nanoparticles in a Heterogeneous Tumor
e47453. doi: 10.1371/[Link].0047453 Microenvironment. J Magn Magn Mater (2020) 498:166089. doi: 10.1016/
25. Zhan W, Xu XY. A Mathematical Model for Thermosensitive Liposomal [Link].2019.166089
Delivery of Doxorubicin to Solid Tumour. J Drug Deliv (2013) 2013:172529. 44. Stillman NR, Kovacevic M, Balaz I, Hauert S. In Silico Modelling of Cancer
doi: 10.1155/2013/172529 Nanomedicine, Across Scales and Transport Barriers. NPJ Comput Mater
26. Stylianopoulos T, Soteriou K, Fukumura D, Jain RK. Cationic Nanoparticles (2020) 6:92. doi: 10.1038/s41524-020-00366-8
Have Superior Transvascular Flux Into Solid Tumors: Insights From a 45. Moradi Kashkooli F, Soltani M, Souri M, Meaney C, Kohandel M. Nexus
Mathematical Model. Ann Biomed Eng (2013) 41:68–77. doi: 10.1007/ Between in Silico and In Vivo Models to Enhance Clinical Translation of
s10439-012-0630-4 Nanomedicine. Nano Today (2021) 36:101057. doi: 10.1016/[Link].
27. Kim MJ, Gillies RJ, Rejniak KA. Current Advances in Mathematical Modeling 2020.101057
of Anti-Cancer Drug Penetration Into Tumor Tissues. Front Oncol (2013) 46. Soltani M, Sefidgar M, Bazmara H, Casey ME, Subramaniam RM, Wahl RL,
3:278. doi: 10.3389/fonc.2013.00278 et al. Spatiotemporal Distribution Modeling of PET Tracer Uptake in Solid
28. Stylianopoulos T, Economides E-A, Baish JW, Fukumura D, Jain RK. Towards Tumors. Ann Nucl Med (2017) 31:109–24. doi: 10.1007/s12149-016-1141-4
Optimal Design of Cancer Nanomedicines: Multi-stage Nanoparticles for the 47. Mascheroni P, Schrefler BA. In Silico Models for Nanomedicine: Recent
Treatment of Solid Tumors. Ann Biomed Eng (2015) 43:2291–300. doi: Developments. Curr med Chem (2018) 25:4192–207. doi: 10.2174/
10.1007/s10439-015-1276-9 0929867324666170417120725
29. Stylianopoulos T, Jain RK. Design Considerations for Nanotherapeutics in 48. Peer D, Karp JM, Hong S, Farokhzad OC, Margalit R, Langer R. Nanocarriers
Oncology. Nanomed: Nanotechnol Biol Med (2015) 11:1893–907. doi: as An Emerging Platform for Cancer Therapy. Nat Nanotechnol (2007) 2
10.1016/[Link].2015.07.015 (12):751–60. doi: 10.1038/nnano.2007.387
30. Chou C-Y, Chang W-I, Horng T-L, Lin W-L. Numerical Modeling of 49. Truskey GA, Yuan F, Katz DF. Transport Phenomena in Biological Systems.
Nanodrug Distribution in Tumors With Heterogeneous Vasculature. PloS Upper Saddle River, New Jersey: Pearson Education, Inc. (2004).
One (2017) 12:e0189802. doi: 10.1371/[Link].0189802 50. Sefidgar M, Soltani M, Raahemifar K, Bazmara H, Nayinian SMM, Bazargan
31. Zhan W, Wang C-H. Convection Enhanced Delivery of Liposome M. Effect of Tumor Shape, Size, and Tissue Transport Properties on Drug
Encapsulated Doxorubicin for Brain Tumour Therapy. J Controlled Release Delivery to Solid Tumors. J Biol Eng (2014) 8:12. doi: 10.1186/1754-1611-8-12
(2018) 285:212–29. doi: 10.1016/[Link].2018.07.006 51. Mpekris F, Angeli S, Pirentis AP, Stylianopoulos T. Stress-Mediated
32. Shamsi M, Sedaghatkish A, Dejam M, Saghafian M, Mohammadi M, Sanati- Progression of Solid Tumors: Effect of Mechanical Stress on Tissue
Nezhad A. Magnetically Assisted Intraperitoneal Drug Delivery for Cancer Oxygenation, Cancer Cell Proliferation, and Drug Delivery. Biomech Model
Chemotherapy. Drug Deliv (2018) 25:846–61. doi: 10.1080/10717544. Mechanobiol (2015) 14:1391–402. doi: 10.1007/s10237-015-0682-0
2018.1455764 52. Eikenberry S. A Tumor Cord Model for Doxorubicin Delivery and Dose
33. Stylianopoulos T, Munn LL, Jain RK. Reengineering the Physical Optimization in Solid Tumors. Theoritical Biol Med Model (2009) 6:16. doi:
Microenvironment of Tumors to Improve Drug Delivery and Efficacy: 10.1186/1742-4682-6-16
From Mathematical Modeling to Bench to Bedside. Trends Cancer (2018) 4 53. Mpekris F, Baish JW, Stylianopoulos T, Jain RK. Role of Vascular
(4):292–319. doi: 10.1016/[Link].2018.02.005 Normalization in Benefit From Metronomic Chemotherapy. Proc Natl Acad
34. Huang Y, Gu B, Liu C, Stebbing J, Gedroyc W, Thanou M, et al. Sci (2017) 114:1994–9. doi: 10.1073/pnas.1700340114
Thermosensitive Liposome-Mediated Drug Delivery in Chemotherapy: 54. Kerr D, Kerr A, Freshney R, Kaye S. Comparative Intracellular Uptake of
Mathematical Modelling for Spatio-Temporal Drug Distribution and Adriamycin and 4’-Deoxydoxorubicin by Non-Small Cell Lung Tumor Cells
Model-Based Optimisation. Pharmaceutics (2019) 11(12):637. doi: 10.3390/ in Culture and Its Relationship to Cell Survival. Biochem Pharmacol (1986) 35
pharmaceutics11120637 (16):2817–23. doi: 10.1016/0006-2952(86)90195-4
35. Rezaeian M, Sedaghatkish A, Soltani M. Numerical Modeling of High- 55. Yorke E, Fuks Z, Norton L, Whitmore W, Ling C. Modeling the Development
Intensity Focused Ultrasound-Mediated Intraperitoneal Delivery of of Metastases From Primary and Locally Recurrent Tumors: Comparison
Thermosensitive Liposomal Doxorubicin for Cancer Chemotherapy. Drug With a Clinical Data Base for Prostatic Cancer. Cancer Res (1993) 53:2987–93.
Delivery (2019) 26:898–917. doi: 10.1080/10717544.2019.1660435 56. Jung KY, Cho SW, Kim YA, Kim D, Oh B-C, Park DJ, et al. Cancers With
36. Shamsi M, Mohammadi A, Manshadi MK, Sanati-Nezhad A. Mathematical Higher Density of Tumor-Associated Macrophages Were Associated With Poor
and Computational Modeling of Nano-Engineered Drug Delivery Systems. Survival Rates. J Pathol Trans Med (2015) 49:318. doi: 10.4132/jptm.2015.06.01
J Controlled Release (2019) 307:150–65. doi: 10.1016/[Link].2019.06.014 57. Moradi Kashkooli F, Soltani M, Rezaeian M, Meaney C, Hamedi MH,
37. He H, Liu C, Liu Y, Liu X, Wu Y, Fan J, et al. Mathematical Modeling of the Kohandel M. Effect of Vascular Normalization on Drug Delivery to
Heterogeneous Distributions of Nanomedicines in Solid Tumors. Eur J Pharm Different Stages of Tumor Progression: In-Silico Analysis. J Drug Delivery
Biopharm (2019) 142:153–64. doi: 10.1016/[Link].2019.06.005 Sci Technol (2020) 60:101989. doi: 10.1016/[Link].2020.101989
38. Dogra P, Butner JD, Chuang Y, Caserta S, Goel S, Brinker CJ, et al. 58. Moradi Kashkooli F, Soltani M, Momeni MM. Computational Modeling of
Mathematical Modeling in Cancer Nanomedicine: A Review. Biomed Drug Delivery to Solid Tumors: A Pilot Study Based on a Real Image. J Drug
Microdevices (2019) 21:40. doi: 10.1007/s10544-019-0380-2 Deliv Sci Technol (2021) 62:102347. doi: 10.1016/[Link].2021.102347

Frontiers in Oncology | [Link] 18 June 2021 | Volume 11 | Article 655781

Moradi Kashkooli et al. Delivery of Drug-Loaded Nanocarriers to Tumors

59. Deen W. Hindered Transport of Large Molecules in Liquid-Filled Pores. Conflict of Interest: The authors declare that the research was conducted in the
AIChE J (1987) 33:1409–25. doi: 10.1002/aic.690330902 absence of any commercial or financial relationships that could be construed as a
60. Roudnicky F, Yoon SY, Poghosyan S, Schwager S, Poyet C, Vella G, et al. potential conflict of interest.
Alternative Transcription of a Shorter, Non-Anti-Angiogenic Thrombospondin-
2 Variant in Cancer-Associated Blood Vessels. Oncogene (2018) 37:2573–85. doi: Copyright © 2021 Moradi Kashkooli, Soltani, Momeni and Rahmim. This is an open-
10.1038/s41388-018-0129-z access article distributed under the terms of the Creative Commons Attribution
61. McDougall SR, Anderson AR, Chaplain MA. Mathematical Modelling of Dynamic License (CC BY). The use, distribution or reproduction in other forums is permitted,
Adaptive Tumour-Induced Angiogenesis: Clinical Implications and Therapeutic provided the original author(s) and the copyright owner(s) are credited and that the
Targeting Strategies. J Theor Biol (2006) 241:564–89. doi: 10.1016/[Link].2005.12.022 original publication in this journal is cited, in accordance with accepted academic
62. Minchinton AI, Tannock IF. Drug Penetration in Solid Tumours. Nat Rev practice. No use, distribution or reproduction is permitted which does not comply with
Cancer (2006) 6:583–92. doi: 10.1038/nrc1893 these terms.

Frontiers in Oncology | [Link] 19 June 2021 | Volume 11 | Article 655781