STAR★Methods

Key resources table

REAGENT or RESOURCE | SOURCE | IDENTIFIER
Software and algorithm

GROMACS (v5.1.4) | Abraham et al.34 | https://www.gromacs.org/
SWISS-MODEL | Waterhouse et al.35 | https://swissmodel.expasy.org/
NAPS webserver | Chakrabarty et al.36,37 | http://bioinf.iiit.ac.in/NAPS/
PDB | PDB ID: 3WPC, 5ZLN | https://www.rcsb.org/
Uniprot | Q9NR96 | https://www.uniprot.org/
Figshare database | Initial and representative structures extracted in this study | https://doi.org/10.6084/m9.figshare.31643326

Method details

TLR9 complex preparation

Although a few structures of horse and mouse TLR9 are available, the structure of human TLR9 is not solved. Therefore, we first retrieved the sequence of TLR9 ectodomain (ECD) from UniProt database (Q9NR96). The sequence identity between horse TLR9 and human TLR9 ECD was 83.6%, whereas that of mouse TLR9 was 71.3%. Hence, we generated 3D models of TLR9 ECD based on the horse TLR9 structure (PDB ID: 3WPC) using SWISS-MODEL webserver35 as a homodimer (hereafter, TLR9 and TLR9∗ monomers) and this structure was considered TLR9_apo. The horse TLR9 structure is available in a complex with CpG DNA; however, no IM is present. Therefore, to construct both TLR9_CpG and TLR9_CpG_IM complexes, we used a mouse crystal structure (PDB ID: 5ZLN). The root-mean-square deviation (RMSD) between the TLR9_apo homodimer and 5ZLN is 0.55 Å (for 1316 Cα atoms). We superimposed the TLR9_apo homodimer on to the mouse TLR9 structure (PDB ID: 5ZLN) and extracted the CpG ssDNA (single stranded) (nucleotide sequence: 5-AGGCGTTTTT-3′) (2:2 stoichiometry).21,24 Moreover, to construct TLR9_CpG_IM, we superimposed the TLR9_CpG DNA complex onto the mouse structure and extracted a single stranded immunomodulator (IM) (5′-xCx nucleotide sequence: 5′-TCGC-3′) at 2:2:2 stoichiometry. Subsequently, the three complexes were analyzed and prepared for simulations. Previous reports have suggested that binding of CpG DNA and IM to TLR9 directly promotes dimerization and augments activation24 and therefore, two monomers of TLR9 were maintained as modeled. Although homology modeling introduces inherent uncertainty, the high sequence identity between human and horse TLR9 and the low RMSD relative to the mouse TLR9 structure support the structural reliability of the model for studying ectodomain-level dynamics.

Molecular dynamics simulations of apo-, CpG ssDNA-, and IM-bound TLR9 complexes

TLR9 is localized in the intracellular compartments of the endosomes after signaling, and an acidic pH (5–6.5) is crucial for its interaction with CpG ssDNA.38 Therefore, we analyzed and modified the protonation of charged residues at pH 5.5, particularly the histidine residues, using H++ webserver.39 Subsequently, all TLR9 complexes were subjected to molecular dynamics simulations using GROMACS 5.1.434 with the AMBER ff99SB-ILDN force field for both proteins and DNA. Furthermore, TLR9 complexes were placed in a dodecahedral box containing TIP3P water molecules. The periodic boundary distance was maintained at 14 Å, and the sufficient counter ions were added to neutralize the system. Simulations were performed in three steps. First, the complexes were subjected to energy minimization with a maximum tolerance of 1000 kJ/mol using the steepest descent gradient algorithm. Subsequently, the system was equilibrated using NVT and NPT ensemble for 500 ps and 1 ns, respectively. Finally, production simulations were performed using a 2 fs time step, integrated with leap-frog algorithm, for 300 ns with three independent replicates (3 × 300 ns). The simulation parameters are similar to our previous studies.40 In brief, particle-mesh Ewald (grid spacing of 0.16 nm) was used for fast Fourier transform and applied for long-range electrostatics with 1.0 cutoff for short-range and van der Waals interactions. The hydrogen bonds were constrained using Lincs algorithm. Temperature (300 K) and pressure (0.1 bar) were applied using modified Berendsen thermostat (v-rescale) and Parrinello-Rahman barostat. Finally, coordinate data were saved every 2 ps for trajectory analysis. Most trajectory analyses were performed using built-in GROMACS tools, such as RMSD, root-mean-square fluctuation (RMSF), hydrogen bonds, and density distributions. To calculate the interaction energy between components during the simulations, we used rerun module integrated into GROMACS. The top principal components (PCs) were extracted using gmx covar and gmx anaeig from a concatenated trajectory comprising the last 250ns from each replicate. The gmx covar calculates and parameterizes the matrix, and the eigenvectors were analyzed using gmx anaeig. The free-energy landscapes (FELs) of the top PCs were determined using gmx Sham module. A summary of the simulation parameters is given in Table S1.

Residue interaction networks

We constructed a non-weighted residue-residue interaction network using the NAPS web server to gain deeper insight into the intramolecular signal transduction mechanisms of the protein.37 The analysis was carried out for the lowest energy structures of apo-, ssDNA- and IM-bound complexes of TLR9 dimer. In brief, each amino acid residue is represented as a node and the edges present the inter-residue contacts. The Cα-Cα distance cutoff was set to 7 Å with non-weighted edges.41 The influence of side chains was accounted for by defining two coarse grained centers per residue: one representing the Cα backbone and the other representing the side chain heavy atoms distant from the Cα atom. The betweenness centrality (CB), a parameter that quantifies signal propagation within the protein, was computed to identify key residues acting as central nodes.36,42 The absolute value was calculated as the difference between the per residue CB values of apo and respective complex forms.