As a general remark, we may see the "Exponential to Ordinary" transformation of generating functions, $$ f(x):=\sum_{r=0}^\infty a_rx^r/r!\mapsto \tilde f(t):=\sum_{r=0}^\infty a_rt^r, $$ as an operator $\mathbb{C}[[x]]\to \mathbb{C}[[t]]$ . Since $r!t^r=\int_0^\infty (tx)^r e^{-x}dx$, the transformation can be computed analytically as $$\tilde f(t)=\int_0^\infty f(tx)e^{-x}dx,$$ at least for suitably convergent $f(x)$, and for special values of $t$. If the RHS is a convergent series $\sum_{r=0}^\infty b_rt^r$, and the equality holds, for a set of values $t$ which accumulates within the disk of convergence, the identity of series $ \sum_{r=0}^\infty a_rt^r=\sum_{r=0}^\infty b_rt^r$ is then established by the principle of isolated zeros.

For instance, we may compute the transform of $f_r(x):=(e^x-1)^r$ for real negative values of $t$ in terms of the Euler's Beta function by a change of variable in the integral:
$$\tilde f_r(t)=\int_0^\infty (e^{tx}-1)^re^{-x}dx=(-1)^{r+1}t^{-1}\int_0^1(1-u)^r u^{-1/t-1}du=$$ $$=(-1)^{r+1}t^{-1}\frac{\Gamma(r+1)\Gamma(-1/t)}{\Gamma(r+1-1/t)}=\frac{r!t^r}{(1-t)\dots(1-rt)}\ .$$ This computation gives your identity, since for the egf $f(x)$ of the $B_r$'s we have $f(x)=e^{e^x-1}=\sum_{r=0}^\infty\frac{1}{r!}f_r(x)$ (in the sense of formal power series) so that $\tilde f(t)=\sum_{r=0}^\infty \frac{t^r}{(1-t)\dots(1-rt)}\ .$

Incidentally, note that by an analogous computation you may, more generally, compute the ogf of the Stirling polynomials of the second kind, $B_r(z)$, starting by their egf $e^{z(e^x-1)}$.

As a general remark, we may see the "Exponential to Ordinary" transformation of generating functions, $$ f(x):=\sum_{r=0}^\infty a_rx^r/r!\mapsto \tilde f(t):=\sum_{r=0}^\infty a_rt^r, $$ as an operator $\mathbb{C}[[x]]\to \mathbb{C}[[t]]$ . Since $r!t^r=\int_0^\infty (tx)^r e^{-x}dx$, the transformation can be computed analytically as $$\tilde f(t)=\int_0^\infty f(tx)e^{-x}dx,$$ at least for suitably convergent $f(x)$, and for special values of $t$. If the RHS is a convergent series $\sum_{r=0}^\infty b_rt^r$, and the equality holds, for a set of values $t$ which accumulates within the disk of convergence, the identity of series $ \sum_{r=0}^\infty a_rt^r=\sum_{r=0}^\infty b_rt^r$ is then established by the principle of isolated zeros.

For instance, we may compute the transform of $f_r(x):=(e^x-1)^r$ for real negative values of $t$ in terms of the Euler's Beta function by a change of variable in the integral:
$$\tilde f_r(t)=\int_0^\infty (e^{tx}-1)^re^{-x}dx=(-1)^{r+1}t^{-1}\int_0^1(1-u)^r u^{-1/t-1}du=$$ $$=(-1)^{r+1}t^{-1}\frac{\Gamma(r+1)\Gamma(-1/t)}{\Gamma(r+1-1/t)}=\frac{r!t^r}{(1-t)\dots(1-rt)}\ .$$ This computation gives your identity, since for the egf $f(x)$ of the $B_r$'s we have $f(x)=e^{e^x-1}=\sum_{r=0}^\infty\frac{1}{r!}f_r(x)$ (in the sense of formal power series) so that $\tilde f(t)=\sum_{r=0}^\infty \frac{t^r}{(1-t)\dots(1-rt)}\ .$

Incidentally, note that by an analogous computation you may, more generally, compute the ordinary gf of the Stirling polynomials of the second kind, $B_r(z) :=\sum_{r=0}^nS(n,r)z^r$, starting by their egf $e^{z(e^x-1)}$.

rmk. I think a more natural proof is, computing directly the ogf of the Stirling polynomials $B_r(z):=\sum_{r=0}^nS(n,r)z^r$, which consists just in translating the inductive relation $S(n+1,k)=kS(n,k)+S(n,k−1)$ in the language of formal power series, and then putting $z=1$, as done in other answers. The present answer is meant as an example of a general alternative approach.

As a general remark, we may see the "Exponential to Ordinary" transformation of generating functions, $$ f(x):=\sum_{r=0}^\infty a_kx^k/k!\mapsto \tilde f(t):=\sum_{r=0}^\infty a_kt^k, $$ as an operator $\mathbb{C}[[x]]\to \mathbb{C}[[t]]$ . Since $n!t^n=\int_0^\infty (tx)^n e^{-x}dx$, the transformation can be computed analytically as $$\tilde f(t)=\int_0^\infty f(tx)e^{-x}dx,$$ at least for suitably convergent $f(x)$, and for special values of $t$. If the RHS is a convergent series $\sum_{r=0}^\infty b_kt^k$, and the equality holds, for a set of values $t$ which accumulates within the disk of convergence, the identity of series $ \sum_{r=0}^\infty a_kt^k=\sum_{r=0}^\infty b_kt^k$ is then established by the principle of isolated zeros.

For instance, we may compute the transform of $f_r(x):=(e^x-1)^r$ for real negative values of $t$ in terms of the Euler's Beta function by a change of variable in the integral:
$$\tilde f_r(t)=\int_0^\infty (e^{tx}-1)^re^{-x}dx=(-1)^{r+1}t^{-1}\int_0^1(1-u)^n u^{-1/t-1}du=$$ $$=(-1)^{r+1}t^{-1}\frac{\Gamma(r+1)\Gamma(-1/t)}{\Gamma(r+1-1/t)}=\frac{r!}{(1-t)\dots(1-rt)}\ .$$ This computation gives your identity, since for the egf $f(x)$ of the $B_r$'s we have $f(x)=e^{e^x-1}=\sum_{r=0}^\infty\frac{1}{r!}f_r(x)$ (in the sense of formal power series) so that $\tilde f(t)=\sum_{r=0}^\infty \frac{t^n}{(1-t)\dots(1-rt)}\ .$

Incidentally, note that by an analogous computation you may, more generally, compute the ogf of the Stirling polynomials of the second kind, $B_r(z)$, starting by their egf $e^{z(e^x-1)}$.

As a general remark, we may see the "Exponential to Ordinary" transformation of generating functions, $$ f(x):=\sum_{r=0}^\infty a_rx^r/r!\mapsto \tilde f(t):=\sum_{r=0}^\infty a_rt^r, $$ as an operator $\mathbb{C}[[x]]\to \mathbb{C}[[t]]$ . Since $r!t^r=\int_0^\infty (tx)^r e^{-x}dx$, the transformation can be computed analytically as $$\tilde f(t)=\int_0^\infty f(tx)e^{-x}dx,$$ at least for suitably convergent $f(x)$, and for special values of $t$. If the RHS is a convergent series $\sum_{r=0}^\infty b_rt^r$, and the equality holds, for a set of values $t$ which accumulates within the disk of convergence, the identity of series $ \sum_{r=0}^\infty a_rt^r=\sum_{r=0}^\infty b_rt^r$ is then established by the principle of isolated zeros.

For instance, we may compute the transform of $f_r(x):=(e^x-1)^r$ for real negative values of $t$ in terms of the Euler's Beta function by a change of variable in the integral:
$$\tilde f_r(t)=\int_0^\infty (e^{tx}-1)^re^{-x}dx=(-1)^{r+1}t^{-1}\int_0^1(1-u)^r u^{-1/t-1}du=$$ $$=(-1)^{r+1}t^{-1}\frac{\Gamma(r+1)\Gamma(-1/t)}{\Gamma(r+1-1/t)}=\frac{r!t^r}{(1-t)\dots(1-rt)}\ .$$ This computation gives your identity, since for the egf $f(x)$ of the $B_r$'s we have $f(x)=e^{e^x-1}=\sum_{r=0}^\infty\frac{1}{r!}f_r(x)$ (in the sense of formal power series) so that $\tilde f(t)=\sum_{r=0}^\infty \frac{t^r}{(1-t)\dots(1-rt)}\ .$

Incidentally, note that by an analogous computation you may, more generally, compute the ogf of the Stirling polynomials of the second kind, $B_r(z)$, starting by their egf $e^{z(e^x-1)}$.