Let $X=(X_1,\dots,X_n)\sim N(\mu,\sigma)$ be a Gaussian random vector of $d$ dimensions with mean $\mu$ and diagonal variance $\sigma^2$. Thus, both $\mu,\sigma\in\mathbb{R}^d$.
Let $I=\arg\max_{i\in[1,n]}X_i$. I would like to get an analytical expression of the distribution of $I$ or an analytical approximation.
Since the variables are independent we can obtain an expression by integrating w.r.t. the maximum $x$ and adding up the probability of variable $i$ being the argmax at exactly that maximum: $$\mathbb{E}_{x\sim N(\mu_i,\sigma_i)}\left[\prod_{j\neq i} cdf(\mu_j,\sigma_j)(x)\right] = \int_{x=-\infty}^{\infty} N(\mu_i,\sigma_i)(x)\left(\prod_{j\neq i} cdf(\mu_j,\sigma_j)(x)\right) dx$$ where the expression within the integral is the probability of variable $i$ being exactly $x$ and all others being less than $x$, multiplying because of independence.
I don't know how to obtain a formula for this integral. I'm searching for an expression without integrals, but which can use the gaussian cdf.
I've found related questions[1,2,3,4], but they either talk about numerical algorithms or show that the question is very hard in the general case when $\Sigma$ is a matrix. I would need an analytical expression, but I'm happy with approximations, diagonal variance, or even $\sigma=I$.
