title: Least Median of Squares robust regression C++
body: <p>I have a set of data <code>z(0), z(1), z(2)...,z(n)</code> that I am currently fitting with a 2 variables polynomial of the kind <code>p(x,y) = a(1)*x^2+a(2)*y^2+a(3)*x*y+a(4)</code>. I have <code>i=1,...,n</code> <code>(x(i),y(i))</code> coordinates that I impose to be <code>p(x(i),y(i))=z(i)</code>. In this way I have a Overdetermined System that I can solve  using <a href="https://eigen.tuxfamily.org/dox-devel/group__LeastSquares.html" rel="nofollow noreferrer">Eigen SVD</a> . I am looking for a more sophisticated method that can take care of outliers, like a Least Median of Squares robust regression (as described <a href="http://web.ipac.caltech.edu/staff/fmasci/home/astro_refs/LeastMedianOfSquares.pdf" rel="nofollow noreferrer">here</a>) but I haven't found a C++ implementation for 2 variables. I looked in GSL but it seems there is nothing for 2 variable functions. The only other solution I can think of is using a TGraph2D in <a href="https://root.cern.ch/doc/master/classTGraph2D.html#a52f2eb4c4cf50b7d086e3e7abe25d647" rel="nofollow noreferrer">ROOT</a>. Do you know any other solution? Numerical recipes maybe? Since I am writing C++ code I would prefer C or C++ implementations.</p>