Working in dimensional regularization, in the MS bar scheme consider the renormalization group equation for the strong coupling
$$\mu\frac{d\alpha}{d\mu}=-\beta_0\alpha(\mu)^2-\beta_1\alpha(\mu)^3-\beta_2\alpha(\mu)^4+\ldots$$
Reordering terms we get an expression we can integrate
$$\frac{d\mu}{\mu}=\frac{d\alpha}{-\beta_0\alpha(\mu)^2-\beta_1\alpha(\mu)^3-\beta_2\alpha(\mu)^4+\ldots}$$
for concreteness' sake let's consider the one loop case (my reasonnig will apply in general but the formulas get cumbersome) where $\beta_1=\beta_2=\ldots=0$
$$\int{}\frac{d\mu}{\mu}=\int\frac{d\alpha}{-\beta_0\alpha(\mu)^2}$$
integration yealds
$$\ln{\frac{\mu}{r}+C_1=\frac{1}{\alpha\beta_0}}+C_2$$
$C_1$ and $C_2$ are integration constants. We can simplify matters merging them simply doing $C\equiv{}C_1-C_2$ which leaves us with
$$\ln{\frac{\mu}{r}+C=\frac{1}{\alpha\beta_0}}$$
You might be wondering about the $r$. We have integrated a first order differential equation, which should leve us only with one integration constant. The $r$ comes from the fact that $\mu$ is a dimensionful parameter, so in order to write the logarithm that appears in the integral we must introduce a dimensionful parameter. It is important to notice that $r$ can be ANY dimensionful parameter (with the same sign as $\mu$ for the log to make sense) since
$$\frac{d}{d\mu}\ln{\frac{\mu}{r}}=\frac{r}{\mu}\frac{1}{r}=\frac{1}{\mu}$$
for ANY value of $r$. $r$ just dissapears after the derivative. Thus, while $C$ will be specified by some intial conditions on $\alpha$ on some point $\mu_0$ with $r$ we have complete freedom to choose whatever we want (as long as it has the same sign as $\mu$).
Of course, different $r$-s will give us different $C$-s. To simplify matters we could wonder if we can always choose an appropriate $r$ such that $C=0$. In order to see if this is indeed possible, assume that $C'$ and $r'$ are such that the initial conditions are satisfied. Then if we wish $C$ to be zero the following equation must be satisfied
$$\ln{\frac{\mu}{r'}}+C'=\ln{\frac{\mu}{r}}$$
whose solution is
$$r=r'e^{-C'}$$
The moral of the story is that for any $r'$ and $C'$ we can always pick $r$ such that $C=0$. This way the solution of the differential equation in the first equation of this answer can be written (at one loop in the beta function)
$$\ln{\frac{\mu}{r}} = \frac{1}{\alpha\beta_0} \quad \Rightarrow \quad \frac{\mu}{r} = e^{1 / \alpha\beta_0} $$
This $r$ is $\Lambda_{QCD}$. Changing names and rearranging terms we get the famous formula (you might see some $2\pi$-s floating around in some places, it's just a matter of how you define the $\beta_i$-s in the first equation of the answer)
$$\Lambda_{QCD}= \mu e^{\frac{-1}{\beta_0\alpha(\mu)}}$$
One point I would like to make. You may have heard sometimes that $\Lambda_{QCD}$ is renormalization group invariant. The analysis above should make this statement crystal clear, after all we have seen that $\Lambda_{QCD}$ is, by construction, a mere constant!
Solving for $\alpha$
$$\alpha(\mu)=\frac{1}{\beta_0\ln\frac{\mu}{\Lambda_{QCD}}}$$
Notice that $\mu=\Lambda_{QCD}$ makes the coupling become singular. Therefore we see that $\Lambda_{QCD}$ is the scale at which nonperturbative effects take over.
In order to establish its numerical value, from the equation above we readily see that $\alpha$ defines $\Lambda_{QCD}$. Therefore we can get it using the experimental value for $\alpha$.
