Numerical derivative at order 3 I'd like to find the formula for the 3rd order numerical derivative in order to further implement a Python function for a time series (Python only have a function from scipy which takes an unknown variable). In order to demonstrate and calculate the general formula I just took the 2 main formulas :
$$f'(x) = \frac{f(x+h) - f(x-h)}{2h}$$
$$f''(x) = \frac{f(x+h) - 2f(x) + f(x+h)}{h^2}$$
then I made the following hypothesis that : 
$$ f'''(x) = \frac{f''(x+h) - f''(x-h)}{2h}$$
where,
$$ f''(x+h) = \frac{f(x+h+h)-2f(x+h)+f(x+h-h)}{2h}$$
$$ = \frac{f(x+2h)-2f(x+h)+f(x)}{2h}$$
and same calculation for the other,
$$ f''(x-h) = \frac{f(x-2h)-2f(x-h) + f(x)}{2h}$$
which gives me at the end :
$$f'''(x) = \frac{f(x+2h) + 2h(f(x-h) - f(x+h))-f(x-2h)}{4h^2}$$
But I just can't find anywhere if my calculation is right or not,


*

*Is my result good ? 

*What about the calculation method ? 


I'm doubting as I'm not studying at school anymore I don't have a teacher to correct me.
 A: If you want to approximate $f'''(x)$ in terms of $f$ evaluated at $x-2h, x-h, x+h,x+2h$, i.e., 
$$f'''(x)\approx af(x-2h)+bf(x-h)+cf(x+h)+df(x+2h), $$
you better make sure that the result for $f'''(0)$ is zero for $f(x)=1$, $f(x)=x$, or $f(x)=x^2$, and is $=6$ for $f(x)=x^3$. This leads us to the conditions
$$ a+b+c+d=0$$
$$ 2a+b-c-2d=0$$
$$ 4a+b+c+4d=0$$
and
$$-8a-b+c+8d=\frac 6{h^3}. $$
The first three quickly lead to $d=-a$, $c=-b$, $b=-2a$, so that the last becomes $-8a+2a+2a-8a=\frac6{h^3}$ and ultimately
$$\tag1 f'''(x)\approx \frac{-f(x-2h)+2f(x-h)-2f(x+h)+f(x+2h)}{2h^3}$$
The interesting question is whether it would be preferable to use other points of evaluation; four points seems adequate because we want to impose four conditions; symmetric with respect to $x$ seems adequate as well; however, evenly spaced points  $x-3h,x-h,x+h,x+3h$ look attractive as well.
By the same method as above, we can determine the best coefficients and find
$$\tag2 f'''(x)\approx\frac{-f(x-3h)+3f(x-h)-3f(x+h)+f(x+3h)}{8h^3}$$
Which of these fomulae is better?
Incidentally, both give us $f'''(0)=0$ for $f(x)=x^4$ (or any even function). For $f(x)=x^5$, we obtain
$$ f'''(x)\approx 30h^2$$
from $(1)$ and
$$ f'''(x)\approx 60h^2$$
from $(2)$. In this respect, it seems that $(1)$ is better, but take into account that by the wider spread of points, $(2)$ effectively uses a $\frac32$ times bigger $h$, and $(3/2)^2>\frac{60}{30}$.
A: Look at the divided difference on equidistant nodes: $$f[x,x+h,x+2h,x+3h]=\frac{f(x+3h)-3f(x+2h)+3f(x+h)-f(x)}{6h^3}.$$ If $f$ is 3-times continuously differentiable at $x$, then this term tends to $\dfrac{f^{\prime\prime\prime}(x)}{6}$ as $h\to 0.$
