docs: def autocorrelation(x, lag):
      """
      Calculates the autocorrelation of the specified lag, according to the formula [1]
  
      .. math::
  
          \\frac{1}{(n-l)\sigma^{2}} \\sum_{t=1}^{n-l}(X_{t}-\\mu )(X_{t+l}-\\mu)
  
      where :math:`n` is the length of the time series :math:`X_i`, :math:`\sigma^2` its variance and :math:`\mu` its
      mean. `l` denotes the lag.
  
      .. rubric:: References
  
      [1] https://en.wikipedia.org/wiki/Autocorrelation#Estimation
  
      :param x: the time series to calculate the feature of
      :type x: numpy.ndarray
      :param lag: the lag
      :type lag: int
      :return: the value of this feature
      :return type: float
      """
      # This is important: If a series is passed, the product below is calculated
      # based on the index, which corresponds to squaring the series.
      if type(x) is pd.Series:
          x = x.values
      if len(x) < lag:
          return np.nan
      # Slice the relevant subseries based on the lag
      y1 = x[:(len(x)-lag)]
      y2 = x[lag:]
      # Subtract the mean of the whole series x
      x_mean = np.mean(x)
      # The result is sometimes referred to as "covariation"
      sum_product = np.sum((y1 - x_mean) * (y2 - x_mean))
      # Return the normalized unbiased covariance
      v = np.var(x)
      if np.isclose(v, 0):
          return np.NaN
      else:
          return sum_product / ((len(x) - lag) * v)