To get the maximum error we need to find the maximum of
$$\left|\prod_{i=0}^n\left(x-x_i\right)\right|=\left|\prod_{i=0}^n\left(x-x_0-ih\right)\right|=\left|h^{n+1}\prod_{i=0}^n\left(\frac{x-x_0}h-i\right)\right|=\left|h^{n+1}\prod_{i=0}^n\left(t-i\right)\right|$$
Where $t=\frac{x-x_0}h$ and $0\le t\le n$. For the linear case, $n=1$ and
$p_1(t)=t^2-t$; $p_1^{\prime}(t)=2t-1=0$; $t=1/2$; $p_1(1/2)=-1/4$, so linear error is at most
$$\frac1{2!}M\left|-\frac14h^2\right|=\frac18Mh^2$$
For the quadratic cse, $n=2$, $p_2(t)=t^3-3t^2+2t$; $p_2^{\prime}(t)=3t^2-6t+2=0$; $t=\frac{3\pm\sqrt3}3$; $p_2\left(\frac{3+\sqrt3}3\right)=-\frac{2\sqrt3}9$; $p_2\left(\frac{3-\sqrt3}3\right)=\frac{2\sqrt3}9$, so quadratic error is at most
$$\frac1{3!}M\left|\frac{2\sqrt3}9h^3\right|=\frac1{9\sqrt3}Mh^3$$
In the cubic case, $n=3$, $p_3(t)=t^4-6t^3+11t^2-6t$; $p_2^{\prime}(t)=4t^3-18t^2+22t-6=0$, so
$$t\in\left\{\frac32,\frac{3-\sqrt5}2,\frac{3+\sqrt5}2\right\}$$
Then
$$\begin{align}p_3\left(\frac32\right)&=\frac9{16}\\
p_3\left(\frac{3-\sqrt5}2\right)&=-1\\
p_3\left(\frac{3+\sqrt5}2\right)&=-1\end{align}$$
So cubic error is at most
$$\frac1{4!}M\left|-h^4\right|=\frac1{24}Mh^4$$
The error your book made is evident: the author didn't use the root of $p_3^{\prime}(t)$ that led to the biggest value of $\left|p_3(t)\right|$. To make this more clear, let's look at the graph of $p_3(x)=x(x-1)(x-2)(x-3)$:
As can be seen, the biggest values of $\left|p_3(x)\right|$ happen at the two outer local extrema, leading to the largest estimate for the error. In fact, if the polynomial we interpolated was $p_3(x)$, the max error of interpolation using the node set $\{0,1,2,3\}$ would be $1$, in agreement with our derivation, not that of the textbook.
This is an aspect of the Runge phenomenon: this polynomial oscillates wildly, with the biggest swings at the outermost local extrema. See what it looks like for $p_7(x)$:
As polynomials of degrees up to $n$ are interpolated exactly, you know from interpolating $p(x)=1$ that
$$
1=\sum_{k=0}^nL_k(x)~\text{ for all }x\in\Bbb R.
$$
Now integrate this over $[a,b]$.
Best Answer
Use the interpolation polynomial on $f(x)=1$