切り捨て正規分布から特定の分散の乱数をサンプルする方法はありますか？

Question

私は、結果の数値の特定の分散と範囲を考慮して、切り捨て正規分布から数値をサンプリングしようとしています。平均値0と単位分散を持つ数値が必要ですが、例えば[-2、2]のような範囲内でなければなりません。

分散を維持しながら分布を切り捨てる方法はわかりません。

import math 
import numpy as np 
import scipy.stats as stats 


truncation = 2 
lower, upper = -truncation, truncation 
mu, sigma = 0, 1 
num_samples = 1000 
if truncation: 
    n = stats.truncnorm((lower - mu)/sigma, (upper - mu)/sigma, loc=mu, scale=sigma) 
    samples = n.rvs(num_samples) 
    std_trunc = np.std(samples) 

    n = stats.norm(loc=mu, scale=sigma) 
    samples = n.rvs(num_samples) 
    std_simple = np.std(samples) 

print(std_trunc, std_simple, sep='\n') 

# outputs 
# 0.859167285015 # I need number close to 1 here 
# 1.01735583631 # like here, but here it's not truncated 

Answer 1

Wikipediaのページがobserved mean and varianceのための式を与え、私たちは私たちが望む結果を与えることをtruncnormに渡すべき値を見つけるために反転するためにそれを使用することができます。

私は標準法の作業に基づいた単純化を利用しません。部分的には一般的であり、部分的には私はまだ朝食がないので一部の算術演算をしたくありません。最小化全体を簡単な計算で置き換えます。私を与え

import numpy as np 
import scipy.stats as stats 
import scipy.optimize 

def truncated_mean_std(mu, sigma, lower, upper): 
    # N.B. lower/upper are the actual values, not Z-scaled 
    alpha = (lower - mu)/sigma 
    beta = (upper - mu)/sigma 
    d_pdf = (stats.norm.pdf(alpha) - stats.norm.pdf(beta)) 
    wd_pdf = (alpha * stats.norm.pdf(alpha) - beta * stats.norm.pdf(beta)) 
    d_cdf = stats.norm.cdf(beta) - stats.norm.cdf(alpha) 
    mu_trunc = mu + sigma * (d_pdf/d_cdf) 
    var_trunc = sigma**2 * (1 + wd_pdf/d_cdf - (d_pdf/d_cdf)**2) 
    std_trunc = var_trunc**0.5 
    return mu_trunc, std_trunc 

def trunc_samples(mu, sigma, lower, upper, num_samples=1000): 
    n = stats.truncnorm((lower - mu)/sigma, (upper - mu)/sigma, loc=mu, scale=sigma) 
    samples = n.rvs(num_samples) 
    return samples 

def corrector(mu, sigma, lower, upper): 
    target = np.array([mu, sigma]) 
    result = scipy.optimize.minimize(
     lambda x: ((target - truncated_mean_std(x[0], x[1], lower, upper))**2).sum(), 
     x0=[mu, sigma]) 
    return result.x 

：

In [79]: s = trunc_samples(mu=0, sigma=1, lower=-2, upper=2, num_samples=10**7) 

In [80]: s.mean(), s.std() 
Out[80]: (-9.8821067931585576e-05, 0.87951241887015619) 

In [81]: mu_to_use, sigma_to_use = corrector(0, 1, -2, 2) 

In [82]: mu_to_use, sigma_to_use 
Out[82]: (-7.4553057719882245e-09, 1.3778928137492246) 

In [83]: s = trunc_samples(mu=mu_to_use, sigma=sigma_to_use, lower=-2, upper=2, num_samples=10**7) 

In [84]: s.mean(), s.std() 
Out[84]: (0.0004091647648333381, 0.99991490259048865) 

In [85]: s.min(), s.max() 
Out[85]: (-1.9999995310631815, 1.9999997070340947)