How to rigorously understand the idea behind Laplace method? Let us consider the integral $$I(\lambda)=\int_a^b f(x) e^{-\lambda g(x)} dx\tag{1}$$ where $f$ and $g$ are smooth functions. The Laplace method (sketched, e.g. here) allows us to take $\lambda \to \infty$ in this expression. The idea seems to be that in the limit $\lambda\to \infty$ only neighborhoods of minima of $g$ will contribute. This is the idea I want to understand from a rigorous point of view.
In other words, suppose that $g$ has just one minimum $c\in[a,b]$ with $g'(c)=0$ and $g''(c)>0$. In that case one approximates $$I(\lambda)\approx \int_{c-\epsilon}^{c+\epsilon} f(x) e^{-\lambda g(x)}dx\tag{2}$$
and since this is one small interval around $c$ one Taylor expands $f$ and $g$. 
Now I want to rigorously understand (2). I mean, one is clearly separating $$I(\lambda)=\int_{(a,b)\setminus (c-\epsilon,c+\epsilon)}f(x)e^{-\lambda g(x)}dx+\int_{c-\epsilon}^{c+\epsilon}f(x)e^{-\lambda g(x)}dx\tag{3}$$
and saying that we can throw the first term away.
But how to understand the accuracy of the approximation with respect to $\epsilon$? Moreover, why the $\lambda \to +\infty$ is needed to do this, and how this limit relates to the choice of $\epsilon$? 
In other words what is the rigorous way to make sense of the asymptotic approximation (2)?
 A: *

*$\lambda \rightarrow \infty $. That is how Laplace Transform is defined. The Laplace Transform of $f\left ( t \right )$, denoted by $\textit{F}\left ( s \right )=\mathfrak{L}\left \{ f\left ( t \right ) \right \}$, is defined by the integral $F\left ( s \right )=\int_{0}^{\infty }e^{-st}f\left ( t \right )dt$, where $f\left ( t \right )$ is defined for $t\geq 0$. This unilateral Laplace Transform is what is meant by "the" Laplace Transform most of the time, although sometimes the bilateral form is also used, $F\left ( s \right )=\int_{-\infty }^{\infty }e^{-st}f\left ( t \right )dt$.

*$I\left ( \lambda  \right )=\int_{a}^{b}f\left ( t \right )e^{-\lambda g\left ( t \right )}dt\approx e^{-\lambda g\left ( c \right )}\int_{c-\epsilon }^{c+\epsilon }f\left ( t \right )e^{-\lambda \left [ g\left ( t \right )-g\left ( c \right ) \right ]}dt$ for large $\lambda$.
If the real continuous function $g\left ( t \right )$ assumes strict minimum on the interval $\left [ a,b \right ]$ at $c \in \left [ a,b \right ]$, then it is only the immediate neighborhood of $c$ that contributes to the asymptotic expansion of $I\left ( \lambda  \right )$ for large $\lambda$. Furthermore, as $g\left ( t \right )$ assumes its strict minimum on the interval, $-g\left ( t \right )$ assumes its strict maximum on the interval.
$I\left ( \lambda  \right )=\int_{a}^{b}f\left ( t \right )e^{-\lambda g\left ( t \right )}dt=\int_{a}^{b}f\left ( t \right )e^{\left(\phi \left ( t \right )\right)\lambda }dt$ where $\phi \left ( t \right )=-g\left ( t \right )$
Write $e^{\phi \left ( t \right )}$ as $e^{\phi \left ( c \right )}e^{\psi \left ( t \right )}$, where $e^{\phi \left ( c \right )}$ is a multiplication constant that can be taken out of the integral. Immediately, it follows that $\psi \left ( t \right )=\phi \left ( t \right )-\phi \left ( c \right )$. The maximum value of $ \psi \left ( t \right )$ is $0$, that is because the largest $\phi \left ( t \right )$ can be is $\phi \left ( c \right )$. The maximum value for $e^{\psi \left ( t \right )}$ is therefore $1$.
Following is a sketch for the typical behavior of $e^{\psi \left ( t \right )}$. Recall that $$\left.\begin{matrix}\phi \left ( t \right )=-g\left ( t \right )\\ \psi \left ( t \right )=\phi \left ( t \right )-\phi \left ( c \right )\end{matrix}\right\}\rightarrow \psi \left ( t \right )=-g\left ( t \right )-\left ( -g\left ( c \right ) \right )=-g\left ( t \right )+g\left ( c \right )\\ \psi \left ( t \right )=-\left [ g\left ( t \right )-g\left ( c \right ) \right ]$$
Hence, this graph is also the typical behavior for $e^{-\left [ g\left ( t \right )-g\left ( c \right ) \right ]}$

As $\lambda$ increases, the maximum for $e^{\psi \left ( t \right )}$ stays at $\left ( c,1 \right )$. However, the wings of the graphs move closer and closer to the vertical line $t=c$.Thus, $\phi \left ( t \right )$ (or $-g\left ( t \right )$) can be replaced with its approximation that only needs to be good in the immediate vicinity of $c$, where $\phi \left ( t \right )$ assumes it maximum (or equivalently, where $g\left ( t \right )$ assumes its minimum.)
Put everything together, we therefore have Expression (2) from the lecture notes you included in your question.
$$I\left ( \lambda  \right )=\int_{a}^{b}f\left ( t \right )e^{\left ( \phi \left ( t \right ) \right )\lambda }\approx e^{\phi \left ( c \right )}\int_{c-\epsilon }^{c+\epsilon }f\left ( t \right )e^{\left ( \phi \left ( t \right ) \right )\lambda }dt=e^{-\lambda g\left ( c \right )}\int_{c-\epsilon }^{c+\epsilon }f\left ( t \right )e^{-\lambda \left [ g\left ( t \right )-g\left ( c \right ) \right ]}dt$$
