docs: def plot_CI(
      ax,
      sampler,
      modelidx=0,
      sed=True,
      confs=[3, 1, 0.5],
      e_unit=u.eV,
      label=None,
      e_range=None,
      e_npoints=100,
      threads=None,
      last_step=False,
  ):
      """Plot confidence interval.
  
      Parameters
      ----------
      ax : `matplotlib.Axes`
          Axes to plot on.
      sampler : `emcee.EnsembleSampler`
          Sampler
      modelidx : int, optional
          Model index. Default is 0
      sed : bool, optional
          Whether to plot SED or differential spectrum. If `None`, the units of
          the observed spectrum will be used.
      confs : list, optional
          List of confidence levels (in sigma) to use for generating the
          confidence intervals. Default is `[3,1,0.5]`
      e_unit : :class:`~astropy.units.Unit` or str parseable to unit
          Unit in which to plot energy axis.
      e_npoints : int, optional
          How many points to compute for the model samples and ML model if
          `e_range` is set.
      threads : int, optional
          How many parallel processing threads to use when computing the samples.
          Defaults to the number of available cores.
      last_step : bool, optional
          Whether to only use the positions in the final step of the run (True,
          default) or the whole chain (False).
      """
      confs.sort(reverse=True)
  
      modelx, CI = _calc_CI(
          sampler,
          modelidx=modelidx,
          confs=confs,
          e_range=e_range,
          e_npoints=e_npoints,
          last_step=last_step,
          threads=threads,
      )
      # pick first confidence interval curve for units
      f_unit, sedf = sed_conversion(modelx, CI[0][0].unit, sed)
  
      for (ymin, ymax), conf in zip(CI, confs):
          color = np.log(conf) / np.log(20) + 0.4
          ax.fill_between(
              modelx.to(e_unit).value,
              (ymax * sedf).to(f_unit).value,
              (ymin * sedf).to(f_unit).value,
              lw=0.001,
              color=(color,) * 3,
              alpha=0.6,
              zorder=-10,
          )
  
      _plot_MLmodel(ax, sampler, modelidx, e_range, e_npoints, e_unit, sed)
  
      if label is not None:
          ax.set_ylabel(
              "{0} [{1}]".format(label, f_unit.to_string("latex_inline"))
          )