Cognate antigen engagement on parenchymal cells stimulates CD8+ T cell proliferation in situ

T-cell responses are initiated upon cognate presentation by professional antigen presenting cells in lymphoid tissue. T cells then migrate to inflamed tissues, but further T-cell stimulation in these parenchymal target sites is not well understood. Here we show that T-cell expansion within inflamed tissues is a distinct phase that is neither a classical primary nor classical secondary response. This response, which we term ‘the mezzanine response', commences within days after initial antigen encounter, unlike the secondary response that usually occurs weeks after priming. A further distinction of this response is that T-cell proliferation is driven by parenchymal cell antigen presentation, without requiring professional antigen presenting cells, but with increased dependence on IL-2. The mezzanine response might, therefore, be a new target for inhibiting T-cell responses in allograft rejection and autoimmunity or for enhancing T-cell responses in the context of microbial or tumour immunity.


Supplementary Note 1 Equivalence to the longitudinal Lipkin-Meshkov-Glick model
We study quantum signatures of synchronization for the many-body Hamiltonian ] .
For L = 2 oscillators, this Hamiltonian describes the Lipkin-Meshkov-Glick (LMG) model with a longitudinal external field as we show in the following. The LMG model was originally introduced to analyze shape transitions of atomic nuclei [1]. By now, it is a standard model in statistical physics for the study of quantum phase transitions and the quantum-to-classical correspondence [2]. The LMG model describes N indistinguishable particles with 2 internal states, e.g. N hadronic spins or N two-level atoms. The particles interact symmetrically with each other and with an external magnetic field h. The Hamiltonian in its most general form is then given byĤ TheŜ x,y,z are the collective spin operators, which form an angular momentum operator with quantum number N/2, where ϵ abc is the totally anti-symmetric tensor.
Traditionally, the case of a purely transversal field (h x,y = 0) has received the most interest. However, a longitudinal field arises naturally in many quantum-optical realizations of the LMG model [3][4][5]. We consider the case that λ x and λ y have opposite signs and rotate the coordinate system according tô with tan(θ) = √ −λ y /λ x . The LMG Hamiltonian then reads up to a constant. Tuning the external field such that and introducing the mapping the LMG Hamiltonian reads Now we can use the fact that the collective spin operators can be represented in terms of bosonic creation and annihilation operators aŝ Using this mapping and the relation (3), it is then easy to see that the LMG Hamiltonian equals the Hamiltonian (2) with L = 2 oscillators up to a constant

Synchronization for cold atoms in tilted optical lattices
We demonstrate how the quantum synchronization Hamiltonian arises for cold bosonic atoms in tilted optical lattices. The quantum dynamics of the atoms in the lattice potential is described by the second-quantized Hamiltonian whereψ(x, t) is the bosonic field operator and we use the s-wave scattering approximation [6]. We consider a tilted or accelerated one-dimensional optical lattice, such that the singleparticle Hamiltonian readsĤ with where d is the lattice period. The one-dimensional limit can be realized by a tight confinement of the atoms in the radial directions, while the external field F can be realized by accelerating the entire lattice (see, e.g. [7]), by gravity or by magnetic gradient fields (see, e.g. [8,9]). Alternatively, such a system is realized for the propagation of light beams in modulated photonic lattices, where the index of refraction is periodically modulated [10,11].
The dynamics can be understood by expanding the field operators into the eigenstates of the single-particle HamiltonianĤ 1 , the so-called Wannier-Stark resonance states [12,13].
The linear system has a fundamental symmetry, it is invariant under a combined translation by the lattice period d and a shift of the energy by dF . The eigenstates are therefore organized in ladders, where Ψ is the wave function and E the energy of the eigenstates. Each ladder α = 0, 1, 2, . . . roughly corresponds to one Bloch band in the field-free case and n ∈ Z labels the rung of the ladder. The fundamental time scale of the system is determined by the energy distance of the rungs E α,n+1 − E α,n = dF . In the non-interacting limit this gives rise to a periodic motion, the Bloch oscillations [7,[12][13][14][15][16][17][18] with the frequency ω B = dF/ .
If the external field is not too strong, Landau-Zener tunneling between the bands can be neglected and the dynamics takes place in the ground ladder α = 0 only [19]. Plugging the expansion of the field operator into the Wannier-Stark resonance stateŝ into equation (12) then yields the Hamiltonian with the overlap integrals We note that this expansion is different from the one based on the Wannier functions which is commonly used in the analysis of the field-free case F = 0 [6]. Wannier functions are strongly localized in one well of the lattice such that nonlinear coupling terms between the different wells can be neglected. However, Wannier functions are no eigenstates such that there is a linear coupling term describing single-particle tunneling between adjacent wells. Wannier-Stark functions are much better suited for tilted or accelerated lattices: The periodicity of the single-particle dynamics becomes apparent and the effects of the nonlinearity can be analyzed rather directly [9,20,21].
The many-body-quantum dynamics described by the Hamiltonian (17) depends crucially on the overlap of the Wannier-Stark resonance states (18). A remarkable property of these states is their localization in real space, which increases with the lattice depth and the field strength. Therefore, several approximations can be made which greatly reduce the complexity of the interaction terms in the Hamiltonian (17). By definition, we have the permutational invariance χ k,ℓ;m,n = χ ℓ,k;m,n = χ k,ℓ;n,m , Furthermore, we note that due to the translational symmetry of the Wannier-Stark states. For a sufficiently deep lattice, only the on-site coupling and the coupling between nearest neighbor states are strong due to the strong localization of the Wannier-Stark states (cf. Supplementary Fig. 2 and [20]).
Using these approximation and neglecting all other interaction terms, the Hamiltonian (17) reduces toĤ Shifting the phase of each mode by a constant factor, b n = e inπ/2â n (22) we finally obtain the quantum synchronization Hamiltonian as discussed in our paper, with U = U 0 × χ 0,0;0,0 and K = 4U 0 × χ 0,0;0,1 . We note that the phase shift (22) has no implications for the physics described by the model. It has been introduced only for the sake of notational convenience, to allow for a better comparability with the established Kuramoto model.
In the paper we consider the dynamics of an initially pure Bose-Einstein condensate, i.e.
a product state where all bosonic atoms are in the same single-partice quantum state In particular, we consider a condensate which is initally at rest, i.e. has homogeneous phases, and is widely distributed in space such that c j = 1/ √ L for j = 1, . . . , L.
In experiments with ultracold atoms one commonly measures the momentum space density ρ(k, t) = ⟨ψ † (k)ψ(k)⟩ using the time-of-flight technique. Using the expansion (16) and the translational symmetry of the Wannier-Stark states the momentum density reads The long-time evolution of the quantum many-body system is thus determined by the coherences ⟨â † ℓâ n ⟩ t . As shown in our paper, synchronization can induce long-time coherence. To understand how this manifests in laboratory experiments, we briefly review the dynamics for K = 0 in the following.
In the linear case U = K = 0 the coherences evolve as which can be directly seen from integrating the Heisenberg equation, Hence, the momentum space density of the atoms is given by In this formulation one can directly see that the quantum dynamics is fully periodic, which is commonly referred to as Bloch oscillations [7,[12][13][14][15][16][17][18]. If t is an integer multiple of the For an initially pure BEC (24), the coherences of this state are simply given by ⟨b † ℓb n ⟩ 0 = N c * ℓ c n such that the momentum density (29) can be further simplified to where C denotes the discrete Fourier transform of the amplitudes c n : If the BEC is initially distributed widely distributed in space, then the Fourier transform | C| 2 is a 'comb' function, i.e. a periodic arrangement of sharp peaks. Hence we obtain a particular simple picture of Bloch oscillations in momentum space: A comb function is moving under the envelope given by |Ψ 0,0 (k)| 2 (cf. Supplementary Fig. 3).
In the regime of strong on-site interaction U and no coupling K = 0, all non-diagonal coherences tend to zero, ⟨b † ℓb n ⟩ → 0 for ℓ ̸ = n [20,[22][23][24], and the momentum density reads Revivals are possible in finite systems for longer time scales. The Wannier-Stark resonance state Ψ 0,0 (k) roughly extends over the first Brillouin zone k ∈ [−π/d, +π/d] (cf. Supplementary Fig. 3). Strong on-site interactions thus lead to a delocalization of ρ(k, t) over the entire first Brillouin zone. A localization in momentum space thus indicates the presence of persistent quantum correlations, which emerge due to classical synchronization, as discussed in the paper.

Supplementary Note 3 The mean-field limit
The mean-field equations of motion for the amplitudes c ℓ = ⟨â ℓ ⟩ are most easily derived from the Heisenberg equation with the HamiltonianĤ ] .
Throughout this paper we use scaled units with = such that we obtain Taking the expectation value on both sides of the equation then yields This equation includes the three point functions of the form ⟨â † jâ kâℓ ⟩ on the right-hand side. This is a general feature of quantum many-body system: The Heisenberg equations induce an infinite hierarchy of coupled equations for the expectation values. To obtain a closed set of equations we must truncate this hierarchy at some point by approximating the three-point functions in terms of one-point functions where the asterisk denotes complex conjugation. The truncation is generally valid for Bose-Einstein condensates with a high number of particles, as the error vanishes as 1/N . A discussion of the validity of this truncation and other approaches can be found in [25] and the references therein. We finally arrive at the mean-field equations of motion, For further analysis we decompose the complex numbers c n into amplitude and phase that is, Within the mean-field approximation, the squared modulus |c n | 2 corresponds to the number of excitations or atoms per mode, respectively: The equations of motion for the amplitudes and phases are then derived from equation (43), which leads to The first line shows that the amplitudes remains constant if all have the same value Denoting the total number of excitations as we thus find that the toric manifolds are invariant under the flow generated by the equations (43). On a torus with a given number of excitations N only the phases evolve according to equation where K nj = N K nj /L. This exactly constitutes the celebrated Kuramoto model [26][27][28].

Supplementary Note 4 Beyond mean-field
The dynamics of quantum fluctuations beyond mean-field and the depletion of the condensate mode can be analyzed in a Bogoliubov-de Gennes approach. To this end we decompose the bosonic annihilation operators into the condensate mode c n and the remaining quantum fluctuations, described by the operatorb n . Together, this gives the ansatz: a n = c n +b n .
Assuming an almost pure condensate with N atoms, c n is of order √ N while the fluctuationŝ b n are of order 1. In the mean-field limit one furthermore assumes that N → ∞, while the macroscopic interaction strength U N/L = g and NK jℓ /L = K jℓ tend to a constant value.
We can thus insert the ansatz (53) into the Heisenberg equations of motion (41) and sort the equation according to the order of the atom number N . To leading order N 1/2 one recovers the mean-field equations of motion (43). In next to leading order one obtains the Bogoliubov-de Gennes equations with the matrix elements On the Kuramoto manifold we have c n = √ N/Le −iϕn and the matrix elements simplify to whereT denotes time ordering. We are particularly interested in the case where the meanfield dynamics is periodic with period T . Then we can define a time-averaged Bogoliubov-de Gennes operatorL as Using this definition we obtain the same condition for dynamical instability as in the timeindependent case: The condensate is stable if all eigenvalues ofL are real, it becomes unstable if an eigenvalue becomes complex with Im(λ n ) > 0. We note that in this case the real parts of the eigenvalues are defined only modulo 2π/T . The Kuramoto model was originally introduced and solved for the case of many globally connected oscillators [26]. In this case we can analytically estimate the eigenvalues of the Bogoliubov-de Gennes operator and thus derive criteria for the growth of quantum fluctuations in the quantum Kuramoto model. So we consider the globally connected case K nj = K/L on the Kuramoto manifold in the limit of many coupled oscillators L → ∞ with a fixed density ν = N/L and define the Kuramoto phase order parameter In the generic case the oscillators relax to a steady state with a fixed value of r. The global phase γ is also constant in time if we transform to an appropriate co-rotating frame of reference, where Ω is given by the mean frequency of the oscillators, The mean-field equations of motion on the Kuramoto manifold the read and the coefficients of the Bogoliubov-de Gennes operator (56) can be simplified to We observe that in the limit L → ∞ all off-diagonal elements of η vanish as L −1 such that both ζ and η become diagonal. Thus the the Bogoliubov-de Gennes equation for the individual oscillators decouple, and we can analyze the stability of all oscillators separately. To evaluate the stability, we have to distinguish between phase-locked and drifting oscillators. All oscillators with |ω n + U ν − Ω| ≤ Kr are phase-locked, i.e. the phase ϕ n assumes a fixed value given by the equation ω n + U ν − Ω + Kr sin(γ − ϕ n ) = 0.
The mean-field dynamics is periodic such that we have to analyze the eigenvalues of the averaged operator L nn as discussed above. In this case the time-ordered integral in equation (60) is evaluated numerically to obtain the eigenvalues λ n,± . We observe a non-zero growth rate in the immediate vicinity of the phase-locking region while λ n,± becomes purely real for all other drifting oscillators (cf. Supplementary Fig. 4).
We stress that a dynamical instability is possible only if the phase order parameter r does The off-diagonal contributions average out such that the eigenvalues ofL nn are purely real.