How can I calculate the α and β parameters for a Beta distribution if I know the mean and 95% upper boundary I want the distribution to have? As an input I have the mean (expected value) of a probability and the 95% upper bound. Now, I would like to define the α and β parameters for a Beta distribution such that I can do Monte Carlo simulations. I know the expressions for these parameters in terms of E[X] and Var[X], but do not have Var[X] directly available. Is it possible to find Var[X] with the given upper bound, or are there other ways to find α and β?

Here's an R solution using a root solver that exploits the relationship between the mean and the parameters. Because $\bar{x}=\frac{\alpha}{\alpha + \beta}$, once we know $\alpha$, we also know $\beta$ through $\beta = \frac{\alpha - \bar{x}\alpha}{\bar{x}}$.
# Setting the true parameters for illustration
alpha <- 2
beta <- 46
q <- 0.975 # The quantile

# The mean and quantile
bmean <- alpha/(alpha + beta)
qval <- qbeta(q, alpha, beta)

# Function to find alpha using root solver
find_alpha <- function(alpha, mean, qval, q) {
  beta <- (alpha - alpha*mean)/mean # Beta can be calculated from the known mean and the estimated alpha
  qbeta(q, alpha, beta) - qval
}

# Estimate alpha using a root solver
alphaest <- uniroot(find_alpha, interval = c(0, 1000), mean = bmean, qval = qval, q = q, extendInt = "upX")$root

# Calculate beta exploiting the relationship between alpha and the known mean
betaest <- (alphaest - alphaest*bmean)/bmean

# The final estimates
c(alphaest, betaest)
[1]  2.00000 46.00001

