Fast approximation to inverse Beta CDF I am looking for a fast approximation to the inverse CDF of the Beta distribution. The approximation need not be precise, but more stress is on simplicity (I'm thinking Taylor expansion of the first 1 or 2 terms). Any ideas?
A few more details: in my Beta distribution, $a,b$ are integers. Also the value I am evaluating the inverse of is always 0.05. 
 A: For applications like this, you have to read the literature carefully, because somebody else's approximation might not reproduce the characteristic you are interested in.
There are simple solutions based on precomputing values within the range where a simple approximation (like a Normal approximation) might not do a good job.  In this case, you would worry most about where skewness is high: that would be where the Beta parameters are very unequal and relatively small.
One solution is to create an interpolation function or lookup table based on this precomputation.  To find the inverse CDF at $0.95$, first look up its value.  If the value is found in the table, use that; otherwise, use the generic approximation.  Given the ease with which millions of values could be precomputed and stored, this is the way to go for intensive simulations on a PC>
If RAM is really tight or, for some reason, you really want a formula, then the lookup table needs to be replaced by a formula.  A good tool for fitting formulas for this purpose is Eureqa Formulize, a free download.  This software identifies and fits functions to tabulated data (using a genetic algorithm).  It is extremely fast, easy to learn, and fun to watch in action.
Using a table of the $0.95$ quantile of Beta distributions for integer values $(a,b)$ in the range $1\ldots 10$ and minimizing worst case error, I have found--by running the software while writing this paragraph--a large number of approximations.  The base case is the normal approximation itself, $z = \text{mean} + \Phi^{-1}(0.95)\times\text{SD}$.  Its worst case error is around $0.085$ (which is not terribly good).  The formulas I searched for are in terms of $a$, $b$, $z$, the skewness (whose formula in terms of $a$ and $b$ you can look up on Wikipedia), and the central third moment$s_3$ (obtained by multiplying the skewness by the cube of the SD).  Among the simpler formulas are
$$0.0012526 + z + 8.90369 s_3 + 0.0344753 z\times \text{skewness}$$
with a worst-case error of $0.011$--pretty good.  This is recognizable as the Normal approximation $z$ plus corrections for (a) using the worst-case error as a criterion, which tends to introduce a small bias ($0.0013$) and (b) the skewness (as expected).  This formula would be used whenever both $a$ and $b$ are $10$ or less; otherwise, you would revert to the Normal approximation ($z$ itself).
Here is a plot comparing my tabulated values (points) against the formula, taken directly from the software:

(The values jump around because this table of course is two-dimensional: it was originally organized by $a$ and $b$, but then flattened into a traditional spreadsheet format for this analysis.)
By applying this approach to a tabulated set of values you are most interested in, you will obtain different formulas.  Pick one that best balances the desired accuracy against the complexity of the formula.
A: whuber's comment is specially true if you want to compute beta inequalities (e.g. the probability that a beta distributed r.v. exceeds another beta distributed r.v.). See here for more details.
Another alternative, pointed out by J. Cook in his blog is to use the Kumaraswamy distribution to approximate the Beta distribution.
To adhere by the cross-validated rule, I'll briefly summarize the main elements here but go read that post. The CDF of a $K(a, b)$ variable is:
$$F(x|a,b)=1–(1–x^a)^b$$ 
This CDF is easy to inverse and to draw from so you can easily generate a draw from $K(a,b)$ by generating a $u$ uniform on $[0,1]$ and:
$$F^{-1}(u)=(1 – (1 – u)^{1/b})^{1/a}$$ 
The idea is that you can pick up the parameters $a$ and $b$ to match most $\text{beta}(\alpha,\beta)$ reasonably well: pick $a=\alpha$ and use a numerical approach to find the optimal $b|a$. Again, more details are given on the original post. Alternatively, you could use a pre-computed table of equivalent parameter sets for Beta and Kumaraswamy distributions (see here) or compute one yourself. 
