0

I transformed the coordinates of two atoms into an array: coord

And I must calculate the root mean squared deviation (RMSD) between the two sets of these coordinates.

For this I have:

def cal_rmsd_numpy(coord_1, coord_2):
    rmsd = np.sqrt(((coord_1 - coord_2) ** 2).mean())    ## this would be the formula
    return rmsd

rmsd = cal_rmsd_numpy(coord_1, coord_2)

print(rmsd)

But the result does not give me the correct number. I think the error is in the formula. Could someone please help me correct what I have wrong?

The formule for RMSD is:

enter image description here

Angellys Correa
  • 51
  • 1
  • 6
  • what are the shapes of `coord_1` and `coord_2`? This is the correct formula if all of the `coord_1[i]` and `coord_2[i]` are scalars (so they have shape like `(n,)`), but if each element is itself a vector, as the formula suggests, you need to multiply by a factor of 3 inside the `sqrt` – AJ Biffl Nov 30 '21 at 02:42
  • Thank you so much. That was what was lacking, multiply by 3 – Angellys Correa Nov 30 '21 at 17:43
  • FYI, there are also more dependable/above-board methods for fixing this (that will also work for vectors of any dimension) that involve summing along a specific axis before taking the mean – AJ Biffl Dec 01 '21 at 01:59

3 Answers3

1

Solution found:

rmsd = np.sqrt(((((coordenadas_1 - coordenadas_2)** 2))*3).mean())
Angellys Correa
  • 51
  • 1
  • 6
0

You have to add up the subtraction of the coordinates firs:

rmsd = np.sqrt((((coordenadas_1 - coordenadas_2)**2).sum()).mean())

Or even more intuitive looking to the formula:

rmsd = np.sqrt((((coordenadas_1 - coordenadas_2)**2).sum())/len(coordenadas_1))

Chema_arguez
  • 57
  • 8
0
np.sqrt(((a - b)**2).sum(-1).mean())
Fairy Dance
  • 1
  • Thank you for contributing to the Stack Overflow community. This may be a correct answer, but it’d be really useful to provide additional explanation of your code so developers can understand your reasoning. This is especially useful for new developers who aren’t as familiar with the syntax or struggling to understand the concepts. **Would you kindly [edit] your answer to include additional details for the benefit of the community?** – Jeremy Caney Aug 16 '23 at 00:31