samples.py (1501B)
1 # imports 2 import numpy as np 3 rng = np.random.default_rng(123) 4 DEFAULT_N = 1000000 5 6 # three simple functions 7 8 9 def normal(mean, std, n=DEFAULT_N): 10 return np.array(rng.normal(mean, std, n)) 11 12 13 def lognormal(mean, std, n=DEFAULT_N): 14 return np.array(rng.lognormal(mean, std, n)) 15 16 17 def to(low, high, n=DEFAULT_N): 18 normal95confidencePoint = 1.6448536269514722 19 logLow = np.log(low) 20 logHigh = np.log(high) 21 meanlog = (logLow + logHigh)/2 22 sdlog = (logHigh - logLow) / (2.0 * normal95confidencePoint) 23 return lognormal(meanlog, sdlog, n) 24 25 26 def mixture(samples_list, weights_array, n=DEFAULT_N): 27 normalized_weights = weights_array/sum(weights_array) 28 cummulative_sums = np.cumsum(normalized_weights) 29 helper_probs = rng.random(n) 30 results = np.empty(n) 31 for i in range(n): 32 helper_list = [j for j in range( 33 len(cummulative_sums)) if cummulative_sums[j] > helper_probs[i]] 34 if len(helper_list) == 0: 35 helper_loc = 0 # continue 36 print("This should never happen") 37 else: 38 helper_loc = helper_list[0] 39 target_samples = samples_list[helper_loc] 40 result = rng.choice(target_samples, 1)[0] 41 results[i] = result 42 return (results) 43 44 45 # Example 46 p_a = 0.8 47 p_b = 0.5 48 p_c = p_a * p_b 49 50 dists = [[0], [1], to(1, 3), to(2, 10)] 51 # print(dists) 52 weights = np.array([1 - p_c, p_c/2, p_c/4, p_c/4]) 53 # print(weights) 54 result = mixture(dists, weights) 55 mean_result = np.mean(result) 56 print(mean_result)