0
# Topology Handling
1
2
MDAnalysis provides comprehensive tools for managing molecular topology information, including atomic connectivity, bonds, angles, dihedrals, and topology attributes. The topology system allows both reading existing topology data and dynamically modifying topology information.
3
4
## Overview
5
6
Topology in MDAnalysis encompasses:
7
- **Static Structure**: Atom types, names, masses, charges
8
- **Connectivity**: Bonds, angles, dihedrals, impropers
9
- **Hierarchy**: Atoms belong to residues, residues belong to segments
10
- **Dynamic Attributes**: User-defined properties that can be added/modified
11
12
## Topology Attributes
13
14
### Core Topology Attributes
15
16
```python { .api }
17
# Access topology attributes from atoms
18
atoms = u.atoms
19
20
# Basic atomic properties
21
names = atoms.names # Atom names (e.g., 'CA', 'CB', 'N')
22
types = atoms.types # Atom types (e.g., 'CT', 'NH1', 'O')
23
masses = atoms.masses # Atomic masses in amu
24
charges = atoms.charges # Partial charges in elementary charge units
25
radii = atoms.radii # Van der Waals radii in Angstrom
26
27
# Residue-level properties
28
resnames = atoms.resnames # Residue names (e.g., 'ALA', 'GLY', 'PRO')
29
resids = atoms.resids # Residue IDs/numbers
30
resnums = atoms.resnums # Sequential residue numbers
31
32
# Segment-level properties
33
segids = atoms.segids # Segment identifiers
34
```
35
36
### Adding Topology Attributes
37
38
```python { .api }
39
def add_TopologyAttr(universe, topologyattr, values=None):
40
"""
41
Add a topology attribute to the Universe.
42
43
Parameters
44
----------
45
universe : Universe
46
Universe to modify.
47
topologyattr : str or TopologyAttr
48
Name of topology attribute or TopologyAttr class instance.
49
values : array-like, optional
50
Values for the new attribute. If None, default values are used.
51
52
Examples
53
--------
54
>>> # Add masses for all atoms
55
>>> u.add_TopologyAttr('masses', [12.0, 14.0, 16.0] * (u.atoms.n_atoms // 3))
56
57
>>> # Add custom attribute
58
>>> u.add_TopologyAttr('occupancy', np.ones(u.atoms.n_atoms))
59
"""
60
61
def del_TopologyAttr(universe, topologyattr):
62
"""
63
Remove a topology attribute from the Universe.
64
65
Parameters
66
----------
67
universe : Universe
68
Universe to modify.
69
topologyattr : str
70
Name of topology attribute to remove.
71
72
Examples
73
--------
74
>>> u.del_TopologyAttr('masses')
75
>>> # Now u.atoms.masses will raise AttributeError
76
"""
77
78
# Available topology attributes
79
class AtomAttr:
80
"""Base class for atom-level topology attributes."""
81
82
# Standard atom attributes
83
names = "Atom names"
84
types = "Atom types"
85
masses = "Atomic masses"
86
charges = "Partial charges"
87
radii = "Van der Waals radii"
88
bfactors = "B-factors from crystallography"
89
occupancies = "Occupancies from crystallography"
90
altlocs = "Alternate locations from crystallography"
91
elements = "Chemical elements"
92
ids = "Atom IDs/indices"
93
94
class ResidueAttr:
95
"""Base class for residue-level topology attributes."""
96
97
resnames = "Residue names"
98
resids = "Residue IDs"
99
resnums = "Sequential residue numbers"
100
resnames = "Residue names"
101
icodes = "Insertion codes"
102
103
class SegmentAttr:
104
"""Base class for segment-level topology attributes."""
105
106
segids = "Segment identifiers"
107
```
108
109
## Connectivity Information
110
111
### Bonds
112
113
```python { .api }
114
def add_bonds(universe, values, types=None, guessed=False, order=None):
115
"""
116
Add bond information to the Universe.
117
118
Parameters
119
----------
120
universe : Universe
121
Universe to modify.
122
values : array-like
123
Bond connectivity as array of shape (n_bonds, 2) containing
124
atom indices for each bond.
125
types : array-like, optional
126
Bond types for each bond (e.g., 'single', 'double', 'aromatic').
127
guessed : bool, optional
128
Whether bonds were guessed from distances (default False).
129
order : array-like, optional
130
Bond orders as numerical values (1.0, 2.0, 1.5 for aromatic).
131
132
Examples
133
--------
134
>>> # Add bonds between consecutive atoms
135
>>> bond_list = [(i, i+1) for i in range(u.atoms.n_atoms-1)]
136
>>> u.add_bonds(bond_list)
137
138
>>> # Add typed bonds
139
>>> u.add_bonds([(0, 1), (1, 2)], types=['single', 'double'])
140
141
>>> # Add bonds with orders
142
>>> u.add_bonds([(0, 1), (2, 3)], order=[1.0, 2.0])
143
"""
144
145
def delete_bonds(universe, values):
146
"""
147
Remove specific bonds from the Universe.
148
149
Parameters
150
----------
151
universe : Universe
152
Universe to modify.
153
values : array-like
154
Bonds to remove as array of shape (n_bonds, 2).
155
156
Examples
157
--------
158
>>> # Remove specific bonds
159
>>> bonds_to_remove = [(10, 11), (15, 16)]
160
>>> u.delete_bonds(bonds_to_remove)
161
"""
162
163
def guess_bonds(atomgroup, vdwradii=None, fudge_factor=0.55, lower_bound=0.1):
164
"""
165
Guess bonds based on inter-atomic distances.
166
167
Parameters
168
----------
169
atomgroup : AtomGroup
170
Atoms for bond guessing.
171
vdwradii : dict, optional
172
Van der Waals radii for elements. If None, uses default values.
173
fudge_factor : float, optional
174
Tolerance factor for bond distance (default 0.55).
175
lower_bound : float, optional
176
Minimum distance for bond formation (default 0.1).
177
178
Returns
179
-------
180
numpy.ndarray
181
Array of shape (n_bonds, 2) containing guessed bonds.
182
183
Examples
184
--------
185
>>> # Guess all bonds in system
186
>>> bonds = u.atoms.guess_bonds()
187
>>> u.add_bonds(bonds, guessed=True)
188
189
>>> # Custom radii for specific elements
190
>>> custom_radii = {'C': 0.77, 'N': 0.70, 'O': 0.66}
191
>>> bonds = u.atoms.guess_bonds(vdwradii=custom_radii)
192
"""
193
194
# Access existing bonds
195
bonds = u.atoms.bonds # TopologyGroup containing all bonds
196
bond_indices = bonds.indices # Array of shape (n_bonds, 2)
197
bond_types = bonds.types # Bond types if available
198
```
199
200
### Angles and Dihedrals
201
202
```python { .api }
203
def add_angles(universe, values, types=None, guessed=False):
204
"""
205
Add angle information to the Universe.
206
207
Parameters
208
----------
209
universe : Universe
210
Universe to modify.
211
values : array-like
212
Angle connectivity as array of shape (n_angles, 3) containing
213
atom indices (i, j, k) where j is the central atom.
214
types : array-like, optional
215
Angle types for each angle.
216
guessed : bool, optional
217
Whether angles were guessed from connectivity (default False).
218
219
Examples
220
--------
221
>>> # Add angles from bond connectivity
222
>>> angles = [(0, 1, 2), (1, 2, 3), (2, 3, 4)]
223
>>> u.add_angles(angles)
224
"""
225
226
def add_dihedrals(universe, values, types=None, guessed=False):
227
"""
228
Add dihedral angle information to the Universe.
229
230
Parameters
231
----------
232
universe : Universe
233
Universe to modify.
234
values : array-like
235
Dihedral connectivity as array of shape (n_dihedrals, 4) containing
236
atom indices (i, j, k, l) defining the dihedral.
237
types : array-like, optional
238
Dihedral types for each dihedral.
239
guessed : bool, optional
240
Whether dihedrals were guessed from connectivity (default False).
241
242
Examples
243
--------
244
>>> # Add backbone dihedrals
245
>>> dihedrals = [(0, 1, 2, 3), (1, 2, 3, 4)]
246
>>> u.add_dihedrals(dihedrals)
247
"""
248
249
def add_impropers(universe, values, types=None, guessed=False):
250
"""
251
Add improper dihedral information to the Universe.
252
253
Parameters
254
----------
255
universe : Universe
256
Universe to modify.
257
values : array-like
258
Improper connectivity as array of shape (n_impropers, 4).
259
types : array-like, optional
260
Improper types for each improper.
261
guessed : bool, optional
262
Whether impropers were guessed (default False).
263
264
Examples
265
--------
266
>>> # Add improper dihedrals for planarity
267
>>> impropers = [(0, 1, 2, 3), (4, 5, 6, 7)]
268
>>> u.add_impropers(impropers)
269
"""
270
271
# Access connectivity information
272
angles = u.atoms.angles # All angles
273
dihedrals = u.atoms.dihedrals # All dihedrals
274
impropers = u.atoms.impropers # All impropers
275
276
# Get indices for connectivity
277
angle_indices = angles.indices # Shape: (n_angles, 3)
278
dihedral_indices = dihedrals.indices # Shape: (n_dihedrals, 4)
279
```
280
281
## Topology Modification
282
283
### Adding Atoms, Residues, and Segments
284
285
```python { .api }
286
def add_Residue(universe, segment=None, **attrs):
287
"""
288
Add a new residue to the Universe.
289
290
Parameters
291
----------
292
universe : Universe
293
Universe to modify.
294
segment : Segment, optional
295
Segment to add residue to. If None, added to first segment.
296
**attrs
297
Attributes for the new residue (resname, resid, etc.).
298
299
Returns
300
-------
301
Residue
302
New residue object.
303
304
Examples
305
--------
306
>>> # Add new residue
307
>>> new_res = u.add_Residue(resname='LIG', resid=999)
308
309
>>> # Add to specific segment
310
>>> seg = u.segments[0]
311
>>> new_res = u.add_Residue(segment=seg, resname='WAT', resid=1000)
312
"""
313
314
def add_Segment(universe, **attrs):
315
"""
316
Add a new segment to the Universe.
317
318
Parameters
319
----------
320
universe : Universe
321
Universe to modify.
322
**attrs
323
Attributes for the new segment (segid, etc.).
324
325
Returns
326
-------
327
Segment
328
New segment object.
329
330
Examples
331
--------
332
>>> # Add new segment
333
>>> new_seg = u.add_Segment(segid='LIG')
334
"""
335
336
# Adding atoms requires more complex topology manipulation
337
# Usually done through specialized functions or parsers
338
```
339
340
### Topology Guessing
341
342
```python { .api }
343
def guess_TopologyAttrs(universe, to_guess=None, skip_fields=None):
344
"""
345
Guess missing topology attributes from available information.
346
347
Parameters
348
----------
349
universe : Universe
350
Universe to modify.
351
to_guess : list, optional
352
List of attributes to guess. If None, guesses all possible.
353
skip_fields : list, optional
354
List of attributes to skip during guessing.
355
356
Examples
357
--------
358
>>> # Guess all missing attributes
359
>>> u.guess_TopologyAttrs()
360
361
>>> # Guess only specific attributes
362
>>> u.guess_TopologyAttrs(to_guess=['masses', 'bonds'])
363
364
>>> # Skip certain fields
365
>>> u.guess_TopologyAttrs(skip_fields=['charges'])
366
"""
367
368
# Available guessers
369
from MDAnalysis.guesser import (
370
guess_bonds, # Guess bonds from distances
371
guess_angles, # Guess angles from bonds
372
guess_dihedrals, # Guess dihedrals from angles
373
guess_atom_type, # Guess atom types from names
374
guess_atom_mass, # Guess masses from types/elements
375
guess_atom_charge # Guess charges from types
376
)
377
```
378
379
## Topology Groups
380
381
### Bond Groups
382
383
```python { .api }
384
class TopologyGroup:
385
"""
386
Base class for topology objects like bonds, angles, dihedrals.
387
388
Provides array-like access to connectivity information with
389
additional methods for analysis.
390
"""
391
392
@property
393
def indices(self):
394
"""
395
Connectivity indices for topology objects.
396
397
Returns
398
-------
399
numpy.ndarray
400
Array of atom indices defining connectivity.
401
Shape depends on topology type:
402
- Bonds: (n_bonds, 2)
403
- Angles: (n_angles, 3)
404
- Dihedrals: (n_dihedrals, 4)
405
"""
406
407
@property
408
def types(self):
409
"""
410
Types for each topology object.
411
412
Returns
413
-------
414
numpy.ndarray or None
415
Array of type strings, None if no type information.
416
"""
417
418
def __len__(self):
419
"""Number of topology objects."""
420
421
def __getitem__(self, key):
422
"""Access specific topology objects by index."""
423
424
# Specific topology groups
425
class BondGroup(TopologyGroup):
426
"""Group of bonds with bond-specific methods."""
427
428
def length(self, pbc=False):
429
"""
430
Calculate bond lengths.
431
432
Parameters
433
----------
434
pbc : bool, optional
435
Whether to consider periodic boundary conditions.
436
437
Returns
438
-------
439
numpy.ndarray
440
Array of bond lengths in Angstrom.
441
"""
442
443
def partner(self, atom):
444
"""
445
Get bonding partners for specified atom.
446
447
Parameters
448
----------
449
atom : Atom or int
450
Atom or atom index.
451
452
Returns
453
-------
454
AtomGroup
455
Atoms bonded to specified atom.
456
"""
457
458
class AngleGroup(TopologyGroup):
459
"""Group of angles with angle-specific methods."""
460
461
def angle(self, pbc=False):
462
"""
463
Calculate angle values.
464
465
Parameters
466
----------
467
pbc : bool, optional
468
Whether to consider periodic boundary conditions.
469
470
Returns
471
-------
472
numpy.ndarray
473
Array of angle values in degrees.
474
"""
475
476
class DihedralGroup(TopologyGroup):
477
"""Group of dihedrals with dihedral-specific methods."""
478
479
def dihedral(self, pbc=False):
480
"""
481
Calculate dihedral angle values.
482
483
Parameters
484
----------
485
pbc : bool, optional
486
Whether to consider periodic boundary conditions.
487
488
Returns
489
-------
490
numpy.ndarray
491
Array of dihedral angle values in degrees.
492
"""
493
```
494
495
## Topology Parsers
496
497
### Parser Classes
498
499
```python { .api }
500
# Available parsers for different formats
501
from MDAnalysis.topology import (
502
PSFParser, # CHARMM PSF format
503
PDBParser, # Protein Data Bank format
504
GROParser, # GROMACS GRO format
505
TPRParser, # GROMACS TPR binary format
506
PrmtopParser, # AMBER topology format
507
MOL2Parser, # Tripos MOL2 format
508
PQRParser, # PQR format (PDB with charges/radii)
509
CRDParser, # CHARMM CRD format
510
XYZParser # Generic XYZ format
511
)
512
513
class TopologyReaderBase:
514
"""
515
Base class for topology parsers.
516
517
All topology parsers inherit from this class and implement
518
the parse() method to extract topology information.
519
"""
520
521
def parse(self, **kwargs):
522
"""
523
Parse topology information from file.
524
525
Returns
526
-------
527
dict
528
Dictionary containing parsed topology data with keys:
529
- 'atoms': atom information
530
- 'bonds': bond connectivity
531
- 'residues': residue information
532
- 'segments': segment information
533
- Additional format-specific data
534
"""
535
536
# Example: Custom topology parsing
537
class CustomParser(TopologyReaderBase):
538
"""Example custom topology parser."""
539
540
format = 'CUSTOM'
541
542
def parse(self, **kwargs):
543
"""Parse custom topology format."""
544
# Implementation details
545
atoms = self._parse_atoms()
546
bonds = self._parse_bonds()
547
residues = self._parse_residues()
548
549
return {
550
'atoms': atoms,
551
'bonds': bonds,
552
'residues': residues
553
}
554
```
555
556
### Format-Specific Features
557
558
#### CHARMM PSF
559
560
```python { .api }
561
# PSF files contain complete topology information
562
u = mda.Universe("system.psf")
563
564
# Access PSF-specific information
565
print(u.atoms.types) # CHARMM atom types
566
print(u.atoms.charges) # Partial charges
567
print(u.atoms.bonds) # Complete bond connectivity
568
print(u.atoms.angles) # All angles
569
print(u.atoms.dihedrals) # All dihedrals
570
print(u.atoms.impropers) # Improper dihedrals
571
```
572
573
#### GROMACS TPR
574
575
```python { .api }
576
# TPR files are binary and contain full force field information
577
u = mda.Universe("topol.tpr")
578
579
# TPR provides extensive topology data
580
print(u.atoms.masses) # Atomic masses
581
print(u.atoms.charges) # Partial charges
582
print(u.atoms.types) # GROMACS atom types
583
print(u.atoms.bonds) # Bond connectivity
584
```
585
586
#### AMBER PRMTOP
587
588
```python { .api }
589
# PRMTOP files contain AMBER force field topology
590
u = mda.Universe("system.prmtop")
591
592
# AMBER-specific topology information
593
print(u.atoms.types) # AMBER atom types
594
print(u.atoms.charges) # Partial charges
595
print(u.atoms.masses) # Atomic masses
596
print(u.atoms.bonds) # Bond information
597
```
598
599
## Usage Examples
600
601
### Basic Topology Operations
602
603
```python { .api }
604
# Create universe and examine topology
605
u = mda.Universe("topology.psf", "trajectory.dcd")
606
607
# Check available topology attributes
608
print("Available attributes:")
609
for attr in ['names', 'types', 'masses', 'charges']:
610
if hasattr(u.atoms, attr):
611
print(f" {attr}: {getattr(u.atoms, attr)[:5]}...") # First 5 values
612
613
# Connectivity information
614
print(f"Number of bonds: {len(u.atoms.bonds)}")
615
print(f"Number of angles: {len(u.atoms.angles)}")
616
print(f"Number of dihedrals: {len(u.atoms.dihedrals)}")
617
618
# Bond analysis
619
bond_lengths = u.atoms.bonds.length()
620
print(f"Average bond length: {np.mean(bond_lengths):.2f} Å")
621
```
622
623
### Adding Missing Topology
624
625
```python { .api }
626
# Start with minimal topology (e.g., from PDB)
627
u = mda.Universe("structure.pdb")
628
629
# Add missing information
630
if not hasattr(u.atoms, 'masses'):
631
print("Guessing atomic masses...")
632
u.guess_TopologyAttrs(to_guess=['masses'])
633
634
if not hasattr(u.atoms, 'bonds'):
635
print("Guessing bond connectivity...")
636
u.guess_TopologyAttrs(to_guess=['bonds'])
637
638
# Add custom attributes
639
occupancies = np.random.random(u.atoms.n_atoms)
640
u.add_TopologyAttr('occupancy', occupancies)
641
642
# Now access new attributes
643
print(f"First 5 masses: {u.atoms.masses[:5]}")
644
print(f"Number of bonds: {len(u.atoms.bonds)}")
645
print(f"Average occupancy: {np.mean(u.atoms.occupancy):.2f}")
646
```
647
648
### Connectivity Analysis
649
650
```python { .api }
651
def analyze_connectivity(universe):
652
"""
653
Analyze connectivity patterns in molecular system.
654
"""
655
atoms = universe.atoms
656
657
# Bond statistics
658
if hasattr(atoms, 'bonds'):
659
bonds = atoms.bonds
660
bond_lengths = bonds.length()
661
662
print("Bond Analysis:")
663
print(f" Total bonds: {len(bonds)}")
664
print(f" Average length: {np.mean(bond_lengths):.2f} ± {np.std(bond_lengths):.2f} Å")
665
print(f" Length range: {np.min(bond_lengths):.2f} - {np.max(bond_lengths):.2f} Å")
666
667
# Coordination numbers
668
coordination = {}
669
for atom in atoms:
670
n_bonds = len(atom.bonds)
671
coordination[n_bonds] = coordination.get(n_bonds, 0) + 1
672
673
print("\nCoordination Analysis:")
674
for coord, count in sorted(coordination.items()):
675
print(f" {count} atoms with {coord} bonds")
676
677
# Fragment analysis
678
if hasattr(atoms, 'bonds'):
679
fragments = atoms.split('fragments')
680
print(f"\nFragment Analysis:")
681
print(f" Number of fragments: {len(fragments)}")
682
fragment_sizes = [len(frag) for frag in fragments]
683
print(f" Largest fragment: {max(fragment_sizes)} atoms")
684
print(f" Average fragment size: {np.mean(fragment_sizes):.1f} atoms")
685
686
# Run connectivity analysis
687
analyze_connectivity(u)
688
```
689
690
### Dynamic Topology Modification
691
692
```python { .api }
693
def modify_topology_example(universe):
694
"""
695
Example of dynamic topology modification during simulation.
696
"""
697
# Monitor and modify bonds based on distance criteria
698
cutoff = 1.8 # Bond formation cutoff in Angstrom
699
700
for ts in universe.trajectory:
701
# Get current positions
702
positions = universe.atoms.positions
703
704
# Calculate all pairwise distances
705
from MDAnalysis.lib.distances import distance_array
706
dist_matrix = distance_array(positions, positions,
707
box=universe.dimensions)
708
709
# Find potential new bonds (close atoms not already bonded)
710
existing_bonds = set()
711
if hasattr(universe.atoms, 'bonds'):
712
for bond in universe.atoms.bonds:
713
i, j = bond.indices[0]
714
existing_bonds.add((min(i, j), max(i, j)))
715
716
# Identify new bonds
717
new_bonds = []
718
for i in range(len(positions)):
719
for j in range(i + 1, len(positions)):
720
if dist_matrix[i, j] < cutoff and (i, j) not in existing_bonds:
721
new_bonds.append((i, j))
722
723
# Add new bonds if found
724
if new_bonds:
725
print(f"Frame {ts.frame}: Adding {len(new_bonds)} new bonds")
726
universe.add_bonds(new_bonds, guessed=True)
727
728
# Break bonds that are too long
729
if hasattr(universe.atoms, 'bonds'):
730
bonds_to_remove = []
731
for bond in universe.atoms.bonds:
732
i, j = bond.indices[0]
733
if dist_matrix[i, j] > 2.5: # Breaking threshold
734
bonds_to_remove.append((i, j))
735
736
if bonds_to_remove:
737
print(f"Frame {ts.frame}: Removing {len(bonds_to_remove)} broken bonds")
738
universe.delete_bonds(bonds_to_remove)
739
740
# Example usage (be careful with trajectory modification!)
741
# modify_topology_example(u)
742
```
743
744
### Topology Validation
745
746
```python { .api }
747
def validate_topology(universe):
748
"""
749
Validate topology consistency and completeness.
750
"""
751
atoms = universe.atoms
752
issues = []
753
754
# Check for missing critical attributes
755
required_attrs = ['names', 'resnames', 'segids']
756
for attr in required_attrs:
757
if not hasattr(atoms, attr):
758
issues.append(f"Missing required attribute: {attr}")
759
760
# Check for reasonable masses
761
if hasattr(atoms, 'masses'):
762
masses = atoms.masses
763
if np.any(masses <= 0):
764
issues.append("Found zero or negative masses")
765
if np.any(masses > 300): # Unusually heavy atoms
766
issues.append("Found unusually heavy atoms (>300 amu)")
767
768
# Check bond lengths
769
if hasattr(atoms, 'bonds'):
770
bond_lengths = atoms.bonds.length()
771
if np.any(bond_lengths < 0.5):
772
issues.append("Found unusually short bonds (<0.5 Å)")
773
if np.any(bond_lengths > 3.0):
774
issues.append("Found unusually long bonds (>3.0 Å)")
775
776
# Check for disconnected residues
777
if hasattr(atoms, 'bonds'):
778
for residue in atoms.residues:
779
res_atoms = residue.atoms
780
res_bonds = [bond for bond in atoms.bonds
781
if all(idx in res_atoms.indices for idx in bond.indices)]
782
if len(res_bonds) == 0 and len(res_atoms) > 1:
783
issues.append(f"Residue {residue.resname}{residue.resid} has no internal bonds")
784
785
# Report results
786
if issues:
787
print("Topology validation found issues:")
788
for issue in issues:
789
print(f" - {issue}")
790
else:
791
print("Topology validation passed - no issues found")
792
793
return len(issues) == 0
794
795
# Validate topology
796
is_valid = validate_topology(u)
797
```