title: Algorithm for superimposition of 3d points
body: <p>I need to superimpose two groups of 3D points on top of each other; i.e. find rotation and translation matrices to minimize the RMSD (root mean square deviation) between their coordinates.</p>
  
  <p>I currently use <a href="http://en.wikipedia.org/wiki/Kabsch_algorithm" rel="nofollow">Kabsch algorithm</a>, which is not very useful for many of the cases I need to deal with. Kabsch requires equal number of points in both data sets, plus, it needs to know which point is going to be aligned with which one beforehand. For my case, the number of points will be different, and I don't care which point corresponds to which in the final alignment, as long as the RMSD is minimized.</p>
  
  <p>So, the algorithm will (presumably) find a 1-1 mapping between the subsets of two point sets such that AFTER rotation&amp;translation, the RMSD is minimized.</p>
  
  <p>I know some algorithms that deal with different number of points, however they all are protein-based, that is, they try to align the backbones together (some continuous segment is aligned with another continuous segment etc), which is not useful for points floating in space, without any connections. (OK, to be clear, some points are connected; but there are points without any connections which I don't want to ignore during superimposition.)</p>
  
  <p>Only algorithm that I found is <a href="http://www.bioinformatics.org/strap/superposition.html" rel="nofollow">DIP-OVL</a>, found in STRAP software module (open source). I tried the code, but the behaviour seems erratic; sometimes it finds good alignments, sometimes it can't align a set of few points with itself after a simple X translation. </p>
  
  <p>Anyone know of an algorithm that deals with such limitations? I'll have at most ~10^2 to ~10^3 points if the performance is an issue.</p>
  
  <hr>
  
  <p>To be honest, the objective function to use is not very clear. RMSD is defined as the RMS of the distance between the <em>corresponding points</em>. If I have two sets with 50 and 100 points, and the algorithm matches 1 or few points within the sets, the resulting RMSD between those few points will be zero, while the overall superposition may not be so great. RMSD between <em>all</em> pairs of points is not a better solution (I think). </p>
  
  <p>Only thing I can think of is to find the closest point in set X for each point in set Y (so there will be exactly min(|X|,|Y|) matches, e.g. 50 in that case) and calculate RMSD from those matches. But the distance calculation and bipartite matching portion seems too computationally complex to call in a batch fashion. Any help in that area will help as well.</p>
  
  <p>Thanks!</p>