Linear Algebra - SVD Fit

Fit a line with Singular Value Decomposition

## imports
import numpy as np
import matplotlib.pyplot as plt
import plotly.graph_objects as go

# fixing numpy random state for reproducibility
np.random.seed(0)

# input data
n = 25
points = np.array( [5*np.arange(n),5+5*np.arange(n),5*np.ones(n)] ).T + 0.5-np.random.rand(n,3)
points

array([[ -0.43239394,   5.2846018 ,   4.64166236],
       [  4.69710663,  10.34085376,   4.89428804],
       [ 10.38433813,  14.77211184,   4.86253772],
       [ 14.68806144,  20.02061545,   4.58513691],
       [ 20.45065105,  25.20711143,   4.7849474 ],
       [ 25.08189079,  30.32704865,   5.39278925],
       [ 29.68266089,  35.02685702,   4.61771633],
       [ 34.76671087,  40.09027379,   5.12648899],
       [ 39.98436165,  44.61094005,   4.76272142],
       [ 45.49484704,  49.80584215,   4.58049259],
       [ 49.78954424,  55.32299422,   5.01648187],
       [ 55.35968398,  60.14100472,   4.56288296],
       [ 59.57669469,  65.21716315,   5.16036896],
       [ 64.89978713,  69.5368027 ,   5.35219867],
       [ 70.24308336,  74.62644317,   5.00810777],
       [ 74.60103891,  80.3144821 ,   4.96733141],
       [ 80.17373037,  85.18345744,   5.05312304],
       [ 85.06692255,  90.14265312,   4.58502923],
       [ 89.76825581,  94.77245301,   5.21008655],
       [ 94.92229058,  99.72082057,   4.70440963],
       [100.15546954, 104.72912724,   4.7641061 ],
       [105.35849351, 109.63405453,   5.05867853],
       [110.01358955, 115.05163082,   4.932154  ],
       [114.87883075, 120.00182043,   4.63321146],
       [119.87226524, 125.09857205,   5.08330824]])

# calculating the mean of the points
avg = np.mean(points, axis=0)

# subtracting the mean from all points
subtracted = points - avg

# performing SVD
_, _, V = np.linalg.svd(subtracted)

# find the direction vector (which is the right singular vector corresponding to the largest singular value)
direction = V[0, :]

# A line is defined by the average and its direction
p0 = avg
d = direction
print(d)

[7.08153575e-01 7.06058100e-01 6.87968476e-04]

d0 = np.array([0, 0, 1])
angle = np.arccos(np.dot(d0,d)/(np.linalg.norm(d0) * np.linalg.norm(d)))
print(angle*180/np.pi)

89.96058230676098

pa = p0 + (-100)*d
pb = p0 + 100*d
print(pa, pb)

[-10.83624093  -5.60662059   4.82481353] [130.79447411 135.60499941   4.96240723]

## plotly
trace1 = go.Scatter3d(
    x=[pa[0],pb[0]],
    y=[pa[1],pb[1]],
    z=[pa[2],pb[2]],
    mode='lines',
    name='3D fitted line',
    line=go.scatter3d.Line(color='rgb(255,0,0)', width=10),
    hoverinfo='none')

labels = []
for i in range(n): labels += [str(i)]  
trace2 = go.Scatter3d(
    x=points[:,0],
    y=points[:,1],
    z=points[:,2],
    mode='markers',
    name='Points',
    marker=go.scatter3d.Marker(
        symbol='cross',
        opacity=1,
        color='rgb(0,200,127)'),
    text=labels,
    hoverinfo='text')

layout = go.Layout(
            title="3D line fitting",
            scene=go.layout.Scene(
                    xaxis_title="x",
                    yaxis_title="y",
                    zaxis_title="z",
                    camera=dict(
                           up=dict(x=0, y=0, z=1),
                           center=dict(x=0, y=0, z=0),
                           eye=dict(x=0, y=2.5, z=0))))     


fig=go.Figure(data=[trace1, trace2], layout=layout)
fig.show(renderer="iframe_connected", config={'showLink': False})

---

source

fitline

 fitline (points)