Article Information

• PDF (875 kb)
• Full Text with Large Figures
• Supporting Material
• Cited by in Scopus (4)

PubMed

• Articles by Anh-Tuan Dinh
• Articles by Samir Mitragotri

Related Articles

• Affected-Sib-Pair Data Can Be Used to Distinguish Two-Locus ...

  Affected-Sib-Pair Data Can Be Used to Distinguish Two-Locus Heterogeneity from Two-Locus Epistasis The American Journal of Human Genetics, Volume 73, Issue 6, 1 December 2003, Pages 1468-1470 Heather J. Cordell Full Text | PDF (54 kb)

• Unbiased Rotational Moves for Rigid-Body Dynamics

  Unbiased Rotational Moves for Rigid-Body Dynamics Biophysical Journal, Volume 85, Issue 5, 1 November 2003, Pages 2973-2976 Daniel A. Beard and Tamar Schlick Abstract We introduce an unbiased protocol for performing rotational moves in rigid-body dynamics simulations. This approach—based on the analytic solution for the rotational equations of motion for an orthogonal coordinate system at constant angular velocity—removes deficiencies that have been largely ignored in Brownian dynamics simulations, namely errors for finite rotations that result from applying the noncommuting rotational matrices in an arbitrary order. Our algorithm should thus replace standard approaches to rotate local coordinate frames in Langevin and Brownian dynamics simulations. Abstract | Full Text | PDF (100 kb)

• Evaluating Intramural Virtual Electrodes in the Myocardial W...

  Evaluating Intramural Virtual Electrodes in the Myocardial Wedge Preparation: Simulations of Experimental Conditions Biophysical Journal, Volume 94, Issue 5, 1 March 2008, Pages 1904-1915 G. Plank, A. Prassl, E. Hofer and N.A. Trayanova Abstract While defibrillation is the only means for prevention of sudden cardiac death, key aspects of the process, such as the intramural virtual electrodes (VEs), remain controversial. Experimental studies had attempted to assess intramural VEs by using wedge preparations and recording activity from the cut surface; however, applicability of this approach remains unclear. These studies found, surprisingly, that for strong shocks, the entire cut surface was negatively polarized, regardless of boundary conditions. The goal of this study is to examine, by means of bidomain simulations, whether VEs on the cut surface represent a good approximation to VEs in depth of the intact wall. Furthermore, we aim to explore mechanisms that could give rise to negative polarization on the cut surface. A model of wedge preparation was used, in which fiber orientation could be changed, and where the cut surface was subjected to permeable and impermeable boundary conditions. Small-scale mechanisms for polarization were also considered. To determine whether any distortions in the recorded VEs arise from averaging during optical mapping, a model of fluorescent recording was employed. The results indicate that, when an applied field is spatially uniform and impermeable boundary conditions are enforced, regardless of the fiber orientation VEs on the cut surface faithfully represent those intramurally, provided tissue properties are not altered by dissection. Results also demonstrate that VEs are sensitive to the conductive layer thickness above the cut surface. Finally, averaging during fluorescent recordings results in large negative VEs on the cut surface, but these do not arise from small-scale heterogeneities. Abstract | Full Text | PDF (1431 kb)

• …more

Copyright © 2007 The Biophysical Society. All rights reserved.
Biophysical Journal, Volume 92, Issue 3, 831-846, 1 February 2007

doi:10.1529/biophysj.106.095521

Biophysical Theory and Modeling

Understanding Intracellular Transport Processes Pertinent to Synthetic Gene Delivery via Stochastic Simulations and Sensitivity Analyses

Anh-Tuan Dinh¹, Chinmay Pangarkar¹, Theo Theofanous and Samir Mitragotri^,  

Department of Chemical Engineering, University of California, Santa Barbara, California

Address reprint requests to Prof. Samir Mitragotri, Dept. of Chemical Engineering, University of California, Santa Barbara, CA 93106. Tel.: 805-893-7532.

¹ Anh-Tuan Dinh and Chinmay Pangarkar contributed equally to this work.

Polymer-based systems are being extensively studied as carriers for gene therapy 1, and have also been used in clinical trials 2. Despite the obvious advantages of safety and malleability over viral vectors, polymer-based vectors (polyplexes) have not been very successful at the clinical level owing to their poor delivery efficiency. Even the best polyplexes are 1000-fold less effective than typical viral vectors such as adenoviruses. This poor efficiency stems from the numerous intracellular barriers that retain or destroy a majority of the gene dose before it can reach the host nucleus 3. A central challenge in the field is identification and synthesis of a polymer that can enhance translocation of the polyplex across these barriers. To this end, a large number (literally thousands) of polymers have been evaluated—poly-imines 4, dendrimers 4, polyamino esters 5, chitosans 6, and cyclodextrins 7, to name a few. However, only a few of these have been found to be significantly better than polyethylenimine (PEI^25kDa), the accepted standard in polymer-based gene delivery. Further, it is not clear whether these polymers, many of which have been developed and optimized for use with cultured cells, will be effective in clinical applications.

The low transfection efficiency of synthetic vectors leads to a key question: Is there an inherent limitation to polyplexes? Or have we not found the magic polymer yet? Two challenges need to be addressed before this question can be answered. The first challenge, biological in nature, involves developing the design criteria of polyplexes that will lead to maximum efficiency. The second challenge, chemical in nature, involves design and synthesis of polymers that will form such polyplexes.

In this article, we exclusively address the first challenge. There are two major issues that make this task difficult. The first issue stems from fragmented understanding of intracellular trafficking of polyplexes, mainly due to the strong dichotomy in the experimental approaches used. One approach is based on macroscopic measurements of the delivery efficiency using in vitro transfection assays. The other approach is based on measurements at the single-cell single-particle level. Macroscopic experiments do not provide in-depth mechanistic understanding of the events leading to the final effect. Single particle data, on the other hand, aim to quantify individual transport steps 8,9, many of which, such as escape from endosomes, remain tremendously difficult to visualize due to the rare and random nature of these events 8. The link between these two distinct scales (i.e., macroscopic and microscopic) is still missing. The second outstanding issue is the limited knowledge of the structure-activity relationships that link vector’s physico-chemistry to the efficiency of individual trafficking steps. This limitation confounds systematic experimental optimization of vector properties.

We address the problem of understanding and optimizing intracellular transport by developing a detailed mathematical model of the gene delivery pathway. Previous models of gene delivery have essentially followed a pharmacokinetic approach 10,11,12,13, which treats the cell as a well-mixed compartment. The major shortcoming of kinetic models is that they approximate all spatial and transport processes by kinetic equations. The kinetic rate constants have to be obtained by fitting the model to experimental data. This makes it difficult to extrapolate the results of the model predictions. Most importantly, kinetic models, due to their simple approximations of intracellular transport, fail to take advantage of the wealth of data made available by single-particle tracking experiments.

We realize that much of the complexity of the vector-cell system arises from the spatio-temporal variation in the rates of various intracellular processes and can be captured only by considering a spatial view of the cell. To address this challenge, we develop an advanced stochastic simulation framework that effectively describes the dynamic and spatial nature of intracellular trafficking of nanoscale carriers. By providing a realistic representation of topology and dynamic organization of the cell interior, the modeling framework serves as a three-dimensional, “live,” computer-constructed cellular map for navigation of nanoscale carriers. The model contains information from the shortest to the longest relevant length and timescales, and wherever possible, is rigorously validated against experimental data. Wherever mechanistic or quantitative information in the literature is incomplete, we perform systematic experiments to measure the required processes. We use the model to calculate the delivery efficiency—the probability that a polyplex, upon internalization, successfully delivers DNA to the host nucleus, and to predict the optimal intracellular itinerary. We show that the optimal pathway is controlled by several cell-specific properties such as cell morphology, endocytic trafficking, and microtubule-dependent transport.

Polyethylenimine (25kDa) was obtained from Sigma (St. Louis, MI). All fluorescent dyes (Oregon green, TMR-Dextran, 70kDa) were obtained from Molecular Probes (Eugene, OR). Amine groups on polyethylenimine (25kDa) were labeled using the succinimidyl ester of Oregon green as per the manufacturer’s protocol. The complexes were always prepared at a Nitrogen/Phosphorus ratio of 9, in 150mM NaCl solution. Complexes were used between 30 and 45min post-synthesis. They were verified to transfect skin fibroblasts.

Human Skin Normal fibroblasts (ATCC, TE.353.Sk, Manassas, VA) were cultured as per the standard ATCC protocols. Briefly, cells were maintained in Dulbecco’s modified eagle’s medium (ATCC), which had been supplemented with 4mM L-glutamine and 1.5g/L sodium bicarbonate. Cells were grown at 37°C and 5% CO₂; and maintained in the log phase by consistent subculture. On the day before the experiment, cells were seeded onto a glass bottom petri dish at a density of ∼7000–10,000cells/cm². The next day, complexes were added to serum-free medium at a concentration of 2μg DNA/ml. The cells were incubated with this medium for 15–30min at 37°C. Post-incubation, the unbound particles were washed away 5× by phosphate-buffered saline and the cells were replenished with fresh serum containing medium, which was supplemented with 25mM HEPES to maintain carbonate-bicarbonate balance in the absence of CO₂.

The cells were placed on the stage of a Zeiss (Thornwood, NY) Axiovert-25 inverted microscope, which was fitted with a Bioptechs Delta T temperature controller (Bioptechs, Butler, PA). The cells were observed using an oil immersion 100× objective. Oregon green was excited and observed using HQ480/40× and HQ 555/40M FITC compatible filters from Chroma (Rockingham, VT). For every cell, bright-field and fluorescence images were acquired using a cooled CCD-camera (CoolSNAP_HQ, Roper Scientific, Duluth, GA) which was controlled by Metamorph imaging software (Universal Imaging, Downington, PA). The images were acquired for 1–2min at 1 frame per second.

Lysosomes in the cells were immunolabeled using standard techniques. Briefly, the cells were fixed for 30min using freshly prepared 2% paraformaldehyde. This was followed by permeabilization with 0.01% Tween-20 in phosphate-buffered saline and blocking with 2% BSA for 2h. The cells were then labeled with mouse monoclonal anti-Lamp1 (primary antibody, BD Biosciences) and subsequently with Rhodamine red-X-Goat anti-mouse IgG (secondary antibody, R6393, Molecular Probes). The cells were then mounted using VectaShield (Vector Labs, Burlingame, CA). Rhodamine-X labeled antibodies were observed using HQ 545/30X and HQ 640/50M filters from Chroma.

All image acquisition and processing steps were performed within the Metamorph software. Phase contrast images of cells were used to extract coordinates of the cell membrane and nuclear membrane. The raw fluorescence images show bright punctuate structures corresponding to diffraction-limited images of fluorescent vectors. Images were processed with a custom FFT-based protocol to filter out very long-range and very short-range correlations. This was followed by a series of Gaussian blurring, convolution, and unsharp masking steps until the resultant image could be thresholded such that most of the bright structures could be individually measured. All particle coordinates were exported to MS Excel (Microsoft, Redmond, WA). All further calculations involving cell shape, vector distribution, etc., were performed in MatLab (The MathWorks, Natick, MA). For single particle tracking, stacks of fluorescent images of cells were filtered to remove high frequency noise. Other filtering steps were avoided to preserve the intensity structure of the raw image. Trajectories of all the particles were individually verified by replaying the video with the trajectory superimposed on the image. The coordinates of a vector in all the frames were exported to MS Excel, and were used for all further calculations. Colocalization of Oregon green labeled PEI25-DNA vectors and anti-Lamp-1 labeled lysosomes was measured by overlaying the Red and Green channel images. The images were self-normalized using standard histogram-stretching methods. Colocalization was considered significant if the Red and Green structures had >85% of their combined pixels in common.

The objective of stochastic simulations is to bring-out the spatial and physical aspects of intracellular trafficking of polyplexes by providing a realistic representation of cell geometry and intracellular organization and a discrete and mechanistic description of transport processes. As in single-particle tracking experiments, the proposed model follows the trajectories of polyplexes inside cells. The computational domain is defined by the cell geometry, which is reconstructed from experimental visualization of cells (see below). We distinguish between in vitro and in vivo experiments. Many cells, especially those considered in this analysis (human skin fibroblasts) spread on the substrate under in vitro conditions, and the cell thickness is often much smaller than other dimensions due to adhesion to the surface. Therefore, it is sufficient to formulate the simulations in two dimensions for in vitro applications (Figure 1a). In vivo applications, however, require considerations of three-dimensional geometry (Figure 1bd).

Display large version of this figure

Figure 1

Computer representation of human skin fibroblasts (a) in vitro, where the cell is represented as a two-dimensional object and (b) in vivo, where the cell is represented by a three-dimensional structure. Green and red objects represent endosomes and lysosomes respectively. Thin yellow lines represent microtubules (MTs). The microtubule-organizing center (MTOC) is randomly located in the vicinity of the nucleus. For clarity, we use fewer endosomes, lysosomes, and MTs to render the reconstructed cells than in simulations. Two regions of the three-dimensional cell are magnified and shown in panels c and d.

The positions and the biophysical states of polyplexes are updated using a stochastic algorithm. Fundamentally, all stochastic events simulated here can be classified into two categories: 1), transport events; and 2), reaction events. Transport events (diffusion or motor-driven movements) are captured by equations of motion. Reaction events, involving association, dissociation or rupturing of entities, are approximated as first-order Markov processes. Definitions of variables and symbols are summarized in Table 1.

Table 1 Nomenclature

	Acell	Area of a cell.
	cE,cL	Concentration of endosomes and lysosomes (# per unit area).
	c(r,t)	Normalized local concentration of polyplexes in the cell at distance r from the nuclear boundary at time t post-transfection.
	cperi	Normalized local concentration of polyplexes in the perinuclear region, which is defined as a layer of width 2μm around the nucleus.
	csupra	Normalized local concentration of polyplexes in the supranuclear region.
	DDNA	Diffusivity of DNA in cytoplasm, μm2/s.
	DE	Diffusivity of endosomes, μm2/s.
	DL,f, DL,c	Diffusivities of free and clustered lysosomes, μm2/s.
	Φ(t)	Spatially averaged success probability of a vector that escapes at any location in the cytoplasm at time t post-transfection.
	fcat	Catastrophe frequency of microtubules, s−1.
		Rate constant for a particle binding to and detachment from MTs, s−1.
	kcluster, kde-cluster	Rates of formation and breakup of lysosomal clusters, s−1.
	kescape	Rate constant for escape of vectors from endosomes or lysosomes, min−1.
	ktr	Rate constant for transfer of vectors from endosomes to lysosomes, min−1.
	kunpack	Rate constant for unpacking of vectors in cytoplasm or lysosomes, min−1.
	kint	Rate constant for internalization of vectors via endocytosis, min−1.
	kdeg-cyto,DNA	Rate constant for degradation of DNA in the cytoplasm, min−1.
	knuc	Rate constant for nuclear entry, min−1.
	kdegL	Rate constant for degradation of vector in lysosomes, min−1.
	Ks	Equilibrium constant for particle binding and unbinding with microtubules (s=±1)
	Φ(r)	Probability that a vector, which escapes from endo/lysosomes at distance r from the nuclear membrane, can successfully deliver DNA into the nucleus.
	Θ(t)	Probability that a vector, which escapes from endo/lysosomes at time t post-transfection, can successfully deliver DNA into the nucleus.
	Pe(t)	Probability distribution of escape time.
	r	Distance of a particle from the nuclear boundary, μm.
	Ψ	Overall delivery efficiency.
	s	Transport substate (0: diffusion, and ±1: toward MT plus and minus ends).
	SV	Biophysical state of a vector (M: cell membrane, E: inside endosomes, L: inside lysosomes, C: inside cytoplasm, and P: unpacked plasmids).
	t	Time.
	Tf	Time at which the experiment/simulation is terminated, in hours.
	tgrow, tshrink	Mean growth and shrinkage times, s.
	vgrow, vshrink	Rates of MT growth and shrinkage, m/s.
		MT-dependent velocity of a particle toward plus or minus ends of MTs, μm/s.

For in vitro applications, we reconstruct cell geometries from phase contrast images of cultured human fibroblasts. Viewing from the top, the cell interior is divided into two regions, cytoplasmic and supranuclear, separated by the nuclear boundary (see Figure 1a). We use an image analysis program (Metamorph, Universal Imaging) to extract the coordinates of the nuclear boundary and the cell boundary. The high-resolution images captured at 100× magnification allow for precise determination of the shape of cell and nucleus. For in vivo applications, we reconstruct the cells from histological and electron microscopy images of fibroblasts embedded in tissues. The procedure for in vivo cells is more complicated and described in detail in the Supplementary Material . The cell geometry is represented by two surfaces describing the cell membrane and the nuclear envelope (see Figure 1b). Our reconstructed three-dimensional cells are in good agreement with past descriptions of fibroblasts under in vivo conditions. All cell geometries remain unchanged during simulations.

Microtubules (MTs) exhibit dynamic instability, i.e., stochastic switching between prolonged phases of assembly and rapid phases of disassembly. The frequency of transition from slow assembly to rapid disassembly (referred to as “catastrophe”) determines the stability of the MT population 24. In this model, MTs are modeled as straight lines, which grow radially in random directions from the microtubule-organizing center (MTOC) toward the cell periphery. MTOC is placed randomly in the vicinity of the nucleus. To describe microtubule polymerization dynamics, we use stochastic equations that follow the individual MT length history. A microtubule switches intermittently between three phases (growth, shrinkage, and collapse) as follows 14:

(1)

where q and m are random numbers between 0 and 1. The values k_grow, k_shrink, and k_cat are the reaction rate constants for MT growth, shrinkage, and collapse, respectively.

We introduce two state variables to characterize intracellular trafficking of polyplexes, namely, biophysical state S, which indicates the biological compartmentalization and transport state s, which represents the type of movements of polyplexes.

After binding to the cell membrane, the intracellular pathway of a polyplexes can be represented by five distinct biophysical states S: bound to cell membrane (M); inside endosomes (E); inside lysosomes (L); inside cytoplasm (C); and unpacked plasmids (P) (Fig. 2). The states correspond to the stages of the gene delivery pathway (see Fig. 2 legend and Supplementary Material , section 1.3).

Display large version of this figure

Figure 2

Gene delivery pathway of polyplexes. Polyplexes are internalized by endocytosis. They are immediately delivered to endosomes and subsequently to lysosomes. To gain entry to the nucleus, the polyplexes must escape from endosomes or lysosomes. Once inside the cytoplasm, the polyplexes unpack to release the DNA for successful nuclear entry. Ultimately, the exogenous DNA is transcribed, and gene products are synthesized. The biophysical states S and the rates of transition between them are shown in the figure.

Transitions from one state to another depend on various biological and physical factors, such as location, interactions with cellular organelles, and physiochemical properties of polyplexes. In the present model, we approximate the transition between two biophysical states of the vectors as a first-order reaction. For instance, we use k_escape (i.e., the rate of escape from endocytic vesicles) to characterize the transition from S=E to S=C (Fig. 2) The probability that a transition from state E to state C occurs in the interval [t, t+Δt] is given by

(2)

Each biophysical state has a distinct transport pattern. Membrane-bound and cytoplasmic polyplexes move only via diffusion. On the other hand, polyplexes inside endosomes switch intermittently between diffusion and directional transport. To account for these distinct transport modes, we introduce the concept of transport state s. Each biophysical state is characterized by a set of transport states and a diagram depicting the transitions between these states. At each point in time, a polyplex occupies one of the movement states. For instance, transport of polyplexes inside endosomes is represented by three transport states, s=±1, −1, and 0 (+1, plus-end directional transport; −1, minus-end directional transport; and 0, diffusion). We define and as velocities of directional movements toward the plus-end and the minus-end of microtubules. The values and are first-order rate constants characterizing particle’s transition from s=0 to s=+1 and s=−1, respectively. The values and are the rate constants of corresponding reverse processes. Thus, the cytoplasmic transport pattern of endosomes is represented by seven parameters, (free diffusion coefficient), , , , , and . These transport coefficients are independent of location inside the cytoplasmic region. Quantitative measurements of movements of PEI-DNA complexes in human fibroblast support this assumption. Table 2 summarizes the transition diagram and the associated transport parameters for all five biophysical states.

Table 2 Transition diagrams and transport parameters for biophysical states

	Biophysical state	Transition map	Parameters	Reference
	Membrane		D0,M∼0.001μm2/s	This work.
	Endosome		D0,E∼0.0005–0.001μm2/s	This work.
	Lysosome		D0,L∼0.0001μm2/s	This work.
			kde-cluster=0.125s−1
			kcluster=0.5s−1
	Cytoplasmic		D0,C∼0.001μm2/s	Estimated based on Stokes-Einstein relationship.
	Plasmid		D0,P∼0.001μm2/s	17

The position of a polyplex occupying a diffusive transport state is updated using two- or three-dimensional random walk. For MT-bound entities, the equations of motion are applied,

(3)

where is the position of polyplex i and characterizes the microtubule j that the polyplex i is associated with.

To simulate transitions between biophysical states, we employ the same principle introduced above (see Eq. (2)). However, the occurrence of several transition events requires specific physical conditions. For instance, binding of a diffusing endosome to MT only occurs if the endosome is within a radius r_e from a neighboring MT, and the probability that it will bind to and travel toward the plus- and minus-ends of the MT in an interval Δt are and , respectively. Similarly, termination of microtubule-dependent runs depend not only on the rates of unbinding, and , but also on the physical layout of the microtubule tracks. The endosome will fall off the MT, that is, become diffusive in the next time step, once reaching the ends of the MT or hitting the nuclear boundary or the cell boundary, which act like solid walls.

The most interesting case is the transition between the clustered state and the free state of polyplexes inside lysosomes. Lysosomes tend to form clusters in the perinuclear region (Figure 4a, inset). When a lysosome-carrying polyplex aggregates with other lysosomes, it exhibits restricted diffusion. When it is free, it can travel bidirectionally along microtubules just like endosomes. The declustering rate, , is presumably space-dependent, while the rate of transition from the free to the clustered state is a function of the normalized clustering rate and the local concentration of lysosomes, C_L(r), which is obtained a priori by analysis of lysosome images after immunolabeling.

Here, we describe the simulation scheme for in vitro applications. The simulation scheme for in vivo applications is similar. For each simulation, an arbitrary cell configuration is chosen from the library of different cell geometries, which defines the computational domain. During one time-step, the cellular environments and the states/positions of polyplexes are updated according to the following scheme:

To guarantee accuracy and stability of the numerical scheme, we use a small time step, Δt=0.2s. The time step is much smaller than the governing timescales of the transport processes. Varying Δt between 0.1s and 2s does not significantly influence the model predictions. For each parameter set, we sample 10,000 trajectories from 200 cells. The delivery efficiency is calculated based on the number of DNA molecules that gain entry to the nucleus at a final time T_f, when we terminate the simulations. The simulation software is available by request.

The model contains a large number of parameters, which describe different aspects of the problem such as cell geometry, microtubules, trafficking of endocytic vesicles, trafficking of vectors in the endocytic pathway, and trafficking of vectors after endosomal/lysosomal escape. We measured parameters describing the trafficking of vectors present within endocytic vesicles using skin-fibroblasts as the model cell-line and PEI^25kDa-DNA as the vector. However, for the steps after endosomal escape, experimental measurements are not possible at the single particle level. In such cases, data from bulk measurements was used. For example, we adapted the rate of plasmid degradation in the cytoplasm from Lechardeur et al. 15, who measured degradation of injected plasmids using fluorescence in situ hybridization.

Based on immunolabeled images of fibroblasts, we estimate that the number of microtubules in the plane of movements under in vitro conditions to be ∼1000–1500. To our knowledge, microtubule polymerization dynamics in human skin fibroblasts have not been quantified yet. Accordingly, we postulate that the parameters that characterize microtubule dynamic instability in human fibroblasts are similar to those of CHO fibroblasts. Komarova et al. 16 investigated the lifecycle of MTs in CHO cells and reported a growth rate (v_grow) of 0.25–0.35μm/s, a shortening rate (v_shrink) of 0.3–0.4μm/s, and a catastrophe frequency (f_cat) of 0.005s⁻¹. Also, the mean growth time t_grow is ∼40–60s and the mean shrinkage time t_shrink is 5–15s. Variations of v_grow, v_shrink, t_grow, and t_shrink within the estimated intervals do not significantly affect the predicted particle distributions. The most sensitive parameter is f_cat, which determines the stability of the microtubule tracks.

The cytoplasmic transport properties of endosomes are represented by seven parameters, D_0,E (free diffusion coefficient), , , , , , and . These parameters are not available in the literature and are therefore measured from our own experiments (see Supplementary Material , Section 1.2.1).

1. Diffusivity, D_0,E. Trajectories of ∼50 fluorescently labeled, almost immobile vectors were selected from experimentally measured movements of PEI^25kDa-DNA vectors (assuming that vectors are present within endosomes). These trajectories were used to calculate the mean-squared displacement (MSD) as a function of the timescale t. Based on the MSD, we estimate D_0,E to be 0.0005–0.001μm²/s. Similar values have been reported previously in the literature. Immunolabeling of endosomes (using anti-rab5, red) showed that all observable green fluorescence of PEI-DNA vectors colocalizes with anti-rab5 at the time of measurements (1h post-transfection). This confirms that, within experimental errors, the imaged vectors were indeed present within endosomes.
2. Rate of unbinding, . The rate of unbinding is the inverse of the mean duration of “spurts” or “runs” in a certain direction. Figure 3a shows the probability distribution of duration of a run in the plus and minus directions, calculated from ∼250 active trajectories. The exponential nature of this distribution also supports the assumption of a first-order unbinding process. The run lengths in either direction show identical distributions and thus yield .

Display large version of this figure

Figure 3

Transport properties of endosomes. (a) Distribution of the duration of runs of endosomes in the plus (solid circles) and minus (open circles) directions and (b) distribution of the frame-to-frame velocities in plus (solid circles) and minus (open circles) directions. The vertical coordinate indicates the frequency of occurrence of an event, i.e., the probability of observing a certain run-length (a) or a certain velocity (b).

3. Rate of binding, . At equilibrium, a certain fraction f of particles must exist in bound states 1 and −1, and . Here and ; and for , we get . Thus, knowing , and measuring f from time-lapse videos of cells at different times, the rate of binding can be calculated. For PEI-DNA containing endosomes, the rate of binding calculated from videos over 20–120min post-transfection is ∼1/8s⁻¹.
4. Velocity distribution, . Endosomes do not exhibit a constant velocity, but a distribution of velocities that ranges from <0.05μm/s to 3μm/s. Figure 3b shows the distribution of frame-to-frame speeds in both the plus and minus directions. It can be seen that the distributions are identical exponential distributions, with a mean velocity of 0.33μm/s.

Lysosomes form clusters and accumulate in the perinuclear region, as illustrated in Figure 4ab. Lysosomes switch stochastically between free and clustered states, and the properties in each of these states are measured as follows:

1. Transport of free lysosomes. We use the same method as in case of endosomes to obtain , , , and D_f,L for free lysosomes. Lysosomes are labeled by fluorescent dextran (10kDa), which has been chased into the lysosomes by incubating the cells for 24h. We confirmed this experimentally by labeling endosomes with lysosome-associated-membrane-protein (lamp)-1. The analysis of the experimental data on lysosome transport yields , , , and D_f,L=10⁻⁴μm²/s.
2. Transport of clustered lysosomes. Clustered perinuclear lysosomes rarely exhibit directed motion. In the event that directed motion is seen, it is always followed by dissociation of the lysosomes from the cluster into the free phase. Thus, transport of a lysosome in the clustered phase is described by diffusion (D_c,L=10⁻⁴μm²/s).
3. Transition between free and clustered states. To find the rate at which a free lysosome becomes part of a cluster and vice versa, we examined videos of lysosomal clusters. In every video, a cluster of 4–5 lysosomes was randomly chosen and tracked for a period of 2min. The cluster is often highly dynamic, with lysosomes constantly leaving the cluster. At the same time, other lysosomes come to join the cluster. Very often, the departed lysosomes rejoin the original cluster. Averaging over ∼75 clusters, we obtain the rate of clustering, i.e., the rate of forming physical bonds between two lysosomes, to be k_cluster=0.5s⁻¹. Based on the global state of the cell (the number of free and the number of clustered lysosomes) and assuming equilibrium, we can find the equilibrium constant for the clustering-declustering processes, Ω=k_cluster/k_de-cluster. Here Ω=4, and hence the rate of declustering is k_de-cluster=0.125s⁻¹.

Display large version of this figure

Figure 4

Intracellular spatial distribution of lysosomes and transfer of PEI-DNA from endosomes to lysosomes. (a) Fluorescent micrograph of lysosomes labeled by incubation of human dermal fibroblasts with fluorescent dextran for 2h followed by a 22-h chase in culture medium. The nuclear and cell membrane are traced for clarity. Inset shows typical clusters of perinuclear lysosomes. (b) Spatial distribution of lysosomes represented by p(r), the probability of finding a lysosome at a distance r from the nucleus. (c) Gradual delivery of PEI²⁵-DNA vectors from endosomes to lysosomes. The ordinate is the fraction of total intracellular PEI²⁵-DNA vectors that colocalize with LAMP1 positive structures (lysosomes) at different times post-transfection. Lysosomal delivery is seen to be substantially complete at 24h post-transfection. The shaded line is a first-order kinetic fit on the data yielding k_transfer=1/500min⁻¹.

When polyplexes are inside endosomes and lysosomes, they display transport patterns according to the parameters reported above. To determine the rate of transfer from endosomes to lysosomes, k_tr, we examine the colocalization of Oregon green labeled PEI-DNA (in endosomes) with rhodamine-immunolabeled lysosomes as a function of time, post-transfection. Quantitative measurements allow us to capture the increasing fraction of PEI-DNA present in lysosomes (Figure 4c). The data are averaged over each cell, and over multiple cells at each time point. The curve is modeled with a first-order reaction to yield k_tr=1/500min⁻¹. This is the global transfer rate. The local transfer rate will depend on the local concentration of lysosomes, c_L.

The parameters which describe the processes after endosomal escape are adapted from the literature, including the rate of DNA degradation 15, the diffusivity of vector in the cytoplasm 17, the rate of unpacking 18, and the rate of binding of the vector to the nucleus 19 (see Supplementary Material , Table S4).

The stochastic model spans several length and timescales, as discussed above. The model is therefore validated at three separate levels against: 1), exact solution of single vesicle transport on single filament at timescales of 100s; 2), experimental measurements of mean-squared displacements of endosomes in a live cell at timescales of 100s; and 3), whole-cell-level spatial distribution of PEI^25kDa-DNA vectors over 24h. Comparisons of model predictions to gene expression data are given in Supplementary Material , section 2.5.

To verify that the stochastic algorithm was implemented properly, we simulate movements of endosomes along a single, fixed filament. Transport of endosome obeys the transition map depicted earlier (see Table 2). The endosome is initially located at x=0. The spatial movements of the endosomes in time are characterized by p(x,t), the probability that an endosome is located at position x at time t. This probability can be theoretically predicted by a system of one-dimensional diffusion-advection-reaction equations 20. Readers are referred to Supplementary Material in section 1.7 for more information on the diffusion-advection-reaction equations. In Figure 5a, we compare p(x,t) predicted by stochastic simulation (sampling of positions of 5000 particles) to the exact solution provided by diffusion-reaction-advection equations, for t=50s and 100s. The good agreement between simulation and analytical solution confirms that the stochastic algorithm is capable of accurately predicting transport mediated by molecular motors.

Display large version of this figure

Figure 5

Model validation. (a) Displacement probability p(x,t) of a particle engaging in motor-assisted transport along a filament at t=25s (solid line, exact solution; open circles, simulation) and t=100s (dashed line, exact solution; solid circles, simulation). (b) Average mean-squared displacement for stochastically generated trajectories (solid line) and the experimentally measured trajectories (triangle). A crossover from subballistic to diffusive regime occurs between t=10–100s.

We simulate movements of endosomes in two-dimensional cultured fibroblasts. The cell geometry is accurately reproduced from phase contrast images, and the MTOC is assigned a random location near the nucleus. The microtubules grow radially from the MTOC and exhibit dynamic instability, thus mimicking the real cells. A comparison between the predicted MSD and the measured MSD is shown in Figure 5b. At short times (t<1s), simulations show that vesicle motion is ballistic, as also seen in experiments (t^γ,γ≈2). However, at long times (t>10s), particle motion exhibits diffusive scaling (t¹) and may be described by an effective diffusive coefficient, D^eff. Kulkarni et al. reported similar observations 21. Over short timescales, vesicles experience limited runs and the mean velocity is nonzero and directional. However, at long times, vesicle motion is averaged over several random (in magnitude and direction) ballistic runs and pauses, and exhibits a diffusionlike behavior. The simulation captures the essential physics of microtubule-dependent transport.

We observe that D^eff depends strongly on the density and the stability of microtubules. Decrease in the number of microtubules, N_MT, and increase in the frequency of catastrophe, f_cat, substantially reduce D^eff and cause accumulation of endosomes in the periphery. The effects of microtubule organization on gene delivery have been discussed in detail elsewhere 22. It should be noted that the effects of microtubule density and stability on the overall delivery efficiency are not significant for the ranges of parameters under consideration.

The model was applied to predict the spatio-temporal distribution of PEI^25kDa-DNA complexes that are present within endocytic vesicles. Prediction of the spatial distribution of vectors in dermal fibroblasts at 30min, 4h, and 11h post-transfection (Fig. 6, lower panel) shows excellent qualitative agreement with experimental images (Fig. 6, upper panel). Looking from the top, the cell cytoplasm was divided into two regions—1), supranuclear; and 2), cytoplasmic—which are separated by the nuclear boundary. The supranuclear region is a thin layer of 2–4μm in width directly above the nucleus. Polyplex distributions in these two regions were treated separately.

Display large version of this figure

Figure 6

Comparison of model predictions (lower panel) of spatiotemporal distribution of PEI-DNA complexes to experimental data (upper panel). PEI-DNA complexes are labeled with Oregon green (green) and lysosomes are immunolabeled with rhodamine-X labeled secondary antibody (red). The yellow regions represent PEI-DNA vectors that colocalize with lysosomes. The same cells are reconstructed using simulations and presented in lower panel. The three cells correspond to (a, b) 30min, (c, d) 4h, and (e, f) 11h post-transfection.

We also performed a quantitative comparison of model predictions with experimental data by calculating normalized polyplex concentration c(r,t) (the number of polyplexes per unit area at a distance r from the nuclear membrane at time t post-transfection). The predicted polyplex concentration is within ±10% of the experimental data (Fig. 7). This is encouraging since we do not use any fitting and all parameters are either measured or adapted from literature.

Display large version of this figure

Figure 7

Validation of the model against experimental measurement of spatio-temporal distribution of PEI^25kDa-vectors. (a) Normalized local concentration c(r,t) of vectors measured experimentally (symbols) and predicted by the model (lines) at 10min (circles), 2h (squares), and 24h (triangles). The abscissa is r, the distance from the nucleus. The distributions are sampled from 20 cells. (b) The temporal evolution of vector concentration in the supranuclear region (squares) and the perinuclear region (circles). The perinuclear region is defined as a region of width 2μm around the nucleus. The lines are the predicted distributions.

Initially, the polyplexes are dispersed uniformly in the cell due to random sites of endocytosis on the cell membrane (c(r,t=0)=1, Figure 7a) and accumulate near the nucleus with time. Perinuclear accumulation of polyplexes occurs in two distinct phases (Figure 7b, c_peri(t)). During the first phase, the average concentration of complexes in the perinuclear region quickly doubles within 1–2h. This rapid accumulation is due to facilitated diffusion of endosomes containing PEI-DNA complexes on MTs as discussed above. Consequently, the complexes gradually disperse along the MTs. Since MTs are organized in an asterlike fashion, a uniform distribution on MTs leads to a higher concentration of polyplexes in the perinuclear region 23. The timescale over which the first phase of accumulation occurs depends upon the size of the cell and the MT-based effective diffusivity. The second phase of perinuclear accumulation is much slower, with a timescale of 6–10h, and arises from gradual transfer of complexes from endosomes to lysosomes (Figure 4c). Delivery to lysosomes immobilizes (and hence concentrates) the polyplexes in perinuclear clusters (Figure 7b, c_peri(t)). Thus, perinuclear accumulation of complexes arises not from any active targeting, but from the transport processes inherent to endocytic trafficking. The model, based on fundamental interactions between vesicles and microtubules, is able to predict the spatio-temporal distribution of polyplexes, which evolve at much longer timescales.

In the supranuclear region, the density of MTs is low and most polyplexes do not exhibit MT-based transport. Supranuclear endosomal polyplexes diffuse slowly toward the perinuclear region, are transferred to lysosomes, and become immobilized in the perinuclear clusters. Effectively, polyplexes are transported out of the supranuclear region (Figure 7b, c_supra(t)). The model accurately predicts the temporal evolution of vector concentration in the supranuclear region. The excellent agreement between the measured and calculated values of c(r,t), c_peri(t), and c_supra(t) is a strong indication of the model’s capability to predict transport of polyplexes inside cells.

The model was first used to examine the individual steps of the gene delivery pathway. The entire pathway was divided into two stages: 1), before endosomal/lysosomal escape, when the transport of vectors is determined by the endocytic vesicles; and 2), after endosomal escape, when the transport is determined by diffusion and stability of the vector. We first calculate the optimal location and time of endosomal escape that maximizes the probability of completing DNA delivery for in vitro cells. We then use this knowledge to predict the optimal pathways for in vitro and in vivo cells. Due to the complexity and multidimensionality of the parameter space, a global optimization of all parameters involved is not possible. Instead, we develop a more intuitive procedure to find the optimal set of parameters.

After escape from endocytic vesicles, polyplexes must unpack to release the plasmid, which can then enter the nucleus. Diffusion coefficients of polyplexes and plasmids in the cytoplasm are small, 10⁻⁴μm²/s and 10⁻³μm²/s, respectively 17. While the polyplex is chemically stable, the plasmid is susceptible to cytoplasmic nucleases and has a short half-life (∼50–90min) 15. Thus, only those plasmids that are close to the nucleus have a reasonable chance of success. To quantify this effect, we calculate the probability, Φ(r), that a polyplex, having escaped from an endocytic vesicle at distance r from the nuclear boundary at time t, will successfully deliver DNA to the cell nucleus (Figure 8a) (see Supplementary Material for mathematical definition of Φ(r)). The success probability Φ(r) decreases as vectors escape away from the cell center, which is located at r=−5μm. The value Φ(r) is generally high for supranuclear escape (r<0) and declines only slightly toward the nuclear boundary (r=0μm). On the contrary, for vectors that escape within the cytoplasmic region (r>0), the probability decreases exponentially with r, and at r =1μm, Φ(r) is already reduced to mere 0.01 (that is, 1%). The prediction that success probability is low for vectors that escape far away from the nucleus is intuitive. However, predictions in Figure 8a provide significant new information. First, they provide a quantitative determination of how success probability decreases with distance. Second, they provide a direct means of comparing the contribution from supranuclear and cytoplasmic vectors. This will become important in discussing in vitro transfection experiments in which the supranuclear region often occupies a significant fraction of the cell (discussed later).

Display large version of this figure

Figure 8

Likelihood of successful DNA delivery. (a) The value Φ(r,t), the probability that a vector, which escapes at a distance r from the nucleus at time t, will reach the nucleus. The probability is evaluated at 24h after the occurrence of the escape event (r>0, cytoplasmic region; r<0, supranuclear region). (b) The value Θ(t_e), the probability that a vector, initially located in the supranuclear (solid) or cytoplasmic (dotted) region, completes DNA delivery to the nucleus at 24h after escape from endosome or lysosome at time t.

A similar analysis was performed with respect to time of escape, t_e. Specifically, we calculated the probability Θ(t_e) that a vector escaping at time t_e post-transfection reaches the nucleus within 24h (see Supplementary Material for mathematical definition of Θ(t_e)). Θ(t_e) is essentially a spatially averaged value of Φ(r) for all vectors released in the cytoplasm at time t_e. As seen in Figure 8b, Θ(t_e) is very different for supranuclear and cytoplasmic vectors. For supranuclear vectors, Θ(t_e) decreases monotonically with t_e as the vectors are gradually transported out of the supranuclear region due to endosome-lysosome fusion (Figure 7b, c_supra). For cytoplasmic vectors, Θ(t_e) exhibits an optimum in t_e at ∼12–13h post-transfection. If endosomal escape occurs too early, when most of the cytoplasmic vectors are far away from the nucleus, the released vectors will degrade before reaching the nucleus. On the other hand, late escape also results in low delivery efficiency. First, the longer the vectors stay inside endocytic vesicles, the more likely they will unpack and degrade inside lysosomes. Second, the released vectors may not have sufficient time to find and to enter the nucleus before the simulation is terminated. The precise contribution of each effect to the existence and the actual value of optimal escape time depend on several parameters.

The data in Figure 8b provides several new insights into transport of gene vectors. First, it clearly establishes that the best strategy for supranuclear vectors is to escape as early as possible. On the other hand, there exists a clear optimum escape time for cytoplasmic vectors. But most importantly, the data reveals that due to the stark differences between the behavior of supranuclear and cytoplasmic vectors, the geometry of the cell must be taken into account while determining the optimal escape time. Different cell lines have different sizes of nucleus and nucleus area to cell area ratios. These factors must be taken into account while comparing data from different cell lines.

We next combined the data from the previous two sections and determined an overall efficiency Ψ, the probability that a vector released anywhere within the cell at any time delivers DNA to the nucleus. We then sought to determine vector parameters that optimize Ψ. We first defined the possible ranges of the parameter values that polyplexes are expected to possess. The parameter space is highly multidimensional, and consequently full optimization is computationally demanding. Instead, we followed a more intuitive path to reach the global maximum. The delivery efficiency exhibits a simple monotonous relationship with many parameters, e.g., rates of internalization and DNA degradation. We maximized these parameters to reach a semioptimal state, around which we performed local optimization within a lower dimensional parameter space. We focused on a three-dimensional region of the parameter space represented by k_escape (the rate of vector escape from endosomes or lysosomes), k_tr (the rate of vector transfer from endosomes to lysosomes), and k_unpack (the rate of vector unpacking, and consequently plasmid degradation either in lysosomes or in the cytoplasm). These parameters account for the delicate balance between endocytic trafficking, endosomal escape, and vector stability. All other vector parameters were held at their respective optima. It is important to note that the parameters under consideration are both vector-specific and cell-specific. For instance, our experiments show that the rate of transfer from endosomes to lysosomes in human fibroblasts, k_tr, depends strongly on the nature of the internalized cargo: k_tr for dextran, poly-lysine-DNA vectors, and PEI^25kDa-DNA vectors is 1/45min⁻¹, 1/60min⁻¹, and 1/600min⁻¹, respectively. Molecular conjugation of gene vectors with proper ligands will significantly alter the value of k_tr. The other parameters, e.g., k_degL, k_escape, and k_unpack, show the same dependency on the physicochemical properties of the vectors.

Under in vitro conditions (flat, adherent, two-dimensional geometry), optimization yielded relatively straightforward results, that is, vectors that escape early lead to high delivery. Supranuclear vectors have a much greater chance of success by virtue of their advantageous location. The relative number of supranuclear and cytoplasmic vectors depend on the ratio of nuclear area to the cell area. For cells with this ratio >0.06 (which is true for most cells), the contribution of supranuclear vectors Ψ is dominant. Hence, high rate of escape and unpacking maximizes efficiency (Fig. 9). Standard PEI^25kDa polyplexes, however, are far away from this maximal delivery state. The calculations show that the ratio of the maximum efficiency to the base state () is ∼80. Note that this enhancement ratio is based only on optimization of the intracellular pathway. Other factors such as toxicity of the polymer are not considered here.

Display large version of this figure

Figure 9

Delivery efficiency in two-dimensional cells Ψ_{in vitro} as a function of k_escape, k_tr, and k_unpack. The value Ψ_{in vitro} measures the fraction of vectors that successfully deliver DNA to the cell nucleus within 24h post-transfection. The minimum and maximum values of Ψ_{in vitro} are 0.8×10⁻⁷ (blue) and 2×10⁻² (red), respectively. The values k_escape, k_tr, and k_unpack are the rates of escape from endocytic vesicles, transfer from endosomes to lysosomes, and vector unpacking, respectively. Slices of the function at different values of (a) k_unpack and (b) k_tr.

In contrast to in vitro situations, cells under in vivo conditions are highly three-dimensional 24 and of a stellate shape, with a large central spherical portion (Figure 1b). Under these conditions, no region of the cell membrane is particularly close to the nuclear membrane. Due to this major difference in cellular morphology, the spatiotemporal distribution of vectors and hence optimal polyplex parameters are radically different from that in vitro. Fig. 10 depicts the predicted the spatial distribution of polyplexes in three-dimensional cells at 10min, 4h, and 24h post-transfection, calculated using the assumption that the principles and the parameters underlying endocytic trafficking are the same as those in vitro. At short times, there are no vectors within 2μm of the nucleus. All vectors are located close to the cell membrane, far-away from the nucleus. They are then gradually delivered to the perinuclear region by means of MT-dependent transport of endosomes, and endosome-to-lysosome transfer.

Display large version of this figure

Figure 10

Spatiotemporal distribution of vectors in three-dimensional cells. Predicted spatial distribution of PEI^25kDa-DNA complexes in a three-dimensional cell at 10min (solid), 4h (shaded), and 24h (dotted) post-transfection. is the probability that a vector is located at distance r from the nucleus at time t in a three-dimensional cell.

Fig. 11 shows optimization of vector properties for in vivo gene delivery. Optimal efficiency is observed at intermediate values of transfer and escape rate. The ratio of the maximum efficiency to the base state () is ≈40. Early escape of vectors from endosomes leads to underutilization of the endocytic trafficking and consequently, reduces perinuclear accumulation of the vectors, whereas a prolonged delay of vector escape promotes lysosomal degradation of the vector. Similarly, fast unpacking causes premature degradation of DNA inside lysosomes and cytoplasm, whereas slow unpacking limits the amount of DNA available for nuclear translocation. Interestingly, the optimal values of k_escape and k_unpack are functions of k_tr. Thus, to arrive at the optimum pathway, the rates of unpacking and escape cannot be tuned independently but are constrained by k_tr.

Display large version of this figure

Figure 11

Delivery efficiency in three-dimensional cells Ψ_{in vivo} as a function of k_escape, k_tr, and k_unpack. The color indicates Ψ_{in vivo}. Ψ_{in vivo} denotes the fraction of vectors that successfully deliver DNA to the cell nucleus within 24h post-transfection. The minimum and maximum values of Ψ_{in vivo} are 0.8×10⁻⁷ (blue) and 3.7×10⁻³ (red), respectively. The values k_escape, k_tr, and k_unpack are the rates of escape from endocytic vesicles, transfer from endosomes to lysosomes, and vector unpacking, respectively. Figure shows slices of the function at different values of (a) k_unpack and (b) k_tr.

Another important factor in the overall delivery efficiency is the stability of the polymer itself. The foregoing analysis is applicable to nonbiodegradable vectors like PEI^25kDa, which are immune to hydrolases, proteases, and nucleases, and offer a great degree of protection to their genetic cargo. For these vectors, it is reasonable to assume that the degradation of DNA is possible only after the vector unpacks. The stability of DNA in both the lysosomes and the cytoplasm is then intimately coupled to the rate of vector unpacking. On the other hand, biodegradable polymers such as poly(β-amino esters) can themselves undergo hydrolysis in the lysosomes or the cytoplasm. The DNA can thus be degraded without the need to unpack. In other words, there is no direct coupling between unpacking and DNA lysosomal degradation. Surprisingly, this latter case, which makes the DNA more susceptible to enzymatic attacks, leads to an optimal pathway slightly more efficient than that for nonbiodegradable polymers (0.6% vs. 0.4% for in vivo applications). This higher efficiency stems from the inherent flexibility afforded to design of biodegradable polyplexes. Since DNA degradation in lysosomes is no longer conditional upon unpacking of the polyplex, it can be independently optimized as another parameter. Under these circumstances, the efficiency varies linearly with the rate of unpacking k_unpack and rate of degradation in lysosomes k_degL, and no optima exist with respect to these parameters (Supplementary Material , section 2.5).

Fig. 12 shows the effects of the shape and the size of in vitro cells on Ψ. The shape of the cell is characterized by its circularity and size. The circularity C is defined as C=4π(A/P²), where A is the area of the cell and P is the perimeter. A small value of C indicates that the cell is elongated. The size is defined as the average distance between the cell membrane and the nuclear membrane. Elongated (less circular) and smaller cells have a larger fraction of the cell area nearer to the nucleus, leading to higher perinuclear concentration and hence a greater delivery efficiency. It is interesting to note that delivery efficiency can vary by almost an order of magnitude, only due to variation of shape and size of the cell.

Display large version of this figure

Figure 12

Effect of cell shape and size on the delivery efficiency Ψ_{in vitro} for dermal fibroblasts in culture. The color indicates Ψ_{in vitro}. The ratio between the maximum and minimum delivery efficiency is 6.2. The delivery efficiency increases with decreasing size and circularity. Low circularity often provides higher access to nucleus, hence the higher transfection efficiency. The circles denoted by A–C show three extremes in the circularity-size space that are captured by the cell population used. Note that k_escape=0.15h⁻¹, k_tr=0.1h⁻¹, k_degL=0.042h⁻¹, D_DNA=10⁻³μm²/s, k_deg-cyto,DNA=1h⁻¹, k_nuc=1/4h⁻¹, and k_unpack=0.25h⁻¹.

Traditionally, gene delivery by synthetic vectors has been studied and modeled as biochemical reactions between the vectors and cellular components. With recent studies based on single-particle tracking, it has become evident that spatially heterogeneous transport processes are intimately involved in every step along the delivery pathway 21,25. However, existing mathematical models of synthetic gene delivery approximate all transport processes by first-order kinetics 11,12,19. Such simple approximations not only undermine the predictability of the model but also fail to take advantage of data provided by microscopic experiments.

To overcome these issues, we introduced three new features in the present model. First, we introduced spatial coordinates to represent the cell geometry. The cell geometry is reconstructed directly from experimental images. Second, we follow the movements of single particles, just like in single-particle tracking experiments. Thus, intracellular trafficking of polyplexes is no longer represented as discrete jumps from one compartment to another, but as continuous movements of single particles. Third, we use stochastic algorithms to update the positions and the states of polyplexes and cellular organization. This is to reflect the inherent randomness of intracellular processes.

This approach has several advantages. First, it allows us to directly compare model predictions to experimental data. This proves to be crucial in not only model validation but also in verification of our understanding of underlying processes. Second, it enables us to develop a realistic and mechanistic description of intracellular transport phenomena. For instance, it accurately captures the basic physics of microtubule-dependent transport. The particles are allowed to switch intermittently between diffusion and motor-driven directional movements on MTs, just like what they do in experiments. Although the model requires several more parameters to characterize binding and detachment of particles to MTs, these parameters are real, have a physical meaning, and can be estimated directly from microscopic single-particle tracking experiments. Third, the stochastic algorithms grant us a great degree of flexibility in model implementation and solution methods, thus facilitating simultaneous representation of diverse physical and chemical processes that occur at multiple length and timescales. PDE-based models can describe spatially variable processes but are only restricted to simple and fixed cellular configurations. To obtain the solutions to the transport equations for complex geometry and dynamic environment, stochastic simulations are far more superior. Fourth, the model provides an integrated framework to systematically study the influence of all model parameters and extract optimal parameter sets.

The model relates molecular-level binding and trafficking events to whole-cell distribution of polyplexes, and eventually to gene delivery efficiency. The mechanistic treatment of transport allows us to study the effects of cell geometry, motor-assisted transport and endocytic trafficking. For the first time, the spatiotemporal distribution of polyplexes before escape from endocytic vesicles can be predicted with great accuracy. This information is crucial to the design of better vectors.

Perhaps the most important feature of the model is that it represents a step toward a more realistic description of intracellular transport in general and gene delivery in particular. Many processes in gene delivery pathway—for example, endosomal escape and DNA degradation—are extremely difficult to measure experimentally in situ. Accordingly, mechanistically sound models provide a valuable tool in studying the impact of these processes on gene delivery. Such models also provide a logical platform for synthesis and integration of information from diverse sources, bridging scales, cell types, and operating conditions, and allow predictions of cellular and subcellular processes that are inaccessible by current experimental tools.

The model presented here brings out several novel insights into mechanisms of gene delivery by polyplexes.

The model reveals, quantitatively, that the morphology of cells (e.g., size, shape, and dimension), which has been largely ignored in gene delivery research, strongly influences the delivery efficiency. This is significant since human cells exhibit a wide range of morphologies—muscle cells in the body are very flat and elongated, whereas dermal fibroblasts are stellate and round. It can be inferred from the model results that the strategies to be used for these two applications must be very different. In addition to the morphology of the cell, other cell-specific properties such as endocytic trafficking and microtubule-dependent transport can also affect design decisions.

The most conspicuous example of impact of geometry is comparison of in vitro and in vivo cells. Traditionally, synthetic vectors are optimized in vitro, based on the assumption that they will also perform optimally under in vivo conditions. It is apparent from our study that such assumption can be erroneous, as not only the optimal parameters but also the pathway to reach the optima are remarkably different for those two conditions. Figure 13a shows the ratio of different vector-dependent parameters at the optimal states and the base states of PEI^25kDa-DNA for gene delivery to dermal fibroblasts under in vitro and in vivo conditions. Under in vitro conditions, it is possible to have an unmitigated increase in the rates of escape and unpacking, while at the same time enhancing the delivery efficiency. For in vivo applications, the extent of increase of k_escape and k_unpack must be carefully controlled, lest the optima be overshot. Figure 13b shows the position of the base state (PEI^25kDa) and the conceptual route to the optimal pathway for in vivo applications. Overshooting the optimum in any direction can cause a significant reduction in the delivery efficiency (up to 1000-fold). Another interesting feature of Figure 13a is the difference in the optimal rate of transfer from endosomes to lysosomes, k_tr under in vitro and in vivo conditions. While this rate does not play a major role in vitro, it plays an important role in vivo and must be reduced to maximize the delivery efficiency. Also, the optimal values of k_escape and k_unpack, as discussed earlier, are strongly dependent upon k_tr. It should be noted that k_tr depends upon the direct interaction of the vector with the endocytic machinery, mechanistic details of which are not well understood. In summary, synergistic considerations of these seemingly independent processes (i.e., endolysosomal escape, transfer to lysosomes, and vector unpacking) are necessary for improving delivery efficiency.

Display large version of this figure

Figure 13

Design strategies to improve PEI^25kDa-DNA polyplexes. (a) Ratio of different polymer-dependent parameters of the optimal solution to those of PEI^25kDa for in vitro (solid bars) and in vivo (hatched bars) and (b) Ψ_{in vivo} as a function of k_unpack and k_escape for k_tr=1/500min⁻¹ (for PEI^25kDa). The position of PEI^25kDa is also indicated (yellow star). The arrow shows the trajectory to be followed to arrive at the optimum, assuming k_tr is held constant.

The model shows that under in vitro conditions, successful gene delivery is dominated by supranuclear vectors for a typical cell. This is interesting because supranuclear population of gene vectors is not as well studied as their cytoplasmic counterpart. The model implies that further detailed experimental studies are warranted on intracellular trafficking and the overall behavior of supranuclear vectors.

The model also makes a quantitative estimate of the extent to which synthetic vectors can be improved. As such, it complements a large number of experimental efforts focused on design and synthesis of new polyplexes. The model shows that PEI^25kDa-DNA, the gold standard of synthetic vectors, is highly suboptimal, with delivery efficiency 20–100-fold less than optimal configurations. However, more interestingly, it predicts that even under optimal conditions, the calculated values of and are only 0.4