bayes.py

#!/usr/bin/python

p = 0.5
sample_size = [10 ** x for x in range(1,6)]
replicates = 10000
biases = []

for n in sample_size:
	bias = np.empty(replicates)
	for i in range(replicates):
		
		true_sample = np.random.normal(size=n)
		negative_values = true_sample<0
		missing = np.random.binomial(1, p ,n).astype(bool)
		obseved_sample = true_sample[-(negative_values & missing)]
		
		bias[i] = obseved_sample.mean()
	
	biases.append(bias)