Area per lipid — lipyphilic.lib.area_per_lipid
¶
- Author
Paul Smith
- Year
2021
- Copyright
GNU Public License v2
This module provides methods for calculating the area per lipid in a bilayer.
The class lipyphilic.lib.area_per_lipid.AreaPerLipid
calculates the
area of each lipid via a 2D Voronoi tessellation of atomic positions. See
Lukat et al. (2013) for
a description of calculating the area per lipid via Voronoi tessellations.
This class uses Freud for performing the Voronoi tessellations from which the area per lipid is calculated.
Input¶
- Required:
universe : an MDAnalysis Universe object.
lipid_sel : atom selection for lipids in the bilayer. These atoms will be used to perform the Voronoi tessellation.
leaflets : leaflet membership (-1: lower leaflet, 0: midplane, 1: upper leaflet) of each lipid in the membrane.
Output¶
area : area per lipid of each lipid as each frame
Area data are returned in a numpy.ndarray
, where each row corresponds
to an individual lipid and each column corresponds to an individual frame, i.e.
areas[i, j] refers to the area of lipid i at frame j. The results are
accessible via the AreaPerLipid.areas
attribute.
Note
No area can be calculated for molecules that are in the midplane, i.e. those for which leaflets==0. These molecules will have NaN values in the results array for the frames at which they are in the midplane.
Example usage of AreaPerLipid
¶
An MDAnalysis Universe must first be created before using AreaPerLipid:
import MDAnalysis as mda
from lipyphilic.lib.assign_leaflets import AssignLeaflets
u = mda.Universe(tpr, trajectory)
Then we need to know which leaflet each lipid is in at each frame. This may be done using
lipyphilic.lib.assign_leaflets.AssignLeaflets
:
leaflets = AssignLeaflets(
universe=u,
lipid_sel="name GL1 GL2 ROH" # assuming we are using the MARTINI forcefield
)
leaflets.run()
The leaflet data are stored in the leaflets.leaflets
attribute. We can now create our
AreaPerLipid object:
areas = AreaPerLipid(
universe=u,
lipid_sel="name GL1 GL2 ROH",
leaflets=leaflets.leaflets
)
The above will use GL1 and GL2 beads to calculate the area of each phospholipid, and the ROH bead to calculate the area of each sterol. Two Voronoi tessellations will be performed at each frame — one for the upper leaflet and one for the lower leaflet.
We then select which frames of the trajectory to analyse (None will use every frame) and choose to display a progress bar (verbose=True):
areas.run(
start=None,
stop=None,
step=None,
verbose=True
)
The results are then available in the areas.areas
attribute as a
numpy.ndarray
. Each row corresponds to an individual lipid and each column
to an individual frame, i.e areas.areas[i, j] contains the area of lipid i
at frame j.
The class and its methods¶
- class lipyphilic.lib.area_per_lipid.AreaPerLipid(universe, lipid_sel, leaflets)[source]¶
Calculate the area of lipids in each leaflet of a bilayer.
Set up parameters for calculating areas.
- Parameters
universe (Universe) – MDAnalysis Universe object
lipid_sel (str) – Selection string for lipids in the bilayer. Typically, in all-atom simulations, one atom per sterol and three atoms per non-sterol lipid would be used. In coarse-grained simulations, one bead per sterol and two beads per non-sterol lipid would typically be used.
leaflets (numpy.ndarray (n_lipids,)) – An array of leaflet membership in which: -1 corresponds to the lower leaflet; 1 corresponds to the upper leaflet; and 0 corresponds to the midplane. If the array is 1D and of shape (n_lipids), each lipid is taken to remain in the same leaflet over the trajectory. If the array is 2D and of shape (n_lipids, n_frames), the leaflet to which each lipid is assisgned at each frame will be taken into account when calculating the area per lipid.
Tip
Leaflet membership can be determined using
lipyphilic.lib.assign_leaflets.AssignLeaflets
.- project_area(lipid_sel=None, start=None, stop=None, step=None, filter_by=None, bins=None, ax=None, cmap=None, vmin=None, vmax=None, cbar=True, cbar_kws={}, imshow_kws={})[source]¶
Project the area per lipid onto the xy plane of the membrane.
The areas per lipid, averaged over a selected range of frames, are projected onto the xy plane based on the center of mass of each lipid. The atoms to be used in calculating the center of mass of the lipids can be specified using the lipid_sel arugment.
This method creates an instance of lipyphilic.lib.plotting.ProjectionPlot with the projected areas interpolated across periodic boundaries. The plot is returned so further modification can be performed if needed.
Note
The lipid positions are taken from the middle frame of the selected range.
- Parameters
lipid_sel (MDAnalysis atom selection, optional) – The center of mass of each lipid will be determined via this selection. The default is None, in which case every atom of a lipid is used to determine its center of mass.
start (int, optional) – Start frame for averaging the area per lipid results.
stop (int, optional) – Final frame for averaging the area per lipid results.
step (int, optional) – Number of frames to skip
filter_by (array-like, optional) – A Boolean mask for selecting a subset of lipids. It may take the following shapes:
(n_lipids)
The mask is used to select a subset of lipids for projecting the areas onto the membrane plane.(n_lipids, n_frames)
This is the same shape as the NumPy array created by the lipyphilic.lib.AreaPerLipid.run() method. Boolean values are used only from the column corresponding to the middle frame of the range selected by start, stop, and step.The default is None, in which case no filtering is applied.
bins (int or array_like or [int, int] or [array, array]) – The bin specification:
int
If int, the number of bins for the two dimensions (nx=ny=bins).
array-like
If array_like, the bin edges for the two dimensions (x_edges=y_edges=bins).
[int, int]
If [int, int], the number of bins in each dimension (nx, ny = bins).
[array, array]
If [array, array], the bin edges in each dimension (x_edges, y_edges = bins).
combination
A combination [int, array] or [array, int], where int is the number of bins and array is the bin edges.
The default is None, in which case a grid with 1 x 1 Angstrom resolution is created.
ax (Axes, optional) – Matplotlib Axes on which to plot the projection. The default is None, in which case a new figure and axes will be created.
cmap (str or ~matplotlib.colors.Colormap, optional) – The Colormap instance or registered colormap name used to map scalar data to colors.
vmin, vmax (float, optional) – Define the data range that the colormap covers. By default, the colormap covers the complete value range of the supplied data.
cbar (bool, optional) – Whether or not to add a colorbar to the plot.
cbar_kws (dict, optional) – A dictionary of keyword options to pass to matplotlib.pyplot.colorbar.
imshow_kws (dict, optional) – A dictionary of keyword options to pass to matplotlib.pyplot.imshow, which is used to plot the 2D density map.
- Returns
area_projection (ProjectionPlot) – The ProjectionPlot object containing the area per lipid data and the matplotlob.pyplot.imshow plot of the projection.