你得到的错误是因为check_gradient只接受平面点数组。如果你使用一个形状为(20,)的数组x0而不是(4, 5)，它应该可以工作。但事实并非如此！
下面是我的安装中approx_fprime的实现（scipy.__version__ = '0.9.0'）：

def approx_fprime(xk,f,epsilon,*args):
    f0 = f(*((xk,)+args))
    grad = numpy.zeros((len(xk),), float)
    ei = numpy.zeros((len(xk),), float)
    for k in range(len(xk)):
        ei[k] = epsilon
        grad[k] = (f(*((xk+ei,)+args)) - f0)/epsilon
        ei[k] = 0.0
    return grad

我已经看过它好几次了，很难相信如此糟糕的代码会出现在scipy发行版中，我相信我一定错过了什么。但恐怕这是不对的。如果将其替换为：

def approx_fprime(xk,f,epsilon,*args):
    return (f(*((xk + epsilon,) + args)) - f(*((xk,) + args))) / epsilon

现在它为我工作。使用x0.shape = (20,)，可以得到：

In [2]: error
Out[2]: 1.746097524556073e-08

关于x0.shape = (4, 5)：

In [4]: error
Out[4]: 
array([  1.03560895e-08,   1.45994321e-08,   8.54143390e-09,
         1.09225833e-08,   9.85988655e-09])

因此，它似乎真的没有准备好在其他地方的非平面阵列了。但无论哪种方式，实现都非常破碎：你应该提交一份bug报告。

关于你的代码（需要 * 矢量化 * -如果你确实需要-做 flatten（） 或 ravel（））：

import numpy as np
from scipy import optimize

def sigmoid_obj(x, gradient=False):
    y = 1 / (1 + np.exp(-1*x))
    res= np.multiply(y, 1 - y) if gradient else y
    return sum(res)  # must return scalar !

def num_gradient(x):
    return optimize.approx_fprime(x, sigmoid_obj, eps)

x = np.array(np.random.uniform(-30, 30, (4, 5)))     # shape(4,5)
##x0= [0,]        # guess !! var(s) used in objective_func

eps = np.sqrt(np.finfo(float).eps)

error = optimize.check_grad(sigmoid_obj, num_gradient, x.flatten() )
print(error)

但通常情况下，objective_func只用于y_target的一个x_point投影，从x0=guess开始优化。在最简单的代码版本中：

import numpy as np
from scipy import optimize

def sigmoid_obj(x, gradient=False):
    y = 1 / (1 + np.exp(-1*x))
    res= np.multiply(y, 1 - y) if gradient else y
    return res

def num_gradient(x):
    return optimize.approx_fprime(x, sigmoid_obj, eps)

x = np.array(np.random.uniform(-30, 30, (4, 5)))
x0= [0,]        # guess !! init_var(s) used in objective_func

eps = np.sqrt(np.finfo(float).eps)

error = optimize.check_grad(sigmoid_obj, num_gradient, x0 )
print(error)

？如果假设initial_guess只有一个x，则guess=x0必须是[0,]
P.S. optimize.check_grad()在内部使用optimize.approx_fprime()，最新的计算是从一个点投影到y [本地]. OR可以分析计算（使用SymPy或手动计算，如docs的示例）：将numpy作为np从scipy import optimize导入

def func(x):
    return x[0]**2 - 0.5 * x[1]**3
def grad(x):
    return [2 * x[0], -1.5 * x[1]**2]   # here can be numerical_gradient=optimize.approx_fprime(...)
from scipy.optimize import check_grad
print(check_grad(func, grad, x0= [1.5, -1.5]))