I have an analytic function defined by a power series $$ f(x):=\sum_{l=0}^\infty f_l x^l, $$ where $f_l \ge 0$ for all $l \in \mathbb N$.
Assume that $f(x)$ has radius of convergence $x_0 > 0$. I want to claim that, for all fixed $r \ge 1$, there exists $C_r>0$ such that $$ R_r(x):=\left|f(x)-\sum_{l=0}^{r-1} f_l x^l\right| \le C_r |x^r|, \qquad (|x| < x_0/2). $$
My argument is that $f^{(r)}(x)$, the $r$-th derivative of $f(x)$, has the same radius of convergence as $f(x)$. So $f^{(r)}(x)$ is bounded when $|x| < x_0/2$, so we can apply Taylor's theorem to bound the remainder. Is this correct?
I faced this problem many years ago, while working to my "Laurea" degree thesis: the solution I found is reported below.
The radius of convergence of $f$ is $x_0$: by Cauchy-Hadamard theorem we have that the series $f(x)$ converges absolutely for every $x\in \mathbb{C}$ such that $|x|< x_0$, i.e.
$$ |f|(x)=\sum_{l=0}^\infty |f_l||x|^l<\infty\quad\forall x\in\mathbb{D}(x_0)\tag{1}\label{1} $$ where $\mathbb{D}(x_0)=\{x\in\mathbb{C}||x|<x_0\}$. This implies that $$ \lim_{l\to\infty}|f_l||x|^l=0 \implies \exists\max_{l\in\mathbb{N}}|f_l||x|^l<\infty\quad\forall x\in\mathbb{D}(x_0) $$ since all the terms of \eqref{1} are $\geq 0$ and their $\limsup$ is monotonically decreasing to $0$. Therefore we can put $$ \max_{l\in\mathbb{N}}|f_l||x|^l=|f_k||x|^k\quad x\in\mathbb{D},\; k=k(x)\in\mathbb{N}\tag{2}\label{2} $$ Now, for each real $\delta$ such that $0<\delta<x_0/2$ we have $x_0/2<x_0-\delta<x_0$: this implies that we can use \eqref{2} for $x=x_0-\delta$ (and correspondingly $k=k(x_0-\delta$) to estimate each term of $f$ for each $x\in\mathbb{D}(x_0/2)$. By using this estimate and proceeding in the evaluation of $R_r(x)$, we obtain $$ \begin{split} R_r(x)=\left|\,f(x)-\sum_{l=0}^{r-1} f_l x^l\right|&=\left|\sum_{l=r}^\infty f_l x^l\right|\leq\sum_{l=r}^\infty |f_l| |x|^l\\ &=\sum_{l=r}^\infty |f_l| |x_0-\delta|^l\left|\frac{x}{x_0-\delta}\right|^l\\ &\leq |f_k|(x_0-\delta)^k\sum_{l=r}^\infty \left|\frac{x}{x_0-\delta}\right|^l=\frac{|f_k|(x_0-\delta)^{k-r+1}}{(x_0-\delta-x)}|x|^{r}\\ &\leq 2\frac{|f_k|(x_0-\delta)^{k-r+1}}{(x_0-2\delta)}|x|^{r}=C_r(x_0-\delta)|x|^{r}\end{split} \tag{3}\label{3} $$ where $C_r$ is a constant with does not depend on $x$.
On the right side you have the sought for polynomial (in fact monomial) estimate of the remainder of power series $f(x)$.
Notes
I changed a bit my answer according to the comments of @ablmf, and changed accordingly even the following comments