code: def approx_fprime(
      x, f, epsilon=None, args=(), kwargs={}, centered=False
  ):
      """
      Gradient of function, or Jacobian if function f returns 1d array
  
      Parameters
      ----------
      x : ndarray
          parameters at which the derivative is evaluated
      f : function
          `f(*((x,)+args), **kwargs)` returning either one value or 1d array
      epsilon : float, optional
          Stepsize, if None, optimal stepsize is used. This is EPS**(1/2)*x for
          `centered` == False and EPS**(1/3)*x for `centered` == True.
      args : tuple
          Tuple of additional arguments for function `f`.
      kwargs : dict
          Dictionary of additional keyword arguments for function `f`.
      centered : bool
          Whether central difference should be returned. If not, does forward
          differencing.
  
      Returns
      -------
      grad : ndarray
          gradient or Jacobian
  
      Notes
      -----
      If f returns a 1d array, it returns a Jacobian. If a 2d array is returned
      by f (e.g., with a value for each observation), it returns a 3d array
      with the Jacobian of each observation with shape xk x nobs x xk. I.e.,
      the Jacobian of the first observation would be [:, 0, :]
      """
      n = len(x)
      f0 = f(*((x,) + args), **kwargs)
      dim = np.atleast_1d(f0).shape  # it could be a scalar
      grad = np.zeros((n,) + dim, np.promote_types(float, x.dtype))
      ei = np.zeros((n,), float)
      if not centered:
          epsilon = _get_epsilon(x, 2, epsilon, n)
          for k in range(n):
              ei[k] = epsilon[k]
              grad[k, :] = (
                  f(*((x + ei,) + args), **kwargs) - f0
              ) / epsilon[k]
              ei[k] = 0.0
      else:
          epsilon = _get_epsilon(x, 3, epsilon, n) / 2.0
          for k in range(n):
              ei[k] = epsilon[k]
              grad[k, :] = (
                  f(*((x + ei,) + args), **kwargs)
                  - f(*((x - ei,) + args), **kwargs)
              ) / (2 * epsilon[k])
              ei[k] = 0.0
  
      if n == 1:
          return grad.T
      else:
          return grad.squeeze().T
  
file_path: statsmodels/statsmodels/tools/numdiff.py
imports: import numpy
  from statsmodels.compat.pandas import Appender, Substitution
label_desc: Use image node_id 2 for calling a global function with example usage: approx_fprime(x, f, epsilon, args, kwargs, centered) and returns: grad, grad