Skip to content

Instantly share code, notes, and snippets.

@BluBb-mADe
Created September 3, 2023 15:09
Show Gist options
  • Save BluBb-mADe/cdd7ab07ff1f682adaf6d657e43407d9 to your computer and use it in GitHub Desktop.
Save BluBb-mADe/cdd7ab07ff1f682adaf6d657e43407d9 to your computer and use it in GitHub Desktop.
Simple Theil-Sen estimator example
#pragma once
template<int Limit, std::floating_point ValueType, std::ranges::sized_range Range, ValueType Epsilon = std::numeric_limits<ValueType>::epsilon() * Limit>
requires std::convertible_to<std::ranges::range_value_t<Range>, std::pair<ValueType, ValueType>>
auto theil_sen_slope(Range value_pairs) {
if (value_pairs.size() <= 1) {
return ValueType{0};
}
std::vector<std::pair<ValueType, ValueType>> sample_holder;
Range sample_range = value_pairs;
if (sample_range.size() > Limit) {
sample_holder.reserve(Limit);
std::ranges::sample(value_pairs, std::back_inserter(sample_holder), Limit, g_gen);
sample_range = sample_holder;
}
std::vector<ValueType> slopes;
slopes.reserve(sample_range.size() * value_pairs.size());
for (auto const& sample : sample_range) {
for (auto const& pair : value_pairs) {
auto [x1, y1] = sample;
auto [x2, y2] = pair;
ValueType den = x1 - x2;
if (std::abs(den) > Epsilon) {
slopes.emplace_back((y1 - y2) / den);
}
}
}
if (slopes.empty()) {
return ValueType{0};
}
// Median
std::nth_element(slopes.begin(), slopes.begin() + slopes.size() / 2, slopes.end());
return slopes[slopes.size() / 2];
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment