SASA Data
GitHub Link to Code.
SASA feature type implementation for molecular dynamics analysis.
SASA feature type implementation for computing Solvent Accessible Surface Area using the Shrake-Rupley algorithm with support for both residue and atom level.
- class mdxplain.feature.feature_type.sasa.sasa.SASA(mode: str = 'residue', probe_radius: float = 0.14)
SASA feature type for computing Solvent Accessible Surface Area.
Computes solvent accessible surface area using the Shrake-Rupley algorithm. Supports both residue-level and atom-level calculations with configurable probe radius. The algorithm determines the surface area accessible to a spherical probe of specified radius.
This is a base feature type with no dependencies that provides insights into protein-solvent interactions and buried surface areas.
Uses mdtraj for SASA calculation under the hood.
Examples
>>> # Residue-level SASA with default water probe >>> sasa = SASA(mode='residue') >>> pipeline.feature.add_feature(sasa)
>>> # Atom-level SASA with custom probe radius >>> sasa = SASA(mode='atom', probe_radius=0.12) >>> pipeline.feature.add_feature(sasa)
>>> # Large probe for detecting cavities >>> sasa = SASA(mode='residue', probe_radius=0.20) >>> pipeline.feature.add_feature(sasa)
- __init__(mode: str = 'residue', probe_radius: float = 0.14) None
Initialize SASA feature type with calculation parameters.
Parameters
- modestr, default=’residue’
Level of SASA calculation:
‘residue’: SASA per residue (sum of constituent atoms)
‘atom’: SASA per individual atom
- probe_radiusfloat, default=0.14
Probe sphere radius in nanometers. Common values:
0.14 nm: Water molecule (default, most common)
0.12 nm: Small probe for tight cavities
0.20 nm: Large probe for detecting large cavities
Returns
None
Examples
>>> # Residue-level SASA (default) >>> sasa = SASA()
>>> # Atom-level SASA >>> sasa = SASA(mode='atom')
>>> # Custom probe radius for cavity detection >>> sasa = SASA(mode='residue', probe_radius=0.20)
>>> # Small probe for detailed surface analysis >>> sasa = SASA(mode='residue', probe_radius=0.12)
- init_calculator(use_memmap: bool = False, cache_path: str = './cache', chunk_size: int = 2000) None
Initialize the SASA calculator with specified configuration.
Parameters
- use_memmapbool, default=False
Whether to use memory mapping for large datasets
- cache_pathstr, optional
Directory path for storing cache files when using memory mapping
- chunk_sizeint, optional
Number of frames to process per chunk for memory-efficient processing
Returns
None
Examples
>>> # Basic initialization >>> sasa.init_calculator()
>>> # With memory mapping for large datasets >>> sasa.init_calculator(use_memmap=True, cache_path='./cache/')
>>> # With custom chunk size >>> sasa.init_calculator(chunk_size=1000)
- compute(input_data: Trajectory, feature_metadata: Dict[str, Any]) Tuple[ndarray, Dict[str, Any]]
Compute SASA from molecular dynamics trajectory.
Parameters
- input_datamdtraj.Trajectory
MD trajectory to compute SASA from
- feature_metadatadict
Residue metadata (passed through for residue-level analysis)
Returns
- tuple[numpy.ndarray, dict]
Tuple containing (sasa_array, feature_metadata) where sasa_array has shape (n_frames, n_residues) or (n_frames, n_atoms) depending on mode, with SASA values in Ų
Examples
>>> # Compute residue-level SASA >>> sasa = SASA(mode='residue') >>> sasa.init_calculator() >>> data, metadata = sasa.compute(trajectory, res_metadata) >>> print(f"SASA array shape: {data.shape}") >>> print(f"Mean SASA per residue: {data.mean(axis=0)}")
>>> # Compute atom-level SASA with memory mapping >>> sasa = SASA(mode='atom') >>> sasa.init_calculator(use_memmap=True) >>> data, metadata = sasa.compute(trajectory, res_metadata)
Raises
- ValueError
If calculator is not initialized
- get_dependencies() List[str]
Get list of feature type dependencies for SASA calculations.
Parameters
None
Returns
- List[str]
Empty list as SASA is a base feature with no dependencies
Examples
>>> sasa = SASA() >>> print(sasa.get_dependencies()) []