diff --git a/.github/workflows/emscripten.yml b/.github/workflows/emscripten.yml
index c499184..d575a71 100644
--- a/.github/workflows/emscripten.yml
+++ b/.github/workflows/emscripten.yml
@@ -9,7 +9,7 @@ jobs:
     steps:
       - uses: actions/checkout@v4
       - name: Build WASM wheel
-        uses: pypa/cibuildwheel@v3.1.2
+        uses: pypa/cibuildwheel@v3.1.3
         env:
           CIBW_PLATFORM: pyodide
       - name: Upload package
diff --git a/.github/workflows/wheel.yml b/.github/workflows/wheel.yml
index 22f4191..be1678a 100644
--- a/.github/workflows/wheel.yml
+++ b/.github/workflows/wheel.yml
@@ -33,7 +33,7 @@ jobs:
     steps:
       - uses: actions/checkout@v4
       - name: Build wheels
-        uses: pypa/cibuildwheel@v3.1.2
+        uses: pypa/cibuildwheel@v3.1.3
         env:
           CIBW_SKIP: "*_i686 *_ppc64le *_s390x *_universal2 *-musllinux_* cp314*"
           CIBW_PROJECT_REQUIRES_PYTHON: ">=3.10"
diff --git a/examples/plot_speed.py b/examples/plot_speed.py
index 2ff7e18..6a168e8 100644
--- a/examples/plot_speed.py
+++ b/examples/plot_speed.py
@@ -169,30 +169,38 @@ def baseline(X, y, t):
 # following iterated selection process, the eta-cosine method will be much faster than
 # the h-correlation method.
 
-X = rng.random((3000, 100))
+from timeit import repeat
+
+X = rng.random((3000, 400))
 y = rng.random((3000, 20))
 
-n_features_max = 30
+feature_num = np.arange(20, 61, step=10, dtype=int)
+
 
-time_h = np.zeros(n_features_max, dtype=float)
-time_eta = np.zeros(n_features_max, dtype=float)
-for i in range(n_features_max):
-    time_h[i] = timeit(
-        f"s = FastCan({i + 1}, verbose=0).fit(X, y)",
-        number=10,
+time_h = np.zeros(len(feature_num), dtype=float)
+time_eta = np.zeros(len(feature_num), dtype=float)
+for i, n_feats in enumerate(feature_num):
+    times_h = repeat(
+        f"s = FastCan({n_feats + 1}, verbose=0).fit(X, y)",
+        number=1,
+        repeat=10,
         globals=globals(),
     )
-    time_eta[i] = timeit(
-        f"s = FastCan({i + 1}, eta=True, verbose=0).fit(X, y)",
-        number=10,
+    time_h[i] = np.median(times_h)
+    times_eta = repeat(
+        f"s = FastCan({n_feats + 1}, eta=True, verbose=0).fit(X, y)",
+        number=1,
+        repeat=10,
         globals=globals(),
     )
+    time_eta[i] = np.median(times_eta)
+
 
-feature_num = np.arange(n_features_max, dtype=int) + 1
 plt.plot(feature_num, time_h, label="h-correlation")
 plt.plot(feature_num, time_eta, label=r"$\eta$-cosine")
 plt.title("Elapsed Time Comparison")
 plt.xlabel("Number of Selected Features")
 plt.ylabel("Elapsed Time (s)")
+plt.xticks(feature_num)
 plt.legend(loc="lower right")
 plt.show()