Let's say I have a sequence $s_n$ of numbers, and I want a series $a_i$ which computes the sequence; that is
$\sum_{i=0}^\infty a_i n^i = s_n$
Clearly $a_0 = s_0$, but after that I am stuck. I need $\sum_{i=0}^\infty a_i = s_1$, so I'm not sure how I should go about picking the $a_1$ -- it seems that I need to choose all the $a_i$ "simultaneously".
Are there any methods for doing this? I am most interested in computable methods, though I suspect it may not be a computable problem in general. I am interested in theoretical results about this as well.
More generally, for any sequences $p_j$ and $s_j$ of complex numbers such that $p_j$ has no limit points, there exist entire functions $f$ such that $f(p_j) = s_j$ for all $j$. See e.g. Rudin, "Real and Complex Analysis", Theorem 15.15. You then want the Maclaurin series of $f$.
EDIT: One way to construct such an $f$ in your case is as follows. Let $ g_n(z)$ be an entire function with $g_n(n) = (-1)^n s_n/\pi$ and $|g_n(z)| \le 1/n$ for $|z| \le n/2$. For example, you might take $g_0(z) = s_1/\pi$ and $g_n(z) = (-1)^n s_n (z/n)^{k_n} \pi$ for $n \ge 1$, where $k_n > \log_2(n |a_n|/\pi)$. Let $g(z) = \sum_{n=0}^\infty \frac{g_n(z)}{z-n}$. This series converges to a meromorphic function with simple poles at the nonnegative integers, having residue $(-1)^n s_n/\pi$ at $z=n$. Then $ f(z) = \sin(\pi z) g(z)$ is entire with $f(n) = s_n$ for all nonnegative integers $n$.
Note that to compute $a_j$ we need only consider the $g_n$ where $k_n \le j$, and we can arrange it so $k_n \ge n$.