I want to point out another interesting solution method, which also generalizes the result to the expectation of $X^{-m}$ for integer $m=1,2,3,\dots$. I will use moment generating functions (mgf) and the results from the paper by N Cressie et.al http://amstat.tandfonline.com/doi/abs/10.1080/00031305.1981.10479334?journalCode=utas20 "The Moment-Generating Function and Negative Integer Moments".
They give the result that when $X$ is a positive random variable with mgf $M_X(t)$ which is defined in an open neighbourhood of the origin, then we have
$$
\DeclareMathOperator{\E}{\mathbb{E}}
\E X^{-m} = \Gamma(m)^{-1} \int_0^\infty t^{m-1} M_X(-t) \; dt
$$
for positive integers $m$.
It is known that for the beta distribution, the mgf is given by a confluent hypergeometric function as
$$
M_X(t) = {}_1F_1(\alpha;\alpha+\beta;t)
$$
so using the result above gives that
$$
\E X^{-m} = \Gamma(m)^{-1} \int_0^\infty t^{m-1} {}_1F_1(\alpha;\alpha+\beta;-t)\; dt
$$
I evaluated that integral with the help of maple:
assume( a>0, b>0 );assume(m-1,posint)
GAMMA(m)^(-1) * int( t^(m-1)*hypergeom([a],[a+b],-t), t=0..infinity )
GAMMA(a + b) GAMMA(a - m)
-------------------------
GAMMA(a) GAMMA(a + b - m)
so finally we can write the result as
$$
\E X^{-m} = \frac{\Gamma(\alpha+\beta)\Gamma(\alpha-m)}{\Gamma(\alpha)\Gamma(\alpha+\beta-m)}
$$
which coinsides with the other answer for $m=1$. Then some human mathematics is needed to conclude that we need the assumption $\alpha > m$ for this to be valid.