Expectation-Maximization Calculator for Two Normal Distributions

EMCalc.py

import math

print("This program assumes a standard deviation of 1 for k=2 normal distributions")
print("Enter your initial values (on separate lines) and a blank line when you've finished.")

# Takes in an unknown number of value arguments.
initialVals = []
while True:
    val = input()
    if val == "":
        break
    initialVals.append(val)

u1, u2 = input('Enter your intial h values: ').split()
u1 = float(u1)
u2 = float(u2)

def probabilityU1(x):
    x = float(x)
    return math.exp( (-1/2) * math.pow((u1 - x), 2) ) / ( math.exp( (-1/2) * math.pow((u1 - x), 2) ) + math.exp( (-1/2) * math.pow((u2 - x), 2) ) )

def probabilityU2(x):
    return 1 - probabilityU1(x)

def newBound(probabilityFunction):
    numerator = []
    denominator = []
    for x in initialVals:
        numerator.append( float(x) * probabilityFunction(x) )
        denominator.append( probabilityFunction(x) )
    return sum(numerator) / sum(denominator)



for index, value in enumerate(initialVals):
    print('x' + str(index + 1) + ': ' + str(value))
    print('E[z' + str(index + 1) + '1]: ' + str(probabilityU1(value)) )
    print('E[z' + str(index + 1) + '2]: ' + str(probabilityU2(value)) )

print("New u1 bound: " + str(newBound(probabilityU1)))
print("New u2 bound: " + str(newBound(probabilityU2)))