[SOLVED] How to plot the free energy landscape of protein structure?

How to plot the free energy landscape of protein structure?

I understand the question is not appropriate for this platform, but I can try if I can get some hints,

I've been trying to plot the free energy landscape of a protein structure ("Chignolin"). I'm completely run out of ideas how to do that!! I've MD simulation trajectory file Trajectory file and using pyemma to plot the energy landscape. But I'm getting the error "" TypeError: plot_free_energy() takes from 2 to 20 positional arguments but 28 were given ""

Could someone figure out where the problem lies? Here is my code


import numpy as np
import matplotlib.pyplot as plt
import mdtraj as md
from itertools import combinations
from simtk.openmm.app.topology import Topology
from simtk.openmm.app.simulation import Simulation
from simtk.openmm.app.dcdreporter import DCDReporter
from simtk.openmm.app.statedatareporter import StateDataReporter
import simtk.unit as u
import simtk.openmm as mm
import simtk.openmm.openmm as openmm
import pyemma.coordinates as coor
import pyemma

pdb = md.load('1uao_Calpha.pdb')
feat = pyemma.coordinates.data.MDFeaturizer(pdb)
feat.add_distances_ca(periodic=False)
files = pyemma.coordinates.load('traj/DESRES/CLN025-0-c-alpha/CLN025-0-c-alpha-005.dcd', features = feat)

pyemma.plots.plot_free_energy(*files.T)
plt.show()

Here is the another pdb file.

Solution

I recommend you start reading the documentation, especially the "learn PyEMMA" section containing Jupyter notebooks teaching you the work-flow to extract properly weighted "pseudo" free-energy surfaces. Usually these surfaces are drawn into the dimensions of the first two slowest dynamical processes, but you can think of any other combination as well. These dimensions are defined by a TICA or VAMP projection, which are basically methods to extract the slow modes from your data, in case of proteins this contains folding and rare events.

As a primer I suggest reading this tutorial first, as it gives you a brief overview how to load and process your data to extract the slow modes. Note that this not yet contain Markov state modelling, so read further in the other examples to learn about that.