A neutron outside the nucleus lives for about 15 minutes and decays mainly through weak decays (beta decay). Many other weakly decaying particles decay with lifetimes between $10^{-10}$ and $10^{-12}$ seconds, which is consistent with $\alpha\_W \simeq 10^{-6}$.

Why does the neutron lives so much longer than the others?
['https://physics.stackexchange.com/questions/31514', 'https://physics.stackexchange.com', 'https://physics.stackexchange.com/users/1502/']
NB: I feel like this is a pretty half-assed job, and I apologize for that but having opened my mouth in the comments I guess I have to write *something* to back it up.

---

We start with Fermi's golden rule for all transitions. The probability of the transition is
$$ P\_{i\to f} = \frac{2\pi}{\hbar} \left|M\_{i,f}\right|^2 \rho $$
where $\rho$ is the density of final states which is proportional to $p^2$ for massive particles. To find the rate1 for all possible final states we sum over these probabilities incoherently. When the mass difference between the initial and final states is much less than the $W$ mass the matrix element $M\_{i,f}$ depends only weakly (hah!) on the particular state and the sum is well approximated by a sum *only* over the density of states:
$$P\_\text{decay} \approx \frac{2\pi}{\hbar} \left|M\_{}\right|^2 \int\_\text{all outcomes} \rho .$$
This sum is collectively called the phase space available to the decay. In these cases the matrix element is also quite small for [the reason that Dr BDO discusses](https://physics.stackexchange.com/a/31526/520).

The phase space computation can be quite complicated as it must be taken over all unconstrained momenta of the products. For decays to two body states it turns out to be easy, there is no freedom in the final states except the $4\pi$ angular distribution in the decay frame (their are eight degrees of freedom in two 4-vectors, but 2 masses and the conservation of four momentum account for all of them except the azimuthal and polar angles of one of the particles).

The decays that you have asked about are to three body states. That gives us twelve degrees of freedom less three constraints from masses, four from conservation of 4-momentum which leaves five. Three of these are the Euler angles describing the orientation of the decay (and a factor of $8\pi^2$ to $\rho$), so our sum is over two non-trival momenta. The integral looks something like
$$
\begin{array}\\
\rho \propto \int p\_1^2 \mathrm{d}p\_1 \int p\_2^2 \mathrm{d}p\_2 \int \mathrm{d}(\cos\theta)
&\delta(m\_0 - E\_1 - E\_2-E\_3 ) \\
&\delta(E\_1^2 - m\_1^2 - p\_1^2) \\
&\delta(E\_2^2 - m\_2^2 - p\_2^2) \\
&\delta(E\_2^2 - m\_2^2 - p\_2^2) \\
&\delta(\vec{p}\_1 + \vec{p}\_2 + \vec{p}\_3)
\end{array}
$$
which is easier to compute in Monte Carlo than by hand. (BTW--the reason for introducing the seemingly redundant integral over the angle $\theta$ between the momenta of particles 1 and 2 will become evident in a little while).

For beta decays the remnant nucleus is very heavy compared to the released energy, [which simplifies the above in one limit](http://hyperphysics.phy-astr.gsu.edu/hbase/quantum/fermi2.html).

In the case of muon decay, it is not unreasonable to treat all the products as ultra-relativistic, and the above reduces to 
$$
\begin{array}\\
\rho \propto \int p\_1^2 \mathrm{d}p\_1 \int p\_2^2 \mathrm{d}p\_2 \int \mathrm{d}(\cos\theta)
&\delta(m\_0 - E\_1 - E\_2 - E\_3 ) \\
&\delta(E\_1 - p\_1) \\
&\delta(E\_2 - p\_2) \\
&\delta(E\_3 - p\_3) \\
&\delta(\vec{p}\_1 + \vec{p}\_2 + \vec{p}\_3) \\
= \int p\_1^2 \mathrm{d}p\_1 \int p\_2^2 \mathrm{d}p\_2 \int \mathrm{d}(\cos\theta)
&\delta(m\_0 - p\_1 - p\_2 - p\_3 ) \\
&\delta(\vec{p}\_1 + \vec{p}\_2 + \vec{p}\_3) \\
= \int p\_1^2 \mathrm{d}p\_1 \int p\_2^2 \mathrm{d}p\_2 \int \mathrm{d}(\cos\theta)
&\delta(m\_0 - p\_1 - p\_2 - \left|\vec{p}\_1 + \vec{p}\_2\right| )\\
= \int p\_1^2 \mathrm{d}p\_1 \int p\_2^2 \mathrm{d}p\_2 \int \mathrm{d}(\cos\theta)
&\delta\left(m\_0 - p\_1 - p\_2 - \sqrt{p\_1^2 + p\_2^2 - p\_1p\_2\cos\theta} \right)
\end{array} 
$$
The integral over the angle will evaluate to one in some regions and zero in others and as such is equivalent to correctly assigning the limits of the other two integrals, so writing $\delta m = m\_0 - m\_1 - m\_2 - m\_3$ we get
$$
\begin{array}
\rho 
& \propto \int\_0^{\delta m/2} p\_1^2 \mathrm{d}p\_1 \int\_0^{\delta m-p\_1} p\_2^2 \mathrm{d}p\_2 \\
& \propto \int\_0^{\delta m/2} p\_1^2 \mathrm{d}p\_1 \left[ \frac{p\_2^3}{3}\right]\_{p\_2=0}^{\delta m-p\_1} \\
& \propto \int\_0^{\delta m/2} p\_1^2 \mathrm{d}p\_1 \frac{(\delta m - p\_1)^3}{3}
\end{array}
$$
which I am not going to bother finishing but shows that that phase space can vary as a high power of the mass difference (up to the sixth power in this case).

---

1 The lifetime of the state is inversely proportional to the probability
As you correctly state, the neutron decay is an decay due to the weak interaction, these are quite a bit slower than other decays due the mass of the
intermediate W boson, 81GeV, which slows the reaction, additionally the neutron decay only liberates a small amount of energy, around 1 MeV, it is the ratio of the liberated energy to the mass of the W which sets the speed of the reaction which is thus much slower than other decays as all the other particle decay release much more energy.