docs: def stationary_distribution_sensitivity(T, j):
      r"""Calculate the sensitivity matrix for entry j the stationary
      distribution vector given transition matrix T.
  
      Parameters
      ----------
      T : numpy.ndarray shape = (n, n)
          Transition matrix
      j : int
          entry of stationary distribution for which the sensitivity is to be computed
  
      Returns
      -------
      x : ndarray, shape=(n, n)
          Sensitivity matrix for entry index around transition matrix T. Reversibility is not assumed.
  
      Remark
      ------
      Note, that this function uses a different normalization convention for the sensitivity compared to
      eigenvector_sensitivity. See there for further information.
      """
  
      n = len(T)
  
      lEV = numpy.ones(n)
      rEV = stationary_distribution(T)
      eVal = 1.0
  
      T = numpy.transpose(T)
  
      vecA = numpy.zeros(n)
      vecA[j] = 1.0
  
      matA = T - eVal * numpy.identity(n)
      # normalize s.t. sum is one using rEV which is constant
      matA = numpy.concatenate((matA, [lEV]))
  
      phi = numpy.linalg.lstsq(numpy.transpose(matA), vecA, rcond=-1)
      phi = numpy.delete(phi[0], -1)
  
      sensitivity = -numpy.outer(rEV, phi) + numpy.dot(phi, rEV) * numpy.outer(rEV, lEV)
  
      return sensitivity