Since you say you have a 2-D array of points, I assume you have a function of two variables f(x, y). That means you don't have a single derivative. Instead you get a set of partial derivatives.
You could approximate the partial derivatives using finite difference formulas.
The partial derivative with respect to x at f(x, y) would be (f(x+h, y) - f(x-h, y))/2h.
The partial derivative with respect to y at f(x, y) would be (f(x, y+h) - f(x, y-h))/2h.
In these formulas, h is the space between nodes on your grid, assuming you have a regularly spaced grid. If the horizontal and vertical spacings are different, use the horizontal spacing for the partial with respect to x and the vertical spacing for the partial with respect to y.
Update: I misunderstood your question. I thought the 2-D array was an array of domain values. If you have a list of x and f(x) values, you can approximate f'(x) as (f(x+h) - f(x-h)) / 2h. This will work everywhere except at the first and last points where one of the terms will be out of range. You can use (f(x + h) - f(x))/h at the left end and (f(x) - f(x-h))/h at the right end. The approximation will be less accurate at the end points but that can't be avoided.