fractional derivative on discrete domain ( Grünwald–Letnikov derivative) In programming, I often use 
$$
\frac{y_{i-1}-y_i}{h}
$$
to represent the derivative of an array of numbers $\{y_i\}$, ($h>0$ is a constant discretization span). Likewise,
$$
\frac{(y_{i-1}-y_{i}) -(y_i-y_{i+1})}{h}=\frac{y_{i-1}-2y_i+y_{i+1}}{h}
$$
represents the $2$nd order derivative.
Now I want to extend this idea to fractional derivative.
The wiki and the paper define the Grünwald–Letnikov derivative:
$$
D^\alpha f(x)= \lim_{h \rightarrow 0} h^{-\alpha}\sum_{n=0}^{\infty} (-1)^n\frac{\Gamma(\alpha +1)}{n! \Gamma(\alpha -n +1)}f(x-nh)
$$
I extract the constant terms into $C$ and have:
$$
D^\alpha f(x) = C \left(\frac{f(x)}{\Gamma(\alpha+1)}-\frac{f(x-h)}{\Gamma(\alpha)} +\frac{f(x-2h)}{2\Gamma(\alpha-1)} -\frac{f(x-3h)}{6\Gamma(\alpha-2)}-\frac{f(x-4h)}{24\Gamma(\alpha-3)}+\dots\right)
$$
where $h$ approaches $0$.
Settings $\alpha = 0.6 $ (between the $1$st and $2$nd order derivative) gives 
$$
D^\alpha f(x) = C \left(\frac{f(x)}{0.98}-\frac{f(x-h)}{1.49} -\frac{f(x-2h)}{2\times 3.72} -\frac{f(x-3h)}{6\times 2.66}-\frac{f(x-4h)}{24\times 1.11}+\dots\right)
$$
Which indicates that we may ignore the terms with large values of $n$ as $n!$ gives rapidly with $n$.
Question: it means that the fractional derivative on a discrete domain is proportional to:
$$
\frac{y_i}{\Gamma(\alpha+1)}-\frac{y_{i+1}}{\Gamma(\alpha)} +\frac{y_{i+2}}{2\Gamma(\alpha-1)} -\frac{y_{i+3}}{6\Gamma(\alpha-3)}+\dots
$$
or
$$
\dots -\frac{y_{i-3}}{6\Gamma(\alpha-3)}+ \frac{y_{i-2}}{2\Gamma(\alpha-1)}-\frac{y_{i-1}}{\Gamma(\alpha)} +\frac{y_i}{\Gamma(\alpha+1)}-\frac{y_{i+1}}{\Gamma(\alpha)} +\frac{y_{i+2}}{2\Gamma(\alpha-1)} -\frac{y_{i+3}}{6\Gamma(\alpha-3)}+\dots
$$
Those terms with $y_{i-3}, y_{i-2}, y_{i-1}$  seem problematic.
 A: The theory of ideal filters
Let $x(t)$ be a smooth fast decreasing signal satisfying the conditions of Shannon's sampling theorem. Let $x[n] = x(n), n \in \mathbb{Z}$ be the sampled signal. Then $$x(t) = \sum_{n=-\infty}^\infty x[n] \frac{\sin( \pi (t-n))}{\pi (t-n)}$$
and
$$x'(t) = \sum_{n=-\infty}^\infty x[n] \frac{d}{dt}\frac{\sin( \pi (t-n))}{\pi (t-n)} = \sum_{n=-\infty}^\infty x'[n] \frac{\sin( \pi (t-n ))}{\pi (t-n)}$$
Where $x'[n]$ is the sampled derivative
$$x'[n] = x'(n) =\sum_{m=-\infty}^\infty x[m] \frac{d}{dt}\frac{\sin( \pi (n-m))}{\pi (n-m)}$$
For the fractional derivative $x^{(a)}(t)$ it works the same way :
$$x^{(a)}(t) = \sum_{n=-\infty}^\infty x[n] \frac{d^a}{dt^a}\frac{\sin( \pi (t-n))}{\pi (t-n)} =\sum_{m=-\infty}^\infty x^{(a)}[m] \frac{\sin( \pi (n-m))}{\pi (n-m)}$$
where $x^{(a)}[n]$ is the sampled fractional derivative
$$x^{(a)}[n] = x^{(a)}(n)=\sum_{n=-\infty}^\infty x[m] \frac{d^a}{dt^a}\frac{\sin( \pi (n-m))}{\pi (n-m)}$$
Note all this works for any convolution $x \ast h(t)$ where $h$ can be any (tempered) distribution. Here $$x'(t) =x \ast \delta'(t), \qquad x^{(a)}(t) = x \ast \delta^{(a)}(t)$$ where $\delta$ is the Dirac delta distribution.
A: Using the Grünwald–Letnikov derivative
$$
D^\alpha f(x)= \lim_{h \rightarrow 0} h^{-\alpha}\sum_{n=0}^{\infty} (-1)^n\frac{\Gamma(\alpha +1)}{n! \Gamma(\alpha -n +1)}f(x-nh)
$$
I derived the coefficients for integer order derivative such as:
$1$st derivative: $ x_{i}-x_{i+1} $
$2$nd derivative: $ x_{i}-2x_{i+1}+x_{i+2} $
and fractional derivative, e.g., $\alpha=1.5$
$$
x_{i}-1.5\;x_{i+1}+0.375\;x_{i+2}-0.0625\;x_{i+3}+\dots
$$
D.A. Benson's paper 

Fractional calculus in hydrologic modeling

shows a comprehensive figure

that facilities numerical implementation of fractional derivative.
