Using linear regression to find the ideal point given a set of trajectory's data

Question

I have a set of points in 2D obtained from a pendular movement with some noise. I want to determine where is equilibrium point ($x_0$, $y_0$) from which the rope is fixed. There are at least two approaches to solve this using linear regression. I'm not sure which one to choose and why the first method gave me no acceptable result (the indicated coordinate was outer the region expected).

Trajectory's data points

Method 1) Each displacement from one experimental point to the next give us a line named "s". Its perpendicular can be obtained through:

$m_{transversal} = - \frac{1}{m_s}$

Now we obtained a collection of lines whose angular and linear coeficients $m_{i}$ and $b_{i}$ are known: $m_{i}x + b_{i} = y$

Ideally, the equilibrium point dwells in every i-line:

\begin{array}{lr} x_{0}m_1 + y_{0} = b_1\\ x_{0}m_2 + y_{0} = b_2\\ ...\\ x_{0}m_i + y_{0} = b_i \end{array}

Using linear regression of N lines one can find $x_0$ and similarly also $y_0$:

$x_{0}$ = $\frac{N \sum (m_ib_i) - (\sum mb_i) (\sum m_i)}{N (\sum m²_i) - (\sum m_i)²}$

My code is shown below:

if((i+1) % 5 == 0){ //i+1 remember: i is initialized as zero
            cout << "Let's sum up all cartesian elements" << endl;
            m_sum += m;
            b_sum += b;
            m2_sum += m*m;
            mb_sum += m*b;
            N = i+1;
        }
        m_sum2 += m_sum*m_sum;
        x0 = (N*mb_sum     - b_sum*m_sum)  / (N*m2_sum - m_sum2);
        y0 = (m2_sum*b_sum - mb_sum*m_sum) / (N*m2_sum - m_sum2);

I implemented it but it returned me (-350, -140). This point dwells outer the region expected.

Method 2) Using a mathematical model that describes the movement (a parabole $y = a(x-x_0)²+b(x-x_0)+c$, for example). And do linear regression dealing with some matrices with which I'm not used to.

Thanks!

Answer 1

For the pendulum, the following equations hold: $$\sin\varphi = \frac{x-x_0}{l}\quad\mbox{and}\quad \cos\varphi=\frac{l-(y-y_0)}{l}$$ where $l$ is the pendulum length, $\varphi$ is the angle of the rope, and $(x_0,y_0)$ is the lowest point of the pendulum movement. With the use of $\cos^2\varphi=1-\sin^2\varphi$, it follows that $$y=y_0+l-\sqrt{l^2-(x-x_0)^2}$$ You can thus use a least squares fit of your points to this function to find estimators for $x_0$, $y_0$ and $l$ and then the fixing point of the rope is $(x_0, y_0+l)$.

Note that when you use the function lm() in R for doing the LSQ fit, you must mask the square in the formula with I(...) (I = capital 'i'). Alternatively you can use minpack.lm (Levenberg-Marquardt algorithm) in R or its equivalent in PythonTo find the LSQ-Minimum.