Distribution of an RBF-transformed normal variable My question might be related to this or this one, but I have reasons to hope my problem is more benign.
Assume I have a normally distributed variable $X \sim  N(0, 1)$. What can be said of $Y = \exp(-\frac{X^2}{2 \sigma^2})$? It is obviously bounded to $(0, 1]$ and it certainly looks like beta-distribution:


However, can this be shown mathematically? How do $\alpha$ and $\beta$ depend on $\sigma$? And, assuming that $Y$ is really beta distributed: What if $X$ is multidimensional? I reckon that in that case we'd be dealing with $Y = \exp(-\frac{Z}{2 \sigma^2})$, with $Z$ having a $\chi^2$ distribution. Would $Y$ still have beta distribution and how would degrees of freedom come into play?

Expanding on whuber's answer (+1), I'd like to cover a more general case, as stated in the question. We have a standard normal random variable $X \sim N(0, 1)$. We transform this variable by passing it through an RBF, $Y = \exp(-\frac{X^2}{2 \sigma^2})$. Notice that $\sigma=1$ for $X$, but may be any positive number for $Y$.
Applying the same reasoning as in the referenced answer, it is easy to see that
$$
f_Y(y) = \frac{\sigma \cdot y^{\sigma^2-1}}{\sqrt{-\pi\log(y)}}.
$$
Of course, by setting $\sigma = 1$ the formula reduced to the one provided by whuber. Neither formula corresponds to beta distribution, but, interestingly, with the right choice of parameters the approximation is quite good. Below are some examples for different values of $\sigma$:





For the reference, here is the code to reproduce the calculations:
import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as st
from scipy.optimize import minimize


# approximate the PDF based on the histogram of X:
def pdf(X):
  h = np.histogram(X, bins=len(X)//100) # the histogram of the data
  x = (h[1][:-1] + h[1][1:])/2          # mid-bin coordinates
  d = (max(h[1]) - min(h[1])) * len(X) / len(h[0])  # dx, for computing the PDF
  return x, h[0]/d
  

# the transform by the Gaussian kernel
def rbf_transform(x, sigma):
  y = np.exp(-x**2/2/sigma**2)
  return y


def transformAndFit(X, sigma, scale='linear', error='abs'):
  Y = rbf_transform(X, sigma)
  y, fY = pdf(Y)
  trueFn = sigma * y**(sigma**2 - 1) / np.sqrt(-np.pi*np.log(y))
  if scale=='linear':
    if error=='abs':
      def objFn(param):
        return np.sum(np.abs(st.beta.pdf(y, param[0], param[1]) - trueFn))
    else:
      def objFn(param):
        return np.sum((st.beta.pdf(y, param[0], param[1]) - trueFn)**2)
  else:
    if error=='abs':
      def objFn(param):
        return np.sum(np.abs(np.log(st.beta.pdf(y, param[0], param[1])) - np.log(trueFn)))
    else:
      def objFn(param):
        return np.sum((np.log(st.beta.pdf(y, param[0], param[1])) - np.log(trueFn))**2)
  
  minFun = np.inf
  alpha = 1.1
  beta  =  .5
  for a in np.arange(.3, 1.3, .1):
    for b in np.arange(.1, .7, .05):
      res = minimize(objFn, [a, b], bounds=((0, None), (0, None)))
      if res.fun < minFun:
        alpha = res.x[0]
        beta  = res.x[1]
  
  plt.plot(y, trueFn, '-', c='steelblue', label='true $f_Y(y)$', lw=2)
  plt.plot(y, fY, '+', c='firebrick', label='simulation')
  plt.plot(y, st.beta.pdf(y, alpha, beta), ':', c='silver', lw=5,
           label=f'B($\\alpha$={alpha:.2f}, $\\beta$={beta:.2f})')
  plt.title(f'RBF($\\sigma$={sigma})-transformed data')
  plt.legend()
  plt.xlabel('y')
  plt.ylabel('$f_Y(y)$')
  if scale == 'log': plt.yscale('log')
  plt.show()


# the starting distribution: standard normal
q = np.linspace(.001, .999, 10000) # 10,000 samples, uniformly over quantiles
X = st.norm.ppf(q)                 # normally distributed samples
x, fX = pdf(X)
plt.plot(x, st.norm.pdf(x), label='theoretical')
plt.plot(x, fX, ':', c='firebrick', label='simulation', lw=3)
plt.title('Original: Standard normal distribution')
plt.legend()
plt.show()

# the rbf-transformed data, the true and beta fit:
for sigma in np.arange(.5, 1.5+.25, .25):
  transformAndFit(X, sigma, scale='log')

