Newton Forward differentiation

function NFD = newton_forward_differentiation(x, y)
    n = length(x);
    diff_table = zeros(n, n);

    for i = 1 : n
        diff_table(i, 1) = y(i);
    end

    for j = 2 : n
        for i = 1 : n - j + 1
            diff_table(i, j) = diff_table(i + 1, j - 1) - diff_table(i, j - 1);
        end
    end

    NFD = 0;
    for i = 2 : n
        odd = mod(i, 2);

        if odd == 1
            NFD = NFD - diff_table(1, i) / (i - 1);
        else
            NFD = NFD + diff_table(1, i) / (i - 1);
        end

    end

    NFD = NFD / (x(2) - x(1));