.. include:: substitutions.txt .. currentmodule:: jscatter Biomacromolecules (bio) ----------------------- Examples using the biomacromolecule module (bio) for protein and DNA scattering. Load PDB structure and calc scattering ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ Load pdb protein structure from the PDB data bank by ID to *scatteringUniverse*. The pdb file is corrected and hydrogen is added automatically. The protein structure including H is saved to 3rn3_h.pdb. .. literalinclude:: ../../examples/example_bio_loadPDB.py :language: python .. image:: ../../examples/images/biosolventmatching.jpg :align: center :width: 50 % :alt: universe formfactor Protein scattering and effective diffusion ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ Protein scattering for different molecular weight proteins at 10 g/l concentration. We assume here rigid protein structures without internal fluctuations. The diffusion corresponds accordingly to the rigid body diffusion of a protein. See `Exploring internal protein dynamics by neutron spin echo spectroscopy `_ for details. .. literalinclude:: ../../examples/example_bio_proteinformfactoranddiffusion.py :lines: 1- The coherent scattering reflects the shape of the protein at low Q (SANS region) while for :math:`Q>3 nm^{-1}` the internal structure is visible (typical backscattering or TOF instruments). For :math:`Q>>3 nm^{-1}` the coherent contribution levels between 10-20% of the incoherent scattering for all proteins. The incoherent is slightly dependent on protein amino acid composition but the coherent/incoherent ratio is independent on concentration. The relative D2O background depends on protein concentration. Above :math:`Q>>10 nm^{-1}` the usage of Babinets principle may be questionable and a different calculation method is needed. :math:`D_{coh}/D_0` reflects the shape and size of the protein like the I(Q) does. The incoherent diffusion equals the coherent at larger Q and :math:`D_{inc}/D_0` depends slightly on shape. Dependent on the used instrument coherent and incoherent diffusion have to be mixed according to the coherent and incoherent contributions to I(Q). E.g for NSE there is a Q where coherent and incoherent (adds negative) compensate and no signal remains. Often for incoherent measurements only incoherent is taken into account. .. image:: ../../examples/images/bio_protein_formfactor+Diffusion.png :align: center :height: 300px :alt: bio_protein_formfactor+Diffusion.png Normal mode relaxation as seen by NSE measuring the intermediate scattering function ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ Alcohol dehydrogenase (yeast) example for normal mode analysis in time domain (NSE) See `Direct Observation of Correlated Interdomain Motion in Alcohol Dehydrogenase `_ for corresponding measurements that show the dynamics. In this example we consecutively - Examine the protein - Create normal modes - Calc effective diffusion - Calc the dynamic mode formfactor from normal modes - Show how the intermediate scattering function (ISF) looks with diffusion and internal normal mode relaxation - Finally, we build from this a model that can be used for fitting including H(Q)/S(Q) correction. .. literalinclude:: ../../examples/example_bio_proteinNormalModesScattering.py :lines: 16- The ISF shows in the initial slope the combined Deff and relaxes for longer times > 30ns to the rigid protein Deff. How strong the change is depends on the mode amplitudes and the relaxation times of the modes. .. image:: ../../examples/images/ADHNM_Deff.jpg :align: center :height: 300px :alt: ADHNM_Deff.jpg .. image:: ../../examples/images/ADHNM_Iqt.jpg :align: center :height: 300px :alt: ADHNM_Iqt.jpg For the model we use the ISF of multiple modes (see :py:func:`~jscatter.bio.scatter.intScatFuncPMode`) and the Rayleigh expansion for diffusing rigid proteins/particles (:py:func:`~jscatter.bio.scatter.intScatFuncYlm`) : Translational and rotational diffusion are corrected for direct interparticle interactions described by the structure factor :math:`S(Q)` and hydrodynamic interactions within the hydrodynamic function :math:`H(Q)` as :math:`D_t(Q) = D_{t,0} H(Q)/S(Q)` and :math:`D_r = D_{r.0}H_r`. The intermediate scattering function :math:`F(Q)` assuming dynamic decoupling with translational/rotational and domain motions is .. math:: F(Q,t) = F_{t,r}(Q,t) * F_a(Q,t) - trans/rot diffusion contribution .. math:: F_{t,r}(Q,t) = e^{-D_{t}Q^2t}\sum_l S_l(Q)e^{-l(l+1)D_{r}t} - normal mode contributions .. math:: F_a(Q,t) &= \frac{F(Q) + \sum_{\alpha} e^{-\lambda_{\alpha} t} P_{\alpha}(Q)}{[F(Q) + \\ \sum_{\alpha} P_{\alpha}(Q)]} &= (1-A(Q)) + A(Q) e^{-\lambda t} with with common relaxation rate :math:`\lambda` and :math:`A(Q) = \frac{\sum_{\alpha} P_{\alpha}(Q)}{[F(Q) + \sum_{\alpha} P_{\alpha}(Q)]}` We can also use the Ornstein-Uhlenbeck like relaxation :py:func:`~jscatter.bio.scatter.intScatFuncOU` for :math:`F_a(Q,t)` that allows description within internal friction in the protein. The above description corresponds to a small displacement approximation of the Ornstein-Uhlenbeck process. Finally NSE data can be read including Q values like :ref:`IN15` or :ref:`JNSE`. Typically we measure 10-15 different Q values between :math:`0.025-0.2 nm^{-1}`. The read dataList should contain for each dataArray an attribute `q` with the scattering vector. For demonstration we just use simulated data here .. literalinclude:: ../../examples/example_bio_proteinNormalModesScattering_2.py :lines: 8- .. image:: ../../examples/images/ADHNM_SqHq.jpg :align: center :height: 300px :alt: ADHNM_Deff.jpg .. image:: ../../examples/images/ADHNM_IQTsim.jpg :align: center :height: 300px :alt: ADHNM_Iqtsim.jpg Load trajectory from MD simulation ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ A *scatteringUniverse* with a complete trajectory from MD simulation is created. The PSF needs atomic types to be guessed from names to identify atoms in the used format. You may need to install MDAnalysisTests to get the files.( `python -m pip install MDAnalysisTests`) It might be necessary to transform the box that the protein is not crossing boundaries of the universe box. .. literalinclude:: ../../examples/example_bio_loadTrajectory.py :language: python .. image:: ../../examples/images/uniformfactorstraj.jpg :align: center :width: 50 % :alt: uniformfactors Compare different resolution options for coarse graining ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ PDB structures without explicit solvent for small angle scattering. The example shows the validity of residue coarse graining up to around 3/nm. With coarse graining in cubes (cubesize) the approximation seems best. This might be useful to speed up computations that take long (e.g. ISF at low Q) There is basically no difference between precalculated and averaged residue formfactors and explicit calculated residue formfactors for each residue (uni.explicitResidueFormFactorAmpl = True) The explicit ones include better deuteration of specific atoms. Cα models loose some precision in volume respectively in forward scattering. C-alpha models need a .calphaCoarseGrainRadius = 0.342 nm because of the missing atoms. In addition, the mean residue position is not the C-alpha position. We use 0.342 nm as a good average to get same forward scattering over a bunch of proteins (see example_bio_proteinCoarseGrainRadius.py). .. literalinclude:: ../../examples/example_bio_comparecoarsegraining.py :language: python .. image:: ../../examples/images/uniformfactors.jpg :align: center :width: 50 % :alt: uniformfactors Selective partial deuteration ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ Partial deuteration of protein complexes is a promising tool for examination of protein structure function relationship. E.g. see an interesting example about the complex in the cyanobacterial circadian clock system Sugiyama et al, Nature Communications Biology (2022). https://doi.org/10.1038/s42003-022-03143-z Here we see how to selectively deuterate domains or specific aminoacids and how scattering is influenced. The deuteration is considered in all bio methods. .. literalinclude:: ../../examples/example_bio_partialdeuteratedprotein.py :lines: 1- .. image:: ../../examples/images/bio_protein_partialmatching.png :align: center :height: 300px :alt: bio_protein_partialmatching.png A protein normal mode animation ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ Creating a trajectory from normal modes and show I(Q) together with the configuration. The calculation can also be done along a MD simulation in ``.trajectory`` . MDAnalysis allows a bunch of trajectory manipulations like center of mass removal or rotations. .. literalinclude:: ../../examples/example_bio_makeNManimation.py :lines: 1- .. image:: ../../examples/images/mode_animation.gif :align: center :width: 50 % :alt: mode_animation Protein density determination ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ We compare measured protein density data with the density calculated by Jscatter to check the SESVolume accuracy in the module :py:mod:`jscatter.bio` . Respective references are given in ``js.examples.datapath+'/proteinSpecificVolumes.txt'`` .. literalinclude:: ../../examples/example_bio_proteinVolume.py :lines: 1- .. image:: ../../examples/images/proteinDensityTest.png :align: center :height: 300px :alt: protein Density Cα coarse grain radius ^^^^^^^^^^^^^^^^^^^^^^ We calculate an appropriate Size of residues dummy atoms for Cα atom models to get a reasonable protein scattering intensity. .. literalinclude:: ../../examples/example_bio_proteinCoarseGrainRadius.py :lines: 1- .. image:: ../../examples/images/proteinCacoarsegrainRadius.png :align: center :height: 300px :alt: proteinCacoarsegrainRadius.png