bayes.py

import random


ACCURACY = 0.75


def simulate(sensor_count, rain_prob):
    real_rain_count, total_count = 0, 0

    for iter in range(100000):
        is_rain = random.random() < rain_prob
        readouts = [
            is_rain if random.random() < ACCURACY else not is_rain
            for _ in range(sensor_count)
        ]
        if all(readouts):
            total_count += 1
            if is_rain:
                real_rain_count += 1

    return real_rain_count / total_count


def solve(sensor_count, rain_prob):
    all_true = rain_prob * ACCURACY ** sensor_count
    all_false = (1.0 - rain_prob) * (1.0 - ACCURACY) ** sensor_count

    return all_true / (all_true + all_false)


if __name__ == '__main__':
    for sensor_count in range(1, 10):
        for rain_prob in [0.1 * i for i in range(1, 10)]:
            print(f'{sensor_count} sensors, {rain_prob} rain probability: {rain_prob}')
            expected = simulate(sensor_count, rain_prob)
            actual = solve(sensor_count, rain_prob)
            print(f'\texpected = {expected:.6f}')
            print(f'\tactual = {actual:.6f}')