proteusPy.Disulfide
This module, Disulfide, is part of the proteusPy package, a Python package for the analysis and modeling of protein structures, with an emphasis on disulfide bonds. It represents the core of the current implementation of proteusPy.
This work is based on the original C/C++ implementation by Eric G. Suchanek.
Author: Eric G. Suchanek, PhD Last revision: 2/17/2024
1""" 2This module, *Disulfide*, is part of the proteusPy package, a Python package for 3the analysis and modeling of protein structures, with an emphasis on disulfide bonds. 4It represents the core of the current implementation of *proteusPy*. 5 6This work is based on the original C/C++ implementation by Eric G. Suchanek. \n 7Author: Eric G. Suchanek, PhD 8Last revision: 2/17/2024 9""" 10 11# Cα N, Cα, Cβ, C', Sγ Å ° ρ 12 13import copy 14import datetime 15import glob 16import math 17import pickle 18import time 19from math import cos 20 21import numpy as np 22import pandas 23import pyvista as pv 24 25np.set_printoptions(suppress=True) 26pv.global_theme.color = "white" 27 28try: 29 # Check if running in Jupyter 30 shell = get_ipython().__class__.__name__ 31 if shell == "ZMQInteractiveShell": 32 from tqdm.notebook import tqdm 33 else: 34 from tqdm import tqdm 35except NameError: 36 from tqdm import tqdm 37 38from Bio.PDB import PDBList, PDBParser, Vector 39from Bio.PDB.vectors import calc_dihedral 40 41import proteusPy 42from proteusPy.atoms import * 43from proteusPy.data import ( 44 PROBLEM_ID_FILE, 45 SS_DICT_PICKLE_FILE, 46 SS_ID_FILE, 47 SS_PICKLE_FILE, 48 SS_TORSIONS_FILE, 49) 50from proteusPy.DisulfideExceptions import * 51from proteusPy.DisulfideList import DisulfideList 52 53# tqdm progress bar width 54from proteusPy.ProteusGlobals import ( 55 _ANG_INIT, 56 _FLOAT_INIT, 57 MODEL_DIR, 58 PBAR_COLS, 59 PDB_DIR, 60 WINSIZE, 61) 62from proteusPy.Residue import build_residue 63from proteusPy.turtle3D import ORIENT_SIDECHAIN, Turtle3D 64from proteusPy.utility import distance3d, prune_extra_ss 65 66# columns for the torsions file dataframe. 67Torsion_DF_Cols = [ 68 "source", 69 "ss_id", 70 "proximal", 71 "distal", 72 "chi1", 73 "chi2", 74 "chi3", 75 "chi4", 76 "chi5", 77 "energy", 78 "ca_distance", 79 "cb_distance", 80 "phi_prox", 81 "psi_prox", 82 "phi_dist", 83 "psi_dist", 84 "torsion_length", 85 "rho", 86] 87 88 89# class for the Disulfide bond 90class Disulfide: 91 """ 92 This class provides a Python object and methods representing a physical disulfide bond 93 either extracted from the RCSB protein databank or built using the 94 [proteusPy.Turtle3D](turtle3D.html) class. The disulfide bond is an important 95 intramolecular stabilizing structural element and is characterized by: 96 97 * Atomic coordinates for the atoms N, Cα, Cβ, C', Sγ for both residues. 98 These are stored as both raw atomic coordinates as read from the RCSB file 99 and internal local coordinates. 100 * The dihedral angles Χ1 - Χ5 for the disulfide bond 101 * A name, by default {pdb_id}{prox_resnumb}{prox_chain}_{distal_resnum}{distal_chain} 102 * Proximal residue number 103 * Distal residue number 104 * Approximate bond torsional energy (kcal/mol): 105 106 $$ 107 E_{kcal/mol} \\approx 2.0 * cos(3.0 * \\chi_{1}) + cos(3.0 * \\chi_{5}) + cos(3.0 * \\chi_{2}) + 108 $$ 109 $$ 110 cos(3.0 * \chi_{4}) + 3.5 * cos(2.0 * \chi_{3}) + 0.6 * cos(3.0 * \chi_{3}) + 10.1 111 $$ 112 113 The equation embodies the typical 3-fold rotation barriers associated with single bonds, 114 (Χ1, Χ5, Χ2, Χ4) and a high 2-fold barrier for Χ3, resulting from the partial double bond 115 character of the S-S bond. This property leads to two major disulfide families, characterized 116 by the sign of Χ3. *Left-handed* disulfides have Χ3 < 0° and *right-handed* disulfides have 117 Χ3 > 0°. Within this breakdown there are numerous subfamilies, broadly known as the *hook*, 118 *spiral* and *staple*. These are under characgterization. 119 120 * Euclidean length of the dihedral angles (degrees) defined as: 121 $$\sqrt(\chi_{1}^{2} + \chi_{2}^{2} + \chi_{3}^{2} + \chi_{4}^{2} + \chi_{5}^{2})$$ 122 * Cα - Cα distance (Å) 123 * Cβ - Cβ distance (Å) 124 * The previous C' and next N for both the proximal and distal residues. These are needed 125 to calculate the backbone dihedral angles Φ and Ψ. 126 * Backbone dihedral angles Φ and Ψ, when possible. Not all structures are complete and 127 in those cases the atoms needed may be undefined. In this case the Φ and Ψ angles are set 128 to -180°. 129 130 The class also provides a rendering capabilities using the excellent [PyVista](https://pyvista.org) 131 library, and can display disulfides interactively in a variety of display styles: 132 * 'sb' - Split Bonds style - bonds colored by their atom type 133 * 'bs' - Ball and Stick style - split bond coloring with small atoms 134 * 'pd' - Proximal/Distal style - bonds colored *Red* for proximal residue and *Green* for 135 the distal residue. 136 * 'cpk' - CPK style rendering, colored by atom type: 137 * Carbon - Grey 138 * Nitrogen - Blue 139 * Sulfur - Yellow 140 * Oxygen - Red 141 * Hydrogen - White 142 143 Individual renderings can be saved to a file, and animations created. 144 """ 145 146 def __init__( 147 self, 148 name: str = "SSBOND", 149 proximal: int = -1, 150 distal: int = -1, 151 proximal_chain: str = "A", 152 distal_chain: str = "A", 153 pdb_id: str = "1egs", 154 quiet: bool = True, 155 permissive: bool = True, 156 ) -> None: 157 """ 158 __init__ Initialize the class to defined internal values. 159 160 :param name: Disulfide name, by default "SSBOND" 161 162 """ 163 self.name = name 164 self.proximal = proximal 165 self.distal = distal 166 self.energy = _FLOAT_INIT 167 self.proximal_chain = proximal_chain 168 self.distal_chain = distal_chain 169 self.pdb_id = pdb_id 170 self.proximal_residue_fullid = str("") 171 self.distal_residue_fullid = str("") 172 self.PERMISSIVE = permissive 173 self.QUIET = quiet 174 self.ca_distance = _FLOAT_INIT 175 self.cb_distance = _FLOAT_INIT 176 self.torsion_array = np.array( 177 (_ANG_INIT, _ANG_INIT, _ANG_INIT, _ANG_INIT, _ANG_INIT) 178 ) 179 self.phiprox = _ANG_INIT 180 self.psiprox = _ANG_INIT 181 self.phidist = _ANG_INIT 182 self.psidist = _ANG_INIT 183 184 # global coordinates for the Disulfide, typically as 185 # returned from the PDB file 186 187 self.n_prox = Vector(0, 0, 0) 188 self.ca_prox = Vector(0, 0, 0) 189 self.c_prox = Vector(0, 0, 0) 190 self.o_prox = Vector(0, 0, 0) 191 self.cb_prox = Vector(0, 0, 0) 192 self.sg_prox = Vector(0, 0, 0) 193 self.sg_dist = Vector(0, 0, 0) 194 self.cb_dist = Vector(0, 0, 0) 195 self.ca_dist = Vector(0, 0, 0) 196 self.n_dist = Vector(0, 0, 0) 197 self.c_dist = Vector(0, 0, 0) 198 self.o_dist = Vector(0, 0, 0) 199 200 # set when we can't find previous or next prox or distal 201 # C' or N atoms. 202 self.missing_atoms = False 203 self.modelled = False 204 self.resolution = -1.0 205 206 # need these to calculate backbone dihedral angles 207 self.c_prev_prox = Vector(0, 0, 0) 208 self.n_next_prox = Vector(0, 0, 0) 209 self.c_prev_dist = Vector(0, 0, 0) 210 self.n_next_dist = Vector(0, 0, 0) 211 212 # local coordinates for the Disulfide, computed using the Turtle3D in 213 # Orientation #1. these are generally private. 214 215 self._n_prox = Vector(0, 0, 0) 216 self._ca_prox = Vector(0, 0, 0) 217 self._c_prox = Vector(0, 0, 0) 218 self._o_prox = Vector(0, 0, 0) 219 self._cb_prox = Vector(0, 0, 0) 220 self._sg_prox = Vector(0, 0, 0) 221 self._sg_dist = Vector(0, 0, 0) 222 self._cb_dist = Vector(0, 0, 0) 223 self._ca_dist = Vector(0, 0, 0) 224 self._n_dist = Vector(0, 0, 0) 225 self._c_dist = Vector(0, 0, 0) 226 self._o_dist = Vector(0, 0, 0) 227 228 # need these to calculate backbone dihedral angles 229 self._c_prev_prox = Vector(0, 0, 0) 230 self._n_next_prox = Vector(0, 0, 0) 231 self._c_prev_dist = Vector(0, 0, 0) 232 self._n_next_dist = Vector(0, 0, 0) 233 234 # Dihedral angles for the disulfide bond itself, set to _ANG_INIT 235 self.chi1 = _ANG_INIT 236 self.chi2 = _ANG_INIT 237 self.chi3 = _ANG_INIT 238 self.chi4 = _ANG_INIT 239 self.chi5 = _ANG_INIT 240 self.rho = _ANG_INIT # new dihedral angle: Nprox - Ca_prox - Ca_dist - N_dist 241 242 self.torsion_length = _FLOAT_INIT 243 244 # comparison operators, used for sorting. keyed to SS bond energy 245 def __lt__(self, other): 246 if isinstance(other, Disulfide): 247 return self.energy < other.energy 248 249 def __le__(self, other): 250 if isinstance(other, Disulfide): 251 return self.energy <= other.energy 252 253 def __gt__(self, other): 254 if isinstance(other, Disulfide): 255 return self.energy > other.energy 256 257 def __ge__(self, other): 258 if isinstance(other, Disulfide): 259 return self.energy >= other.energy 260 261 def __eq__(self, other): 262 if isinstance(other, Disulfide): 263 return self.energy == other.energy 264 265 def __ne__(self, other): 266 if isinstance(other, Disulfide): 267 return self.energy != other.energy 268 269 def __repr__(self): 270 """ 271 Representation for the Disulfide class 272 """ 273 s1 = self.repr_ss_info() 274 res = f"{s1}>" 275 return res 276 277 def _render( 278 self, 279 pvplot: pv.Plotter, 280 style="bs", 281 plain=False, 282 bondcolor=BOND_COLOR, 283 bs_scale=BS_SCALE, 284 spec=SPECULARITY, 285 specpow=SPEC_POWER, 286 translate=True, 287 bond_radius=BOND_RADIUS, 288 res=100, 289 ): 290 """ 291 Update the passed pyVista plotter() object with the mesh data for the 292 input Disulfide Bond. Used internally 293 294 Parameters 295 ---------- 296 pvplot : pv.Plotter 297 pyvista.Plotter object 298 299 style : str, optional 300 Rendering style, by default 'bs'. One of 'bs', 'st', 'cpk', Render as \ 301 CPK, ball-and-stick or stick. Bonds are colored by atom color, unless \ 302 'plain' is specified. 303 304 plain : bool, optional 305 Used internally, by default False 306 307 bondcolor : pyVista color name, optional bond color for simple bonds, by default BOND_COLOR 308 309 bs_scale : float, optional 310 scale factor (0-1) to reduce the atom sizes for ball and stick, by default BS_SCALE 311 312 spec : float, optional 313 specularity (0-1), where 1 is totally smooth and 0 is rough, by default SPECULARITY 314 315 specpow : int, optional 316 exponent used for specularity calculations, by default SPEC_POWER 317 318 translate : bool, optional 319 Flag used internally to indicate if we should translate \ 320 the disulfide to its geometric center of mass, by default True. 321 322 Returns 323 ------- 324 pv.Plotter 325 Updated pv.Plotter object with atoms and bonds. 326 """ 327 328 _bradius = bond_radius 329 coords = self.internal_coords() 330 missing_atoms = self.missing_atoms 331 clen = coords.shape[0] 332 333 model = self.modelled 334 if model: 335 all_atoms = False 336 else: 337 all_atoms = True 338 339 if translate: 340 cofmass = self.cofmass() 341 for i in range(clen): 342 coords[i] = coords[i] - cofmass 343 344 atoms = ( 345 "N", 346 "C", 347 "C", 348 "O", 349 "C", 350 "SG", 351 "N", 352 "C", 353 "C", 354 "O", 355 "C", 356 "SG", 357 "C", 358 "N", 359 "C", 360 "N", 361 ) 362 pvp = pvplot 363 364 # bond connection table with atoms in the specific order shown above: 365 # returned by ss.get_internal_coords() 366 367 def draw_bonds( 368 pvp, 369 bradius=BOND_RADIUS, 370 style="sb", 371 bcolor=BOND_COLOR, 372 missing=True, 373 all_atoms=True, 374 res=100, 375 ): 376 """ 377 Generate the appropriate pyVista cylinder objects to represent 378 a particular disulfide bond. This utilizes a connection table 379 for the starting and ending atoms and a color table for the 380 bond colors. Used internally. 381 382 :param pvp: input plotter object to be updated 383 :param bradius: bond radius 384 :param style: bond style. One of sb, plain, pd 385 :param bcolor: pyvista color 386 :param missing: True if atoms are missing, False othersie 387 :param all_atoms: True if rendering O, False if only backbone rendered 388 389 :return pvp: Updated Plotter object. 390 391 """ 392 _bond_conn = np.array( 393 [ 394 [0, 1], # n-ca 395 [1, 2], # ca-c 396 [2, 3], # c-o 397 [1, 4], # ca-cb 398 [4, 5], # cb-sg 399 [6, 7], # n-ca 400 [7, 8], # ca-c 401 [8, 9], # c-o 402 [7, 10], # ca-cb 403 [10, 11], # cb-sg 404 [5, 11], # sg -sg 405 [12, 0], # cprev_prox-n 406 [2, 13], # c-nnext_prox 407 [14, 6], # cprev_dist-n_dist 408 [8, 15], # c-nnext_dist 409 ] 410 ) 411 412 # modeled disulfides only have backbone atoms since 413 # phi and psi are undefined, which makes the carbonyl 414 # oxygen (O) undefined as well. Their previous and next N 415 # are also undefined. 416 417 _bond_conn_backbone = np.array( 418 [ 419 [0, 1], # n-ca 420 [1, 2], # ca-c 421 [1, 4], # ca-cb 422 [4, 5], # cb-sg 423 [6, 7], # n-ca 424 [7, 8], # ca-c 425 [7, 10], # ca-cb 426 [10, 11], # cb-sg 427 [5, 11], # sg -sg 428 ] 429 ) 430 431 # colors for the bonds. Index into ATOM_COLORS array 432 _bond_split_colors = np.array( 433 [ 434 ("N", "C"), 435 ("C", "C"), 436 ("C", "O"), 437 ("C", "C"), 438 ("C", "SG"), 439 ("N", "C"), 440 ("C", "C"), 441 ("C", "O"), 442 ("C", "C"), 443 ("C", "SG"), 444 ("SG", "SG"), 445 # prev and next C-N bonds - color by atom Z 446 ("C", "N"), 447 ("C", "N"), 448 ("C", "N"), 449 ("C", "N"), 450 ] 451 ) 452 453 _bond_split_colors_backbone = np.array( 454 [ 455 ("N", "C"), 456 ("C", "C"), 457 ("C", "C"), 458 ("C", "SG"), 459 ("N", "C"), 460 ("C", "C"), 461 ("C", "C"), 462 ("C", "SG"), 463 ("SG", "SG"), 464 ] 465 ) 466 # work through connectivity and colors 467 orig_col = dest_col = bcolor 468 469 if all_atoms: 470 bond_conn = _bond_conn 471 bond_split_colors = _bond_split_colors 472 else: 473 bond_conn = _bond_conn_backbone 474 bond_split_colors = _bond_split_colors_backbone 475 476 for i in range(len(bond_conn)): 477 if all_atoms: 478 if i > 10 and missing_atoms == True: # skip missing atoms 479 continue 480 481 bond = bond_conn[i] 482 483 # get the indices for the origin and destination atoms 484 orig = bond[0] 485 dest = bond[1] 486 487 col = bond_split_colors[i] 488 489 # get the coords 490 prox_pos = coords[orig] 491 distal_pos = coords[dest] 492 493 # compute a direction vector 494 direction = distal_pos - prox_pos 495 496 # compute vector length. divide by 2 since split bond 497 height = math.dist(prox_pos, distal_pos) / 2.0 498 499 # the cylinder origins are actually in the 500 # middle so we translate 501 502 origin = prox_pos + 0.5 * direction # for a single plain bond 503 origin1 = prox_pos + 0.25 * direction 504 origin2 = prox_pos + 0.75 * direction 505 506 bradius = _bradius 507 508 if style == "plain": 509 orig_col = dest_col = bcolor 510 511 # proximal-distal red/green coloring 512 elif style == "pd": 513 if i <= 4 or i == 11 or i == 12: 514 orig_col = dest_col = "red" 515 else: 516 orig_col = dest_col = "green" 517 if i == 10: 518 orig_col = dest_col = "yellow" 519 else: 520 orig_col = ATOM_COLORS[col[0]] 521 dest_col = ATOM_COLORS[col[1]] 522 523 if i >= 11: # prev and next residue atoms for phi/psi calcs 524 bradius = _bradius * 0.5 # make smaller to distinguish 525 526 cap1 = pv.Sphere(center=prox_pos, radius=bradius) 527 cap2 = pv.Sphere(center=distal_pos, radius=bradius) 528 529 if style == "plain": 530 cyl = pv.Cylinder( 531 origin, direction, radius=bradius, height=height * 2.0 532 ) 533 pvp.add_mesh(cyl, color=orig_col) 534 else: 535 cyl1 = pv.Cylinder( 536 origin1, 537 direction, 538 radius=bradius, 539 height=height, 540 capping=False, 541 resolution=res, 542 ) 543 cyl2 = pv.Cylinder( 544 origin2, 545 direction, 546 radius=bradius, 547 height=height, 548 capping=False, 549 resolution=res, 550 ) 551 pvp.add_mesh(cyl1, color=orig_col) 552 pvp.add_mesh(cyl2, color=dest_col) 553 554 pvp.add_mesh(cap1, color=orig_col) 555 pvp.add_mesh(cap2, color=dest_col) 556 557 return pvp # end draw_bonds 558 559 if style == "cpk": 560 i = 0 561 for atom in atoms: 562 rad = ATOM_RADII_CPK[atom] 563 pvp.add_mesh( 564 pv.Sphere(center=coords[i], radius=rad), 565 color=ATOM_COLORS[atom], 566 smooth_shading=True, 567 specular=spec, 568 specular_power=specpow, 569 ) 570 i += 1 571 572 elif style == "cov": 573 i = 0 574 for atom in atoms: 575 rad = ATOM_RADII_COVALENT[atom] 576 pvp.add_mesh( 577 pv.Sphere(center=coords[i], radius=rad), 578 color=ATOM_COLORS[atom], 579 smooth_shading=True, 580 specular=spec, 581 specular_power=specpow, 582 ) 583 i += 1 584 585 elif style == "bs": # ball and stick 586 i = 0 587 for atom in atoms: 588 rad = ATOM_RADII_CPK[atom] * bs_scale 589 if i > 11: 590 rad = rad * 0.75 591 592 pvp.add_mesh( 593 pv.Sphere(center=coords[i], radius=rad), 594 color=ATOM_COLORS[atom], 595 smooth_shading=True, 596 specular=spec, 597 specular_power=specpow, 598 ) 599 i += 1 600 pvp = draw_bonds( 601 pvp, style="bs", missing=missing_atoms, all_atoms=all_atoms 602 ) 603 604 elif style == "sb": # splitbonds 605 pvp = draw_bonds( 606 pvp, style="sb", missing=missing_atoms, all_atoms=all_atoms 607 ) 608 609 elif style == "pd": # proximal-distal 610 pvp = draw_bonds( 611 pvp, style="pd", missing=missing_atoms, all_atoms=all_atoms 612 ) 613 614 else: # plain 615 pvp = draw_bonds( 616 pvp, 617 style="plain", 618 bcolor=bondcolor, 619 missing=missing_atoms, 620 all_atoms=all_atoms, 621 ) 622 623 return pvp 624 625 def _plot( 626 self, 627 pvplot, 628 style="bs", 629 plain=False, 630 bondcolor=BOND_COLOR, 631 bs_scale=BS_SCALE, 632 spec=SPECULARITY, 633 specpow=SPEC_POWER, 634 translate=True, 635 bond_radius=BOND_RADIUS, 636 res=100, 637 ): 638 """ 639 Update the passed pyVista plotter() object with the mesh data for the 640 input Disulfide Bond. Used internally 641 642 Parameters 643 ---------- 644 pvplot : pv.Plotter 645 pyvista.Plotter object 646 647 style : str, optional 648 Rendering style, by default 'bs'. One of 'bs', 'st', 'cpk', Render as \ 649 CPK, ball-and-stick or stick. Bonds are colored by atom color, unless \ 650 'plain' is specified. 651 652 plain : bool, optional 653 Used internally, by default False 654 655 bondcolor : pyVista color name, optional bond color for simple bonds, by default BOND_COLOR 656 657 bs_scale : float, optional 658 scale factor (0-1) to reduce the atom sizes for ball and stick, by default BS_SCALE 659 660 spec : float, optional 661 specularity (0-1), where 1 is totally smooth and 0 is rough, by default SPECULARITY 662 663 specpow : int, optional 664 exponent used for specularity calculations, by default SPEC_POWER 665 666 translate : bool, optional 667 Flag used internally to indicate if we should translate \ 668 the disulfide to its geometric center of mass, by default True. 669 670 Returns 671 ------- 672 pv.Plotter 673 Updated pv.Plotter object with atoms and bonds. 674 """ 675 676 _bradius = bond_radius 677 coords = self.internal_coords() 678 missing_atoms = self.missing_atoms 679 clen = coords.shape[0] 680 681 model = self.modelled 682 if model: 683 all_atoms = False 684 else: 685 all_atoms = True 686 687 if translate: 688 cofmass = self.cofmass() 689 for i in range(clen): 690 coords[i] = coords[i] - cofmass 691 692 atoms = ( 693 "N", 694 "C", 695 "C", 696 "O", 697 "C", 698 "SG", 699 "N", 700 "C", 701 "C", 702 "O", 703 "C", 704 "SG", 705 "C", 706 "N", 707 "C", 708 "N", 709 ) 710 pvp = pvplot.copy() 711 712 # bond connection table with atoms in the specific order shown above: 713 # returned by ss.get_internal_coords() 714 715 def draw_bonds( 716 pvp, 717 bradius=BOND_RADIUS, 718 style="sb", 719 bcolor=BOND_COLOR, 720 missing=True, 721 all_atoms=True, 722 res=100, 723 ): 724 """ 725 Generate the appropriate pyVista cylinder objects to represent 726 a particular disulfide bond. This utilizes a connection table 727 for the starting and ending atoms and a color table for the 728 bond colors. Used internally. 729 730 :param pvp: input plotter object to be updated 731 :param bradius: bond radius 732 :param style: bond style. One of sb, plain, pd 733 :param bcolor: pyvista color 734 :param missing: True if atoms are missing, False othersie 735 :param all_atoms: True if rendering O, False if only backbone rendered 736 737 :return pvp: Updated Plotter object. 738 739 """ 740 _bond_conn = np.array( 741 [ 742 [0, 1], # n-ca 743 [1, 2], # ca-c 744 [2, 3], # c-o 745 [1, 4], # ca-cb 746 [4, 5], # cb-sg 747 [6, 7], # n-ca 748 [7, 8], # ca-c 749 [8, 9], # c-o 750 [7, 10], # ca-cb 751 [10, 11], # cb-sg 752 [5, 11], # sg -sg 753 [12, 0], # cprev_prox-n 754 [2, 13], # c-nnext_prox 755 [14, 6], # cprev_dist-n_dist 756 [8, 15], # c-nnext_dist 757 ] 758 ) 759 760 # modeled disulfides only have backbone atoms since 761 # phi and psi are undefined, which makes the carbonyl 762 # oxygen (O) undefined as well. Their previous and next N 763 # are also undefined. 764 765 _bond_conn_backbone = np.array( 766 [ 767 [0, 1], # n-ca 768 [1, 2], # ca-c 769 [1, 4], # ca-cb 770 [4, 5], # cb-sg 771 [6, 7], # n-ca 772 [7, 8], # ca-c 773 [7, 10], # ca-cb 774 [10, 11], # cb-sg 775 [5, 11], # sg -sg 776 ] 777 ) 778 779 # colors for the bonds. Index into ATOM_COLORS array 780 _bond_split_colors = np.array( 781 [ 782 ("N", "C"), 783 ("C", "C"), 784 ("C", "O"), 785 ("C", "C"), 786 ("C", "SG"), 787 ("N", "C"), 788 ("C", "C"), 789 ("C", "O"), 790 ("C", "C"), 791 ("C", "SG"), 792 ("SG", "SG"), 793 # prev and next C-N bonds - color by atom Z 794 ("C", "N"), 795 ("C", "N"), 796 ("C", "N"), 797 ("C", "N"), 798 ] 799 ) 800 801 _bond_split_colors_backbone = np.array( 802 [ 803 ("N", "C"), 804 ("C", "C"), 805 ("C", "C"), 806 ("C", "SG"), 807 ("N", "C"), 808 ("C", "C"), 809 ("C", "C"), 810 ("C", "SG"), 811 ("SG", "SG"), 812 ] 813 ) 814 # work through connectivity and colors 815 orig_col = dest_col = bcolor 816 817 if all_atoms: 818 bond_conn = _bond_conn 819 bond_split_colors = _bond_split_colors 820 else: 821 bond_conn = _bond_conn_backbone 822 bond_split_colors = _bond_split_colors_backbone 823 824 for i in range(len(bond_conn)): 825 if all_atoms: 826 if i > 10 and missing_atoms == True: # skip missing atoms 827 continue 828 829 bond = bond_conn[i] 830 831 # get the indices for the origin and destination atoms 832 orig = bond[0] 833 dest = bond[1] 834 835 col = bond_split_colors[i] 836 837 # get the coords 838 prox_pos = coords[orig] 839 distal_pos = coords[dest] 840 841 # compute a direction vector 842 direction = distal_pos - prox_pos 843 844 # compute vector length. divide by 2 since split bond 845 height = math.dist(prox_pos, distal_pos) / 2.0 846 847 # the cylinder origins are actually in the 848 # middle so we translate 849 850 origin = prox_pos + 0.5 * direction # for a single plain bond 851 origin1 = prox_pos + 0.25 * direction 852 origin2 = prox_pos + 0.75 * direction 853 854 bradius = _bradius 855 856 if style == "plain": 857 orig_col = dest_col = bcolor 858 859 # proximal-distal red/green coloring 860 elif style == "pd": 861 if i <= 4 or i == 11 or i == 12: 862 orig_col = dest_col = "red" 863 else: 864 orig_col = dest_col = "green" 865 if i == 10: 866 orig_col = dest_col = "yellow" 867 else: 868 orig_col = ATOM_COLORS[col[0]] 869 dest_col = ATOM_COLORS[col[1]] 870 871 if i >= 11: # prev and next residue atoms for phi/psi calcs 872 bradius = _bradius * 0.5 # make smaller to distinguish 873 874 cap1 = pv.Sphere(center=prox_pos, radius=bradius) 875 cap2 = pv.Sphere(center=distal_pos, radius=bradius) 876 877 if style == "plain": 878 cyl = pv.Cylinder( 879 origin, direction, radius=bradius, height=height * 2.0 880 ) 881 # pvp.add_mesh(cyl, color=orig_col) 882 pvp.append(cyl) 883 else: 884 cyl1 = pv.Cylinder( 885 origin1, 886 direction, 887 radius=bradius, 888 height=height, 889 capping=False, 890 resolution=res, 891 ) 892 cyl2 = pv.Cylinder( 893 origin2, 894 direction, 895 radius=bradius, 896 height=height, 897 capping=False, 898 resolution=res, 899 ) 900 # pvp.add_mesh(cyl1, color=orig_col) 901 # pvp.add_mesh(cyl2, color=dest_col) 902 pvp.append(cyl1) 903 pvp.append(cyl2) 904 905 # pvp.add_mesh(cap1, color=orig_col) 906 # pvp.add_mesh(cap2, color=dest_col) 907 pvp.append(cap1) 908 pvp.append(cap2) 909 910 return pvp.copy() # end draw_bonds 911 912 if style == "cpk": 913 i = 0 914 for atom in atoms: 915 rad = ATOM_RADII_CPK[atom] 916 pvp.append(pv.Sphere(center=coords[i], radius=rad)) 917 i += 1 918 919 elif style == "cov": 920 i = 0 921 for atom in atoms: 922 rad = ATOM_RADII_COVALENT[atom] 923 pvp.append(pv.Sphere(center=coords[i], radius=rad)) 924 i += 1 925 926 elif style == "bs": # ball and stick 927 i = 0 928 for atom in atoms: 929 rad = ATOM_RADII_CPK[atom] * bs_scale 930 if i > 11: 931 rad = rad * 0.75 932 933 pvp.append(pv.Sphere(center=coords[i])) 934 i += 1 935 pvp = draw_bonds( 936 pvp, style="bs", missing=missing_atoms, all_atoms=all_atoms 937 ) 938 939 elif style == "sb": # splitbonds 940 pvp = draw_bonds( 941 pvp, style="sb", missing=missing_atoms, all_atoms=all_atoms 942 ) 943 944 elif style == "pd": # proximal-distal 945 pvp = draw_bonds( 946 pvp, style="pd", missing=missing_atoms, all_atoms=all_atoms 947 ) 948 949 else: # plain 950 pvp = draw_bonds( 951 pvp, 952 style="plain", 953 bcolor=bondcolor, 954 missing=missing_atoms, 955 all_atoms=all_atoms, 956 ) 957 958 return 959 960 def _handle_SS_exception(self, message: str): 961 """ 962 This method catches an exception that occurs in the Disulfide 963 object (if PERMISSIVE), or raises it again, this time adding the 964 PDB line number to the error message. (private). 965 966 :param message: Error message 967 :raises DisulfideConstructionException: Fatal construction exception. 968 969 """ 970 # message = "%s at line %i." % (message) 971 message = f"{message}" 972 973 if self.PERMISSIVE: 974 # just print a warning - some residues/atoms may be missing 975 warnings.warn( 976 "DisulfideConstructionException: %s\n" 977 "Exception ignored.\n" 978 "Some atoms may be missing in the data structure." % message, 979 DisulfideConstructionWarning, 980 ) 981 else: 982 # exceptions are fatal - raise again with new message (including line nr) 983 raise DisulfideConstructionException(message) from None 984 985 @property 986 def dihedrals(self) -> list: 987 """ 988 Return a list containing the dihedral angles for the disulfide. 989 990 """ 991 return [self.chi1, self.chi2, self.chi3, self.chi4, self.chi5] 992 993 @dihedrals.setter 994 def dihedrals(self, dihedrals: list) -> None: 995 """ 996 Sets the disulfide dihedral angles to the inputs specified in the list. 997 998 :param dihedrals: list of dihedral angles. 999 """ 1000 self.chi1 = dihedrals[0] 1001 self.chi2 = dihedrals[1] 1002 self.chi3 = dihedrals[2] 1003 self.chi4 = dihedrals[3] 1004 self.chi5 = dihedrals[4] 1005 self.compute_torsional_energy() 1006 self.compute_torsion_length() 1007 1008 def bounding_box(self) -> np.array: 1009 """ 1010 Return the bounding box array for the given disulfide 1011 1012 Returns 1013 ------- 1014 :return: np.Array(3,2): Array containing the min, max for X, Y, and Z respectively. 1015 Does not currently take the atom's radius into account. 1016 1017 """ 1018 res = np.zeros(shape=(3, 2)) 1019 xmin, xmax = self.compute_extents("x") 1020 ymin, ymax = self.compute_extents("y") 1021 zmin, zmax = self.compute_extents("z") 1022 1023 res[0] = [xmin, xmax] 1024 res[1] = [ymin, ymax] 1025 res[2] = [zmin, zmax] 1026 1027 return res 1028 1029 def build_yourself(self) -> None: 1030 """ 1031 Build a model Disulfide based its internal dihedral state 1032 Routine assumes turtle is in orientation #1 (at Ca, headed toward 1033 Cb, with N on left), builds disulfide, and updates the object's internal 1034 coordinates. It also adds the distal protein backbone, 1035 and computes the disulfide conformational energy. 1036 """ 1037 chi1 = self.chi1 1038 chi2 = self.chi2 1039 chi3 = self.chi3 1040 chi4 = self.chi4 1041 chi5 = self.chi5 1042 self.build_model(chi1, chi2, chi3, chi4, chi5) 1043 1044 def build_model( 1045 self, chi1: float, chi2: float, chi3: float, chi4: float, chi5: float 1046 ) -> None: 1047 """ 1048 Build a model Disulfide based on the input dihedral angles. 1049 Routine assumes turtle is in orientation #1 (at Ca, headed toward 1050 Cb, with N on left), builds disulfide, and updates the object's internal 1051 coordinates. It also adds the distal protein backbone, 1052 and computes the disulfide conformational energy. 1053 1054 :param chi1: Chi1 (degrees) 1055 :param chi2: Chi2 (degrees) 1056 :param chi3: Chi3 (degrees) 1057 :param chi4: Chi4 (degrees) 1058 :param chi5: Chi5 (degrees) 1059 1060 Example: 1061 >>> from proteusPy.Disulfide import Disulfide 1062 >>> modss = Disulfide('model') 1063 >>> modss.build_model(-60, -60, -90, -60, -60) 1064 >>> modss.display(style='sb') 1065 """ 1066 1067 self.set_dihedrals(chi1, chi2, chi3, chi4, chi5) 1068 self.proximal = 1 1069 self.distal = 2 1070 1071 tmp = Turtle3D("tmp") 1072 tmp.Orientation = 1 1073 1074 n = Vector(0, 0, 0) 1075 ca = Vector(0, 0, 0) 1076 cb = Vector(0, 0, 0) 1077 c = Vector(0, 0, 0) 1078 1079 self.ca_prox = tmp._position 1080 tmp.schain_to_bbone() 1081 n, ca, cb, c = build_residue(tmp) 1082 1083 self.n_prox = n 1084 self.ca_prox = ca 1085 self.c_prox = c 1086 self.cb_prox = cb 1087 1088 tmp.bbone_to_schain() 1089 tmp.move(1.53) 1090 tmp.roll(self.chi1) 1091 tmp.yaw(112.8) 1092 self.cb_prox = Vector(tmp._position) 1093 1094 tmp.move(1.86) 1095 tmp.roll(self.chi2) 1096 tmp.yaw(103.8) 1097 self.sg_prox = Vector(tmp._position) 1098 1099 tmp.move(2.044) 1100 tmp.roll(self.chi3) 1101 tmp.yaw(103.8) 1102 self.sg_dist = Vector(tmp._position) 1103 1104 tmp.move(1.86) 1105 tmp.roll(self.chi4) 1106 tmp.yaw(112.8) 1107 self.cb_dist = Vector(tmp._position) 1108 1109 tmp.move(1.53) 1110 tmp.roll(self.chi5) 1111 tmp.pitch(180.0) 1112 1113 tmp.schain_to_bbone() 1114 1115 n, ca, cb, c = build_residue(tmp) 1116 1117 self.n_dist = n 1118 self.ca_dist = ca 1119 self.c_dist = c 1120 self.compute_torsional_energy() 1121 self.compute_local_coords() 1122 self.ca_distance = distance3d(self.ca_prox, self.ca_dist) 1123 self.cb_distance = distance3d(self.cb_prox, self.cb_dist) 1124 self.torsion_array = np.array( 1125 (self.chi1, self.chi2, self.chi3, self.chi4, self.chi5) 1126 ) 1127 self.compute_torsion_length() 1128 self.compute_rho() 1129 self.missing_atoms = True 1130 self.modelled = True 1131 1132 def cofmass(self) -> np.array: 1133 """ 1134 Return the geometric center of mass for the internal coordinates of 1135 the given Disulfide. Missing atoms are not included. 1136 1137 :return: 3D array for the geometric center of mass 1138 """ 1139 1140 res = self.internal_coords() 1141 return res.mean(axis=0) 1142 1143 def copy(self): 1144 """ 1145 Copy the Disulfide. 1146 1147 :return: A copy of self. 1148 """ 1149 return copy.deepcopy(self) 1150 1151 def compute_extents(self, dim="z"): 1152 """ 1153 Calculate the internal coordinate extents for the input axis. 1154 1155 :param dim: Axis, one of 'x', 'y', 'z', by default 'z' 1156 :return: min, max 1157 """ 1158 1159 ic = self.internal_coords() 1160 # set default index to 'z' 1161 idx = 2 1162 1163 if dim == "x": 1164 idx = 0 1165 elif dim == "y": 1166 idx = 1 1167 elif dim == "z": 1168 idx = 2 1169 1170 _min = ic[:, idx].min() 1171 _max = ic[:, idx].max() 1172 return _min, _max 1173 1174 def compute_local_coords(self) -> None: 1175 """ 1176 Compute the internal coordinates for a properly initialized Disulfide Object. 1177 1178 :param self: SS initialized Disulfide object 1179 :returns: None, modifies internal state of the input 1180 """ 1181 1182 turt = Turtle3D("tmp") 1183 # get the coordinates as np.array for Turtle3D use. 1184 cpp = self.c_prev_prox.get_array() 1185 nnp = self.n_next_prox.get_array() 1186 1187 n = self.n_prox.get_array() 1188 ca = self.ca_prox.get_array() 1189 c = self.c_prox.get_array() 1190 cb = self.cb_prox.get_array() 1191 o = self.o_prox.get_array() 1192 sg = self.sg_prox.get_array() 1193 1194 sg2 = self.sg_dist.get_array() 1195 cb2 = self.cb_dist.get_array() 1196 ca2 = self.ca_dist.get_array() 1197 c2 = self.c_dist.get_array() 1198 n2 = self.n_dist.get_array() 1199 o2 = self.o_dist.get_array() 1200 1201 cpd = self.c_prev_dist.get_array() 1202 nnd = self.n_next_dist.get_array() 1203 1204 turt.orient_from_backbone(n, ca, c, cb, ORIENT_SIDECHAIN) 1205 1206 # internal (local) coordinates, stored as Vector objects 1207 # to_local returns np.array objects 1208 1209 self._n_prox = Vector(turt.to_local(n)) 1210 self._ca_prox = Vector(turt.to_local(ca)) 1211 self._c_prox = Vector(turt.to_local(c)) 1212 self._o_prox = Vector(turt.to_local(o)) 1213 self._cb_prox = Vector(turt.to_local(cb)) 1214 self._sg_prox = Vector(turt.to_local(sg)) 1215 1216 self._c_prev_prox = Vector(turt.to_local(cpp)) 1217 self._n_next_prox = Vector(turt.to_local(nnp)) 1218 self._c_prev_dist = Vector(turt.to_local(cpd)) 1219 self._n_next_dist = Vector(turt.to_local(nnd)) 1220 1221 self._n_dist = Vector(turt.to_local(n2)) 1222 self._ca_dist = Vector(turt.to_local(ca2)) 1223 self._c_dist = Vector(turt.to_local(c2)) 1224 self._o_dist = Vector(turt.to_local(o2)) 1225 self._cb_dist = Vector(turt.to_local(cb2)) 1226 self._sg_dist = Vector(turt.to_local(sg2)) 1227 1228 def compute_torsional_energy(self) -> float: 1229 """ 1230 Compute the approximate torsional energy for the Disulfide's 1231 conformation and sets its internal state. 1232 1233 :return: Energy (kcal/mol) 1234 """ 1235 # @TODO find citation for the ss bond energy calculation 1236 1237 def torad(deg): 1238 return np.radians(deg) 1239 1240 chi1 = self.chi1 1241 chi2 = self.chi2 1242 chi3 = self.chi3 1243 chi4 = self.chi4 1244 chi5 = self.chi5 1245 1246 energy = 2.0 * (cos(torad(3.0 * chi1)) + cos(torad(3.0 * chi5))) 1247 energy += cos(torad(3.0 * chi2)) + cos(torad(3.0 * chi4)) 1248 energy += 3.5 * cos(torad(2.0 * chi3)) + 0.6 * cos(torad(3.0 * chi3)) + 10.1 1249 1250 self.energy = energy 1251 return energy 1252 1253 def display(self, single=True, style="sb", light=True, shadows=False) -> None: 1254 """ 1255 Display the Disulfide bond in the specific rendering style. 1256 1257 :param single: Display the bond in a single panel in the specific style. 1258 :param style: Rendering style: One of: 1259 * 'sb' - split bonds 1260 * 'bs' - ball and stick 1261 * 'cpk' - CPK style 1262 * 'pd' - Proximal/Distal style - Red=proximal, Green=Distal 1263 * 'plain' - boring single color 1264 :param light: If True, light background, if False, dark 1265 1266 Example: 1267 >>> import proteusPy 1268 >>> from proteusPy.Disulfide import Disulfide 1269 >>> from proteusPy.DisulfideLoader import DisulfideLoader, Load_PDB_SS 1270 1271 >>> PDB_SS = Load_PDB_SS(verbose=False, subset=True) 1272 >>> ss = PDB_SS[0] 1273 >>> ss.display(style='cpk') 1274 >>> ss.screenshot(style='bs', fname='proteus_logo_sb.png') 1275 """ 1276 src = self.pdb_id 1277 enrg = self.energy 1278 1279 title = f"{src}: {self.proximal}{self.proximal_chain}-{self.distal}{self.distal_chain}: {enrg:.2f} kcal/mol. Cα: {self.ca_distance:.2f} Å Cβ: {self.cb_distance:.2f} Å Tors: {self.torsion_length:.2f}°" 1280 1281 if light: 1282 pv.set_plot_theme("document") 1283 else: 1284 pv.set_plot_theme("dark") 1285 1286 if single == True: 1287 _pl = pv.Plotter(window_size=WINSIZE) 1288 _pl.add_title(title=title, font_size=FONTSIZE) 1289 _pl.enable_anti_aliasing("msaa") 1290 # _pl.add_camera_orientation_widget() 1291 1292 self._render( 1293 _pl, 1294 style=style, 1295 bs_scale=BS_SCALE, 1296 spec=SPECULARITY, 1297 specpow=SPEC_POWER, 1298 ) 1299 _pl.reset_camera() 1300 if shadows == True: 1301 _pl.enable_shadows() 1302 _pl.show() 1303 1304 else: 1305 pl = pv.Plotter(window_size=WINSIZE, shape=(2, 2)) 1306 pl.subplot(0, 0) 1307 1308 pl.add_title(title=title, font_size=FONTSIZE) 1309 pl.enable_anti_aliasing("msaa") 1310 1311 # pl.add_camera_orientation_widget() 1312 1313 self._render( 1314 pl, 1315 style="cpk", 1316 bondcolor=BOND_COLOR, 1317 bs_scale=BS_SCALE, 1318 spec=SPECULARITY, 1319 specpow=SPEC_POWER, 1320 ) 1321 1322 pl.subplot(0, 1) 1323 pl.add_title(title=title, font_size=FONTSIZE) 1324 1325 self._render( 1326 pl, 1327 style="bs", 1328 bondcolor=BOND_COLOR, 1329 bs_scale=BS_SCALE, 1330 spec=SPECULARITY, 1331 specpow=SPEC_POWER, 1332 ) 1333 1334 pl.subplot(1, 0) 1335 pl.add_title(title=title, font_size=FONTSIZE) 1336 1337 self._render( 1338 pl, 1339 style="sb", 1340 bondcolor=BOND_COLOR, 1341 bs_scale=BS_SCALE, 1342 spec=SPECULARITY, 1343 specpow=SPEC_POWER, 1344 ) 1345 1346 pl.subplot(1, 1) 1347 pl.add_title(title=title, font_size=FONTSIZE) 1348 1349 self._render( 1350 pl, 1351 style="pd", 1352 bondcolor=BOND_COLOR, 1353 bs_scale=BS_SCALE, 1354 spec=SPECULARITY, 1355 specpow=SPEC_POWER, 1356 ) 1357 1358 pl.link_views() 1359 pl.reset_camera() 1360 if shadows == True: 1361 pl.enable_shadows() 1362 pl.show() 1363 return 1364 1365 def plot( 1366 self, pl, single=True, style="sb", light=True, shadows=False 1367 ) -> pv.Plotter: 1368 """ 1369 Return the pyVista Plotter object for the Disulfide bond in the specific rendering style. 1370 1371 :param single: Display the bond in a single panel in the specific style. 1372 :param style: Rendering style: One of: 1373 * 'sb' - split bonds 1374 * 'bs' - ball and stick 1375 * 'cpk' - CPK style 1376 * 'pd' - Proximal/Distal style - Red=proximal, Green=Distal 1377 * 'plain' - boring single color 1378 :param light: If True, light background, if False, dark 1379 """ 1380 src = self.pdb_id 1381 enrg = self.energy 1382 1383 title = f"{src}: {self.proximal}{self.proximal_chain}-{self.distal}{self.distal_chain}: {enrg:.2f} kcal/mol. Cα: {self.ca_distance:.2f} Å Cβ: {self.cb_distance:.2f} Å Tors: {self.torsion_length:.2f}°" 1384 1385 if light: 1386 pv.set_plot_theme("document") 1387 else: 1388 pv.set_plot_theme("dark") 1389 1390 if single == True: 1391 # _pl = pv.Plotter(window_size=WINSIZE) 1392 # _pl.add_title(title=title, font_size=FONTSIZE) 1393 pl.clear() 1394 pl.enable_anti_aliasing("msaa") 1395 # pl.add_camera_orientation_widget() 1396 1397 self._render( 1398 pl, 1399 style=style, 1400 bs_scale=BS_SCALE, 1401 spec=SPECULARITY, 1402 specpow=SPEC_POWER, 1403 ) 1404 pl.reset_camera() 1405 if shadows == True: 1406 pl.enable_shadows() 1407 else: 1408 pl = pv.Plotter(shape=(2, 2)) 1409 pl.subplot(0, 0) 1410 1411 # pl.add_title(title=title, font_size=FONTSIZE) 1412 pl.enable_anti_aliasing("msaa") 1413 1414 # pl.add_camera_orientation_widget() 1415 1416 self._render( 1417 pl, 1418 style="cpk", 1419 bondcolor=BOND_COLOR, 1420 bs_scale=BS_SCALE, 1421 spec=SPECULARITY, 1422 specpow=SPEC_POWER, 1423 ) 1424 1425 pl.subplot(0, 1) 1426 1427 self._render( 1428 pl, 1429 style="bs", 1430 bondcolor=BOND_COLOR, 1431 bs_scale=BS_SCALE, 1432 spec=SPECULARITY, 1433 specpow=SPEC_POWER, 1434 ) 1435 1436 pl.subplot(1, 0) 1437 1438 self._render( 1439 pl, 1440 style="sb", 1441 bondcolor=BOND_COLOR, 1442 bs_scale=BS_SCALE, 1443 spec=SPECULARITY, 1444 specpow=SPEC_POWER, 1445 ) 1446 1447 pl.subplot(1, 1) 1448 self._render( 1449 pl, 1450 style="pd", 1451 bondcolor=BOND_COLOR, 1452 bs_scale=BS_SCALE, 1453 spec=SPECULARITY, 1454 specpow=SPEC_POWER, 1455 ) 1456 1457 pl.link_views() 1458 pl.reset_camera() 1459 if shadows == True: 1460 pl.enable_shadows() 1461 return pl 1462 1463 def Distance_neighbors(self, others: DisulfideList, cutoff: float) -> DisulfideList: 1464 """ 1465 Return list of Disulfides whose RMS atomic distance is within 1466 the cutoff (Å) in the others list. 1467 1468 :param others: DisulfideList to search 1469 :param cutoff: Distance cutoff (Å) 1470 :return: DisulfideList within the cutoff 1471 """ 1472 1473 res = [ss.copy() for ss in others if self.Distance_RMS(ss) < cutoff] 1474 return DisulfideList(res, "neighbors") 1475 1476 def Distance_RMS(self, other) -> float: 1477 """ 1478 Calculate the RMS distance between the internal coordinates of self and another Disulfide. 1479 :param other: Comparison Disulfide 1480 :return: RMS distance (Å) 1481 """ 1482 1483 # Get internal coordinates of both objects 1484 ic1 = self.internal_coords() 1485 ic2 = other.internal_coords() 1486 1487 # Compute the sum of squared differences between corresponding internal coordinates 1488 totsq = sum(math.dist(p1, p2) ** 2 for p1, p2 in zip(ic1, ic2)) 1489 1490 # Compute the mean of the squared distances 1491 totsq /= len(ic1) 1492 1493 # Take the square root of the mean to get the RMS distance 1494 return math.sqrt(totsq) 1495 1496 def Torsion_RMS(self, other) -> float: 1497 """ 1498 Calculate the RMS distance between the dihedral angles of self and another Disulfide. 1499 :param other: Comparison Disulfide 1500 :return: RMS distance (deg). 1501 """ 1502 import math 1503 1504 # Get internal coordinates of both objects 1505 ic1 = self.torsion_array 1506 ic2 = other.torsion_array 1507 1508 # Compute the sum of squared differences between corresponding internal coordinates 1509 totsq = sum((p1 - p2) ** 2 for p1, p2 in zip(ic1, ic2)) 1510 # Compute the mean of the squared distances 1511 totsq /= len(ic1) 1512 1513 # Take the square root of the mean to get the RMS distance 1514 return math.sqrt(totsq) 1515 1516 def get_chains(self) -> tuple: 1517 """ 1518 Return the proximal and distal chain IDs for the Disulfide. 1519 1520 :return: tuple (proximal, distal) chain IDs 1521 """ 1522 prox = self.proximal_chain 1523 dist = self.distal_chain 1524 return tuple(prox, dist) 1525 1526 def get_permissive(self) -> bool: 1527 """ 1528 Return the Permissive flag state. (Used in PDB parsing) 1529 1530 :return: Permissive state 1531 """ 1532 return self.PERMISIVE 1533 1534 def get_full_id(self) -> tuple: 1535 """ 1536 Return the Disulfide full IDs (Used with BIO.PDB) 1537 1538 :return: Disulfide full IDs 1539 """ 1540 return (self.proximal_residue_fullid, self.distal_residue_fullid) 1541 1542 def initialize_disulfide_from_chain( 1543 self, chain1, chain2, proximal, distal, resolution, quiet=True 1544 ) -> None: 1545 """ 1546 Initialize a new Disulfide object with atomic coordinates from 1547 the proximal and distal coordinates, typically taken from a PDB file. 1548 This routine is primarily used internally when building the compressed 1549 database. 1550 1551 :param chain1: list of Residues in the model, eg: chain = model['A'] 1552 :param chain2: list of Residues in the model, eg: chain = model['A'] 1553 :param proximal: proximal residue sequence ID 1554 :param distal: distal residue sequence ID 1555 :param resolution: structure resolution 1556 :param quiet: Quiet or noisy parsing, defaults to True 1557 :raises DisulfideConstructionWarning: Raised when not parsed correctly 1558 """ 1559 id = chain1.get_full_id()[0] 1560 self.pdb_id = id 1561 1562 chi1 = chi2 = chi3 = chi4 = chi5 = _ANG_INIT 1563 1564 prox = int(proximal) 1565 dist = int(distal) 1566 1567 prox_residue = chain1[prox] 1568 dist_residue = chain2[dist] 1569 1570 if prox_residue.get_resname() != "CYS" or dist_residue.get_resname() != "CYS": 1571 print( 1572 f"build_disulfide() requires CYS at both residues: {prox} {prox_residue.get_resname()} {dist} {dist_residue.get_resname()} Chain: {prox_residue.get_segid()}" 1573 ) 1574 1575 # set the objects proximal and distal values 1576 self.set_resnum(proximal, distal) 1577 1578 if resolution is not None: 1579 self.resolution = resolution 1580 1581 self.proximal_chain = chain1.get_id() 1582 self.distal_chain = chain2.get_id() 1583 1584 self.proximal_residue_fullid = prox_residue.get_full_id() 1585 self.distal_residue_fullid = dist_residue.get_full_id() 1586 1587 if quiet: 1588 warnings.filterwarnings("ignore", category=DisulfideConstructionWarning) 1589 else: 1590 warnings.simplefilter("always") 1591 1592 # grab the coordinates for the proximal and distal residues as vectors 1593 # so we can do math on them later 1594 # proximal residue 1595 1596 try: 1597 n1 = prox_residue["N"].get_vector() 1598 ca1 = prox_residue["CA"].get_vector() 1599 c1 = prox_residue["C"].get_vector() 1600 o1 = prox_residue["O"].get_vector() 1601 cb1 = prox_residue["CB"].get_vector() 1602 sg1 = prox_residue["SG"].get_vector() 1603 1604 except Exception: 1605 raise DisulfideConstructionWarning( 1606 f"Invalid or missing coordinates for proximal residue {proximal}" 1607 ) from None 1608 1609 # distal residue 1610 try: 1611 n2 = dist_residue["N"].get_vector() 1612 ca2 = dist_residue["CA"].get_vector() 1613 c2 = dist_residue["C"].get_vector() 1614 o2 = dist_residue["O"].get_vector() 1615 cb2 = dist_residue["CB"].get_vector() 1616 sg2 = dist_residue["SG"].get_vector() 1617 1618 except Exception: 1619 raise DisulfideConstructionWarning( 1620 f"Invalid or missing coordinates for distal residue {distal}" 1621 ) from None 1622 1623 # previous residue and next residue - optional, used for phi, psi calculations 1624 try: 1625 prevprox = chain1[prox - 1] 1626 nextprox = chain1[prox + 1] 1627 1628 prevdist = chain2[dist - 1] 1629 nextdist = chain2[dist + 1] 1630 1631 cprev_prox = prevprox["C"].get_vector() 1632 nnext_prox = nextprox["N"].get_vector() 1633 1634 cprev_dist = prevdist["C"].get_vector() 1635 nnext_dist = nextdist["N"].get_vector() 1636 1637 # compute phi, psi for prox and distal 1638 self.phiprox = np.degrees(calc_dihedral(cprev_prox, n1, ca1, c1)) 1639 self.psiprox = np.degrees(calc_dihedral(n1, ca1, c1, nnext_prox)) 1640 self.phidist = np.degrees(calc_dihedral(cprev_dist, n2, ca2, c2)) 1641 self.psidist = np.degrees(calc_dihedral(n2, ca2, c2, nnext_dist)) 1642 1643 except Exception: 1644 mess = f"Missing coords for: {id} {prox-1} or {dist+1} for SS {proximal}-{distal}" 1645 cprev_prox = nnext_prox = cprev_dist = nnext_dist = Vector(-1.0, -1.0, -1.0) 1646 self.missing_atoms = True 1647 warnings.warn(mess, DisulfideConstructionWarning) 1648 1649 # update the positions and conformation 1650 self.set_positions( 1651 n1, 1652 ca1, 1653 c1, 1654 o1, 1655 cb1, 1656 sg1, 1657 n2, 1658 ca2, 1659 c2, 1660 o2, 1661 cb2, 1662 sg2, 1663 cprev_prox, 1664 nnext_prox, 1665 cprev_dist, 1666 nnext_dist, 1667 ) 1668 1669 # calculate and set the disulfide dihedral angles 1670 self.chi1 = np.degrees(calc_dihedral(n1, ca1, cb1, sg1)) 1671 self.chi2 = np.degrees(calc_dihedral(ca1, cb1, sg1, sg2)) 1672 self.chi3 = np.degrees(calc_dihedral(cb1, sg1, sg2, cb2)) 1673 self.chi4 = np.degrees(calc_dihedral(sg1, sg2, cb2, ca2)) 1674 self.chi5 = np.degrees(calc_dihedral(sg2, cb2, ca2, n2)) 1675 self.rho = np.degrees(calc_dihedral(n1, ca1, ca2, n2)) 1676 1677 self.ca_distance = distance3d(self.ca_prox, self.ca_dist) 1678 self.cb_distance = distance3d(self.cb_prox, self.cb_dist) 1679 self.torsion_array = np.array( 1680 (self.chi1, self.chi2, self.chi3, self.chi4, self.chi5) 1681 ) 1682 self.compute_torsion_length() 1683 1684 # calculate and set the SS bond torsional energy 1685 self.compute_torsional_energy() 1686 1687 # compute and set the local coordinates 1688 self.compute_local_coords() 1689 1690 def internal_coords(self) -> np.array: 1691 """ 1692 Return the internal coordinates for the Disulfide. 1693 If there are missing atoms the extra atoms for the proximal 1694 and distal N and C are set to [0,0,0]. This is needed for the center of 1695 mass calculations, used when rendering. 1696 1697 :return: Array containing the coordinates, [16][3]. 1698 """ 1699 1700 # if we don't have the prior and next atoms we initialize those 1701 # atoms to the origin so as to not effect the center of mass calculations 1702 if self.missing_atoms: 1703 res_array = np.array( 1704 ( 1705 self._n_prox.get_array(), 1706 self._ca_prox.get_array(), 1707 self._c_prox.get_array(), 1708 self._o_prox.get_array(), 1709 self._cb_prox.get_array(), 1710 self._sg_prox.get_array(), 1711 self._n_dist.get_array(), 1712 self._ca_dist.get_array(), 1713 self._c_dist.get_array(), 1714 self._o_dist.get_array(), 1715 self._cb_dist.get_array(), 1716 self._sg_dist.get_array(), 1717 [0, 0, 0], 1718 [0, 0, 0], 1719 [0, 0, 0], 1720 [0, 0, 0], 1721 ) 1722 ) 1723 else: 1724 res_array = np.array( 1725 ( 1726 self._n_prox.get_array(), 1727 self._ca_prox.get_array(), 1728 self._c_prox.get_array(), 1729 self._o_prox.get_array(), 1730 self._cb_prox.get_array(), 1731 self._sg_prox.get_array(), 1732 self._n_dist.get_array(), 1733 self._ca_dist.get_array(), 1734 self._c_dist.get_array(), 1735 self._o_dist.get_array(), 1736 self._cb_dist.get_array(), 1737 self._sg_dist.get_array(), 1738 self._c_prev_prox.get_array(), 1739 self._n_next_prox.get_array(), 1740 self._c_prev_dist.get_array(), 1741 self._n_next_dist.get_array(), 1742 ) 1743 ) 1744 return res_array 1745 1746 def internal_coords_res(self, resnumb) -> np.array: 1747 """ 1748 Return the internal coordinates for the internal coordinates of 1749 the given Disulfide. Missing atoms are not included. 1750 1751 :param resnumb: Residue number for disulfide 1752 :raises DisulfideConstructionWarning: Warning raised if the residue number is invalid 1753 :return: Array containing the internal coordinates for the disulfide 1754 """ 1755 res_array = np.zeros(shape=(6, 3)) 1756 1757 if resnumb == self.proximal: 1758 res_array = np.array( 1759 ( 1760 self._n_prox.get_array(), 1761 self._ca_prox.get_array(), 1762 self._c_prox.get_array(), 1763 self._o_prox.get_array(), 1764 self._cb_prox.get_array(), 1765 self._sg_prox.get_array(), 1766 ) 1767 ) 1768 return res_array 1769 1770 elif resnumb == self.distal: 1771 res_array = np.array( 1772 ( 1773 self._n_dist.get_array(), 1774 self._ca_dist.get_array(), 1775 self._c_dist.get_array(), 1776 self._o_dist.get_array(), 1777 self._cb_dist.get_array(), 1778 self._sg_dist.get_array(), 1779 ) 1780 ) 1781 return res_array 1782 else: 1783 mess = f"-> Disulfide.internal_coords(): Invalid argument. \ 1784 Unable to find residue: {resnumb} " 1785 raise DisulfideConstructionWarning(mess) 1786 1787 def make_movie( 1788 self, style="sb", fname="ssbond.mp4", verbose=False, steps=360 1789 ) -> None: 1790 """ 1791 Create an animation for ```self``` rotating one revolution about the Y axis, 1792 in the given ```style```, saving to ```filename```. 1793 1794 :param style: Rendering style, defaults to 'sb', one of: 1795 * 'sb' - split bonds 1796 * 'bs' - ball and stick 1797 * 'cpk' - CPK style 1798 * 'pd' - Proximal/Distal style - Red=proximal, Green=Distal 1799 * 'plain' - boring single color 1800 1801 :param fname: Output filename, defaults to ```ssbond.mp4``` 1802 :param verbose: Verbosity, defaults to False 1803 :param steps: Number of steps for one complete rotation, defaults to 360. 1804 """ 1805 src = self.pdb_id 1806 name = self.name 1807 enrg = self.energy 1808 1809 title = f"{src} {name}: {self.proximal}{self.proximal_chain}-{self.distal}{self.distal_chain}: {enrg:.2f} kcal/mol, Cα: {self.ca_distance:.2f} Å, Tors: {self.torsion_length:.2f}" 1810 1811 if verbose: 1812 print(f"Rendering animation to {fname}...") 1813 1814 pl = pv.Plotter(window_size=WINSIZE, off_screen=True) 1815 pl.open_movie(fname) 1816 path = pl.generate_orbital_path(n_points=steps) 1817 1818 # 1819 # pl.add_title(title=title, font_size=FONTSIZE) 1820 pl.enable_anti_aliasing("msaa") 1821 # pl.add_camera_orientation_widget() 1822 pl = self._render( 1823 pl, 1824 style=style, 1825 bondcolor=BOND_COLOR, 1826 bs_scale=BS_SCALE, 1827 spec=SPECULARITY, 1828 specpow=SPEC_POWER, 1829 ) 1830 pl.reset_camera() 1831 pl.orbit_on_path(path, write_frames=True) 1832 pl.close() 1833 1834 if verbose: 1835 print(f"Saved mp4 animation to: {fname}") 1836 1837 def pprint(self) -> None: 1838 """ 1839 Pretty print general info for the Disulfide 1840 """ 1841 s1 = self.repr_ss_info() 1842 s2 = self.repr_ss_ca_dist() 1843 s3 = self.repr_ss_conformation() 1844 s4 = self.repr_ss_torsion_length() 1845 res = f"{s1} \n{s3} \n{s2} \n{s4}>" 1846 print(res) 1847 1848 def pprint_all(self) -> None: 1849 """ 1850 Pretty print all info for the Disulfide 1851 """ 1852 s1 = self.repr_ss_info() + "\n" 1853 s2 = self.repr_ss_coords() 1854 s3 = self.repr_ss_local_coords() 1855 s4 = self.repr_ss_conformation() 1856 s5 = self.repr_chain_ids() 1857 s6 = self.repr_ss_ca_dist() 1858 s7 = self.repr_ss_torsion_length() 1859 1860 res = f"{s1} {s5} {s2} {s3} {s4}\n {s6}\n {s7}>" 1861 1862 print(res) 1863 1864 # repr functions. The class is large, so I split it up into sections 1865 def repr_ss_info(self) -> str: 1866 """ 1867 Representation for the Disulfide class 1868 """ 1869 s1 = f"<Disulfide {self.name}, Source: {self.pdb_id}, Resolution: {self.resolution} Å" 1870 return s1 1871 1872 def repr_ss_coords(self) -> str: 1873 """ 1874 Representation for Disulfide coordinates 1875 """ 1876 s2 = f"\nProximal Coordinates:\n N: {self.n_prox}\n Cα: {self.ca_prox}\n C: {self.c_prox}\n O: {self.o_prox}\n Cβ: {self.cb_prox}\n Sγ: {self.sg_prox}\n Cprev {self.c_prev_prox}\n Nnext: {self.n_next_prox}\n" 1877 s3 = f"Distal Coordinates:\n N: {self.n_dist}\n Cα: {self.ca_dist}\n C: {self.c_dist}\n O: {self.o_dist}\n Cβ: {self.cb_dist}\n Sγ: {self.sg_dist}\n Cprev {self.c_prev_dist}\n Nnext: {self.n_next_dist}\n\n" 1878 stot = f"{s2} {s3}" 1879 return stot 1880 1881 def repr_ss_conformation(self) -> str: 1882 """ 1883 Representation for Disulfide conformation 1884 """ 1885 s4 = f"Χ1-Χ5: {self.chi1:.2f}°, {self.chi2:.2f}°, {self.chi3:.2f}°, {self.chi4:.2f}° {self.chi5:.2f}°, {self.rho:.2f}°, {self.energy:.2f} kcal/mol" 1886 stot = f"{s4}" 1887 return stot 1888 1889 def repr_ss_local_coords(self) -> str: 1890 """ 1891 Representation for the Disulfide internal coordinates. 1892 """ 1893 s2i = f"Proximal Internal Coords:\n N: {self._n_prox}\n Cα: {self._ca_prox}\n C: {self._c_prox}\n O: {self._o_prox}\n Cβ: {self._cb_prox}\n Sγ: {self._sg_prox}\n Cprev {self.c_prev_prox}\n Nnext: {self.n_next_prox}\n" 1894 s3i = f"Distal Internal Coords:\n N: {self._n_dist}\n Cα: {self._ca_dist}\n C: {self._c_dist}\n O: {self._o_dist}\n Cβ: {self._cb_dist}\n Sγ: {self._sg_dist}\n Cprev {self.c_prev_dist}\n Nnext: {self.n_next_dist}\n" 1895 stot = f"{s2i}{s3i}" 1896 return stot 1897 1898 def repr_ss_chain_ids(self) -> str: 1899 """ 1900 Representation for Disulfide chain IDs 1901 """ 1902 return f"Proximal Chain fullID: <{self.proximal_residue_fullid}> Distal Chain fullID: <{self.distal_residue_fullid}>" 1903 1904 def repr_ss_ca_dist(self) -> str: 1905 """ 1906 Representation for Disulfide Ca distance 1907 """ 1908 s1 = f"Cα Distance: {self.ca_distance:.2f} Å" 1909 return s1 1910 1911 def repr_ss_cb_dist(self) -> str: 1912 """ 1913 Representation for Disulfide Ca distance 1914 """ 1915 s1 = f"Cβ Distance: {self.cb_distance:.2f} Å" 1916 return s1 1917 1918 def repr_ss_torsion_length(self) -> str: 1919 """ 1920 Representation for Disulfide torsion length 1921 """ 1922 s1 = f"Torsion length: {self.torsion_length:.2f} deg" 1923 return s1 1924 1925 def repr_all(self) -> str: 1926 """ 1927 Return a string representation for all Disulfide information 1928 contained in self. 1929 """ 1930 1931 s1 = self.repr_ss_info() + "\n" 1932 s2 = self.repr_ss_coords() 1933 s3 = self.repr_ss_local_coords() 1934 s4 = self.repr_ss_conformation() 1935 s5 = self.repr_chain_ids() 1936 s6 = self.repr_ss_ca_dist() 1937 s8 = self.repr_ss_cb_dist() 1938 s7 = self.repr_ss_torsion_length() 1939 1940 res = f"{s1} {s5} {s2} {s3} {s4} {s6} {s7} {s8}>" 1941 return res 1942 1943 def repr_compact(self) -> str: 1944 """ 1945 Return a compact representation of the Disulfide object 1946 :return: string 1947 """ 1948 return f"{self.repr_ss_info()} {self.repr_ss_conformation()}" 1949 1950 def repr_conformation(self) -> str: 1951 """ 1952 Return a string representation of the Disulfide object's conformation. 1953 :return: string 1954 """ 1955 return f"{self.repr_ss_conformation()}" 1956 1957 def repr_coords(self) -> str: 1958 """ 1959 Return a string representation of the Disulfide object's coordinates. 1960 :return: string 1961 """ 1962 return f"{self.repr_ss_coords()}" 1963 1964 def repr_internal_coords(self) -> str: 1965 """ 1966 Return a string representation of the Disulfide object's internal coordinaes. 1967 :return: string 1968 """ 1969 return f"{self.repr_ss_local_coords()}" 1970 1971 def repr_chain_ids(self) -> str: 1972 """ 1973 Return a string representation of the Disulfide object's chain ids. 1974 :return: string 1975 """ 1976 return f"{self.repr_ss_chain_ids()}" 1977 1978 def compute_rho(self) -> float: 1979 self.rho = calc_dihedral(self.n_prox, self.ca_prox, self.ca_dist, self.n_dist) 1980 return self.rho 1981 1982 def reset(self) -> None: 1983 """ 1984 Resets the disulfide object to its initial state. All distances, 1985 angles and positions are reset. The name is unchanged. 1986 """ 1987 self.__init__(self) 1988 1989 def same_chains(self) -> bool: 1990 """ 1991 Function checks if the Disulfide is cross-chain or not. 1992 1993 Returns 1994 ------- 1995 bool \n 1996 True if the proximal and distal residues are on the same chains, 1997 False otherwise. 1998 """ 1999 2000 (prox, dist) = self.get_chains() 2001 return prox == dist 2002 2003 def screenshot( 2004 self, 2005 single=True, 2006 style="sb", 2007 fname="ssbond.png", 2008 verbose=False, 2009 shadows=False, 2010 light=True, 2011 ) -> None: 2012 """ 2013 Create and save a screenshot of the Disulfide in the given style 2014 and filename 2015 2016 :param single: Display a single vs panel view, defaults to True 2017 :param style: Rendering style, one of: 2018 * 'sb' - split bonds 2019 * 'bs' - ball and stick 2020 * 'cpk' - CPK style 2021 * 'pd' - Proximal/Distal style - Red=proximal, Green=Distal 2022 * 'plain' - boring single color, 2023 :param fname: output filename,, defaults to 'ssbond.png' 2024 :param verbose: Verbosit, defaults to False 2025 :param shadows: Enable shadows, defaults to False 2026 """ 2027 src = self.pdb_id 2028 name = self.name 2029 enrg = self.energy 2030 2031 title = f"{src} {name}: {self.proximal}{self.proximal_chain}-{self.distal}{self.distal_chain}: {enrg:.2f} kcal/mol, Cα: {self.ca_distance:.2f} Å, Cβ: {self.cb_distance:.2f} Å, Tors: {self.torsion_length:.2f}" 2032 2033 if light: 2034 pv.set_plot_theme("document") 2035 else: 2036 pv.set_plot_theme("dark") 2037 2038 if verbose: 2039 print(f"-> screenshot(): Rendering screenshot to file {fname}") 2040 2041 if single: 2042 pl = pv.Plotter(window_size=WINSIZE) 2043 # pl.add_title(title=title, font_size=FONTSIZE) 2044 pl.enable_anti_aliasing("msaa") 2045 # pl.add_camera_orientation_widget() 2046 self._render( 2047 pl, 2048 style=style, 2049 bondcolor=BOND_COLOR, 2050 bs_scale=BS_SCALE, 2051 spec=SPECULARITY, 2052 specpow=SPEC_POWER, 2053 ) 2054 pl.reset_camera() 2055 if shadows: 2056 pl.enable_shadows() 2057 2058 pl.show(auto_close=False) 2059 pl.screenshot(fname) 2060 pl.clear() 2061 2062 else: 2063 pl = pv.Plotter(window_size=WINSIZE, shape=(2, 2)) 2064 pl.subplot(0, 0) 2065 2066 # pl.add_title(title=title, font_size=FONTSIZE) 2067 pl.enable_anti_aliasing("msaa") 2068 2069 # pl.add_camera_orientation_widget() 2070 self._render( 2071 pl, 2072 style="cpk", 2073 bondcolor=BOND_COLOR, 2074 bs_scale=BS_SCALE, 2075 spec=SPECULARITY, 2076 specpow=SPEC_POWER, 2077 ) 2078 2079 pl.subplot(0, 1) 2080 # pl.add_title(title=title, font_size=FONTSIZE) 2081 self._render( 2082 pl, 2083 style="pd", 2084 bondcolor=BOND_COLOR, 2085 bs_scale=BS_SCALE, 2086 spec=SPECULARITY, 2087 specpow=SPEC_POWER, 2088 ) 2089 2090 pl.subplot(1, 0) 2091 # pl.add_title(title=title, font_size=FONTSIZE) 2092 self._render( 2093 pl, 2094 style="bs", 2095 bondcolor=BOND_COLOR, 2096 bs_scale=BS_SCALE, 2097 spec=SPECULARITY, 2098 specpow=SPEC_POWER, 2099 ) 2100 2101 pl.subplot(1, 1) 2102 # pl.add_title(title=title, font_size=FONTSIZE) 2103 self._render( 2104 pl, 2105 style="sb", 2106 bondcolor=BOND_COLOR, 2107 bs_scale=BS_SCALE, 2108 spec=SPECULARITY, 2109 specpow=SPEC_POWER, 2110 ) 2111 2112 pl.link_views() 2113 pl.reset_camera() 2114 if shadows: 2115 pl.enable_shadows() 2116 2117 pl.show(auto_close=False) 2118 pl.screenshot(fname) 2119 2120 if verbose: 2121 print(f"Saved: {fname}") 2122 2123 def save_meshes_as_stl(self, meshes, filename) -> None: 2124 """Save a list of meshes as a single STL file. 2125 2126 Args: 2127 meshes (list): List of pyvista mesh objects to save. 2128 filename (str): Path to save the STL file to. 2129 """ 2130 merged_mesh = pv.UnstructuredGrid() 2131 for mesh in meshes: 2132 merged_mesh += mesh 2133 merged_mesh.save(filename) 2134 2135 def export(self, style="sb", verbose=True, fname="ssbond_plt") -> None: 2136 """ 2137 Create and save a screenshot of the Disulfide in the given style and filename. 2138 2139 :param single: Display a single vs panel view, defaults to True 2140 :param style: Rendering style, one of: 2141 * 'sb' - split bonds 2142 * 'bs' - ball and stick 2143 * 'cpk' - CPK style 2144 * 'pd' - Proximal/Distal style - Red=proximal, Green=Distal 2145 * 'plain' - boring single color, 2146 2147 :param fname: output filename,, defaults to 'ssbond.stl' 2148 :param verbose: Verbosit, defaults to False 2149 """ 2150 2151 if verbose: 2152 print(f"-> screenshot(): Rendering screenshot to file {fname}") 2153 2154 pl = pv.PolyData() 2155 2156 self._plot( 2157 pl, 2158 style=style, 2159 bondcolor=BOND_COLOR, 2160 bs_scale=BS_SCALE, 2161 spec=SPECULARITY, 2162 specpow=SPEC_POWER, 2163 ) 2164 2165 self.save_meshes_as_stl(pl, fname) 2166 2167 return 2168 2169 def set_permissive(self, perm: bool) -> None: 2170 """ 2171 Set PERMISSIVE flag for Disulfide parsing. 2172 2173 :return: None 2174 """ 2175 2176 self.PERMISSIVE = perm 2177 2178 def set_positions( 2179 self, 2180 n_prox: Vector, 2181 ca_prox: Vector, 2182 c_prox: Vector, 2183 o_prox: Vector, 2184 cb_prox: Vector, 2185 sg_prox: Vector, 2186 n_dist: Vector, 2187 ca_dist: Vector, 2188 c_dist: Vector, 2189 o_dist: Vector, 2190 cb_dist: Vector, 2191 sg_dist: Vector, 2192 c_prev_prox: Vector, 2193 n_next_prox: Vector, 2194 c_prev_dist: Vector, 2195 n_next_dist: Vector, 2196 ) -> None: 2197 """ 2198 Set the atomic coordinates for all atoms in the Disulfide object. 2199 2200 :param n_prox: Proximal N position 2201 :param ca_prox: Proximal Cα position 2202 :param c_prox: Proximal C' position 2203 :param o_prox: Proximal O position 2204 :param cb_prox: Proximal Cβ position 2205 :param sg_prox: Proximal Sγ position 2206 :param n_dist: Distal N position 2207 :param ca_dist: Distal Cα position 2208 :param c_dist: Distal C' position 2209 :param o_dist: Distal O position 2210 :param cb_dist: Distal Cβ position 2211 :param sg_dist: Distal Sγ position 2212 :param c_prev_prox: Proximal previous C' 2213 :param n_next_prox: Proximal next N 2214 :param c_prev_dist: Distal previous C' 2215 :param n_next_dist: Distal next N 2216 """ 2217 2218 # deep copy 2219 self.n_prox = n_prox.copy() 2220 self.ca_prox = ca_prox.copy() 2221 self.c_prox = c_prox.copy() 2222 self.o_prox = o_prox.copy() 2223 self.cb_prox = cb_prox.copy() 2224 self.sg_prox = sg_prox.copy() 2225 self.sg_dist = sg_dist.copy() 2226 self.cb_dist = cb_dist.copy() 2227 self.ca_dist = ca_dist.copy() 2228 self.n_dist = n_dist.copy() 2229 self.c_dist = c_dist.copy() 2230 self.o_dist = o_dist.copy() 2231 2232 self.c_prev_prox = c_prev_prox.copy() 2233 self.n_next_prox = n_next_prox.copy() 2234 self.c_prev_dist = c_prev_dist.copy() 2235 self.n_next_dist = n_next_dist.copy() 2236 2237 def set_dihedrals( 2238 self, chi1: float, chi2: float, chi3: float, chi4: float, chi5: float 2239 ) -> None: 2240 """ 2241 Set the disulfide's dihedral angles, Chi1-Chi5. -180 - 180 degrees. 2242 2243 :param chi1: Chi1 2244 :param chi2: Chi2 2245 :param chi3: Chi3 2246 :param chi4: Chi4 2247 :param chi5: Chi5 2248 """ 2249 self.chi1 = chi1 2250 self.chi2 = chi2 2251 self.chi3 = chi3 2252 self.chi4 = chi4 2253 self.chi5 = chi5 2254 self.torsion_array = np.array([chi1, chi2, chi3, chi4, chi5]) 2255 self.compute_torsional_energy() 2256 self.compute_torsion_length() 2257 2258 def set_name(self, namestr="Disulfide") -> None: 2259 """ 2260 Set the Disulfide's name. 2261 2262 :param namestr: Name, by default "Disulfide" 2263 """ 2264 self.name = namestr 2265 2266 def set_resnum(self, proximal: int, distal: int) -> None: 2267 """ 2268 Set the proximal and residue numbers for the Disulfide. 2269 2270 :param proximal: Proximal residue number 2271 :param distal: Distal residue number 2272 """ 2273 self.proximal = proximal 2274 self.distal = distal 2275 2276 def compute_torsion_length(self) -> float: 2277 """ 2278 Compute the 5D Euclidean length of the Disulfide object. Update the disulfide internal state. 2279 2280 :return: Torsion length (Degrees) 2281 """ 2282 # Use numpy array to compute element-wise square 2283 tors2 = np.square(self.torsion_array) 2284 2285 # Compute the sum of squares using numpy's sum function 2286 dist = math.sqrt(np.sum(tors2)) 2287 2288 # Update the internal state 2289 self.torsion_length = dist 2290 2291 return dist 2292 2293 def Torsion_Distance(self, other) -> float: 2294 """ 2295 Calculate the 5D Euclidean distance between ```self``` and another Disulfide 2296 object. This is used to compare Disulfide Bond torsion angles to 2297 determine their torsional similarity via a 5-Dimensional Euclidean distance metric. 2298 2299 :param other: Comparison Disulfide 2300 :raises ProteusPyWarning: Warning if ```other``` is not a Disulfide object 2301 :return: Euclidean distance (Degrees) between ```self``` and ```other```. 2302 """ 2303 2304 from proteusPy.ProteusPyWarning import ProteusPyWarning 2305 2306 # Check length of torsion arrays 2307 if len(self.torsion_array) != 5 or len(other.torsion_array) != 5: 2308 raise ProteusPyWarning( 2309 "--> Torsion_Distance() requires vectors of length 5!" 2310 ) 2311 2312 # Convert to numpy arrays and add 180 to each element 2313 p1 = np.array(self.torsion_array) + 180.0 2314 p2 = np.array(other.torsion_array) + 180.0 2315 2316 # Compute the 5D Euclidean distance using numpy's linalg.norm function 2317 dist = np.linalg.norm(p1 - p2) 2318 2319 return dist 2320 2321 def Torsion_neighbors(self, others, cutoff) -> DisulfideList: 2322 """ 2323 Return a list of Disulfides within the angular cutoff in the others list. 2324 This routine is used to find Disulfides having the same torsion length 2325 within the others list. This is used to find families of Disulfides with 2326 similar conformations. Assumes self is properly initialized. 2327 2328 *NB* The routine will not distinguish between +/- 2329 dihedral angles. *i.e.* [-60, -60, -90, -60, -60] would have the same 2330 torsion length as [60, 60, 90, 60, 60], two clearly different structures. 2331 2332 :param others: ```DisulfideList``` to search 2333 :param cutoff: Dihedral angle degree cutoff 2334 :return: DisulfideList within the cutoff 2335 2336 Example: 2337 In this example we load the disulfide database subset, find the disulfides with 2338 the lowest and highest energies, and then find the nearest conformational neighbors. 2339 Finally, we display the neighbors overlaid against a common reference frame. 2340 2341 >>> from proteusPy import * 2342 >>> from proteusPy.DisulfideLoader import Load_PDB_SS 2343 >>> from proteusPy.DisulfideList import DisulfideList 2344 >>> from proteusPy.Disulfide import Disulfide 2345 >>> PDB_SS = Load_PDB_SS(verbose=False, subset=True) 2346 >>> ss_list = DisulfideList([], 'tmp') 2347 2348 We point to the complete list to search for lowest and highest energies. 2349 >>> sslist = PDB_SS.SSList 2350 >>> ssmin_enrg, ssmax_enrg = PDB_SS.SSList.minmax_energy 2351 2352 Make an empty list and find the nearest neighbors within 10 degrees avg RMS in 2353 sidechain dihedral angle space. 2354 2355 >>> low_energy_neighbors = DisulfideList([],'Neighbors') 2356 >>> low_energy_neighbors = ssmin_enrg.Torsion_neighbors(sslist, 10) 2357 2358 Display the number found, and then display them overlaid onto their common reference frame. 2359 2360 >>> tot = low_energy_neighbors.length 2361 >>> print(f'Neighbors: {tot}') 2362 Neighbors: 2 2363 >>> low_energy_neighbors.display_overlay() 2364 2365 """ 2366 res = [ss for ss in others if self.Torsion_Distance(ss) <= cutoff] 2367 return DisulfideList(res, "neighbors") 2368 2369 def torsion_to_sixclass(self) -> str: 2370 """ 2371 Return the sextant class string for ``self``. 2372 2373 :raises DisulfideIOException: _description_ 2374 :return: Sextant string 2375 """ 2376 from proteusPy.DisulfideClasses import get_sixth_quadrant 2377 2378 tors = self.torsion_array 2379 res = [get_sixth_quadrant(x) for x in tors] 2380 return "".join([str(r) for r in res]) 2381 2382 2383# Class defination ends 2384 2385 2386def parse_ssbond_header_rec(ssbond_dict: dict) -> list: 2387 """ 2388 Parse the SSBOND dict returned by parse_pdb_header. 2389 NB: Requires EGS-Modified BIO.parse_pdb_header.py. 2390 This is used internally. 2391 2392 :param ssbond_dict: the input SSBOND dict 2393 :return: A list of tuples representing the proximal, 2394 distal residue ids for the Disulfide. 2395 2396 """ 2397 disulfide_list = [] 2398 for ssb in ssbond_dict.items(): 2399 disulfide_list.append(ssb[1]) 2400 2401 return disulfide_list 2402 2403 2404# 2405# Function reads a comma separated list of PDB IDs and download the corresponding 2406# .ent files to the PDB_DIR global. 2407# Used to download the list of proteins containing at least one SS bond 2408# with the ID list generated from: http://www.rcsb.org/ 2409# 2410 2411 2412def Download_Disulfides( 2413 pdb_home=PDB_DIR, model_home=MODEL_DIR, verbose=False, reset=False 2414) -> None: 2415 """ 2416 Read a comma separated list of PDB IDs and download them 2417 to the pdb_home path. 2418 2419 This utility function is used to download proteins containing at least 2420 one SS bond with the ID list generated from: http://www.rcsb.org/. 2421 2422 This is the primary data loader for the proteusPy Disulfide 2423 analysis package. The list of IDs represents files in the 2424 RCSB containing > 1 disulfide bond, and it contains 2425 over 39000 structures. The total download takes about 12 hours. The 2426 function keeps track of downloaded files so it's possible to interrupt and 2427 restart the download without duplicating effort. 2428 2429 :param pdb_home: Path for downloaded files, defaults to PDB_DIR 2430 :param model_home: Path for extracted data, defaults to MODEL_DIR 2431 :param verbose: Verbosity, defaults to False 2432 :param reset: Reset the downloaded files index. Used to restart the download. 2433 :raises DisulfideIOException: I/O error raised when the PDB file is not found. 2434 """ 2435 import os 2436 import time 2437 2438 start = time.time() 2439 donelines = [] 2440 SS_done = [] 2441 ssfile = None 2442 2443 cwd = os.getcwd() 2444 os.chdir(pdb_home) 2445 2446 pdblist = PDBList(pdb=pdb_home, verbose=verbose) 2447 ssfilename = f"{model_home}{SS_ID_FILE}" 2448 print(ssfilename) 2449 2450 # list of IDs containing >1 SSBond record 2451 try: 2452 ssfile = open(ssfilename) 2453 Line = ssfile.readlines() 2454 except Exception: 2455 raise DisulfideIOException(f"Cannot open file: {ssfile}") 2456 2457 for line in Line: 2458 entries = line.split(",") 2459 2460 print(f"Found: {len(entries)} entries") 2461 completed = {"xxx"} # set to keep track of downloaded 2462 2463 # file to track already downloaded entries. 2464 if reset is True: 2465 completed_file = open(f"{model_home}ss_completed.txt", "w") 2466 donelines = [] 2467 SS_DONE = [] 2468 else: 2469 completed_file = open(f"{model_home}ss_completed.txt", "w+") 2470 donelines = completed_file.readlines() 2471 2472 if len(donelines) > 0: 2473 for dl in donelines[0]: 2474 # create a list of pdb id already downloaded 2475 SS_done = dl.split(",") 2476 2477 count = len(SS_done) - 1 2478 completed.update(SS_done) # update the completed set with what's downloaded 2479 2480 # Loop over all entries, 2481 pbar = tqdm(entries, ncols=PBAR_COLS) 2482 for entry in pbar: 2483 pbar.set_postfix({"Entry": entry}) 2484 if entry not in completed: 2485 if pdblist.retrieve_pdb_file(entry, file_format="pdb", pdir=pdb_home): 2486 completed.update(entry) 2487 completed_file.write(f"{entry},") 2488 count += 1 2489 2490 completed_file.close() 2491 2492 end = time.time() 2493 elapsed = end - start 2494 2495 print(f"Overall files processed: {count}") 2496 print(f"Complete. Elapsed time: {datetime.timedelta(seconds=elapsed)} (h:m:s)") 2497 os.chdir(cwd) 2498 return 2499 2500 2501# Function extracts the disulfide bonds from the PDB files and creates the .pkl files 2502# needed for the proteusPy.DisulfideLoader.DisulfideLoader class. 2503 2504 2505def Extract_Disulfides( 2506 numb=-1, 2507 verbose=False, 2508 quiet=True, 2509 pdbdir=PDB_DIR, 2510 datadir=MODEL_DIR, 2511 picklefile=SS_PICKLE_FILE, 2512 torsionfile=SS_TORSIONS_FILE, 2513 problemfile=PROBLEM_ID_FILE, 2514 dictfile=SS_DICT_PICKLE_FILE, 2515 dist_cutoff=-1.0, 2516) -> None: 2517 """ 2518 Read the PDB files contained in ``pdbdir`` and create the .pkl files needed for the 2519 proteusPy.DisulfideLoader.DisulfideLoader class. 2520 The ```Disulfide``` objects are contained in a ```DisulfideList``` object and 2521 ```Dict``` within these files. In addition, .csv files containing all of 2522 the torsions for the disulfides and problem IDs are written. The optional 2523 ```dist_cutoff``` allows for removal of Disufides whose Cα-Cα distance is > 2524 than the cutoff value. If it's -1.0 then the function keeps all Disulfides. 2525 2526 :param numb: number of entries to process, defaults to all 2527 :param verbose: more messages 2528 :param quiet: turns off DisulfideConstruction warnings 2529 :param pdbdir: path to PDB files 2530 :param datadir: path to resulting .pkl files 2531 :param picklefile: name of the disulfide .pkl file 2532 :param torsionfile: name of the disulfide torsion file .csv created 2533 :param problemfile: name of the .csv file containing problem ids 2534 :param dictfile: name of the .pkl file 2535 :param dist_cutoff: Ca distance cutoff to reject a Disulfide. 2536 """ 2537 2538 def name_to_id(fname: str) -> str: 2539 """ 2540 Returns the PDB ID from the filename. 2541 2542 :param fname: Complete PDB filename 2543 :return: PDB ID 2544 """ 2545 ent = fname[3:-4] 2546 return ent 2547 2548 import os 2549 2550 entrylist = [] 2551 problem_ids = [] 2552 bad = bad_dist = 0 2553 2554 # we use the specialized list class DisulfideList to contain our disulfides 2555 # we'll use a dict to store DisulfideList objects, indexed by the structure ID 2556 All_ss_dict = {} 2557 All_ss_list = DisulfideList([], "PDB_SS") 2558 All_ss_dict2 = {} # new dict of pointers to indices 2559 2560 start = time.time() 2561 cwd = os.getcwd() 2562 2563 # Build a list of PDB files in PDB_DIR that are readable. These files were downloaded 2564 # via the RCSB web query interface for structures containing >= 1 SS Bond. 2565 2566 os.chdir(pdbdir) 2567 2568 ss_filelist = glob.glob("*.ent") 2569 tot = len(ss_filelist) 2570 2571 if verbose: 2572 print(f"PDB Directory {pdbdir} contains: {tot} files") 2573 2574 # the filenames are in the form pdb{entry}.ent, I loop through them and extract 2575 # the PDB ID, with Disulfide.name_to_id(), then add to entrylist. 2576 2577 for entry in ss_filelist: 2578 entrylist.append(name_to_id(entry)) 2579 2580 # create a dataframe with the following columns for the disulfide conformations 2581 # extracted from the structure 2582 2583 SS_df = pandas.DataFrame(columns=Torsion_DF_Cols) 2584 2585 # define a tqdm progressbar using the fully loaded entrylist list. 2586 # If numb is passed then 2587 # only do the last numb entries. 2588 2589 if numb > 0: 2590 pbar = tqdm(entrylist[:numb], ncols=PBAR_COLS) 2591 else: 2592 pbar = tqdm(entrylist, ncols=PBAR_COLS) 2593 2594 tot = 0 2595 cnt = 0 2596 # loop over ss_filelist, create disulfides and initialize them 2597 for entry in pbar: 2598 pbar.set_postfix( 2599 {"ID": entry, "Bad": bad, "Ca": bad_dist, "Cnt": tot} 2600 ) # update the progress bar 2601 2602 # returns an empty list if none are found. 2603 _sslist = DisulfideList([], entry) 2604 _sslist = proteusPy.DisulfideList.load_disulfides_from_id( 2605 entry, model_numb=0, verbose=verbose, quiet=quiet, pdb_dir=pdbdir 2606 ) 2607 sslist, xchain = prune_extra_ss(_sslist) 2608 2609 if len(sslist) > 0: 2610 sslist2 = [] # list to hold indices for ss_dict2 2611 for ss in sslist: 2612 # Ca distance cutoff 2613 dist = ss.ca_distance 2614 if dist >= dist_cutoff and dist_cutoff != -1.0: 2615 bad_dist += 1 2616 continue 2617 2618 All_ss_list.append(ss) 2619 new_row = [ 2620 ss.pdb_id, 2621 ss.name, 2622 ss.proximal, 2623 ss.distal, 2624 ss.chi1, 2625 ss.chi2, 2626 ss.chi3, 2627 ss.chi4, 2628 ss.chi5, 2629 ss.energy, 2630 ss.ca_distance, 2631 ss.cb_distance, 2632 ss.phiprox, 2633 ss.psiprox, 2634 ss.phidist, 2635 ss.psidist, 2636 ss.torsion_length, 2637 ss.rho, 2638 ] 2639 2640 # add the row to the end of the dataframe 2641 SS_df.loc[len(SS_df.index)] = new_row.copy() # deep copy 2642 sslist2.append(cnt) 2643 cnt += 1 2644 tot += 1 2645 2646 # All_ss_dict[entry] = sslist 2647 # print(f'Entry: {entry}. Dict indices: {sslist2}') 2648 All_ss_dict2[entry] = sslist2 2649 # print(f'{entry} ss dict adding: {sslist2}') 2650 2651 else: 2652 # at this point I really shouldn't have any bad non-parsible file 2653 bad += 1 2654 problem_ids.append(entry) 2655 # os.remove(f'pdb{entry}.ent') 2656 2657 if bad > 0: 2658 prob_cols = ["id"] 2659 problem_df = pandas.DataFrame(columns=prob_cols) 2660 problem_df["id"] = problem_ids 2661 2662 print( 2663 f"-> Extract_Disulfides(): Found and removed: {len(problem_ids)} non-parsable structures." 2664 ) 2665 print( 2666 f"-> Extract_Disulfides(): Saving problem IDs to file: {datadir}{problemfile}" 2667 ) 2668 2669 problem_df.to_csv(f"{datadir}{problemfile}") 2670 else: 2671 if verbose: 2672 print("-> Extract_Disulfides(): No non-parsable structures found.") 2673 2674 if bad_dist > 0: 2675 print(f"-> Extract_Disulfides(): Found and ignored: {bad_dist} long SS bonds.") 2676 else: 2677 if verbose: 2678 print("No problems found.") 2679 2680 # dump the all_ss list of disulfides to a .pkl file. ~520 MB. 2681 fname = f"{datadir}{picklefile}" 2682 print( 2683 f"-> Extract_Disulfides(): Saving {len(All_ss_list)} Disulfides to file: {fname}" 2684 ) 2685 2686 with open(fname, "wb+") as f: 2687 pickle.dump(All_ss_list, f) 2688 2689 # dump the dict2 disulfides to a .pkl file. ~520 MB. 2690 dict_len = len(All_ss_dict2) 2691 fname = f"{datadir}{dictfile}" 2692 print( 2693 f"-> Extract_Disulfides(): Saving indices of {dict_len} Disulfide-containing PDB IDs to file: {fname}" 2694 ) 2695 2696 with open(fname, "wb+") as f: 2697 pickle.dump(All_ss_dict2, f) 2698 2699 # save the torsions 2700 fname = f"{datadir}{torsionfile}" 2701 print(f"-> Extract_Disulfides(): Saving torsions to file: {fname}") 2702 SS_df.to_csv(fname) 2703 2704 end = time.time() 2705 elapsed = end - start 2706 2707 print( 2708 f"-> Extract_Disulfides(): Disulfide Extraction complete! Elapsed time:\ 2709 {datetime.timedelta(seconds=elapsed)} (h:m:s)" 2710 ) 2711 2712 # return to original directory 2713 os.chdir(cwd) 2714 return 2715 2716 2717def check_header_from_file( 2718 filename: str, model_numb=0, verbose=False, dbg=False 2719) -> bool: 2720 """ 2721 Check the Disulfides in the PDB file. 2722 2723 NB: Requires EGS-Modified BIO.parse_pdb_header.py from https://github.com/suchanek/biopython/ 2724 2725 :param filename: the name of the PDB entry. 2726 :param model_numb: model number to use, defaults to 0 for single structure files. 2727 :param verbose: print info while parsing 2728 :return: True if parsable 2729 2730 Example: 2731 Assuming ```MODEL_DIR``` has the pdb5rsa.ent file we can load the disulfides 2732 with the following: 2733 2734 >>> from proteusPy.Disulfide import Disulfide, check_header_from_file 2735 >>> MODEL_DIR = '../data/' 2736 >>> OK = False 2737 >>> OK = check_header_from_file(f'{MODEL_DIR}pdb5rsa.ent', verbose=True) 2738 -> check_header_from_file() - Parsing file: ../data/pdb5rsa.ent: 2739 -> SSBond: 1: tmp: 26A - 84A 2740 -> SSBond: 2: tmp: 40A - 95A 2741 -> SSBond: 3: tmp: 58A - 110A 2742 -> SSBond: 4: tmp: 65A - 72A 2743 >>> OK 2744 True 2745 """ 2746 import os 2747 2748 i = 1 2749 proximal = distal = -1 2750 _chaina = None 2751 _chainb = None 2752 2753 parser = PDBParser(PERMISSIVE=True) 2754 2755 # Biopython uses the Structure -> Model -> Chain hierarchy to organize 2756 # structures. All are iterable. 2757 2758 structure = parser.get_structure("tmp", file=filename) 2759 struct_name = structure.get_id() 2760 2761 model = structure[model_numb] 2762 2763 if verbose: 2764 print(f"-> check_header_from_file() - Parsing file: {filename}:") 2765 2766 ssbond_dict = structure.header["ssbond"] # NB: this requires the modified code 2767 2768 # list of tuples with (proximal distal chaina chainb) 2769 ssbonds = parse_ssbond_header_rec(ssbond_dict) 2770 2771 for pair in ssbonds: 2772 # in the form (proximal, distal, chain) 2773 proximal = pair[0] 2774 distal = pair[1] 2775 2776 if not proximal.isnumeric() or not distal.isnumeric(): 2777 if verbose: 2778 mess = f" ! Cannot parse SSBond record (non-numeric IDs):\ 2779 {struct_name} Prox: {proximal} {chain1_id} Dist: {distal} {chain2_id}" 2780 warnings.warn(mess, DisulfideParseWarning) 2781 continue # was pass 2782 else: 2783 proximal = int(proximal) 2784 distal = int(distal) 2785 2786 chain1_id = pair[2] 2787 chain2_id = pair[3] 2788 2789 _chaina = model[chain1_id] 2790 _chainb = model[chain2_id] 2791 2792 if chain1_id != chain2_id: 2793 if verbose: 2794 mess = f" -> Cross Chain SS for: Prox: {proximal}{chain1_id} Dist: {distal}{chain2_id}" 2795 warnings.warn(mess, DisulfideParseWarning) 2796 pass # was break 2797 2798 try: 2799 prox_res = _chaina[proximal] 2800 dist_res = _chainb[distal] 2801 except KeyError: 2802 print( 2803 f" ! Cannot parse SSBond record (KeyError): {struct_name} Prox: <{proximal}> {chain1_id} Dist: <{distal}> {chain2_id}" 2804 ) 2805 return False 2806 2807 if (_chaina is not None) and (_chainb is not None): 2808 if _chaina[proximal].is_disordered() or _chainb[distal].is_disordered(): 2809 continue 2810 else: 2811 if verbose: 2812 print( 2813 f" -> SSBond: {i}: {struct_name}: {proximal}{chain1_id} - {distal}{chain2_id}" 2814 ) 2815 else: 2816 if dbg: 2817 print( 2818 f" -> NULL chain(s): {struct_name}: {proximal}{chain1_id} - {distal}{chain2_id}" 2819 ) 2820 i += 1 2821 return True 2822 2823 2824def check_header_from_id( 2825 struct_name: str, pdb_dir=".", model_numb=0, verbose=False, dbg=False 2826) -> bool: 2827 """ 2828 Check parsability PDB ID and initializes the Disulfide objects. 2829 Assumes the file is downloaded in ```MODEL_DIR``` path. 2830 2831 NB: Requires EGS-Modified BIO.parse_pdb_header.py from https://github.com/suchanek/biopython/ 2832 2833 :param struct_name: the name of the PDB entry. 2834 :param pdb_dir: path to the PDB files, defaults to PDB_DIR 2835 :param model_numb: model number to use, defaults to 0 for single structure files. 2836 :param verbose: print info while parsing 2837 :param dbg: Debugging Flag 2838 :return: True if OK, False otherwise 2839 2840 Example: 2841 Assuming the PDB_DIR has the pdb5rsa.ent file we can check the file thusly: 2842 2843 >>> from proteusPy.Disulfide import Disulfide, check_header_from_id 2844 >>> MODEL_DIR = '/Users/egs/PDB/good/' 2845 >>> OK = False 2846 >>> OK = check_header_from_id('5rsa', pdb_dir=MODEL_DIR, verbose=True) 2847 -> SSBond: 1: 5rsa: 26A - 84A 2848 -> SSBond: 2: 5rsa: 40A - 95A 2849 -> SSBond: 3: 5rsa: 58A - 110A 2850 -> SSBond: 4: 5rsa: 65A - 72A 2851 >>> OK 2852 True 2853 """ 2854 parser = PDBParser(PERMISSIVE=True, QUIET=True) 2855 structure = parser.get_structure(struct_name, file=f"{pdb_dir}pdb{struct_name}.ent") 2856 model = structure[0] 2857 2858 ssbond_dict = structure.header["ssbond"] # NB: this requires the modified code 2859 2860 bondlist = [] 2861 i = 0 2862 2863 # get a list of tuples containing the proximal, distal residue IDs for 2864 # all SSBonds in the chain. 2865 bondlist = parse_ssbond_header_rec(ssbond_dict) 2866 2867 if len(bondlist) == 0: 2868 if verbose: 2869 print("-> check_header_from_id(): no bonds found in bondlist.") 2870 return False 2871 2872 for pair in bondlist: 2873 # in the form (proximal, distal, chain) 2874 proximal = pair[0] 2875 distal = pair[1] 2876 chain1 = pair[2] 2877 chain2 = pair[3] 2878 2879 chaina = model[chain1] 2880 chainb = model[chain2] 2881 2882 try: 2883 prox_residue = chaina[proximal] 2884 dist_residue = chainb[distal] 2885 2886 prox_residue.disordered_select("CYS") 2887 dist_residue.disordered_select("CYS") 2888 2889 if ( 2890 prox_residue.get_resname() != "CYS" 2891 or dist_residue.get_resname() != "CYS" 2892 ): 2893 if verbose: 2894 print( 2895 f"build_disulfide() requires CYS at both residues:\ 2896 {prox_residue.get_resname()} {dist_residue.get_resname()}" 2897 ) 2898 return False 2899 except KeyError: 2900 if dbg: 2901 print( 2902 f"Keyerror: {struct_name}: {proximal} {chain1} - {distal} {chain2}" 2903 ) 2904 return False 2905 2906 if verbose: 2907 print( 2908 f" -> SSBond: {i+1}: {struct_name}: {proximal}{chain1} - {distal}{chain2}" 2909 ) 2910 2911 i += 1 2912 return True 2913 2914 2915def Disulfide_Energy_Function(x: list) -> float: 2916 """ 2917 Compute the approximate torsional energy (kcal/mpl) for the input dihedral angles. 2918 2919 :param x: A list of dihedral angles: [chi1, chi2, chi3, chi4, chi5] 2920 :return: Energy in kcal/mol 2921 2922 Example: 2923 >>> from proteusPy.Disulfide import Disulfide_Energy_Function 2924 >>> dihed = [-60.0, -60.0, -90.0, -60.0, -90.0] 2925 >>> res = Disulfide_Energy_Function(dihed) 2926 >>> res 2927 2.5999999999999996 2928 """ 2929 import numpy as np 2930 2931 chi1, chi2, chi3, chi4, chi5 = x 2932 energy = 2.0 * (np.cos(np.deg2rad(3.0 * chi1)) + np.cos(np.deg2rad(3.0 * chi5))) 2933 energy += np.cos(np.deg2rad(3.0 * chi2)) + np.cos(np.deg2rad(3.0 * chi4)) 2934 energy += ( 2935 3.5 * np.cos(np.deg2rad(2.0 * chi3)) 2936 + 0.6 * np.cos(np.deg2rad(3.0 * chi3)) 2937 + 10.1 2938 ) 2939 return energy 2940 2941 2942def Minimize(inputSS: Disulfide) -> Disulfide: 2943 """ 2944 Minimizes the energy of a Disulfide object using the Nelder-Mead optimization method. 2945 2946 Parameters: 2947 inputSS (Disulfide): The Disulfide object to be minimized. 2948 2949 Returns: 2950 Disulfide: The minimized Disulfide object. 2951 2952 """ 2953 from scipy.optimize import minimize 2954 2955 from proteusPy import Disulfide, Disulfide_Energy_Function 2956 2957 initial_guess = inputSS.torsion_array 2958 result = minimize(Disulfide_Energy_Function, initial_guess, method="Nelder-Mead") 2959 minimum_conformation = result.x 2960 modelled_min = Disulfide("minimized") 2961 modelled_min.dihedrals = minimum_conformation 2962 modelled_min.build_yourself() 2963 return modelled_min 2964 2965 2966if __name__ == "__main__": 2967 import doctest 2968 2969 doctest.testmod() 2970 2971# End of file
Torsion_DF_Cols =
['source', 'ss_id', 'proximal', 'distal', 'chi1', 'chi2', 'chi3', 'chi4', 'chi5', 'energy', 'ca_distance', 'cb_distance', 'phi_prox', 'psi_prox', 'phi_dist', 'psi_dist', 'torsion_length', 'rho']
class
Disulfide:
91class Disulfide: 92 """ 93 This class provides a Python object and methods representing a physical disulfide bond 94 either extracted from the RCSB protein databank or built using the 95 [proteusPy.Turtle3D](turtle3D.html) class. The disulfide bond is an important 96 intramolecular stabilizing structural element and is characterized by: 97 98 * Atomic coordinates for the atoms N, Cα, Cβ, C', Sγ for both residues. 99 These are stored as both raw atomic coordinates as read from the RCSB file 100 and internal local coordinates. 101 * The dihedral angles Χ1 - Χ5 for the disulfide bond 102 * A name, by default {pdb_id}{prox_resnumb}{prox_chain}_{distal_resnum}{distal_chain} 103 * Proximal residue number 104 * Distal residue number 105 * Approximate bond torsional energy (kcal/mol): 106 107 $$ 108 E_{kcal/mol} \\approx 2.0 * cos(3.0 * \\chi_{1}) + cos(3.0 * \\chi_{5}) + cos(3.0 * \\chi_{2}) + 109 $$ 110 $$ 111 cos(3.0 * \chi_{4}) + 3.5 * cos(2.0 * \chi_{3}) + 0.6 * cos(3.0 * \chi_{3}) + 10.1 112 $$ 113 114 The equation embodies the typical 3-fold rotation barriers associated with single bonds, 115 (Χ1, Χ5, Χ2, Χ4) and a high 2-fold barrier for Χ3, resulting from the partial double bond 116 character of the S-S bond. This property leads to two major disulfide families, characterized 117 by the sign of Χ3. *Left-handed* disulfides have Χ3 < 0° and *right-handed* disulfides have 118 Χ3 > 0°. Within this breakdown there are numerous subfamilies, broadly known as the *hook*, 119 *spiral* and *staple*. These are under characgterization. 120 121 * Euclidean length of the dihedral angles (degrees) defined as: 122 $$\sqrt(\chi_{1}^{2} + \chi_{2}^{2} + \chi_{3}^{2} + \chi_{4}^{2} + \chi_{5}^{2})$$ 123 * Cα - Cα distance (Å) 124 * Cβ - Cβ distance (Å) 125 * The previous C' and next N for both the proximal and distal residues. These are needed 126 to calculate the backbone dihedral angles Φ and Ψ. 127 * Backbone dihedral angles Φ and Ψ, when possible. Not all structures are complete and 128 in those cases the atoms needed may be undefined. In this case the Φ and Ψ angles are set 129 to -180°. 130 131 The class also provides a rendering capabilities using the excellent [PyVista](https://pyvista.org) 132 library, and can display disulfides interactively in a variety of display styles: 133 * 'sb' - Split Bonds style - bonds colored by their atom type 134 * 'bs' - Ball and Stick style - split bond coloring with small atoms 135 * 'pd' - Proximal/Distal style - bonds colored *Red* for proximal residue and *Green* for 136 the distal residue. 137 * 'cpk' - CPK style rendering, colored by atom type: 138 * Carbon - Grey 139 * Nitrogen - Blue 140 * Sulfur - Yellow 141 * Oxygen - Red 142 * Hydrogen - White 143 144 Individual renderings can be saved to a file, and animations created. 145 """ 146 147 def __init__( 148 self, 149 name: str = "SSBOND", 150 proximal: int = -1, 151 distal: int = -1, 152 proximal_chain: str = "A", 153 distal_chain: str = "A", 154 pdb_id: str = "1egs", 155 quiet: bool = True, 156 permissive: bool = True, 157 ) -> None: 158 """ 159 __init__ Initialize the class to defined internal values. 160 161 :param name: Disulfide name, by default "SSBOND" 162 163 """ 164 self.name = name 165 self.proximal = proximal 166 self.distal = distal 167 self.energy = _FLOAT_INIT 168 self.proximal_chain = proximal_chain 169 self.distal_chain = distal_chain 170 self.pdb_id = pdb_id 171 self.proximal_residue_fullid = str("") 172 self.distal_residue_fullid = str("") 173 self.PERMISSIVE = permissive 174 self.QUIET = quiet 175 self.ca_distance = _FLOAT_INIT 176 self.cb_distance = _FLOAT_INIT 177 self.torsion_array = np.array( 178 (_ANG_INIT, _ANG_INIT, _ANG_INIT, _ANG_INIT, _ANG_INIT) 179 ) 180 self.phiprox = _ANG_INIT 181 self.psiprox = _ANG_INIT 182 self.phidist = _ANG_INIT 183 self.psidist = _ANG_INIT 184 185 # global coordinates for the Disulfide, typically as 186 # returned from the PDB file 187 188 self.n_prox = Vector(0, 0, 0) 189 self.ca_prox = Vector(0, 0, 0) 190 self.c_prox = Vector(0, 0, 0) 191 self.o_prox = Vector(0, 0, 0) 192 self.cb_prox = Vector(0, 0, 0) 193 self.sg_prox = Vector(0, 0, 0) 194 self.sg_dist = Vector(0, 0, 0) 195 self.cb_dist = Vector(0, 0, 0) 196 self.ca_dist = Vector(0, 0, 0) 197 self.n_dist = Vector(0, 0, 0) 198 self.c_dist = Vector(0, 0, 0) 199 self.o_dist = Vector(0, 0, 0) 200 201 # set when we can't find previous or next prox or distal 202 # C' or N atoms. 203 self.missing_atoms = False 204 self.modelled = False 205 self.resolution = -1.0 206 207 # need these to calculate backbone dihedral angles 208 self.c_prev_prox = Vector(0, 0, 0) 209 self.n_next_prox = Vector(0, 0, 0) 210 self.c_prev_dist = Vector(0, 0, 0) 211 self.n_next_dist = Vector(0, 0, 0) 212 213 # local coordinates for the Disulfide, computed using the Turtle3D in 214 # Orientation #1. these are generally private. 215 216 self._n_prox = Vector(0, 0, 0) 217 self._ca_prox = Vector(0, 0, 0) 218 self._c_prox = Vector(0, 0, 0) 219 self._o_prox = Vector(0, 0, 0) 220 self._cb_prox = Vector(0, 0, 0) 221 self._sg_prox = Vector(0, 0, 0) 222 self._sg_dist = Vector(0, 0, 0) 223 self._cb_dist = Vector(0, 0, 0) 224 self._ca_dist = Vector(0, 0, 0) 225 self._n_dist = Vector(0, 0, 0) 226 self._c_dist = Vector(0, 0, 0) 227 self._o_dist = Vector(0, 0, 0) 228 229 # need these to calculate backbone dihedral angles 230 self._c_prev_prox = Vector(0, 0, 0) 231 self._n_next_prox = Vector(0, 0, 0) 232 self._c_prev_dist = Vector(0, 0, 0) 233 self._n_next_dist = Vector(0, 0, 0) 234 235 # Dihedral angles for the disulfide bond itself, set to _ANG_INIT 236 self.chi1 = _ANG_INIT 237 self.chi2 = _ANG_INIT 238 self.chi3 = _ANG_INIT 239 self.chi4 = _ANG_INIT 240 self.chi5 = _ANG_INIT 241 self.rho = _ANG_INIT # new dihedral angle: Nprox - Ca_prox - Ca_dist - N_dist 242 243 self.torsion_length = _FLOAT_INIT 244 245 # comparison operators, used for sorting. keyed to SS bond energy 246 def __lt__(self, other): 247 if isinstance(other, Disulfide): 248 return self.energy < other.energy 249 250 def __le__(self, other): 251 if isinstance(other, Disulfide): 252 return self.energy <= other.energy 253 254 def __gt__(self, other): 255 if isinstance(other, Disulfide): 256 return self.energy > other.energy 257 258 def __ge__(self, other): 259 if isinstance(other, Disulfide): 260 return self.energy >= other.energy 261 262 def __eq__(self, other): 263 if isinstance(other, Disulfide): 264 return self.energy == other.energy 265 266 def __ne__(self, other): 267 if isinstance(other, Disulfide): 268 return self.energy != other.energy 269 270 def __repr__(self): 271 """ 272 Representation for the Disulfide class 273 """ 274 s1 = self.repr_ss_info() 275 res = f"{s1}>" 276 return res 277 278 def _render( 279 self, 280 pvplot: pv.Plotter, 281 style="bs", 282 plain=False, 283 bondcolor=BOND_COLOR, 284 bs_scale=BS_SCALE, 285 spec=SPECULARITY, 286 specpow=SPEC_POWER, 287 translate=True, 288 bond_radius=BOND_RADIUS, 289 res=100, 290 ): 291 """ 292 Update the passed pyVista plotter() object with the mesh data for the 293 input Disulfide Bond. Used internally 294 295 Parameters 296 ---------- 297 pvplot : pv.Plotter 298 pyvista.Plotter object 299 300 style : str, optional 301 Rendering style, by default 'bs'. One of 'bs', 'st', 'cpk', Render as \ 302 CPK, ball-and-stick or stick. Bonds are colored by atom color, unless \ 303 'plain' is specified. 304 305 plain : bool, optional 306 Used internally, by default False 307 308 bondcolor : pyVista color name, optional bond color for simple bonds, by default BOND_COLOR 309 310 bs_scale : float, optional 311 scale factor (0-1) to reduce the atom sizes for ball and stick, by default BS_SCALE 312 313 spec : float, optional 314 specularity (0-1), where 1 is totally smooth and 0 is rough, by default SPECULARITY 315 316 specpow : int, optional 317 exponent used for specularity calculations, by default SPEC_POWER 318 319 translate : bool, optional 320 Flag used internally to indicate if we should translate \ 321 the disulfide to its geometric center of mass, by default True. 322 323 Returns 324 ------- 325 pv.Plotter 326 Updated pv.Plotter object with atoms and bonds. 327 """ 328 329 _bradius = bond_radius 330 coords = self.internal_coords() 331 missing_atoms = self.missing_atoms 332 clen = coords.shape[0] 333 334 model = self.modelled 335 if model: 336 all_atoms = False 337 else: 338 all_atoms = True 339 340 if translate: 341 cofmass = self.cofmass() 342 for i in range(clen): 343 coords[i] = coords[i] - cofmass 344 345 atoms = ( 346 "N", 347 "C", 348 "C", 349 "O", 350 "C", 351 "SG", 352 "N", 353 "C", 354 "C", 355 "O", 356 "C", 357 "SG", 358 "C", 359 "N", 360 "C", 361 "N", 362 ) 363 pvp = pvplot 364 365 # bond connection table with atoms in the specific order shown above: 366 # returned by ss.get_internal_coords() 367 368 def draw_bonds( 369 pvp, 370 bradius=BOND_RADIUS, 371 style="sb", 372 bcolor=BOND_COLOR, 373 missing=True, 374 all_atoms=True, 375 res=100, 376 ): 377 """ 378 Generate the appropriate pyVista cylinder objects to represent 379 a particular disulfide bond. This utilizes a connection table 380 for the starting and ending atoms and a color table for the 381 bond colors. Used internally. 382 383 :param pvp: input plotter object to be updated 384 :param bradius: bond radius 385 :param style: bond style. One of sb, plain, pd 386 :param bcolor: pyvista color 387 :param missing: True if atoms are missing, False othersie 388 :param all_atoms: True if rendering O, False if only backbone rendered 389 390 :return pvp: Updated Plotter object. 391 392 """ 393 _bond_conn = np.array( 394 [ 395 [0, 1], # n-ca 396 [1, 2], # ca-c 397 [2, 3], # c-o 398 [1, 4], # ca-cb 399 [4, 5], # cb-sg 400 [6, 7], # n-ca 401 [7, 8], # ca-c 402 [8, 9], # c-o 403 [7, 10], # ca-cb 404 [10, 11], # cb-sg 405 [5, 11], # sg -sg 406 [12, 0], # cprev_prox-n 407 [2, 13], # c-nnext_prox 408 [14, 6], # cprev_dist-n_dist 409 [8, 15], # c-nnext_dist 410 ] 411 ) 412 413 # modeled disulfides only have backbone atoms since 414 # phi and psi are undefined, which makes the carbonyl 415 # oxygen (O) undefined as well. Their previous and next N 416 # are also undefined. 417 418 _bond_conn_backbone = np.array( 419 [ 420 [0, 1], # n-ca 421 [1, 2], # ca-c 422 [1, 4], # ca-cb 423 [4, 5], # cb-sg 424 [6, 7], # n-ca 425 [7, 8], # ca-c 426 [7, 10], # ca-cb 427 [10, 11], # cb-sg 428 [5, 11], # sg -sg 429 ] 430 ) 431 432 # colors for the bonds. Index into ATOM_COLORS array 433 _bond_split_colors = np.array( 434 [ 435 ("N", "C"), 436 ("C", "C"), 437 ("C", "O"), 438 ("C", "C"), 439 ("C", "SG"), 440 ("N", "C"), 441 ("C", "C"), 442 ("C", "O"), 443 ("C", "C"), 444 ("C", "SG"), 445 ("SG", "SG"), 446 # prev and next C-N bonds - color by atom Z 447 ("C", "N"), 448 ("C", "N"), 449 ("C", "N"), 450 ("C", "N"), 451 ] 452 ) 453 454 _bond_split_colors_backbone = np.array( 455 [ 456 ("N", "C"), 457 ("C", "C"), 458 ("C", "C"), 459 ("C", "SG"), 460 ("N", "C"), 461 ("C", "C"), 462 ("C", "C"), 463 ("C", "SG"), 464 ("SG", "SG"), 465 ] 466 ) 467 # work through connectivity and colors 468 orig_col = dest_col = bcolor 469 470 if all_atoms: 471 bond_conn = _bond_conn 472 bond_split_colors = _bond_split_colors 473 else: 474 bond_conn = _bond_conn_backbone 475 bond_split_colors = _bond_split_colors_backbone 476 477 for i in range(len(bond_conn)): 478 if all_atoms: 479 if i > 10 and missing_atoms == True: # skip missing atoms 480 continue 481 482 bond = bond_conn[i] 483 484 # get the indices for the origin and destination atoms 485 orig = bond[0] 486 dest = bond[1] 487 488 col = bond_split_colors[i] 489 490 # get the coords 491 prox_pos = coords[orig] 492 distal_pos = coords[dest] 493 494 # compute a direction vector 495 direction = distal_pos - prox_pos 496 497 # compute vector length. divide by 2 since split bond 498 height = math.dist(prox_pos, distal_pos) / 2.0 499 500 # the cylinder origins are actually in the 501 # middle so we translate 502 503 origin = prox_pos + 0.5 * direction # for a single plain bond 504 origin1 = prox_pos + 0.25 * direction 505 origin2 = prox_pos + 0.75 * direction 506 507 bradius = _bradius 508 509 if style == "plain": 510 orig_col = dest_col = bcolor 511 512 # proximal-distal red/green coloring 513 elif style == "pd": 514 if i <= 4 or i == 11 or i == 12: 515 orig_col = dest_col = "red" 516 else: 517 orig_col = dest_col = "green" 518 if i == 10: 519 orig_col = dest_col = "yellow" 520 else: 521 orig_col = ATOM_COLORS[col[0]] 522 dest_col = ATOM_COLORS[col[1]] 523 524 if i >= 11: # prev and next residue atoms for phi/psi calcs 525 bradius = _bradius * 0.5 # make smaller to distinguish 526 527 cap1 = pv.Sphere(center=prox_pos, radius=bradius) 528 cap2 = pv.Sphere(center=distal_pos, radius=bradius) 529 530 if style == "plain": 531 cyl = pv.Cylinder( 532 origin, direction, radius=bradius, height=height * 2.0 533 ) 534 pvp.add_mesh(cyl, color=orig_col) 535 else: 536 cyl1 = pv.Cylinder( 537 origin1, 538 direction, 539 radius=bradius, 540 height=height, 541 capping=False, 542 resolution=res, 543 ) 544 cyl2 = pv.Cylinder( 545 origin2, 546 direction, 547 radius=bradius, 548 height=height, 549 capping=False, 550 resolution=res, 551 ) 552 pvp.add_mesh(cyl1, color=orig_col) 553 pvp.add_mesh(cyl2, color=dest_col) 554 555 pvp.add_mesh(cap1, color=orig_col) 556 pvp.add_mesh(cap2, color=dest_col) 557 558 return pvp # end draw_bonds 559 560 if style == "cpk": 561 i = 0 562 for atom in atoms: 563 rad = ATOM_RADII_CPK[atom] 564 pvp.add_mesh( 565 pv.Sphere(center=coords[i], radius=rad), 566 color=ATOM_COLORS[atom], 567 smooth_shading=True, 568 specular=spec, 569 specular_power=specpow, 570 ) 571 i += 1 572 573 elif style == "cov": 574 i = 0 575 for atom in atoms: 576 rad = ATOM_RADII_COVALENT[atom] 577 pvp.add_mesh( 578 pv.Sphere(center=coords[i], radius=rad), 579 color=ATOM_COLORS[atom], 580 smooth_shading=True, 581 specular=spec, 582 specular_power=specpow, 583 ) 584 i += 1 585 586 elif style == "bs": # ball and stick 587 i = 0 588 for atom in atoms: 589 rad = ATOM_RADII_CPK[atom] * bs_scale 590 if i > 11: 591 rad = rad * 0.75 592 593 pvp.add_mesh( 594 pv.Sphere(center=coords[i], radius=rad), 595 color=ATOM_COLORS[atom], 596 smooth_shading=True, 597 specular=spec, 598 specular_power=specpow, 599 ) 600 i += 1 601 pvp = draw_bonds( 602 pvp, style="bs", missing=missing_atoms, all_atoms=all_atoms 603 ) 604 605 elif style == "sb": # splitbonds 606 pvp = draw_bonds( 607 pvp, style="sb", missing=missing_atoms, all_atoms=all_atoms 608 ) 609 610 elif style == "pd": # proximal-distal 611 pvp = draw_bonds( 612 pvp, style="pd", missing=missing_atoms, all_atoms=all_atoms 613 ) 614 615 else: # plain 616 pvp = draw_bonds( 617 pvp, 618 style="plain", 619 bcolor=bondcolor, 620 missing=missing_atoms, 621 all_atoms=all_atoms, 622 ) 623 624 return pvp 625 626 def _plot( 627 self, 628 pvplot, 629 style="bs", 630 plain=False, 631 bondcolor=BOND_COLOR, 632 bs_scale=BS_SCALE, 633 spec=SPECULARITY, 634 specpow=SPEC_POWER, 635 translate=True, 636 bond_radius=BOND_RADIUS, 637 res=100, 638 ): 639 """ 640 Update the passed pyVista plotter() object with the mesh data for the 641 input Disulfide Bond. Used internally 642 643 Parameters 644 ---------- 645 pvplot : pv.Plotter 646 pyvista.Plotter object 647 648 style : str, optional 649 Rendering style, by default 'bs'. One of 'bs', 'st', 'cpk', Render as \ 650 CPK, ball-and-stick or stick. Bonds are colored by atom color, unless \ 651 'plain' is specified. 652 653 plain : bool, optional 654 Used internally, by default False 655 656 bondcolor : pyVista color name, optional bond color for simple bonds, by default BOND_COLOR 657 658 bs_scale : float, optional 659 scale factor (0-1) to reduce the atom sizes for ball and stick, by default BS_SCALE 660 661 spec : float, optional 662 specularity (0-1), where 1 is totally smooth and 0 is rough, by default SPECULARITY 663 664 specpow : int, optional 665 exponent used for specularity calculations, by default SPEC_POWER 666 667 translate : bool, optional 668 Flag used internally to indicate if we should translate \ 669 the disulfide to its geometric center of mass, by default True. 670 671 Returns 672 ------- 673 pv.Plotter 674 Updated pv.Plotter object with atoms and bonds. 675 """ 676 677 _bradius = bond_radius 678 coords = self.internal_coords() 679 missing_atoms = self.missing_atoms 680 clen = coords.shape[0] 681 682 model = self.modelled 683 if model: 684 all_atoms = False 685 else: 686 all_atoms = True 687 688 if translate: 689 cofmass = self.cofmass() 690 for i in range(clen): 691 coords[i] = coords[i] - cofmass 692 693 atoms = ( 694 "N", 695 "C", 696 "C", 697 "O", 698 "C", 699 "SG", 700 "N", 701 "C", 702 "C", 703 "O", 704 "C", 705 "SG", 706 "C", 707 "N", 708 "C", 709 "N", 710 ) 711 pvp = pvplot.copy() 712 713 # bond connection table with atoms in the specific order shown above: 714 # returned by ss.get_internal_coords() 715 716 def draw_bonds( 717 pvp, 718 bradius=BOND_RADIUS, 719 style="sb", 720 bcolor=BOND_COLOR, 721 missing=True, 722 all_atoms=True, 723 res=100, 724 ): 725 """ 726 Generate the appropriate pyVista cylinder objects to represent 727 a particular disulfide bond. This utilizes a connection table 728 for the starting and ending atoms and a color table for the 729 bond colors. Used internally. 730 731 :param pvp: input plotter object to be updated 732 :param bradius: bond radius 733 :param style: bond style. One of sb, plain, pd 734 :param bcolor: pyvista color 735 :param missing: True if atoms are missing, False othersie 736 :param all_atoms: True if rendering O, False if only backbone rendered 737 738 :return pvp: Updated Plotter object. 739 740 """ 741 _bond_conn = np.array( 742 [ 743 [0, 1], # n-ca 744 [1, 2], # ca-c 745 [2, 3], # c-o 746 [1, 4], # ca-cb 747 [4, 5], # cb-sg 748 [6, 7], # n-ca 749 [7, 8], # ca-c 750 [8, 9], # c-o 751 [7, 10], # ca-cb 752 [10, 11], # cb-sg 753 [5, 11], # sg -sg 754 [12, 0], # cprev_prox-n 755 [2, 13], # c-nnext_prox 756 [14, 6], # cprev_dist-n_dist 757 [8, 15], # c-nnext_dist 758 ] 759 ) 760 761 # modeled disulfides only have backbone atoms since 762 # phi and psi are undefined, which makes the carbonyl 763 # oxygen (O) undefined as well. Their previous and next N 764 # are also undefined. 765 766 _bond_conn_backbone = np.array( 767 [ 768 [0, 1], # n-ca 769 [1, 2], # ca-c 770 [1, 4], # ca-cb 771 [4, 5], # cb-sg 772 [6, 7], # n-ca 773 [7, 8], # ca-c 774 [7, 10], # ca-cb 775 [10, 11], # cb-sg 776 [5, 11], # sg -sg 777 ] 778 ) 779 780 # colors for the bonds. Index into ATOM_COLORS array 781 _bond_split_colors = np.array( 782 [ 783 ("N", "C"), 784 ("C", "C"), 785 ("C", "O"), 786 ("C", "C"), 787 ("C", "SG"), 788 ("N", "C"), 789 ("C", "C"), 790 ("C", "O"), 791 ("C", "C"), 792 ("C", "SG"), 793 ("SG", "SG"), 794 # prev and next C-N bonds - color by atom Z 795 ("C", "N"), 796 ("C", "N"), 797 ("C", "N"), 798 ("C", "N"), 799 ] 800 ) 801 802 _bond_split_colors_backbone = np.array( 803 [ 804 ("N", "C"), 805 ("C", "C"), 806 ("C", "C"), 807 ("C", "SG"), 808 ("N", "C"), 809 ("C", "C"), 810 ("C", "C"), 811 ("C", "SG"), 812 ("SG", "SG"), 813 ] 814 ) 815 # work through connectivity and colors 816 orig_col = dest_col = bcolor 817 818 if all_atoms: 819 bond_conn = _bond_conn 820 bond_split_colors = _bond_split_colors 821 else: 822 bond_conn = _bond_conn_backbone 823 bond_split_colors = _bond_split_colors_backbone 824 825 for i in range(len(bond_conn)): 826 if all_atoms: 827 if i > 10 and missing_atoms == True: # skip missing atoms 828 continue 829 830 bond = bond_conn[i] 831 832 # get the indices for the origin and destination atoms 833 orig = bond[0] 834 dest = bond[1] 835 836 col = bond_split_colors[i] 837 838 # get the coords 839 prox_pos = coords[orig] 840 distal_pos = coords[dest] 841 842 # compute a direction vector 843 direction = distal_pos - prox_pos 844 845 # compute vector length. divide by 2 since split bond 846 height = math.dist(prox_pos, distal_pos) / 2.0 847 848 # the cylinder origins are actually in the 849 # middle so we translate 850 851 origin = prox_pos + 0.5 * direction # for a single plain bond 852 origin1 = prox_pos + 0.25 * direction 853 origin2 = prox_pos + 0.75 * direction 854 855 bradius = _bradius 856 857 if style == "plain": 858 orig_col = dest_col = bcolor 859 860 # proximal-distal red/green coloring 861 elif style == "pd": 862 if i <= 4 or i == 11 or i == 12: 863 orig_col = dest_col = "red" 864 else: 865 orig_col = dest_col = "green" 866 if i == 10: 867 orig_col = dest_col = "yellow" 868 else: 869 orig_col = ATOM_COLORS[col[0]] 870 dest_col = ATOM_COLORS[col[1]] 871 872 if i >= 11: # prev and next residue atoms for phi/psi calcs 873 bradius = _bradius * 0.5 # make smaller to distinguish 874 875 cap1 = pv.Sphere(center=prox_pos, radius=bradius) 876 cap2 = pv.Sphere(center=distal_pos, radius=bradius) 877 878 if style == "plain": 879 cyl = pv.Cylinder( 880 origin, direction, radius=bradius, height=height * 2.0 881 ) 882 # pvp.add_mesh(cyl, color=orig_col) 883 pvp.append(cyl) 884 else: 885 cyl1 = pv.Cylinder( 886 origin1, 887 direction, 888 radius=bradius, 889 height=height, 890 capping=False, 891 resolution=res, 892 ) 893 cyl2 = pv.Cylinder( 894 origin2, 895 direction, 896 radius=bradius, 897 height=height, 898 capping=False, 899 resolution=res, 900 ) 901 # pvp.add_mesh(cyl1, color=orig_col) 902 # pvp.add_mesh(cyl2, color=dest_col) 903 pvp.append(cyl1) 904 pvp.append(cyl2) 905 906 # pvp.add_mesh(cap1, color=orig_col) 907 # pvp.add_mesh(cap2, color=dest_col) 908 pvp.append(cap1) 909 pvp.append(cap2) 910 911 return pvp.copy() # end draw_bonds 912 913 if style == "cpk": 914 i = 0 915 for atom in atoms: 916 rad = ATOM_RADII_CPK[atom] 917 pvp.append(pv.Sphere(center=coords[i], radius=rad)) 918 i += 1 919 920 elif style == "cov": 921 i = 0 922 for atom in atoms: 923 rad = ATOM_RADII_COVALENT[atom] 924 pvp.append(pv.Sphere(center=coords[i], radius=rad)) 925 i += 1 926 927 elif style == "bs": # ball and stick 928 i = 0 929 for atom in atoms: 930 rad = ATOM_RADII_CPK[atom] * bs_scale 931 if i > 11: 932 rad = rad * 0.75 933 934 pvp.append(pv.Sphere(center=coords[i])) 935 i += 1 936 pvp = draw_bonds( 937 pvp, style="bs", missing=missing_atoms, all_atoms=all_atoms 938 ) 939 940 elif style == "sb": # splitbonds 941 pvp = draw_bonds( 942 pvp, style="sb", missing=missing_atoms, all_atoms=all_atoms 943 ) 944 945 elif style == "pd": # proximal-distal 946 pvp = draw_bonds( 947 pvp, style="pd", missing=missing_atoms, all_atoms=all_atoms 948 ) 949 950 else: # plain 951 pvp = draw_bonds( 952 pvp, 953 style="plain", 954 bcolor=bondcolor, 955 missing=missing_atoms, 956 all_atoms=all_atoms, 957 ) 958 959 return 960 961 def _handle_SS_exception(self, message: str): 962 """ 963 This method catches an exception that occurs in the Disulfide 964 object (if PERMISSIVE), or raises it again, this time adding the 965 PDB line number to the error message. (private). 966 967 :param message: Error message 968 :raises DisulfideConstructionException: Fatal construction exception. 969 970 """ 971 # message = "%s at line %i." % (message) 972 message = f"{message}" 973 974 if self.PERMISSIVE: 975 # just print a warning - some residues/atoms may be missing 976 warnings.warn( 977 "DisulfideConstructionException: %s\n" 978 "Exception ignored.\n" 979 "Some atoms may be missing in the data structure." % message, 980 DisulfideConstructionWarning, 981 ) 982 else: 983 # exceptions are fatal - raise again with new message (including line nr) 984 raise DisulfideConstructionException(message) from None 985 986 @property 987 def dihedrals(self) -> list: 988 """ 989 Return a list containing the dihedral angles for the disulfide. 990 991 """ 992 return [self.chi1, self.chi2, self.chi3, self.chi4, self.chi5] 993 994 @dihedrals.setter 995 def dihedrals(self, dihedrals: list) -> None: 996 """ 997 Sets the disulfide dihedral angles to the inputs specified in the list. 998 999 :param dihedrals: list of dihedral angles. 1000 """ 1001 self.chi1 = dihedrals[0] 1002 self.chi2 = dihedrals[1] 1003 self.chi3 = dihedrals[2] 1004 self.chi4 = dihedrals[3] 1005 self.chi5 = dihedrals[4] 1006 self.compute_torsional_energy() 1007 self.compute_torsion_length() 1008 1009 def bounding_box(self) -> np.array: 1010 """ 1011 Return the bounding box array for the given disulfide 1012 1013 Returns 1014 ------- 1015 :return: np.Array(3,2): Array containing the min, max for X, Y, and Z respectively. 1016 Does not currently take the atom's radius into account. 1017 1018 """ 1019 res = np.zeros(shape=(3, 2)) 1020 xmin, xmax = self.compute_extents("x") 1021 ymin, ymax = self.compute_extents("y") 1022 zmin, zmax = self.compute_extents("z") 1023 1024 res[0] = [xmin, xmax] 1025 res[1] = [ymin, ymax] 1026 res[2] = [zmin, zmax] 1027 1028 return res 1029 1030 def build_yourself(self) -> None: 1031 """ 1032 Build a model Disulfide based its internal dihedral state 1033 Routine assumes turtle is in orientation #1 (at Ca, headed toward 1034 Cb, with N on left), builds disulfide, and updates the object's internal 1035 coordinates. It also adds the distal protein backbone, 1036 and computes the disulfide conformational energy. 1037 """ 1038 chi1 = self.chi1 1039 chi2 = self.chi2 1040 chi3 = self.chi3 1041 chi4 = self.chi4 1042 chi5 = self.chi5 1043 self.build_model(chi1, chi2, chi3, chi4, chi5) 1044 1045 def build_model( 1046 self, chi1: float, chi2: float, chi3: float, chi4: float, chi5: float 1047 ) -> None: 1048 """ 1049 Build a model Disulfide based on the input dihedral angles. 1050 Routine assumes turtle is in orientation #1 (at Ca, headed toward 1051 Cb, with N on left), builds disulfide, and updates the object's internal 1052 coordinates. It also adds the distal protein backbone, 1053 and computes the disulfide conformational energy. 1054 1055 :param chi1: Chi1 (degrees) 1056 :param chi2: Chi2 (degrees) 1057 :param chi3: Chi3 (degrees) 1058 :param chi4: Chi4 (degrees) 1059 :param chi5: Chi5 (degrees) 1060 1061 Example: 1062 >>> from proteusPy.Disulfide import Disulfide 1063 >>> modss = Disulfide('model') 1064 >>> modss.build_model(-60, -60, -90, -60, -60) 1065 >>> modss.display(style='sb') 1066 """ 1067 1068 self.set_dihedrals(chi1, chi2, chi3, chi4, chi5) 1069 self.proximal = 1 1070 self.distal = 2 1071 1072 tmp = Turtle3D("tmp") 1073 tmp.Orientation = 1 1074 1075 n = Vector(0, 0, 0) 1076 ca = Vector(0, 0, 0) 1077 cb = Vector(0, 0, 0) 1078 c = Vector(0, 0, 0) 1079 1080 self.ca_prox = tmp._position 1081 tmp.schain_to_bbone() 1082 n, ca, cb, c = build_residue(tmp) 1083 1084 self.n_prox = n 1085 self.ca_prox = ca 1086 self.c_prox = c 1087 self.cb_prox = cb 1088 1089 tmp.bbone_to_schain() 1090 tmp.move(1.53) 1091 tmp.roll(self.chi1) 1092 tmp.yaw(112.8) 1093 self.cb_prox = Vector(tmp._position) 1094 1095 tmp.move(1.86) 1096 tmp.roll(self.chi2) 1097 tmp.yaw(103.8) 1098 self.sg_prox = Vector(tmp._position) 1099 1100 tmp.move(2.044) 1101 tmp.roll(self.chi3) 1102 tmp.yaw(103.8) 1103 self.sg_dist = Vector(tmp._position) 1104 1105 tmp.move(1.86) 1106 tmp.roll(self.chi4) 1107 tmp.yaw(112.8) 1108 self.cb_dist = Vector(tmp._position) 1109 1110 tmp.move(1.53) 1111 tmp.roll(self.chi5) 1112 tmp.pitch(180.0) 1113 1114 tmp.schain_to_bbone() 1115 1116 n, ca, cb, c = build_residue(tmp) 1117 1118 self.n_dist = n 1119 self.ca_dist = ca 1120 self.c_dist = c 1121 self.compute_torsional_energy() 1122 self.compute_local_coords() 1123 self.ca_distance = distance3d(self.ca_prox, self.ca_dist) 1124 self.cb_distance = distance3d(self.cb_prox, self.cb_dist) 1125 self.torsion_array = np.array( 1126 (self.chi1, self.chi2, self.chi3, self.chi4, self.chi5) 1127 ) 1128 self.compute_torsion_length() 1129 self.compute_rho() 1130 self.missing_atoms = True 1131 self.modelled = True 1132 1133 def cofmass(self) -> np.array: 1134 """ 1135 Return the geometric center of mass for the internal coordinates of 1136 the given Disulfide. Missing atoms are not included. 1137 1138 :return: 3D array for the geometric center of mass 1139 """ 1140 1141 res = self.internal_coords() 1142 return res.mean(axis=0) 1143 1144 def copy(self): 1145 """ 1146 Copy the Disulfide. 1147 1148 :return: A copy of self. 1149 """ 1150 return copy.deepcopy(self) 1151 1152 def compute_extents(self, dim="z"): 1153 """ 1154 Calculate the internal coordinate extents for the input axis. 1155 1156 :param dim: Axis, one of 'x', 'y', 'z', by default 'z' 1157 :return: min, max 1158 """ 1159 1160 ic = self.internal_coords() 1161 # set default index to 'z' 1162 idx = 2 1163 1164 if dim == "x": 1165 idx = 0 1166 elif dim == "y": 1167 idx = 1 1168 elif dim == "z": 1169 idx = 2 1170 1171 _min = ic[:, idx].min() 1172 _max = ic[:, idx].max() 1173 return _min, _max 1174 1175 def compute_local_coords(self) -> None: 1176 """ 1177 Compute the internal coordinates for a properly initialized Disulfide Object. 1178 1179 :param self: SS initialized Disulfide object 1180 :returns: None, modifies internal state of the input 1181 """ 1182 1183 turt = Turtle3D("tmp") 1184 # get the coordinates as np.array for Turtle3D use. 1185 cpp = self.c_prev_prox.get_array() 1186 nnp = self.n_next_prox.get_array() 1187 1188 n = self.n_prox.get_array() 1189 ca = self.ca_prox.get_array() 1190 c = self.c_prox.get_array() 1191 cb = self.cb_prox.get_array() 1192 o = self.o_prox.get_array() 1193 sg = self.sg_prox.get_array() 1194 1195 sg2 = self.sg_dist.get_array() 1196 cb2 = self.cb_dist.get_array() 1197 ca2 = self.ca_dist.get_array() 1198 c2 = self.c_dist.get_array() 1199 n2 = self.n_dist.get_array() 1200 o2 = self.o_dist.get_array() 1201 1202 cpd = self.c_prev_dist.get_array() 1203 nnd = self.n_next_dist.get_array() 1204 1205 turt.orient_from_backbone(n, ca, c, cb, ORIENT_SIDECHAIN) 1206 1207 # internal (local) coordinates, stored as Vector objects 1208 # to_local returns np.array objects 1209 1210 self._n_prox = Vector(turt.to_local(n)) 1211 self._ca_prox = Vector(turt.to_local(ca)) 1212 self._c_prox = Vector(turt.to_local(c)) 1213 self._o_prox = Vector(turt.to_local(o)) 1214 self._cb_prox = Vector(turt.to_local(cb)) 1215 self._sg_prox = Vector(turt.to_local(sg)) 1216 1217 self._c_prev_prox = Vector(turt.to_local(cpp)) 1218 self._n_next_prox = Vector(turt.to_local(nnp)) 1219 self._c_prev_dist = Vector(turt.to_local(cpd)) 1220 self._n_next_dist = Vector(turt.to_local(nnd)) 1221 1222 self._n_dist = Vector(turt.to_local(n2)) 1223 self._ca_dist = Vector(turt.to_local(ca2)) 1224 self._c_dist = Vector(turt.to_local(c2)) 1225 self._o_dist = Vector(turt.to_local(o2)) 1226 self._cb_dist = Vector(turt.to_local(cb2)) 1227 self._sg_dist = Vector(turt.to_local(sg2)) 1228 1229 def compute_torsional_energy(self) -> float: 1230 """ 1231 Compute the approximate torsional energy for the Disulfide's 1232 conformation and sets its internal state. 1233 1234 :return: Energy (kcal/mol) 1235 """ 1236 # @TODO find citation for the ss bond energy calculation 1237 1238 def torad(deg): 1239 return np.radians(deg) 1240 1241 chi1 = self.chi1 1242 chi2 = self.chi2 1243 chi3 = self.chi3 1244 chi4 = self.chi4 1245 chi5 = self.chi5 1246 1247 energy = 2.0 * (cos(torad(3.0 * chi1)) + cos(torad(3.0 * chi5))) 1248 energy += cos(torad(3.0 * chi2)) + cos(torad(3.0 * chi4)) 1249 energy += 3.5 * cos(torad(2.0 * chi3)) + 0.6 * cos(torad(3.0 * chi3)) + 10.1 1250 1251 self.energy = energy 1252 return energy 1253 1254 def display(self, single=True, style="sb", light=True, shadows=False) -> None: 1255 """ 1256 Display the Disulfide bond in the specific rendering style. 1257 1258 :param single: Display the bond in a single panel in the specific style. 1259 :param style: Rendering style: One of: 1260 * 'sb' - split bonds 1261 * 'bs' - ball and stick 1262 * 'cpk' - CPK style 1263 * 'pd' - Proximal/Distal style - Red=proximal, Green=Distal 1264 * 'plain' - boring single color 1265 :param light: If True, light background, if False, dark 1266 1267 Example: 1268 >>> import proteusPy 1269 >>> from proteusPy.Disulfide import Disulfide 1270 >>> from proteusPy.DisulfideLoader import DisulfideLoader, Load_PDB_SS 1271 1272 >>> PDB_SS = Load_PDB_SS(verbose=False, subset=True) 1273 >>> ss = PDB_SS[0] 1274 >>> ss.display(style='cpk') 1275 >>> ss.screenshot(style='bs', fname='proteus_logo_sb.png') 1276 """ 1277 src = self.pdb_id 1278 enrg = self.energy 1279 1280 title = f"{src}: {self.proximal}{self.proximal_chain}-{self.distal}{self.distal_chain}: {enrg:.2f} kcal/mol. Cα: {self.ca_distance:.2f} Å Cβ: {self.cb_distance:.2f} Å Tors: {self.torsion_length:.2f}°" 1281 1282 if light: 1283 pv.set_plot_theme("document") 1284 else: 1285 pv.set_plot_theme("dark") 1286 1287 if single == True: 1288 _pl = pv.Plotter(window_size=WINSIZE) 1289 _pl.add_title(title=title, font_size=FONTSIZE) 1290 _pl.enable_anti_aliasing("msaa") 1291 # _pl.add_camera_orientation_widget() 1292 1293 self._render( 1294 _pl, 1295 style=style, 1296 bs_scale=BS_SCALE, 1297 spec=SPECULARITY, 1298 specpow=SPEC_POWER, 1299 ) 1300 _pl.reset_camera() 1301 if shadows == True: 1302 _pl.enable_shadows() 1303 _pl.show() 1304 1305 else: 1306 pl = pv.Plotter(window_size=WINSIZE, shape=(2, 2)) 1307 pl.subplot(0, 0) 1308 1309 pl.add_title(title=title, font_size=FONTSIZE) 1310 pl.enable_anti_aliasing("msaa") 1311 1312 # pl.add_camera_orientation_widget() 1313 1314 self._render( 1315 pl, 1316 style="cpk", 1317 bondcolor=BOND_COLOR, 1318 bs_scale=BS_SCALE, 1319 spec=SPECULARITY, 1320 specpow=SPEC_POWER, 1321 ) 1322 1323 pl.subplot(0, 1) 1324 pl.add_title(title=title, font_size=FONTSIZE) 1325 1326 self._render( 1327 pl, 1328 style="bs", 1329 bondcolor=BOND_COLOR, 1330 bs_scale=BS_SCALE, 1331 spec=SPECULARITY, 1332 specpow=SPEC_POWER, 1333 ) 1334 1335 pl.subplot(1, 0) 1336 pl.add_title(title=title, font_size=FONTSIZE) 1337 1338 self._render( 1339 pl, 1340 style="sb", 1341 bondcolor=BOND_COLOR, 1342 bs_scale=BS_SCALE, 1343 spec=SPECULARITY, 1344 specpow=SPEC_POWER, 1345 ) 1346 1347 pl.subplot(1, 1) 1348 pl.add_title(title=title, font_size=FONTSIZE) 1349 1350 self._render( 1351 pl, 1352 style="pd", 1353 bondcolor=BOND_COLOR, 1354 bs_scale=BS_SCALE, 1355 spec=SPECULARITY, 1356 specpow=SPEC_POWER, 1357 ) 1358 1359 pl.link_views() 1360 pl.reset_camera() 1361 if shadows == True: 1362 pl.enable_shadows() 1363 pl.show() 1364 return 1365 1366 def plot( 1367 self, pl, single=True, style="sb", light=True, shadows=False 1368 ) -> pv.Plotter: 1369 """ 1370 Return the pyVista Plotter object for the Disulfide bond in the specific rendering style. 1371 1372 :param single: Display the bond in a single panel in the specific style. 1373 :param style: Rendering style: One of: 1374 * 'sb' - split bonds 1375 * 'bs' - ball and stick 1376 * 'cpk' - CPK style 1377 * 'pd' - Proximal/Distal style - Red=proximal, Green=Distal 1378 * 'plain' - boring single color 1379 :param light: If True, light background, if False, dark 1380 """ 1381 src = self.pdb_id 1382 enrg = self.energy 1383 1384 title = f"{src}: {self.proximal}{self.proximal_chain}-{self.distal}{self.distal_chain}: {enrg:.2f} kcal/mol. Cα: {self.ca_distance:.2f} Å Cβ: {self.cb_distance:.2f} Å Tors: {self.torsion_length:.2f}°" 1385 1386 if light: 1387 pv.set_plot_theme("document") 1388 else: 1389 pv.set_plot_theme("dark") 1390 1391 if single == True: 1392 # _pl = pv.Plotter(window_size=WINSIZE) 1393 # _pl.add_title(title=title, font_size=FONTSIZE) 1394 pl.clear() 1395 pl.enable_anti_aliasing("msaa") 1396 # pl.add_camera_orientation_widget() 1397 1398 self._render( 1399 pl, 1400 style=style, 1401 bs_scale=BS_SCALE, 1402 spec=SPECULARITY, 1403 specpow=SPEC_POWER, 1404 ) 1405 pl.reset_camera() 1406 if shadows == True: 1407 pl.enable_shadows() 1408 else: 1409 pl = pv.Plotter(shape=(2, 2)) 1410 pl.subplot(0, 0) 1411 1412 # pl.add_title(title=title, font_size=FONTSIZE) 1413 pl.enable_anti_aliasing("msaa") 1414 1415 # pl.add_camera_orientation_widget() 1416 1417 self._render( 1418 pl, 1419 style="cpk", 1420 bondcolor=BOND_COLOR, 1421 bs_scale=BS_SCALE, 1422 spec=SPECULARITY, 1423 specpow=SPEC_POWER, 1424 ) 1425 1426 pl.subplot(0, 1) 1427 1428 self._render( 1429 pl, 1430 style="bs", 1431 bondcolor=BOND_COLOR, 1432 bs_scale=BS_SCALE, 1433 spec=SPECULARITY, 1434 specpow=SPEC_POWER, 1435 ) 1436 1437 pl.subplot(1, 0) 1438 1439 self._render( 1440 pl, 1441 style="sb", 1442 bondcolor=BOND_COLOR, 1443 bs_scale=BS_SCALE, 1444 spec=SPECULARITY, 1445 specpow=SPEC_POWER, 1446 ) 1447 1448 pl.subplot(1, 1) 1449 self._render( 1450 pl, 1451 style="pd", 1452 bondcolor=BOND_COLOR, 1453 bs_scale=BS_SCALE, 1454 spec=SPECULARITY, 1455 specpow=SPEC_POWER, 1456 ) 1457 1458 pl.link_views() 1459 pl.reset_camera() 1460 if shadows == True: 1461 pl.enable_shadows() 1462 return pl 1463 1464 def Distance_neighbors(self, others: DisulfideList, cutoff: float) -> DisulfideList: 1465 """ 1466 Return list of Disulfides whose RMS atomic distance is within 1467 the cutoff (Å) in the others list. 1468 1469 :param others: DisulfideList to search 1470 :param cutoff: Distance cutoff (Å) 1471 :return: DisulfideList within the cutoff 1472 """ 1473 1474 res = [ss.copy() for ss in others if self.Distance_RMS(ss) < cutoff] 1475 return DisulfideList(res, "neighbors") 1476 1477 def Distance_RMS(self, other) -> float: 1478 """ 1479 Calculate the RMS distance between the internal coordinates of self and another Disulfide. 1480 :param other: Comparison Disulfide 1481 :return: RMS distance (Å) 1482 """ 1483 1484 # Get internal coordinates of both objects 1485 ic1 = self.internal_coords() 1486 ic2 = other.internal_coords() 1487 1488 # Compute the sum of squared differences between corresponding internal coordinates 1489 totsq = sum(math.dist(p1, p2) ** 2 for p1, p2 in zip(ic1, ic2)) 1490 1491 # Compute the mean of the squared distances 1492 totsq /= len(ic1) 1493 1494 # Take the square root of the mean to get the RMS distance 1495 return math.sqrt(totsq) 1496 1497 def Torsion_RMS(self, other) -> float: 1498 """ 1499 Calculate the RMS distance between the dihedral angles of self and another Disulfide. 1500 :param other: Comparison Disulfide 1501 :return: RMS distance (deg). 1502 """ 1503 import math 1504 1505 # Get internal coordinates of both objects 1506 ic1 = self.torsion_array 1507 ic2 = other.torsion_array 1508 1509 # Compute the sum of squared differences between corresponding internal coordinates 1510 totsq = sum((p1 - p2) ** 2 for p1, p2 in zip(ic1, ic2)) 1511 # Compute the mean of the squared distances 1512 totsq /= len(ic1) 1513 1514 # Take the square root of the mean to get the RMS distance 1515 return math.sqrt(totsq) 1516 1517 def get_chains(self) -> tuple: 1518 """ 1519 Return the proximal and distal chain IDs for the Disulfide. 1520 1521 :return: tuple (proximal, distal) chain IDs 1522 """ 1523 prox = self.proximal_chain 1524 dist = self.distal_chain 1525 return tuple(prox, dist) 1526 1527 def get_permissive(self) -> bool: 1528 """ 1529 Return the Permissive flag state. (Used in PDB parsing) 1530 1531 :return: Permissive state 1532 """ 1533 return self.PERMISIVE 1534 1535 def get_full_id(self) -> tuple: 1536 """ 1537 Return the Disulfide full IDs (Used with BIO.PDB) 1538 1539 :return: Disulfide full IDs 1540 """ 1541 return (self.proximal_residue_fullid, self.distal_residue_fullid) 1542 1543 def initialize_disulfide_from_chain( 1544 self, chain1, chain2, proximal, distal, resolution, quiet=True 1545 ) -> None: 1546 """ 1547 Initialize a new Disulfide object with atomic coordinates from 1548 the proximal and distal coordinates, typically taken from a PDB file. 1549 This routine is primarily used internally when building the compressed 1550 database. 1551 1552 :param chain1: list of Residues in the model, eg: chain = model['A'] 1553 :param chain2: list of Residues in the model, eg: chain = model['A'] 1554 :param proximal: proximal residue sequence ID 1555 :param distal: distal residue sequence ID 1556 :param resolution: structure resolution 1557 :param quiet: Quiet or noisy parsing, defaults to True 1558 :raises DisulfideConstructionWarning: Raised when not parsed correctly 1559 """ 1560 id = chain1.get_full_id()[0] 1561 self.pdb_id = id 1562 1563 chi1 = chi2 = chi3 = chi4 = chi5 = _ANG_INIT 1564 1565 prox = int(proximal) 1566 dist = int(distal) 1567 1568 prox_residue = chain1[prox] 1569 dist_residue = chain2[dist] 1570 1571 if prox_residue.get_resname() != "CYS" or dist_residue.get_resname() != "CYS": 1572 print( 1573 f"build_disulfide() requires CYS at both residues: {prox} {prox_residue.get_resname()} {dist} {dist_residue.get_resname()} Chain: {prox_residue.get_segid()}" 1574 ) 1575 1576 # set the objects proximal and distal values 1577 self.set_resnum(proximal, distal) 1578 1579 if resolution is not None: 1580 self.resolution = resolution 1581 1582 self.proximal_chain = chain1.get_id() 1583 self.distal_chain = chain2.get_id() 1584 1585 self.proximal_residue_fullid = prox_residue.get_full_id() 1586 self.distal_residue_fullid = dist_residue.get_full_id() 1587 1588 if quiet: 1589 warnings.filterwarnings("ignore", category=DisulfideConstructionWarning) 1590 else: 1591 warnings.simplefilter("always") 1592 1593 # grab the coordinates for the proximal and distal residues as vectors 1594 # so we can do math on them later 1595 # proximal residue 1596 1597 try: 1598 n1 = prox_residue["N"].get_vector() 1599 ca1 = prox_residue["CA"].get_vector() 1600 c1 = prox_residue["C"].get_vector() 1601 o1 = prox_residue["O"].get_vector() 1602 cb1 = prox_residue["CB"].get_vector() 1603 sg1 = prox_residue["SG"].get_vector() 1604 1605 except Exception: 1606 raise DisulfideConstructionWarning( 1607 f"Invalid or missing coordinates for proximal residue {proximal}" 1608 ) from None 1609 1610 # distal residue 1611 try: 1612 n2 = dist_residue["N"].get_vector() 1613 ca2 = dist_residue["CA"].get_vector() 1614 c2 = dist_residue["C"].get_vector() 1615 o2 = dist_residue["O"].get_vector() 1616 cb2 = dist_residue["CB"].get_vector() 1617 sg2 = dist_residue["SG"].get_vector() 1618 1619 except Exception: 1620 raise DisulfideConstructionWarning( 1621 f"Invalid or missing coordinates for distal residue {distal}" 1622 ) from None 1623 1624 # previous residue and next residue - optional, used for phi, psi calculations 1625 try: 1626 prevprox = chain1[prox - 1] 1627 nextprox = chain1[prox + 1] 1628 1629 prevdist = chain2[dist - 1] 1630 nextdist = chain2[dist + 1] 1631 1632 cprev_prox = prevprox["C"].get_vector() 1633 nnext_prox = nextprox["N"].get_vector() 1634 1635 cprev_dist = prevdist["C"].get_vector() 1636 nnext_dist = nextdist["N"].get_vector() 1637 1638 # compute phi, psi for prox and distal 1639 self.phiprox = np.degrees(calc_dihedral(cprev_prox, n1, ca1, c1)) 1640 self.psiprox = np.degrees(calc_dihedral(n1, ca1, c1, nnext_prox)) 1641 self.phidist = np.degrees(calc_dihedral(cprev_dist, n2, ca2, c2)) 1642 self.psidist = np.degrees(calc_dihedral(n2, ca2, c2, nnext_dist)) 1643 1644 except Exception: 1645 mess = f"Missing coords for: {id} {prox-1} or {dist+1} for SS {proximal}-{distal}" 1646 cprev_prox = nnext_prox = cprev_dist = nnext_dist = Vector(-1.0, -1.0, -1.0) 1647 self.missing_atoms = True 1648 warnings.warn(mess, DisulfideConstructionWarning) 1649 1650 # update the positions and conformation 1651 self.set_positions( 1652 n1, 1653 ca1, 1654 c1, 1655 o1, 1656 cb1, 1657 sg1, 1658 n2, 1659 ca2, 1660 c2, 1661 o2, 1662 cb2, 1663 sg2, 1664 cprev_prox, 1665 nnext_prox, 1666 cprev_dist, 1667 nnext_dist, 1668 ) 1669 1670 # calculate and set the disulfide dihedral angles 1671 self.chi1 = np.degrees(calc_dihedral(n1, ca1, cb1, sg1)) 1672 self.chi2 = np.degrees(calc_dihedral(ca1, cb1, sg1, sg2)) 1673 self.chi3 = np.degrees(calc_dihedral(cb1, sg1, sg2, cb2)) 1674 self.chi4 = np.degrees(calc_dihedral(sg1, sg2, cb2, ca2)) 1675 self.chi5 = np.degrees(calc_dihedral(sg2, cb2, ca2, n2)) 1676 self.rho = np.degrees(calc_dihedral(n1, ca1, ca2, n2)) 1677 1678 self.ca_distance = distance3d(self.ca_prox, self.ca_dist) 1679 self.cb_distance = distance3d(self.cb_prox, self.cb_dist) 1680 self.torsion_array = np.array( 1681 (self.chi1, self.chi2, self.chi3, self.chi4, self.chi5) 1682 ) 1683 self.compute_torsion_length() 1684 1685 # calculate and set the SS bond torsional energy 1686 self.compute_torsional_energy() 1687 1688 # compute and set the local coordinates 1689 self.compute_local_coords() 1690 1691 def internal_coords(self) -> np.array: 1692 """ 1693 Return the internal coordinates for the Disulfide. 1694 If there are missing atoms the extra atoms for the proximal 1695 and distal N and C are set to [0,0,0]. This is needed for the center of 1696 mass calculations, used when rendering. 1697 1698 :return: Array containing the coordinates, [16][3]. 1699 """ 1700 1701 # if we don't have the prior and next atoms we initialize those 1702 # atoms to the origin so as to not effect the center of mass calculations 1703 if self.missing_atoms: 1704 res_array = np.array( 1705 ( 1706 self._n_prox.get_array(), 1707 self._ca_prox.get_array(), 1708 self._c_prox.get_array(), 1709 self._o_prox.get_array(), 1710 self._cb_prox.get_array(), 1711 self._sg_prox.get_array(), 1712 self._n_dist.get_array(), 1713 self._ca_dist.get_array(), 1714 self._c_dist.get_array(), 1715 self._o_dist.get_array(), 1716 self._cb_dist.get_array(), 1717 self._sg_dist.get_array(), 1718 [0, 0, 0], 1719 [0, 0, 0], 1720 [0, 0, 0], 1721 [0, 0, 0], 1722 ) 1723 ) 1724 else: 1725 res_array = np.array( 1726 ( 1727 self._n_prox.get_array(), 1728 self._ca_prox.get_array(), 1729 self._c_prox.get_array(), 1730 self._o_prox.get_array(), 1731 self._cb_prox.get_array(), 1732 self._sg_prox.get_array(), 1733 self._n_dist.get_array(), 1734 self._ca_dist.get_array(), 1735 self._c_dist.get_array(), 1736 self._o_dist.get_array(), 1737 self._cb_dist.get_array(), 1738 self._sg_dist.get_array(), 1739 self._c_prev_prox.get_array(), 1740 self._n_next_prox.get_array(), 1741 self._c_prev_dist.get_array(), 1742 self._n_next_dist.get_array(), 1743 ) 1744 ) 1745 return res_array 1746 1747 def internal_coords_res(self, resnumb) -> np.array: 1748 """ 1749 Return the internal coordinates for the internal coordinates of 1750 the given Disulfide. Missing atoms are not included. 1751 1752 :param resnumb: Residue number for disulfide 1753 :raises DisulfideConstructionWarning: Warning raised if the residue number is invalid 1754 :return: Array containing the internal coordinates for the disulfide 1755 """ 1756 res_array = np.zeros(shape=(6, 3)) 1757 1758 if resnumb == self.proximal: 1759 res_array = np.array( 1760 ( 1761 self._n_prox.get_array(), 1762 self._ca_prox.get_array(), 1763 self._c_prox.get_array(), 1764 self._o_prox.get_array(), 1765 self._cb_prox.get_array(), 1766 self._sg_prox.get_array(), 1767 ) 1768 ) 1769 return res_array 1770 1771 elif resnumb == self.distal: 1772 res_array = np.array( 1773 ( 1774 self._n_dist.get_array(), 1775 self._ca_dist.get_array(), 1776 self._c_dist.get_array(), 1777 self._o_dist.get_array(), 1778 self._cb_dist.get_array(), 1779 self._sg_dist.get_array(), 1780 ) 1781 ) 1782 return res_array 1783 else: 1784 mess = f"-> Disulfide.internal_coords(): Invalid argument. \ 1785 Unable to find residue: {resnumb} " 1786 raise DisulfideConstructionWarning(mess) 1787 1788 def make_movie( 1789 self, style="sb", fname="ssbond.mp4", verbose=False, steps=360 1790 ) -> None: 1791 """ 1792 Create an animation for ```self``` rotating one revolution about the Y axis, 1793 in the given ```style```, saving to ```filename```. 1794 1795 :param style: Rendering style, defaults to 'sb', one of: 1796 * 'sb' - split bonds 1797 * 'bs' - ball and stick 1798 * 'cpk' - CPK style 1799 * 'pd' - Proximal/Distal style - Red=proximal, Green=Distal 1800 * 'plain' - boring single color 1801 1802 :param fname: Output filename, defaults to ```ssbond.mp4``` 1803 :param verbose: Verbosity, defaults to False 1804 :param steps: Number of steps for one complete rotation, defaults to 360. 1805 """ 1806 src = self.pdb_id 1807 name = self.name 1808 enrg = self.energy 1809 1810 title = f"{src} {name}: {self.proximal}{self.proximal_chain}-{self.distal}{self.distal_chain}: {enrg:.2f} kcal/mol, Cα: {self.ca_distance:.2f} Å, Tors: {self.torsion_length:.2f}" 1811 1812 if verbose: 1813 print(f"Rendering animation to {fname}...") 1814 1815 pl = pv.Plotter(window_size=WINSIZE, off_screen=True) 1816 pl.open_movie(fname) 1817 path = pl.generate_orbital_path(n_points=steps) 1818 1819 # 1820 # pl.add_title(title=title, font_size=FONTSIZE) 1821 pl.enable_anti_aliasing("msaa") 1822 # pl.add_camera_orientation_widget() 1823 pl = self._render( 1824 pl, 1825 style=style, 1826 bondcolor=BOND_COLOR, 1827 bs_scale=BS_SCALE, 1828 spec=SPECULARITY, 1829 specpow=SPEC_POWER, 1830 ) 1831 pl.reset_camera() 1832 pl.orbit_on_path(path, write_frames=True) 1833 pl.close() 1834 1835 if verbose: 1836 print(f"Saved mp4 animation to: {fname}") 1837 1838 def pprint(self) -> None: 1839 """ 1840 Pretty print general info for the Disulfide 1841 """ 1842 s1 = self.repr_ss_info() 1843 s2 = self.repr_ss_ca_dist() 1844 s3 = self.repr_ss_conformation() 1845 s4 = self.repr_ss_torsion_length() 1846 res = f"{s1} \n{s3} \n{s2} \n{s4}>" 1847 print(res) 1848 1849 def pprint_all(self) -> None: 1850 """ 1851 Pretty print all info for the Disulfide 1852 """ 1853 s1 = self.repr_ss_info() + "\n" 1854 s2 = self.repr_ss_coords() 1855 s3 = self.repr_ss_local_coords() 1856 s4 = self.repr_ss_conformation() 1857 s5 = self.repr_chain_ids() 1858 s6 = self.repr_ss_ca_dist() 1859 s7 = self.repr_ss_torsion_length() 1860 1861 res = f"{s1} {s5} {s2} {s3} {s4}\n {s6}\n {s7}>" 1862 1863 print(res) 1864 1865 # repr functions. The class is large, so I split it up into sections 1866 def repr_ss_info(self) -> str: 1867 """ 1868 Representation for the Disulfide class 1869 """ 1870 s1 = f"<Disulfide {self.name}, Source: {self.pdb_id}, Resolution: {self.resolution} Å" 1871 return s1 1872 1873 def repr_ss_coords(self) -> str: 1874 """ 1875 Representation for Disulfide coordinates 1876 """ 1877 s2 = f"\nProximal Coordinates:\n N: {self.n_prox}\n Cα: {self.ca_prox}\n C: {self.c_prox}\n O: {self.o_prox}\n Cβ: {self.cb_prox}\n Sγ: {self.sg_prox}\n Cprev {self.c_prev_prox}\n Nnext: {self.n_next_prox}\n" 1878 s3 = f"Distal Coordinates:\n N: {self.n_dist}\n Cα: {self.ca_dist}\n C: {self.c_dist}\n O: {self.o_dist}\n Cβ: {self.cb_dist}\n Sγ: {self.sg_dist}\n Cprev {self.c_prev_dist}\n Nnext: {self.n_next_dist}\n\n" 1879 stot = f"{s2} {s3}" 1880 return stot 1881 1882 def repr_ss_conformation(self) -> str: 1883 """ 1884 Representation for Disulfide conformation 1885 """ 1886 s4 = f"Χ1-Χ5: {self.chi1:.2f}°, {self.chi2:.2f}°, {self.chi3:.2f}°, {self.chi4:.2f}° {self.chi5:.2f}°, {self.rho:.2f}°, {self.energy:.2f} kcal/mol" 1887 stot = f"{s4}" 1888 return stot 1889 1890 def repr_ss_local_coords(self) -> str: 1891 """ 1892 Representation for the Disulfide internal coordinates. 1893 """ 1894 s2i = f"Proximal Internal Coords:\n N: {self._n_prox}\n Cα: {self._ca_prox}\n C: {self._c_prox}\n O: {self._o_prox}\n Cβ: {self._cb_prox}\n Sγ: {self._sg_prox}\n Cprev {self.c_prev_prox}\n Nnext: {self.n_next_prox}\n" 1895 s3i = f"Distal Internal Coords:\n N: {self._n_dist}\n Cα: {self._ca_dist}\n C: {self._c_dist}\n O: {self._o_dist}\n Cβ: {self._cb_dist}\n Sγ: {self._sg_dist}\n Cprev {self.c_prev_dist}\n Nnext: {self.n_next_dist}\n" 1896 stot = f"{s2i}{s3i}" 1897 return stot 1898 1899 def repr_ss_chain_ids(self) -> str: 1900 """ 1901 Representation for Disulfide chain IDs 1902 """ 1903 return f"Proximal Chain fullID: <{self.proximal_residue_fullid}> Distal Chain fullID: <{self.distal_residue_fullid}>" 1904 1905 def repr_ss_ca_dist(self) -> str: 1906 """ 1907 Representation for Disulfide Ca distance 1908 """ 1909 s1 = f"Cα Distance: {self.ca_distance:.2f} Å" 1910 return s1 1911 1912 def repr_ss_cb_dist(self) -> str: 1913 """ 1914 Representation for Disulfide Ca distance 1915 """ 1916 s1 = f"Cβ Distance: {self.cb_distance:.2f} Å" 1917 return s1 1918 1919 def repr_ss_torsion_length(self) -> str: 1920 """ 1921 Representation for Disulfide torsion length 1922 """ 1923 s1 = f"Torsion length: {self.torsion_length:.2f} deg" 1924 return s1 1925 1926 def repr_all(self) -> str: 1927 """ 1928 Return a string representation for all Disulfide information 1929 contained in self. 1930 """ 1931 1932 s1 = self.repr_ss_info() + "\n" 1933 s2 = self.repr_ss_coords() 1934 s3 = self.repr_ss_local_coords() 1935 s4 = self.repr_ss_conformation() 1936 s5 = self.repr_chain_ids() 1937 s6 = self.repr_ss_ca_dist() 1938 s8 = self.repr_ss_cb_dist() 1939 s7 = self.repr_ss_torsion_length() 1940 1941 res = f"{s1} {s5} {s2} {s3} {s4} {s6} {s7} {s8}>" 1942 return res 1943 1944 def repr_compact(self) -> str: 1945 """ 1946 Return a compact representation of the Disulfide object 1947 :return: string 1948 """ 1949 return f"{self.repr_ss_info()} {self.repr_ss_conformation()}" 1950 1951 def repr_conformation(self) -> str: 1952 """ 1953 Return a string representation of the Disulfide object's conformation. 1954 :return: string 1955 """ 1956 return f"{self.repr_ss_conformation()}" 1957 1958 def repr_coords(self) -> str: 1959 """ 1960 Return a string representation of the Disulfide object's coordinates. 1961 :return: string 1962 """ 1963 return f"{self.repr_ss_coords()}" 1964 1965 def repr_internal_coords(self) -> str: 1966 """ 1967 Return a string representation of the Disulfide object's internal coordinaes. 1968 :return: string 1969 """ 1970 return f"{self.repr_ss_local_coords()}" 1971 1972 def repr_chain_ids(self) -> str: 1973 """ 1974 Return a string representation of the Disulfide object's chain ids. 1975 :return: string 1976 """ 1977 return f"{self.repr_ss_chain_ids()}" 1978 1979 def compute_rho(self) -> float: 1980 self.rho = calc_dihedral(self.n_prox, self.ca_prox, self.ca_dist, self.n_dist) 1981 return self.rho 1982 1983 def reset(self) -> None: 1984 """ 1985 Resets the disulfide object to its initial state. All distances, 1986 angles and positions are reset. The name is unchanged. 1987 """ 1988 self.__init__(self) 1989 1990 def same_chains(self) -> bool: 1991 """ 1992 Function checks if the Disulfide is cross-chain or not. 1993 1994 Returns 1995 ------- 1996 bool \n 1997 True if the proximal and distal residues are on the same chains, 1998 False otherwise. 1999 """ 2000 2001 (prox, dist) = self.get_chains() 2002 return prox == dist 2003 2004 def screenshot( 2005 self, 2006 single=True, 2007 style="sb", 2008 fname="ssbond.png", 2009 verbose=False, 2010 shadows=False, 2011 light=True, 2012 ) -> None: 2013 """ 2014 Create and save a screenshot of the Disulfide in the given style 2015 and filename 2016 2017 :param single: Display a single vs panel view, defaults to True 2018 :param style: Rendering style, one of: 2019 * 'sb' - split bonds 2020 * 'bs' - ball and stick 2021 * 'cpk' - CPK style 2022 * 'pd' - Proximal/Distal style - Red=proximal, Green=Distal 2023 * 'plain' - boring single color, 2024 :param fname: output filename,, defaults to 'ssbond.png' 2025 :param verbose: Verbosit, defaults to False 2026 :param shadows: Enable shadows, defaults to False 2027 """ 2028 src = self.pdb_id 2029 name = self.name 2030 enrg = self.energy 2031 2032 title = f"{src} {name}: {self.proximal}{self.proximal_chain}-{self.distal}{self.distal_chain}: {enrg:.2f} kcal/mol, Cα: {self.ca_distance:.2f} Å, Cβ: {self.cb_distance:.2f} Å, Tors: {self.torsion_length:.2f}" 2033 2034 if light: 2035 pv.set_plot_theme("document") 2036 else: 2037 pv.set_plot_theme("dark") 2038 2039 if verbose: 2040 print(f"-> screenshot(): Rendering screenshot to file {fname}") 2041 2042 if single: 2043 pl = pv.Plotter(window_size=WINSIZE) 2044 # pl.add_title(title=title, font_size=FONTSIZE) 2045 pl.enable_anti_aliasing("msaa") 2046 # pl.add_camera_orientation_widget() 2047 self._render( 2048 pl, 2049 style=style, 2050 bondcolor=BOND_COLOR, 2051 bs_scale=BS_SCALE, 2052 spec=SPECULARITY, 2053 specpow=SPEC_POWER, 2054 ) 2055 pl.reset_camera() 2056 if shadows: 2057 pl.enable_shadows() 2058 2059 pl.show(auto_close=False) 2060 pl.screenshot(fname) 2061 pl.clear() 2062 2063 else: 2064 pl = pv.Plotter(window_size=WINSIZE, shape=(2, 2)) 2065 pl.subplot(0, 0) 2066 2067 # pl.add_title(title=title, font_size=FONTSIZE) 2068 pl.enable_anti_aliasing("msaa") 2069 2070 # pl.add_camera_orientation_widget() 2071 self._render( 2072 pl, 2073 style="cpk", 2074 bondcolor=BOND_COLOR, 2075 bs_scale=BS_SCALE, 2076 spec=SPECULARITY, 2077 specpow=SPEC_POWER, 2078 ) 2079 2080 pl.subplot(0, 1) 2081 # pl.add_title(title=title, font_size=FONTSIZE) 2082 self._render( 2083 pl, 2084 style="pd", 2085 bondcolor=BOND_COLOR, 2086 bs_scale=BS_SCALE, 2087 spec=SPECULARITY, 2088 specpow=SPEC_POWER, 2089 ) 2090 2091 pl.subplot(1, 0) 2092 # pl.add_title(title=title, font_size=FONTSIZE) 2093 self._render( 2094 pl, 2095 style="bs", 2096 bondcolor=BOND_COLOR, 2097 bs_scale=BS_SCALE, 2098 spec=SPECULARITY, 2099 specpow=SPEC_POWER, 2100 ) 2101 2102 pl.subplot(1, 1) 2103 # pl.add_title(title=title, font_size=FONTSIZE) 2104 self._render( 2105 pl, 2106 style="sb", 2107 bondcolor=BOND_COLOR, 2108 bs_scale=BS_SCALE, 2109 spec=SPECULARITY, 2110 specpow=SPEC_POWER, 2111 ) 2112 2113 pl.link_views() 2114 pl.reset_camera() 2115 if shadows: 2116 pl.enable_shadows() 2117 2118 pl.show(auto_close=False) 2119 pl.screenshot(fname) 2120 2121 if verbose: 2122 print(f"Saved: {fname}") 2123 2124 def save_meshes_as_stl(self, meshes, filename) -> None: 2125 """Save a list of meshes as a single STL file. 2126 2127 Args: 2128 meshes (list): List of pyvista mesh objects to save. 2129 filename (str): Path to save the STL file to. 2130 """ 2131 merged_mesh = pv.UnstructuredGrid() 2132 for mesh in meshes: 2133 merged_mesh += mesh 2134 merged_mesh.save(filename) 2135 2136 def export(self, style="sb", verbose=True, fname="ssbond_plt") -> None: 2137 """ 2138 Create and save a screenshot of the Disulfide in the given style and filename. 2139 2140 :param single: Display a single vs panel view, defaults to True 2141 :param style: Rendering style, one of: 2142 * 'sb' - split bonds 2143 * 'bs' - ball and stick 2144 * 'cpk' - CPK style 2145 * 'pd' - Proximal/Distal style - Red=proximal, Green=Distal 2146 * 'plain' - boring single color, 2147 2148 :param fname: output filename,, defaults to 'ssbond.stl' 2149 :param verbose: Verbosit, defaults to False 2150 """ 2151 2152 if verbose: 2153 print(f"-> screenshot(): Rendering screenshot to file {fname}") 2154 2155 pl = pv.PolyData() 2156 2157 self._plot( 2158 pl, 2159 style=style, 2160 bondcolor=BOND_COLOR, 2161 bs_scale=BS_SCALE, 2162 spec=SPECULARITY, 2163 specpow=SPEC_POWER, 2164 ) 2165 2166 self.save_meshes_as_stl(pl, fname) 2167 2168 return 2169 2170 def set_permissive(self, perm: bool) -> None: 2171 """ 2172 Set PERMISSIVE flag for Disulfide parsing. 2173 2174 :return: None 2175 """ 2176 2177 self.PERMISSIVE = perm 2178 2179 def set_positions( 2180 self, 2181 n_prox: Vector, 2182 ca_prox: Vector, 2183 c_prox: Vector, 2184 o_prox: Vector, 2185 cb_prox: Vector, 2186 sg_prox: Vector, 2187 n_dist: Vector, 2188 ca_dist: Vector, 2189 c_dist: Vector, 2190 o_dist: Vector, 2191 cb_dist: Vector, 2192 sg_dist: Vector, 2193 c_prev_prox: Vector, 2194 n_next_prox: Vector, 2195 c_prev_dist: Vector, 2196 n_next_dist: Vector, 2197 ) -> None: 2198 """ 2199 Set the atomic coordinates for all atoms in the Disulfide object. 2200 2201 :param n_prox: Proximal N position 2202 :param ca_prox: Proximal Cα position 2203 :param c_prox: Proximal C' position 2204 :param o_prox: Proximal O position 2205 :param cb_prox: Proximal Cβ position 2206 :param sg_prox: Proximal Sγ position 2207 :param n_dist: Distal N position 2208 :param ca_dist: Distal Cα position 2209 :param c_dist: Distal C' position 2210 :param o_dist: Distal O position 2211 :param cb_dist: Distal Cβ position 2212 :param sg_dist: Distal Sγ position 2213 :param c_prev_prox: Proximal previous C' 2214 :param n_next_prox: Proximal next N 2215 :param c_prev_dist: Distal previous C' 2216 :param n_next_dist: Distal next N 2217 """ 2218 2219 # deep copy 2220 self.n_prox = n_prox.copy() 2221 self.ca_prox = ca_prox.copy() 2222 self.c_prox = c_prox.copy() 2223 self.o_prox = o_prox.copy() 2224 self.cb_prox = cb_prox.copy() 2225 self.sg_prox = sg_prox.copy() 2226 self.sg_dist = sg_dist.copy() 2227 self.cb_dist = cb_dist.copy() 2228 self.ca_dist = ca_dist.copy() 2229 self.n_dist = n_dist.copy() 2230 self.c_dist = c_dist.copy() 2231 self.o_dist = o_dist.copy() 2232 2233 self.c_prev_prox = c_prev_prox.copy() 2234 self.n_next_prox = n_next_prox.copy() 2235 self.c_prev_dist = c_prev_dist.copy() 2236 self.n_next_dist = n_next_dist.copy() 2237 2238 def set_dihedrals( 2239 self, chi1: float, chi2: float, chi3: float, chi4: float, chi5: float 2240 ) -> None: 2241 """ 2242 Set the disulfide's dihedral angles, Chi1-Chi5. -180 - 180 degrees. 2243 2244 :param chi1: Chi1 2245 :param chi2: Chi2 2246 :param chi3: Chi3 2247 :param chi4: Chi4 2248 :param chi5: Chi5 2249 """ 2250 self.chi1 = chi1 2251 self.chi2 = chi2 2252 self.chi3 = chi3 2253 self.chi4 = chi4 2254 self.chi5 = chi5 2255 self.torsion_array = np.array([chi1, chi2, chi3, chi4, chi5]) 2256 self.compute_torsional_energy() 2257 self.compute_torsion_length() 2258 2259 def set_name(self, namestr="Disulfide") -> None: 2260 """ 2261 Set the Disulfide's name. 2262 2263 :param namestr: Name, by default "Disulfide" 2264 """ 2265 self.name = namestr 2266 2267 def set_resnum(self, proximal: int, distal: int) -> None: 2268 """ 2269 Set the proximal and residue numbers for the Disulfide. 2270 2271 :param proximal: Proximal residue number 2272 :param distal: Distal residue number 2273 """ 2274 self.proximal = proximal 2275 self.distal = distal 2276 2277 def compute_torsion_length(self) -> float: 2278 """ 2279 Compute the 5D Euclidean length of the Disulfide object. Update the disulfide internal state. 2280 2281 :return: Torsion length (Degrees) 2282 """ 2283 # Use numpy array to compute element-wise square 2284 tors2 = np.square(self.torsion_array) 2285 2286 # Compute the sum of squares using numpy's sum function 2287 dist = math.sqrt(np.sum(tors2)) 2288 2289 # Update the internal state 2290 self.torsion_length = dist 2291 2292 return dist 2293 2294 def Torsion_Distance(self, other) -> float: 2295 """ 2296 Calculate the 5D Euclidean distance between ```self``` and another Disulfide 2297 object. This is used to compare Disulfide Bond torsion angles to 2298 determine their torsional similarity via a 5-Dimensional Euclidean distance metric. 2299 2300 :param other: Comparison Disulfide 2301 :raises ProteusPyWarning: Warning if ```other``` is not a Disulfide object 2302 :return: Euclidean distance (Degrees) between ```self``` and ```other```. 2303 """ 2304 2305 from proteusPy.ProteusPyWarning import ProteusPyWarning 2306 2307 # Check length of torsion arrays 2308 if len(self.torsion_array) != 5 or len(other.torsion_array) != 5: 2309 raise ProteusPyWarning( 2310 "--> Torsion_Distance() requires vectors of length 5!" 2311 ) 2312 2313 # Convert to numpy arrays and add 180 to each element 2314 p1 = np.array(self.torsion_array) + 180.0 2315 p2 = np.array(other.torsion_array) + 180.0 2316 2317 # Compute the 5D Euclidean distance using numpy's linalg.norm function 2318 dist = np.linalg.norm(p1 - p2) 2319 2320 return dist 2321 2322 def Torsion_neighbors(self, others, cutoff) -> DisulfideList: 2323 """ 2324 Return a list of Disulfides within the angular cutoff in the others list. 2325 This routine is used to find Disulfides having the same torsion length 2326 within the others list. This is used to find families of Disulfides with 2327 similar conformations. Assumes self is properly initialized. 2328 2329 *NB* The routine will not distinguish between +/- 2330 dihedral angles. *i.e.* [-60, -60, -90, -60, -60] would have the same 2331 torsion length as [60, 60, 90, 60, 60], two clearly different structures. 2332 2333 :param others: ```DisulfideList``` to search 2334 :param cutoff: Dihedral angle degree cutoff 2335 :return: DisulfideList within the cutoff 2336 2337 Example: 2338 In this example we load the disulfide database subset, find the disulfides with 2339 the lowest and highest energies, and then find the nearest conformational neighbors. 2340 Finally, we display the neighbors overlaid against a common reference frame. 2341 2342 >>> from proteusPy import * 2343 >>> from proteusPy.DisulfideLoader import Load_PDB_SS 2344 >>> from proteusPy.DisulfideList import DisulfideList 2345 >>> from proteusPy.Disulfide import Disulfide 2346 >>> PDB_SS = Load_PDB_SS(verbose=False, subset=True) 2347 >>> ss_list = DisulfideList([], 'tmp') 2348 2349 We point to the complete list to search for lowest and highest energies. 2350 >>> sslist = PDB_SS.SSList 2351 >>> ssmin_enrg, ssmax_enrg = PDB_SS.SSList.minmax_energy 2352 2353 Make an empty list and find the nearest neighbors within 10 degrees avg RMS in 2354 sidechain dihedral angle space. 2355 2356 >>> low_energy_neighbors = DisulfideList([],'Neighbors') 2357 >>> low_energy_neighbors = ssmin_enrg.Torsion_neighbors(sslist, 10) 2358 2359 Display the number found, and then display them overlaid onto their common reference frame. 2360 2361 >>> tot = low_energy_neighbors.length 2362 >>> print(f'Neighbors: {tot}') 2363 Neighbors: 2 2364 >>> low_energy_neighbors.display_overlay() 2365 2366 """ 2367 res = [ss for ss in others if self.Torsion_Distance(ss) <= cutoff] 2368 return DisulfideList(res, "neighbors") 2369 2370 def torsion_to_sixclass(self) -> str: 2371 """ 2372 Return the sextant class string for ``self``. 2373 2374 :raises DisulfideIOException: _description_ 2375 :return: Sextant string 2376 """ 2377 from proteusPy.DisulfideClasses import get_sixth_quadrant 2378 2379 tors = self.torsion_array 2380 res = [get_sixth_quadrant(x) for x in tors] 2381 return "".join([str(r) for r in res])
This class provides a Python object and methods representing a physical disulfide bond either extracted from the RCSB protein databank or built using the proteusPy.Turtle3D class. The disulfide bond is an important intramolecular stabilizing structural element and is characterized by:
- Atomic coordinates for the atoms N, Cα, Cβ, C', Sγ for both residues. These are stored as both raw atomic coordinates as read from the RCSB file and internal local coordinates.
- The dihedral angles Χ1 - Χ5 for the disulfide bond
- A name, by default {pdb_id}{prox_resnumb}{prox_chain}_{distal_resnum}{distal_chain}
- Proximal residue number
- Distal residue number
- Approximate bond torsional energy (kcal/mol):
$$ E_{kcal/mol} \approx 2.0 * cos(3.0 * \chi_{1}) + cos(3.0 * \chi_{5}) + cos(3.0 * \chi_{2}) + $$ $$ cos(3.0 * \chi_{4}) + 3.5 * cos(2.0 * \chi_{3}) + 0.6 * cos(3.0 * \chi_{3}) + 10.1 $$
The equation embodies the typical 3-fold rotation barriers associated with single bonds, (Χ1, Χ5, Χ2, Χ4) and a high 2-fold barrier for Χ3, resulting from the partial double bond character of the S-S bond. This property leads to two major disulfide families, characterized by the sign of Χ3. Left-handed disulfides have Χ3 < 0° and right-handed disulfides have Χ3 > 0°. Within this breakdown there are numerous subfamilies, broadly known as the hook, spiral and staple. These are under characgterization.
- Euclidean length of the dihedral angles (degrees) defined as: $$\sqrt(\chi_{1}^{2} + \chi_{2}^{2} + \chi_{3}^{2} + \chi_{4}^{2} + \chi_{5}^{2})$$
- Cα - Cα distance (Å)
- Cβ - Cβ distance (Å)
- The previous C' and next N for both the proximal and distal residues. These are needed to calculate the backbone dihedral angles Φ and Ψ.
- Backbone dihedral angles Φ and Ψ, when possible. Not all structures are complete and in those cases the atoms needed may be undefined. In this case the Φ and Ψ angles are set to -180°.
The class also provides a rendering capabilities using the excellent PyVista library, and can display disulfides interactively in a variety of display styles:
- 'sb' - Split Bonds style - bonds colored by their atom type
- 'bs' - Ball and Stick style - split bond coloring with small atoms
- 'pd' - Proximal/Distal style - bonds colored Red for proximal residue and Green for the distal residue.
- 'cpk' - CPK style rendering, colored by atom type:
- Carbon - Grey
- Nitrogen - Blue
- Sulfur - Yellow
- Oxygen - Red
- Hydrogen - White
Individual renderings can be saved to a file, and animations created.
Disulfide( name: str = 'SSBOND', proximal: int = -1, distal: int = -1, proximal_chain: str = 'A', distal_chain: str = 'A', pdb_id: str = '1egs', quiet: bool = True, permissive: bool = True)
147 def __init__( 148 self, 149 name: str = "SSBOND", 150 proximal: int = -1, 151 distal: int = -1, 152 proximal_chain: str = "A", 153 distal_chain: str = "A", 154 pdb_id: str = "1egs", 155 quiet: bool = True, 156 permissive: bool = True, 157 ) -> None: 158 """ 159 __init__ Initialize the class to defined internal values. 160 161 :param name: Disulfide name, by default "SSBOND" 162 163 """ 164 self.name = name 165 self.proximal = proximal 166 self.distal = distal 167 self.energy = _FLOAT_INIT 168 self.proximal_chain = proximal_chain 169 self.distal_chain = distal_chain 170 self.pdb_id = pdb_id 171 self.proximal_residue_fullid = str("") 172 self.distal_residue_fullid = str("") 173 self.PERMISSIVE = permissive 174 self.QUIET = quiet 175 self.ca_distance = _FLOAT_INIT 176 self.cb_distance = _FLOAT_INIT 177 self.torsion_array = np.array( 178 (_ANG_INIT, _ANG_INIT, _ANG_INIT, _ANG_INIT, _ANG_INIT) 179 ) 180 self.phiprox = _ANG_INIT 181 self.psiprox = _ANG_INIT 182 self.phidist = _ANG_INIT 183 self.psidist = _ANG_INIT 184 185 # global coordinates for the Disulfide, typically as 186 # returned from the PDB file 187 188 self.n_prox = Vector(0, 0, 0) 189 self.ca_prox = Vector(0, 0, 0) 190 self.c_prox = Vector(0, 0, 0) 191 self.o_prox = Vector(0, 0, 0) 192 self.cb_prox = Vector(0, 0, 0) 193 self.sg_prox = Vector(0, 0, 0) 194 self.sg_dist = Vector(0, 0, 0) 195 self.cb_dist = Vector(0, 0, 0) 196 self.ca_dist = Vector(0, 0, 0) 197 self.n_dist = Vector(0, 0, 0) 198 self.c_dist = Vector(0, 0, 0) 199 self.o_dist = Vector(0, 0, 0) 200 201 # set when we can't find previous or next prox or distal 202 # C' or N atoms. 203 self.missing_atoms = False 204 self.modelled = False 205 self.resolution = -1.0 206 207 # need these to calculate backbone dihedral angles 208 self.c_prev_prox = Vector(0, 0, 0) 209 self.n_next_prox = Vector(0, 0, 0) 210 self.c_prev_dist = Vector(0, 0, 0) 211 self.n_next_dist = Vector(0, 0, 0) 212 213 # local coordinates for the Disulfide, computed using the Turtle3D in 214 # Orientation #1. these are generally private. 215 216 self._n_prox = Vector(0, 0, 0) 217 self._ca_prox = Vector(0, 0, 0) 218 self._c_prox = Vector(0, 0, 0) 219 self._o_prox = Vector(0, 0, 0) 220 self._cb_prox = Vector(0, 0, 0) 221 self._sg_prox = Vector(0, 0, 0) 222 self._sg_dist = Vector(0, 0, 0) 223 self._cb_dist = Vector(0, 0, 0) 224 self._ca_dist = Vector(0, 0, 0) 225 self._n_dist = Vector(0, 0, 0) 226 self._c_dist = Vector(0, 0, 0) 227 self._o_dist = Vector(0, 0, 0) 228 229 # need these to calculate backbone dihedral angles 230 self._c_prev_prox = Vector(0, 0, 0) 231 self._n_next_prox = Vector(0, 0, 0) 232 self._c_prev_dist = Vector(0, 0, 0) 233 self._n_next_dist = Vector(0, 0, 0) 234 235 # Dihedral angles for the disulfide bond itself, set to _ANG_INIT 236 self.chi1 = _ANG_INIT 237 self.chi2 = _ANG_INIT 238 self.chi3 = _ANG_INIT 239 self.chi4 = _ANG_INIT 240 self.chi5 = _ANG_INIT 241 self.rho = _ANG_INIT # new dihedral angle: Nprox - Ca_prox - Ca_dist - N_dist 242 243 self.torsion_length = _FLOAT_INIT
__init__ Initialize the class to defined internal values.
Parameters
- name: Disulfide name, by default "SSBOND"
dihedrals: list
986 @property 987 def dihedrals(self) -> list: 988 """ 989 Return a list containing the dihedral angles for the disulfide. 990 991 """ 992 return [self.chi1, self.chi2, self.chi3, self.chi4, self.chi5]
Return a list containing the dihedral angles for the disulfide.
def
bounding_box(self) -> <built-in function array>:
1009 def bounding_box(self) -> np.array: 1010 """ 1011 Return the bounding box array for the given disulfide 1012 1013 Returns 1014 ------- 1015 :return: np.Array(3,2): Array containing the min, max for X, Y, and Z respectively. 1016 Does not currently take the atom's radius into account. 1017 1018 """ 1019 res = np.zeros(shape=(3, 2)) 1020 xmin, xmax = self.compute_extents("x") 1021 ymin, ymax = self.compute_extents("y") 1022 zmin, zmax = self.compute_extents("z") 1023 1024 res[0] = [xmin, xmax] 1025 res[1] = [ymin, ymax] 1026 res[2] = [zmin, zmax] 1027 1028 return res
Return the bounding box array for the given disulfide
Returns
Returns
np.Array(3,2): Array containing the min, max for X, Y, and Z respectively. Does not currently take the atom's radius into account.
def
build_yourself(self) -> None:
1030 def build_yourself(self) -> None: 1031 """ 1032 Build a model Disulfide based its internal dihedral state 1033 Routine assumes turtle is in orientation #1 (at Ca, headed toward 1034 Cb, with N on left), builds disulfide, and updates the object's internal 1035 coordinates. It also adds the distal protein backbone, 1036 and computes the disulfide conformational energy. 1037 """ 1038 chi1 = self.chi1 1039 chi2 = self.chi2 1040 chi3 = self.chi3 1041 chi4 = self.chi4 1042 chi5 = self.chi5 1043 self.build_model(chi1, chi2, chi3, chi4, chi5)
Build a model Disulfide based its internal dihedral state Routine assumes turtle is in orientation #1 (at Ca, headed toward Cb, with N on left), builds disulfide, and updates the object's internal coordinates. It also adds the distal protein backbone, and computes the disulfide conformational energy.
def
build_model( self, chi1: float, chi2: float, chi3: float, chi4: float, chi5: float) -> None:
1045 def build_model( 1046 self, chi1: float, chi2: float, chi3: float, chi4: float, chi5: float 1047 ) -> None: 1048 """ 1049 Build a model Disulfide based on the input dihedral angles. 1050 Routine assumes turtle is in orientation #1 (at Ca, headed toward 1051 Cb, with N on left), builds disulfide, and updates the object's internal 1052 coordinates. It also adds the distal protein backbone, 1053 and computes the disulfide conformational energy. 1054 1055 :param chi1: Chi1 (degrees) 1056 :param chi2: Chi2 (degrees) 1057 :param chi3: Chi3 (degrees) 1058 :param chi4: Chi4 (degrees) 1059 :param chi5: Chi5 (degrees) 1060 1061 Example: 1062 >>> from proteusPy.Disulfide import Disulfide 1063 >>> modss = Disulfide('model') 1064 >>> modss.build_model(-60, -60, -90, -60, -60) 1065 >>> modss.display(style='sb') 1066 """ 1067 1068 self.set_dihedrals(chi1, chi2, chi3, chi4, chi5) 1069 self.proximal = 1 1070 self.distal = 2 1071 1072 tmp = Turtle3D("tmp") 1073 tmp.Orientation = 1 1074 1075 n = Vector(0, 0, 0) 1076 ca = Vector(0, 0, 0) 1077 cb = Vector(0, 0, 0) 1078 c = Vector(0, 0, 0) 1079 1080 self.ca_prox = tmp._position 1081 tmp.schain_to_bbone() 1082 n, ca, cb, c = build_residue(tmp) 1083 1084 self.n_prox = n 1085 self.ca_prox = ca 1086 self.c_prox = c 1087 self.cb_prox = cb 1088 1089 tmp.bbone_to_schain() 1090 tmp.move(1.53) 1091 tmp.roll(self.chi1) 1092 tmp.yaw(112.8) 1093 self.cb_prox = Vector(tmp._position) 1094 1095 tmp.move(1.86) 1096 tmp.roll(self.chi2) 1097 tmp.yaw(103.8) 1098 self.sg_prox = Vector(tmp._position) 1099 1100 tmp.move(2.044) 1101 tmp.roll(self.chi3) 1102 tmp.yaw(103.8) 1103 self.sg_dist = Vector(tmp._position) 1104 1105 tmp.move(1.86) 1106 tmp.roll(self.chi4) 1107 tmp.yaw(112.8) 1108 self.cb_dist = Vector(tmp._position) 1109 1110 tmp.move(1.53) 1111 tmp.roll(self.chi5) 1112 tmp.pitch(180.0) 1113 1114 tmp.schain_to_bbone() 1115 1116 n, ca, cb, c = build_residue(tmp) 1117 1118 self.n_dist = n 1119 self.ca_dist = ca 1120 self.c_dist = c 1121 self.compute_torsional_energy() 1122 self.compute_local_coords() 1123 self.ca_distance = distance3d(self.ca_prox, self.ca_dist) 1124 self.cb_distance = distance3d(self.cb_prox, self.cb_dist) 1125 self.torsion_array = np.array( 1126 (self.chi1, self.chi2, self.chi3, self.chi4, self.chi5) 1127 ) 1128 self.compute_torsion_length() 1129 self.compute_rho() 1130 self.missing_atoms = True 1131 self.modelled = True
Build a model Disulfide based on the input dihedral angles. Routine assumes turtle is in orientation #1 (at Ca, headed toward Cb, with N on left), builds disulfide, and updates the object's internal coordinates. It also adds the distal protein backbone, and computes the disulfide conformational energy.
Parameters
- chi1: Chi1 (degrees)
- chi2: Chi2 (degrees)
- chi3: Chi3 (degrees)
- chi4: Chi4 (degrees)
- chi5: Chi5 (degrees)
Example:
>>> from proteusPy.Disulfide import Disulfide
>>> modss = Disulfide('model')
>>> modss.build_model(-60, -60, -90, -60, -60)
>>> modss.display(style='sb')
def
cofmass(self) -> <built-in function array>:
1133 def cofmass(self) -> np.array: 1134 """ 1135 Return the geometric center of mass for the internal coordinates of 1136 the given Disulfide. Missing atoms are not included. 1137 1138 :return: 3D array for the geometric center of mass 1139 """ 1140 1141 res = self.internal_coords() 1142 return res.mean(axis=0)
Return the geometric center of mass for the internal coordinates of the given Disulfide. Missing atoms are not included.
Returns
3D array for the geometric center of mass
def
copy(self):
1144 def copy(self): 1145 """ 1146 Copy the Disulfide. 1147 1148 :return: A copy of self. 1149 """ 1150 return copy.deepcopy(self)
Copy the Disulfide.
Returns
A copy of self.
def
compute_extents(self, dim='z'):
1152 def compute_extents(self, dim="z"): 1153 """ 1154 Calculate the internal coordinate extents for the input axis. 1155 1156 :param dim: Axis, one of 'x', 'y', 'z', by default 'z' 1157 :return: min, max 1158 """ 1159 1160 ic = self.internal_coords() 1161 # set default index to 'z' 1162 idx = 2 1163 1164 if dim == "x": 1165 idx = 0 1166 elif dim == "y": 1167 idx = 1 1168 elif dim == "z": 1169 idx = 2 1170 1171 _min = ic[:, idx].min() 1172 _max = ic[:, idx].max() 1173 return _min, _max
Calculate the internal coordinate extents for the input axis.
Parameters
- dim: Axis, one of 'x', 'y', 'z', by default 'z'
Returns
min, max
def
compute_local_coords(self) -> None:
1175 def compute_local_coords(self) -> None: 1176 """ 1177 Compute the internal coordinates for a properly initialized Disulfide Object. 1178 1179 :param self: SS initialized Disulfide object 1180 :returns: None, modifies internal state of the input 1181 """ 1182 1183 turt = Turtle3D("tmp") 1184 # get the coordinates as np.array for Turtle3D use. 1185 cpp = self.c_prev_prox.get_array() 1186 nnp = self.n_next_prox.get_array() 1187 1188 n = self.n_prox.get_array() 1189 ca = self.ca_prox.get_array() 1190 c = self.c_prox.get_array() 1191 cb = self.cb_prox.get_array() 1192 o = self.o_prox.get_array() 1193 sg = self.sg_prox.get_array() 1194 1195 sg2 = self.sg_dist.get_array() 1196 cb2 = self.cb_dist.get_array() 1197 ca2 = self.ca_dist.get_array() 1198 c2 = self.c_dist.get_array() 1199 n2 = self.n_dist.get_array() 1200 o2 = self.o_dist.get_array() 1201 1202 cpd = self.c_prev_dist.get_array() 1203 nnd = self.n_next_dist.get_array() 1204 1205 turt.orient_from_backbone(n, ca, c, cb, ORIENT_SIDECHAIN) 1206 1207 # internal (local) coordinates, stored as Vector objects 1208 # to_local returns np.array objects 1209 1210 self._n_prox = Vector(turt.to_local(n)) 1211 self._ca_prox = Vector(turt.to_local(ca)) 1212 self._c_prox = Vector(turt.to_local(c)) 1213 self._o_prox = Vector(turt.to_local(o)) 1214 self._cb_prox = Vector(turt.to_local(cb)) 1215 self._sg_prox = Vector(turt.to_local(sg)) 1216 1217 self._c_prev_prox = Vector(turt.to_local(cpp)) 1218 self._n_next_prox = Vector(turt.to_local(nnp)) 1219 self._c_prev_dist = Vector(turt.to_local(cpd)) 1220 self._n_next_dist = Vector(turt.to_local(nnd)) 1221 1222 self._n_dist = Vector(turt.to_local(n2)) 1223 self._ca_dist = Vector(turt.to_local(ca2)) 1224 self._c_dist = Vector(turt.to_local(c2)) 1225 self._o_dist = Vector(turt.to_local(o2)) 1226 self._cb_dist = Vector(turt.to_local(cb2)) 1227 self._sg_dist = Vector(turt.to_local(sg2))
Compute the internal coordinates for a properly initialized Disulfide Object.
Parameters
- self: SS initialized Disulfide object :returns: None, modifies internal state of the input
def
compute_torsional_energy(self) -> float:
1229 def compute_torsional_energy(self) -> float: 1230 """ 1231 Compute the approximate torsional energy for the Disulfide's 1232 conformation and sets its internal state. 1233 1234 :return: Energy (kcal/mol) 1235 """ 1236 # @TODO find citation for the ss bond energy calculation 1237 1238 def torad(deg): 1239 return np.radians(deg) 1240 1241 chi1 = self.chi1 1242 chi2 = self.chi2 1243 chi3 = self.chi3 1244 chi4 = self.chi4 1245 chi5 = self.chi5 1246 1247 energy = 2.0 * (cos(torad(3.0 * chi1)) + cos(torad(3.0 * chi5))) 1248 energy += cos(torad(3.0 * chi2)) + cos(torad(3.0 * chi4)) 1249 energy += 3.5 * cos(torad(2.0 * chi3)) + 0.6 * cos(torad(3.0 * chi3)) + 10.1 1250 1251 self.energy = energy 1252 return energy
Compute the approximate torsional energy for the Disulfide's conformation and sets its internal state.
Returns
Energy (kcal/mol)
def
display(self, single=True, style='sb', light=True, shadows=False) -> None:
1254 def display(self, single=True, style="sb", light=True, shadows=False) -> None: 1255 """ 1256 Display the Disulfide bond in the specific rendering style. 1257 1258 :param single: Display the bond in a single panel in the specific style. 1259 :param style: Rendering style: One of: 1260 * 'sb' - split bonds 1261 * 'bs' - ball and stick 1262 * 'cpk' - CPK style 1263 * 'pd' - Proximal/Distal style - Red=proximal, Green=Distal 1264 * 'plain' - boring single color 1265 :param light: If True, light background, if False, dark 1266 1267 Example: 1268 >>> import proteusPy 1269 >>> from proteusPy.Disulfide import Disulfide 1270 >>> from proteusPy.DisulfideLoader import DisulfideLoader, Load_PDB_SS 1271 1272 >>> PDB_SS = Load_PDB_SS(verbose=False, subset=True) 1273 >>> ss = PDB_SS[0] 1274 >>> ss.display(style='cpk') 1275 >>> ss.screenshot(style='bs', fname='proteus_logo_sb.png') 1276 """ 1277 src = self.pdb_id 1278 enrg = self.energy 1279 1280 title = f"{src}: {self.proximal}{self.proximal_chain}-{self.distal}{self.distal_chain}: {enrg:.2f} kcal/mol. Cα: {self.ca_distance:.2f} Å Cβ: {self.cb_distance:.2f} Å Tors: {self.torsion_length:.2f}°" 1281 1282 if light: 1283 pv.set_plot_theme("document") 1284 else: 1285 pv.set_plot_theme("dark") 1286 1287 if single == True: 1288 _pl = pv.Plotter(window_size=WINSIZE) 1289 _pl.add_title(title=title, font_size=FONTSIZE) 1290 _pl.enable_anti_aliasing("msaa") 1291 # _pl.add_camera_orientation_widget() 1292 1293 self._render( 1294 _pl, 1295 style=style, 1296 bs_scale=BS_SCALE, 1297 spec=SPECULARITY, 1298 specpow=SPEC_POWER, 1299 ) 1300 _pl.reset_camera() 1301 if shadows == True: 1302 _pl.enable_shadows() 1303 _pl.show() 1304 1305 else: 1306 pl = pv.Plotter(window_size=WINSIZE, shape=(2, 2)) 1307 pl.subplot(0, 0) 1308 1309 pl.add_title(title=title, font_size=FONTSIZE) 1310 pl.enable_anti_aliasing("msaa") 1311 1312 # pl.add_camera_orientation_widget() 1313 1314 self._render( 1315 pl, 1316 style="cpk", 1317 bondcolor=BOND_COLOR, 1318 bs_scale=BS_SCALE, 1319 spec=SPECULARITY, 1320 specpow=SPEC_POWER, 1321 ) 1322 1323 pl.subplot(0, 1) 1324 pl.add_title(title=title, font_size=FONTSIZE) 1325 1326 self._render( 1327 pl, 1328 style="bs", 1329 bondcolor=BOND_COLOR, 1330 bs_scale=BS_SCALE, 1331 spec=SPECULARITY, 1332 specpow=SPEC_POWER, 1333 ) 1334 1335 pl.subplot(1, 0) 1336 pl.add_title(title=title, font_size=FONTSIZE) 1337 1338 self._render( 1339 pl, 1340 style="sb", 1341 bondcolor=BOND_COLOR, 1342 bs_scale=BS_SCALE, 1343 spec=SPECULARITY, 1344 specpow=SPEC_POWER, 1345 ) 1346 1347 pl.subplot(1, 1) 1348 pl.add_title(title=title, font_size=FONTSIZE) 1349 1350 self._render( 1351 pl, 1352 style="pd", 1353 bondcolor=BOND_COLOR, 1354 bs_scale=BS_SCALE, 1355 spec=SPECULARITY, 1356 specpow=SPEC_POWER, 1357 ) 1358 1359 pl.link_views() 1360 pl.reset_camera() 1361 if shadows == True: 1362 pl.enable_shadows() 1363 pl.show() 1364 return
Display the Disulfide bond in the specific rendering style.
Parameters
- single: Display the bond in a single panel in the specific style.
- style: Rendering style: One of:
- 'sb' - split bonds
- 'bs' - ball and stick
- 'cpk' - CPK style
- 'pd' - Proximal/Distal style - Red=proximal, Green=Distal
- 'plain' - boring single color
- light: If True, light background, if False, dark
Example:
>>> import proteusPy
>>> from proteusPy.Disulfide import Disulfide
>>> from proteusPy.DisulfideLoader import DisulfideLoader, Load_PDB_SS
>>> PDB_SS = Load_PDB_SS(verbose=False, subset=True)
>>> ss = PDB_SS[0]
>>> ss.display(style='cpk')
>>> ss.screenshot(style='bs', fname='proteus_logo_sb.png')
def
plot( self, pl, single=True, style='sb', light=True, shadows=False) -> pyvista.plotting.plotter.Plotter:
1366 def plot( 1367 self, pl, single=True, style="sb", light=True, shadows=False 1368 ) -> pv.Plotter: 1369 """ 1370 Return the pyVista Plotter object for the Disulfide bond in the specific rendering style. 1371 1372 :param single: Display the bond in a single panel in the specific style. 1373 :param style: Rendering style: One of: 1374 * 'sb' - split bonds 1375 * 'bs' - ball and stick 1376 * 'cpk' - CPK style 1377 * 'pd' - Proximal/Distal style - Red=proximal, Green=Distal 1378 * 'plain' - boring single color 1379 :param light: If True, light background, if False, dark 1380 """ 1381 src = self.pdb_id 1382 enrg = self.energy 1383 1384 title = f"{src}: {self.proximal}{self.proximal_chain}-{self.distal}{self.distal_chain}: {enrg:.2f} kcal/mol. Cα: {self.ca_distance:.2f} Å Cβ: {self.cb_distance:.2f} Å Tors: {self.torsion_length:.2f}°" 1385 1386 if light: 1387 pv.set_plot_theme("document") 1388 else: 1389 pv.set_plot_theme("dark") 1390 1391 if single == True: 1392 # _pl = pv.Plotter(window_size=WINSIZE) 1393 # _pl.add_title(title=title, font_size=FONTSIZE) 1394 pl.clear() 1395 pl.enable_anti_aliasing("msaa") 1396 # pl.add_camera_orientation_widget() 1397 1398 self._render( 1399 pl, 1400 style=style, 1401 bs_scale=BS_SCALE, 1402 spec=SPECULARITY, 1403 specpow=SPEC_POWER, 1404 ) 1405 pl.reset_camera() 1406 if shadows == True: 1407 pl.enable_shadows() 1408 else: 1409 pl = pv.Plotter(shape=(2, 2)) 1410 pl.subplot(0, 0) 1411 1412 # pl.add_title(title=title, font_size=FONTSIZE) 1413 pl.enable_anti_aliasing("msaa") 1414 1415 # pl.add_camera_orientation_widget() 1416 1417 self._render( 1418 pl, 1419 style="cpk", 1420 bondcolor=BOND_COLOR, 1421 bs_scale=BS_SCALE, 1422 spec=SPECULARITY, 1423 specpow=SPEC_POWER, 1424 ) 1425 1426 pl.subplot(0, 1) 1427 1428 self._render( 1429 pl, 1430 style="bs", 1431 bondcolor=BOND_COLOR, 1432 bs_scale=BS_SCALE, 1433 spec=SPECULARITY, 1434 specpow=SPEC_POWER, 1435 ) 1436 1437 pl.subplot(1, 0) 1438 1439 self._render( 1440 pl, 1441 style="sb", 1442 bondcolor=BOND_COLOR, 1443 bs_scale=BS_SCALE, 1444 spec=SPECULARITY, 1445 specpow=SPEC_POWER, 1446 ) 1447 1448 pl.subplot(1, 1) 1449 self._render( 1450 pl, 1451 style="pd", 1452 bondcolor=BOND_COLOR, 1453 bs_scale=BS_SCALE, 1454 spec=SPECULARITY, 1455 specpow=SPEC_POWER, 1456 ) 1457 1458 pl.link_views() 1459 pl.reset_camera() 1460 if shadows == True: 1461 pl.enable_shadows() 1462 return pl
Return the pyVista Plotter object for the Disulfide bond in the specific rendering style.
Parameters
- single: Display the bond in a single panel in the specific style.
- style: Rendering style: One of:
- 'sb' - split bonds
- 'bs' - ball and stick
- 'cpk' - CPK style
- 'pd' - Proximal/Distal style - Red=proximal, Green=Distal
- 'plain' - boring single color
- light: If True, light background, if False, dark
def
Distance_neighbors( self, others: proteusPy.DisulfideList.DisulfideList, cutoff: float) -> proteusPy.DisulfideList.DisulfideList:
1464 def Distance_neighbors(self, others: DisulfideList, cutoff: float) -> DisulfideList: 1465 """ 1466 Return list of Disulfides whose RMS atomic distance is within 1467 the cutoff (Å) in the others list. 1468 1469 :param others: DisulfideList to search 1470 :param cutoff: Distance cutoff (Å) 1471 :return: DisulfideList within the cutoff 1472 """ 1473 1474 res = [ss.copy() for ss in others if self.Distance_RMS(ss) < cutoff] 1475 return DisulfideList(res, "neighbors")
Return list of Disulfides whose RMS atomic distance is within the cutoff (Å) in the others list.
Parameters
- others: DisulfideList to search
- cutoff: Distance cutoff (Å)
Returns
DisulfideList within the cutoff
def
Distance_RMS(self, other) -> float:
1477 def Distance_RMS(self, other) -> float: 1478 """ 1479 Calculate the RMS distance between the internal coordinates of self and another Disulfide. 1480 :param other: Comparison Disulfide 1481 :return: RMS distance (Å) 1482 """ 1483 1484 # Get internal coordinates of both objects 1485 ic1 = self.internal_coords() 1486 ic2 = other.internal_coords() 1487 1488 # Compute the sum of squared differences between corresponding internal coordinates 1489 totsq = sum(math.dist(p1, p2) ** 2 for p1, p2 in zip(ic1, ic2)) 1490 1491 # Compute the mean of the squared distances 1492 totsq /= len(ic1) 1493 1494 # Take the square root of the mean to get the RMS distance 1495 return math.sqrt(totsq)
Calculate the RMS distance between the internal coordinates of self and another Disulfide.
Parameters
- other: Comparison Disulfide
Returns
RMS distance (Å)
def
Torsion_RMS(self, other) -> float:
1497 def Torsion_RMS(self, other) -> float: 1498 """ 1499 Calculate the RMS distance between the dihedral angles of self and another Disulfide. 1500 :param other: Comparison Disulfide 1501 :return: RMS distance (deg). 1502 """ 1503 import math 1504 1505 # Get internal coordinates of both objects 1506 ic1 = self.torsion_array 1507 ic2 = other.torsion_array 1508 1509 # Compute the sum of squared differences between corresponding internal coordinates 1510 totsq = sum((p1 - p2) ** 2 for p1, p2 in zip(ic1, ic2)) 1511 # Compute the mean of the squared distances 1512 totsq /= len(ic1) 1513 1514 # Take the square root of the mean to get the RMS distance 1515 return math.sqrt(totsq)
Calculate the RMS distance between the dihedral angles of self and another Disulfide.
Parameters
- other: Comparison Disulfide
Returns
RMS distance (deg).
def
get_chains(self) -> tuple:
1517 def get_chains(self) -> tuple: 1518 """ 1519 Return the proximal and distal chain IDs for the Disulfide. 1520 1521 :return: tuple (proximal, distal) chain IDs 1522 """ 1523 prox = self.proximal_chain 1524 dist = self.distal_chain 1525 return tuple(prox, dist)
Return the proximal and distal chain IDs for the Disulfide.
Returns
tuple (proximal, distal) chain IDs
def
get_permissive(self) -> bool:
1527 def get_permissive(self) -> bool: 1528 """ 1529 Return the Permissive flag state. (Used in PDB parsing) 1530 1531 :return: Permissive state 1532 """ 1533 return self.PERMISIVE
Return the Permissive flag state. (Used in PDB parsing)
Returns
Permissive state
def
get_full_id(self) -> tuple:
1535 def get_full_id(self) -> tuple: 1536 """ 1537 Return the Disulfide full IDs (Used with BIO.PDB) 1538 1539 :return: Disulfide full IDs 1540 """ 1541 return (self.proximal_residue_fullid, self.distal_residue_fullid)
Return the Disulfide full IDs (Used with BIO.PDB)
Returns
Disulfide full IDs
def
initialize_disulfide_from_chain(self, chain1, chain2, proximal, distal, resolution, quiet=True) -> None:
1543 def initialize_disulfide_from_chain( 1544 self, chain1, chain2, proximal, distal, resolution, quiet=True 1545 ) -> None: 1546 """ 1547 Initialize a new Disulfide object with atomic coordinates from 1548 the proximal and distal coordinates, typically taken from a PDB file. 1549 This routine is primarily used internally when building the compressed 1550 database. 1551 1552 :param chain1: list of Residues in the model, eg: chain = model['A'] 1553 :param chain2: list of Residues in the model, eg: chain = model['A'] 1554 :param proximal: proximal residue sequence ID 1555 :param distal: distal residue sequence ID 1556 :param resolution: structure resolution 1557 :param quiet: Quiet or noisy parsing, defaults to True 1558 :raises DisulfideConstructionWarning: Raised when not parsed correctly 1559 """ 1560 id = chain1.get_full_id()[0] 1561 self.pdb_id = id 1562 1563 chi1 = chi2 = chi3 = chi4 = chi5 = _ANG_INIT 1564 1565 prox = int(proximal) 1566 dist = int(distal) 1567 1568 prox_residue = chain1[prox] 1569 dist_residue = chain2[dist] 1570 1571 if prox_residue.get_resname() != "CYS" or dist_residue.get_resname() != "CYS": 1572 print( 1573 f"build_disulfide() requires CYS at both residues: {prox} {prox_residue.get_resname()} {dist} {dist_residue.get_resname()} Chain: {prox_residue.get_segid()}" 1574 ) 1575 1576 # set the objects proximal and distal values 1577 self.set_resnum(proximal, distal) 1578 1579 if resolution is not None: 1580 self.resolution = resolution 1581 1582 self.proximal_chain = chain1.get_id() 1583 self.distal_chain = chain2.get_id() 1584 1585 self.proximal_residue_fullid = prox_residue.get_full_id() 1586 self.distal_residue_fullid = dist_residue.get_full_id() 1587 1588 if quiet: 1589 warnings.filterwarnings("ignore", category=DisulfideConstructionWarning) 1590 else: 1591 warnings.simplefilter("always") 1592 1593 # grab the coordinates for the proximal and distal residues as vectors 1594 # so we can do math on them later 1595 # proximal residue 1596 1597 try: 1598 n1 = prox_residue["N"].get_vector() 1599 ca1 = prox_residue["CA"].get_vector() 1600 c1 = prox_residue["C"].get_vector() 1601 o1 = prox_residue["O"].get_vector() 1602 cb1 = prox_residue["CB"].get_vector() 1603 sg1 = prox_residue["SG"].get_vector() 1604 1605 except Exception: 1606 raise DisulfideConstructionWarning( 1607 f"Invalid or missing coordinates for proximal residue {proximal}" 1608 ) from None 1609 1610 # distal residue 1611 try: 1612 n2 = dist_residue["N"].get_vector() 1613 ca2 = dist_residue["CA"].get_vector() 1614 c2 = dist_residue["C"].get_vector() 1615 o2 = dist_residue["O"].get_vector() 1616 cb2 = dist_residue["CB"].get_vector() 1617 sg2 = dist_residue["SG"].get_vector() 1618 1619 except Exception: 1620 raise DisulfideConstructionWarning( 1621 f"Invalid or missing coordinates for distal residue {distal}" 1622 ) from None 1623 1624 # previous residue and next residue - optional, used for phi, psi calculations 1625 try: 1626 prevprox = chain1[prox - 1] 1627 nextprox = chain1[prox + 1] 1628 1629 prevdist = chain2[dist - 1] 1630 nextdist = chain2[dist + 1] 1631 1632 cprev_prox = prevprox["C"].get_vector() 1633 nnext_prox = nextprox["N"].get_vector() 1634 1635 cprev_dist = prevdist["C"].get_vector() 1636 nnext_dist = nextdist["N"].get_vector() 1637 1638 # compute phi, psi for prox and distal 1639 self.phiprox = np.degrees(calc_dihedral(cprev_prox, n1, ca1, c1)) 1640 self.psiprox = np.degrees(calc_dihedral(n1, ca1, c1, nnext_prox)) 1641 self.phidist = np.degrees(calc_dihedral(cprev_dist, n2, ca2, c2)) 1642 self.psidist = np.degrees(calc_dihedral(n2, ca2, c2, nnext_dist)) 1643 1644 except Exception: 1645 mess = f"Missing coords for: {id} {prox-1} or {dist+1} for SS {proximal}-{distal}" 1646 cprev_prox = nnext_prox = cprev_dist = nnext_dist = Vector(-1.0, -1.0, -1.0) 1647 self.missing_atoms = True 1648 warnings.warn(mess, DisulfideConstructionWarning) 1649 1650 # update the positions and conformation 1651 self.set_positions( 1652 n1, 1653 ca1, 1654 c1, 1655 o1, 1656 cb1, 1657 sg1, 1658 n2, 1659 ca2, 1660 c2, 1661 o2, 1662 cb2, 1663 sg2, 1664 cprev_prox, 1665 nnext_prox, 1666 cprev_dist, 1667 nnext_dist, 1668 ) 1669 1670 # calculate and set the disulfide dihedral angles 1671 self.chi1 = np.degrees(calc_dihedral(n1, ca1, cb1, sg1)) 1672 self.chi2 = np.degrees(calc_dihedral(ca1, cb1, sg1, sg2)) 1673 self.chi3 = np.degrees(calc_dihedral(cb1, sg1, sg2, cb2)) 1674 self.chi4 = np.degrees(calc_dihedral(sg1, sg2, cb2, ca2)) 1675 self.chi5 = np.degrees(calc_dihedral(sg2, cb2, ca2, n2)) 1676 self.rho = np.degrees(calc_dihedral(n1, ca1, ca2, n2)) 1677 1678 self.ca_distance = distance3d(self.ca_prox, self.ca_dist) 1679 self.cb_distance = distance3d(self.cb_prox, self.cb_dist) 1680 self.torsion_array = np.array( 1681 (self.chi1, self.chi2, self.chi3, self.chi4, self.chi5) 1682 ) 1683 self.compute_torsion_length() 1684 1685 # calculate and set the SS bond torsional energy 1686 self.compute_torsional_energy() 1687 1688 # compute and set the local coordinates 1689 self.compute_local_coords()
Initialize a new Disulfide object with atomic coordinates from the proximal and distal coordinates, typically taken from a PDB file. This routine is primarily used internally when building the compressed database.
Parameters
- chain1: list of Residues in the model, eg: chain = model['A']
- chain2: list of Residues in the model, eg: chain = model['A']
- proximal: proximal residue sequence ID
- distal: distal residue sequence ID
- resolution: structure resolution
- quiet: Quiet or noisy parsing, defaults to True
Raises
- DisulfideConstructionWarning: Raised when not parsed correctly
def
internal_coords(self) -> <built-in function array>:
1691 def internal_coords(self) -> np.array: 1692 """ 1693 Return the internal coordinates for the Disulfide. 1694 If there are missing atoms the extra atoms for the proximal 1695 and distal N and C are set to [0,0,0]. This is needed for the center of 1696 mass calculations, used when rendering. 1697 1698 :return: Array containing the coordinates, [16][3]. 1699 """ 1700 1701 # if we don't have the prior and next atoms we initialize those 1702 # atoms to the origin so as to not effect the center of mass calculations 1703 if self.missing_atoms: 1704 res_array = np.array( 1705 ( 1706 self._n_prox.get_array(), 1707 self._ca_prox.get_array(), 1708 self._c_prox.get_array(), 1709 self._o_prox.get_array(), 1710 self._cb_prox.get_array(), 1711 self._sg_prox.get_array(), 1712 self._n_dist.get_array(), 1713 self._ca_dist.get_array(), 1714 self._c_dist.get_array(), 1715 self._o_dist.get_array(), 1716 self._cb_dist.get_array(), 1717 self._sg_dist.get_array(), 1718 [0, 0, 0], 1719 [0, 0, 0], 1720 [0, 0, 0], 1721 [0, 0, 0], 1722 ) 1723 ) 1724 else: 1725 res_array = np.array( 1726 ( 1727 self._n_prox.get_array(), 1728 self._ca_prox.get_array(), 1729 self._c_prox.get_array(), 1730 self._o_prox.get_array(), 1731 self._cb_prox.get_array(), 1732 self._sg_prox.get_array(), 1733 self._n_dist.get_array(), 1734 self._ca_dist.get_array(), 1735 self._c_dist.get_array(), 1736 self._o_dist.get_array(), 1737 self._cb_dist.get_array(), 1738 self._sg_dist.get_array(), 1739 self._c_prev_prox.get_array(), 1740 self._n_next_prox.get_array(), 1741 self._c_prev_dist.get_array(), 1742 self._n_next_dist.get_array(), 1743 ) 1744 ) 1745 return res_array
Return the internal coordinates for the Disulfide. If there are missing atoms the extra atoms for the proximal and distal N and C are set to [0,0,0]. This is needed for the center of mass calculations, used when rendering.
Returns
Array containing the coordinates, [16][3].
def
internal_coords_res(self, resnumb) -> <built-in function array>:
1747 def internal_coords_res(self, resnumb) -> np.array: 1748 """ 1749 Return the internal coordinates for the internal coordinates of 1750 the given Disulfide. Missing atoms are not included. 1751 1752 :param resnumb: Residue number for disulfide 1753 :raises DisulfideConstructionWarning: Warning raised if the residue number is invalid 1754 :return: Array containing the internal coordinates for the disulfide 1755 """ 1756 res_array = np.zeros(shape=(6, 3)) 1757 1758 if resnumb == self.proximal: 1759 res_array = np.array( 1760 ( 1761 self._n_prox.get_array(), 1762 self._ca_prox.get_array(), 1763 self._c_prox.get_array(), 1764 self._o_prox.get_array(), 1765 self._cb_prox.get_array(), 1766 self._sg_prox.get_array(), 1767 ) 1768 ) 1769 return res_array 1770 1771 elif resnumb == self.distal: 1772 res_array = np.array( 1773 ( 1774 self._n_dist.get_array(), 1775 self._ca_dist.get_array(), 1776 self._c_dist.get_array(), 1777 self._o_dist.get_array(), 1778 self._cb_dist.get_array(), 1779 self._sg_dist.get_array(), 1780 ) 1781 ) 1782 return res_array 1783 else: 1784 mess = f"-> Disulfide.internal_coords(): Invalid argument. \ 1785 Unable to find residue: {resnumb} " 1786 raise DisulfideConstructionWarning(mess)
Return the internal coordinates for the internal coordinates of the given Disulfide. Missing atoms are not included.
Parameters
- resnumb: Residue number for disulfide
Raises
- DisulfideConstructionWarning: Warning raised if the residue number is invalid
Returns
Array containing the internal coordinates for the disulfide
def
make_movie(self, style='sb', fname='ssbond.mp4', verbose=False, steps=360) -> None:
1788 def make_movie( 1789 self, style="sb", fname="ssbond.mp4", verbose=False, steps=360 1790 ) -> None: 1791 """ 1792 Create an animation for ```self``` rotating one revolution about the Y axis, 1793 in the given ```style```, saving to ```filename```. 1794 1795 :param style: Rendering style, defaults to 'sb', one of: 1796 * 'sb' - split bonds 1797 * 'bs' - ball and stick 1798 * 'cpk' - CPK style 1799 * 'pd' - Proximal/Distal style - Red=proximal, Green=Distal 1800 * 'plain' - boring single color 1801 1802 :param fname: Output filename, defaults to ```ssbond.mp4``` 1803 :param verbose: Verbosity, defaults to False 1804 :param steps: Number of steps for one complete rotation, defaults to 360. 1805 """ 1806 src = self.pdb_id 1807 name = self.name 1808 enrg = self.energy 1809 1810 title = f"{src} {name}: {self.proximal}{self.proximal_chain}-{self.distal}{self.distal_chain}: {enrg:.2f} kcal/mol, Cα: {self.ca_distance:.2f} Å, Tors: {self.torsion_length:.2f}" 1811 1812 if verbose: 1813 print(f"Rendering animation to {fname}...") 1814 1815 pl = pv.Plotter(window_size=WINSIZE, off_screen=True) 1816 pl.open_movie(fname) 1817 path = pl.generate_orbital_path(n_points=steps) 1818 1819 # 1820 # pl.add_title(title=title, font_size=FONTSIZE) 1821 pl.enable_anti_aliasing("msaa") 1822 # pl.add_camera_orientation_widget() 1823 pl = self._render( 1824 pl, 1825 style=style, 1826 bondcolor=BOND_COLOR, 1827 bs_scale=BS_SCALE, 1828 spec=SPECULARITY, 1829 specpow=SPEC_POWER, 1830 ) 1831 pl.reset_camera() 1832 pl.orbit_on_path(path, write_frames=True) 1833 pl.close() 1834 1835 if verbose: 1836 print(f"Saved mp4 animation to: {fname}")
Create an animation for self rotating one revolution about the Y axis,
in the given style, saving to filename.
Parameters
style: Rendering style, defaults to 'sb', one of:
- 'sb' - split bonds
- 'bs' - ball and stick
- 'cpk' - CPK style
- 'pd' - Proximal/Distal style - Red=proximal, Green=Distal
- 'plain' - boring single color
fname: Output filename, defaults to
ssbond.mp4- verbose: Verbosity, defaults to False
- steps: Number of steps for one complete rotation, defaults to 360.
def
pprint(self) -> None:
1838 def pprint(self) -> None: 1839 """ 1840 Pretty print general info for the Disulfide 1841 """ 1842 s1 = self.repr_ss_info() 1843 s2 = self.repr_ss_ca_dist() 1844 s3 = self.repr_ss_conformation() 1845 s4 = self.repr_ss_torsion_length() 1846 res = f"{s1} \n{s3} \n{s2} \n{s4}>" 1847 print(res)
Pretty print general info for the Disulfide
def
pprint_all(self) -> None:
1849 def pprint_all(self) -> None: 1850 """ 1851 Pretty print all info for the Disulfide 1852 """ 1853 s1 = self.repr_ss_info() + "\n" 1854 s2 = self.repr_ss_coords() 1855 s3 = self.repr_ss_local_coords() 1856 s4 = self.repr_ss_conformation() 1857 s5 = self.repr_chain_ids() 1858 s6 = self.repr_ss_ca_dist() 1859 s7 = self.repr_ss_torsion_length() 1860 1861 res = f"{s1} {s5} {s2} {s3} {s4}\n {s6}\n {s7}>" 1862 1863 print(res)
Pretty print all info for the Disulfide
def
repr_ss_info(self) -> str:
1866 def repr_ss_info(self) -> str: 1867 """ 1868 Representation for the Disulfide class 1869 """ 1870 s1 = f"<Disulfide {self.name}, Source: {self.pdb_id}, Resolution: {self.resolution} Å" 1871 return s1
Representation for the Disulfide class
def
repr_ss_coords(self) -> str:
1873 def repr_ss_coords(self) -> str: 1874 """ 1875 Representation for Disulfide coordinates 1876 """ 1877 s2 = f"\nProximal Coordinates:\n N: {self.n_prox}\n Cα: {self.ca_prox}\n C: {self.c_prox}\n O: {self.o_prox}\n Cβ: {self.cb_prox}\n Sγ: {self.sg_prox}\n Cprev {self.c_prev_prox}\n Nnext: {self.n_next_prox}\n" 1878 s3 = f"Distal Coordinates:\n N: {self.n_dist}\n Cα: {self.ca_dist}\n C: {self.c_dist}\n O: {self.o_dist}\n Cβ: {self.cb_dist}\n Sγ: {self.sg_dist}\n Cprev {self.c_prev_dist}\n Nnext: {self.n_next_dist}\n\n" 1879 stot = f"{s2} {s3}" 1880 return stot
Representation for Disulfide coordinates
def
repr_ss_conformation(self) -> str:
1882 def repr_ss_conformation(self) -> str: 1883 """ 1884 Representation for Disulfide conformation 1885 """ 1886 s4 = f"Χ1-Χ5: {self.chi1:.2f}°, {self.chi2:.2f}°, {self.chi3:.2f}°, {self.chi4:.2f}° {self.chi5:.2f}°, {self.rho:.2f}°, {self.energy:.2f} kcal/mol" 1887 stot = f"{s4}" 1888 return stot
Representation for Disulfide conformation
def
repr_ss_local_coords(self) -> str:
1890 def repr_ss_local_coords(self) -> str: 1891 """ 1892 Representation for the Disulfide internal coordinates. 1893 """ 1894 s2i = f"Proximal Internal Coords:\n N: {self._n_prox}\n Cα: {self._ca_prox}\n C: {self._c_prox}\n O: {self._o_prox}\n Cβ: {self._cb_prox}\n Sγ: {self._sg_prox}\n Cprev {self.c_prev_prox}\n Nnext: {self.n_next_prox}\n" 1895 s3i = f"Distal Internal Coords:\n N: {self._n_dist}\n Cα: {self._ca_dist}\n C: {self._c_dist}\n O: {self._o_dist}\n Cβ: {self._cb_dist}\n Sγ: {self._sg_dist}\n Cprev {self.c_prev_dist}\n Nnext: {self.n_next_dist}\n" 1896 stot = f"{s2i}{s3i}" 1897 return stot
Representation for the Disulfide internal coordinates.
def
repr_ss_chain_ids(self) -> str:
1899 def repr_ss_chain_ids(self) -> str: 1900 """ 1901 Representation for Disulfide chain IDs 1902 """ 1903 return f"Proximal Chain fullID: <{self.proximal_residue_fullid}> Distal Chain fullID: <{self.distal_residue_fullid}>"
Representation for Disulfide chain IDs
def
repr_ss_ca_dist(self) -> str:
1905 def repr_ss_ca_dist(self) -> str: 1906 """ 1907 Representation for Disulfide Ca distance 1908 """ 1909 s1 = f"Cα Distance: {self.ca_distance:.2f} Å" 1910 return s1
Representation for Disulfide Ca distance
def
repr_ss_cb_dist(self) -> str:
1912 def repr_ss_cb_dist(self) -> str: 1913 """ 1914 Representation for Disulfide Ca distance 1915 """ 1916 s1 = f"Cβ Distance: {self.cb_distance:.2f} Å" 1917 return s1
Representation for Disulfide Ca distance
def
repr_ss_torsion_length(self) -> str:
1919 def repr_ss_torsion_length(self) -> str: 1920 """ 1921 Representation for Disulfide torsion length 1922 """ 1923 s1 = f"Torsion length: {self.torsion_length:.2f} deg" 1924 return s1
Representation for Disulfide torsion length
def
repr_all(self) -> str:
1926 def repr_all(self) -> str: 1927 """ 1928 Return a string representation for all Disulfide information 1929 contained in self. 1930 """ 1931 1932 s1 = self.repr_ss_info() + "\n" 1933 s2 = self.repr_ss_coords() 1934 s3 = self.repr_ss_local_coords() 1935 s4 = self.repr_ss_conformation() 1936 s5 = self.repr_chain_ids() 1937 s6 = self.repr_ss_ca_dist() 1938 s8 = self.repr_ss_cb_dist() 1939 s7 = self.repr_ss_torsion_length() 1940 1941 res = f"{s1} {s5} {s2} {s3} {s4} {s6} {s7} {s8}>" 1942 return res
Return a string representation for all Disulfide information contained in self.
def
repr_compact(self) -> str:
1944 def repr_compact(self) -> str: 1945 """ 1946 Return a compact representation of the Disulfide object 1947 :return: string 1948 """ 1949 return f"{self.repr_ss_info()} {self.repr_ss_conformation()}"
Return a compact representation of the Disulfide object
Returns
string
def
repr_conformation(self) -> str:
1951 def repr_conformation(self) -> str: 1952 """ 1953 Return a string representation of the Disulfide object's conformation. 1954 :return: string 1955 """ 1956 return f"{self.repr_ss_conformation()}"
Return a string representation of the Disulfide object's conformation.
Returns
string
def
repr_coords(self) -> str:
1958 def repr_coords(self) -> str: 1959 """ 1960 Return a string representation of the Disulfide object's coordinates. 1961 :return: string 1962 """ 1963 return f"{self.repr_ss_coords()}"
Return a string representation of the Disulfide object's coordinates.
Returns
string
def
repr_internal_coords(self) -> str:
1965 def repr_internal_coords(self) -> str: 1966 """ 1967 Return a string representation of the Disulfide object's internal coordinaes. 1968 :return: string 1969 """ 1970 return f"{self.repr_ss_local_coords()}"
Return a string representation of the Disulfide object's internal coordinaes.
Returns
string
def
repr_chain_ids(self) -> str:
1972 def repr_chain_ids(self) -> str: 1973 """ 1974 Return a string representation of the Disulfide object's chain ids. 1975 :return: string 1976 """ 1977 return f"{self.repr_ss_chain_ids()}"
Return a string representation of the Disulfide object's chain ids.
Returns
string
def
reset(self) -> None:
1983 def reset(self) -> None: 1984 """ 1985 Resets the disulfide object to its initial state. All distances, 1986 angles and positions are reset. The name is unchanged. 1987 """ 1988 self.__init__(self)
Resets the disulfide object to its initial state. All distances, angles and positions are reset. The name is unchanged.
def
same_chains(self) -> bool:
1990 def same_chains(self) -> bool: 1991 """ 1992 Function checks if the Disulfide is cross-chain or not. 1993 1994 Returns 1995 ------- 1996 bool \n 1997 True if the proximal and distal residues are on the same chains, 1998 False otherwise. 1999 """ 2000 2001 (prox, dist) = self.get_chains() 2002 return prox == dist
Function checks if the Disulfide is cross-chain or not.
Returns
bool
True if the proximal and distal residues are on the same chains,
False otherwise.
def
screenshot( self, single=True, style='sb', fname='ssbond.png', verbose=False, shadows=False, light=True) -> None:
2004 def screenshot( 2005 self, 2006 single=True, 2007 style="sb", 2008 fname="ssbond.png", 2009 verbose=False, 2010 shadows=False, 2011 light=True, 2012 ) -> None: 2013 """ 2014 Create and save a screenshot of the Disulfide in the given style 2015 and filename 2016 2017 :param single: Display a single vs panel view, defaults to True 2018 :param style: Rendering style, one of: 2019 * 'sb' - split bonds 2020 * 'bs' - ball and stick 2021 * 'cpk' - CPK style 2022 * 'pd' - Proximal/Distal style - Red=proximal, Green=Distal 2023 * 'plain' - boring single color, 2024 :param fname: output filename,, defaults to 'ssbond.png' 2025 :param verbose: Verbosit, defaults to False 2026 :param shadows: Enable shadows, defaults to False 2027 """ 2028 src = self.pdb_id 2029 name = self.name 2030 enrg = self.energy 2031 2032 title = f"{src} {name}: {self.proximal}{self.proximal_chain}-{self.distal}{self.distal_chain}: {enrg:.2f} kcal/mol, Cα: {self.ca_distance:.2f} Å, Cβ: {self.cb_distance:.2f} Å, Tors: {self.torsion_length:.2f}" 2033 2034 if light: 2035 pv.set_plot_theme("document") 2036 else: 2037 pv.set_plot_theme("dark") 2038 2039 if verbose: 2040 print(f"-> screenshot(): Rendering screenshot to file {fname}") 2041 2042 if single: 2043 pl = pv.Plotter(window_size=WINSIZE) 2044 # pl.add_title(title=title, font_size=FONTSIZE) 2045 pl.enable_anti_aliasing("msaa") 2046 # pl.add_camera_orientation_widget() 2047 self._render( 2048 pl, 2049 style=style, 2050 bondcolor=BOND_COLOR, 2051 bs_scale=BS_SCALE, 2052 spec=SPECULARITY, 2053 specpow=SPEC_POWER, 2054 ) 2055 pl.reset_camera() 2056 if shadows: 2057 pl.enable_shadows() 2058 2059 pl.show(auto_close=False) 2060 pl.screenshot(fname) 2061 pl.clear() 2062 2063 else: 2064 pl = pv.Plotter(window_size=WINSIZE, shape=(2, 2)) 2065 pl.subplot(0, 0) 2066 2067 # pl.add_title(title=title, font_size=FONTSIZE) 2068 pl.enable_anti_aliasing("msaa") 2069 2070 # pl.add_camera_orientation_widget() 2071 self._render( 2072 pl, 2073 style="cpk", 2074 bondcolor=BOND_COLOR, 2075 bs_scale=BS_SCALE, 2076 spec=SPECULARITY, 2077 specpow=SPEC_POWER, 2078 ) 2079 2080 pl.subplot(0, 1) 2081 # pl.add_title(title=title, font_size=FONTSIZE) 2082 self._render( 2083 pl, 2084 style="pd", 2085 bondcolor=BOND_COLOR, 2086 bs_scale=BS_SCALE, 2087 spec=SPECULARITY, 2088 specpow=SPEC_POWER, 2089 ) 2090 2091 pl.subplot(1, 0) 2092 # pl.add_title(title=title, font_size=FONTSIZE) 2093 self._render( 2094 pl, 2095 style="bs", 2096 bondcolor=BOND_COLOR, 2097 bs_scale=BS_SCALE, 2098 spec=SPECULARITY, 2099 specpow=SPEC_POWER, 2100 ) 2101 2102 pl.subplot(1, 1) 2103 # pl.add_title(title=title, font_size=FONTSIZE) 2104 self._render( 2105 pl, 2106 style="sb", 2107 bondcolor=BOND_COLOR, 2108 bs_scale=BS_SCALE, 2109 spec=SPECULARITY, 2110 specpow=SPEC_POWER, 2111 ) 2112 2113 pl.link_views() 2114 pl.reset_camera() 2115 if shadows: 2116 pl.enable_shadows() 2117 2118 pl.show(auto_close=False) 2119 pl.screenshot(fname) 2120 2121 if verbose: 2122 print(f"Saved: {fname}")
Create and save a screenshot of the Disulfide in the given style and filename
Parameters
- single: Display a single vs panel view, defaults to True
- style: Rendering style, one of:
- 'sb' - split bonds
- 'bs' - ball and stick
- 'cpk' - CPK style
- 'pd' - Proximal/Distal style - Red=proximal, Green=Distal
- 'plain' - boring single color,
- fname: output filename,, defaults to 'ssbond.png'
- verbose: Verbosit, defaults to False
- shadows: Enable shadows, defaults to False
def
save_meshes_as_stl(self, meshes, filename) -> None:
2124 def save_meshes_as_stl(self, meshes, filename) -> None: 2125 """Save a list of meshes as a single STL file. 2126 2127 Args: 2128 meshes (list): List of pyvista mesh objects to save. 2129 filename (str): Path to save the STL file to. 2130 """ 2131 merged_mesh = pv.UnstructuredGrid() 2132 for mesh in meshes: 2133 merged_mesh += mesh 2134 merged_mesh.save(filename)
Save a list of meshes as a single STL file.
Args: meshes (list): List of pyvista mesh objects to save. filename (str): Path to save the STL file to.
def
export(self, style='sb', verbose=True, fname='ssbond_plt') -> None:
2136 def export(self, style="sb", verbose=True, fname="ssbond_plt") -> None: 2137 """ 2138 Create and save a screenshot of the Disulfide in the given style and filename. 2139 2140 :param single: Display a single vs panel view, defaults to True 2141 :param style: Rendering style, one of: 2142 * 'sb' - split bonds 2143 * 'bs' - ball and stick 2144 * 'cpk' - CPK style 2145 * 'pd' - Proximal/Distal style - Red=proximal, Green=Distal 2146 * 'plain' - boring single color, 2147 2148 :param fname: output filename,, defaults to 'ssbond.stl' 2149 :param verbose: Verbosit, defaults to False 2150 """ 2151 2152 if verbose: 2153 print(f"-> screenshot(): Rendering screenshot to file {fname}") 2154 2155 pl = pv.PolyData() 2156 2157 self._plot( 2158 pl, 2159 style=style, 2160 bondcolor=BOND_COLOR, 2161 bs_scale=BS_SCALE, 2162 spec=SPECULARITY, 2163 specpow=SPEC_POWER, 2164 ) 2165 2166 self.save_meshes_as_stl(pl, fname) 2167 2168 return
Create and save a screenshot of the Disulfide in the given style and filename.
Parameters
- single: Display a single vs panel view, defaults to True
style: Rendering style, one of:
- 'sb' - split bonds
- 'bs' - ball and stick
- 'cpk' - CPK style
- 'pd' - Proximal/Distal style - Red=proximal, Green=Distal
- 'plain' - boring single color,
fname: output filename,, defaults to 'ssbond.stl'
- verbose: Verbosit, defaults to False
def
set_permissive(self, perm: bool) -> None:
2170 def set_permissive(self, perm: bool) -> None: 2171 """ 2172 Set PERMISSIVE flag for Disulfide parsing. 2173 2174 :return: None 2175 """ 2176 2177 self.PERMISSIVE = perm
Set PERMISSIVE flag for Disulfide parsing.
Returns
None
def
set_positions( self, n_prox: Bio.PDB.vectors.Vector, ca_prox: Bio.PDB.vectors.Vector, c_prox: Bio.PDB.vectors.Vector, o_prox: Bio.PDB.vectors.Vector, cb_prox: Bio.PDB.vectors.Vector, sg_prox: Bio.PDB.vectors.Vector, n_dist: Bio.PDB.vectors.Vector, ca_dist: Bio.PDB.vectors.Vector, c_dist: Bio.PDB.vectors.Vector, o_dist: Bio.PDB.vectors.Vector, cb_dist: Bio.PDB.vectors.Vector, sg_dist: Bio.PDB.vectors.Vector, c_prev_prox: Bio.PDB.vectors.Vector, n_next_prox: Bio.PDB.vectors.Vector, c_prev_dist: Bio.PDB.vectors.Vector, n_next_dist: Bio.PDB.vectors.Vector) -> None:
2179 def set_positions( 2180 self, 2181 n_prox: Vector, 2182 ca_prox: Vector, 2183 c_prox: Vector, 2184 o_prox: Vector, 2185 cb_prox: Vector, 2186 sg_prox: Vector, 2187 n_dist: Vector, 2188 ca_dist: Vector, 2189 c_dist: Vector, 2190 o_dist: Vector, 2191 cb_dist: Vector, 2192 sg_dist: Vector, 2193 c_prev_prox: Vector, 2194 n_next_prox: Vector, 2195 c_prev_dist: Vector, 2196 n_next_dist: Vector, 2197 ) -> None: 2198 """ 2199 Set the atomic coordinates for all atoms in the Disulfide object. 2200 2201 :param n_prox: Proximal N position 2202 :param ca_prox: Proximal Cα position 2203 :param c_prox: Proximal C' position 2204 :param o_prox: Proximal O position 2205 :param cb_prox: Proximal Cβ position 2206 :param sg_prox: Proximal Sγ position 2207 :param n_dist: Distal N position 2208 :param ca_dist: Distal Cα position 2209 :param c_dist: Distal C' position 2210 :param o_dist: Distal O position 2211 :param cb_dist: Distal Cβ position 2212 :param sg_dist: Distal Sγ position 2213 :param c_prev_prox: Proximal previous C' 2214 :param n_next_prox: Proximal next N 2215 :param c_prev_dist: Distal previous C' 2216 :param n_next_dist: Distal next N 2217 """ 2218 2219 # deep copy 2220 self.n_prox = n_prox.copy() 2221 self.ca_prox = ca_prox.copy() 2222 self.c_prox = c_prox.copy() 2223 self.o_prox = o_prox.copy() 2224 self.cb_prox = cb_prox.copy() 2225 self.sg_prox = sg_prox.copy() 2226 self.sg_dist = sg_dist.copy() 2227 self.cb_dist = cb_dist.copy() 2228 self.ca_dist = ca_dist.copy() 2229 self.n_dist = n_dist.copy() 2230 self.c_dist = c_dist.copy() 2231 self.o_dist = o_dist.copy() 2232 2233 self.c_prev_prox = c_prev_prox.copy() 2234 self.n_next_prox = n_next_prox.copy() 2235 self.c_prev_dist = c_prev_dist.copy() 2236 self.n_next_dist = n_next_dist.copy()
Set the atomic coordinates for all atoms in the Disulfide object.
Parameters
- n_prox: Proximal N position
- ca_prox: Proximal Cα position
- c_prox: Proximal C' position
- o_prox: Proximal O position
- cb_prox: Proximal Cβ position
- sg_prox: Proximal Sγ position
- n_dist: Distal N position
- ca_dist: Distal Cα position
- c_dist: Distal C' position
- o_dist: Distal O position
- cb_dist: Distal Cβ position
- sg_dist: Distal Sγ position
- c_prev_prox: Proximal previous C'
- n_next_prox: Proximal next N
- c_prev_dist: Distal previous C'
- n_next_dist: Distal next N
def
set_dihedrals( self, chi1: float, chi2: float, chi3: float, chi4: float, chi5: float) -> None:
2238 def set_dihedrals( 2239 self, chi1: float, chi2: float, chi3: float, chi4: float, chi5: float 2240 ) -> None: 2241 """ 2242 Set the disulfide's dihedral angles, Chi1-Chi5. -180 - 180 degrees. 2243 2244 :param chi1: Chi1 2245 :param chi2: Chi2 2246 :param chi3: Chi3 2247 :param chi4: Chi4 2248 :param chi5: Chi5 2249 """ 2250 self.chi1 = chi1 2251 self.chi2 = chi2 2252 self.chi3 = chi3 2253 self.chi4 = chi4 2254 self.chi5 = chi5 2255 self.torsion_array = np.array([chi1, chi2, chi3, chi4, chi5]) 2256 self.compute_torsional_energy() 2257 self.compute_torsion_length()
Set the disulfide's dihedral angles, Chi1-Chi5. -180 - 180 degrees.
Parameters
- chi1: Chi1
- chi2: Chi2
- chi3: Chi3
- chi4: Chi4
- chi5: Chi5
def
set_name(self, namestr='Disulfide') -> None:
2259 def set_name(self, namestr="Disulfide") -> None: 2260 """ 2261 Set the Disulfide's name. 2262 2263 :param namestr: Name, by default "Disulfide" 2264 """ 2265 self.name = namestr
Set the Disulfide's name.
Parameters
- namestr: Name, by default "Disulfide"
def
set_resnum(self, proximal: int, distal: int) -> None:
2267 def set_resnum(self, proximal: int, distal: int) -> None: 2268 """ 2269 Set the proximal and residue numbers for the Disulfide. 2270 2271 :param proximal: Proximal residue number 2272 :param distal: Distal residue number 2273 """ 2274 self.proximal = proximal 2275 self.distal = distal
Set the proximal and residue numbers for the Disulfide.
Parameters
- proximal: Proximal residue number
- distal: Distal residue number
def
compute_torsion_length(self) -> float:
2277 def compute_torsion_length(self) -> float: 2278 """ 2279 Compute the 5D Euclidean length of the Disulfide object. Update the disulfide internal state. 2280 2281 :return: Torsion length (Degrees) 2282 """ 2283 # Use numpy array to compute element-wise square 2284 tors2 = np.square(self.torsion_array) 2285 2286 # Compute the sum of squares using numpy's sum function 2287 dist = math.sqrt(np.sum(tors2)) 2288 2289 # Update the internal state 2290 self.torsion_length = dist 2291 2292 return dist
Compute the 5D Euclidean length of the Disulfide object. Update the disulfide internal state.
Returns
Torsion length (Degrees)
def
Torsion_Distance(self, other) -> float:
2294 def Torsion_Distance(self, other) -> float: 2295 """ 2296 Calculate the 5D Euclidean distance between ```self``` and another Disulfide 2297 object. This is used to compare Disulfide Bond torsion angles to 2298 determine their torsional similarity via a 5-Dimensional Euclidean distance metric. 2299 2300 :param other: Comparison Disulfide 2301 :raises ProteusPyWarning: Warning if ```other``` is not a Disulfide object 2302 :return: Euclidean distance (Degrees) between ```self``` and ```other```. 2303 """ 2304 2305 from proteusPy.ProteusPyWarning import ProteusPyWarning 2306 2307 # Check length of torsion arrays 2308 if len(self.torsion_array) != 5 or len(other.torsion_array) != 5: 2309 raise ProteusPyWarning( 2310 "--> Torsion_Distance() requires vectors of length 5!" 2311 ) 2312 2313 # Convert to numpy arrays and add 180 to each element 2314 p1 = np.array(self.torsion_array) + 180.0 2315 p2 = np.array(other.torsion_array) + 180.0 2316 2317 # Compute the 5D Euclidean distance using numpy's linalg.norm function 2318 dist = np.linalg.norm(p1 - p2) 2319 2320 return dist
Calculate the 5D Euclidean distance between self and another Disulfide
object. This is used to compare Disulfide Bond torsion angles to
determine their torsional similarity via a 5-Dimensional Euclidean distance metric.
Parameters
- other: Comparison Disulfide
Raises
- ProteusPyWarning: Warning if
otheris not a Disulfide object
Returns
Euclidean distance (Degrees) between
selfandother.
2322 def Torsion_neighbors(self, others, cutoff) -> DisulfideList: 2323 """ 2324 Return a list of Disulfides within the angular cutoff in the others list. 2325 This routine is used to find Disulfides having the same torsion length 2326 within the others list. This is used to find families of Disulfides with 2327 similar conformations. Assumes self is properly initialized. 2328 2329 *NB* The routine will not distinguish between +/- 2330 dihedral angles. *i.e.* [-60, -60, -90, -60, -60] would have the same 2331 torsion length as [60, 60, 90, 60, 60], two clearly different structures. 2332 2333 :param others: ```DisulfideList``` to search 2334 :param cutoff: Dihedral angle degree cutoff 2335 :return: DisulfideList within the cutoff 2336 2337 Example: 2338 In this example we load the disulfide database subset, find the disulfides with 2339 the lowest and highest energies, and then find the nearest conformational neighbors. 2340 Finally, we display the neighbors overlaid against a common reference frame. 2341 2342 >>> from proteusPy import * 2343 >>> from proteusPy.DisulfideLoader import Load_PDB_SS 2344 >>> from proteusPy.DisulfideList import DisulfideList 2345 >>> from proteusPy.Disulfide import Disulfide 2346 >>> PDB_SS = Load_PDB_SS(verbose=False, subset=True) 2347 >>> ss_list = DisulfideList([], 'tmp') 2348 2349 We point to the complete list to search for lowest and highest energies. 2350 >>> sslist = PDB_SS.SSList 2351 >>> ssmin_enrg, ssmax_enrg = PDB_SS.SSList.minmax_energy 2352 2353 Make an empty list and find the nearest neighbors within 10 degrees avg RMS in 2354 sidechain dihedral angle space. 2355 2356 >>> low_energy_neighbors = DisulfideList([],'Neighbors') 2357 >>> low_energy_neighbors = ssmin_enrg.Torsion_neighbors(sslist, 10) 2358 2359 Display the number found, and then display them overlaid onto their common reference frame. 2360 2361 >>> tot = low_energy_neighbors.length 2362 >>> print(f'Neighbors: {tot}') 2363 Neighbors: 2 2364 >>> low_energy_neighbors.display_overlay() 2365 2366 """ 2367 res = [ss for ss in others if self.Torsion_Distance(ss) <= cutoff] 2368 return DisulfideList(res, "neighbors")
Return a list of Disulfides within the angular cutoff in the others list. This routine is used to find Disulfides having the same torsion length within the others list. This is used to find families of Disulfides with similar conformations. Assumes self is properly initialized.
NB The routine will not distinguish between +/- dihedral angles. i.e. [-60, -60, -90, -60, -60] would have the same torsion length as [60, 60, 90, 60, 60], two clearly different structures.
Parameters
- others:
DisulfideListto search - cutoff: Dihedral angle degree cutoff
Returns
DisulfideList within the cutoff
Example: In this example we load the disulfide database subset, find the disulfides with the lowest and highest energies, and then find the nearest conformational neighbors. Finally, we display the neighbors overlaid against a common reference frame.
>>> from proteusPy import *
>>> from proteusPy.DisulfideLoader import Load_PDB_SS
>>> from proteusPy.DisulfideList import DisulfideList
>>> from proteusPy.Disulfide import Disulfide
>>> PDB_SS = Load_PDB_SS(verbose=False, subset=True)
>>> ss_list = DisulfideList([], 'tmp')
We point to the complete list to search for lowest and highest energies.
>>> sslist = PDB_SS.SSList
>>> ssmin_enrg, ssmax_enrg = PDB_SS.SSList.minmax_energy
Make an empty list and find the nearest neighbors within 10 degrees avg RMS in sidechain dihedral angle space.
>>> low_energy_neighbors = DisulfideList([],'Neighbors')
>>> low_energy_neighbors = ssmin_enrg.Torsion_neighbors(sslist, 10)
Display the number found, and then display them overlaid onto their common reference frame.
>>> tot = low_energy_neighbors.length
>>> print(f'Neighbors: {tot}')
Neighbors: 2
>>> low_energy_neighbors.display_overlay()
def
torsion_to_sixclass(self) -> str:
2370 def torsion_to_sixclass(self) -> str: 2371 """ 2372 Return the sextant class string for ``self``. 2373 2374 :raises DisulfideIOException: _description_ 2375 :return: Sextant string 2376 """ 2377 from proteusPy.DisulfideClasses import get_sixth_quadrant 2378 2379 tors = self.torsion_array 2380 res = [get_sixth_quadrant(x) for x in tors] 2381 return "".join([str(r) for r in res])
Return the sextant class string for self.
Raises
- DisulfideIOException: _description_
Returns
Sextant string
def
parse_ssbond_header_rec(ssbond_dict: dict) -> list:
2387def parse_ssbond_header_rec(ssbond_dict: dict) -> list: 2388 """ 2389 Parse the SSBOND dict returned by parse_pdb_header. 2390 NB: Requires EGS-Modified BIO.parse_pdb_header.py. 2391 This is used internally. 2392 2393 :param ssbond_dict: the input SSBOND dict 2394 :return: A list of tuples representing the proximal, 2395 distal residue ids for the Disulfide. 2396 2397 """ 2398 disulfide_list = [] 2399 for ssb in ssbond_dict.items(): 2400 disulfide_list.append(ssb[1]) 2401 2402 return disulfide_list
Parse the SSBOND dict returned by parse_pdb_header. NB: Requires EGS-Modified BIO.parse_pdb_header.py. This is used internally.
Parameters
- ssbond_dict: the input SSBOND dict
Returns
A list of tuples representing the proximal, distal residue ids for the Disulfide.
def
Download_Disulfides( pdb_home='/Users/egs/PDB/', model_home='/Users/egs/PDB/good/', verbose=False, reset=False) -> None:
2413def Download_Disulfides( 2414 pdb_home=PDB_DIR, model_home=MODEL_DIR, verbose=False, reset=False 2415) -> None: 2416 """ 2417 Read a comma separated list of PDB IDs and download them 2418 to the pdb_home path. 2419 2420 This utility function is used to download proteins containing at least 2421 one SS bond with the ID list generated from: http://www.rcsb.org/. 2422 2423 This is the primary data loader for the proteusPy Disulfide 2424 analysis package. The list of IDs represents files in the 2425 RCSB containing > 1 disulfide bond, and it contains 2426 over 39000 structures. The total download takes about 12 hours. The 2427 function keeps track of downloaded files so it's possible to interrupt and 2428 restart the download without duplicating effort. 2429 2430 :param pdb_home: Path for downloaded files, defaults to PDB_DIR 2431 :param model_home: Path for extracted data, defaults to MODEL_DIR 2432 :param verbose: Verbosity, defaults to False 2433 :param reset: Reset the downloaded files index. Used to restart the download. 2434 :raises DisulfideIOException: I/O error raised when the PDB file is not found. 2435 """ 2436 import os 2437 import time 2438 2439 start = time.time() 2440 donelines = [] 2441 SS_done = [] 2442 ssfile = None 2443 2444 cwd = os.getcwd() 2445 os.chdir(pdb_home) 2446 2447 pdblist = PDBList(pdb=pdb_home, verbose=verbose) 2448 ssfilename = f"{model_home}{SS_ID_FILE}" 2449 print(ssfilename) 2450 2451 # list of IDs containing >1 SSBond record 2452 try: 2453 ssfile = open(ssfilename) 2454 Line = ssfile.readlines() 2455 except Exception: 2456 raise DisulfideIOException(f"Cannot open file: {ssfile}") 2457 2458 for line in Line: 2459 entries = line.split(",") 2460 2461 print(f"Found: {len(entries)} entries") 2462 completed = {"xxx"} # set to keep track of downloaded 2463 2464 # file to track already downloaded entries. 2465 if reset is True: 2466 completed_file = open(f"{model_home}ss_completed.txt", "w") 2467 donelines = [] 2468 SS_DONE = [] 2469 else: 2470 completed_file = open(f"{model_home}ss_completed.txt", "w+") 2471 donelines = completed_file.readlines() 2472 2473 if len(donelines) > 0: 2474 for dl in donelines[0]: 2475 # create a list of pdb id already downloaded 2476 SS_done = dl.split(",") 2477 2478 count = len(SS_done) - 1 2479 completed.update(SS_done) # update the completed set with what's downloaded 2480 2481 # Loop over all entries, 2482 pbar = tqdm(entries, ncols=PBAR_COLS) 2483 for entry in pbar: 2484 pbar.set_postfix({"Entry": entry}) 2485 if entry not in completed: 2486 if pdblist.retrieve_pdb_file(entry, file_format="pdb", pdir=pdb_home): 2487 completed.update(entry) 2488 completed_file.write(f"{entry},") 2489 count += 1 2490 2491 completed_file.close() 2492 2493 end = time.time() 2494 elapsed = end - start 2495 2496 print(f"Overall files processed: {count}") 2497 print(f"Complete. Elapsed time: {datetime.timedelta(seconds=elapsed)} (h:m:s)") 2498 os.chdir(cwd) 2499 return
Read a comma separated list of PDB IDs and download them to the pdb_home path.
This utility function is used to download proteins containing at least one SS bond with the ID list generated from: http://www.rcsb.org/.
This is the primary data loader for the proteusPy Disulfide analysis package. The list of IDs represents files in the RCSB containing > 1 disulfide bond, and it contains over 39000 structures. The total download takes about 12 hours. The function keeps track of downloaded files so it's possible to interrupt and restart the download without duplicating effort.
Parameters
- pdb_home: Path for downloaded files, defaults to PDB_DIR
- model_home: Path for extracted data, defaults to MODEL_DIR
- verbose: Verbosity, defaults to False
- reset: Reset the downloaded files index. Used to restart the download.
Raises
- DisulfideIOException: I/O error raised when the PDB file is not found.
def
Extract_Disulfides( numb=-1, verbose=False, quiet=True, pdbdir='/Users/egs/PDB/', datadir='/Users/egs/PDB/good/', picklefile='PDB_all_ss.pkl', torsionfile='PDB_all_ss_torsions.csv', problemfile='PDB_all_SS_problems.csv', dictfile='PDB_all_ss_dict.pkl', dist_cutoff=-1.0) -> None:
2506def Extract_Disulfides( 2507 numb=-1, 2508 verbose=False, 2509 quiet=True, 2510 pdbdir=PDB_DIR, 2511 datadir=MODEL_DIR, 2512 picklefile=SS_PICKLE_FILE, 2513 torsionfile=SS_TORSIONS_FILE, 2514 problemfile=PROBLEM_ID_FILE, 2515 dictfile=SS_DICT_PICKLE_FILE, 2516 dist_cutoff=-1.0, 2517) -> None: 2518 """ 2519 Read the PDB files contained in ``pdbdir`` and create the .pkl files needed for the 2520 proteusPy.DisulfideLoader.DisulfideLoader class. 2521 The ```Disulfide``` objects are contained in a ```DisulfideList``` object and 2522 ```Dict``` within these files. In addition, .csv files containing all of 2523 the torsions for the disulfides and problem IDs are written. The optional 2524 ```dist_cutoff``` allows for removal of Disufides whose Cα-Cα distance is > 2525 than the cutoff value. If it's -1.0 then the function keeps all Disulfides. 2526 2527 :param numb: number of entries to process, defaults to all 2528 :param verbose: more messages 2529 :param quiet: turns off DisulfideConstruction warnings 2530 :param pdbdir: path to PDB files 2531 :param datadir: path to resulting .pkl files 2532 :param picklefile: name of the disulfide .pkl file 2533 :param torsionfile: name of the disulfide torsion file .csv created 2534 :param problemfile: name of the .csv file containing problem ids 2535 :param dictfile: name of the .pkl file 2536 :param dist_cutoff: Ca distance cutoff to reject a Disulfide. 2537 """ 2538 2539 def name_to_id(fname: str) -> str: 2540 """ 2541 Returns the PDB ID from the filename. 2542 2543 :param fname: Complete PDB filename 2544 :return: PDB ID 2545 """ 2546 ent = fname[3:-4] 2547 return ent 2548 2549 import os 2550 2551 entrylist = [] 2552 problem_ids = [] 2553 bad = bad_dist = 0 2554 2555 # we use the specialized list class DisulfideList to contain our disulfides 2556 # we'll use a dict to store DisulfideList objects, indexed by the structure ID 2557 All_ss_dict = {} 2558 All_ss_list = DisulfideList([], "PDB_SS") 2559 All_ss_dict2 = {} # new dict of pointers to indices 2560 2561 start = time.time() 2562 cwd = os.getcwd() 2563 2564 # Build a list of PDB files in PDB_DIR that are readable. These files were downloaded 2565 # via the RCSB web query interface for structures containing >= 1 SS Bond. 2566 2567 os.chdir(pdbdir) 2568 2569 ss_filelist = glob.glob("*.ent") 2570 tot = len(ss_filelist) 2571 2572 if verbose: 2573 print(f"PDB Directory {pdbdir} contains: {tot} files") 2574 2575 # the filenames are in the form pdb{entry}.ent, I loop through them and extract 2576 # the PDB ID, with Disulfide.name_to_id(), then add to entrylist. 2577 2578 for entry in ss_filelist: 2579 entrylist.append(name_to_id(entry)) 2580 2581 # create a dataframe with the following columns for the disulfide conformations 2582 # extracted from the structure 2583 2584 SS_df = pandas.DataFrame(columns=Torsion_DF_Cols) 2585 2586 # define a tqdm progressbar using the fully loaded entrylist list. 2587 # If numb is passed then 2588 # only do the last numb entries. 2589 2590 if numb > 0: 2591 pbar = tqdm(entrylist[:numb], ncols=PBAR_COLS) 2592 else: 2593 pbar = tqdm(entrylist, ncols=PBAR_COLS) 2594 2595 tot = 0 2596 cnt = 0 2597 # loop over ss_filelist, create disulfides and initialize them 2598 for entry in pbar: 2599 pbar.set_postfix( 2600 {"ID": entry, "Bad": bad, "Ca": bad_dist, "Cnt": tot} 2601 ) # update the progress bar 2602 2603 # returns an empty list if none are found. 2604 _sslist = DisulfideList([], entry) 2605 _sslist = proteusPy.DisulfideList.load_disulfides_from_id( 2606 entry, model_numb=0, verbose=verbose, quiet=quiet, pdb_dir=pdbdir 2607 ) 2608 sslist, xchain = prune_extra_ss(_sslist) 2609 2610 if len(sslist) > 0: 2611 sslist2 = [] # list to hold indices for ss_dict2 2612 for ss in sslist: 2613 # Ca distance cutoff 2614 dist = ss.ca_distance 2615 if dist >= dist_cutoff and dist_cutoff != -1.0: 2616 bad_dist += 1 2617 continue 2618 2619 All_ss_list.append(ss) 2620 new_row = [ 2621 ss.pdb_id, 2622 ss.name, 2623 ss.proximal, 2624 ss.distal, 2625 ss.chi1, 2626 ss.chi2, 2627 ss.chi3, 2628 ss.chi4, 2629 ss.chi5, 2630 ss.energy, 2631 ss.ca_distance, 2632 ss.cb_distance, 2633 ss.phiprox, 2634 ss.psiprox, 2635 ss.phidist, 2636 ss.psidist, 2637 ss.torsion_length, 2638 ss.rho, 2639 ] 2640 2641 # add the row to the end of the dataframe 2642 SS_df.loc[len(SS_df.index)] = new_row.copy() # deep copy 2643 sslist2.append(cnt) 2644 cnt += 1 2645 tot += 1 2646 2647 # All_ss_dict[entry] = sslist 2648 # print(f'Entry: {entry}. Dict indices: {sslist2}') 2649 All_ss_dict2[entry] = sslist2 2650 # print(f'{entry} ss dict adding: {sslist2}') 2651 2652 else: 2653 # at this point I really shouldn't have any bad non-parsible file 2654 bad += 1 2655 problem_ids.append(entry) 2656 # os.remove(f'pdb{entry}.ent') 2657 2658 if bad > 0: 2659 prob_cols = ["id"] 2660 problem_df = pandas.DataFrame(columns=prob_cols) 2661 problem_df["id"] = problem_ids 2662 2663 print( 2664 f"-> Extract_Disulfides(): Found and removed: {len(problem_ids)} non-parsable structures." 2665 ) 2666 print( 2667 f"-> Extract_Disulfides(): Saving problem IDs to file: {datadir}{problemfile}" 2668 ) 2669 2670 problem_df.to_csv(f"{datadir}{problemfile}") 2671 else: 2672 if verbose: 2673 print("-> Extract_Disulfides(): No non-parsable structures found.") 2674 2675 if bad_dist > 0: 2676 print(f"-> Extract_Disulfides(): Found and ignored: {bad_dist} long SS bonds.") 2677 else: 2678 if verbose: 2679 print("No problems found.") 2680 2681 # dump the all_ss list of disulfides to a .pkl file. ~520 MB. 2682 fname = f"{datadir}{picklefile}" 2683 print( 2684 f"-> Extract_Disulfides(): Saving {len(All_ss_list)} Disulfides to file: {fname}" 2685 ) 2686 2687 with open(fname, "wb+") as f: 2688 pickle.dump(All_ss_list, f) 2689 2690 # dump the dict2 disulfides to a .pkl file. ~520 MB. 2691 dict_len = len(All_ss_dict2) 2692 fname = f"{datadir}{dictfile}" 2693 print( 2694 f"-> Extract_Disulfides(): Saving indices of {dict_len} Disulfide-containing PDB IDs to file: {fname}" 2695 ) 2696 2697 with open(fname, "wb+") as f: 2698 pickle.dump(All_ss_dict2, f) 2699 2700 # save the torsions 2701 fname = f"{datadir}{torsionfile}" 2702 print(f"-> Extract_Disulfides(): Saving torsions to file: {fname}") 2703 SS_df.to_csv(fname) 2704 2705 end = time.time() 2706 elapsed = end - start 2707 2708 print( 2709 f"-> Extract_Disulfides(): Disulfide Extraction complete! Elapsed time:\ 2710 {datetime.timedelta(seconds=elapsed)} (h:m:s)" 2711 ) 2712 2713 # return to original directory 2714 os.chdir(cwd) 2715 return
Read the PDB files contained in pdbdir and create the .pkl files needed for the
proteusPy.DisulfideLoader.DisulfideLoader class.
The Disulfide objects are contained in a DisulfideList object and
Dict within these files. In addition, .csv files containing all of
the torsions for the disulfides and problem IDs are written. The optional
dist_cutoff allows for removal of Disufides whose Cα-Cα distance is >
than the cutoff value. If it's -1.0 then the function keeps all Disulfides.
Parameters
- numb: number of entries to process, defaults to all
- verbose: more messages
- quiet: turns off DisulfideConstruction warnings
- pdbdir: path to PDB files
- datadir: path to resulting .pkl files
- picklefile: name of the disulfide .pkl file
- torsionfile: name of the disulfide torsion file .csv created
- problemfile: name of the .csv file containing problem ids
- dictfile: name of the .pkl file
- dist_cutoff: Ca distance cutoff to reject a Disulfide.
def
check_header_from_file(filename: str, model_numb=0, verbose=False, dbg=False) -> bool:
2718def check_header_from_file( 2719 filename: str, model_numb=0, verbose=False, dbg=False 2720) -> bool: 2721 """ 2722 Check the Disulfides in the PDB file. 2723 2724 NB: Requires EGS-Modified BIO.parse_pdb_header.py from https://github.com/suchanek/biopython/ 2725 2726 :param filename: the name of the PDB entry. 2727 :param model_numb: model number to use, defaults to 0 for single structure files. 2728 :param verbose: print info while parsing 2729 :return: True if parsable 2730 2731 Example: 2732 Assuming ```MODEL_DIR``` has the pdb5rsa.ent file we can load the disulfides 2733 with the following: 2734 2735 >>> from proteusPy.Disulfide import Disulfide, check_header_from_file 2736 >>> MODEL_DIR = '../data/' 2737 >>> OK = False 2738 >>> OK = check_header_from_file(f'{MODEL_DIR}pdb5rsa.ent', verbose=True) 2739 -> check_header_from_file() - Parsing file: ../data/pdb5rsa.ent: 2740 -> SSBond: 1: tmp: 26A - 84A 2741 -> SSBond: 2: tmp: 40A - 95A 2742 -> SSBond: 3: tmp: 58A - 110A 2743 -> SSBond: 4: tmp: 65A - 72A 2744 >>> OK 2745 True 2746 """ 2747 import os 2748 2749 i = 1 2750 proximal = distal = -1 2751 _chaina = None 2752 _chainb = None 2753 2754 parser = PDBParser(PERMISSIVE=True) 2755 2756 # Biopython uses the Structure -> Model -> Chain hierarchy to organize 2757 # structures. All are iterable. 2758 2759 structure = parser.get_structure("tmp", file=filename) 2760 struct_name = structure.get_id() 2761 2762 model = structure[model_numb] 2763 2764 if verbose: 2765 print(f"-> check_header_from_file() - Parsing file: {filename}:") 2766 2767 ssbond_dict = structure.header["ssbond"] # NB: this requires the modified code 2768 2769 # list of tuples with (proximal distal chaina chainb) 2770 ssbonds = parse_ssbond_header_rec(ssbond_dict) 2771 2772 for pair in ssbonds: 2773 # in the form (proximal, distal, chain) 2774 proximal = pair[0] 2775 distal = pair[1] 2776 2777 if not proximal.isnumeric() or not distal.isnumeric(): 2778 if verbose: 2779 mess = f" ! Cannot parse SSBond record (non-numeric IDs):\ 2780 {struct_name} Prox: {proximal} {chain1_id} Dist: {distal} {chain2_id}" 2781 warnings.warn(mess, DisulfideParseWarning) 2782 continue # was pass 2783 else: 2784 proximal = int(proximal) 2785 distal = int(distal) 2786 2787 chain1_id = pair[2] 2788 chain2_id = pair[3] 2789 2790 _chaina = model[chain1_id] 2791 _chainb = model[chain2_id] 2792 2793 if chain1_id != chain2_id: 2794 if verbose: 2795 mess = f" -> Cross Chain SS for: Prox: {proximal}{chain1_id} Dist: {distal}{chain2_id}" 2796 warnings.warn(mess, DisulfideParseWarning) 2797 pass # was break 2798 2799 try: 2800 prox_res = _chaina[proximal] 2801 dist_res = _chainb[distal] 2802 except KeyError: 2803 print( 2804 f" ! Cannot parse SSBond record (KeyError): {struct_name} Prox: <{proximal}> {chain1_id} Dist: <{distal}> {chain2_id}" 2805 ) 2806 return False 2807 2808 if (_chaina is not None) and (_chainb is not None): 2809 if _chaina[proximal].is_disordered() or _chainb[distal].is_disordered(): 2810 continue 2811 else: 2812 if verbose: 2813 print( 2814 f" -> SSBond: {i}: {struct_name}: {proximal}{chain1_id} - {distal}{chain2_id}" 2815 ) 2816 else: 2817 if dbg: 2818 print( 2819 f" -> NULL chain(s): {struct_name}: {proximal}{chain1_id} - {distal}{chain2_id}" 2820 ) 2821 i += 1 2822 return True
Check the Disulfides in the PDB file.
NB: Requires EGS-Modified BIO.parse_pdb_header.py from https://github.com/suchanek/biopython/
Parameters
- filename: the name of the PDB entry.
- model_numb: model number to use, defaults to 0 for single structure files.
- verbose: print info while parsing
Returns
True if parsable
Example:
Assuming MODEL_DIR has the pdb5rsa.ent file we can load the disulfides
with the following:
>>> from proteusPy.Disulfide import Disulfide, check_header_from_file
>>> MODEL_DIR = '../data/'
>>> OK = False
>>> OK = check_header_from_file(f'{MODEL_DIR}pdb5rsa.ent', verbose=True)
-> check_header_from_file() - Parsing file: ../data/pdb5rsa.ent:
-> SSBond: 1: tmp: 26A - 84A
-> SSBond: 2: tmp: 40A - 95A
-> SSBond: 3: tmp: 58A - 110A
-> SSBond: 4: tmp: 65A - 72A
>>> OK
True
def
check_header_from_id( struct_name: str, pdb_dir='.', model_numb=0, verbose=False, dbg=False) -> bool:
2825def check_header_from_id( 2826 struct_name: str, pdb_dir=".", model_numb=0, verbose=False, dbg=False 2827) -> bool: 2828 """ 2829 Check parsability PDB ID and initializes the Disulfide objects. 2830 Assumes the file is downloaded in ```MODEL_DIR``` path. 2831 2832 NB: Requires EGS-Modified BIO.parse_pdb_header.py from https://github.com/suchanek/biopython/ 2833 2834 :param struct_name: the name of the PDB entry. 2835 :param pdb_dir: path to the PDB files, defaults to PDB_DIR 2836 :param model_numb: model number to use, defaults to 0 for single structure files. 2837 :param verbose: print info while parsing 2838 :param dbg: Debugging Flag 2839 :return: True if OK, False otherwise 2840 2841 Example: 2842 Assuming the PDB_DIR has the pdb5rsa.ent file we can check the file thusly: 2843 2844 >>> from proteusPy.Disulfide import Disulfide, check_header_from_id 2845 >>> MODEL_DIR = '/Users/egs/PDB/good/' 2846 >>> OK = False 2847 >>> OK = check_header_from_id('5rsa', pdb_dir=MODEL_DIR, verbose=True) 2848 -> SSBond: 1: 5rsa: 26A - 84A 2849 -> SSBond: 2: 5rsa: 40A - 95A 2850 -> SSBond: 3: 5rsa: 58A - 110A 2851 -> SSBond: 4: 5rsa: 65A - 72A 2852 >>> OK 2853 True 2854 """ 2855 parser = PDBParser(PERMISSIVE=True, QUIET=True) 2856 structure = parser.get_structure(struct_name, file=f"{pdb_dir}pdb{struct_name}.ent") 2857 model = structure[0] 2858 2859 ssbond_dict = structure.header["ssbond"] # NB: this requires the modified code 2860 2861 bondlist = [] 2862 i = 0 2863 2864 # get a list of tuples containing the proximal, distal residue IDs for 2865 # all SSBonds in the chain. 2866 bondlist = parse_ssbond_header_rec(ssbond_dict) 2867 2868 if len(bondlist) == 0: 2869 if verbose: 2870 print("-> check_header_from_id(): no bonds found in bondlist.") 2871 return False 2872 2873 for pair in bondlist: 2874 # in the form (proximal, distal, chain) 2875 proximal = pair[0] 2876 distal = pair[1] 2877 chain1 = pair[2] 2878 chain2 = pair[3] 2879 2880 chaina = model[chain1] 2881 chainb = model[chain2] 2882 2883 try: 2884 prox_residue = chaina[proximal] 2885 dist_residue = chainb[distal] 2886 2887 prox_residue.disordered_select("CYS") 2888 dist_residue.disordered_select("CYS") 2889 2890 if ( 2891 prox_residue.get_resname() != "CYS" 2892 or dist_residue.get_resname() != "CYS" 2893 ): 2894 if verbose: 2895 print( 2896 f"build_disulfide() requires CYS at both residues:\ 2897 {prox_residue.get_resname()} {dist_residue.get_resname()}" 2898 ) 2899 return False 2900 except KeyError: 2901 if dbg: 2902 print( 2903 f"Keyerror: {struct_name}: {proximal} {chain1} - {distal} {chain2}" 2904 ) 2905 return False 2906 2907 if verbose: 2908 print( 2909 f" -> SSBond: {i+1}: {struct_name}: {proximal}{chain1} - {distal}{chain2}" 2910 ) 2911 2912 i += 1 2913 return True
Check parsability PDB ID and initializes the Disulfide objects.
Assumes the file is downloaded in MODEL_DIR path.
NB: Requires EGS-Modified BIO.parse_pdb_header.py from https://github.com/suchanek/biopython/
Parameters
- struct_name: the name of the PDB entry.
- pdb_dir: path to the PDB files, defaults to PDB_DIR
- model_numb: model number to use, defaults to 0 for single structure files.
- verbose: print info while parsing
- dbg: Debugging Flag
Returns
True if OK, False otherwise
Example: Assuming the PDB_DIR has the pdb5rsa.ent file we can check the file thusly:
>>> from proteusPy.Disulfide import Disulfide, check_header_from_id
>>> MODEL_DIR = '/Users/egs/PDB/good/'
>>> OK = False
>>> OK = check_header_from_id('5rsa', pdb_dir=MODEL_DIR, verbose=True)
-> SSBond: 1: 5rsa: 26A - 84A
-> SSBond: 2: 5rsa: 40A - 95A
-> SSBond: 3: 5rsa: 58A - 110A
-> SSBond: 4: 5rsa: 65A - 72A
>>> OK
True
def
Disulfide_Energy_Function(x: list) -> float:
2916def Disulfide_Energy_Function(x: list) -> float: 2917 """ 2918 Compute the approximate torsional energy (kcal/mpl) for the input dihedral angles. 2919 2920 :param x: A list of dihedral angles: [chi1, chi2, chi3, chi4, chi5] 2921 :return: Energy in kcal/mol 2922 2923 Example: 2924 >>> from proteusPy.Disulfide import Disulfide_Energy_Function 2925 >>> dihed = [-60.0, -60.0, -90.0, -60.0, -90.0] 2926 >>> res = Disulfide_Energy_Function(dihed) 2927 >>> res 2928 2.5999999999999996 2929 """ 2930 import numpy as np 2931 2932 chi1, chi2, chi3, chi4, chi5 = x 2933 energy = 2.0 * (np.cos(np.deg2rad(3.0 * chi1)) + np.cos(np.deg2rad(3.0 * chi5))) 2934 energy += np.cos(np.deg2rad(3.0 * chi2)) + np.cos(np.deg2rad(3.0 * chi4)) 2935 energy += ( 2936 3.5 * np.cos(np.deg2rad(2.0 * chi3)) 2937 + 0.6 * np.cos(np.deg2rad(3.0 * chi3)) 2938 + 10.1 2939 ) 2940 return energy
Compute the approximate torsional energy (kcal/mpl) for the input dihedral angles.
Parameters
- x: A list of dihedral angles: [chi1, chi2, chi3, chi4, chi5]
Returns
Energy in kcal/mol
Example:
>>> from proteusPy.Disulfide import Disulfide_Energy_Function
>>> dihed = [-60.0, -60.0, -90.0, -60.0, -90.0]
>>> res = Disulfide_Energy_Function(dihed)
>>> res
2.5999999999999996
2943def Minimize(inputSS: Disulfide) -> Disulfide: 2944 """ 2945 Minimizes the energy of a Disulfide object using the Nelder-Mead optimization method. 2946 2947 Parameters: 2948 inputSS (Disulfide): The Disulfide object to be minimized. 2949 2950 Returns: 2951 Disulfide: The minimized Disulfide object. 2952 2953 """ 2954 from scipy.optimize import minimize 2955 2956 from proteusPy import Disulfide, Disulfide_Energy_Function 2957 2958 initial_guess = inputSS.torsion_array 2959 result = minimize(Disulfide_Energy_Function, initial_guess, method="Nelder-Mead") 2960 minimum_conformation = result.x 2961 modelled_min = Disulfide("minimized") 2962 modelled_min.dihedrals = minimum_conformation 2963 modelled_min.build_yourself() 2964 return modelled_min
Minimizes the energy of a Disulfide object using the Nelder-Mead optimization method.
Parameters: inputSS (Disulfide): The Disulfide object to be minimized.
Returns: Disulfide: The minimized Disulfide object.