diff --git a/lib/iris/analysis/__init__.py b/lib/iris/analysis/__init__.py
index 41b6ec389c..25697266fa 100644
--- a/lib/iris/analysis/__init__.py
+++ b/lib/iris/analysis/__init__.py
@@ -1309,10 +1309,13 @@ def _calc_percentile(data, percent, fast_percentile_method=False, **kwargs):
                 "ignore",
                 "Warning: 'partition' will ignore the 'mask' of the MaskedArray.",
             )
-            result = np.percentile(data, percent, axis=-1)
+            result = np.percentile(data, percent, axis=-1, **kwargs)
+
         result = result.T
     else:
         quantiles = percent / 100.0
+        for key in ["alphap", "betap"]:
+            kwargs.setdefault(key, 1)
         result = scipy.stats.mstats.mquantiles(
             data, quantiles, axis=-1, **kwargs
         )
@@ -1344,9 +1347,9 @@ def _percentile(data, percent, fast_percentile_method=False, **kwargs):
         alternative to the scipy.mstats.mquantiles method. Does not handle
         masked arrays.
 
-    **kwargs
+    **kwargs : dict, optional
         passed to scipy.stats.mstats.mquantiles if fast_percentile_method is
-        False
+        False.  Otherwise passed to numpy.percentile.
 
     """
     if not isinstance(percent, Iterable):
@@ -1967,7 +1970,7 @@ def interp_order(length):
 """
 
 
-PERCENTILE = PercentileAggregator(alphap=1, betap=1)
+PERCENTILE = PercentileAggregator()
 """
 A :class:`~iris.analysis.PercentileAggregator` instance that calculates the
 percentile over a :class:`~iris.cube.Cube`, as computed by
@@ -1976,23 +1979,25 @@ def interp_order(length):
 
 **Required** kwargs associated with the use of this aggregator:
 
-* percent (float or sequence of floats):
+percent : float or sequence of floats
     Percentile rank/s at which to extract value/s.
 
 Additional kwargs associated with the use of this aggregator:
 
-* alphap (float):
+alphap : float
     Plotting positions parameter, see :func:`scipy.stats.mstats.mquantiles`.
     Defaults to 1.
-* betap (float):
+betap : float
     Plotting positions parameter, see :func:`scipy.stats.mstats.mquantiles`.
     Defaults to 1.
-* fast_percentile_method (boolean):
+fast_percentile_method : bool
     When set to True, uses :func:`numpy.percentile` method as a faster
-    alternative to the :func:`scipy.stats.mstats.mquantiles` method.  alphap and
-    betap are ignored. An exception is raised if the data are masked and the
-    missing data tolerance is not 0.
-    Defaults to False.
+    alternative to the :func:`scipy.stats.mstats.mquantiles` method.  An
+    exception is raised if the data are masked and the missing data tolerance
+    is not 0.  Defaults to False.
+
+kwargs : dict, optional
+    Passed to :func:`scipy.stats.mstats.mquantiles` or :func:`numpy.percentile`.
 
 **For example**:
 
diff --git a/lib/iris/tests/unit/analysis/test_PERCENTILE.py b/lib/iris/tests/unit/analysis/test_PERCENTILE.py
index a74dd891ba..bfd3234d26 100644
--- a/lib/iris/tests/unit/analysis/test_PERCENTILE.py
+++ b/lib/iris/tests/unit/analysis/test_PERCENTILE.py
@@ -196,6 +196,13 @@ def test_missing_mandatory_kwarg(self):
         with self.assertRaisesRegex(ValueError, emsg):
             PERCENTILE.aggregate("dummy", axis=0)
 
+    def test_wrong_kwarg(self):
+        # Test we get an error out of scipy if we pass the numpy keyword.
+        data = range(5)
+        emsg = "unexpected keyword argument"
+        with self.assertRaisesRegex(TypeError, emsg):
+            PERCENTILE.aggregate(data, percent=50, axis=0, method="nearest")
+
 
 class Test_fast_aggregate(tests.IrisTest, AggregateMixin):
     """Tests for fast percentile method on real data."""
@@ -239,6 +246,26 @@ def test_numpy_percentile_called(self, mocked_percentile):
         self.agg_method(data, axis=0, percent=42, fast_percentile_method=True)
         mocked_percentile.assert_called_once()
 
+        # Check that we left "method" keyword to numpy's default.
+        self.assertNotIn("method", mocked_percentile.call_args.kwargs)
+
+    @mock.patch("numpy.percentile")
+    def test_chosen_kwarg_passed(self, mocked_percentile):
+        data = np.arange(5)
+        percent = [42, 75]
+        axis = 0
+
+        self.agg_method(
+            data,
+            axis=axis,
+            percent=percent,
+            fast_percentile_method=True,
+            method="nearest",
+        )
+        self.assertEqual(
+            mocked_percentile.call_args.kwargs["method"], "nearest"
+        )
+
 
 class MultiAxisMixin:
     """
@@ -336,6 +363,29 @@ def test_numpy_percentile_called(self, mocked_percentile):
         as_concrete_data(result)
         mocked_percentile.assert_called()
 
+        # Check we have left "method" keyword to numpy's default.
+        self.assertNotIn("method", mocked_percentile.call_args.kwargs)
+
+    @mock.patch("numpy.percentile")
+    def test_chosen_method_kwarg_passed(self, mocked_percentile):
+        data = da.arange(5)
+        percent = [42, 75]
+        axis = 0
+
+        result = self.agg_method(
+            data,
+            axis=axis,
+            percent=percent,
+            fast_percentile_method=True,
+            method="nearest",
+        )
+
+        self.assertTrue(is_lazy_data(result))
+        as_concrete_data(result)
+        self.assertEqual(
+            mocked_percentile.call_args.kwargs["method"], "nearest"
+        )
+
 
 class Test_lazy_aggregate(
     tests.IrisTest, AggregateMixin, ScipyAggregateMixin, MultiAxisMixin