Lipid-rafts are functional platforms in the plasma membrane that consist of $L_o$ lipids, and they act as sites for the aggregation of membrane-signalling proteins. Lipid-rafts are thought be very small, on the order of nanometers. Some pure lipid mixtures for so-called nanodomains - small regions of $L_o$ lipids in an otherwise $L_d$ bulk phase. These pure lipid mixtures are frequently used as approximations to lipid-rafts as they are much easier to study.
There is much interest in understanding under what conditions nanodomains in apposing leaflets are spatially aligned. Such alignment is known as interleaflet registration. The class lipyphilic.lib.registration.Registration
can be used to quantify the degree of interleaflet registration in a planar bilayer.
Registration
is an implementation of the registration analysis described by Thallmair et al. The degree of registration is calculated as the Pearson correlation coefficient of molecular densities in the upper and lower leaflets. First, the two-dimensional density of each leaflet is calculated:
where the $(x, y)$ positions of lipid atoms in leaflet $L$ are binned into two-dimensional histograms with bin lengths of 1 Å. $L$ is either the upper ($u$) or lower ($l$) leaflet. The two-dimensional density is then convolved with a circular Gaussian density of standard deviation $\sigma$.
The registration between the two leaflets, $r_{u/l}$, is then calculated as the Pearson correlation coefficient between $\rho(x, y)_{u}$ and $\rho(x, y)_{l}$. Values of $r_{u/l}=1$ correspond to perfectly registered domains and values of $r_{u/l}=-1$ correspond to perfectly anti-registered domains.
Below, we will look at the registration of $L_o$ lipids in a DPPC, DOPC, and cholesterol membrane simulated by Smith et al..
import pathlib
import pickle
import numpy as np
import matplotlib.pyplot as plt
import MDAnalysis as mda
from lipyphilic.lib.registration import Registration
u = mda.Universe("../datafiles/dppc-dopc-chol.tpr", "../datafiles/dppc-dopc-chol.xtc")
To calculate interleaflet registration, we first need to know which leaflet each lipid is in at each frame. We will use the results from the notebook on assigning lipids to leaflets.
with open("../results/leaflets/leaflets-dppc-dopc-chol.pkl", 'rb') as f:
leaflets = pickle.load(f)
To calculate interleaflet registration of $L_o$ lipids, we first need to know whether each lipid at each frame is $L_d$ or $L_o$.
In another tutorial, we construct a hidden Markov model based on lipid thicknesses in order to determine whether each lipid in each frame is liquid-ordered ($L_o$) or liquid-disordered ($L_d$). We will use this information here in order to quantify to local lipid environments.
# Load the lipid order data
filename = "../results/HMM/lipid-order.npy"
lipid_order = np.load(filename)
lipid_order
is a two-dimensional NumPy. Like many analyses in lipyphilic, the array is of shape ($N_{\rm lipids}$, $N_{\rm lipids}$) The:
In the array, the ordered state of each lipid is defined as follows:
Let's look at the array:
lipid_order
array([[ 1, 1, 1, ..., 1, 1, 1], [ 1, 1, 1, ..., -1, -1, -1], [ 1, 1, 1, ..., 1, 1, 1], ..., [ 1, 1, 1, ..., 1, 1, 1], [-1, -1, -1, ..., -1, -1, -1], [ 1, 1, 1, ..., 1, 1, 1]], dtype=int8)
We see the first lipid (first row) is $L_o$ both at the beginning and end of the frames used in the analysis.
The second lipid (second row) is $L_o$ at the beginning of the analysis but $L_d$ at the end.
To calculate the interleaflet registration, we will used the class Registration
.
We need to pass a lipid selection for each leaflet (upper_sel
and lower_sel
), along with the leaflet membership (leaflets
). As we wish to calculate the registration of $L_o$ lipids, we will pass a boolean mask to the filter_by
keyword of Registration
. A boolean mask is an array of True
and False
values, and in our case it will tell Registration
whether or not to consider each lipid.
registration = Registration(
universe=u,
upper_sel="name GL1 GL2 ROH", # select all lipids in the membrane
lower_sel="name GL1 GL2 ROH",
leaflets=leaflets.leaflets, # pass the leaflet information
filter_by=lipid_order == 1 # consider only L_o lipids
)
We then select which frames of the trajectory to analyse (None
will use every frame) and choose to display a progress bar (verbose=True
):
registration.run(
start=None,
stop=None,
step=None,
verbose=True
)
0%| | 0/51 [00:00<?, ?it/s]
<lipyphilic.lib.registration.Registration at 0x15ee49d10>
The results are available in the registration.registration
attribute as a NumPy array, in which:
# There is one value for each frame
registration.registration.size
51
And we see the domains are highly registered
registration.registration
array([0.97071996, 0.96829909, 0.9694673 , 0.96837016, 0.97048573, 0.97032585, 0.96955231, 0.96604412, 0.96699038, 0.96914026, 0.96991722, 0.97008281, 0.97089963, 0.97107136, 0.97198107, 0.97105928, 0.97273829, 0.97266528, 0.97164116, 0.96786753, 0.96917645, 0.97076293, 0.97269181, 0.97321275, 0.97250479, 0.97076373, 0.96878143, 0.97098904, 0.97179541, 0.97321968, 0.97401061, 0.97211395, 0.97215704, 0.97248677, 0.97265972, 0.97194608, 0.97395345, 0.97145309, 0.97340571, 0.97201339, 0.97229573, 0.97255324, 0.97186511, 0.97338536, 0.97579757, 0.97749242, 0.97452066, 0.97218495, 0.97316014, 0.97361326, 0.97548453])
This is because the frame we are analysing are taken from the end of a simulation. At this time, the phase-separation between $L_o$ and $L_d$ regions is complete and the domains are fully registered.
We can plot the results to check whether the domains are becoming more registered over time:
plt.plot(registration.times, registration.registration)
plt.ylabel("Registration")
plt.xlabel(r"Time ($\rm \mu s$)")
Text(0.5, 0, 'Time ($\\rm \\mu s$)')