ORIGINAL RESEARCH published: 11 June 2019 doi: 10.3389/fimmu.2019.01213 Spatial Lymphocyte Dynamics in Lymph Nodes Predicts the Cytotoxic T Cell Frequency Needed for HIV Infection Control Dmitry Grebennikov 1,2,3*, Anass Bouchnita 4 , Vitaly Volpert 3,5,6 , Nikolay Bessonov 7 , Andreas Meyerhans 8,9 and Gennady Bocharov 2,10* 1 Moscow Institute of Physics and Technology, National Research University, Dolgoprudny, Russia, 2 Marchuk Institute of Numerical Mathematics, Russian Academy of Sciences, Moscow, Russia, 3 Peoples’ Friendship University of Russia (RUDN University), Moscow, Russia, 4 Division of Scientific Computing, Department of Information Technology, Uppsala University, Uppsala, Sweden, 5 Institut Camille Jordan, UMR 5208 CNRS, University Lyon 1, Villeurbanne, France, 6 INRIA Team Dracula, INRIA Lyon La Doua, Villeurbanne, France, 7 Institute of Problems of Mechanical Engineering, Russian Academy of Sciences, Saint Petersburg, Russia, 8 Infection Biology Laboratory, Department of Experimental and Health Sciences, Universitat Pompeu Fabra, Barcelona, Spain, 9 Institució Catalana de Recerca i Estudis Avançats (ICREA), Barcelona, Spain, 10 Sechenov First Moscow State Medical University, Moscow, Russia Edited by: Michael Loran Dustin, The surveillance of host body tissues by immune cells is central for mediating their University of Oxford, United Kingdom defense function. In vivo imaging technologies have been used to quantitatively Reviewed by: Jens Volker Stein, characterize target cell scanning and migration of lymphocytes within lymph nodes Université de Fribourg, Switzerland (LNs). The translation of these quantitative insights into a predictive understanding of Martin Meier-Schellersheim, immune system functioning in response to various perturbations critically depends on National Institutes of Health (NIH), United States computational tools linking the individual immune cell properties with the emergent Frederik Graw, behavior of the immune system. By choosing the Newtonian second law for the governing Universität Heidelberg, Germany equations, we developed a broadly applicable mathematical model linking individual and *Correspondence: Dmitry Grebennikov coordinated T-cell behaviors. The spatial cell dynamics is described by a superposition dmitry.ew@gmail.com of autonomous locomotion, intercellular interaction, and viscous damping processes. Gennady Bocharov The model is calibrated using in vivo data on T-cell motility metrics in LNs such as the bocharov@m.inm.ras.ru translational speeds, turning angle speeds, and meandering indices. The model is applied Specialty section: to predict the impact of T-cell motility on protection against HIV infection, i.e., to estimate This article was submitted to the threshold frequency of HIV-specific cytotoxic T cells (CTLs) that is required to detect T Cell Biology, a section of the journal productively infected cells before the release of viral particles starts. With this, it provides Frontiers in Immunology guidance for HIV vaccine studies allowing for the migration of cells in fibrotic LNs. Received: 08 February 2019 Keywords: lymphoid tissue, cell motility, HIV infection, cytotoxic T cell scanning, multicellular dynamics, Accepted: 13 May 2019 dissipative particle dynamics, stochastic differential equation Published: 11 June 2019 Citation: Grebennikov D, Bouchnita A, INTRODUCTION Volpert V, Bessonov N, Meyerhans A and Bocharov G (2019) Spatial The surveillance of host body tissues by cells of the immune system is central for mediating Lymphocyte Dynamics in Lymph Nodes Predicts the Cytotoxic T Cell defense functions against invading pathogens and tumor cells (1, 2). The initial recognition of Frequency Needed for HIV Infection foreign antigens that leads to the induction of adaptive immune responses takes place in lymph Control. Front. Immunol. 10:1213. nodes (LNs), which, by virtue of their location and structure, facilitate the interactions between doi: 10.3389/fimmu.2019.01213 immune cells (3). The motility of pathogen spread and immune cells represents relevant parameters Frontiers in Immunology | www.frontiersin.org 1 June 2019 | Volume 10 | Article 1213 Grebennikov et al. T Cell Motility in HIV Infection Control controlling the fate of the pathogen–host interaction. In vivo dynamics in LNs. By choosing a first principles approach imaging technologies have been used to quantitatively in formulating the governing equations in conjunction with characterize target cell scanning and migration dynamics published experimental data on T-cell motility in lymphoid of lymphocytes within LNs (4, 5). The translation of these tissues, we offer a broadly applicable generic mathematical tool quantitative insights into a predictive understanding of immune linking individual and coordinated cell behaviors. The potential system functioning in response to various perturbations critically of the model is illustrated by an analysis of the combined effects depends on the availability of computational tools linking the of antigen-specific T-cell numbers and intrinsic T-cell motility individual immune cell properties with the systems response as a parameters in LNs on the time needed to locate both mobile whole (6). and non-motile HIV-infected target cells. Computed predictions Multiscale models of the immune system provide the in of the ratio of effector CTLs to infected T cells in the LN silico tool to embed immune processes into their spatial context paracortex needed for a timely detection of infected cells within (7–9). A core module of the models is the mathematical 18 h postinfection, i.e., before the release of viral particles starts framework used to describe individual cell migration in complex (20), provide a novel quantitative guide for an informed design multicellular environments. One can distinguish two general of HIV vaccines. types of modeling approaches, cellular automata-based models (CAMs), and physical models (PMs). CAMs consider a regular MATERIALS AND METHODS grid with cells that change their state in time and space according to some rules (functions of the system state). The respective Programming Languages and computational algorithms can take the form of random walks Computing Resources (10) or cellular Potts models (11). Although CAMs incorporate All algorithms were written in C++ and compiled using G++ experimentally defined characteristics of cell motion and, thus, (version 5.4.0). Pseudorandom numbers were generated using simulate cell dynamics based on actual data, they lack quantifiable the PCG random library (version 0.98) and the PCG64-XSL- links to the underlying biophysical interactions between cells RR algorithm (21). The seed was either specified manually (for in multicellular environments and to intrinsic cell motility code development) or set based on the system’s random device parameters (12). PMs of lymphocyte migration dynamics derived (for computational experiments). Simulations were run on a 2- from the Newtonian second law offer the possibility to define cell core Xeon E3-1220 v5 @3.0 GHz × 4 processor. The wxWidgets motions in terms of the forces generated by the environment and library (version 2.8.12) was used for visualization purposes. the cell itself. Using the experimental data on cell movement, The processing of the simulation results (i.e., calculating the potential functions underlying cell-to-cell interactions and statistical motility profiles, comparing CDFs, and plotting) was intrinsic cell motility can be identified and can provide a deeper implemented in Python and R scripts. insight into the mechanical properties of cells. Thus, PMs of individual cells and coordinated cell migration represent a Model Equations of Multicellular Dynamics general and generic way to describe and predict the multicellular According to a basic mechanics view, a system consisting of N system dynamics for a broad range of cell numbers and external cells of some mass located in a liquid milieu, interacting with each conditions (13, 14). other and affected by some external field, is uniquely determined It is widely accepted in immunology that the physiological by their coordinates and velocities and is governed by the classical function of cytotoxic T cell (CTL) motility is to search mechanics motion equations. In our model, each cell i, i = 1, N, for target cells, i.e., for virus-infected cells or cancer cells is represented as the circle with certain mass mi , radius ri , and (15). Computational modeling studies have revealed that the position of its center xi . The fundamental equation governing search efficiency depends on the organization of the stromal locomotion of cells is Newton’s second law of motion. It can be environment of a tissue (16). In addition, the spatial behavior, expressed as follows: for example, of HIV-infected target cells scanned for foreign X antigens by CTLs strongly impacts the elimination efficiency mi ẍ i = Fi = fijint + fimot + fidis , i = 1, N, (1) of the infected targets (17, 18). Experimental investigation of j6=i live attenuated SIV vaccines clearly suggested that a robust protection against intravenous wild-type SIVmac239 challenge where the first term on the right side specifies the net effect of strongly correlates with the number and function of antigen- the pairwise interaction forces with contacting neighbor cells, the specific effector CTLs in LN rather than the responses of such second term stands for the cell intrinsic locomotion force, by cells in the blood (19). However, the quantitative effects of T- which the cell establishes motility within the extracellular matrix cell migration parameters in LNs on the efficiency of antiviral (ECM) of the LN reticular network, and the last one takes into immune responses in vivo remain unknown. account the action of a dissipative force, taken to be proportional In the current study, we have developed a physics-based to the cell velocity fidis = −µẋi . We neglect the impact of gravity. description of spatial T-lymphocyte dynamics in the multicellular environment of LNs. A fundamental relationship between a cell Random Motility Force Sampling motion and the forces acting on it is provided by Newton’s The random motility force fimot for the ith cell is modeled as a second law. It is used to formulate, calibrate, and apply a stochastic vector fi sampled every 30 s from certain probability generic mathematical model of coordinated T-cell migration distributions analogously to the inverse homogeneous correlated Frontiers in Immunology | www.frontiersin.org 2 June 2019 | Volume 10 | Article 1213 Grebennikov et al. T Cell Motility in HIV Infection Control random walk (IHomoCRW) model (22). The motility magnitude neighboring cells in contact (such that the distance between cell fi = ηi · Ki is sampled from the following Gaussian distribution: centers hij ≤ ri + rj ), and R̂i determines the net shift of the Ki ∈ N µ (K) , σ 2 (K) . To obtain the motility magnitude fi ,  inner motility direction away from the neighboring cells, cinh is the sampled value Ki is multiplied by the arresting coefficient the weighting coefficient varying the level of CIL. The hat above ηi . The arresting coefficients are increased for both T cells and the vector indicates that it is normalized. DCs if they establish a sufficiently long contact to temporarily arrest their inner motility as follows: (1) ηi : = 10ηi default for T Numerical Integration of the Equations of cells and DCs when the duration of an uninterrupted contact Cell Motion default To numerically integrate the equations of motion (Equation 1), exceeded 30 s, and (2) ηi : = 100ηi if the contact outlasted 20 min. The cell inner motility is restored back to a default value we used the first-order semi-implicit (i.e., the cell coordinate at if the contact lasted for a time longer than the sampled value time t n+1 is computed using the velocity vector vin+1 rather than tcontact ∈ N (2, 0.4) hours. The parameter ηi is also used to vin ) Euler method: decrease intrinsic motility when performing in silico simulations mi vin + h · Fiint t n , xin + fimot t n , xin   to study the effect of decreased T-cell motility on target cell vin+1 = (4) location efficiency (see details in Supplementary Text). mi + h · µ The motility direction f̂i is turned from the previous direction xin+1 = xin + h · vin+1 (5) on the angle θi :  ! in which xin and vin are the coordinate and velocity of cell i at the Ki β  2  αi ∈ N 0, σ (α) , θi = αi · 1 − , time t n after n steps t n = t 0 +h·n. We note that the second-order Kmax generalization of this method, i.e., the Störmer–Verlet method, Kmax = µ (K) + 3σ (K) . (2) could be developed. However, it will be computationally more demanding as the cell acceleration depending on velocity due Here, N 0, σ 2 (·) denotes a Gaussian distribution, and β to the presence of dissipative velocity-damping viscosity forces  is a scalar coefficient. The angle sampled from the normal needs to be reevaluated at each time step t n+1 . We verified that distribution is multiplied by a factor depending on the sampled the time step h = 0.02 min used in the simulations is sufficient motility magnitude to reproduce the experimentally observed for a stable integration of the initial value problem with the negative correlation between cell translational and turning angle semi-implicit Euler method. To efficiently locate the neighboring speeds. Indeed, the cells do not simultaneously perform fast cells (which is needed for intercellular force calculations and for translational movements and large reorientations (22). Note determining the effect of CIL), we use a simple uniform-grid- that a similar feature was named “directional propensity” and based spatial neighbor search, which performs well for a densely modeled with trigonometric parameterization in a cellular Potts packed multicellular environment. Note that the convergence of model to describe the motion of T cells (11). The Gaussian the integration scheme was verified by repeating simulations for distribution for the motility magnitude is set so that the a smaller time step. (µ − 3σ , µ + 3σ ) range is positive. The absolute value is taken to ensure that the magnitude is non-negative. The parameter Boundary Conditions Kmax provides an upper boundary for sampled values Ki During the model calibration process, we used periodic boundary (approximately 1 of 370 cases falls outside of the three-sigma conditions for all boundaries of a square domain. To perform in interval). The hat above the vector denotes the normalized silico simulations in a closed ellipse-shaped domain representing unit vector. a LN, we implemented a biologically based boundary condition of cell repolarization. We do not model explicitly the interaction Implementation of Contact Inhibition forces between cells and the boundary (i.e., the subcapsular sinus of Locomotion wall). At the stage of coordinate updates (in accordance with the After the stochastic vector fi is sampled, it is modified in numerical scheme specified in the section Numerical Integration accordance with the contact inhibition of locomotion (CIL) of the Equations of Cell Motion), if the proposed coordinate of model, as described (23). The resultant vector fimot is then used cell xin+1 is outside the boundary, the current coordinate of the in the right-hand side of Equation (1). The modification consists cell is preserved (xin+1 = xin ), while the direction of the motility of shifting the direction of vector fi away from the neighboring vector Efimot is changed to be the opposite direction of vector vin+1 , cells and decreasing the magnitude of vector fi proportionally to thus resulting in cell repolarization. the number of neighboring cells:   Generating the Initial Spatial Configuration fi · cinh f̂i + R̂i X xi − xj for Simulations Within a LN fimot = , R̂i = , (3) To generate the initial spatial configuration of the immune cinh + n hij j,hij ≤ri +rj cells within a LN, we followed the descriptions from a LN imaging study (24). The following cell subsets were considered: in which fi is the magnitude and f̂i is the direction of the inner CD4+ T cells, CD8+ T cells, and cross-presenting migratory motility as it would be if unaffected by CIL, n is the number of CD8αint CD103hi DCs. Both T-cell subsets are distributed Frontiers in Immunology | www.frontiersin.org 3 June 2019 | Volume 10 | Article 1213 Grebennikov et al. T Cell Motility in HIV Infection Control uniformly through the whole LN, while migratory DCs are found are expected to be present between all cells according to the mainly deep in the paracortex area. To arrange cells in agreement model of Bell (28). We calibrate this parameter by the typical with the experimental data, we approximated the DC-rich area value of low-adhesive forces, with which integrins present on T- α=0.99 . The spatial positions for DC locations are as an ellipse ΩDC cell membrane bind to their ligands present on the other cells iteratively sampled from the 2D Gaussian distribution with a 99- (29). For cognate T cell/APC interactions the attraction force α=0.99 and accepted if the DC with sampled percentile ellipse ΩDC is much stronger as it is determined by a broad spectrum of coordinates lies within the LN domain ΩLN and does not overlap various adhesion molecules involved in T-cell activation clusters, with the other seeded DCs. After DCs are placed, the T cells are i.e., the immunological synapse (30). The estimated values of positioned uniformly through the remaining non-occupied space the intercellular interaction forces are given in Table 1. For of ΩLN . details on the data-based T-cell motility model calibration, see Supplementary Text. RESULTS The dissipative (friction) force acting on T cells describes the effect of viscous damping, which reduces the velocity Biophysical Parametrization of the Spatial of the cell. It is assumed to be proportional to the cell Multicellular Dynamics velocity fidis = −µẋi . The dissipative force acts along the Multicellular systems dynamics can be accurately described by line of the cell center and in opposite direction to the cell biophysical models as reviewed recently (14, 25). Here, we displacement. Consideration of viscous damping is appropriate develop a physics-based mathematical model of coordinated for the highly viscous low-Reynolds-number environment of immune cell motion that belongs to the class of self-propelled LNs (40). The viscous damping parameter estimate is listed particle models (14) and, more generally, to the dissipative in Table 1. particle dynamics (26, 27) framework. The random motility force fimot determines the traction of Immune cells in LNs are continuously interacting with each self-propelled lymphocytes. It represents a stochastic process other and with stromal cells via forces of different origin, i.e., of receptor-mediated cell–ECM interactions regulated by either elastic (membranes), chemical (receptors), and electric. The cytoskeletal or membrane reorganizations and governed by respective forces in combination with cell intrinsic locomotion biomechanical and intracellular molecular mechanisms (4, 13). events act in concert to determine the basal intranodal motility Basically, cells establish directed caterpillar-like movement of T cells. Figure 1A presents the overall summary of physical by polarizing, forming contacts between their leading edge forces included in the model with some implementation details. and collagen fibers of ECM, detaching their trailing edge The scheme of the forces exerted on cell i interacting with cells from ECM, and contracting. However, T lymphocytes and p and k is shown in Figure 1B. The quantitative features of dendritic cells (DCs) are characterized by low-adhesive integrin the force functions are detailed in Figure 1C. Here, fijint is the interactions with the microenvironment. This allows them to intercellular force acting on cell i due to interaction with cell j. adapt their direction and morphology with no need to reorganize The pairwise cell-to-cell interactions are assumed to have a finite microstructure while effectively sliding along the stromal cutoff distance and are considered to be elastic acting along the network of fibroblastic reticular cells (41, 42). As we do not line of cell centers. The intercellular forces fijint can be considered model the reticular network and the ECM microstructure in this as the gradients of pairwise potentials, which are repulsive at study explicitly, this motility behavior is considered implicitly short distances and attractive at larger distances, thus accounting in the stochastic nature of fimot . Note that the autonomous cell for volume exclusion at the cell body and cell–cell adhesion motility can also be affected by external signaling, e.g., through near membranes. We consider the following cubic polynomial chemotaxis, CIL, or immunological synapse formation. The cell function to model the force exerted by cell j on cell i: trajectory in the model is characterized by three quantifiable  3 values, i.e., the translational speed, the turning angle speed, xi − xj −a · f adh · rj −x + b · f adh · rj −x , hij < ri + rj ,  int r r and the meandering index as explained in Figure 1D and as fij = · j j (6) hij  0, hij ≥ ri + rj , described in Read et al. (22). The corresponding experimental data are shown in Figure 1E. To capture the experimentally observed patterns of T-lymphocyte migration in lymphoid tissues where ri is a radius of the ith cell membrane, hij is the distance (see Figure 1E), the T-cell motility is modeled using a random between cell centers (see Figure 1C), and x = hij − ri is the variable fi with its magnitude and angle values updated every distance between the center of cell j and the surface membrane 1t seconds according to the IHomoCRW recently suggested r −x of cell i. The function a · f adh · j rj describes the attraction force and validated (22). The IHomoCRW model was shown to r −x 3 reproduce the experimentally measured statistical profiles of   between two cells, and the function b · f adh · j rj corresponds T-cell locomotion (22). In the present model, the magnitude to a repulsive force, both calibrated as shown in Figure 1C. and direction of the random vector fimot are sampled from The coefficients a and b are set such that the minimum of distributions provided by the experimental data (the specific rules function fijint is equal to f adh . Thus, the only remaining free are defined in section Materials and Methods). The key difference parameter is f adh , the adhesive interaction strength. In the from the original IHomoCRW model is that it is the cell motility case of T cell/T cell interaction, it corresponds to weak non- inducing force fimot rather than the cell velocity ẋi that is sampled specific electrical forces (electrostatic and electrodynamic) that and then substituted into equation (Equation 1). In addition, Frontiers in Immunology | www.frontiersin.org 4 June 2019 | Volume 10 | Article 1213 Grebennikov et al. T Cell Motility in HIV Infection Control FIGURE 1 | Physics-based model of multicellular system dynamics reproduces experimental data on T-cell locomotion. (A) The set of forces considered in the model with description of their features and implementation details. (B) The fundamental equation governing locomotion of cells determined by the forces exerted on cell i, (Continued) Frontiers in Immunology | www.frontiersin.org 5 June 2019 | Volume 10 | Article 1213 Grebennikov et al. T Cell Motility in HIV Infection Control FIGURE 1 | including the repulsive–attractive interaction with neighbor cells p and k, respectively. (C) The parameterization of intercellular interaction force fijint and formula definition. The calibrated force for non-specific interaction of two T cells with a radius of 3 µm is depicted. By simulation, the parameters a and b are calculated at each time step depending on the radii ri , rj and the distances hij , x, so that the condition fij λrj = fij rj = 0, minfij (x) = −fiadh is satisfied. The   parameter λ determines the relative deformation of the cells that separates the repulsive and attractive interactions between them. Parameter fiadh represents the adhesive strength between the membranes of cells i and j. (D) Schematic illustration and definition of the metrics characterizing T-cell motility: translational speed, turning angle speed, and meandering index. All metrics are measured for each cell every 1t seconds and pooled together to form statistical distributions. (E) Statistical profiles characterizing the T-cell locomotion consists of distribution histograms of translational speeds, turning angle speeds, and meandering indices. The histograms are derived from the corresponding empirical cumulative distribution functions (CDFs) available in Figure S17 from Read et al. (22), in which original in vivo data are presented. (F) The details of the 2D geometric setup for simulations used in the model calibration: spatial configuration, initial and boundary conditions, and the experimental protocol used to sample the statistical profile. (G) The statistical characteristics of T-cell motility coming from simulations of the calibrated model plotted against the in vivo histogram data (22). The statistical distributions of each metric are depicted as CDFs. The Kolmogorov–Smirnov statistics comparing the model and target CDFs are indicated with their respective p-values. TABLE 1 | Set of calibrated model parameters used as a baseline for all simulations. Parameter Description Value References mTC T-cell mass 215 pg (31–34)†* rTC T-cell radius 3 µm (35)† µ Viscous damping coefficient 0.2 nN·min/µm (= 12 g/s) (36, 37)*t adh fT−T Adhesive strength between T cells 0.01 nN (28–30)† λTC The normalized distance between the cell-to-cell interaction synapse and the cell 0.83 (with respect to r TC ) (38)*t center separating repulsive and attractive modes of T-cell interaction (model analog of nuclear-to-cytoplasmic ratio in experimental cell biology) 1t Time step for inner motility mTC update 30 s (22) ‡ (22, 39)‡*t  µ mTC Mean of the inner motility force magnitude distribution 3 nN (22, 39) ‡*t  σ mTC SD of the inner motility force magnitude distribution 0.3 nN 60◦ (22)‡t  σ αTC SD of the inner motility turning angle distribution cinh Scalar coefficient varying the level of CIL 1.0 (23)‡t β Scalar coefficient in parameterization of negative correlation between magnitude 2.0 (23)‡t and turning angle of sampled inner motility † Parameters obtained directly from experimental measurements. *Parameters estimated indirectly from experimental measurements. ‡ Parameters derived from underlying computational models. t Parameters tuned to fit cell motility profiles within the model calibration. the random vector fimot can be influenced by contact effects acting on a randomly selected T cell in a simulation of from neighboring cells, resulting in (1) a shift of the vector fimot multicellular dynamics are shown in Figure S1. away from neighboring cells and (2) a decrease of its magnitude proportionally to their number, similar to the CIL model (23) (see details in section Materials and Methods). By default, the Calibration of T-Cell Motility arresting coefficient for T cells is equal to one. For DCs, its value Our model mostly operates with biophysical parameters that is estimated so that the resultant DC velocities do not exceed 5 are either directly measurable or can be estimated indirectly µm/min (the estimated value is specified in Table 2). such as the mass m (wet weight) and the radius r of a cell, the Overall, the mechanistic model of the spatial multicellular adhesive strength between T-cell membranes fijadh (measured by dynamics is formulated as a system of N random ordinary single cell force spectroscopy), the viscous damping coefficient differential equations (44) represented by Equation (1) and µ, typical forces and velocities of T cells, and the location embedded into the 2D geometric domain as detailed in of demarcation between repulsive and attractive areas of a Figure 1F. Essentially, the system is a deterministic system cell λ (nuclear-to-cytoplasma ratio). The other parameters that of ordinary differential equations on each interval of 1t describe the random motility force or the contact inhibition seconds, until the force fimot becomes updated. The quantitative of locomotion are derived using the information presented consistency of the computational model of multicellular in the original IHomoCRW model (22) and the CIL model dynamics with experimental data on translation speed, (23) with the underlying experimental data. To calibrate our turning angle speed, and the meandering index is illustrated model, we evaluated admissible ranges of parameters and tuned in Figure 1G. The relevant components of the numerical them manually to match the statistical characteristics of T- implementation of the model (computational domain, boundary cell locomotion (22). The baseline sets of the estimated model conditions, integration algorithm) are described in Materials and parameters are presented in Tables 1, 2. For details of the Methods. The dynamics of the net forces and their contributions parameter estimation, see Supplementary Text. Frontiers in Immunology | www.frontiersin.org 6 June 2019 | Volume 10 | Article 1213 Grebennikov et al. T Cell Motility in HIV Infection Control TABLE 2 | Extra parameters which are used for LN simulations. Parameter Description Value/range References mCD4 CD4+ T-cell mass 215 ± 28 pg (31–34) mCD8 CD8+ T-cell mass 290 ± 28 pg (31–34) mDC Dendritic cell mass 350 ± 28 pg (32–34) rDC Dendritic cell radius 6.5 µm (38, 43) λDC The normalized distance between the cell-to-cell interaction synapse and the cell 0.5 (38) center separating repulsive and attractive modes of DC interaction (model analog of nuclear-to-cytoplasmic ratio in experimental cell biology) adh fDC Adhesive strength between specific T cells and DCs 1 nN (28, 30, 36) ηDC Default value for coefficient arresting inner motility of DCs 3 Tuned so that DC velocities are <5 µm/min Computational Domain, Immune Cell Subsets, and motility force was generated randomly for all cells. The positions Initial Configuration of cells were saved every 30 s during 10 numerical experiments of The computational domain was implemented as an ellipse- 30-min simulation time after a 30-min pre-run to randomly mix shaped 2D approximation of the bean-like cross section of a the cells. Cells with total displacements <27 µm were excluded murine skin-draining LN (see Figure 2A). At the beginning, as was done in the original experimental protocol. Likewise, cells both CD4+ T cells (green) and CD8+ T cells (blue) are evenly that passed through the boundary and left the imaging volume distributed throughout the domain. Some randomly chosen T were also excluded. The saved cell positions were post-processed cells are considered to be antigen specific and marked in light to calculate the target metrics (defined in Figure 1D), which green and blue, respectively (their numbers are specified below). were pooled together to form three separate distributions. The The antigen-presenting cells considered in this study represent pooled cell motility distributions were calibrated with in vivo the subset of cross-presenting migratory CD8αint CD103hi DCs, data. The simultaneous adjustment of all distributions was which are mainly involved in CD8+ T-cell immune responses computationally challenging due to the different uncorrelated and which immigrate into LNs from the periphery (24). They aspects of cell migration captured in each of the motility metric as are normally localized in the deep parts of the T-cell zone previously outlined (22). Figure 1G shows the best-fit cumulative and leave LNs slowly with a turnover rate of 6 days. For distribution functions (CDFs) of the calibrated model with the initial configuration, these DCs are spatially placed according baseline parameter set from Table 1 and the target experimentally to a Gaussian distribution with 99-percentile ellipse ΩDC α=0.99 observed distributions with Kolmogorov–Smirnov statistics and representing the T-cell zone (see Figure 2A and section Materials p-values describing the discrepancy between CDFs. and Methods). The evolution of the above multicellular system was The numbers of antigen-specific DC and T-cell subsets simulated over a 12-h period. The visualization of the systems are estimated using published data (24), which were rescaled spatiotemporal dynamics is presented in Movie S1. Figures 2B,C according to the size of the computational domain. A total population of 12,469 immune cells was considered. The total shows the kinetics of median velocities of antigen-specific CD4+ number of non-antigen-specific T cells was estimated so that T and CD8+ T cells, and the median distances between the about 80% of the computation domain was filled up. The T cells and the centroid of their cognate antigen-presenting precursor frequency of antigen-specific T cells, that is, their DCs throughout 12 h of an in silico experiment. The model proportion in the total amount of T cells, was set to be about demonstrates that antigen-specific CD8+ T cells that interact 1%. We consider the inflow and outflow of immune cells to with their cognate CD8αint DCs but not the CD4+ T cells the region of interest to be negligible because of the short decrease their velocities, move closer to the area of DCs in simulation time of 12 h. The closed boundary conditions used in the first 4–6 h, and remain there with low velocities afterward. the simulations are specified in section Materials and Methods. Figures 2B,C is quantitatively consistent with experimental data The overall geometrical scheme of the computational domain shown in Figures 1E,F, and in Figure 2B from Kitano et al. (24). and the initial configuration of the multicellular system generated for simulations are presented in Figure 2A. Quantitation of the DC and T-Cell Contact Interactions Data Assimilation and Model Validation The calibrated mathematical model of T-cell locomotion was To assimilate the statistical data on the three T-cell locomotion validated by confronting its predictions with data from the measures (i.e., the translational speed, turning angle speed, and intranodal spatiotemporal dynamics of different immune cell meandering index), the following numerical simulation protocol subsets after soluble antigen immunization presented in a recent was used, which is close to the original experimental protocol experimental study (24). The data specify the evolution of the (22). First, the same 2D 412 × 412 µm2 domain was used, in distances between the centroid of the migratory DC area and which we initialized 4,489 squarely tiled T cells with 3-µm radii individual CD4+ T and CD8+ T cells. The model was adjusted and η ≈ 80% packing density. The initial direction of the intrinsic to the functional configuration of skin-draining LNs specified Frontiers in Immunology | www.frontiersin.org 7 June 2019 | Volume 10 | Article 1213 Grebennikov et al. T Cell Motility in HIV Infection Control FIGURE 2 | Heterogeneous dynamics of T cells in LNs. (A) The scheme of a LN and illustration of the initial configuration generated for simulations. DCs, CD4+ T cells, and CD8+ T cells are placed within a LN as described in the Supplementary Text with total cellularity of 12,469 cells, ≈80% packing density and ≈1% precursor frequency. (B) Twelve-hour kinetics of median velocities of antigen-specific CD8+ T and CD4+ T cells, and their distributions at the start and at the end of a 12-h simulation. (C) Twelve-hour kinetics of median distances from T cells to the centroid of DCs, measured for antigen-specific CD8+ T and CD4+ T cells, and their distributions at the start and at the end of a 12-h simulation. TC, T cell; DC, dendritic cell. in the above study. A representative example of the numerical antigen-presenting DCs and antigen-specific T cells, and (3) simulation of individual cell trajectories is shown in Figure 3A. antigen-presenting DCs and polyclonal T cells. The values of the An example of multicellular dynamics in a LN during 12 h is corresponding parameters are presented in Table 2. The physical shown in Figure 3B. forces driving the dynamics of individual cells in the LN and the respective velocities of the cells predicted by the model are Quantitation of the Forces Determining T-Cell and DC shown in Figure 3C. We assume that the intrinsic motility of Motility and Their Interaction DCs can be represented by the same type of force function as To consider DCs in multicellular system simulations, we carried that for T cells (Figure 1C); however, due to their much smaller out a parameterization of their intrinsic motility forces and average velocity, the respective DC force function value was the intercellular forces for contacts between (1) two DCs, (2) of small magnitude. The adhesive force for cognate contacts Frontiers in Immunology | www.frontiersin.org 8 June 2019 | Volume 10 | Article 1213 Grebennikov et al. T Cell Motility in HIV Infection Control FIGURE 3 | Quantitation of immune cell motility, driving forces, and contacts. (A) Representative example of individual cell trajectories obtained with numerical realization of the calibrated model. The trajectories illustrate the 5-h dynamics of 15 cells randomly chosen from 4,489 cells presented in the 412 × 412 µm2 domain with periodic boundary conditions. (B) Twelve-hour multicellular dynamics of T-cell trajectories in a lymph node obtained by numerical simulation with an initial configuration specified in Figure 2A. Only cells with total displacement longer than 27 µm are shown. (C) Values of forces and cell velocities driving the multicellular system dynamics in a square subdomain of a LN. In a center pane, the velocity field is represented as a contour plot of the field of cell velocity magnitudes linearly interpolated at uniform grid, as well as detected streamlines of possible cell flow patterns. (D) Kinetics of the numbers of cognate DC–T cell contacts at different stages of the simulation and distribution of durations of all cognate contact durations occurring within a 12-h simulation. DC, dendritic cell. (i.e., of antigen-specific T cells with antigen-presenting DC) is computational procedure to temporarily arrest the motility for T around 100 times higher (∼1 nN) than the non-specific adhesion cells in a sufficiently long cognate contact (see section Materials force for T cell/T cell contacts (36). We also implemented a and Methods). Frontiers in Immunology | www.frontiersin.org 9 June 2019 | Volume 10 | Article 1213 Grebennikov et al. T Cell Motility in HIV Infection Control Ag-Specific CD8+ T Cells Migrate Toward target cells before they can release viral progeny. We consider Cross-Presenting DCs and Form Cognate Contacts only one HIV-infected cell in the computational domain, which With Them is consistent with the frequency of productively infected CD4+ T Figure 3D presents the model prediction for the kinetics of cells of about 0.0001–0.001 (54). Specifically, the newly infected the number of cognate DC–CD8+ T-cell contacts occurring target cell should be located by the nearby effector cells before it at different time intervals during the in silico simulation. can release viral progeny, i.e., before completion of the 18–24 h Antigen-specific CD8+ T cells robustly increase the number life cycle of HIV (20). of contacts with DCs over time in the process of T-cell zone The overall simulation setup is the same as described in the scanning for antigen-presenting target cells. Although most of model validation subsection above. Randomly chosen cells in the cognate contacts are of short duration, i.e., they last for the stochastically generated multicellular system configurations <5 min, the distribution has a heavy tail of stable more than 1-h representing the LN cortex zone were marked as infected in length contacts. These predictions are in agreement with previous yellow (see Figure 4A). Both the motile CD4+ T cells and data (45). the non-motile DCs were considered as HIV-infected targets. In simulations, we varied the frequency of HIV-specific CD8+ T cells and the intrinsic motility of T cells (searching CD8+ CTL Frequency Needed to Locate T cells, infected- and uninfected CD4+ T cells) (Figure 4A) HIV-Infected Target Cells Before to analyze the effect of variations on the target cell detection Viral Release time. A 10-fold range of HIV-specific CD8+ T-cell frequencies During viral infections, the induction of cellular immune typical for HIV infection, i.e., from 0.4 to 5%, was examined. responses takes place in secondary lymphoid organs such as The intrinsic motility of T cell was varied within 100 and LNs and spleen. Antigen-presenting cells such as DCs take up 50% relative to the calibrated baseline parameters of average antigens and migrate to LNs to encounter specific lymphocytes, T-cell velocity (see details in the Supplementary Text). A e.g., CD4+ T and CD8+ T cells, to induce their activation and decreased intranodal T-cell motility (of searching CD8+ T differentiation into effector and memory cell subtypes (46). The cells, infected- and uninfected CD4+ T cells) is expected to low frequency of antigen-specific T cells in unprimed hosts turns take place during the chronic stage of an HIV infection when the scanning of cognate DCs by specific T cells in a highly LN tissues become fibrotic, i.e., when collagen formation in crowded cellular LN environment into a “needle-in-a-haystack” T-cell zones takes place (55). Then, T cells have to move problem (47). It was revealed that optimal LN scanning depends through increased collagen deposition with major consequences on a combination of intrinsic T-cell motility, the chemokine for search patterns (56). In addition, CD4+ T-cell migration milieu, and the microarchitecture of the LN (1). When virus- is also inhibited by the HIV-1 Nef protein as shown in infected DCs reach the LN, the less the time needed to locate chemotaxis assays (57). In our study, the motility of all virus-specific T cells and to form stable DC–T cell contacts, the considered types of T-cell subsets, i.e., the searching CD8+ T more likely is that the precursor CTL activation will happen cells, uninfected CD4+ T cells, and infected CD4+ T cells, is before the viruses will be released from infected cells, therefore decreased uniformly. making the elimination of local clusters of infection spread more Figures 4B,C illustrate the model predictions for the decrease probable. This aspect of CD8+ T-cell activity is crucial for a of time to locate HIV-infected target cells with the increase prompt activation of specific CTL immune responses and the of HIV-specific CD8+ T-cell frequency. The modeling results elimination of viruses. The precursor frequency in blood can be imply that 5% is a sufficient effector CTL frequency for a as small as 0.0001% (48), reaching about 5–10% in the chronic timely detection of both types of target cells within 18 h stage of an HIV infection (49). The here-developed physics- post-infection, i.e., before the beginning of HIV release from based model of T-cell dynamics can be directly used to study productively infected cells (20). A stepwise increase by five- the efficiency of scanning the paracortical T-cell zone of the LN fold of the HIV-specific CTL frequency from 0.04 to 5% for target cells expressing cognate antigen as a function of the increases the probability of detection of HIV-infected cell within frequency of CTL and their motility. 24 h from 0.07 to 0.34 to 0.84 and to 1, respectively. In Development of an effective AIDS vaccine remains a addition, the model shows that infected motile CD4+ T cells global priority, and there is a need for a vaccine to induce are located faster than non-motile DCs with the probability of cellular immune responses capable of eradicating or efficiently detecting them within 24 h increasing from 0.35 to 0.86 and containing virus replication (50). Experimental studies with to 1 with a CTL frequency rising from 0.04 to 0.2% and to attenuated SIV vaccines indicated that SIV-specific CTLs, if 1%, respectively. present in sufficient frequencies, can completely control and even Figure 4D shows the increase of time to locate infected non- clear an infection (19). Similar to SIV, HIV infection is sustained motile DCs for a decrease of the T-cell motility from a basal by the activation of CD4+ T cells, which occurs in the form level by 10 and 50% considering an HIV-specific CD8+ T-cell of transient bursts in the local microenvironment of lymphoid frequency of 1%. If the average T-cell velocity is decreased by tissues (51, 52). The proximal activation and transmission 50%, then the probability to locate DCs within 24 h is <0.5. involving latently infected cells represent locally propagating Figure 4E depicts a similar dynamics for locating motile infected events (53). Therefore, we applied our calibrated model of spatial CD4+ T cells. Note that the motile targets were located within immune cell dynamics in LNs to study the necessary conditions 24 h even with a 50% decrease of the average CD8+ T-cell velocity for effector HIV-specific CTLs to promptly locate HIV-infected in all performed simulations. Frontiers in Immunology | www.frontiersin.org 10 June 2019 | Volume 10 | Article 1213 Grebennikov et al. T Cell Motility in HIV Infection Control FIGURE 4 | Conditions to locate HIV-infected target cells within a LN before viral release. (A) General scheme of in silico simulations. Time since the HIV-infected target cell was introduced until it was located by effector HIV-specific CTLs was measured in 24-h simulations. The infected cell was either non-motile DC (B,D) or motile CD4+ T cell (C,E). In (B,C), the precursor frequency, i.e., the frequency of effector T cells, was varied, from 0.04 to 5%. In (D,E), the effect of decreased intrinsic motility of T cells was studied. The average T-cell velocity was decreased up to 50%. In all plots, the fraction of cases with location time >24 h is indicated, thus providing the estimates for probability to locate target cells within 24 h. The time range between the start and the peak HIV release from infected T cells (20) is shown in pink. It is used to estimate the probability of a virus burst to escape effector CTLs and, thus, to contribute to the spread of HIV-infected cells within a LN. TC, T cell; DC, dendritic cell. DISCUSSION three spatial metrics of multicellular systems behavior, i.e., translational speed, turning angle speeds, and meandering index. We have developed a biophysics-based computational model The model provides the tool to quantify the velocity and of T-lymphocyte motility that is calibrated using empirical the driving force fields in the LN. It enabled us to predict in vivo data on T-cell migration in LN tissue representing frequency and motility parameters that are required for a Frontiers in Immunology | www.frontiersin.org 11 June 2019 | Volume 10 | Article 1213 Grebennikov et al. T Cell Motility in HIV Infection Control timely detection of productively HIV-infected cells within LNs range of modeling tools available for implementing integrative before they release viral progeny. As such, our study provides approaches to the exploration of the immune system. a quantitative guide for an informed design of HIV vaccines. CPMs have also been applied previously to study intranodal Furthermore, as the immunological principles of antigen-specific T-cell migration (58). The CPM framework is a valuable tool T-cell activation and immune surveillance imbedded in our for a phenomenological description of multicellular patterning, model also apply to other infections and cancers, our findings providing realistic simulations of morphological changes for may be used to define the general requirements for any efficient various cell types. The strength of this approach stems from immunotherapeutic intervention against pathogens or cancers its flexible energetic formalism that allows for extensions to in relation to disease-specific parameters and states of lymphoid incorporate various biological processes (62). Although the tissue and T cells. Thus, our model has a significant potential to CPM framework has a richer potential for describing individual guide the search for better and more efficient immunotherapies cell dynamics, including the cell shape, this comes at the in the near future. expense of (i) a higher-dimensional representation of the cell Other processes, e.g., chemotaxis, haptotaxis, and others, configuration (e.g., the number of voxels or pixels), (ii) the can influence the efficacy of target cell search by CTL. The use of phenomenologically rather than biophysically defined impact of chemotactic migration of T cells toward DCs has parameters, and (iii) a much higher computational cost to been computationally analyzed using a cellular Potts Model (58), perform simulations required to explore T-cell search strategies. an agent-based model (59), and a multicompartmental spatially Besides, there is no direct correspondence of most of the CPM resolved stochastic model of T-cell circulation (60). The results parameters with biophysical properties of cells, and the meaning suggest that the chemoattraction toward target cells modestly of some CPM parameters is still under debate (62, 63). Moreover, speeds up the search process for T cells that successfully find CPM temporal kinetics obtained with the modified Metropolis the chemokine-producing DCs. However, a qualitative model algorithm does not preserve the detailed balance condition for presented in (59) suggested that with even weak chemotaxis, the underlying stochastic process. This implies that the exact substantially lower numbers of CTL are required for sterilizing relation between forces of cell interactions and energy terms of immunity. Further data-based model-driven research is needed CPMs cannot be obtained even for the overdamped dynamics to clarify the contribution of chemotaxis to T-cell migration approximation (62, 63). under normal conditions and during inflammation (61). Computational modeling of multicellular dynamics in Phenomenological Ordinary Differential Equation (ODE) lymphoid tissues provides a theoretical tool to be used for a models may also be developed to simulate the interactions better understanding of the determinants of efficient immune between cell populations in the LN. However, these models are responses against pathogens with a final aim of an optimal not suitable for the present study for three reasons. First, data on manipulation of the immune systems performance (2, 15). T-cell motility in the LN cannot be directly used to calibrate such Given that the quest for an effective HIV vaccine remains models, thereby limiting the validity of their predictions. Second, a global priority (64) and that the localization, migration, the objective of our study, which is the early detection of HIV- and frequency of CTLs in LNs determine the extent of virus infected T cells and DCs, requires the monitoring of the spatial elimination (17, 19, 36, 56, 65), we sought to use our modeling density of T cells in the LN rather than the total number of T cells. approach to define threshold frequencies of CTLs in LNs for Changes in the spatial distribution of T cells in the LN can be protection against HIV. Since an HIV infection can influence related to spatial mechanisms such as chemotaxis and migration. (i) CD4+ T-cell motility by a direct mechanism involving Therefore, it is crucial to consider spatial aspects in the model. the HIV Nef protein and (ii) CTL locomotion via an indirect Finally, ODE models based on “mass action”- or “predator– mechanism related to the induction of lymphatic tissue fibrosis, prey”-type parameterizations would require the parameter values we considered both phenomena to predict the effect of reduction specifying a per capita killing rate of target cells. The respective of T-cell motility. We estimated that the frequency of antigen- parameter can be determined by the mean time needed for a specific CTL should be about 5% to timely detect and completely migrating CTL to locate infected cells. A priori estimates of this eliminate productively infected DCs within 18 h. The time parameter are not available. It is the spatially resolved model- reduces to 4 h for productively infected CD4+ T cells, which are based simulation that needs to be implemented in order to motile. For an HIV-specific T-cell frequency of 1%, we computed quantify the killing rate coefficient of the ODE model. that the inhibition of CTL locomotion by two-fold would Moving from phenomenological models of spatiotemporal reduce the probability of detection of infected target cells within dynamics of immune processes (e.g., the compartmental 24 h post-infection from 0.84 to 0.42. Thus, the requirements models, CAMs) to a physics-based description of immune cell for a prophylactic vaccine for seronegative individuals and migration in complex multicellular tissue environments presents an immunotherapeutic intervention of already HIV-infected a challenge to mathematical immunology. Advances in the direct individuals may differ significantly and are influenced by the visualization of antigen-specific T-cell mobility during their state of the lymphatic tissue structure. search for and their interaction with antigen-presenting cells Understanding the spatiotemporal dynamics of immune cells within LNs set the basis for diverse modeling approaches (7, 10, globally in the lymphatic system and locally in LNs is considered 11), which have been so far based on ad hoc postulated rules to be a prerequisite for the development of novel immune of cell behaviors. Our study gives a biophysics perspective on interventions in the context of HIV cure strategies (15, 56). To coordinated cell motility in lymphoid tissues, thus extending the this end, mathematical tools are being increasingly applied to Frontiers in Immunology | www.frontiersin.org 12 June 2019 | Volume 10 | Article 1213 Grebennikov et al. T Cell Motility in HIV Infection Control predict the impacts of trafficking and motility parameters on applicable generic mathematical tool for linking individual and the efficiency of immune surveillance in health and disease. For coordinated cell behaviors that can be used for in silico studies example, an optimal surveillance strategy for T cells was analyzed to embed immune processes into their spatial context. The by compartmental modeling of their systemic recirculation and physics-based computational model of multicellular dynamics LN transit times using a multicompartmental consideration of the immune response in lymphoid tissues provides a solid (66). The protective effect of increased CD4+ T-cell trafficking module that can be universally used in systems immunology on the dynamics of HIV infection has been recently shown studies (2, 6) for the benefit of patients suffering from chronic using another compartmental model (67), thus providing a basis virus diseases. for considering cell trafficking as an adjunct therapy option. A multiscale model of Mycobacterium tuberculosis infection AUTHOR CONTRIBUTIONS including an agent-based description of the cellular movement in a two-dimensional simulation grid representing the granuloma DG, AM, VV, and GB conceived and designed the study. DG, VV, was developed and calibrated using non-human primates to AM, and GB wrote the paper. DG, AB, and NB performed the derive the prediction of parameters underlying granuloma computational implementation of the model. DG performed the sterilization (8). However, such modeling attempts are still rare. model calibration. In conclusion, the large number of existing mathematical models based on low-resolution descriptions of immune FUNDING functions has to be further extended and embedded into physiologically distinct compartments and 3D morphological This work is supported by the Russian Science Foundation constraints inherent to cells, tissues, and the whole organism. (grant 18-11-00171). VV was supported by a grant from the This will then allow the research community not only to RUDN University Program 5-100. AM is also supported get a better quantitative understanding of immune system by a grant from the Spanish Ministry of Economy, functioning in infections such as HIV but also enable to Industry and Competitiveness and FEDER grant no. build integrative models for antiviral and immunomodulatory SAF2016-75505-R (AEI/MINEICO/FEDER, UE) and the drugs of various physical and chemical nature as well as the María de Maeztu Programme for Units of Excellence in effects of adoptive cell transfer therapies. We believe that a R&D (MDM-2014-0370). comprehensive approach to combination therapies based on ART and immunomodulatory drugs affecting a range of processes, SUPPLEMENTARY MATERIAL including LN fibrosis, the exhaustion of CTLs, and T-cell motility, should rely on formulation and implementation of The Supplementary Material for this article can be found hybrid spatially resolved multiscale mathematical models of virus online at: https://www.frontiersin.org/articles/10.3389/fimmu. infections (8, 9, 68). The here-developed model offers a broadly 2019.01213/full#supplementary-material REFERENCES 8. Pienaar E, Matern WM, Linderman JJ, Bader JS, Kirschner DE. Multi-scale model of Mycobacterium tuberculosis infection maps metabolite and gene 1. Munoz MA, Biro M, Weninger W. T cell migration in intact lymph nodes perturbations to granuloma sterilization predictions. Infect Immun. (2016) in vivo. Curr Opin Cell Biol. (2014) 30:17–24. doi: 10.1016/j.ceb.2014.05.002 84:1650–69. doi: 10.1128/IAI.01438-15 2. Wong HS, Germain RN. Robust control of the adaptive immune system. 9. Bouchnita A, Bocharov G, Meyerhans A, Volpert V. Hybrid approach to Semin Immunol. (2017) 36:17–27. doi: 10.1016/j.smim.2017.12.009 model the spatial regulation of T cell responses. BMC Immunol. (2017) 18:29. 3. Ludewig B, Stein JV, Sharpe J, Cervantes-Barragan L, Thiel V, Bocharov G. A doi: 10.1186/s12865-017-0205-0 global “imaging” view on systems approaches in immunology. Eur J Immunol. 10. Beauchemin C, Dixit NM, Perelson AS. Characterizing T cell movement (2012) 42:3116–25. doi: 10.1002/eji.201242508 within lymph nodes in the absence of antigen. J Immunol. (2007) 178:5505–12. 4. Worbs T, Förster R. T cell migration dynamics within lymph nodes during doi: 10.4049/jimmunol.178.9.5505 steady state: an overview of extracellular and intracellular factors influencing 11. Beltman JB, Marée AFM, Lynch JN, Miller MJ. Boer RJ de. Lymph node the basal intranodal T cell motility. In: Dustin M, McGavern D, editors. topology dictates T cell migration behavior. J Exp Med. (2007) 204:771–80. Visualizing Immunity Current Topics in Microbiology and Immunology. Berlin; doi: 10.1084/jem.20061278 Heidelberg: Springer (2009). p. 71–105. doi: 10.1007/978-3-540-93864-4_4 12. Drasdo D, Hoehme S, Block M. On the role of physics in the 5. Mandl JN, Liou R, Klauschen F, Vrisekoop N, Monteiro JP, Yates AJ, et al. growth and pattern formation of multi-cellular systems: what can we Quantification of lymph node transit times reveals differences in antigen learn from individual-cell based models? J Stat Phys. (2007) 128:287. surveillance strategies of naïve CD4+ and CD8+ T cells. Proc Natl Acad Sci doi: 10.1007/s10955-007-9289-x USA. (2012) 109:18036–41. doi: 10.1073/pnas.1211717109 13. Risler T. Cytoskeleton and cell motility. In: Meyers RA, editor. Encyclopedia of 6. Germain RN. Will systems biology deliver its promise and contribute to Complexity and Systems Science. New York, NY: Springer (2013). p. 1–53. the development of new or improved vaccines? What really constitutes doi: 10.1007/978-3-642-27737-5_112-3 the study of “systems biology” and how might such an approach 14. Camley BA, Rappel W-J. Physical models of collective cell motility: from facilitate vaccine design. Cold Spring Harb Perspect Biol. (2018) 10:a033308. cell to tissue. J Phys Appl Phys. (2017) 50:113002. doi: 10.1088/1361-6463/ doi: 10.1101/cshperspect.a033308 aa56fe 7. Gong C, Mattila JT, Miller M, Flynn JL, Linderman JJ, Kirschner D. Predicting 15. Krummel MF, Bartumeus F, Gérard A. T cell migration, search strategies lymph node output efficiency using systems biology. J Theor Biol. (2013) and mechanisms. Nat Rev Immunol. (2016) 16:193–201. doi: 10.1038/nri. 335:169–84. doi: 10.1016/j.jtbi.2013.06.016 2015.16 Frontiers in Immunology | www.frontiersin.org 13 June 2019 | Volume 10 | Article 1213 Grebennikov et al. T Cell Motility in HIV Infection Control 16. Gadhamsetty S, Marée AFM, Boer DJR, Beltman JB. Tissue dimensionality 39. Basu R, Whitlock BM, Husson J, Le Floc’h A, Jin W, Oyler-Yaniv A, et al. influences the functional response of cytotoxic T lymphocyte-mediated killing Cytotoxic T cells use mechanical force to potentiate target cell killing. Cell. of targets. Front Immunol. (2017) 7:668. doi: 10.3389/fimmu.2016.00668 (2016) 165:100–10. doi: 10.1016/j.cell.2016.01.021 17. Foley MH, Forcier T, McAndrew E, Gonzalez M, Chen H, Juelg B, et al. High 40. Purcell EM. Life at low Reynolds number. Am J Phys. (1977) 45:3–11. avidity CD8+ T cells efficiently eliminate motile HIV-infected targets and doi: 10.1119/1.10903 execute a locally focused program of anti-viral function. PLoS ONE. (2014) 41. Friedl P, Bröcker E-B. The biology of cell locomotion within three- 9:e87873. doi: 10.1371/journal.pone.0087873 dimensional extracellular matrix. Cell Mol Life Sci. (2000) 57:41–64. 18. Murooka TT, Deruaz M, Marangoni F, Vrbanac VD, Seung E, von Andrian doi: 10.1007/s000180050498 UH, et al. HIV-infected T cells are migratory vehicles for viral dissemination. 42. Hons M, Kopf A, Hauschild R, Leithner A, Gaertner F, Abe J, et al. Nature. (2012) 490:283–7. doi: 10.1038/nature11398 Chemokines and integrins independently tune actin flow and substrate 19. Fukazawa Y, Park H, Cameron MJ, Lefebvre F, Lum R, Coombes N, et al. friction during intranodal migration of T cells. Nat Immunol. (2018) 19:606– Lymph node T cell responses predict the efficacy of live attenuated SIV 16. doi: 10.1038/s41590-018-0109-z vaccines. Nat Med. (2012) 18:1673–81. doi: 10.1038/nm.2934 43. Goya GF, Marcos-Campos I, Fernández-Pacheco R, Sáez B, Godino J, Asín L, 20. Mohammadi P, Desfarges S, Bartha I, Joos B, Zangger N, Muñoz M, et al. 24 et al. Dendritic cell uptake of iron-based magnetic nanoparticles. Cell Biol Int. hours in the life of HIV-1 in a T cell line. PLoS Pathog. (2013) 9:e1003161. (2008) 32:1001–5. doi: 10.1016/j.cellbi.2008.04.001 doi: 10.1371/journal.ppat.1003161 44. Han X, Kloeden PE. Random Ordinary Differential Equations 21. O’Neill ME. PCG: A Family of Simple Fast Space-Efficient Statistically Good and Their Numerical Solution. Singapore: Springer (2017). Algorithms for Random Number Generation. Claremont, CA: Harvey Mudd doi: 10.1007/978-981-10-6265-0 College (2014). 45. Rothoeft T, Balkow S, Krummen M, Beissert S, Varga G, Loser K, et al. 22. Read MN, Bailey J, Timmis J, Chtanova T. Leukocyte motility models assessed Structure and duration of contact between dendritic cells and T cells are through simulation and multi-objective optimization-based model selection. controlled by T cell activation state. Eur J Immunol. (2006) 36:3105–17. PLoS Comput Biol. (2016) 12:e1005082. doi: 10.1371/journal.pcbi.1005082 doi: 10.1002/eji.200636145 23. Zimmermann J, Camley BA, Rappel W-J, Levine H. Contact inhibition of 46. Mempel TR, Henrickson SE. Andrian UH von. T-cell priming by dendritic locomotion determines cell–cell and cell–substrate forces in tissues. Proc Natl cells in lymph nodes occurs in three distinct phases. Nature. (2004) 427:154. Acad Sci USA. (2016) 113:2660–5. doi: 10.1073/pnas.1522330113 doi: 10.1038/nature02238 24. Kitano M, Yamazaki C, Takumi A, Ikeno T, Hemmi H, Takahashi 47. Miller MJ, Hejazi AS, Wei SH, Cahalan MD, Parker I. T cell repertoire N, et al. Imaging of the cross-presenting dendritic cell subsets in the scanning is promoted by dynamic dendritic cell behavior and random T cell skin-draining lymph node. Proc Natl Acad Sci USA. (2016) 113:1044–9. motility in the lymph node. Proc Natl Acad Sci USA. (2004) 101:998–1003. doi: 10.1073/pnas.1513607113 doi: 10.1073/pnas.0306407101 25. Hakim V, Silberzan P. Collective cell migration: a physics perspective. Rep Prog 48. Blattman JN, Antia R, Sourdive DJD, Wang X, Kaech SM, Murali-Krishna Phys. (2017) 80:076601. doi: 10.1088/1361-6633/aa65ef K, et al. Estimating the precursor frequency of naive antigen-specific CD8 T 26. Basan M, Prost J, Joanny J-F, Elgeti J. Dissipative particle dynamics Cells. J Exp Med. (2002) 195:657–64. doi: 10.1084/jem.20001021 simulations for biological tissues: rheology and competition. Phys Biol. (2011) 49. Bastidas S, Graw F, Smith MZ, Kuster H, Günthard HF, Oxenius A. CD8+ 8:026014. doi: 10.1088/1478-3975/8/2/026014 T cells are activated in an antigen-independent manner in HIV-infected 27. Liu MB, Liu GR, Zhou LW, Chang JZ. Dissipative particle dynamics (DPD): individuals. J Immunol. (2014) 192:1732–44. doi: 10.4049/jimmunol.1302027 an overview and recent developments. Arch Comput Methods Eng. (2015) 50. Haynes BF, Shaw GM, Korber B, Kelsoe G, Sodroski J, Hahn BH, et al. HIV– 22:529–56. doi: 10.1007/s11831-014-9124-x host interactions: implications for vaccine design. Cell Host Microbe. (2016) 28. Bell GI. Models for the specific adhesion of cells to cells. Science. (1978) 19:292–303. doi: 10.1016/j.chom.2016.02.002 200:618–27. doi: 10.1126/science.347575 51. Grossman Z, Feinberg MB, Paul WE. Multiple modes of cellular activation 29. Kong F, García AJ, Mould AP, Humphries MJ, Zhu C. Demonstration of catch and virus transmission in HIV infection: a role for chronically and latently bonds between an integrin and its ligand. J Cell Biol. (2009) 185:1275–84. infected cells in sustaining viral replication. Proc Natl Acad Sci USA. (1998) doi: 10.1083/jcb.200810002 95:6314–9. doi: 10.1073/pnas.95.11.6314 30. Basu R, Huse M. Mechanical communication at the immunological synapse. 52. Grossman Z, Meier-Schellersheim M, Paul WE, Picker LJ. Pathogenesis of Trends Cell Biol. (2017) 27:241–54. doi: 10.1016/j.tcb.2016.10.005 HIV infection: what the virus spares is as important as what it destroys. Nat 31. Segel GB, Cokelet GR, Lichtman MA. The measurement of lymphocyte Med. (2006) 12:289–95. doi: 10.1038/nm1380 volume: importance of reference particle deformability and counting solution 53. Cheynier R, Henrichwark S, Hadida F, Pelletier E, Oksenhendler E, Autran tonicity. Blood. (1981) 57:894–9. B, et al. HIV and T cell expansion in splenic white pulps is accompanied by 32. Wiśniewski JR, Hein MY, Cox J, Mann M. A “proteomic ruler” for protein infiltration of HIV-specific cytotoxic T lymphocytes. Cell. (1994) 78:373–87. copy number and concentration estimation without spike-in standards. Mol doi: 10.1016/0092-8674(94)90417-0 Cell Proteomics. (2014) 13:3497–506. doi: 10.1074/mcp.M113.037309 54. Haase AT. Population biology of HIV-1 infection: viral and CD4+ T cell 33. Illmer P, Erlebach C, Schinner F. A practicable and accurate method to demographics and dynamics in lymphatic tissues. Annu Rev Immunol. (1999) differentiate between intra- and extracellular water of microbial cells. FEMS 17:625–56. doi: 10.1146/annurev.immunol.17.1.625 Microbiol Lett. (1999) 178:135–9. doi: 10.1111/j.1574-6968.1999.tb13769.x 55. Kityo C, Makamdop KN, Rothenberger M, Chipman JG, Hoskuldsson T, 34. Brown GC. Total cell protein concentration as an evolutionary constraint on Beilman GJ, et al. Lymphoid tissue fibrosis is associated with impaired vaccine the metabolic control distribution in cells. J Theor Biol. (1991) 153:195–203. responses. J Clin Invest. (2018) 128:2763–73. doi: 10.1172/JCI97377 doi: 10.1016/S0022-5193(05)80422-9 56. Dimopoulos Y, Moysi E, Petrovas C. The lymph node in HIV pathogenesis. 35. Turgeon ML. Clinical Hematology: Theory and Procedures. Philadelphia, PA: Curr HIV/AIDS Rep. (2017) 14:133–40. doi: 10.1007/s11904-017- Lippincott Williams & Wilkins (2005). 0359-7 36. Lim TS, Goh JKH, Mortellaro A, Lim CT, Hämmerling GJ, Ricciardi- 57. Fackler OT, Murooka TT, Imle A, Mempel TR. Adding new dimensions: Castagnoli P. CD80 and CD86 differentially regulate mechanical interactions towards an integrative understanding of HIV-1 spread. Nat Rev Microbiol. of T-cells with antigen-presenting dendritic cells and B-cells. PLoS ONE. (2014) 12:563–74. doi: 10.1038/nrmicro3309 (2012) 7:e45185. doi: 10.1371/journal.pone.0045185 58. Vroomans RMA, Marée AFM, de Boer RJ, Beltman JB. Chemotactic migration 37. Miller MJ, Wei SH, Cahalan MD, Parker I. Autonomous T cell trafficking of T cells towards dendritic cells promotes the detection of rare antigens. PLoS examined in vivo with intravital two-photon microscopy. Proc Natl Acad Sci Comput Biol. (2012) 8:e1002763. doi: 10.1371/journal.pcbi.1002763 USA. (2003) 100:2604–9. doi: 10.1073/pnas.2628040100 59. Kadolsky UD, Yates AJ. How is the effectiveness of immune surveillance 38. Zhang ML, Guo AX, VandenBussche CJ. Morphologists overestimate impacted by the spatial distribution of spreading infections? Philos the nuclear-to-cytoplasmic ratio. Cancer Cytopathol. (2016) 124:669–77. Trans R Soc B Biol Sci. (2015) 370:20140289. doi: 10.1098/rstb. doi: 10.1002/cncy.21735 2014.0289 Frontiers in Immunology | www.frontiersin.org 14 June 2019 | Volume 10 | Article 1213 Grebennikov et al. T Cell Motility in HIV Infection Control 60. Textor J, Henrickson SE, Mandl JN, von Andrian UH, Westermann J, de 66. Lee M, Mandl JN, Germain RN, Yates AJ. The race for the prize: T cell Boer RJ, et al. Random migration and signal integration promote rapid trafficking strategies for optimal surveillance. Blood. (2012) 120:1432–8. and robust T cell recruitment. PLoS Comput Biol. (2014) 10:e1003752. doi: 10.1182/blood-2012-04-424655 doi: 10.1371/journal.pcbi.1003752 67. Fryer HR, Wolinsky SM, McLean AR. Increased T cell trafficking as 61. Banigan EJ, Harris TH, Christian DA, Hunter CA, Liu AJ. Heterogeneous adjunct therapy for HIV-1. PLoS Comput. Biol. (2018) 14:e1006028. CD8+ T cell migration in the lymph node in the absence of inflammation doi: 10.1371/journal.pcbi.1006028 revealed by quantitative migration analysis. PLoS Comput Biol. (2015) 68. Bocharov G, Volpert V, Ludewig B, Meyerhans A. Mathematical Immunology 11:e1004058. doi: 10.1371/journal.pcbi.1004058 of Virus Infections. Cham: Springer (2018). doi: 10.1007/978-3-319-72317-4 62. Scianna M, Preziosi L. Cellular Potts Models: Multiscale Extensions and Biological Applications. New York, NY: CRC Press (2013). doi: 10.1201/b14075 Conflict of Interest Statement: The authors declare that the research was 63. Anderson ARA, Chaplain MAJ, Rejniak KA (eds.). Single-Cell-Based conducted in the absence of any commercial or financial relationships that could Models in Biology and Medicine. Basel: Birkhäuser Basel (2007). be construed as a potential conflict of interest. doi: 10.1007/978-3-7643-8123-3 64. Haynes BF, Burton DR. Developing an HIV vaccine. Science. (2017) 355:1129– Copyright © 2019 Grebennikov, Bouchnita, Volpert, Bessonov, Meyerhans and 30. doi: 10.1126/science.aan0662 Bocharov. This is an open-access article distributed under the terms of the Creative 65. Connick E, Folkvord JM, Lind KT, Rakasz EG, Miles B, Wilson NA, et al. Commons Attribution License (CC BY). The use, distribution or reproduction in Compartmentalization of simian immunodeficiency virus replication within other forums is permitted, provided the original author(s) and the copyright owner(s) secondary lymphoid tissues of rhesus macaques is linked to disease stage are credited and that the original publication in this journal is cited, in accordance and inversely related to localization of virus-specific CTL. J Immunol. (2014) with accepted academic practice. No use, distribution or reproduction is permitted 193:5613–25. doi: 10.4049/jimmunol.1401161 which does not comply with these terms. Frontiers in Immunology | www.frontiersin.org 15 June 2019 | Volume 10 | Article 1213