Merge two lipid bilayers

Dear all,

I have two separated lipid bilayers (two gro files), and I want to combine them together in a way they don’t mix (please see the picture). Is it possible using GROMACS? pic

My GROMACS version: 2020.3

Thank you in advance
Mehrnoosh

A brute force approach would be,

Say you want to join along x-axis
Step-1: Shift gro2 by box length (x length) of gro1 one by gmx trjconv

gmx trjconv -f gro2 -trans xlength1 0.0 0.0 -o new_gro2  

** Use slightly large x-length to avoid overlap

Step-2: Now you need to put the contents of two files together. Remove the last line of gro1, and first two lines of gro2. append gro2 into gro1 (e.g. cat gro1 new_gro2 > new_gro). And adjust total number of atoms

Step-3: Fix serial numbers with gmx editconf

gmx editconf -f new_gro -o new_gro_fixed 

adjust new box length in x-direction in the new_gro_fixed. Which will be xlength1 + xlength2 + offset

There might be a more smarter way to do this, above one is a dirty way to accomplish the job.

Cheers,
Masrul

Hi Masrul,

Thank you for your reply. Well, your method actually worked, but as you said in a dirty way :) I think I should look for a cleaner way to do this. But I really liked your idea. Thanks for sharing it.

Cheers,
Mehrnoosh

Hello Mehrnoosh,

Long back, I wrote a general purpose python script, which facilitates different transformation of coordinates. I think my script can be handy for your purpose too. I pushed that script to my github. You can use my script to do the same thing if you are comfortable with python. Download MolAligner.py and import it.

Example:

# merge.py 

from MolAligner import Aligner 

# load layers as object 
layer1= Aligner("layer1.gro")
layer2= Aligner("layer2.gro") 

# tranlate 
offset = 0.20 
transVector = [layer1.box[0]+ offset, 0.0, 0.0]
layer2.moveBy(transVector)  

#merge 
layer1.merge(layer2) 

# write 
lx = layer1.box[0] + layer2.box[0] + offset   # box length along x 
ly = layer1.box[1]
lz = layer1.box[2] 
boxLen = [lx, ly, lz]
 
layer1.write(outFile="merged.gro",boxLen=boxLen) 

Don’t edit MolAligner.py, but import it into a new python, per se merge.py. Adjust example script accordingly to meet your need. Script has some other functionalities, you might want to play with.

Let me know If you have any questions.

Regards,
Masrul

Hi Masrul,

Thank you very much for your help. I merged the bilayers with the first method you told me. It was a bit tricky but worked.
Thank you for the Python code, I am just learning Python, so not too familiar with the code :)

Best wishes,
Mehrnoosh

Hi,

Actually, now the only problem is that when I solvated the system, there are some holes between the two bilayers that water molecules can enter. If I reduce the distance between the two lipids, some atoms would overlap.
Minimization of the system didn’t help.

I appreciate any suggestion.

Best,
Mehrnoosh

How did you solvate the system? I mean which tool?

I use gmx solvate in GROMACS.
I also changed vdwradii.dat for C atoms to 0.5, but water molecules can still enter the holes between the two lipids.

One easy way is to select the index number of the water molecules accommodating out of the gap between the bilayers using gmx select and then gmx trjconv would give you the water molecules of your interest.

The steps are as follows:

  1. gmx grompp -f tst.mdp -c merged-solvated.gro -o select.tpr

  2. echo “resname SOL and x > x2 and x < x1” | gmx select -f merged-solvated.gro -s select.tpr -selrpos mol_cog -on sol.ndx
    (# x2 and x1 are the approximate borders of the gap along x direction)

  3. gmx trjconv -f merged-solvated.gro -s select.tpr -n sol.ndx -o sol.gro

The sol.gro contains the water molecules of your desire which you can cat(add) it to the initial merged.gro file.

Regards,
Salman

Hi Salman,

Thank you for your reply. Actually, I could solve the problem by a longer minimization and reducing dt.

Best wishes,
Mehrnoosh