82.193.241.6 at 17:15, 11 November 2024

2024-11-11T17:15:55Z

New page

{{Short description|Type of algorithm}}
The '''Kabsch algorithm''', also known as the '''Kabsch-Umeyama algorithm''',<ref>{{Cite journal |last1=Lawrence |first1=Jim |last2=Bernal |first2=Javier |last3=Witzgall |first3=Christoph |date=2019-10-09 |title=A Purely Algebraic Justification of the Kabsch-Umeyama Algorithm |url=https://nvlpubs.nist.gov/nistpubs/jres/124/jres.124.028.pdf |journal=Journal of Research of the National Institute of Standards and Technology |language=en |volume=124 |pages=124028 |doi=10.6028/jres.124.028 |issn=2165-7254 |pmc=7340555 |pmid=34877177}}</ref> named after Wolfgang Kabsch and Shinji Umeyama, is a method for calculating the optimal [[rotation matrix]] that minimizes the [[RMSD]] ([[root mean square]]d deviation) between two paired sets of points. It is useful for [[point-set registration]] in [[computer graphics]], and in [[cheminformatics]] and [[bioinformatics]] to compare molecular and [[protein]] structures (in particular, see [[root-mean-square deviation (bioinformatics)]]).

The algorithm only computes the rotation matrix, but it also requires the computation of a translation vector. When both the translation and rotation are actually performed, the algorithm is sometimes called partial [[Procrustes superimposition]] (see also [[orthogonal Procrustes problem]]).

== Description ==

Let {{mvar|P}} and {{mvar|Q}} be two sets, each containing {{mvar|N}} points in <math>\mathbb{R}^n</math>. We want to find the transformation from {{mvar|Q}} to {{mvar|P}}. For simplicity, we will consider the three-dimensional case (<math>n = 3</math>).
The sets {{mvar|P}} and {{mvar|Q}} can each be represented by {{math|''N'' × 3}} [[matrix (mathematics)|matrices]] with the first row containing the coordinates of the first point, the second row containing the coordinates of the second point, and so on, as shown in this matrix:

<math display="block">\begin{pmatrix}
x_1 & y_1 & z_1 \\
x_2 & y_2 & z_2 \\
\vdots & \vdots & \vdots \\
x_N & y_N & z_N \end{pmatrix}</math>

The algorithm works in three steps: a '''translation,''' the '''computation of a covariance matrix''', and '''the computation''' of the optimal rotation matrix.

=== Translation ===
Both sets of coordinates must be translated first, so that their [[centroid]] coincides with the origin of the [[coordinate system]]. This is done by subtracting the centroid coordinates from the point coordinates.

=== Computation of the covariance matrix ===
The second step consists of calculating a matrix {{mvar|H}}. In matrix notation,

:<math> H = P^\mathsf{T}Q \, </math>

or, using summation notation,

:<math> H_{ij} = \sum_{k = 1}^N P_{ki} Q_{kj}, </math>

which is a [[cross-covariance|cross-covariance matrix]] when {{mvar|P}} and {{mvar|Q}} are seen as [[design_matrix|data matrices]].

=== Computation of the optimal rotation matrix ===
It is possible to calculate the optimal rotation {{mvar|R}} based on the matrix formula

:<math> R = \left(H^\mathsf{T} H\right)^\frac12 H^{-1}, </math>

but implementing a numerical solution to this formula becomes complicated when all special cases are accounted for (for example, the case of {{mvar|H}} not having an inverse).

If [[singular value decomposition]] (SVD) routines are available the optimal rotation, {{mvar|R}}, can be calculated using the following algorithm.

First, calculate the SVD of the covariance matrix {{mvar|H}},

:<math> H = U \Sigma V^\mathsf{T} </math>

where {{mvar|U}} and {{mvar|V}} are orthogonal and <math>\Sigma</math> is diagonal. Next, record if the orthogonal matrices contain a reflection,

:<math> d = \det\left(U V^\mathsf{T}\right) = \det(U) \det(V).</math>

Finally, calculate our optimal rotation matrix {{mvar|R}} as

:<math> R = U \begin{pmatrix}
1 & 0 & 0 \\
0 & 1 & 0 \\
0 & 0 & d \end{pmatrix} V^\mathsf{T}. </math>

This {{mvar|R}} minimizes <math>\sum_{k = 1}^N|R q_k - p_k|</math>, where <math>q_k</math> and <math>p_k</math> are rows in {{mvar|Q}} and {{mvar|P}} respectively.

Alternatively, optimal rotation matrix can also be directly evaluated as [[quaternion]].<ref>{{Cite journal|last=Horn|first=Berthold K. P.|authorlink=Berthold K.P. Horn|date=1987-04-01|title=Closed-form solution of absolute orientation using unit quaternions|journal=Journal of the Optical Society of America A|language=EN|volume=4|issue=4|pages=629|doi=10.1364/josaa.4.000629|bibcode=1987JOSAA...4..629H|issn=1520-8532|citeseerx=10.1.1.68.7320|s2cid=11038004 }}</ref><ref>{{Cite journal|last=Kneller|first=Gerald R.|date=1991-05-01|title=Superposition of Molecular Structures using Quaternions|journal=Molecular Simulation|volume=7|issue=1–2|pages=113–119|doi=10.1080/08927029108022453|issn=0892-7022}}</ref><ref name="Coutsias2004">{{cite journal |last1=Coutsias |first1=E. A. |last2=Seok |first2=C. |last3=Dill |first3=K. A. | title = Using quaternions to calculate RMSD | journal = J. Comput. Chem. | volume = 25 | issue = 15 | pages = 1849–1857 | year = 2004 | pmid = 15376254 | doi = 10.1002/jcc.20110|s2cid=18224579 }}</ref><ref name="Petitjean1999">{{cite journal | last = Petitjean | first = M. | title = On the root mean square quantitative chirality and quantitative symmetry measures | journal = J. Math. Phys. | volume = 40 | issue = 9 | pages = 4587–4595 | year = 1999 | doi = 10.1063/1.532988| bibcode = 1999JMP....40.4587P | url = https://hal.archives-ouvertes.fr/hal-02122820/file/PMP.JMP_1999.pdf }}</ref> This alternative description has been used in the development of a rigorous method for removing rigid-body motions from [[molecular dynamics]] trajectories of flexible molecules.<ref>{{Cite journal|date=2011-08-24|title=Least constraint approach to the extraction of internal motions from molecular dynamics trajectories of flexible macromolecules|journal=J. Chem. Phys.|volume=135|issue=8|pages=084110|doi=10.1063/1.3626275|pmid=21895162|issn=0021-9606|last1=Chevrot|first1=Guillaume|last2=Calligari|first2=Paolo|last3=Hinsen|first3=Konrad|last4=Kneller|first4=Gerald R.|bibcode=2011JChPh.135h4110C}}</ref> In 2002 a generalization for the application to probability distributions (continuous or not) was also proposed.<ref name="Petitjean2002">{{cite journal | last = Petitjean | first = M. | title = Chiral mixtures | journal = J. Math. Phys. | volume = 43 | issue = 8 | pages = 4147–4157 | year = 2002 | doi = 10.1063/1.1484559| bibcode = 2002JMP....43.4147P | s2cid = 85454709 | url = https://hal.archives-ouvertes.fr/hal-02122882/file/PMP.JMP_2002.pdf }}</ref>

=== Generalizations ===

The algorithm was described for points in a three-dimensional space. The generalization to {{mvar|D}} dimensions is immediate.

== External links ==
This SVD algorithm is described in more detail at https://web.archive.org/web/20140225050055/http://cnx.org/content/m11608/latest/

A [[Matlab]] function is available at http://www.mathworks.com/matlabcentral/fileexchange/25746-kabsch-algorithm

A [https://github.com/oleg-alexandrov/projects/blob/master/eigen/Kabsch.cpp C++ implementation] (and unit test) using [[Eigen (C++ library)|Eigen]]

A [[Python (programming language)|Python]] script is available at https://github.com/charnley/rmsd. Another implementation can be found in
[https://docs.scipy.org/doc/scipy/reference/generated/scipy.spatial.transform.Rotation.align_vectors.html SciPy].

A free [[PyMol]] plugin easily implementing Kabsch is [https://www.pymolwiki.org/index.php/Kabsch]. (This previously linked to CEalign [https://wiki.pymol.org/index.php/Cealign], but this uses the Combinatorial Extension (CE) algorithm.) [[Visual Molecular Dynamics|VMD]] uses the Kabsch algorithm for its alignment.

The [[FoldX]] modeling toolsuite incorporates the Kabsch algorithm to measure RMSD between Wild Type and Mutated protein structures.

== See also ==

* [[Wahba's problem|Wahba's Problem]]
* [[Orthogonal Procrustes problem]]

== References ==
{{reflist}}

* {{cite journal|last=Kabsch|first=Wolfgang|date=1976|title=A solution for the best rotation to relate two sets of vectors|journal=Acta Crystallographica|volume=A32|issue=5|page=922|doi=10.1107/S0567739476001873|bibcode=1976AcCrA..32..922K}}
** With a correction in {{cite journal|last=Kabsch|first=Wolfgang|date=1978|title=A discussion of the solution for the best rotation to relate two sets of vectors|journal=Acta Crystallographica|volume=A34|issue=5|pages=827–828|doi=10.1107/S0567739478001680|bibcode=1978AcCrA..34..827K|doi-access=free}}
* {{cite conference|last1=Lin|first1=Ying-Hung|last2=Chang|first2=Hsun-Chang|last3=Lin|first3=Yaw-Ling|date=December 15–17, 2004|title=A Study on Tools and Algorithms for 3-D Protein Structures Alignment and Comparison|conference=International Computer Symposium|location=Taipei, Taiwan}}
* {{cite journal|last=Umeyama|first=Shinji|date=1991|title=Least-Squares Estimation of Transformation Parameters Between Two Point Patterns|journal=IEEE Trans. Pattern Anal. Mach. Intell.|volume=13|issue=4|pages=376–380|doi=10.1109/34.88573}}

[[Category:Bioinformatics algorithms]]

Kabsch algorithm - Revision history

82.193.241.6 at 17:15, 11 November 2024