Commit 59879f11 authored by Ninjani's avatar Ninjani
Browse files

fixed prody superposition

parent ef979c90
......@@ -3,6 +3,7 @@ import typing
from dataclasses import dataclass, field
from pathlib import Path
import matplotlib.pyplot as plt
import numba as nb
import numpy as np
import prody as pd
......@@ -547,22 +548,33 @@ class StructureMultiple:
if not output_pdb_folder.exists():
output_pdb_folder.mkdir()
reference_name = self.structures[0].name
reference_pdb = pd.parsePDB(self.structures[0].pdb_file)
reference_pdb = pd.parsePDB(str(self.structures[0].pdb_file))
core_indices = np.array([i for i in range(len(alignments[reference_name])) if '-' not in [alignments[n][i] for n in alignments]])
aln_ref = helper.aligned_string_to_array(alignments[reference_name])
ref_coords_core = reference_pdb[helper.get_alpha_indices(reference_pdb)].getCoords().astype(np.float64)[
np.array([aln_ref[c] for c in core_indices])]
ref_centroid = rmsd_calculations.nb_mean_axis_0(ref_coords_core)
ref_coords_core -= ref_centroid
transformation = pd.Transformation(np.eye(3), -ref_centroid)
reference_pdb = pd.applyTransformation(transformation, reference_pdb)
# reference_pdb.setCoords(reference_pdb.getCoords().astype(np.float64) - ref_centroid)
pd.writePDB(str(output_pdb_folder / f"{reference_name}.pdb"), reference_pdb)
pdb_ref_2 = pd.parsePDB(str(output_pdb_folder / f"{reference_name}.pdb"))
plt.scatter(pdb_ref_2.getCoords()[:, 0], pdb_ref_2.getCoords()[:, 1], alpha=0.3)
# plt.show()
for i in range(1, len(self.structures)):
name = self.structures[i].name
pdb = pd.parsePDB(self.structures[i].pdb_file)
pdb = pd.parsePDB(str(self.structures[i].pdb_file))
aln_name = helper.aligned_string_to_array(alignments[name])
common_coords_2 = pdb[helper.get_alpha_indices(pdb)].getCoords().astype(np.float64)[np.array([aln_name[c] for c in core_indices])]
rotation_matrix, translation_matrix = rmsd_calculations.svd_superimpose(ref_coords_core, common_coords_2)
transformation = pd.Transformation(rotation_matrix, translation_matrix)
# coords = rmsd_calculations.apply_rotran(pdb.getCoords().astype(np.float64), rotation_matrix, translation_matrix)
# pdb.setCoords(coords)
transformation = pd.Transformation(rotation_matrix.T, translation_matrix)
pdb = pd.applyTransformation(transformation, pdb)
pd.writePDB(str(output_pdb_folder / f"{name}.pdb"), pdb)
pdb_2 = pd.parsePDB(str(output_pdb_folder / f"{name}.pdb"))
plt.scatter(pdb_2.getCoords()[:, 0], pdb_2.getCoords()[:, 1], alpha=0.3)
# plt.show()
......@@ -10,3 +10,4 @@ setup(name='caretta',
"dash-core-components==1.2.1", "dash-html-components==1.0.1", "dash-renderer==1.1.0",
"dash-table==4.3.0", "plotly==3.7.1", "flask"]
})
q
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment