Normal approximation of tail probability in binomial distribution My problem: From the Berry--Esseen theorem I know, that $$\sup_{x\in\mathbb R}|P(B_n \le x)-\Phi(x)|=O\left(\frac 1{\sqrt n}\right),$$ where $B_n$ has the standardized binomial distribution and $\Phi$ is the standard normal distribution function. I can prove this for $x\approx 0$ with Stirling's formula and in a way similar to this proof.
Unfortunately Stirling's approximation becomes worse the bigger $|x|$ is, so that I can only show that
$$\sup_{|x| \le c}|P(B_n \le x)-\Phi(x)|=O\left(\frac 1{\sqrt n}\right)$$ for a fixed $c > 0$.
My question: Is there a direct proof for $\sup_{|x| > c}|P(B_n \le x)-\Phi(x)|=O\left(\frac 1{\sqrt n}\right)$, based on Stirling's formula? How can I estimate the error of the difference in the tail probability of the binomial distribution to the normal distribution?
Note: Chebyshev's inequality does not provide the right convergence speed. As zhoraster shows in his answer to Finding an error estimation for the De Moivre–Laplace theorem, also Chernoff's inequality does not provide the right convergence speed, too. That's why I am looking for another method. I also want to use more easy and direct approximations than Stein's method which is used in the proof of the Berry--Esseen theorem.
I already asked this question on MSE (see this thread), but I haven't got an answer there (even after adding a bounty). That's why I want to ask this question here. 
 A: The needed observation here is that, whereas the relative error of the Stirling approximation is worse for large $|x|$, it gets multiplied by a fast decreasing normal density function, and then the product gets integrated to produce the desired result. Here are the details. 
Let $X$ be a binomial random variable with parameters $n\in\mathbb{N}$ and $p\in(0,1)$, where $n\to\infty$ and $p$ is fixed. For any integer $k$, let 
$$z_k:=\frac{k-np}{\sqrt{npq}},\quad q:=1-p. 
$$
Let 
$$K:=\{k=0,\dots,n\colon|z_k|\le n^{1/4}\},$$
so that $k/n\to p$ uniformly over $k\in K$. 
Let 
$$\hat p:=k/n,\quad\hat q:=1-\hat p=(n-k)/n,$$
$$\epsilon:=\epsilon_k:=z_k/\sqrt n,\quad\lambda:=\sqrt{q/p}, 
$$
whence 
$$\frac{\hat p}{p}=1+\epsilon\lambda,\quad\frac{\hat q}{q}=1-\epsilon/\lambda.$$
Uniformly over $k\in K$ one has $k\sim np$ and $n-k\sim nq$, so that $\frac1k=O(\frac1n)$ and $\frac1{n-k}=O(\frac1n)$. 
Here and in what follows, the constant in $O(\cdot)$ depends only on $p$. So, by Stirling's formula, over all $k\in K$ 
$$P(X=k)=\frac{n!}{k!(n-k)!} p^kq^{n-k}
=\frac1{\sqrt{2\pi n\hat{p}\hat{q}}}\,\frac1R\,\Big(1+O\Big(\frac1n\Big)\Big)$$
$$
=\frac1{\sqrt{2\pi npq}}\,\frac1R\,\Big(1+O\Big(|\epsilon|+\frac1n\Big)\Big)
=\frac\delta{\sqrt{2\pi}}\,e^{-nH}\Big(1+O\Big(|\epsilon|+\frac1n\Big)\Big)
$$
where
$$R:=\Big(\frac{\hat p}{p}\Big)^k \Big(\frac{\hat q}{q}\Big)^{n-k},  
$$
$$\delta:=\frac1{\sqrt{npq}},\quad H:=H_k:=\ln(R^{1/n})=\hat p\ln\frac{\hat p}p+\hat q\ln\frac{\hat q}q
:=p\psi(\epsilon\lambda)+q\psi(-\epsilon/\lambda),$$
$$\psi(s):=(1+s)\ln(1+s).
$$ 
Since $\psi(s)=s + s^2/2 +O(|s|^3)$ for $|s|\le1/2$ and 
$\frac{|z_k|^3}{\sqrt n}+|\epsilon|=\frac{|z_k|^3}{\sqrt n}+\frac{|z_k|}{\sqrt n}=O\big(\frac{|z_k|^3+1}{\sqrt n}\big)$, we find that
$$(1)\qquad P(X=k)=\delta\,\phi(z_k)\Big(1+O\Big(\frac{|z_k|^3}{\sqrt n}\Big)\Big)\Big(1+O\Big(|\epsilon|+\frac1n\Big)\Big)$$
$$
=\delta\,\phi(z_k)\Big(1+O\Big(\frac{|z_k|^3+1}{\sqrt n}\Big)\Big)  
$$
over all $k\in K$, 
where $\phi$ is the standard normal density. 
Next, 
for $u=O(n^{1/4})$ and $|h|\le\delta/2$, 
$$\frac{\phi(u+h)+\phi(u-h)}{2}=\phi(u)\cosh(hu)e^{-h^2/2}
=\phi(u)(1+O(h^2u^2+h^2))=\phi(u)\Big(1+O\Big(\frac1{\sqrt n}\Big)\Big). 
$$
So, for $k\in K$
$$\delta\,\phi(z_k)=\frac12\,\int_{-\delta/2}^{\delta/2}[\phi(z_k+h)+\phi(z_k-h)]\,dh\,\Big(1+O\Big(\frac1{\sqrt n}\Big)\Big)$$
$$
=\int_{z_k-\delta/2}^{z_k+\delta/2}\phi(z)\,dz\,\Big(1+O\Big(\frac1{\sqrt n}\Big)\Big). 
$$
Also, 
for any $k\in K$ and any $z\in[z_k-\delta/2,z_k+\delta/2]$ one has 
$|z_k|^3-|z|^3\le3\delta(|z_k|+\delta/2)^2=O(1)$, whence $|z_k|^3/\sqrt n=O((|z|^3+1)/\sqrt n)$. 
So, by $(1)$, for some positive real $c$ depending only on $p$ and for all $k\in K$ 
$${J_k\!}^-\le P(X=k)\le J_k^+,$$
where
$$J_k^\pm:=\int_{z_k-\delta/2}^{z_k+\delta/2}\phi(z)
\Big(1\pm\frac{c(|z|^3+1)}{\sqrt n}\Big)\,dz.  
$$
Hence, for any $x\in[0,n^{1/4}]$, letting $K_x:=\{k\in\mathbb{Z}\colon0<z_k\le x\}$, one has 
$$P(0<B_n\le x)=\sum_{k\in K_x}P(X=k)
\le\sum_{k\in K_x}J_k^+$$
$$
\le\int_{-\delta/2}^{x+\delta/2}\phi(z)\Big(1+\frac{c(|z|^3+1)}{\sqrt n}\Big)\,dz$$
$$
\le\int_0^x\phi(z)\,dz+A\delta+\frac B{\sqrt n}
\le\int_0^x\phi(z)\,dz+\frac C{\sqrt n}, 
$$ 
where $A,B,C$ depend only on $p$. Quite similarly one obtains the corresponding lower bound on $P(0<B_n\le x)$, so that 
$$(2)\qquad |P(0<B_n\le x)-(\Phi(x)-1/2)|\le\frac C{\sqrt n}  
$$
for $x\in[0,n^{1/4}]$, where $\Phi$ is the standard normal cumulative distribution function. 
On the other hand, for real $x>n^{1/4}$ one has $1-\Phi(x)=O(1/\sqrt n)$ and, 
by Bernstein's inequality, 
$$P(B_n>x)\le\exp\Big(-\frac{x^2}{2+x/\sqrt{npq}}\Big)
\le e^{-x^2/4}+e^{-x\sqrt{npq}/2}=O(1/\sqrt n). 
$$
It follows that $(2)$ holds for all real $x\ge0$ (perhaps with a greater constant $C$). Similarly, $|P(B_n\le0)-1/2)|\le\frac C{\sqrt n}$. 
Thus, for all real $x$
$$|P(B_n\le x)-\Phi(x)|\le\frac{2C}{\sqrt n},   
$$
as desired. 
