None of the ones you mentioned.
Use modified
min-heap,
of size k
.
Sketch:
- Any correct algorithm will necessarily examine every element, with O(n) cost.
- k ≪ n
- Typically the heap will be full, containing k entries.
- Inserts are cheap, O(1) time, but they're paired with removes that have O(log k) cost.
- During most of the scan, we can almost always avoid those costs.
Procedure:
- At all times let
threshold
be the minimum heap value; we read it without removing it.
- Read
k
elements to populate the heap. Now it is full.
- Scan over each remaining element
e
.
- 3.1. If
e <= threshold
then loop back to (3.) to keep scanning.
- 3.2. Do
heap.remove_min()
.
- 3.3. Do
heap.insert(e)
- 3.4. Loop back to (3.) to continue scanning.
The threshold will ratchet up, up, ever more slowly but always upward.
At the end the heap contains the k
largest values, in sorted order.
Worst case is O(n + n log k), meaning O(n log k),
when input is already sorted.
In which case prefer to scan in reverse direction,
for cost of O(n).
Making an initial sample of k random locations
is also a valid approach for filling the heap,
which is robust against adversarial input patterns.
For a random permutation of distinct values
the cost is O(n + m log k), where m relates to
the distribution of input values (slow ratcheting behavior)
and we expect that m ≪ n.
The key here is our threshold comparison has constant cost,
unrelated to magnitude of k.
And for each given element it nearly always says
that we needn't bother with heap operations for it.
Here is a graphical depiction of performance.
Number of heap push() operations looks logarithmic in n
, to me.
Here is x="ratio"
, on a log-log scale.
Notice that for each value of k
we run
through four decades of n
, 1e4 up to 1e7.
from array import array
from collections.abc import Generator
from dataclasses import asdict, dataclass
from heapq import heapify, heappop, heappushpop
from pathlib import Path
from random import shuffle
from time import time
from typing import Any
from matplotlib import pyplot as plt
from seaborn.palettes import SEABORN_PALETTES as palette
from tqdm import tqdm
import pandas as pd
import seaborn as sns
def ratchet(ranks: array[int]) -> Generator[dict[str, int], None, None]:
"""Demonstrates the speed at which we slowly ratchet up a min bound."""
bound = ranks[0]
prev = bound
for i, rank in enumerate(ranks):
if bound < rank:
bound = rank
yield {"idx": i, "bound": bound, "delta": bound - prev}
prev = bound
def _get_xs(n: int) -> array[int]:
xs = array("L", range(n))
shuffle(xs)
return xs
@dataclass
class Result:
pushes: int
elapsed: float
def experiment(n: int, k: int) -> Result:
xs = list(_get_xs(n))
t0 = time()
h, xs = xs[:k], xs[k:]
heapify(h)
i = 0
for x in xs:
if x > h[0]:
heappushpop(h, x)
i += 1
elapsed = round(time() - t0, 6)
for x in range(n - k, n):
assert x == heappop(h)
return Result(i, elapsed)
def run_experiments() -> Generator[dict[str, float], None, None]:
"""Sweeps through a few decades of problem sizes."""
for trial in tqdm(range(12)):
n = 10_000
while n <= 10_000_000:
for k in [1, 10, 100, 1_000]:
param = {"n": n, "k": k}
result = experiment(**param)
yield {**param, **asdict(result)}
n *= 10
def main(csv: Path = Path("/tmp/heap.csv")) -> None:
ratchet_df = pd.DataFrame(ratchet(_get_xs(int(1e7))))
no_scientific_notation: dict[str, Any] = {"disable_numparse": True}
print(ratchet_df.to_markdown(index=False, **no_scientific_notation))
if not csv.exists():
df = pd.DataFrame(run_experiments())
df.to_csv(csv, index=False)
df = pd.read_csv(csv)
df["ratio"] = df.n / df.k
print(df)
sns.catplot(data=df, x="n", y="pushes", hue="k", palette=palette["bright"])
# plt.gca().set(xscale="log", yscale="log")
plt.show()
if __name__ == "__main__":
main()
n
?