| Title: | Occlusion Surface Using the Occluded Surface and Fibonacci Occluded Surface |
|---|---|
| Description: | The Occluded Surface (OS) algorithm is a widely used approach for analyzing atomic packing in biomolecules as described by Pattabiraman N, Ward KB, Fleming PJ (1995) <doi:10.1002/jmr.300080603>. Here, we introduce 'fibos', an 'R' and 'Python' package that extends the 'OS' methodology, as presented in Soares HHM, Romanelli JPR, Fleming PJ, da Silveira CH (2024) <doi:10.1101/2024.11.01.621530>. |
| Authors: | Herson Soares [cre, aut], Carlos Silveira [aut], João Romanelli [aut], Patrick Fleming [aut], Posit Software, PBC [cph, fnd] |
| Maintainer: | Herson Soares <[email protected]> |
| License: | GPL-3 |
| Version: | 2.0.1 |
| Built: | 2026-06-01 08:57:02 UTC |
| Source: | https://github.com/cran/fibos |
This function creates a 'Python' virtual environment and installs the 'Python' module 'fibos' required for the full functionality of this package. It handles different system configurations and ensures that the correct compiler paths are set.
fibos_config()fibos_config()
This function will install external software (a 'Python' package) on your system. Administrator/sudo privileges might be required on some systems.
# Set up the 'Python' environment and install the required module fibos_config()# Set up the 'Python' environment and install the required module fibos_config()
The 'get_radii' function is responsible for loading the atomic radii values used for surface-occlusion calculations. The values it returns are those currently employed in those calculations.
get_radii()get_radii()
A data frame containing the radii values.
Carlos Henrique da Silveira ([email protected])
Herson Hebert Mendes Soares ([email protected])
Joao Paulo Roquim Romanelli ([email protected])
Patrick Fleming ([email protected])
library(fibos) fibos_config() #Loads the radius values that have been configured for code execution. radii = get_radii() #Displays the first three lines. radii |> utils::head(3) |> print()library(fibos) fibos_config() #Loads the radius values that have been configured for code execution. radii = get_radii() #Displays the first three lines. radii |> utils::head(3) |> print()
The 'Occluded Surface (OS)' algorithm is a widely used approach for analyzing atomic packing in biomolecules. Here, we introduce 'FIBOS', an 'R' and 'Python' package that extends the 'OS' methodology with enhancements. The homonymous function 'occluded_surface' calculates 'OS' per atom.
occluded_surface(pdb, method = "FIBOS", density_dots = 5)occluded_surface(pdb, method = "FIBOS", density_dots = 5)
pdb |
4-digit PDB id (will fetch it from the RCSB repository) or the path to a PDB local file. |
method |
Method to be used: 'OS' (classic) or 'FIBOS'(default).The classic 'OS' covers the surface radially with one of the axes as a reference when allocating the dots. In 'FIBOS', Fibonacci spirals were used to allocate the dots, which is known to produce lower axial anisotropy as well as more evenly spaced points on a sphere. |
density_dots |
Distribution density of atomic dots for surface occlusion calculation. |
'Occluded Surface (OS)' (Pattabiraman et al. 1995) method distributes dots (representing patches of area) across the atom surfaces. Each dot has a normal that extends until it reaches either a van der Waals surface of a neighboring atom (the dot is considered occluded) or covers a distance greater than the diameter of a water molecule (the dot is considered non-occluded and disregarded). Thus, with the summed areas of dots and the lengths of normals, it is possible to compose robust metrics capable of inferring the average packing density of atoms, residues, proteins, as well as any other group of biomolecules.
For more details, see (Fleming et al, 2000) and (Soares, et al, 2024)
A table containing:
ATOMthe atomic contacts for each atom.
NUMBER OF POINTSthe number of dots (patches of area) on atomic surface.
AREAthe summed areas of dots.
RAYLENGTHthe average lengths of normals normalized by 2.8 (water diameter). So, raylen is a value between 0 and 1. A raylen close to 1 indicates worse packaging.
DISTANCEthe average distances of contacts in ().
Herson Soares, Joao Romanelli, Patrick Fleming, Carlos Silveira.
Fleming PJ, Richards FM (2000). "Protein packing: Dependence on protein size, secondary structure and amino acid composition." https://doi.org/10.1006/jmbi.2000.3750
Pattabiraman N, Ward KB, Fleming PJ (1995). "Occluded molecular surface: Analysis of protein packing." https://doi.org/10.1002/jmr.300080603
Soares HHM, Romanelli JPR, Fleming PJ, da Silveira CH (2024). "bioRxiv, 2024.11.01.621530." https://doi.org/10.1101/2024.11.01.621530
library(fibos) #Configure the environment fibos_config() # Calculate FIBOS per atom and create .srf files in fibos_files folder pdb_fibos <- occluded_surface("1ptx", method = "FIBOS", density_dots = 5.0) # Calculate OSP metric per residue from .srf file in fibos_files folder pdb_osp <- osp(fs::path("fibos_files","prot_1ptx.srf"))library(fibos) #Configure the environment fibos_config() # Calculate FIBOS per atom and create .srf files in fibos_files folder pdb_fibos <- occluded_surface("1ptx", method = "FIBOS", density_dots = 5.0) # Calculate OSP metric per residue from .srf file in fibos_files folder pdb_osp <- osp(fs::path("fibos_files","prot_1ptx.srf"))
Implements the 'occluded surface' packing density metric (OSP) averaged by residue, as described in (Fleming and Richards 2000).
osp(file)osp(file)
file |
a SRF File (.srf) generated by 'occluded_surface' in fibos_files folder. |
A table containing:
Resnumresidue id.
Resnameresidue name.
OSthe summed areas of dots in residue.
`os*[1-raylen]`'OS' areas weighted by (1-raylen). Raylen is the average lengths of normals normalized by 2.8 (water diameter). So, raylen is a value between 0 and 1. A raylen close to 1 indicates worse packaging, and the 'OS' will be reduced.
OSPaverage occluded surface packing value (OSP) by residue.
Herson Soares
Joao Romanelli
Patrick Fleming
Carlos Silveira.
Fleming PJ, Richards FM (2000). "Protein packing: Dependence on protein size, secondary structure and amino acid composition." https://doi.org/10.1006/jmbi.2000.3750
Pattabiraman N, Ward KB, Fleming PJ (1995). "Occluded molecular surface: Analysis of protein packing." https://doi.org/10.1002/jmr.300080603
Soares HHM, Romanelli JPR, Fleming PJ, da Silveira CH (2024). "bioRxiv, 2024.11.01.621530." https://doi.org/10.1101/2024.11.01.621530
library(fibos) #Configure the Environment fibos_config() # Calculate FIBOS per atom and create .srf files in fibos_files folder pdb_fibos <- occluded_surface("1ptx", method = "FIBOS", density_dots = 5.0) # Calculate OSP metric per residue from .srf file in fibos_files folder pdb_osp <- osp(fs::path("fibos_files","prot_1ptx.srf"))library(fibos) #Configure the Environment fibos_config() # Calculate FIBOS per atom and create .srf files in fibos_files folder pdb_fibos <- occluded_surface("1ptx", method = "FIBOS", density_dots = 5.0) # Calculate OSP metric per residue from .srf file in fibos_files folder pdb_osp <- osp(fs::path("fibos_files","prot_1ptx.srf"))
This function reloads the 'OS' default atomic radii values.
reset_radii()reset_radii()
Carlos Henrique da Silveira ([email protected])
Herson Hebert Mendes Soares ([email protected])
Joao Paulo Roquim Romanelli ([email protected])
Patrick Fleming ([email protected])
library(fibos) fibos_config() #Loads the radius values that have been configured for code execution. radii = get_radii() #Displays the first three lines. radii |> utils::head(3) |> print() #Modifies the value of a specific radius. radii$RAY[1] = 2.15 #Sets the radius value from a tibble. set_radii(radii) #Displays the first three lines. get_radii() |> utils::head(3) |> print() #Loads the default radius values. reset_radii() #Displays the first three lines. get_radii() |> utils::head(3) |> print()library(fibos) fibos_config() #Loads the radius values that have been configured for code execution. radii = get_radii() #Displays the first three lines. radii |> utils::head(3) |> print() #Modifies the value of a specific radius. radii$RAY[1] = 2.15 #Sets the radius value from a tibble. set_radii(radii) #Displays the first three lines. get_radii() |> utils::head(3) |> print() #Loads the default radius values. reset_radii() #Displays the first three lines. get_radii() |> utils::head(3) |> print()
This function enables modification of the radius values by passing a 'data.frame' as an argument.
set_radii(radii_values)set_radii(radii_values)
radii_values |
A 'data.frame' containing atomic radii values. |
Carlos Henrique da Silveira ([email protected])
Herson Hebert Mendes Soares ([email protected])
Joao Paulo Roquim Romanelli ([email protected])
Patrick Fleming ([email protected])
library(fibos) fibos_config() #Loads the radius values that have been configured for code execution. radii = get_radii() #Displays the first three lines. radii |> utils::head(3) |> print() #Modifies the value of a specific radius. radii$RAY[1] = 2.15 #Sets the radius value from a tibble. set_radii(radii) #Displays the first three lines. get_radii() |> utils::head(3) |> print()library(fibos) fibos_config() #Loads the radius values that have been configured for code execution. radii = get_radii() #Displays the first three lines. radii |> utils::head(3) |> print() #Modifies the value of a specific radius. radii$RAY[1] = 2.15 #Sets the radius value from a tibble. set_radii(radii) #Displays the first three lines. get_radii() |> utils::head(3) |> print()