13
$\begingroup$

I have a simple setting: I simulate demand patterns that are distributed according to a truncated normal distribution with a given mean and variance after truncation. The truncation is from the left at zero (as demand is always positive).

However, all functions that draw from truncated normal distributions require me to specify the mean and variance of the normal distribution before truncation as e.g. the truncnorm function from the truncnorm package in R.

Example:

val = rtruncnorm(10000, a=0,  mean = 100, sd = 240)
print(mean(val))
[1] 232.2385
print(sd(val))
[1] 162.853

I would like the resulting variable val to have mean 100 and standard deviation 240. Is there any function in R to do this or is there a formula to compute the mean and variance before truncation given the truncation limits and the post-mean and post-variance?

$\endgroup$
6
  • 3
    $\begingroup$ Formulas for post mean and variance in terms of pre mean and pre variance are in wikipedia. This gives you two non-linear equations that can be solved with rootSolve::multiroot. $\endgroup$ May 13, 2019 at 20:39
  • $\begingroup$ Thanks for the suggestion. I checked wikipedia. The equations involve $\alpha$, which depends on the untruncated mean and standard deviation as well. Together with this, I can't see a way to use these equations to get values for the truncated mean and sd. $\endgroup$
    – chris12
    May 13, 2019 at 21:20
  • 1
    $\begingroup$ That comment appears to contradict your post. The post strongly implies there is truncation to the left of $a=0$ only. Is that how we should continue understanding the question? $\endgroup$
    – whuber
    May 13, 2019 at 21:58
  • 2
    $\begingroup$ I wonder whether mean 100 and standard deviation 240 is actually possible for the truncated distribution? Is there a theoretical result that all combinations of (positive) mean and variance can be attained? $\endgroup$
    – jochen
    May 24, 2019 at 8:44
  • 1
    $\begingroup$ The question you are asking is solved in this paper. $\endgroup$ Aug 22, 2020 at 5:12

2 Answers 2

10
$\begingroup$

What Jarle Tufto points out is that, if $X\sim\mathcal N^+(\mu,\sigma^2)$, then, defining $$\alpha=-\mu/\sigma\quad\text{and}\quad\beta=\infty$$ Wikipedia states that $$\Bbb E_{\mu,\sigma}[X]= \mu + \frac{\varphi(\alpha)-\varphi(\beta)}{\Phi(\beta)-\Phi(\alpha)}\sigma$$ and $$\text{var}_{\mu,\sigma}(X)=\sigma^2\left[1+\frac{\alpha\varphi(\alpha)-\beta\varphi(\beta)}{\Phi(\beta)-\Phi(\alpha)} -\left(\frac{\varphi(\alpha)-\varphi(\beta)}{\Phi(\beta)-\Phi(\alpha)}\right)^2\right]$$ That is, $$\Bbb E_{\mu,\sigma}[X]= \mu + \frac{\varphi(\mu/\sigma)}{1-\Phi(-\mu/\sigma)}\sigma\tag{1}$$ and $$\text{var}_{\mu,\sigma}(X)=\sigma^2\left[1-\frac{\mu\varphi(\mu/\sigma)/\sigma}{1-\Phi(-\mu/\sigma)} -\left(\frac{\varphi(\mu/\sigma)}{1-\Phi(-\mu/\sigma)}\right)^2\right]\tag{2}$$ Given the numerical values of the truncated moments $(\Bbb E_{\mu,\sigma}[X],\text{var}_{\mu,\sigma}(X))$, one can then solve numerically (1) and (2) as a system of two equations in $(\mu,\sigma)$, assuming $(\Bbb E_{\mu,\sigma}[X],\text{var}_{\mu,\sigma}(X))$ is a possible value for a truncated Normal $\mathcal N^+(\mu,\sigma^2)$.

$\endgroup$
1
$\begingroup$

To complement Xi'an's answer, one simple approach to solve numerically such equations would be by using a Taylor series approximation.

If we denote $m,s^2$ the target mean and variance of the truncated distributrion, and $\mu,\sigma^2$ the mean and variance of the underlying Gaussian distribution then: \begin{align} m &= \mu+z(\mu/\sigma)\,\sigma\\ s^2 &= \left[1-z(\mu/\sigma)\,\mu/\sigma-z(\mu/\sigma)^2 \right]\sigma^2\\ z(\mu/\sigma) &= \phi(\mu/\sigma)/\Phi(\mu/\sigma) \end{align}

The function $z(\mu/\sigma)$ can be approximated with a Taylor's series expansion around $m/s$. Denote $\hat{z}=z(m/s)$, then: \begin{align} z(\mu/\sigma) &\approx \hat{z}+(\mu/\sigma-m/s)\frac{\phi'(m/s)\Phi(m/s)-\phi(m/s)^2}{\Phi(m/s)^2}\\ &= \hat{z}+(\mu/\sigma-m/s)\,(-\hat{z} m/s-\hat{z}^2)\\ &= [1+(m/s)^2+\hat{z} m/s]\hat{z}-(\mu/\sigma)(m/s+\hat{z})\hat{z}\\ &= A - B\mu/\sigma \end{align}

Where $A=[1+(m/s)^2+\hat{z} m/s]\hat{z}>0$ and $B=(m/s+\hat{z})\hat{z}>0$. Using this approximation into the first equation we get: \begin{align} m &\approx \mu + (A-B\mu/\sigma)\sigma =(1-B)\mu+A\sigma\\ \Rightarrow \mu &\approx (1-B)^{-1}(m-A\sigma) \end{align}

Replacing both approximations into the second equation we get: \begin{align} s^2 &\approx [1-(A-B\mu/\sigma)\mu/\sigma-(A-B\mu/\sigma)^2]\sigma^2\\ &= (1-A^2)\sigma^2-A(1-2B)\mu\sigma+B(1-B)\mu^2\\ &= (1-A^2)\sigma^2-\frac{A(1-2B)}{1-B}(m-A\sigma)\sigma+\frac{B}{1-B}(m-A\sigma)^2\\ \Rightarrow (1-B)s^2 &\approx (1-B)(1-A^2) \sigma^2 -A(1-2B)(m-A\sigma)\sigma +B(m-A\sigma)^2\\ 0 &\approx (1-B)\sigma^2-Am\sigma + Bm^2-(1-B)s^2 \end{align}

Such equation has as an approximate solution for $\sigma^2$: \begin{equation} \hat{\sigma}^2 =\frac{2BD+C^2-2D+C\sqrt{C^2-4D(1-B)}}{2(1-B)^2} \end{equation}

With $C=Am$, $D=Bm^2-(1-B)s^2$, and provided that $C^2-4D(1-B)>0$ and $\hat{\sigma}^2>0$. The approximate solution for $\mu$ is then: $$ \hat{\mu} = \frac{m-A\hat{\sigma}}{1-B} $$

Thus, given $(m,s)$, one approximate solution $(\hat{\mu},\hat{\sigma}^2)$ for $({\mu},{\sigma}^2)$ can be found following these steps:

  1. $\hat{z}=\phi(m/s)/\Phi(m/s)$
  2. $A=[1+(m/s)^2+\hat{z} m/s]\hat{z}$; $B=(m/s+\hat{z})\hat{z}$
  3. $C=Am$; $D=Bm^2-(1-B)s^2$
  4. $\hat{\sigma}^2 =\frac{2BD+C^2-2D+C\sqrt{C^2-4D(1-B)}}{2(1-B)^2}$ (provided $C^2-4D(1-B)>0$ and $\hat{\sigma}^2>0$)
  5. $\hat{\mu}=\frac{m-A\hat{\sigma}}{1-B}$
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.