The Hypervariable Loops of Free TCRs Sample Multiple Distinct Metastable Conformations in Solution

Abstract

CD4⁺ and CD8⁺ αβ T cell antigen recognition is determined by the interaction between the TCR Complementarity Determining Region (CDR) loops and the peptide-presenting MHC complex. These T cells are known for their ability to recognize multiple pMHC complexes, and for a necessary promiscuity that is required for their selection and function in the periphery. Crystallographic studies have previously elucidated the role of structural interactions in TCR engagement, but our understanding of the dynamic process that occurs during TCR binding is limited. To better understand the dynamic states that exist for TCR CDR loops in solution, and how this relates to their states when in complex with pMHC, we simulated the 2C T cell receptor in solution using all-atom molecular dynamics in explicit water and constructed a Markov State Model for each of the CDR3α and CDR3β loops. These models reveal multiple metastable states for the CDR3 loops in solution. Simulation data and metastable states reproduce known CDR3β crystal conformations, and reveal several novel conformations suggesting that CDR3β bound states are the result of search processes from nearby pre-existing equilibrium conformational states. Similar simulations of the invariant, Type I Natural Killer T cell receptor NKT15, which engages the monomorphic, MHC-like CD1d ligand, demonstrate that iNKT TCRs also have distinct states, but comparatively restricted degrees of motion.

Keywords: T cell receptor dynamics; adaptive immunity; independent component (IC) analysis; innate immunity; markov state models; molecular dynamics (MD).

Figure 1

Crystal structures of the 2C and NKT15 TCRs highlight their hypothesized dynamics and the starting points used for simulations. (A) Variable domains (gray) of 2C shown from the perspective of the pMHC surface. CDR3 loops shown for the unbound structure (cyan, PDB: 1TCR), bound to H-2K^b/SIYR (red, PDB: 1G6R), H-2K^b/QL9 (green, PDB:2CKB), and H-2L^d/QL9 (blue, PDB:2OI9). (B) Similar view of NKT15 variable domains (gray) with CDR loops shown for the unbound structure (cyan, PDB:2EYS), and bound to CD1d presenting αGalCer (orange, PDB: 3HUJ) or C20:2 (pink, PDB: 3VWJ).

Figure 2

The first 8 tICs of the tICA-decomposed dihedral angles show relatively simple motions of the CDR3α loop, and more complex behavior for the CDR3β loop for the 2C TCR. (A,B) 2-D Kernel Density Estimate of the simulation data projected onto the first two degrees of freedom discovered by tICA for the CDR3α and CDR3β loops, respectively, using a 2-D Gaussian kernel. The KDEs estimate the probability density function for finding a randomly selected frame in a region of conformational space described by the tICA degrees of freedom. (C,D) 1-D probability density graphs of the first eight tICA degrees of freedom for CDR3α and CDR3β, respectively, using a Gaussian kernel.

Figure 3

The states identified by the k-mediods clustering of the tICA decomposition provide insight in to the complex motions that the 2C CDR3β loop displays. (A) Ball and stick model of the CDR3β loop showing the centroids of the four macrostate clusters determined by the MSM. Centroids were determined by finding the orientation that minimizes the distance to all other members of the cluster under the tICA projection distance. (B) Equilibrium populations of the four clusters, determined by eigenvalue analysis of the macrostate MSM. (C) Projection of the centroids onto the first two tICA dimensions overlaid on the kernel density estimates of the projected data. (D) ϕ/ψ backbone angles of eight residues along the CDR3β loop. Colors are consistent throughout for state 1 (green), 2 (light blue), 3 (purple), and 4 (dark blue).

Figure 4

Macrostate Markov State Model of the 2C CDR3β, showing the complex dynamic behavior the loop adopts in solution. The Markov model was built with a time-lag of 115 ns, distinct from the 5 ns time-lag used for the tICA decompositions. State clusters are represented by their centroids as initially described in Figure 3A, and jump probabilities are described by arrows labeled by the probability of that state transition occurring in a 115 ns time step. Arrow size is proportional to jump probability.

Figure 5

tICA decomposition and k-mediods clustering of the resulting data highlight the significant, but simple, motion of the NKT15 TCR CDR loops. (A) Probability distributions of the NKT15 CDR3α and CDR3β conformations projected onto each of the first eight tICA degrees of freedom, computed by a 1-D kernel density estimate with a Gaussian kernel. (B) 2-D probability distribution of NKT15 CDR3α projected on the first two tICA degrees of freedom; selected conformations from the simulation are shown in orange and gold and overlaid on the probability distribution plot. (C) As in panel (D) for the CDR3β loop with selected conformations shown in light green and ochre. Probability distributions were computed by a 2-D kernel density estimate with a Gaussian kernel over all collected trajectory data.

Figure 6

Inter-loop hydrogen bonding plays a significant role in the simulations of the 2C TCR, comprising one of the 4 stable wells of the CDR3β tICA decomposition. (A) Structural rendering of the 2C hydrogen bond interaction between Ser93 and Gly207 with the Vα domain shown on the left and the Vβ domain shown on the right. Surrounding hydrophobic residues are shown in pink. (B) Projection of the data frame onto the tICA projections of the CDR3α (left) and CDR3β (right) loops overlaid in cyan.

Figure 7

Simulations of the 2C TCR CDR3β loop recapitulate bound crystal structure conformations and support the conformational melding hypothesis based on their location in the tICA projections. (A) tICA projections of the bound 2C CDR3β loop conformations for 2C bound to H-2K^b/SIYR (red), H-2K^b/dEV8 (green), and H-2L^d/QL9 (blue) overlaid on the 2-D probability density. (B) Ball and stick render of the CDR3β bound crystal structures overlaid with the nearest simulation frame by RMSD of the Cα's after aligning the β variable domain. Simulation data is shown in cyan.