Divided differences are defined like this:
\begin{align} [y_n] & = y_n \quad (1) \\ [y_{n}, y_{n+1}] & = \dfrac{[y_{n+1}] - [y_n]}{x_{n+1} - x_n} \quad (2) \\ [y_{n}, y_{n+1}, \dots, y_{n+k}] & = \dfrac{[y_{n+1}, y_{n+2}, \dots, y_{n+k}] - [y_{n}, y_{n+1}, \dots, y_{n+k-1}]}{x_{n+k} - x_n} \quad (3) \end{align}
where $y_n = f(x_n)$.
See this answer for the derivation of the above recursion. The question is how to connect the divided differences with function derivatives. Looking at the above equations, if we say
$x_{n+k} = x_n + kh$,
we can write
\begin{equation} \lim_{h\to 0} [y_{n}, y_{n+1}] = [y_n,y_n] = \lim_{h\to 0} \dfrac{f(x_n + h) - f(x_n)}{h} = f'(x_n) \quad (4) \end{equation}
and from (4) and (3), we get
\begin{equation} \lim_{h\to 0} [y_{n}, y_{n+1}, y_{n+2}] = \lim_{h\to 0} \dfrac{[y_{n+1}, y_{n+2}] - [y_{n}, y_{n+1}]}{x_{n+2} - x_n} = \lim_{h\to 0} \dfrac{f'(x_n+h) - f'(x_n) }{2h} = \frac{ f''(x_n)}{2} \end{equation}
I can go on and write this further and "see" that
\begin{equation} \lim_{h\to 0} [y_{n}, y_{n+1}, y_{n+2}, \dots y_{n+k}] = \dfrac{f^k(x_n)}{k!} \quad (5) \end{equation}
But how to prove this?
Using (3) I got this:
\begin{align} \lim_{h\to 0} [y_{n}, \dots y_{n+k}] & = \lim_{h\to0} \dfrac{[y_{n+1}, \dots, y_{n+k}]-[y_{n}, \dots, y_{n+k-1}]}{kh} \\ & = \lim_{h \to 0}\dfrac{\dfrac{[y_{n+2}, \dots, y_{n+k}]-[y_{n+1}, \dots, y_{n+k-1}]}{(k-1)h} - \dfrac{[y_{n+1}, \dots, y_{n+k-1}]-[y_{n}, \dots, y_{n+k-2}]}{(k-1)h} }{kh} \end{align}
and got stuck here. I can see that if I were to write out recursively the fractions, I would continue to get $(k-l)$ for every level $l$, ultimately leading to (4), that can be back-substituted to give (5).
Another way might be to use the polynomial interpolation error:
$$f(x) - P_n(x) = \dfrac{f^{n+1}(\zeta)}{(n+1)!}\Pi_{i=0}^{n}(x-x_i) \quad (6)$$
where $\zeta \in [x_0, \dots, x_n]$.
We can have $P_n(x)$ interpolate $f(x)$ over $[x_{n_0}, x_{n_1}, \dots x_{n_n}]$ using Hermite interpolation
$$P_n(x_{n_0}) = y_{n_0} $$
$$P_n'(x_{n_1}) = y'_{n_1} $$
$$P_n''(x_{n_2}) = y''_{n_2} $$
$$P_n^3(x_{n_3}) = y^3_{n_3} $$
$$\dots$$
$$P_n^n(x_{n_n}) = y^n_{n_n} $$
If $P_n(x)$ interpolates $f(x)$, then from (6), $\zeta = x_n$, so
$$f(x) - P_n(x) = \dfrac{f^{n+1}(x_n)}{(n+1)!}(x-x_{n})^{n+1} \quad (7)$$
We can extend the support with one more point $(\tilde{x}, \tilde{y}), \tilde{x} \ne x_n$ and get a Newton polynomial extended by a term:
$$P_{n+1}(x) = P_n(x) + [y_{n_0},\dots,y_{n_n}, \tilde{y}](x-x_{n})^{n+1} \quad (8) $$
Inserting (8) into (7) we get
$$[y_{n_0},\dots,y_{n_n}, \tilde{y}] = \dfrac{f^{n+1}(\zeta)}{(n+1)!}, \quad \zeta \in [x_n, \tilde{x}]$$
and from this
$$[y_{n_0},\dots,y_{n_n}] = \dfrac{f^{n}(x_n)}{(n+1)!}$$
Is this fine?
What about the proof by induction?
