1"""Module calculates interactions between two molecules 2(proein-protein, protein-ligand, small-small). 3Currently following interacions are implemented: 4 * hydrogen bonds 5 * halogen bonds 6 * pi stacking (parallel and perpendicular) 7 * salt bridges 8 * hydrophobic contacts 9 * pi-cation 10 * metal coordination 11 * pi-metal 12""" 13 14import numpy as np 15from oddt.spatial import angle, angle_2v, distance 16 17__all__ = ['close_contacts', 18 'hbond_acceptor_donor', 19 'hbonds', 20 'halogenbond_acceptor_halogen', 21 'halogenbonds', 22 'pi_stacking', 23 'salt_bridge_plus_minus', 24 'salt_bridges', 25 'hydrophobic_contacts', 26 'pi_cation', 27 'acceptor_metal', 28 'pi_metal'] 29 30 31def close_contacts(x, y, cutoff, x_column='coords', y_column='coords', 32 cutoff_low=0.): 33 """Returns pairs of atoms which are within close contac distance cutoff. 34 The cutoff is semi-inclusive, i.e (cutoff_low, cutoff]. 35 36 Parameters 37 ---------- 38 x, y : atom_dict-type numpy array 39 Atom dictionaries generated by oddt.toolkit.Molecule objects. 40 41 cutoff : float 42 Cutoff distance for close contacts 43 44 x_column, ycolumn : string, (default='coords') 45 Column containing coordinates of atoms (or pseudo-atoms, 46 i.e. ring centroids) 47 48 cutoff_low : float (default=0.) 49 Lower bound of contacts to find (exclusive). Zero by default. 50 .. versionadded:: 0.6 51 52 Returns 53 ------- 54 x_, y_ : atom_dict-type numpy array 55 Aligned pairs of atoms in close contact for further processing. 56 """ 57 if len(x[x_column]) > 0 and len(x[x_column]) > 0: 58 d = distance(x[x_column], y[y_column]) 59 index = np.argwhere((d > cutoff_low) & (d <= cutoff)) 60 return x[index[:, 0]], y[index[:, 1]] 61 else: 62 return x[[]], y[[]] 63 64 65def hbond_acceptor_donor(mol1, mol2, cutoff=3.5, base_angle=120, tolerance=30): 66 """Returns pairs of acceptor-donor atoms, which meet H-bond criteria 67 68 Parameters 69 ---------- 70 mol1, mol2 : oddt.toolkit.Molecule object 71 Molecules to compute H-bond acceptor and H-bond donor pairs 72 73 cutoff : float, (default=3.5) 74 Distance cutoff for A-D pairs 75 76 base_angle : int, (default=120) 77 Base angle determining allowed direction of hydrogen bond formation, 78 which is devided by the number of neighbors of acceptor atom 79 to establish final directional angle 80 81 tolerance : int, (default=30) 82 Range (+/- tolerance) from perfect direction (base_angle/n_neighbors) 83 in which H-bonds are considered as strict. 84 85 Returns 86 ------- 87 a, d : atom_dict-type numpy array 88 Aligned arrays of atoms forming H-bond, firstly acceptors, 89 secondly donors. 90 91 strict : numpy array, dtype=bool 92 Boolean array align with atom pairs, informing whether atoms 93 form 'strict' H-bond (pass all angular cutoffs). If false, 94 only distance cutoff is met, therefore the bond is 'crude'. 95 """ 96 a, d = close_contacts(mol1.atom_dict[mol1.atom_dict['isacceptor']], 97 mol2.atom_dict[mol2.atom_dict['isdonor']], 98 cutoff) 99 # skip empty values 100 if len(a) > 0 and len(d) > 0: 101 angle1 = angle(d['coords'][:, np.newaxis, :], 102 a['coords'][:, np.newaxis, :], 103 a['neighbors']) 104 a_neighbors_num = np.sum(~np.isnan(a['neighbors'][:, :, 0]), axis=-1)[:, np.newaxis] 105 angle2 = angle(a['coords'][:, np.newaxis, :], 106 d['coords'][:, np.newaxis, :], 107 d['neighbors']) 108 d_neighbors_num = np.sum(~np.isnan(d['neighbors'][:, :, 0]), axis=-1)[:, np.newaxis] 109 strict = (((np.nan_to_num(angle1) > (base_angle / a_neighbors_num - tolerance)) | 110 np.isnan(angle1)) & 111 ((np.nan_to_num(angle2) > (base_angle / d_neighbors_num - tolerance)) | 112 np.isnan(angle2))).all(axis=-1) 113 return a, d, strict 114 else: 115 return a, d, np.array([], dtype=bool) 116 117 118def hbonds(mol1, mol2, *args, **kwargs): 119 """Calculates H-bonds between molecules 120 121 Parameters 122 ---------- 123 mol1, mol2 : oddt.toolkit.Molecule object 124 Molecules to compute H-bond acceptor and H-bond donor pairs 125 126 cutoff : float, (default=3.5) 127 Distance cutoff for A-D pairs 128 129 base_angle : int, (default=120) 130 Base angle determining allowed direction of hydrogen bond formation, 131 which is devided by the number of neighbors of acceptor atom 132 to establish final directional angle 133 134 tolerance : int, (default=30) 135 Range (+/- tolerance) from perfect direction (base_angle/n_neighbors) 136 in which H-bonds are considered as strict. 137 138 Returns 139 ------- 140 mol1_atoms, mol2_atoms : atom_dict-type numpy array 141 Aligned arrays of atoms forming H-bond 142 143 strict : numpy array, dtype=bool 144 Boolean array align with atom pairs, informing whether atoms 145 form 'strict' H-bond (pass all angular cutoffs). If false, 146 only distance cutoff is met, therefore the bond is 'crude'. 147 """ 148 a1, d1, s1 = hbond_acceptor_donor(mol1, mol2, *args, **kwargs) 149 a2, d2, s2 = hbond_acceptor_donor(mol2, mol1, *args, **kwargs) 150 return np.concatenate((a1, d2)), np.concatenate((d1, a2)), np.concatenate((s1, s2)) 151 152 153def halogenbond_acceptor_halogen(mol1, 154 mol2, 155 base_angle_acceptor=120, 156 base_angle_halogen=180, 157 tolerance=30, 158 cutoff=4): 159 """Returns pairs of acceptor-halogen atoms, which meet halogen bond criteria 160 161 Parameters 162 ---------- 163 mol1, mol2 : oddt.toolkit.Molecule object 164 Molecules to compute halogen bond acceptor and halogen pairs 165 166 cutoff : float, (default=4) 167 Distance cutoff for A-H pairs 168 169 base_angle_acceptor : int, (default=120) 170 Base angle determining allowed direction of halogen bond formation, 171 which is devided by the number of neighbors of acceptor atom 172 to establish final directional angle 173 174 base_angle_halogen : int (default=180) 175 Ideal base angle between halogen bond and halogen-neighbor bond 176 177 tolerance : int, (default=30) 178 Range (+/- tolerance) from perfect direction (base_angle/n_neighbors) 179 in which halogen bonds are considered as strict. 180 181 Returns 182 ------- 183 a, h : atom_dict-type numpy array 184 Aligned arrays of atoms forming halogen bond, firstly acceptors, 185 secondly halogens 186 187 strict : numpy array, dtype=bool 188 Boolean array align with atom pairs, informing whether atoms 189 form 'strict' halogen bond (pass all angular cutoffs). If false, 190 only distance cutoff is met, therefore the bond is 'crude'. 191 """ 192 a, h = close_contacts(mol1.atom_dict[mol1.atom_dict['isacceptor']], 193 mol2.atom_dict[mol2.atom_dict['ishalogen']], 194 cutoff) 195 # skip empty values 196 if len(a) > 0 and len(h) > 0: 197 angle1 = angle(h['coords'][:, np.newaxis, :], 198 a['coords'][:, np.newaxis, :], 199 a['neighbors']) 200 angle2 = angle(a['coords'][:, np.newaxis, :], 201 h['coords'][:, np.newaxis, :], 202 h['neighbors']) 203 a_neighbors_num = np.sum(~np.isnan(a['neighbors'][:, :, 0]), axis=-1)[:, np.newaxis] 204 h_neighbors_num = np.sum(~np.isnan(h['neighbors'][:, :, 0]), axis=-1)[:, np.newaxis] 205 strict = (((np.nan_to_num(angle1) > (base_angle_acceptor / a_neighbors_num - tolerance)) | 206 np.isnan(angle1)) & 207 ((np.nan_to_num(angle2) > (base_angle_halogen / h_neighbors_num - tolerance)) | 208 np.isnan(angle2))).all(axis=-1) 209 return a, h, strict 210 else: 211 return a, h, np.array([], dtype=bool) 212 213 214def halogenbonds(mol1, mol2, **kwargs): 215 """Calculates halogen bonds between molecules 216 217 Parameters 218 ---------- 219 mol1, mol2 : oddt.toolkit.Molecule object 220 Molecules to compute halogen bond acceptor and halogen pairs 221 222 cutoff : float, (default=4) 223 Distance cutoff for A-H pairs 224 225 base_angle_acceptor : int, (default=120) 226 Base angle determining allowed direction of halogen bond formation, 227 which is devided by the number of neighbors of acceptor atom 228 to establish final directional angle 229 230 base_angle_halogen : int (default=180) 231 Ideal base angle between halogen bond and halogen-neighbor bond 232 233 tolerance : int, (default=30) 234 Range (+/- tolerance) from perfect direction (base_angle/n_neighbors) 235 in which halogen bonds are considered as strict. 236 237 Returns 238 ------- 239 mol1_atoms, mol2_atoms : atom_dict-type numpy array 240 Aligned arrays of atoms forming halogen bond 241 242 strict : numpy array, dtype=bool 243 Boolean array align with atom pairs, informing whether atoms 244 form 'strict' halogen bond (pass all angular cutoffs). If false, 245 only distance cutoff is met, therefore the bond is 'crude'. 246 """ 247 a1, h1, s1 = halogenbond_acceptor_halogen(mol1, mol2, **kwargs) 248 a2, h2, s2 = halogenbond_acceptor_halogen(mol2, mol1, **kwargs) 249 return np.concatenate((a1, h2)), np.concatenate((h1, a2)), np.concatenate((s1, s2)) 250 251 252def pi_stacking(mol1, mol2, cutoff=5, tolerance=30): 253 """Returns pairs of rings, which meet pi stacking criteria 254 255 Parameters 256 ---------- 257 mol1, mol2 : oddt.toolkit.Molecule object 258 Molecules to compute ring pairs 259 260 cutoff : float, (default=5) 261 Distance cutoff for Pi-stacking pairs 262 263 tolerance : int, (default=30) 264 Range (+/- tolerance) from perfect direction (parallel or 265 perpendicular) in which pi-stackings are considered as strict. 266 267 Returns 268 ------- 269 r1, r2 : ring_dict-type numpy array 270 Aligned arrays of rings forming pi-stacking 271 272 strict_parallel : numpy array, dtype=bool 273 Boolean array align with ring pairs, informing whether rings 274 form 'strict' parallel pi-stacking. If false, only distance cutoff is met, 275 therefore the stacking is 'crude'. 276 277 strict_perpendicular : numpy array, dtype=bool 278 Boolean array align with ring pairs, informing whether rings 279 form 'strict' perpendicular pi-stacking (T-shaped, T-face, etc.). 280 If false, only distance cutoff is met, therefore the stacking is 'crude'. 281 """ 282 r1, r2 = close_contacts(mol1.ring_dict, 283 mol2.ring_dict, 284 cutoff, 285 x_column='centroid', 286 y_column='centroid') 287 if len(r1) > 0 and len(r2) > 0: 288 angle1 = angle_2v(r1['vector'], r2['vector']) 289 angle2 = angle(r1['vector'] + r1['centroid'], 290 r1['centroid'], 291 r2['centroid']) 292 strict_parallel = (((angle1 > 180 - tolerance) | (angle1 < tolerance)) & 293 ((angle2 > 180 - tolerance) | (angle2 < tolerance))) 294 strict_perpendicular = (((angle1 > 90 - tolerance) & (angle1 < 90 + tolerance)) & 295 ((angle2 > 180 - tolerance) | (angle2 < tolerance))) 296 return r1, r2, strict_parallel, strict_perpendicular 297 else: 298 return r1, r2, np.array([], dtype=bool), np.array([], dtype=bool) 299 300 301def salt_bridge_plus_minus(mol1, mol2, cutoff=4): 302 """Returns pairs of plus-mins atoms, which meet salt bridge criteria 303 304 Parameters 305 ---------- 306 mol1, mol2 : oddt.toolkit.Molecule object 307 Molecules to compute plus and minus pairs 308 309 cutoff : float, (default=4) 310 Distance cutoff for A-H pairs 311 312 Returns 313 ------- 314 plus, minus : atom_dict-type numpy array 315 Aligned arrays of atoms forming salt bridge, firstly plus, secondly minus 316 317 """ 318 m1_plus, m2_minus = close_contacts(mol1.atom_dict[mol1.atom_dict['isplus']], 319 mol2.atom_dict[mol2.atom_dict['isminus']], 320 cutoff) 321 return m1_plus, m2_minus 322 323 324def salt_bridges(mol1, mol2, *args, **kwargs): 325 """Calculates salt bridges between molecules 326 327 Parameters 328 ---------- 329 mol1, mol2 : oddt.toolkit.Molecule object 330 Molecules to compute plus and minus pairs 331 332 cutoff : float, (default=4) 333 Distance cutoff for plus-minus pairs 334 335 Returns 336 ------- 337 mol1_atoms, mol2_atoms : atom_dict-type numpy array 338 Aligned arrays of atoms forming salt bridges 339 """ 340 m1_plus, m2_minus = salt_bridge_plus_minus(mol1, mol2, *args, **kwargs) 341 m2_plus, m1_minus = salt_bridge_plus_minus(mol2, mol1, *args, **kwargs) 342 return np.concatenate((m1_plus, m1_minus)), np.concatenate((m2_minus, m2_plus)) 343 344 345def hydrophobic_contacts(mol1, mol2, cutoff=4): 346 """Calculates hydrophobic contacts between molecules 347 348 Parameters 349 ---------- 350 mol1, mol2 : oddt.toolkit.Molecule object 351 Molecules to compute hydrophobe pairs 352 353 cutoff : float, (default=4) 354 Distance cutoff for hydrophobe pairs 355 356 Returns 357 ------- 358 mol1_atoms, mol2_atoms : atom_dict-type numpy array 359 Aligned arrays of atoms forming hydrophobic contacts 360 361 """ 362 h1, h2 = close_contacts(mol1.atom_dict[mol1.atom_dict['ishydrophobe']], 363 mol2.atom_dict[mol2.atom_dict['ishydrophobe']], 364 cutoff) 365 return h1, h2 366 367 368def pi_cation(mol1, mol2, cutoff=5, tolerance=30): 369 """Returns pairs of ring-cation atoms, which meet pi-cation criteria 370 371 Parameters 372 ---------- 373 mol1, mol2 : oddt.toolkit.Molecule object 374 Molecules to compute ring-cation pairs 375 376 cutoff : float, (default=5) 377 Distance cutoff for Pi-cation pairs 378 379 tolerance : int, (default=30) 380 Range (+/- tolerance) from perfect direction (perpendicular) 381 in which pi-cation are considered as strict. 382 383 Returns 384 ------- 385 r1 : ring_dict-type numpy array 386 Aligned rings forming pi-stacking 387 388 plus2 : atom_dict-type numpy array 389 Aligned cations forming pi-cation 390 391 strict_parallel : numpy array, dtype=bool 392 Boolean array align with ring-cation pairs, informing whether 393 they form 'strict' pi-cation. If false, only distance cutoff is met, 394 therefore the interaction is 'crude'. 395 396 """ 397 r1, plus2 = close_contacts(mol1.ring_dict, 398 mol2.atom_dict[mol2.atom_dict['isplus']], 399 cutoff, 400 x_column='centroid') 401 if len(r1) > 0 and len(plus2) > 0: 402 angle1 = angle_2v(r1['vector'], plus2['coords'] - r1['centroid']) 403 strict = (angle1 > 180 - tolerance) | (angle1 < tolerance) 404 return r1, plus2, strict 405 else: 406 return r1, plus2, np.array([], dtype=bool) 407 408 409def acceptor_metal(mol1, mol2, base_angle=120, tolerance=30, cutoff=4): 410 """Returns pairs of acceptor-metal atoms, which meet metal coordination criteria 411 Note: This function is directional (mol1 holds acceptors, mol2 holds metals) 412 413 Parameters 414 ---------- 415 mol1, mol2 : oddt.toolkit.Molecule object 416 Molecules to compute acceptor and metal pairs 417 418 cutoff : float, (default=4) 419 Distance cutoff for A-M pairs 420 421 base_angle : int, (default=120) 422 Base angle determining allowed direction of metal coordination, 423 which is devided by the number of neighbors of acceptor atom 424 to establish final directional angle 425 426 tolerance : int, (default=30) 427 Range (+/- tolerance) from perfect direction (base_angle/n_neighbors) 428 in metal coordination are considered as strict. 429 430 Returns 431 ------- 432 a, d : atom_dict-type numpy array 433 Aligned arrays of atoms forming metal coordination, 434 firstly acceptors, secondly metals. 435 436 strict : numpy array, dtype=bool 437 Boolean array align with atom pairs, informing whether atoms 438 form 'strict' metal coordination (pass all angular cutoffs). 439 If false, only distance cutoff is met, therefore the interaction 440 is 'crude'. 441 """ 442 a, m = close_contacts(mol1.atom_dict[mol1.atom_dict['isacceptor']], 443 mol2.atom_dict[mol2.atom_dict['ismetal']], 444 cutoff) 445 # skip empty values 446 if len(a) > 0 and len(m) > 0: 447 angle1 = angle(m['coords'][:, np.newaxis, :], 448 a['coords'][:, np.newaxis, :], 449 a['neighbors']) 450 a_neighbors_num = np.sum(~np.isnan(a['neighbors'][:, :, 0]), axis=-1)[:, np.newaxis] 451 strict = ((np.nan_to_num(angle1) > (base_angle / a_neighbors_num - tolerance)) | 452 np.isnan(angle1)).all(axis=-1) 453 return a, m, strict 454 else: 455 return a, m, np.array([], dtype=bool) 456 457 458def pi_metal(mol1, mol2, cutoff=5, tolerance=30): 459 """Returns pairs of ring-metal atoms, which meet pi-metal criteria 460 461 Parameters 462 ---------- 463 mol1, mol2 : oddt.toolkit.Molecule object 464 Molecules to compute ring-metal pairs 465 466 cutoff : float, (default=5) 467 Distance cutoff for Pi-metal pairs 468 469 tolerance : int, (default=30) 470 Range (+/- tolerance) from perfect direction (perpendicular) 471 in which pi-metal are considered as strict. 472 473 Returns 474 ------- 475 r1 : ring_dict-type numpy array 476 Aligned rings forming pi-metal 477 478 m : atom_dict-type numpy array 479 Aligned metals forming pi-metal 480 481 strict_parallel : numpy array, dtype=bool 482 Boolean array align with ring-metal pairs, informing whether 483 they form 'strict' pi-metal. If false, only distance cutoff is met, 484 therefore the interaction is 'crude'. 485 486 """ 487 r1, m = close_contacts(mol1.ring_dict, 488 mol2.atom_dict[mol2.atom_dict['ismetal']], 489 cutoff, 490 x_column='centroid') 491 if len(r1) > 0 and len(m) > 0: 492 angle1 = angle_2v(r1['vector'], m['coords'] - r1['centroid']) 493 strict = (angle1 > 180 - tolerance) | (angle1 < tolerance) 494 return r1, m, strict 495 else: 496 return r1, m, np.array([], dtype=bool) 497