Find mean of truncated normal distribution in R

Question

How can I find the mean of a truncated normal distribution in R, where the lower bound a is 1, sigma is 2, and mu is 2.5? I have used the truncnorm library, but it does not have any functions for moments. In Python, I tried the following code:

import numpy as np
from scipy.stats import truncnorm
a, b = 1, np.inf
mean, var, skew, kurt = truncnorm.stats(a, b, moments='mvsk')
print(mean)

which gives mean = 1.52513528. How can I achieve the same result in R?

Answer 1

In your python code, you are not setting location and scale and thus taking the default values for location = 0 and scale = 1. Thats why you get 1.525. You should consider the following:

Your example in python:

import numpy as np
from scipy.stats import truncnorm
a, b = 1, np.inf
mean, var, skew, kurt = truncnorm.stats(a, b, moments='mvsk')
print(mean)
1.525135276160982

In R you could simply do:

a <- 1
b <- Inf
diff(dnorm(c(b, a)))/diff(pnorm(c(a,b)))
[1] 1.525135

truncnorm::etruncnorm(a, b)
[1] 1.525135

To make use of the provided data:

Python

import numpy as np
from scipy.stats import truncnorm
a, b = 1, np.inf
mu, sigma = 2.5, 2
mean, var= truncnorm.stats((a-mu)/sigma, (b - mu)/sigma, mu, sigma)
print(mean)
3.278764113471854

In R you can write a simple code to compute the mean:

truncnorm::etruncnorm(a, b, mu, sigma)
[1] 3.278764

You can always confirm your answer using base R:

qnorm(pnorm(diff(dnorm(c(b, a), mu, sigma))/diff(pnorm(c(a,b), mu, sigma))), mu, sigma^2)
[1] 3.278764