Skip to content

Instantly share code, notes, and snippets.

@geohot
Last active July 31, 2023 06:39
Show Gist options
  • Save geohot/9743ad59598daf61155bf0d43a10838c to your computer and use it in GitHub Desktop.
Save geohot/9743ad59598daf61155bf0d43a10838c to your computer and use it in GitHub Desktop.
RANSAC polyfit. Fit polynomials with RANSAC in Python
def ransac_polyfit(x, y, order=3, n=20, k=100, t=0.1, d=100, f=0.8):
# Thanks https://en.wikipedia.org/wiki/Random_sample_consensus
# n – minimum number of data points required to fit the model
# k – maximum number of iterations allowed in the algorithm
# t – threshold value to determine when a data point fits a model
# d – number of close data points required to assert that a model fits well to data
# f – fraction of close data points required
besterr = np.inf
bestfit = None
for kk in xrange(k):
maybeinliers = np.random.randint(len(x), size=n)
maybemodel = np.polyfit(x[maybeinliers], y[maybeinliers], order)
alsoinliers = np.abs(np.polyval(maybemodel, x)-y) < t
if sum(alsoinliers) > d and sum(alsoinliers) > len(x)*f:
bettermodel = np.polyfit(x[alsoinliers], y[alsoinliers], order)
thiserr = np.sum(np.abs(np.polyval(bettermodel, x[alsoinliers])-y[alsoinliers]))
if thiserr < besterr:
bestfit = bettermodel
besterr = thiserr
return bestfit
@arielisrael
Copy link

does x[maybeinliers] work for you?
I got:
TypeError: only integer scalar arrays can be converted to a scalar index
Needed to create lists of x and y values through list comprehension to use instead of x[maybeinliers] and y[maybeinliers]

@jonathan-cartica
Copy link

instead of np.random.randint(len(x), size=n) -> sample with replacement
it is better to use random.sample(list(range(len(x))), size=n) -> sample without replacement
because randint can choose the same index multiple times

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment