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?
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