pythonrstatisticsnormal-distributiontruncated

Find mean of truncated normal distribution in R


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?


Solution

  • 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