diff --git a/R/utils.R b/R/utils.R
index e668606..f680e4d 100644
--- a/R/utils.R
+++ b/R/utils.R
@@ -8,6 +8,26 @@ power_function <- function(x) {
    10^(1:x)
 }
 
-# power_function <- function(x) {
-#   as.integer(outer(c(1), 10^(1:x)))
-# }
+#' Thin samples
+#' Thin a large sample of sorted z-scores and stages for efficiency.
+#' Retains 5,000 evenly-spaced points across the full range plus the last
+#' 10,000 points for dense coverage at the extreme tail.
+#'
+#' @param n            Total number of samples
+#' @param z_sorted     Sorted z-scores (ascending), length n
+#' @param stage_sorted Sorted stage values (ascending), length n
+#' @return A tibble with columns z_aep and stage
+thin_samples <- function(n, z_sorted, stage_sorted) {
+
+  thin_idx <- unique(c(
+    round(seq(1, n, length.out = 5000)),
+    (n - 10000):n
+  ))
+  thin_idx <- sort(unique(thin_idx))
+  thin_idx <- thin_idx[thin_idx >= 1 & thin_idx <= n]
+
+  data.frame(
+    z_aep = z_sorted[thin_idx],
+    stage = stage_sorted[thin_idx]
+  )
+}