計算生命科学の基礎９（北尾 彰朗）

分子シミュレーションで見る タンパク質や核酸の複合体形成と解離 2022年11月30日（水）神戸大学統合研究拠点２階セミナー室 「計算生命科学の基礎９」（第９回） 北尾 彰朗 東京工業大学生命理工学院

原子解像度で生体分子を観察する方法 すごく冷やして形を観る（低温電子顕微鏡） 規則的に並べて回折像を観る（X線、中性子線、電子線結晶解析） 局所的幾何学情報から構造構築（核磁気共鳴） 10？分子の平均 and/or 時間平均

Protein association ー function ー dissociation Association μs ≥ Diffusion Long-range interactions Function Dehydration Conformational change Complex formation kon 結合速度定数 Dissociation ms ≥ Catalysis Activation/Inhibition Signal transduction … Hydration Conformational change Kd , 解離定数 ∆G°, 標準結合自由エネルギー K i , K cat , koff 解離速度定数

Standard binding free energy, dissociation constant, and association and dissociation rate constants 𝑅𝑅 𝐿𝐿 kon → ← koff + d [ RL ] = kon [ R ][ L ] − koff [ RL ] dt 𝑅𝑅𝑅𝑅 = Kd [ R= ][ L ] [ RL ] koff kon , = ∆G° β −1 ln K d / C0

「分子動力学 Molecular Dynamics (MD)」 コンピュータで原子レベルから機能を観察 ニュートンの運動方程式 質量×加速度＝力 ma = F Takemura et al. Sci. Rep. 2017.

Requirements for molecular simulation 分子シミュレーションに必要なこと • Model for molecular systems (force field) • Simulation methods • Analysis and visualization methods

Force field (分子力場） for molecular dynamics • Typical force field E = Ebond + Eangle + Etorsion + ELJ + Eel ri Bond stretching energy = Ebond θi Bond angle bending energy Torsion (dihedral) energy Lennard-Jones energy Electrostatic energy Etorsion = ϕi Rij i = ELJ j Eel = ∑ ∑ ki θi −θi0  + bond i Eangle = 2 ki  ri − ri0  + anglei   2   Vi 1+ cos(nϕi + δ i ) + ∑  torsioni    σ 4εij  ij R atom pair ij  ij  ∑ qi q j ∑ atom pair ij ε Rij 12      σ −  ij  Rij       6     

Simulation setup Simulated system All atoms in a box Interaction Consider periodic and infinite boxes ※Water and ions are included but not shown for visualizability

Solving equations of motion • Leap-frog scheme 1 1 v(t + ∆t ) = v(t − ∆t ) + a (t )∆t 2 2 1 x(t + ∆= t ) x(t ) + v(t + ∆t )∆t 2 • Velocity Verlet ∆t 1 v(t + ) = v(t ) + a (t )∆t 2 2 Satisfy ‘Simplectic condition’ シンプレクティック条件満たす （Time reversible 時間反転ＯＫ） 1 x(t + ∆= t ) x(t ) + v(t )∆t + a (t )∆t 2 2 ∆t = x(t ) + v(t + )∆t 2 ∆t 1 v(t + ∆t= ) v(t + ) + a (t + ∆t )∆t 2 2

10 Vibration frequencies, Δt, and MD length • Constraining hydrogen related bonds 水素の結合長振動を拘束 1.0 [ fs ] ≤ ∆t ≤ 2.0 [ fs ] , fs =10−15 s Vibrational frequencies of proteins • No constraints 拘束なし 0.5 [ fs ] ≤ ∆t ≤ 1.0 [ fs ] ν = ν c = 3333 cm −1  ⇔ T = 1 ν  10 [ fs ] • MDで観察できる時間 ∆t =2.0 [ fs ]の場合 106 steps = 2.0 ×10−9 s = 2.0 ns 109 steps = 2.0 ×10−6 s = 2.0 µ s 1012 steps = 2.0 ×10−3 s = 2.0 ms N‒H bond stretching

Strategies for binding free energy calculation Direct method ΔGb ΔG1 Alchemic cycle ΔGb = ΔG1 − ΔG2 0 Double Annihilation Method (DAM), etc. ΔG2

Binding free calculation by association/dissociation simulation Association μs ≥ ΔGb Dissociation ms ≥ • Applicable to pairs of large molecules • Additional Information can be obtained. A/dissociation pathways Kinetics • Problem: Real a/dissociation time > Feasible MD time Biased simulations Unbiased simulations

Parallel Cascade Selection MD (PaCS-MD) Ryuhei Harada • Cycle of parallel independent all-atom MDs without bias • Selection of snapshots closet to a target value (Targeted) PaCS-MD R. Harada & A. Kitao, J. Chem. Phys. 139, 035103 (2013). Nontargeted PaCS-MD R. Harada, A. Kitao, J. Chem Theory Comput. 148 105101 (2018). nrep The number of replica PaCS-MD cycle ncycle

Procedure of PaCS-MD Multiple Independent ＭＤ Typically, 0.1 ̶ 2.0 ns for each MD Ranking and Selection of nrep Snapshots Short MD Short MD Short MD ・ ・ ・ × nrep ・ ・ nrep = 10 ̶ 100 Short MD Restart with Random Initial Velocities × nrep

Effects of the selection in PaCS-MD d d Quantity for selection Target Start preliminary MD P(r) select & restart 1st cycle P(r) Selection effectively enhances probability of rare event occurence. select & restart 2nd cycle select & restart P(r) Cycle select & restart Nth cycle P(r)

Effects of velocity reinitialization in PaCS-MD MD j in Cycle (n+1) Route change MD i in Cycle n Selected snapshot MD i starting point

Comparison with related unbiased samplings • Weighted Ensemble (WE) Splitting of trajectory Statistical weights are clearly defined. Predetermined order parameters to describe reaction coordinates Huber GA, Kim S. Biophys J., 70, 97 (1996). Zuckerman, D. M. Chong, L. T., Annu Rev Biophys, 46, 43 (2017). • Forward Flux Sampling (FFS) Reaction coordinate is predetermined. Allen RJ, Warren PB, ten Wolde PR. Phys Rev Lett. 94, 4 (2005). • PaCS-MD Quantity for selection is required but it is independent of reaction coordinates. Statistical weights are directly not defined but can be estimated with Markov state model.

Dissociation PaCS-MD (dPaCS-MD) Selection of snapshots with longer center of mass (COM) distances Ligands bound deep inside of GPCR pocket (A2A/T4E) Protein/intrinsically disordered peptide (MDM2/NTD-p53) Enzyme/inhibitor (barnase/barster) Protein/DNA (p53 core domain/DNA) Real time︓ a few hours Total sim. time︓101 ns dPaCS-MD/MSM successfully reproduced standard binding free energy in >20 cases.

Calculating binding affinities kon 結合速度定数 Kd , 解離定数 ∆G° 標準結合自由エネルギー koff 解離速度定数

Mt Fuji climbing routes

Conformational space sampled by targeted PaCS-MD Better statistics should be achieved by multiple trials (ntrial) No acceleration Other directions Conformational space Selected snapshots Acceleration Target direction

Markov State Model (MSM) Discretization of sampled conformational states into microstates →Calculation of transition matrix→ Stationary probabilities → Free energies, kinetics, etc … i pki pji pij k j … pjk Free energy of i –th state → ∆G, Kd MFPT (Mean fast passage time) → kon, koff MSM can be constructed from many short MD trajectories. J. D. Chodera, et al., J. Chem. Phys. 126, 155101 (2007). J. Prinz, et al., J.Chem. Pyhs. 134, 244108 (2011)

Binding free energy calculation by dPaCS-MD/MSM Dissociation d Standard binding free energy ∆𝑊𝑊 ∆𝐺𝐺 ∘ = −∆𝑊𝑊 + ∆𝐺𝐺𝑣𝑣 1 𝑉𝑉𝑏𝑏 ∆𝐺𝐺𝑣𝑣 = − ln 𝑜𝑜 𝛽𝛽 𝑉𝑉 DP Tran, K Takemura, K Kuwata, A Kitao, J Chem Theory Comput 2018 DP Tran, A Kitao, J. Phys. Chem. B 2019

Biased versus unbiased simulations • Biased: Steered MD • Unbiased : PaCS-MD Selecting snapshots with longer inter-center of mass distance HEW lysozyme/tri-NAG DP Tran, K Takemura, K Kuwata, A Kitao, J Chem Theory Comput 2018

Functions and structure of tumor suppressor p53 • p53 is a transcriptional factor that regulates the cell's response to a variety of stresses, including DNA damage, hypoxia, oncogene activation and other stress signals. 1 Intrinsically disordered TADs 67 Trans activation domains Pro 98 Intrinsically disordered DBD DNA binding domain 303 323 363 392 NLS I Tetra. D CTD Tetrameric domain C-terminal domain

Binding of p53-DBD with DNA • DBD binds to a specific sequence : RRRCWWGYYY (R=A/G, W=A/T, Y=C/T) • About a half of all human tumors contain the p53 mutations. 95% of them occur in DBD. 1 TADs 67 Trans activation domains Pro 98 DBD DNA binding domain 303 323 363 392 NLS I Tetra D CTD

Interactions between p53-DBD monomer and DNA Mohamed Marzouk Sobeh R283 R280 R248 K120 Major groove Minor g. Adapted with permission from Marzouk & Kitao, J Chem Info Modeling 2022. Copyright 2022 American Chemical Society.

Dissociation pathways Average simulation time for dissociation :11.2 ± 2.2 ns Total accumulated MD cost: 8.4 µs Primarily dissociation : +X direction Majority: −Y direction Minority: +Y direction Adapted with permission from Marzouk & Kitao, J Chem Info Modeling 2022. Copyright 2022 American Chemical Society.

Majority: −Y direction Major groove dis. ⇒ Minor groove dis. Minority: +Y direction Minor groove dis. ⇒ Major groove dis. Adapted with permission from Marzouk & Kitao, J Chem Info Modeling 2022. Copyright 2022 American Chemical Society.

Free energy landscape projected onto intermolecular position Free energy change along X Minor pathways Major pathways Why is the steeper direction major? Binding affinity of p53-DBD monomer ΔG° = − 10.9 ± 0.4 kcal/mol (This work) – 11.1 kcal/mol (He et al. PNAS, 2019 by ITC) Adapted with permission from Marzouk & Kitao, J Chem Info Modeling 2022. Copyright 2022 American Chemical Society.

Why is dissociation along +Y minor? Rate limiting transition Rate limiting step Larger free energy gap Low transition probability Smaller free energy gap higher transition probability Marzouk & Kitao, in preparation

Key residues in dissociation and p53 mutations in human cancers • R248 (minor groove) • R248 bridges two DNA strands by interacting 6 nucleotides including TTG in the middle of the consensus sequence. • R248Q and R248W are ranked in top 5 (7.9 %). • R248 was the last residue dissociated from DNA in 75 % of trials. • R280 (major groove) • • • • R280 is important for DNA recognition. R280K impairs DNA transcription. Two R280 mutations are below top 25 (0.8 %). R280 was the last residue dissociated from DNA in 20 % of trials. Adapted with permission from Marzouk & Kitao, J Chem Info Modeling 2022. Copyright 2022 American Chemical Society.

Roles of trans activation domains (TADs) in p53 1 Intrinsically disordered TADs 67 Trans activation domains Pro 98 Intrinsically disordered DBD DNA binding domain 転写活性化ドメイン (TAD1: 残基1-42, TAD2:43-62) 303 323 363 392 NLS I Tetra. D CTD Tetrameric domain C-terminal domain

MDM2/TAD-p53 complex Murine double minute clone 2 protein (MDM2) bound to the transactivation domain 1 of p53 protein (TAD-p53) • TAD-p53 is an intrinsically disordered region that can bind to different proteins. • MDM2 ubiquitinates p53 upon binding, interacts with p53 to block its binding with target DNA, and promotes the export of p53 out of the cell nucleus, preventing it from interacting with target DNA. • ∆Go, koff, and kon were measured experimentally.

Free energy landscape of TAD-p53 monomer Metadynamics Adapted with permission from Tran & Kitao, J Chem Theory Comput 2020. Copyright 2022 American Chemical Society.

Obtaining more information of association processes • dPaCS-MD typically generates “straight” dissociation pathways, which may not sufficiently sample encounter complexes. • If such transient conformations significantly contribute to kinetics, other sampling methods to better sample encounter complexes may be needed.

Association and dissociation simulation by a/dPaCS-MD a/dPaCS-MD Association PaCS-MD (Select shorter distance) Dissociation PaCS-MD (Select longer distance) Adapted with permission from Tran & Kitao, J Chem Theory Comput 2020. Copyright 2022 American Chemical Society.

Docking simulation by association and dissociation PaCS-MD/MSM Tran & Kitao, J Chem Theory Comput 2020 Duy Phuoc Tran Tran & Kitao, J Chem Theory Comput 2020.

Relaxation of obtained bound poses with 60 ns standard MD Adapted with permission from Tran & Kitao, J Chem Theory Comput 2020. Copyright 2022 American Chemical Society.

Relaxation of obtained bound poses with 60 ns standard MD Induced fit and dehydration Adapted with permission from Tran & Kitao, J Chem Theory Comput 2020. Copyright 2022 American Chemical Society.

Obtained minimum free energy structure and ∆Go / koff / kon Total accumulated MD cost: 11.7 µs Calculated ∆Go , koff, and kon agree with experimental values. Orange: Global free energy minimum Green: X-ray ∆Go koff kon (kcal/mol) (s-1) (107M-1s-1) a/dPaCS-MD/Cα-MSM -8.9 ± 0.2 1.1 ± 0.3 0.25 ± 0.01 dPaCS-MD/Cα-MSM -8.4±0.7 3.1±1.3 1.4±0.2 -9.1±0.1 2.06±0.09 0.92±0.004 -9.9±0.1 1.43±0.24 2.2±0.02 SPR -8.8±0.2 0.8±0.3 0.21±0.02 DFS -8.4±0.5 Method ITC and stopped-flow Adapted with permission from Tran & Kitao, J Chem Theory Comput 2020. Copyright 2022 American Chemical Society.

Analysis of association and dissociation pathways Free energy landscape Transition probabilities between states CI3 CI2 CI1 CN CI4 UB 5 complex intermediates (CI1-CI5) close to complex native (CN) in free energy CI5 Adapted with permission from Tran & Kitao, J Chem Theory Comput 2020. Copyright 2022 American Chemical Society.

Analysis of association and dissociation pathways Free energy surface Association pathways and fluxes UB Kinetic selection CI2 Induced fit CN Adapted with permission from Tran & Kitao, J Chem Theory Comput 2020. Copyright 2022 American Chemical Society.

Simulation of a/dissociation processes by PaCS-MD/MSM • MD length 0.1 ns MD is the best choice for efficient PaCS-MD. For MSM, 0.1 ns is the minimum length. • Selection In the cases we tried up to now (> 20), the use of intermolecular distance worked well for enhancing a/dissociations. • Computational time In the cases presented today, only 100 ~ 101 μs MD in accumulated MD time is sufficient. Preceding works for similar cases spent submilliseconds to milliseconds.

Computational requirements for calculating binding affinity with PaCS-MD/MSM • The result of binding free energy is robust, but kinetics depends on the choice of MSM reaction coordinate (RC). RC should be sufficiently finegrained to distinguish important microstates, depending on phenomena. • To calculate binding free energy only, 1D-MSM is sufficient. To construct 1D-MSM, nrep should be at least 30 ~ 50. Alternatively, more MD trajectories should be added later if PaCS-MD is conducted with less nrep. • To better calculate kinetic rates, high dimensional MSM is needed. To construct 3D or higher dimensional MSM, ntrial should be at least 30 ~ 50. • For better computational efficiency, conformational sampling of less sampled areas should be conducted with more sophisticated ways.

Thank you for your attention! Acknowledgements ～ March 2022 Kitao lab member Tran Phuoc Duy Mohamed Marzouk Sobeh Oct 2022 ～