Probability proof of inversion formula for Laplace transform Let $f:[0, \infty[\longrightarrow \mathbb{R}$ be  bounded and continuous and define $L(\lambda)=\int_0^\infty e^{-\lambda x}f(x)dx$.
Let $X_n$ be a sequence of independent random variables with exponential distribution of rate $\lambda$.
Using the fact that the sum $S_n=X_1+...+X_n$ has Gamma distribution we can see that
$$  (-1)^{n-1}\frac{\lambda^nL^{(n-1)}(\lambda)}{(n-1)!} =Ef(S_n)$$
where $L^{(n-1)}$ is the $n-1$ derivative of $L$.
How can we use this to 
prove $$ f(y)=\lim_{n}(-1)^{n-1}\frac{(\frac{n}{y})^nL^{(n-1)}(\frac{n}{y})}{(n-1)!} $$
We know that by the strong law of large numbers $S_n/n$ converges to $E(X_1)$ almost everywhere, so I could consider 
the parameter $\lambda=f(y)$ but that doesn't seem useful..
 A: Try the Weierstrass Approximation Theorem.
Define $B_n(y) := E[f(S_n)](n/y)$, $Y_n := |f(S_n) - f(y)|$ and $Z_n := |S_n - y|$
Since f is bounded and continuous on $[0,\infty)$, it is bounded and continuous on $[0,M] \ \forall M > 0$ and hence uniformly continuous in $[0,M]$.
$\to \forall \epsilon > 0, \exists \delta > 0 \forall x, y \in [0,M]$ s.t.
$$Z_n \le \delta \to Y_n < \epsilon/2$$
Let us try to show that $\forall \epsilon > 0, \exists N > 0$ s.t. $n > N \ \to \ |B_n(y) - f(y)| < \epsilon$:
$$|B_n(y) - f(y)|$$
$$ = |E[f(S_n)](n/y) - f(y)|$$
$$ = |E[f(S_n)](\lambda)|_{\lambda = n/y} - f(y)|$$
$$ = |(E[f(S_n)] - f(y))(\lambda)|_{\lambda = n/y}|$$
$$ = |(E[f(S_n) - f(y)])(\lambda)|_{\lambda = n/y}|$$
$$ = (|E[f(S_n) - f(y)]|)(\lambda)|_{\lambda = n/y}$$
By Jensen's Inequality,
$$ \le (E[|f(S_n) - f(y)|])(\lambda)|_{\lambda = n/y}$$
Omitting $\lambda$ for now
$$ = E[|f(S_n) - f(y)|]$$
$$ = E[Y_n]$$
$$ = E[Y_n 1_{Z_n \le \delta} + Y_n 1_{Z_n > \delta}]$$
$$ = E[Y_n 1_{Z_n \le \delta}] + E[Y_n 1_{Z_n > \delta}]$$
$$ < E[\epsilon/2 1_{Z_n \le \delta}] + E[Y_n 1_{Z_n > \delta}]$$
$$ = \epsilon/2 E[1_{Z_n \le \delta}] + E[Y_n 1_{Z_n > \delta}]$$
$$ = \epsilon/2 P(Z_n \le \delta) + E[Y_n 1_{Z_n > \delta}]$$
$$ = \epsilon/2 (1) + E[Y_n 1_{Z_n > \delta}]$$
$$ = \epsilon/2 + E[Y_n 1_{Z_n > \delta}]$$
$$ \le \epsilon/2 + E[2K 1_{Z_n > \delta}]$$

Note: $\exists K > 0$ s.t. $|f(x)| \le K \forall x \in [0, M]$
$\to \ Y_n = |f(S_n) - f(y)| \le |f(S_n)| + |f(y)| \le K + K = 2K$

$$ = \epsilon/2 + 2K E[ 1_{Z_n > \delta}]$$
$$ = \epsilon/2 + 2K P(Z_n > \delta)$$

By Chebyshev's Inequality, $\forall \delta > 0, P(Z_n > \delta) = P(|S_n -
 \frac{n}{\lambda}| > \delta) = P(|S_n - \frac{n}{\lambda}| > n\delta)$
$= P(|\frac{S_n}{n} - \frac{1}{\lambda}| > \delta) \le \frac{1}{n
 \lambda^2 \delta^2}$

$$ \le \epsilon/2 + 2K \frac{1}{n \lambda^2 \delta^2}$$
$$\le \epsilon/2 + \epsilon/2 = \epsilon$$
if we choose $\color{red}{N}$ as such:
$$2K \frac{1}{n \lambda^2 \delta^2} < \epsilon/2$$
$$\frac{1}{n \lambda^2 \delta^2} < \frac{\epsilon}{4K}$$
$$\frac{1}{n} < \frac{\lambda^2 \delta^2 \epsilon}{4K}$$
$$n > \frac{4K}{\lambda^2 \delta^2 \epsilon} \color{red}{= N}$$
Not sure what to do about the $\lambda$, but it seems that $\forall \epsilon > 0$,
$$n > \frac{4K}{\lambda^2 \delta^2 \epsilon} \to Y_n < \epsilon/2$$
$$\to |B_n(y) - f(y)| < \epsilon \ QED$$
$$\therefore \ \lim_{n \to \infty} B_n(y) = f(y)$$
A: Old question, but I was working on the same problem recently. Here's an idea:
The functions $g_n$ defined by $g(x) = f(x/n)$ are bounded and continuous. Let $L_h$ denote the Laplace transform of a function $h$, and let the independent exponential random variables $X_1,X_2,\dots$ have rate $\frac{1}{y}$ (using $y\neq 0$). By the first part of the problem
$$
\mathbb{E}\left(f\left(\frac{S_n}{n}\right)\right) = \mathbb{E}(g_n(S_n)) = (-1)^{n-1}\frac{\frac{1}{y^n}L^{(n-1)}_{g_n}\left(\frac{1}{y}\right)}{(n-1)!}
$$
Yet the strong law from Chapter 7 gives
$$\lim_{n\to\infty}\frac{S_n}{n} = \mathbb{E}(X_1) = y \quad \text{almost surely}$$
and since $f$ is bounded the expected values converge:
$$
\lim_{n\to\infty} (-1)^{n-1}\frac{\frac{1}{y^n}L^{(n-1)}_{g_n}\left(\frac{1}{y}\right)}{(n-1)!} = \lim_{n\to\infty} \mathbb{E}\left(f\left(\frac{S_n}{n}\right)\right) = f(y)
$$
The rest is just a matter of rewriting the terms $L_{g_n}^{(n-1)}$ in terms of $L_f^{(n-1)}$. For any bounded continuous (say) function $h$ one can compute for each $\lambda > 0$ that
$$
\quad L_h^{(n-1)}(\lambda) = (-1)^{n-1}\int_0^\infty e^{-\lambda x} x^{n-1} h(x)\,dx
$$
In our case, then,
$$
L_{g_n}^{(n-1)}(\lambda) = (-1)^{n-1}\int_0^\infty e^{-\lambda x} x^{n-1} f(x/n) \, \mathrel{\underset{u:=x/n}{=}} (-1)^{n-1}\int_0^\infty e^{-n\lambda u} n^{n-1}u f(u) n \, du = n^n L_f(n\lambda)
$$
Setting $\lambda = \frac{1}{y}$ yields
$$
f(y) = \lim_{n\to\infty} (-1)^{n-1}\frac{\frac{1}{y^n}L^{(n-1)}_{g_n}\left(\frac{1}{y}\right)}{(n-1)!} = \lim_{n\to\infty} (-1)^{n-1}\frac{\frac{n^n}{y^n}L^{(n-1)}_{f}\left(\frac{n}{y}\right)}{(n-1)!}
$$
