Fitting lines through large point clouds

Question

I have a large set of points (order of 10k points) formed by particle tracks (movement in the xy plane in time filmed by a camera, so 3d - 256x256px and ca 3k frames in my example set) and noise. These particles travel on approximately straight lines roughly (but only roughly) in the same direction, and so for the analysis of their trajectories I am trying to fit lines through the points. I tried to use Sequential RANSAC, but can't find a criterion to reliably single out false positives, as well as T- and J-Linkage, which were too slow and also not reliable enough.

Here is an image of a part of the dataset with good and bad fits I got with sequential Ransac: I'm using the centroids of the particle blobs here, blob sizes vary between 1 and about 20 pixels.

I found that subsamples using for example only every 10th frame worked quite well too, so the data size to be processed can be reduced this way.

I read a blog post about all the things neural networks can accomplish, and would like to ask you if this would be a feasible application for one before I start reading (I come from a non-maths background, so I would have to do quite a bit of reading)?

Or could you suggest a different method?

Thanks!

Addendum: Here is code for a Matlab function to generate a sample point cloud containing 30 parallel noisy lines, which I can't distinguish yet:

function coords = generateSampleData()
coords = [];
for i = 1:30
    randOffset = i*2;
    coords = vertcat(coords, makeLine([100+randOffset 100 100], [200+randOffset 200 200], 150, 0.2));
end

figure
scatter3(coords(:,1),coords(:,2),coords(:,3),'.')

function linepts = makeLine(startpt, endpt, numpts, noiseOffset)
    dirvec = endpt - startpt;
    linepts = bsxfun( @plus, startpt, rand(numpts,1)*dirvec); % random points on line
    linepts = linepts + noiseOffset*randn(numpts,3); % add random offsets to points
end

end

Answer 1

Based on the feedback and trying to find more effective approach I developed the following algorithm using a dedicated distance measure.

Following steps are performed:

1) Define a distance metric returning:

zero - if the points do not belong to a line

Euclidian distance of the points - if the points constitute a line according to the defined parameters, i.e.

• their distance is higher or equal than the min_line_length and

• their distance is lower or equal than the max_line_length and

• the line consists of at least min_line_points points with a distance lower that line_width/2 from the line

2) Calculate distance matrix using this distance measure (use sample of the data for large data sets; adjust the line parameters accordingly)

3) Find the points A and B with maximum distance - break to step 5) if the distance is zero

Note that if the distance is higher than zero the points A and B are building a line based on our definition

4) Get all points belonging to the line AB and remove them from the distance matrix. Repeat the step 3) to find another line

5) Check the coverage of the point with the selected lines, if substantial number of points remains uncovered, repeat the whole algorithm with adjusted line parameters.

6) In case that data sample was used - reassign all points to the lines and recalculate the boundary points.

Following parameters are used:

line width - line_width/2 is the allowed distance of the point from the ideal line = r line_width

minimum line length - points with shorter distance are not considered to belong to the same line = r min_line_length

maximum line length - points with longer distance are not considered to belong to the same line = r max_line_length

minimum points on a line - lines with less points are ignored = r min_line_points

With your data (after some fiddling with parameters) I got a good result covering all 30 lines.

More details can be found in the knitr script

Answer 2

I solved similar, though simpler, task with a brute force approach. The simplification was in the assumption, that the line is a linear function (in my case even the coefficients and intercept were in some known range).

This will therefore not solve your problem in general, where a particle can move orthogonal with axis x (i.e. it trace no function), but I post the solution as a possible inspiration.

1) Take all combinations of two points A and B with A(x) > B(x) + constant (to avoid the symmetry and a high error while calculating the coefficient)

2) Calculate the coefficient c and intercept i of the line AB

 A(y) = i + c * A(x)
 B(y) = i + c * B(x)
 A(y) - B(y) = c * (A(x) - B(x))
 c = (A(y) - B(y)) / (A(x) - B(x))
 i = A(y) - c * A(x)

3) Round the coefficient and intercept (this should eliminate / lower the problems with errors caused by the points in a grid)

4) For each intercept and coefficient calculate the number of points on this line

5) Consider only lines with points above some threshold.

Simple example see here