A flip-flop event occurs when a molecule - typically a sterol - moves from one leaflet of a bilayer into the opposing leaflet.
We will identify flip-flop events in a ternary mixture of DPPC, DOPC, and Cholesterol simulated by Smith et al..
import pathlib
import pickle
import numpy as np
import MDAnalysis as mda
from lipyphilic.lib.flip_flop import FlipFlop
u = mda.Universe("../datafiles/dppc-dopc-chol.tpr", "../datafiles/dppc-dopc-chol.xtc")
To identify flip-flop events, 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)
The class lipyphilic.lib.flip_flop.FlipFlop
finds the frames at which a flip-flop event begins and ends, as well as the direction of travel (upper-to-lower or lower-to-upper).
flip_flop = FlipFlop(
universe=u,
lipid_sel="resname CHOL and name ROH", # only find flip-flop events for cholesterol
leaflets=leaflets.filter_leaflets("resname CHOL and name ROH") # pass only leaflet information on cholesterol
)
Select the frames to use in the analysis (None
will use every frame):
flip_flop.run(
start=None,
stop=None,
step=None
)
0%| | 0/4000 [00:00<?, ?it/s]
<lipyphilic.lib.flip_flop.FlipFlop at 0x11189b5d0>
The results are then available in the flipflop.flip_flops
attribute as a NumPy array.
Each row corresponds to an individual flip-flop event. In each row, the four columns correspond to the:
# We see there were 15166 flip-flop events in total
flip_flop.flip_flops.shape
(15166, 4)
Let's look at the first flip-flop event:
first_event = flip_flop.flip_flops[0]
first_event
array([ 1, 11, 13, -1])
In the above case, cholesterol with residue index 1 left it's original leaflet at frame 11, and entered its new leaflet at frame 39. This new leaflet is the lower leaflet (as the fourth column is equal to -1).
With this information, it is possible to, for example, determine the local lipid environemnt immediately before and after the flip-flop event occured. This is useful to know as Gu et al. showed that translocation is highly influenced by the local lipid environment of cholesterol.
Due to thermal fluctuations, cholesterol may briefly go to the midplane before returning to its original lealfet. It may also briefly leave the midplane before returning to the hydrophobic core. As such, it is useful to be able to ignore these small fluctuations.
With FlipFlop
, we can specify the minumum number of frames a molecule must reside in its new leaflet for the flip-flop to be considered successful. We do this using the frame_cutoff
keyword:
flip_flop = FlipFlop(
universe=u,
lipid_sel="name ROH",
leaflets=leaflets.filter_leaflets("name ROH"),
frame_cutoff=2,
)
With frame_cutoff=2, a molecule must remain in its new leaflet for at least 2 consecutive frames for the flip-flop to be considered successful. If this condition is not, the flip-flop event is ignored.
flip_flop.run(
start=None,
stop=None,
step=None
)
0%| | 0/4000 [00:00<?, ?it/s]
<lipyphilic.lib.flip_flop.FlipFlop at 0x111986b10>
Because we have ignored any flip-flop events where cholesterol residues in its new leaflet for one frame only, we see fewer flip-flop events in total:
# We see there are now 4665 flip-flop events in total (compared to 15166 previously)
flip_flop.flip_flops.shape
(4665, 4)
We will store the results in the following directory:
results_directory = pathlib.Path(
"../results/flip-flop"
)
# Create the directory if it doesn't already exist
results_directory.resolve().mkdir(exist_ok=True, parents=True)
# Location to store the results
filename = results_directory.joinpath("flip-flop-dppc-dopc-chol.pkl")
print(filename)
../results/flip-flop/flip-flop-dppc-dopc-chol.pkl
# store the object
with open(filename, 'wb') as f:
pickle.dump(leaflets, f)
FlipFlop
returns information on whether each event was successful (the molecule resides in the opposing leaflet for at least a given length of time), or not (the molecule went to the midplane but returned to its original leaflet).
This information is stored as a NumPy array in the flip_flop.flip_flop_success
attribute:
flip_flop.flip_flop_success
array(['Success', 'Success', 'Success', ..., 'Success', 'Fail', 'Fail'], dtype='<U7')
There is one value ("Success" or "Fail") for each flip-flop event:
flip_flop.flip_flop_success.size
4665
The number of successful flip-flops can be used to calculate the rate of cholesterol translocation. The flip-flop rate, $k$, is given by the number of successful flip-flop events, $N_{\rm Success}$, divided by the product of the number of cholesterol molecles, ($N_{\rm Chol}$), and the total simulation time, $t_{\rm Seconds}$, used in the analysis:
$$ k = \frac{N_{\rm Success}}{N_{\rm Chol} t_{\rm Seconds}} $$cholesterol = u.select_atoms("resname CHOL")
n_cholesterol = cholesterol.n_residues
ps_to_seconds = 1e-12 # convert ps to seconds
total_time = flip_flop.n_frames * flip_flop.step * u.trajectory.dt * ps_to_seconds
number_successful = np.sum(flip_flop.flip_flop_success == "Success")
flip_flop_rate = number_successful / (n_cholesterol * total_time) # per second
# Print the rate per second in scientific notation
print(f"{flip_flop_rate:e}")
2.753922e+06
See Baral et al. for further discussion on flip-flop in lipid bilayers, including the affect on the flip-flop rate of the buffer size used to assign molecules to the midplane of the bilayer.