In [1]:
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
import datetime
sns.set(rc={'figure.figsize':(10,6.180)})
sns.set_style("whitegrid")
%matplotlib inline
In [2]:
plt.rcParams['figure.figsize'] = [16.18033, 10] #golden ratio
plt.rcParams['figure.facecolor'] = 'w'
plt.rcParams['figure.dpi'] = 100
In [31]:
pdb_list = ["2bg9", "1j4n", "1py6_SD", "2bl2", "1rhz", "1iwg", "2ic8", "1pv6", "1occ", "1kpl", "2bs2", "1py6", "1u19"]
simulationType = "openMM_membrane_structure_prediction"
folder = "test_2"
location = f"/Users/weilu/Research/server/jun_2019/{simulationType}/{folder}/energy_new_no_burial.dat"
new = pd.read_csv(location, sep="\s+").assign(Folder="new")
location = f"/Users/weilu/Research/server/jun_2019/{simulationType}/{folder}/energy_no_burial.dat"
original = pd.read_csv(location, sep="\s+").assign(Folder="original")
# d = pd.concat([new, original])
plt.scatter(new["Qc"], new["Contact"])
plt.scatter(original["Qc"], original["Contact"], c="red")
# original.plot.scatter("Qc", "Contact")
Out[31]:
In [32]:
pdb_list = ["2bg9", "1j4n", "1py6_SD", "2bl2", "1rhz", "1iwg", "2ic8", "1pv6", "1occ", "1kpl", "2bs2", "1py6", "1u19"]
simulationType = "openMM_membrane_structure_prediction"
folder = "test_1_1"
location = f"/Users/weilu/Research/server/jun_2019/{simulationType}/{folder}/energy_new_no_burial.dat"
new = pd.read_csv(location, sep="\s+").assign(Folder="new")
location = f"/Users/weilu/Research/server/jun_2019/{simulationType}/{folder}/energy_no_burial.dat"
original = pd.read_csv(location, sep="\s+").assign(Folder="original")
# d = pd.concat([new, original])
plt.scatter(new["Qc"], new["Contact"])
plt.scatter(original["Qc"], original["Contact"], c="red")
# original.plot.scatter("Qc", "Contact")
Out[32]:
In [30]:
Out[30]:
In [ ]:
pdb_list = ["2bg9", "1j4n", "1py6_SD", "2bl2", "1rhz", "1iwg", "2ic8", "1pv6", "1occ", "1kpl", "2bs2", "1py6", "1u19"]
simulationType = "openMM_membrane_structure_prediction"
folder = "test_2"
n_run = 1
all_data = []
for pdb in pdb_list:
for i in range(n_run):
for restart in range(2):
location = f"/Users/weilu/Research/server/jun_2019/{simulationType}/{folder}/{pdb}/{i}_{restart}/info.dat"
try:
tmp = pd.read_csv(location, sep="\s+")
tmp = tmp.assign(Run=i, Protein=pdb, Restart=restart)
all_data.append(tmp)
except:
print(pdb, i, restart)
pass
data = pd.concat(all_data)
today = datetime.today().strftime('%m-%d')
data.reset_index(drop=True).to_csv(f"/Users/weilu/Research/data/openMM/{simulationType}_{folder}_{today}.csv")
In [3]:
def my_reorder(a, first):
# move first to the top. and keep the rest
new_order = first.copy()
for col in a:
if col not in first:
new_order.append(col)
return new_order
def read_pdb(pre, name, run=30, rerun=2):
all_data = []
if run == -1:
run_list = ["native"]
else:
run_list = list(range(run))
for i in run_list:
if rerun == -1:
rerun_list = ["rerun"]
else:
rerun_list = list(range(rerun))
for j in rerun_list:
# pre = "/Users/weilu/Research/server/nov_2018/iterative_optimization_4/all_simulations/"
location = pre + f"{name}/simulation/{i}/{j}/"
try:
wham = pd.read_csv(location+"wham.dat")
except:
print(f"PDB: {name}, Run: {i}, Rerun: {j} not exist")
print(location+"wham.dat")
continue
wham.columns = wham.columns.str.strip()
remove_columns = ['Tc', 'Energy']
wham = wham.drop(remove_columns, axis=1)
energy = pd.read_csv(location+"energy.dat")
energy.columns = energy.columns.str.strip()
remove_columns = ['Steps', 'Shake', 'Excluded', 'Helix', 'AMH-Go', 'Vec_FM', 'SSB']
energy = energy.drop(remove_columns, axis=1)
data = pd.concat([wham, energy], axis=1).assign(Repeat=i, Run=j)
all_data.append(data)
data = pd.concat(all_data).reset_index(drop=True)
data = data.reindex(columns=my_reorder(data.columns, ["Steps", "Qw", "VTotal", "Run", "Repeat"]))
print(name, len(data))
return data
In [5]:
dataset = {"old":"1R69, 1UTG, 3ICB, 256BA, 4CPV, 1CCR, 2MHR, 1MBA, 2FHA".split(", "),
"new":"1FC2C, 1ENH, 2GB1, 2CRO, 1CTF, 4ICB".split(", "),
"test":["t089", "t120", "t251", "top7", "1ubq", "t0766", "t0778", "t0782", "t0792", "t0803", "t0815", "t0833", "t0842", "t0844"]}
dataset["combined"] = dataset["old"] + dataset["new"]
# folder = "first"
# steps = 30
def get_complete_data(pre, folder_list, pdb_list, formatName=True, **kwargs):
complete_all_data = []
for folder in folder_list:
# pre = "/Users/weilu/Research/server/april_2019/iterative_optimization_old_set/"
pre_folder = f"{pre}{folder}/"
all_data = []
for p in pdb_list:
if formatName:
name = p.lower()[:4]
else:
name = p
tmp = read_pdb(pre_folder, name, **kwargs)
all_data.append(tmp.assign(Name=name))
data = pd.concat(all_data)
complete_all_data.append(data.assign(Folder=folder))
data = pd.concat(complete_all_data)
data = data.reindex(columns=my_reorder(data.columns, ["Name", "Folder"]))
return data
In [ ]:
In [7]:
pre = "/Users/weilu/Research/server/jun_2019/membrane_potein_prediction/"
# folder_list = ["multi_iter0", "original"]
folder_list = ["first"]
# pdb_list = ['T0759-D1', 'T0953s2-D1', 'T0943-D1', 'T0773-D1', 'T0816-D1', 'T0854-D2', 'T0767-D1', 'T0853-D1', 'T0958-D1', 'T0834-D2', 'T0960-D3', 'T0862-D1', 'T0912-D3', 'T0898-D1', 'T0824-D1', 'T0782-D1', 'T0830-D2', 'T0761-D2', 'T0968s1-D1', 'T0870-D1', 'T0838-D1', 'T0803-D1']
pdb_list = ["2bg9", "1j4n", "1py6_SD", "2bl2", "1rhz", "1iwg", "2ic8", "1pv6", "1occ", "1kpl", "2bs2", "1py6", "1u19"]
# data = get_complete_data(pre, folder_list, pdb_list, run=30, rerun=-1, formatName=True)
data = get_complete_data(pre, folder_list, pdb_list, run=10, rerun=1, formatName=False)
data.Steps = data.Steps.astype(int)
data["Contact"] = data["Water"] + data["Burial"]
subset_data = data.query("Steps % 80000 == 0 and Steps != 0")
a = data
In [10]:
In [12]:
data = get_complete_data(pre, folder_list, pdb_list, run=-1, rerun=-1, formatName=False)
data.Steps = data.Steps.astype(int)
data["Contact"] = data["Water"] + data["Burial"]
native = data
In [14]:
native_energy = native.groupby(["Name", "Folder"]).head(1)
prediction_energy = data.groupby(["Name", "Repeat", "Folder"]).head(1)
native_energy["Contact"] = native_energy["Water"] + native_energy["Burial"]
In [17]:
native_energy
Out[17]:
In [15]:
y_show = "Contact"
g = sns.FacetGrid(a, col="Name",col_wrap=4, hue="Folder", sharey=False, sharex=False)
g = (g.map(plt.scatter, "Qw", y_show, alpha=0.5).add_legend())
# energy = native_energy.query("Name == 'T0759-D1' and Folder == 'multi_iter0_with_minimization'")["VTotal"][0]
# g.axes[0].axhline(energy, ls="--", color="blue", linewidth=4)
# energy = native_energy.query("Name == 'T0759-D1' and Folder == 'original_with_minimization'")["VTotal"][0]
# g.axes[0].axhline(energy, ls="--", color="orange", linewidth=4)
for ax in g.axes:
name= ax.title.get_text().split(" ")[-1]
# print(name)
energy = native_energy.query(f"Name == '{name}'")[y_show].iloc[0]
ax.axhline(energy, ls="--", color="blue", linewidth=4)
# energy = native_energy.query(f"Name == '{name}'")[y_show].iloc[1]
# ax.axhline(energy, ls="--", color="orange", linewidth=4)
In [ ]:
subset_data.tail()
In [ ]:
In [84]:
a = subset_data.query("Folder != 'multi_iter0_A_norm'")
In [85]:
y_show = "Contact"
g = sns.FacetGrid(a, col="Name",col_wrap=4, hue="Folder", sharey=False, sharex=False)
g = (g.map(plt.scatter, "Qw", y_show, alpha=0.5).add_legend())
# energy = native_energy.query("Name == 'T0759-D1' and Folder == 'multi_iter0_with_minimization'")["VTotal"][0]
# g.axes[0].axhline(energy, ls="--", color="blue", linewidth=4)
# energy = native_energy.query("Name == 'T0759-D1' and Folder == 'original_with_minimization'")["VTotal"][0]
# g.axes[0].axhline(energy, ls="--", color="orange", linewidth=4)
for ax in g.axes:
name= ax.title.get_text().split(" ")[-1]
# print(name)
energy = native_energy.query(f"Name == '{name}'")[y_show].iloc[0]
ax.axhline(energy, ls="--", color="blue", linewidth=4)
energy = native_energy.query(f"Name == '{name}'")[y_show].iloc[1]
ax.axhline(energy, ls="--", color="orange", linewidth=4)
In [40]:
data.to_csv("/Users/weilu/Research/data/optimization/may08.csv")
In [41]:
last_frame = data.groupby(["Name", "Repeat", "Folder"]).tail(1)
In [55]:
def do(cmd):
return subprocess.Popen(cmd, shell=True).wait()
In [51]:
import subprocess
In [ ]:
In [57]:
for i, line in last_frame.iterrows():
frame = int(int(line["Steps"]) / 4000)
folder = line["Folder"]
name = line["Name"]
repeat = int(line["Repeat"])
source = f"/Users/weilu/Research/server/may_2019/database/{folder}_{name}_{repeat}/frame{frame}.pdb"
pre = "/Users/weilu/Research/server/may_2019/single_memory/chosen_structures/"
do(f"mkdir -p {pre}{name}")
target = f"{pre}{name}/{folder}_{repeat}.pdb"
a = do(f"cp {source} {target}")
# print(a)
# print(frame, line)
# break
In [65]:
In [61]:
name_list = [a.lower()[:4] for a in pdb_list]
In [69]:
name = name_list[0]
for name in name_list:
out_file = f'''
load cleaned_pdbs/{name}.pdb
load {name}/original_0.pdb
load {name}/multi_iter0_0.pdb
alignto {name},
orient
util.color_deep("gray80", '{name}', 0)
util.color_deep("green", 'original_0', 0)
util.color_deep("red", 'multi_iter0_0', 0)
'''
with open(f"/Users/weilu/Research/server/may_2019/single_memory/chosen_structures/{name}.pml", "w") as out:
out.write(out_file)
In [25]:
subset_data = subset_data.query("Steps % 80000 == 0 and Steps != 0")
In [30]:
subset_data.head()
Out[30]:
In [31]:
subset_data["Contact"] = subset_data["Water"] + subset_data["Burial"]
In [21]:
data = get_complete_data(pre, folder_list, pdb_list, run=-1, rerun=-1, formatName=True)
data.Steps = data.Steps.astype(int)
data["Contact"] = data["Water"] + data["Burial"]
native = data
In [7]:
subset_data.head()
Out[7]:
In [23]:
native_energy = native.groupby(["Name", "Folder"]).head(1)
prediction_energy = data.groupby(["Name", "Repeat", "Folder"]).head(1)
native_energy["Contact"] = native_energy["Water"] + native_energy["Burial"]
In [11]:
In [70]:
y_show = "VTotal"
g = sns.FacetGrid(subset_data, col="Name",col_wrap=4, hue="Folder", sharey=False, sharex=False)
g = (g.map(plt.scatter, "Qw", y_show, alpha=0.5).add_legend())
# energy = native_energy.query("Name == 'T0759-D1' and Folder == 'multi_iter0_with_minimization'")["VTotal"][0]
# g.axes[0].axhline(energy, ls="--", color="blue", linewidth=4)
# energy = native_energy.query("Name == 'T0759-D1' and Folder == 'original_with_minimization'")["VTotal"][0]
# g.axes[0].axhline(energy, ls="--", color="orange", linewidth=4)
for ax in g.axes:
name= ax.title.get_text().split(" ")[-1]
# print(name)
energy = native_energy.query(f"Name == '{name}'")[y_show].iloc[0]
ax.axhline(energy, ls="--", color="blue", linewidth=4)
energy = native_energy.query(f"Name == '{name}'")[y_show].iloc[1]
ax.axhline(energy, ls="--", color="orange", linewidth=4)
In [35]:
y_show = "Contact"
g = sns.FacetGrid(subset_data, col="Name",col_wrap=4, hue="Folder", sharey=False, sharex=False)
g = (g.map(plt.scatter, "Qw", y_show, alpha=0.5).add_legend())
# energy = native_energy.query("Name == 'T0759-D1' and Folder == 'multi_iter0_with_minimization'")["VTotal"][0]
# g.axes[0].axhline(energy, ls="--", color="blue", linewidth=4)
# energy = native_energy.query("Name == 'T0759-D1' and Folder == 'original_with_minimization'")["VTotal"][0]
# g.axes[0].axhline(energy, ls="--", color="orange", linewidth=4)
for ax in g.axes:
name= ax.title.get_text().split(" ")[-1]
# print(name)
energy = native_energy.query(f"Name == '{name}'")[y_show].iloc[0]
ax.axhline(energy, ls="--", color="blue", linewidth=4)
energy = native_energy.query(f"Name == '{name}'")[y_show].iloc[1]
ax.axhline(energy, ls="--", color="orange", linewidth=4)
In [29]:
y_show = "Water"
g = sns.FacetGrid(subset_data, col="Name",col_wrap=4, hue="Folder", sharey=False, sharex=False)
g = (g.map(plt.scatter, "Qw", y_show, alpha=0.5).add_legend())
# energy = native_energy.query("Name == 'T0759-D1' and Folder == 'multi_iter0_with_minimization'")["VTotal"][0]
# g.axes[0].axhline(energy, ls="--", color="blue", linewidth=4)
# energy = native_energy.query("Name == 'T0759-D1' and Folder == 'original_with_minimization'")["VTotal"][0]
# g.axes[0].axhline(energy, ls="--", color="orange", linewidth=4)
for ax in g.axes:
name= ax.title.get_text().split(" ")[-1]
# print(name)
energy = native_energy.query(f"Name == '{name}'")[y_show].iloc[0]
ax.axhline(energy, ls="--", color="blue", linewidth=4)
energy = native_energy.query(f"Name == '{name}'")[y_show].iloc[1]
ax.axhline(energy, ls="--", color="orange", linewidth=4)
In [36]:
y_show = "Frag_Mem"
g = sns.FacetGrid(subset_data, col="Name",col_wrap=4, hue="Folder", sharey=False, sharex=False)
g = (g.map(plt.scatter, "Qw", y_show, alpha=0.5).add_legend())
# energy = native_energy.query("Name == 'T0759-D1' and Folder == 'multi_iter0_with_minimization'")["VTotal"][0]
# g.axes[0].axhline(energy, ls="--", color="blue", linewidth=4)
# energy = native_energy.query("Name == 'T0759-D1' and Folder == 'original_with_minimization'")["VTotal"][0]
# g.axes[0].axhline(energy, ls="--", color="orange", linewidth=4)
for ax in g.axes:
name= ax.title.get_text().split(" ")[-1]
# print(name)
energy = native_energy.query(f"Name == '{name}'")[y_show].iloc[0]
ax.axhline(energy, ls="--", color="blue", linewidth=4)
energy = native_energy.query(f"Name == '{name}'")[y_show].iloc[1]
ax.axhline(energy, ls="--", color="orange", linewidth=4)
In [ ]: